CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MuonErrorMatrixAdjuster.cc
Go to the documentation of this file.
2 
3 #include "TString.h"
4 #include "TMath.h"
7 
10 
12 
15 
17 {
18  theCategory="MuonErrorMatrixAdjuster";
19  theInstanceName = iConfig.getParameter<std::string>("instanceName");
20  //register your products
21  produces<reco::TrackCollection>( theInstanceName);
22  produces<TrackingRecHitCollection>();
23  produces<reco::TrackExtraCollection>();
24 
25  theTrackLabel = iConfig.getParameter<edm::InputTag>("trackLabel");
26  theRescale = iConfig.getParameter<bool>("rescale");
27 
28  theMatrixProvider_pset = iConfig.getParameter<edm::ParameterSet>("errorMatrix_pset");
29 }
30 
31 
33 {
34 
35  // do anything here that needs to be done at desctruction time
36  // (e.g. close files, deallocate resources etc.)
37 
38 }
39 
40 
41 //
42 // member functions
43 //
44 
45 //take the error matrix and rescale it or just replace it
47 {
48  //CovarianceMatrix is template for SMatrix
49  reco::TrackBase::CovarianceMatrix revised_matrix( theMatrixProvider->get(momentum));
50 
51  if( theRescale){
52  //rescale old error matrix up by a factor
53  multiply(revised_matrix,error_matrix);
54  }
55  return revised_matrix;
56 }
57 
59  //scale term by term the matrix
60  // the true type of the matrix is such that [i][j] is the same memory object as [j][i]: looping i:0-5, j:0-5 double multiply the terms
61  // need to loop only on i:0-5, j:i-5
62  for(int i = 0;i!=5;i++){for(int j = i;j!=5;j++){
63  revised_matrix(i,j)*=scale_matrix(i,j);
64  }}
65 }
67  //divide term by term the matrix
68  // the true type of the matrix is such that [i][j] is the same memory object as [j][i]: looping i:0-5, j:0-5 double multiply the terms
69  // need to loop only on i:0-5, j:i-5
70  for(int i = 0;i!=5;i++){for(int j = i;j!=5;j++){
71  if (denom_matrix(i,j)==0) return false;
72  num_matrix(i,j)/=denom_matrix(i,j);
73  }}
74  return true;
75 }
76 
78  const FreeTrajectoryState & PCAstate){
79  //get the parameters of the track so I can reconstruct it
80  double chi2 = recotrack_orig.chi2();
81  double ndof = recotrack_orig.ndof();
82  math::XYZPoint refpos = recotrack_orig.referencePoint();
83  math::XYZVector mom = recotrack_orig.momentum();
84  int charge = recotrack_orig.charge();
85 
86  reco::TrackBase::CovarianceMatrix covariance_matrix = fix_cov_matrix(recotrack_orig.covariance(),PCAstate.momentum());
87 
88  LogDebug(theCategory)<<"chi2: "<<chi2
89  <<"\n ndof: "<<ndof
90  <<"\n refpos: "<<refpos
91  <<"\n mom: "<<mom
92  <<"\n charge: "<<charge
93  <<"\n covariance:\n"<<recotrack_orig.covariance()
94  <<"\n replaced by:\n"<<covariance_matrix;
95 
96  return reco::Track(chi2,ndof,refpos,mom,charge,covariance_matrix);
97 }
98 
102  //get the 5x5 matrix of recotrack/recotrack_orig
103  reco::TrackBase::CovarianceMatrix scale_matrix(recotrack.covariance());
104  if (!divide(scale_matrix,recotrack_orig.covariance())){
105  edm::LogError( theCategory)<<"original track error matrix has term ==0... skipping.";
106  return 0; }
107 
108  const reco::TrackExtraRef & trackExtra_orig = recotrack_orig.extra();
109  if (trackExtra_orig.isNull()) {
110  edm::LogError( theCategory)<<"original track has no track extra... skipping.";
111  return 0;}
112 
113  //copy the outer state. rescaling the error matrix
114  reco::TrackBase::CovarianceMatrix outerCov(trackExtra_orig->outerStateCovariance());
115  multiply(outerCov,scale_matrix);
116 
117  //copy the inner state, rescaling the error matrix
118  reco::TrackBase::CovarianceMatrix innerCov(trackExtra_orig->innerStateCovariance());
119  multiply(innerCov,scale_matrix);
120 
121  //put the trackExtra
122  TEcol.push_back(reco::TrackExtra (trackExtra_orig->outerPosition(), trackExtra_orig->outerMomentum(), true,
123  trackExtra_orig->innerPosition(), trackExtra_orig->innerMomentum(), true,
124  outerCov, trackExtra_orig->outerDetId(),
125  innerCov, trackExtra_orig->innerDetId(),
126  trackExtra_orig->seedDirection()));
127 
128  //add a reference to the trackextra on the track
130 
131  //return the reference to the last inserted then
132  return &(TEcol.back());
133 }
134 
135 
138  reco::TrackExtra & trackextra,
139  TrackingRecHitCollection& RHcol){
140  //loop over the hits of the original track
141  trackingRecHit_iterator recHit = recotrack_orig.recHitsBegin();
142  unsigned int irh=0;
143  for (; recHit!=recotrack_orig.recHitsEnd();++recHit){
144  //clone it. this is meandatory
145  TrackingRecHit * hit = (*recHit)->clone();
146 
147  //put it on the new track
148  recotrack.setHitPattern(*hit,irh++);
149  //copy them in the new collection
150  RHcol.push_back(hit);
151  //do something with the trackextra
152  trackextra.add(TrackingRecHitRef( theRefprodRH, theRHi++));
153 
154  }//loop over original rechits
155 
156  return true; //if nothing fails
157 }
158 
159 bool MuonErrorMatrixAdjuster::selectTrack(const reco::Track & recotrack_orig){return true;}
160 
161 // ------------ method called to produce the data ------------
162 void
164 {
165  using namespace edm;
166 
167  //open a collection of track
170  LogDebug( theCategory)<<"considering: "<<tracks->size()<<" uncorrected reco::Track from the event.("<< theTrackLabel<<")";
171 
172  //get the mag field
173  iSetup.get<IdealMagneticFieldRecord>().get( theField);
174 
175  //prepare the output collection
176  std::auto_ptr<reco::TrackCollection> Toutput(new reco::TrackCollection());
177  std::auto_ptr<TrackingRecHitCollection> TRHoutput(new TrackingRecHitCollection());
178  std::auto_ptr<reco::TrackExtraCollection> TEoutput(new reco::TrackExtraCollection());
180  theTEi=0;
182  theRHi=0;
183 
184 
185 
186  for(unsigned int it=0;it!=tracks->size();it++)
187  {
188  const reco::Track & recotrack_orig = (*tracks)[it];
190  if (PCAstate.position().mag()==0)
191  {edm::LogError( theCategory)<<"invalid state from track initial state in "<< theTrackLabel <<". skipping.";continue;}
192 
193  //create a reco::Track
194  reco::Track recotrack = makeTrack(recotrack_orig,PCAstate);
195 
196  //make a selection on the create reco::Track
197  if (!selectTrack(recotrack)) continue;
198 
199  Toutput->push_back(recotrack);
200  reco::Track & recotrackref = Toutput->back();
201 
202  //build the track extra
203  reco::TrackExtra * extra = makeTrackExtra(recotrack_orig,recotrackref,*TEoutput);
204  if (!extra){
205  edm::LogError( theCategory)<<"cannot create the track extra for this track.";
206  //pop the inserted track
207  Toutput->pop_back();
208  continue;}
209 
210  //attach the collection of rechits
211  if (!attachRecHits(recotrack_orig,recotrackref,*extra,*TRHoutput)){
212  edm::LogError( theCategory)<<"cannot attach any rechits on this track";
213  //pop the inserted track
214  Toutput->pop_back();
215  //pop the track extra
216  TEoutput->pop_back();
217  theTEi--;
218  continue;}
219 
220  }//loop over the original reco tracks
221 
222  LogDebug( theCategory)<<"writing: "<<Toutput->size()<<" corrected reco::Track to the event.";
223 
224  iEvent.put(Toutput, theInstanceName);
225  iEvent.put(TEoutput);
226  iEvent.put(TRHoutput);
227 
228 }
229 
230 // ------------ method called once each job just before starting event loop ------------
231 void
233 {
235 }
236 
237 // ------------ method called once each job just after ending the event loop ------------
238 void
241 }
#define LogDebug(id)
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
edm::RefProd< reco::TrackExtraCollection > theRefprodTE
get reference before put track extra to the event, in order to create edm::Ref
const Vector & momentum() const
track momentum vector
Definition: TrackBase.h:148
virtual void beginJob()
framework method
const TrackExtraRef & extra() const
reference to &quot;extra&quot; object
Definition: Track.h:97
reco::TrackBase::CovarianceMatrix fix_cov_matrix(const reco::TrackBase::CovarianceMatrix &error_matrix, const GlobalVector &momentum)
return a corrected error matrix
virtual void produce(edm::Event &, const edm::EventSetup &)
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:10
bool attachRecHits(const reco::Track &recotrack_orig, reco::Track &recotrack, reco::TrackExtra &trackextra, TrackingRecHitCollection &RHcol)
attached rechits to the newly created reco::Track and reco::TrackExtra
CurvilinearTrajectoryError get(GlobalVector momentum, bool convolute=true)
main method to be used. Retrieve a 5x5 symetrical matrix according to parametrization of error or sca...
std::string theInstanceName
instrance name of the created track collecion. rechit and trackextra have no instance name ...
edm::InputTag theTrackLabel
input tag of the reco::Track collection to be corrected
const Point & referencePoint() const
Reference point on the track.
Definition: TrackBase.h:151
double charge(const std::vector< uint8_t > &Ampls)
CovarianceMatrix covariance() const
return track covariance matrix
Definition: TrackBase.h:180
std::string theCategory
log category: MuonErrorMatrixAdjuster
void multiply(reco::TrackBase::CovarianceMatrix &revised_matrix, const reco::TrackBase::CovarianceMatrix &scale_matrix)
mutliply revised_matrix (first argument) by second matrix TERM by TERM
edm::ESHandle< MagneticField > theField
hold on to the magnetic field
MuonErrorMatrixAdjuster(const edm::ParameterSet &)
constructor
void push_back(D *&d)
Definition: OwnVector.h:274
MuonErrorMatrix * theMatrixProvider
int iEvent
Definition: GenABIO.cc:230
bool isNull() const
Checks for null.
Definition: Ref.h:247
T mag() const
Definition: PV3DBase.h:67
double chi2() const
chi-squared of the fit
Definition: TrackBase.h:105
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:116
bool divide(reco::TrackBase::CovarianceMatrix &num_matrix, const reco::TrackBase::CovarianceMatrix &denom_matrix)
divide the num_matrix (first argument) by second matrix, TERM by TERM
double ndof() const
number of degrees of freedom of the fit
Definition: TrackBase.h:107
edm::Ref< TrackingRecHitCollection > TrackingRecHitRef
persistent reference to a TrackingRecHit
int j
Definition: DBlmapReader.cc:9
edm::Ref< reco::TrackExtraCollection >::key_type theTEi
edm::RefProd< TrackingRecHitCollection > theRefprodRH
get reference before put rechit to the event, in order to create edm::Ref
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:63
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
GlobalVector momentum() const
RefProd< PROD > getRefBeforePut()
Definition: Event.h:128
edm::ParameterSet theMatrixProvider_pset
holds the error matrix parametrization
GlobalPoint position() const
void setHitPattern(const C &c)
set hit patterns from vector of hit references
Definition: TrackBase.h:244
std::vector< TrackExtra > TrackExtraCollection
collection of TrackExtra objects
Definition: TrackExtraFwd.h:9
edm::OwnVector< TrackingRecHit > TrackingRecHitCollection
collection of TrackingRecHits
void setExtra(const TrackExtraRef &ref)
set reference to &quot;extra&quot; object
Definition: Track.h:95
reco::TrackExtra * makeTrackExtra(const reco::Track &recotrack_orig, reco::Track &recotrack, reco::TrackExtraCollection &TEcol)
create a track extra for the newly created recotrack, scaling the outer/inner measurment error matrix...
tuple tracks
Definition: testEve_cfg.py:39
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
FreeTrajectoryState initialFreeState(const reco::Track &tk, const MagneticField *field)
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:62
void add(const TrackingRecHitRef &r)
add a reference to a RecHit
bool selectTrack(const reco::Track &recotrack_orig)
make a selection on the reco:Track. (dummy for the moment)
reco::Track makeTrack(const reco::Track &recotrack_orig, const FreeTrajectoryState &PCAstate)
create a corrected reco::Track from itself and trajectory state (redundant information) ...
bool theRescale
select the rescaling or replacing method to correct the error matrix
int charge() const
track electric charge
Definition: TrackBase.h:111
math::Error< dimension >::type CovarianceMatrix
5 parameter covariance matrix
Definition: TrackBase.h:70
edm::Ref< TrackingRecHitCollection >::key_type theRHi
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:65