ATLAS Offline Software
Loading...
Searching...
No Matches
ConstrainedTrackProvider.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3*/
4
5
6#include "TrkFitterInterfaces/IGlobalTrackFitter.h"
7#include "TrkFitterInterfaces/ITrackFitter.h"
8#include "TrkAlignGenTools/ConstrainedTrackProvider.h"
9#include "TrkPseudoMeasurementOnTrack/PseudoMeasurementOnTrack.h"
10#include "TrkSurfaces/PerigeeSurface.h"
11#include "TrkParameters/TrackParameters.h"
12#include "TrkMeasurementBase/MeasurementBase.h"
13
14
15#include "TFile.h"
16#include "TH2F.h"
17
18
19namespace Trk {
20
21 //________________________________________________________________________
22 ConstrainedTrackProvider::ConstrainedTrackProvider(const std::string& type,
23 const std::string& name,
24 const IInterface* parent)
25
26 : AthAlgTool(type,name,parent)
27 , m_trackFitter("Trk::GlobalChi2Fitter/InDetTrackFitter")
28 , m_scalepmaptogev(false)
29 , m_doTrackSelection (true)
30 , m_constrainedTracks(0)
31 , m_passconstrainedRequirements(0)
32 , m_unconstrainedTracks(0)
33 , m_inputTrackCollection("Tracks")
34 , m_constraintInputFile_P(nullptr)
35 , m_etaphiMap_P(nullptr)
36 , m_constraintInputFile_d0(nullptr)
37 , m_etaphiMap_d0(nullptr)
38 , m_constraintInputFile_z0(nullptr)
39 , m_etaphiMap_z0(nullptr)
40 , m_CorrectMeanD0()
41 , m_SelectByCharge()
42 , m_SelectPositive(false)
43 {
44 declareInterface<ITrackCollectionProvider>(this);
45 declareProperty("doTrackSelection", m_doTrackSelection );
46 declareProperty("TrackFitter", m_trackFitter );
47 declareProperty("InputTracksCollection", m_inputTrackCollection );
48 declareProperty("RunOutlierRemoval", m_runOutlierRemoval = true );
49 declareProperty("MaxRetrievalErrors", m_maxRetrievalErrors = 10 );
50 declareProperty("UseConstrainedTrkOnly", m_useConstrainedTrkOnly = false );
51 declareProperty("MinPt", m_minPt = 15.0 );
52 declareProperty("MaxPt", m_maxPt = 100.0 );
53 declareProperty("SelectByCharge", m_SelectByCharge = false );
54 declareProperty("SelectPositive", m_SelectPositive = true );
55 declareProperty("MinPIXHits", m_minPIXHits = 0 ); //1
56 declareProperty("MinSCTHits", m_minSCTHits = 0 ); //6
57 declareProperty("MinTRTHits", m_minTRTHits = 0 );
58 declareProperty("Maxd0", m_maxd0 = 500. );
59 declareProperty("Maxz0", m_maxz0 = 500. );
60 declareProperty("CorrectMomentum", m_CorrectMomentum = true );
61 declareProperty("MomentumConstraintFileName", m_constraintFileName_P = "Constraint.root" );
62 declareProperty("MomentumConstraintHistName", m_constraintHistName_P = "EtaPhiMap" );
63 declareProperty("CorrectD0", m_CorrectD0 = false );
64 declareProperty("CorrectMeanD0", m_CorrectMeanD0 = false );
65 declareProperty("d0ConstraintFileName", m_constraintFileName_d0 = "Constraint.root" );
66 declareProperty("d0ConstraintHistName", m_constraintHistName_d0 = "EtaPhiMap" );
67 declareProperty("CorrectZ0", m_CorrectZ0 = false );
68 declareProperty("z0ConstraintFileName", m_constraintFileName_z0 = "Constraint.root" );
69 declareProperty("z0ConstraintHistName", m_constraintHistName_z0 = "EtaPhiMap" );
70 declareProperty("UseConstraintError", m_useConstraintError = true ,"Bla bla " );
71 declareProperty("ReduceConstraintUncertainty",m_reduceConstraintUncertainty = 1., "Reduce the uncertainty on teh track parmater constraint by this amount" );
72 declareProperty("ReduceConstraintUncert_z0" ,m_reduceConstraintUncert_z0 = 1., "Reduce the uncertainty on z0 track parameter constraint by this amount" );
73 declareProperty("DeltaScaling", m_deltaScaling = 1.);
74 declareProperty("ScalePMapToGeV" ,m_scalepmaptogev);
75 }
76
77 //________________________________________________________________________
78 ConstrainedTrackProvider::~ConstrainedTrackProvider()
79 = default;
80
81 //________________________________________________________________________
82 StatusCode ConstrainedTrackProvider::initialize()
83 {
84 // configure main track fitter
85 if(m_trackFitter.retrieve().isFailure()) {
86 msg(MSG::FATAL) << "Could not get " << m_trackFitter << endmsg;
87 return StatusCode::FAILURE;
88 }
89 ATH_MSG_INFO("Retrieved " << m_trackFitter);
90
91
92 ATH_MSG_INFO( "init m_useConstraintError " <<m_useConstraintError );
93
94 ATH_MSG_INFO( "Init m_reduceConstraintUncertainty " << m_reduceConstraintUncertainty );
95
96 ATH_MSG_INFO( "Init m_reduceConstraintUncert_z0 " << m_reduceConstraintUncert_z0 );
97
98 if(m_CorrectMomentum){
99
100 m_constraintInputFile_P = new TFile(m_constraintFileName_P.c_str() ,"read");
101 if ( m_constraintInputFile_P->IsZombie() || !(m_constraintInputFile_P->IsOpen()) ) {
102 ATH_MSG_FATAL( " Problem reading TFile " << m_constraintFileName_P );
103 return StatusCode::FAILURE;
104 }
105 ATH_MSG_INFO("Opened file containing the deltaSagitta constraints" << m_constraintFileName_P);
106 m_etaphiMap_P = (TH2F*)m_constraintInputFile_P->Get(m_constraintHistName_P.c_str());
107 if(!m_etaphiMap_P){
108 ATH_MSG_FATAL( " Problem getting constraints Hist. Name " << m_constraintHistName_P );
109 m_constraintInputFile_P->Close();
110 delete m_constraintInputFile_P;
111 return StatusCode::FAILURE;
112 }
113 ATH_MSG_INFO("Opened constraints histogram " << m_constraintHistName_P);
114
115 }
116
117 if(m_CorrectD0){
118
119 m_constraintInputFile_d0 = new TFile(m_constraintFileName_d0.c_str() ,"read");
120 if ( m_constraintInputFile_d0->IsZombie() || !(m_constraintInputFile_d0->IsOpen()) ) {
121 ATH_MSG_FATAL( " Problem reading TFile " << m_constraintFileName_d0 );
122 return StatusCode::FAILURE;
123 }
124 ATH_MSG_INFO("Opened file containing the d0 constraints" << m_constraintFileName_d0);
125 m_etaphiMap_d0 = (TH2F*)m_constraintInputFile_d0->Get(m_constraintHistName_d0.c_str());
126 if(!m_etaphiMap_d0){
127 ATH_MSG_FATAL( " Problem getting constraints Hist. Name " << m_constraintHistName_d0 );
128 m_constraintInputFile_d0->Close();
129 delete m_constraintInputFile_d0;
130 return StatusCode::FAILURE;
131 }
132 ATH_MSG_INFO("Opened constraints histogram " << m_constraintHistName_d0);
133
134 }
135
136 if(m_CorrectZ0){
137
138 m_constraintInputFile_z0 = new TFile(m_constraintFileName_z0.c_str() ,"read");
139 if ( m_constraintInputFile_z0->IsZombie() || !(m_constraintInputFile_z0->IsOpen()) ) {
140 ATH_MSG_FATAL( " Problem reading TFile " << m_constraintFileName_z0 );
141 return StatusCode::FAILURE;
142 }
143 ATH_MSG_INFO("Opened file containing the z0 constraints" << m_constraintFileName_z0);
144 m_etaphiMap_z0 = (TH2F*)m_constraintInputFile_z0->Get(m_constraintHistName_z0.c_str());
145 if(!m_etaphiMap_z0){
146 ATH_MSG_FATAL( " Problem getting constraints Hist. Name " << m_constraintHistName_z0 );
147 m_constraintInputFile_z0->Close();
148 delete m_constraintInputFile_z0;
149 return StatusCode::FAILURE;
150 }
151 ATH_MSG_INFO("Opened constraints histogram " << m_constraintHistName_z0);
152
153 }
154
155 /* Commented to enable the scaling of the uncertainty
156 *
157 if(m_reduceConstraintUncertainty<=0){
158 ATH_MSG_FATAL( " reduceConstraintUncertainty is <= 0 it must be >0 . Currently " << m_reduceConstraintUncertainty );
159 return StatusCode::FAILURE;
160 }
161 */
162
163 return StatusCode::SUCCESS;
164 }
165
166
167
168
169
170 //________________________________________________________________________
171 StatusCode ConstrainedTrackProvider::finalize()
172 {
173 if(m_constraintInputFile_P){
174 m_constraintInputFile_P->Close();
175 delete m_constraintInputFile_P;
176 }
177
178 if(m_constraintInputFile_d0){
179 m_constraintInputFile_d0->Close();
180 delete m_constraintInputFile_d0;
181 }
182
183 if(m_constraintInputFile_z0){
184 m_constraintInputFile_z0->Close();
185 delete m_constraintInputFile_z0;
186 }
187
188 return StatusCode::SUCCESS;
189 }
190
191
192 StatusCode ConstrainedTrackProvider::trackCollection(const TrackCollection*& finalTracks)
193 {
194
195
196 const TrackCollection* originalTracks = nullptr;
197
198 if ( StatusCode::SUCCESS != evtStore()->retrieve(originalTracks, m_inputTrackCollection) ){
199 ATH_MSG_WARNING(" Can't retrieve " << m_inputTrackCollection << " from the StoreGate ");
200 finalTracks = originalTracks; // return empty collection
201 return StatusCode::SUCCESS;
202 }
203
204 TrackCollection* trackCollection = new TrackCollection;
205
206 ATH_MSG_DEBUG("Input track collection size: " << originalTracks->size());
207
208 TrackCollection::const_iterator trackIt = originalTracks->begin();
209 TrackCollection::const_iterator trackItE = originalTracks->end();
210 int trackCount = 0;
211
212 for ( ; trackIt != trackItE; ++trackIt ) {
213 trackCount++;
214 ATH_MSG_DEBUG("Dealing with track: " << trackCount << " / " << originalTracks->size());
215
216 bool acceptedTrack = true;
217 if (m_doTrackSelection) acceptedTrack = passTrackSelection( *trackIt );
218
219 if( acceptedTrack ){
220 ATH_MSG_DEBUG("Track selection OK for track: " << trackCount);
221 const Trk::Perigee* measuredPerigee = (*trackIt)->perigeeParameters();
222 if(!measuredPerigee){
223 ATH_MSG_DEBUG(" no measuredPerigee");
224 continue;
225 }
226
227 if(!(measuredPerigee->covariance()) ){
228 ATH_MSG_DEBUG(" no measuredPerigee");
229 continue;
230 }
231
232 const Trk::Surface* surf = dynamic_cast<const Trk::Surface*>( &(measuredPerigee->associatedSurface()) ) ;
233 if( !surf ) {
234 ATH_MSG_DEBUG("NO surface of the measuredPerigee");
235 continue;
236 }
237
238 // Get corrections from histograms if they are requested
239
240 double correctedQoverP(0), correctedQoverPError(0);
241 if(m_CorrectMomentum){
242 getCorrectedValues_P( measuredPerigee, correctedQoverP, correctedQoverPError);
243 } else {
244 correctedQoverP = measuredPerigee->parameters()[Trk::qOverP];
245 correctedQoverPError = 1e6;
246 }
247
248 double correctedD0(0), correctedD0Error(0);
249 if(m_CorrectD0){
250 getCorrectedValues_d0( measuredPerigee, correctedD0, correctedD0Error);
251 } else {
252 correctedD0 = measuredPerigee->parameters()[Trk::d0];
253 correctedD0Error = 1e6;
254 }
255
256 double correctedZ0(0), correctedZ0Error(0);
257 if(m_CorrectZ0){
258 getCorrectedValues_z0( measuredPerigee, correctedZ0, correctedZ0Error);
259 } else {
260 correctedZ0 = measuredPerigee->parameters()[Trk::z0];
261 correctedZ0Error = 1e6;
262 }
263
264 // Construct pseudomeasurement
265 Trk::DefinedParameter constrainedD0( correctedD0 , Trk::d0) ;
266 Trk::DefinedParameter constrainedZ0( correctedZ0 , Trk::z0) ;
267 Trk::DefinedParameter constrainedqOverP( correctedQoverP , Trk::qOverP) ;
268
269 std::vector<Trk::DefinedParameter> constrainedPars;
270 constrainedPars.push_back( constrainedD0 ) ;
271 constrainedPars.push_back( constrainedZ0 ) ;
272 constrainedPars.push_back( constrainedqOverP ) ;
273
274 Amg::MatrixX constrainedCov( 3,3 ) ;
275 constrainedCov.setZero();
276 constrainedCov( 0, 0 ) = correctedD0Error;
277 constrainedCov( 1, 1 ) = correctedZ0Error;
278 constrainedCov( 2, 2 ) = correctedQoverPError;
279 // constrainedCov( 1, 0 ) = 0; //initially we will assume the constraints are independant
280 // constrainedCov( 2, 0 ) = 0;
281 // constrainedCov( 2, 1 ) = 0;
282
283 Trk::PseudoMeasurementOnTrack *pmot = new Trk::PseudoMeasurementOnTrack( Trk::LocalParameters( constrainedPars ),
284 std::move(constrainedCov), *surf) ;
285
286
287 // Add pseusdo measurement to a list of hits to be fit;
288 std::vector<const Trk::MeasurementBase*> vecOfMB;
289 if (pmot) vecOfMB.push_back(pmot);
290
291 // Add the remaining hits;
292 DataVector<const Trk::MeasurementBase>::const_iterator it = (*trackIt)->measurementsOnTrack()->begin();
293 DataVector<const Trk::MeasurementBase>::const_iterator itend = (*trackIt)->measurementsOnTrack()->end();
294 for (;it!=itend;++it)
295 vecOfMB.push_back(*it);
296
297 if(m_runOutlierRemoval){
298 it = (*trackIt)->outliersOnTrack()->begin();
299 itend = (*trackIt)->outliersOnTrack()->end();
300 // Add the remaining hits;
301 for (;it!=itend;++it)
302 vecOfMB.push_back(*it);
303 }
304
305 Trk::Track* constrainedFittedTrack = m_trackFitter->fit( Gaudi::Hive::currentContext(),
306 vecOfMB, *measuredPerigee,
307 m_runOutlierRemoval, Trk::pion).release();
308 delete pmot;
309
310 if (constrainedFittedTrack){
311 ATH_MSG_DEBUG("Fit was successful");
312 const Trk::Perigee* constrainedPerigee = constrainedFittedTrack->perigeeParameters();
313 if(!constrainedPerigee){
314 ATH_MSG_DEBUG(" no constrainedmeasuredPerigee");
315 } else {
316 ATH_MSG_DEBUG("Initial track parameters --> momentum: " << 1/measuredPerigee->parameters()[Trk::qOverP] * 1e-3
317 << " pt: " << 1/measuredPerigee->parameters()[Trk::qOverP] * 1e-3 * sin(measuredPerigee->parameters()[Trk::theta])
318 << " d0: " << measuredPerigee->parameters()[Trk::d0]
319 << " z0: " << measuredPerigee->parameters()[Trk::z0] );
320 ATH_MSG_DEBUG("constraint values --> momentum: " << 1/correctedQoverP * 1e-3
321 << " pt: " << 1/correctedQoverP * 1e-3 * sin(measuredPerigee->parameters()[Trk::theta])
322 << " d0: " << correctedD0
323 << " z0: " << correctedZ0 );
324 ATH_MSG_DEBUG("Constrained fit result --> momentum: " << 1/constrainedPerigee->parameters()[Trk::qOverP] * 1e-3
325 << " pt: " << 1/constrainedPerigee->parameters()[Trk::qOverP] * 1e-3 * sin(constrainedPerigee->parameters()[Trk::theta])
326 << " d0: " << constrainedPerigee->parameters()[Trk::d0]
327 << " z0: " << constrainedPerigee->parameters()[Trk::z0] );
328 }
329
330 const Track* track = *trackIt;
331 constrainedFittedTrack->setTrackSummary( std::make_unique<Trk::TrackSummary> (*track->trackSummary()) );
332
333 trackCollection->push_back(constrainedFittedTrack);
334 ++m_constrainedTracks;
335 } else{
336 ATH_MSG_DEBUG("Constrained fit was unsuccessful");
337 }
338 } //--> if ( acceptedTrack)
339 else {
340
341 if(m_useConstrainedTrkOnly) continue;
342
343 //Prob could just copy the track?
344 Trk::Track* unconstrainedFittedTrack =
345 m_trackFitter
346 ->fit(Gaudi::Hive::currentContext(),
347 **trackIt,
348 m_runOutlierRemoval,
349 Trk::pion)
350 .release();
351
352 if(unconstrainedFittedTrack) {
353 ATH_MSG_DEBUG("Unconstrained fit was successful");
354 trackCollection->push_back(unconstrainedFittedTrack);
355 ++m_unconstrainedTracks;
356 } else {
357 ATH_MSG_DEBUG("Unconstrained fit was unsuccessful");
358 }
359 }
360 }
361
362
363 if (StatusCode::SUCCESS != evtStore()->record(trackCollection, "AlignmentConstrainedTracks")){
364 ATH_MSG_WARNING("Problem with recording AlignmentConstrainedTracks to StoreGate!");
365 delete trackCollection;
366 return StatusCode::SUCCESS;
367 }
368
369 ATH_MSG_INFO ("Size of the trackCollection in this event : " << trackCollection->size() );
370 //track Collection cannot be null here; it has already been dereferenced
371 finalTracks = trackCollection;
372
373 return StatusCode::SUCCESS;
374
375 }
376
377
378
379
380 bool ConstrainedTrackProvider::passTrackSelection(const Trk::Track *track)
381 {
382 ATH_MSG_DEBUG("== TrackSelection == STARTED ==");
383
384 const Trk::Perigee* perigee = track->perigeeParameters();
385 if(!perigee){
386 ATH_MSG_DEBUG("== TrackSelection == NO perigee on this track");
387 return false;
388 }
389
390 const Trk::TrackSummary* summary = track->trackSummary();
391
392 if(!summary){
393 ATH_MSG_DEBUG("== TrackSelection == No summary on this track");
394 return false;
395 }
396
397 const int nPixHits = summary->get(Trk::numberOfPixelHits);
398 const int nSCTHits = summary->get(Trk::numberOfSCTHits);
399 const int nTRTHits = summary->get(Trk::numberOfTRTHits);
400
401 const double qoverP = perigee->parameters()[Trk::qOverP] * 1000.;
402 const double z0 = perigee->parameters()[Trk::z0];
403 const double d0 = perigee->parameters()[Trk::d0];
404 double pt = 0.;
405 if ( qoverP != 0 ) pt = fabs(1.0/qoverP)*sin(perigee->parameters()[Trk::theta]);
406
407 ATH_MSG_DEBUG( "== TrackSelection == track pt : "<< pt );
408 double charge(0);
409 if( qoverP > 0)
410 charge = 1.;
411 else
412 charge = -1.;
413 ATH_MSG_DEBUG( " pt : "<< pt << " ; Charge : "<< charge);
414
415 if(pt < m_minPt || pt > m_maxPt ){
416 ATH_MSG_DEBUG("== TrackSelection == Track fails pt cut [ "<<m_minPt<<", "<< m_maxPt << "] " << pt );
417 return false;
418 }
419
420 if (m_SelectByCharge) {
421 if (m_SelectPositive){
422 if (charge<0){
423 ATH_MSG_DEBUG("Track fails charge cut: negative charge " << charge );
424 return false;
425 }
426 }
427 if (not m_SelectPositive){
428 if (charge>0) {
429 ATH_MSG_DEBUG("Track fails charge cut: positive charge " << charge);
430 return false;
431 }
432 }
433 }
434
435 ATH_MSG_VERBOSE("== TrackSelection == Track N PIX hits: " << nPixHits << " requested: "<< m_minPIXHits );
436 ATH_MSG_VERBOSE("== TrackSelection == Track N SCT hits: " << nSCTHits << " requested: "<< m_minSCTHits );
437 ATH_MSG_VERBOSE("== TrackSelection == Track N TRT hits: " << nTRTHits << " requested: "<< m_minTRTHits );
438 ATH_MSG_VERBOSE("== TrackSelection == Track d0: " << d0 << " requested: "<< m_maxd0 );
439 ATH_MSG_VERBOSE("== TrackSelection == Track z0: " << z0 << " requested: "<< m_maxz0 );
440
441 if( nPixHits < m_minPIXHits ||
442 nSCTHits < m_minSCTHits ||
443 nTRTHits < m_minTRTHits ||
444 d0 > m_maxd0 ||
445 z0 > m_maxz0 ) {
446 ATH_MSG_DEBUG("This track did not pass cuts --- nPixHits: " << nPixHits << " nSCTHits: " << nSCTHits <<
447 " nTRTHits: " << nTRTHits << " idd0atIP: " << d0 << " idz0atIP: " << z0 );
448 return false;
449 }
450
451 return true;
452 }
453
454 void ConstrainedTrackProvider::getCorrectedValues_P(const Trk::Perigee* measuredPerigee, double& correctedQoverP,double& correctedQoverPError)
455 {
456
457 // scale q/p not p
458
459 double charge(0);
460 if( measuredPerigee->parameters()[Trk::qOverP] > 0)
461 charge = 1.;
462 else
463 charge = -1.;
464
465 double pt = fabs(1./(measuredPerigee->parameters()[Trk::qOverP]))* sin(measuredPerigee->parameters()[Trk::theta]) *1e-3;
466 double eta = -log(tan(measuredPerigee->parameters()[Trk::theta]/2.));
467 double phi = measuredPerigee->parameters()[Trk::phi];
468 double perr = (*measuredPerigee->covariance())( Trk::qOverP, Trk::qOverP );
469
470
471 int binNumber = m_etaphiMap_P->FindBin(eta, phi);
472 double constraintErr = m_etaphiMap_P->GetBinError(binNumber);
473
474 double delta = m_etaphiMap_P->GetBinContent(binNumber) * m_deltaScaling;
475
476 if (m_scalepmaptogev)
477 delta = delta * 0.001;
478
479 correctedQoverP = measuredPerigee->parameters()[Trk::qOverP] * (1.+ charge *pt *delta );
480 correctedQoverPError = perr;
481
482 if(m_useConstraintError){
483 double qoverpFracError = perr/pow(measuredPerigee->parameters()[Trk::qOverP], 2);
484 qoverpFracError += pow(constraintErr,2);
485 correctedQoverPError = qoverpFracError*correctedQoverP*correctedQoverP;
486 }
487 ATH_MSG_DEBUG("Scaling by 1/m_reduceConstraintUncertainty " << m_reduceConstraintUncertainty << '\t'<< pow( m_reduceConstraintUncertainty,-2)) ;
488 correctedQoverPError = correctedQoverPError * pow( m_reduceConstraintUncertainty,-2);
489 ATH_MSG_DEBUG(" == input pt: " << pt << " deltaSagitta: " << delta << " final pt: " << sin(measuredPerigee->parameters()[Trk::theta]) * 1e-3/correctedQoverP);
490 }
491
492
493 void ConstrainedTrackProvider::getCorrectedValues_d0(const Trk::Perigee* measuredPerigee, double& correctedD0, double& correctedD0Error)
494 {
495 ATH_MSG_DEBUG(" == getCorrectedValues_d0 == STARTED ==");
496 // scale d0
497
498 double charge(0);
499 if( measuredPerigee->parameters()[Trk::qOverP] > 0)
500 charge = 1.;
501 else
502 charge = -1.;
503
504 double d0 = measuredPerigee->parameters()[Trk::d0];
505 double eta = -log(tan(measuredPerigee->parameters()[Trk::theta]/2.));
506 double phi = measuredPerigee->parameters()[Trk::phi];
507 double d0err = (*measuredPerigee->covariance())( Trk::d0, Trk::d0 );
508
509
510 int binNumber = m_etaphiMap_d0->FindBin(eta, phi);
511 double constraintErr = m_etaphiMap_d0->GetBinError(binNumber);
512 double delta = m_etaphiMap_d0->GetBinContent(binNumber) * m_deltaScaling;
513
514 correctedD0 = d0 - charge * delta * 0.5 ; // delta d0 = d0_pos - d0_neg . To make the difference null, substract half to the pos and add half to the neg track.
515 if (m_CorrectMeanD0) {correctedD0 = d0 + delta;}
516 correctedD0Error = d0err;
517 ATH_MSG_DEBUG("Correcting d0 by delta = " << delta << "\t, corrected d0 = "<< correctedD0 << "\t, original d0 = "<< d0 ) ;
518
519 if(m_useConstraintError){
520 correctedD0Error = d0err + pow( constraintErr, 2 );
521 }
522 ATH_MSG_DEBUG("Scaling by 1/m_reduceConstraintUncertainty " << m_reduceConstraintUncertainty << '\t'<< pow( m_reduceConstraintUncertainty,-2)) ;
523 correctedD0Error = correctedD0Error * pow( m_reduceConstraintUncertainty,-2);
524
525 ATH_MSG_DEBUG(" == input d0: " << d0 << " deltad0: " << delta << " final d0: " << correctedD0 << " +- " << correctedD0Error);
526 ATH_MSG_DEBUG(" == getCorrectedValues_d0 == Completed ==");
527 }
528
529 void ConstrainedTrackProvider::getCorrectedValues_z0(const Trk::Perigee* measuredPerigee, double& correctedZ0,double& correctedZ0Error)
530 {
531
532 // scale z0
533
534 double charge(0);
535 if( measuredPerigee->parameters()[Trk::qOverP] > 0)
536 charge = 1.;
537 else
538 charge = -1.;
539
540 double z0 = measuredPerigee->parameters()[Trk::z0];
541 double eta = -log(tan(measuredPerigee->parameters()[Trk::theta]/2.));
542 double phi = measuredPerigee->parameters()[Trk::phi];
543 double z0err = (*measuredPerigee->covariance())( Trk::z0, Trk::z0 );
544
545 int binNumber = m_etaphiMap_z0->FindBin(eta, phi);
546 double constraintErr = m_etaphiMap_z0->GetBinError(binNumber);
547 double delta = m_etaphiMap_z0->GetBinContent(binNumber) * m_deltaScaling;
548
549 correctedZ0 = z0 - charge * delta * 0.5 ; // delta z0 = z0_pos - z0_neg . To make the difference null, substract half to the pos and add half to the neg track.
550 correctedZ0Error = z0err;
551
552 if(m_useConstraintError){
553 correctedZ0Error = z0err + pow( constraintErr, 2 );
554 }
555 ATH_MSG_DEBUG("Scaling by 1/m_reduceConstraintUncertainty " << m_reduceConstraintUncertainty << '\t'<< pow( m_reduceConstraintUncertainty,-2)) ;
556 correctedZ0Error = correctedZ0Error * pow( m_reduceConstraintUncertainty,-2);
557 ATH_MSG_DEBUG("Scaling by 1/m_reduceConstraintUncert_z0 " << m_reduceConstraintUncert_z0 << '\t'<< pow( m_reduceConstraintUncert_z0,-2)) ;
558 correctedZ0Error = correctedZ0Error * pow( m_reduceConstraintUncert_z0,-2);
559 ATH_MSG_DEBUG(" == input z0: " << z0 << " deltaz0: " << delta << " final z0: " << correctedZ0 << " +- " << correctedZ0Error);
560 }
561
562
563
564
565
566 void ConstrainedTrackProvider::printSummary(){
567
568 if(m_logStream) {
569
570 *m_logStream<<"*************************************************************"<<std::endl;
571 *m_logStream<<"****** ConstrainedTrackProvider Summary ******"<<std::endl;
572 *m_logStream<<"*"<<std::endl;
573 *m_logStream<<"* number of combined muons failed in refit: " << m_passconstrainedRequirements << std::endl;
574 *m_logStream<<"* number of combined muons failed in refit: " << m_constrainedTracks << std::endl;
575 *m_logStream<<"* number of combined muons succeeded in refit: " << m_unconstrainedTracks << std::endl;
576
577 *m_logStream<<"*"<<std::endl;
578 }
579 }
580
581
582
583
584
585
586} // end namespace
eta
Scalar eta() const
pseudorapidity method
Definition AmgMatrixBasePlugin.h:83
endmsg
#define endmsg
Definition AnalysisConfig_Ntuple.cxx:63
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition AthMsgStreamMacros.h:34
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition AthMsgStreamMacros.h:31
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition AthMsgStreamMacros.h:28
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition AthMsgStreamMacros.h:32
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition AthMsgStreamMacros.h:29
charge
double charge(const T &p)
Definition AtlasPID.h:997
TrackCollection
DataVector< Trk::Track > TrackCollection
This typedef represents a collection of Trk::Track objects.
Definition TrackCollection.h:19
pow
constexpr int pow(int base, int exp) noexcept
Definition ap_fixedTest.cxx:15
AthAlgTool::AthAlgTool
AthAlgTool(const std::string &type, const std::string &name, const IInterface *parent)
Constructor with parameters:
Definition AthAlgTool.cxx:16
AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
Definition AthCommonDataStore.h:145
AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore
ServiceHandle< StoreGateSvc > & evtStore()
Definition AthCommonDataStore.h:85
AthCommonMsg< AlgTool >::msg
MsgStream & msg() const
Definition AthCommonMsg.h:24
DataVector< Trk::Track >::const_iterator
DataModel_detail::const_iterator< DataVector > const_iterator
Definition DataVector.h:838
DataVector::end
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
DataVector::begin
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
DataVector::size
size_type size() const noexcept
Returns the number of elements in the collection.
Trk::ConstrainedTrackProvider::getCorrectedValues_d0
void getCorrectedValues_d0(const Trk::Perigee *mp, double &corrected_d0, double &corrected_d0Error)
Definition ConstrainedTrackProvider.cxx:493
Trk::ConstrainedTrackProvider::getCorrectedValues_P
void getCorrectedValues_P(const Trk::Perigee *mp, double &correctedQoverP, double &correctedQoverPError)
Definition ConstrainedTrackProvider.cxx:454
Trk::ConstrainedTrackProvider::getCorrectedValues_z0
void getCorrectedValues_z0(const Trk::Perigee *mp, double &corrected_z0, double &corrected_z0Error)
Definition ConstrainedTrackProvider.cxx:529
Trk::ConstrainedTrackProvider::passTrackSelection
bool passTrackSelection(const Trk::Track *track)
Definition ConstrainedTrackProvider.cxx:380
Trk::ConstrainedTrackProvider::printSummary
virtual void printSummary()
Print statistical summary to logfile.
Definition ConstrainedTrackProvider.cxx:566
Trk::ConstrainedTrackProvider::m_maxRetrievalErrors
int m_maxRetrievalErrors
maximum allowed number of retrieval errors at the beginning of the job (-1 for infinite)
Definition ConstrainedTrackProvider.h:56
Trk::ConstrainedTrackProvider::m_trackFitter
ToolHandle< IGlobalTrackFitter > m_trackFitter
normal track fitter
Definition ConstrainedTrackProvider.h:51
Trk::ConstrainedTrackProvider::ConstrainedTrackProvider
ConstrainedTrackProvider(const std::string &type, const std::string &name, const IInterface *parent)
Definition ConstrainedTrackProvider.cxx:22
Trk::ConstrainedTrackProvider::m_runOutlierRemoval
RunOutlierRemoval m_runOutlierRemoval
run outlier removal in the GX2 fitter
Definition ConstrainedTrackProvider.h:53
Trk::ConstrainedTrackProvider::trackCollection
virtual StatusCode trackCollection(const TrackCollection *&tracks)
Definition ConstrainedTrackProvider.cxx:192
Trk::ITrackCollectionProvider::m_logStream
std::ostream * m_logStream
logfile output stream
Definition ITrackCollectionProvider.h:53
Trk::ParametersT::associatedSurface
virtual const S & associatedSurface() const override final
Access to the Surface method.
Trk::PseudoMeasurementOnTrack
Class to handle pseudo-measurements in fitters and on track objects.
Definition PseudoMeasurementOnTrack.h:44
Trk::Surface
Abstract Base Class for tracking surfaces.
Definition Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/Surface.h:79
Trk::TrackSummary
A summary of the information contained by a track.
Definition Tracking/TrkEvent/TrkTrackSummary/TrkTrackSummary/TrackSummary.h:287
Trk::Track::setTrackSummary
void setTrackSummary(std::unique_ptr< Trk::TrackSummary > input)
Set the track summary.
Trk::Track::perigeeParameters
const Perigee * perigeeParameters() const
return Perigee.
Definition Tracking/TrkEvent/TrkTrack/src/Track.cxx:163
Amg::MatrixX
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Definition EventPrimitives.h:27
Trk
Ensure that the ATLAS eigen extensions are properly loaded.
Definition FakeTrackBuilder.h:9
Trk::Perigee
ParametersT< TrackParametersDim, Charged, PerigeeSurface > Perigee
Definition Tracking/TrkEvent/TrkParameters/TrkParameters/TrackParameters.h:33
Trk::theta
@ theta
Definition ParamDefs.h:66
Trk::qOverP
@ qOverP
perigee
Definition ParamDefs.h:67
Trk::phi
@ phi
Definition ParamDefs.h:75
Trk::d0
@ d0
Definition ParamDefs.h:63
Trk::z0
@ z0
Definition ParamDefs.h:64
Trk::DefinedParameter
std::pair< double, ParamDefs > DefinedParameter
Typedef to of a std::pair<double, ParamDefs> to identify a passed-through double as a specific type o...
Definition DefinedParameter.h:27
Trk::numberOfPixelHits
@ numberOfPixelHits
number of pixel layers on track with absence of hits
Definition Tracking/TrkEvent/TrkTrackSummary/TrkTrackSummary/TrackSummary.h:57
type
Generated on for ATLAS Offline Software by doxygen 1.14.0