ATLAS Offline Software
TrigPrimaryVertexFitter.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 // TrigPrimaryVertexFitter tool
7 // -------------------------------
8 // ATLAS Collaboration
9 //
10 // 04.05.2008 Package created
11 //
12 // Author: Dmitry Emeliyanov, RAL
13 // e-mail: D.Emeliyanov@rl.ac.uk
14 //
16 #include <cmath>
17 #include <iostream>
18 #include <algorithm>
19 
20 #include "TrigInDetEvent/TrigL2Vertex.h"
21 #include "TrigInDetEvent/TrigVertex.h"
22 
23 #include "TrkTrack/Track.h"
24 #include "TrkTrack/TrackCollection.h"
25 
26 #include "TrigInDetToolInterfaces/ITrigPrimaryVertexFitter.h"
27 #include "TrigPrimaryVertexFitter.h"
28 #include "TrigPrimaryVertexTrack.h"
29 
30 
31 TrigPrimaryVertexFitter::TrigPrimaryVertexFitter(const std::string& t,
32  const std::string& n,
33  const IInterface* p ): AthAlgTool(t,n,p),
34  m_createTrackLists(false),
35  m_xyVariance(10.0)
36 {
37  declareInterface< ITrigPrimaryVertexFitter >( this );
38  declareProperty( "numberOfIterations", m_numIter=3);
39  declareProperty( "maxChi2Increase", m_maxChi2Increase=1000.0);
40  declareProperty( "zVariance", m_zVariance=1.0);
41  declareProperty( "Chi2Cut", m_chi2cut=15.0);
42  declareProperty( "Chi2CutLoose", m_chi2cut_loose=100.0);
43  declareProperty( "CreateTrackLists", m_createTrackLists=false);
44  declareProperty( "xyVariance", m_xyVariance=10.0);
45 }
46 
47 StatusCode TrigPrimaryVertexFitter::initialize()
48 {
49  ATH_CHECK( AthAlgTool::initialize() );
50 
51  ATH_MSG_INFO( "Number of iterations is set to " << m_numIter );
52  ATH_MSG_INFO( "Accuracy of starting z-position is " << m_zVariance<<" mm" );
53  ATH_MSG_INFO( "Track lists created ? " << std::boolalpha<<m_createTrackLists );
54  return StatusCode::SUCCESS;
55 }
56 
57 StatusCode TrigPrimaryVertexFitter::finalize()
58 {
59  ATH_CHECK( AthAlgTool::finalize() );
60  return StatusCode::SUCCESS;
61 }
62 
63 TrigPrimaryVertexFitter::~TrigPrimaryVertexFitter()
64 {
65 
66 }
67 
68 //TrigVertex has associated TrigInDetTracks, but we want Trk::Tracks
69 TrigVertex* TrigPrimaryVertexFitter::fit(const TrackCollection* tc, TrackCollection& output_tc, double z0) const
70 {
71 
72  TrigVertex* pVertex=NULL;
73 
74  double chi2;
75  int ndof;
76 
77  std::vector<TrigPrimaryVertexTrack*> tracks;
78  tracks.clear();
79  int idx=1;
80  for(TrackCollection::const_iterator ptIt=tc->begin();ptIt!=tc->end();++ptIt)
81  {
82  TrigPrimaryVertexTrack* pT= new TrigPrimaryVertexTrack(*ptIt);
83  if(pT!=NULL) {
84  pT->setIndex(idx++);
85  tracks.push_back(pT);
86  }
87  }
88 
89  TrigL2Vertex* pV = new TrigL2Vertex();
90 
91  pV->getParametersVector()[0]=0.0;
92  pV->getParametersVector()[1]=0.0;
93  pV->getParametersVector()[2]=z0;
94  bool fitFailed=false;
95  chi2=0.0;ndof=-3;
96  for(int nIter=0;nIter<m_numIter;nIter++)
97  {
98  memset(&pV->m_Gk[0][0],0,sizeof(pV->m_Gk));
99  pV->m_Gk[0][0]=m_xyVariance*m_xyVariance;
100  pV->m_Gk[1][1]=m_xyVariance*m_xyVariance;
101  pV->m_Gk[2][2]=m_zVariance*m_zVariance;
102  chi2=0.0;ndof=-3;
103  fitFailed=false;
104 
105  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
106  {
107  if(!(*it)->isActive())
108  (*it)->activate();
109 
110  double dchi2=(*it)->getChi2Distance(pV);
111  ATH_MSG_VERBOSE("Track "<<(*it)->getIndex()<<" dchi2="<<dchi2 );
112  if(std::isnan(dchi2)||(dchi2<0.0)||(dchi2>m_maxChi2Increase))
113  {
114  fitFailed=true;
115  break;
116  }
117  if(dchi2<m_chi2cut)
118  {
119  chi2+=dchi2;ndof+=2;
120  (*it)->updateVertex(pV);
121  }
122  else
123  {
124  (*it)->mask();
125  ATH_MSG_DEBUG("Skipping track "<<(*it)->getIndex()<<"due to large dchi2="<<dchi2 );
126  }
127  }
128  if(fitFailed) break;
129  }
130 
131 
132  if(!fitFailed && (ndof>0))
133  {
134  ATH_MSG_DEBUG("Primary vertex fit OK, chi2 = "<<chi2<<" NDOF = "<<ndof );
135  ATH_MSG_DEBUG("x= "<<pV->getParametersVector()[0]<<
136  " y= "<<pV->getParametersVector()[1]<<" z= "<<pV->getParametersVector()[2] );
137  double cv[6];
138  cv[0]=pV->m_Gk[0][0];cv[1]=pV->m_Gk[0][1];cv[2]=pV->m_Gk[1][1];
139  cv[3]=pV->m_Gk[0][2];cv[4]=pV->m_Gk[1][2];cv[5]=pV->m_Gk[2][2];
140  if(m_createTrackLists) {
141  output_tc.clear();
142  output_tc.reserve(tracks.size());
143 
144  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it) {
145  if ((*it)->isActive())
146  output_tc.push_back(new Trk::Track(*(*it)->getTrkTrack()));
147  }
148  }
149  pVertex = new TrigVertex(pV->getParametersVector()[0],
150  pV->getParametersVector()[1],
151  pV->getParametersVector()[2],cv,chi2,ndof, nullptr);
152  }
153  delete pV;
154  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
155  {
156  delete (*it);
157  }
158  return pVertex;
159 }
160 
161 struct TrigPrimaryTrackSortPredicate
162 {
163  bool operator()(TrigPrimaryVertexTrack* pA, TrigPrimaryVertexTrack* pB)
164  {
165  return (pA->getChi2Contribution() < pB->getChi2Contribution());
166  }
167 };
168 
169 TrigVertex* TrigPrimaryVertexFitter::fit(const TrackCollection* tc, std::vector<double>& Chi2V, double z0) const
170 {
171 
172  TrigVertex* pVertex=NULL;
173 
174  double chi2;
175  int ndof;
176 
177  std::vector<TrigPrimaryVertexTrack*> tracks;
178  tracks.clear();
179  Chi2V.clear();
180  int idx=1;
181  for(TrackCollection::const_iterator ptIt=tc->begin();ptIt!=tc->end();++ptIt)
182  {
183  TrigPrimaryVertexTrack* pT= new TrigPrimaryVertexTrack(*ptIt);
184  if(pT!=NULL)
185  {
186  pT->setIndex(idx++);
187  tracks.push_back(pT);
188  Chi2V.push_back(-100.0);
189  }
190  }
191 
192  TrigL2Vertex* pV = new TrigL2Vertex();
193 
194  // 1. preliminary fit
195 
196  pV->getParametersVector()[0]=0.0;
197  pV->getParametersVector()[1]=0.0;
198  pV->getParametersVector()[2]=z0;
199 
200  memset(&pV->m_Gk[0][0],0,sizeof(pV->m_Gk));
201  pV->m_Gk[0][0]=m_xyVariance*m_xyVariance;
202  pV->m_Gk[1][1]=m_xyVariance*m_xyVariance;
203  pV->m_Gk[2][2]=m_zVariance*m_zVariance;
204 
205  bool fitFailed=false;
206  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
207  {
208  double dchi2=(*it)->getChi2Distance(pV);
209  if(std::isnan(dchi2)||(dchi2<0.0)||(dchi2>m_maxChi2Increase))
210  {
211  fitFailed=true;
212  break;
213  }
214  if(dchi2<m_chi2cut_loose) (*it)->updateVertex(pV);
215  }
216  if(fitFailed)
217  {
218  delete pV;
219  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();
220  it!=tracks.end();++it) delete (*it);
221  return pVertex;
222  }
223 
224  // 2. calculation of "smoothed" chi2s
225 
226  memset(&pV->m_Gk[0][0],0,sizeof(pV->m_Gk));
227 
228  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
229  {
230  double dchi2=(*it)->getChi2Distance(pV);
231  ATH_MSG_DEBUG("Track "<<(*it)->getIndex()<<" dchi2="<<dchi2 );
232  }
233 
234  // 3. sort tracks in accordance with their chi2s
235 
236  std::sort(tracks.begin(),tracks.end(),TrigPrimaryTrackSortPredicate());
237 
238  // 4. Final track fit
239 
240  pV->getParametersVector()[0]=0.0;
241  pV->getParametersVector()[1]=0.0;
242  pV->getParametersVector()[2]=z0;
243  fitFailed=false;
244  chi2=0.0;ndof=-3;
245  for(int nIter=0;nIter<m_numIter;nIter++)
246  {
247  memset(&pV->m_Gk[0][0],0,sizeof(pV->m_Gk));
248  pV->m_Gk[0][0]=100.0;
249  pV->m_Gk[1][1]=100.0;
250  pV->m_Gk[2][2]=m_zVariance*m_zVariance;
251  chi2=0.0;ndof=-3;
252  fitFailed=false;
253 
254  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
255  {
256  double dchi2=(*it)->getChi2Distance(pV);
257  Chi2V[(*it)->getIndex()-1]=dchi2;
258  ATH_MSG_VERBOSE("Track "<<(*it)->getIndex()<<" dchi2="<<dchi2 );
259  if(std::isnan(dchi2)||(dchi2<0.0)||(dchi2>m_maxChi2Increase))
260  {
261  fitFailed=true;
262  break;
263  }
264  if(dchi2<m_chi2cut)
265  {
266  chi2+=dchi2;ndof+=2;
267  (*it)->updateVertex(pV);
268  }
269  else
270  {
271  ATH_MSG_DEBUG("Skipping track "<<(*it)->getIndex()<<"due to large dchi2="<<dchi2 );
272  }
273  }
274  if(fitFailed) break;
275  }
276  for(std::vector<TrigPrimaryVertexTrack*>::iterator it=tracks.begin();it!=tracks.end();++it)
277  {
278  delete (*it);
279  }
280 
281  if(!fitFailed && (ndof>0))
282  {
283  ATH_MSG_DEBUG("Primary vertex fit OK, chi2 = "<<chi2<<" NDOF = "<<ndof );
284  ATH_MSG_DEBUG("x= "<<pV->getParametersVector()[0]<<
285  " y= "<<pV->getParametersVector()[1]<<" z= "<<pV->getParametersVector()[2] );
286  double cv[6];
287  cv[0]=pV->m_Gk[0][0];cv[1]=pV->m_Gk[0][1];cv[2]=pV->m_Gk[1][1];
288  cv[3]=pV->m_Gk[0][2];cv[4]=pV->m_Gk[1][2];cv[5]=pV->m_Gk[2][2];
289  pVertex = new TrigVertex(pV->getParametersVector()[0],
290  pV->getParametersVector()[1],
291  pV->getParametersVector()[2],cv,chi2,ndof,NULL);
292  }
293  delete pV;
294 
295  return pVertex;
296 }
297 
DataVector::reserve
void reserve(size_type n)
Attempt to preallocate enough memory for a specified number of elements.
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
TrigPrimaryVertexFitter::m_createTrackLists
bool m_createTrackLists
Definition: TrigPrimaryVertexFitter.h:24
IDTPM::ndof
float ndof(const U &p)
Definition: TrackParametersHelper.h:134
CalculateHighPtTerm.pT
pT
Definition: ICHEP2016/CalculateHighPtTerm.py:57
python.tests.PyTestsLib.finalize
def finalize(self)
_info( "content of StoreGate..." ) self.sg.dump()
Definition: PyTestsLib.py:50
DataModel_detail::const_iterator
Const iterator class for DataVector/DataList.
Definition: DVLIterator.h:82
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
Trk::Track
The ATLAS Track class.
Definition: Tracking/TrkEvent/TrkTrack/TrkTrack/Track.h:73
AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
TrigPrimaryVertexTrack.h
initialize
void initialize()
Definition: run_EoverP.cxx:894
TrigPrimaryVertexFitter::~TrigPrimaryVertexFitter
virtual ~TrigPrimaryVertexFitter()
Definition: TrigPrimaryVertexFitter.cxx:63
skel.it
it
Definition: skel.GENtoEVGEN.py:396
TrigPrimaryVertexFitter::m_chi2cut
double m_chi2cut
Definition: TrigPrimaryVertexFitter.h:23
TrigPrimaryVertexFitter::m_numIter
int m_numIter
Definition: TrigPrimaryVertexFitter.h:22
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
Track.h
TrigPrimaryVertexFitter::m_xyVariance
double m_xyVariance
Definition: TrigPrimaryVertexFitter.h:25
TrigPrimaryVertexFitter.h
TrigL2Vertex.h
TrigL2Vertex
Definition: TrigL2Vertex.h:146
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:210
TrigPrimaryVertexTrack::getChi2Contribution
double getChi2Contribution()
chi2-contribution to the vertex fit
Definition: TrigPrimaryVertexTrack.cxx:243
beamspotman.n
n
Definition: beamspotman.py:731
TrigPrimaryVertexTrack
Definition: TrigPrimaryVertexTrack.h:15
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
TrackCollection.h
chi2
double chi2(TH1 *h0, TH1 *h1)
Definition: comparitor.cxx:523
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
fitFailed
bool fitFailed
Definition: fbtTestToyMC.cxx:116
TRT::Track::z0
@ z0
Definition: InnerDetector/InDetCalibEvent/TRT_CalibData/TRT_CalibData/TrackInfo.h:63
DataVector< Trk::Track >
DataVector::clear
void clear()
Erase all the elements in the collection.
TrigPrimaryTrackSortPredicate::operator()
bool operator()(TrigPrimaryVertexTrack *pA, TrigPrimaryVertexTrack *pB)
Definition: TrigPrimaryVertexFitter.cxx:163
TrigL2Vertex::m_Gk
double m_Gk[MAX_SIZE_VERT_COVM][MAX_SIZE_VERT_COVM]
Definition: TrigL2Vertex.h:182
DataVector::push_back
value_type push_back(value_type pElem)
Add an element to the end of the collection.
DataVector::end
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
TrigVertex.h
TrigVertex
Definition: TrigVertex.h:28
TrigL2Vertex::getParametersVector
double * getParametersVector()
returns vector of vertex fit parameters: vertex position + refitted track momenta at-perigee (sic !...
Definition: TrigL2Vertex.cxx:651
TrigPrimaryVertexFitter::finalize
virtual StatusCode finalize()
Definition: TrigPrimaryVertexFitter.cxx:57
LArNewCalib_DelayDump_OFC_Cali.idx
idx
Definition: LArNewCalib_DelayDump_OFC_Cali.py:69
TrigPrimaryVertexFitter::m_maxChi2Increase
double m_maxChi2Increase
Definition: TrigPrimaryVertexFitter.h:23
ITrigPrimaryVertexFitter.h
TrigPrimaryVertexFitter::fit
virtual TrigVertex * fit(const TrackCollection *, TrackCollection &, double z=0.0) const
Definition: TrigPrimaryVertexFitter.cxx:69
TrigPrimaryVertexFitter::TrigPrimaryVertexFitter
TrigPrimaryVertexFitter(const std::string &, const std::string &, const IInterface *)
Definition: TrigPrimaryVertexFitter.cxx:31
AthAlgTool
Definition: AthAlgTool.h:26
TrigPrimaryVertexFitter::m_zVariance
double m_zVariance
Definition: TrigPrimaryVertexFitter.h:23
TrigPrimaryTrackSortPredicate
Definition: TrigPrimaryVertexFitter.cxx:162
TrigPrimaryVertexFitter::m_chi2cut_loose
double m_chi2cut_loose
Definition: TrigPrimaryVertexFitter.h:23
DataVector::begin
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
TrigPrimaryVertexFitter::initialize
virtual StatusCode initialize()
Definition: TrigPrimaryVertexFitter.cxx:47
Generated on Mon Dec 23 2024 21:21:32 for ATLAS Offline Software by   doxygen 1.8.18