![Logo](../../ATLAS-Logo-Square-Blue-RGB.png) |
ATLAS Offline Software
|
Go to the documentation of this file.
18 #include "GaudiKernel/NTuple.h"
19 #include "GaudiKernel/INTupleSvc.h"
20 #include "GaudiKernel/SmartDataPtr.h"
22 #include "GaudiKernel/IPartPropSvc.h"
42 #include "CLHEP/GenericFunctions/CumulativeChiSquare.hh"
43 #include "CLHEP/Geometry/Point3D.h"
44 #include "CLHEP/Units/SystemOfUnits.h"
52 static const int maxTracks = 10000;
58 const std::string&
name,
78 m_truthToTrack(
"Trk::TruthToTrack", this),
79 m_extrapolator(
"Trk::Extrapolator/CosmicsExtrapolator", this),
80 m_particleCreator(
"Trk::TrackParticleCreatorTool/TrackParticleCreatorTool", this) {
81 declareInterface<IInDetAlignFillTrack>(
this);
100 "tool to produce perigee track parameters from generated parameters");
102 "tool to extrapolate tracks");
104 "tool to build TrackParticle");
129 IPartPropSvc* p_PartPropSvc;
130 static const bool CREATEIFNOTTHERE(
true);
131 ATH_CHECK(svcLoc()->service(
"PartPropSvc", p_PartPropSvc, CREATEIFNOTTHERE));
139 return StatusCode::SUCCESS;
183 return StatusCode::SUCCESS;
194 const EventContext& ctx = Gaudi::Hive::currentContext();
232 if (StatusCode::SUCCESS !=
dumpMatching(Uptracks, Lowtracks)) {
234 return StatusCode::FAILURE;
247 return StatusCode::FAILURE;
251 return StatusCode::FAILURE;
256 msg(
MSG::DEBUG) <<
"Retrieved " << truthCol->size() <<
" truth tracks from StoreGate" <<
endmsg;
265 for (; trackItr != trackItrE && nTracks < maxTracks; ++trackItr) {
267 if (
track ==
nullptr) {
283 TrackTruthCollection::const_iterator
found = truthCol->find(trackLink2);
285 if (
found != truthCol->end()) {
305 Amg::Vector3D productionVertex(genParticle->production_vertex()->position().x(),
306 genParticle->production_vertex()->position().y(),
307 genParticle->production_vertex()->position().z());
312 <<
", Status " << genParticle->status()
313 <<
", mass " << genParticle->momentum().m() <<
" CLHEP::MeV/c"
317 float genPt = std::sqrt((genParticle->momentum().x()) * (genParticle->momentum().x())
318 + (genParticle->momentum().y()) * (genParticle->momentum().y()));
321 <<
", p " << genParticle->momentum().e() /
CLHEP::GeV <<
" CLHEP::GeV/c"
322 <<
", eta " << genParticle->momentum().eta()
323 <<
", phi " << genParticle->momentum().phi() <<
" CLHEP::rad");
328 float pX = genParticle->momentum().px();
329 float pY = genParticle->momentum().py();
330 float genParticlePt = std::sqrt((pX * pX) + (pY * pY));
335 if (genParticle->pdg_id() == 0)
ATH_MSG_WARNING(
"Particle with PDG 0!");
337 "No GenVertex (generator level) production vertex found!");
344 generatedTrackPerigee =
m_truthToTrack->makePerigeeParameters(genParticle);
346 if (!generatedTrackPerigee) {
349 msg(
MSG::DEBUG) <<
"Trying to extrapolate directly to exclude material effects!" <<
endmsg;
354 sc =
evtStore()->retrieve(recordCollection,
"CaloEntryLayer");
355 if (
sc == StatusCode::FAILURE) {
360 << recordCollection->
size());
365 record != recordCollection->
end(); ++record) {
366 if (std::abs((*record).GetPDGCode()) != 13)
continue;
367 HepGeom::Point3D<double> productionVertex = (*record).GetPosition();
368 double charge = (*record).GetPDGCode() > 0 ? -1.0 : 1.0;
372 (*record).GetMomentum().y(),
373 (*record).GetMomentum().z());
377 double genPar_qOverP = 1. / direction.mag();
378 direction *= genPar_qOverP;
379 if (
charge < 0) genPar_qOverP = -genPar_qOverP;
382 double genPar_phi = direction.phi();
383 if (genPar_phi < 0.) genPar_phi = genPar_phi + 2.0*
M_PI;
386 << productionVertex.x() <<
", "
387 << productionVertex.y() <<
", "
388 << productionVertex.z() <<
")");
390 double genPar_theta = direction.theta();
394 globalSurfaceCentre.setIdentity();
396 productionVertex.y(), productionVertex.z());
400 const Amg::Vector3D productionVertexAsGlobalPosition(productionVertex.x(),
401 productionVertex.y(),
402 productionVertex.z());
406 genPar_phi, genPar_theta, genPar_qOverP,
418 << perigeeGlobalPosition.x() <<
", "
419 << perigeeGlobalPosition.y() <<
", "
420 << perigeeGlobalPosition.z() <<
")");
426 const Amg::Vector3D difference = productionVertexAsGlobalPosition - perigeeGlobalPosition;
428 double distance = std::sqrt(difference.x() * difference.x() + difference.y() * difference.y());
429 ATH_MSG_DEBUG(
"Distance between perigee point and generated vertex: "
436 ATH_MSG_DEBUG(
"Distance between perigee and generated vertex exceeds tolerance ("
437 << 1.
e-4 <<
" mm)... Extrapolating!");
440 *productionVertexTrackParams,
445 if (!generatedTrackPerigee)
continue;
447 ATH_MSG_DEBUG(
"Distance between perigee and generated vertex is less than tolerance ("
448 << 1.
e-4 <<
" CLHEP::mm)... " <<
"No propagation to perigee required");
452 genPar_phi, genPar_theta, genPar_qOverP,
461 delete productionVertexTrackParams;
462 delete generatedTrackPerigee;
468 if (generatedTrackPerigee) {
470 m_nt_mc_Trk_vtxX[nTracks] = genParticle->production_vertex()->position().x();
471 m_nt_mc_Trk_vtxY[nTracks] = genParticle->production_vertex()->position().y();
472 m_nt_mc_Trk_vtxZ[nTracks] = genParticle->production_vertex()->position().z();
474 delete generatedTrackPerigee;
508 std::string nt3id =
m_ntupleName +
"/TrkTrack_Matching";
544 return StatusCode::SUCCESS;
555 std::string
comments =
"Trk::Track Information";
629 else ATH_MSG_DEBUG(
"Ntuple " << nt0id <<
" has been booked successfully! ");
645 std::string
comments =
"Trk::UpTrack Information";
690 else ATH_MSG_DEBUG(
"Ntuple " << nt1id <<
" has been booked successfully! ");
704 std::string
comments =
"Trk::LowTrack Information";
759 ATH_MSG_DEBUG(
"Booking Matching between Up and Low Trk::Track Collections...");
761 std::string nt3id =
m_ntupleName +
"/TrkTrack_Matching";
762 std::string
comments =
"Matching between Up and Low Trk::Track Collections";
810 const std::string& TrkColName) {
819 for (; trackItr != trackItrE && itrk < maxTracks; ++trackItr) {
820 if (*trackItr !=
nullptr)
dumpTrack(itrk, (*trackItr), TrkColName);
832 const std::string& TrkColName) {
837 if (aMeasPer ==
nullptr) {
840 double d0 = aMeasPer->parameters()[
Trk::d0];
841 double z0 = aMeasPer->parameters()[
Trk::z0];
848 <<
" Track Parameters (d0, z0, phi0, theta, q/p)" <<
endmsg;
853 float transverseMomentum = std::sqrt((aMeasPer->momentum().x()) * (aMeasPer->momentum().x())
854 + (aMeasPer->momentum().y()) * (aMeasPer->momentum().y()));
857 <<
" pT = " << transverseMomentum /
CLHEP::GeV <<
" CLHEP::GeV/c");
860 int nHits = (*trk).measurementsOnTrack()->size();
869 int nhitspix = 0, nhitssct = 0, nhitstrt = 0;
870 int nshared = 0, nshpix = 0, nshsct = 0;
871 int nholes = 0, nhpix = 0, nhsct = 0;
876 if (not trackPart)
ATH_MSG_ERROR(
"Could not get xAOD::TrackParticle");
909 nshared = nshpix + nshsct;
911 if (nshpix < 0) nshpix = 0;
912 if (nshsct < 0) nshsct = 0;
913 if (nshared < 0) nshared = 0;
924 nholes = nhpix + nhsct;
926 if (nhpix < 0) nhpix = 0;
927 if (nhsct < 0) nhsct = 0;
928 if (nholes < 0) nholes = 0;
939 double chi2Prob = 0.;
941 if (fitQual ==
nullptr) {
946 Genfun::CumulativeChiSquare probabilityFunction(fitQual->
numberDoF());
947 chi2Prob = 1 - probabilityFunction(fitQual->
chiSquared());
957 if (TrkColName ==
"Up") {
985 }
else if (TrkColName ==
"Low") {
1054 float d0 = generatedTrackPerigee->parameters()[
Trk::d0];
1055 float z0 = generatedTrackPerigee->parameters()[
Trk::z0];
1056 float phi0 = generatedTrackPerigee->parameters()[
Trk::phi0];
1058 float eta = generatedTrackPerigee->
eta();
1060 float qoverp = generatedTrackPerigee->parameters()[
Trk::qOverP];
1062 float pt = (1 / qoverpt) * (
charge);
1066 msg(
MSG::DEBUG) <<
" - Extrapolated genParticle perigee parameters "
1067 <<
"(d0, z0, phi0, theta, q/p)" <<
endmsg;
1095 int nTracksUpper = 0;
1100 for (; trackItrUpper != trackItrUpperE; ++trackItrUpper) {
1101 const Trk::Track* trackUpper = *trackItrUpper;
1102 if (trackUpper ==
nullptr) {
1103 ATH_MSG_DEBUG(
"No associated Trk::Track object found for track " << nTracksUpper);
1110 float d0Up = measUpperPer->parameters()[
Trk::d0];
1111 float phi0Up = measUpperPer->parameters()[
Trk::phi0];
1112 float eta0Up = measUpperPer->eta();
1113 float z0Up = measUpperPer->parameters()[
Trk::z0];
1114 float thetaUp = measUpperPer->parameters()[
Trk::theta];
1116 float chargeUp = measUpperPer->charge();
1117 float ptUp = measUpperPer->pT() / 1000.;
1121 msg(
MSG::DEBUG) <<
" d0, z0, phi0, theta, q/p" << d0Up <<
", " << z0Up <<
", "
1122 << phi0Up <<
", " << thetaUp <<
", " << qOverPtUp <<
endmsg;
1125 int nTracksLower = 0;
1127 bool matchFound =
false;
1128 float Matched_Low_d0 = -999;
1129 float Matched_Low_phi0 = -999;
1131 float Matched_Low_theta = -999;
1132 float Matched_Low_eta0 = -999;
1133 float Matched_Low_z0 = -999;
1134 float Matched_Low_qOverPt = -999;
1135 float Matched_Low_charge = -999;
1136 float Matched_Low_pt = -999;
1141 for (; trackItrLower != trackItrLowerE; ++trackItrLower) {
1142 const Trk::Track* trackLower = *trackItrLower;
1143 if (trackLower ==
nullptr) {
1144 ATH_MSG_DEBUG(
"No associated Trk::Track object found for track " << nTracksLower);
1150 float d0Low = measLowerPer->parameters()[
Trk::d0];
1151 float phi0Low = measLowerPer->parameters()[
Trk::phi0];
1152 float eta0Low = measLowerPer->eta();
1153 float z0Low = measLowerPer->parameters()[
Trk::z0];
1154 float thetaLow = measLowerPer->parameters()[
Trk::theta];
1156 float chargeLow = measLowerPer->charge();
1157 float ptLow = measLowerPer->pT() / 1000.;
1160 float dphi2 = (phi0Up - phi0Low) * (phi0Up - phi0Low);
1164 float deta2 = (eta0Up - eta0Low) * (eta0Up - eta0Low);
1167 float dR = std::sqrt(dphi2 + deta2);
1172 Matched_Low_d0 = d0Low;
1173 Matched_Low_phi0 = phi0Low;
1175 Matched_Low_theta = thetaLow;
1176 Matched_Low_eta0 = eta0Low;
1177 Matched_Low_z0 = z0Low;
1178 Matched_Low_qOverPt = qOverPtLow;
1179 Matched_Low_charge = chargeLow;
1180 Matched_Low_pt = ptLow;
1185 msg(
MSG::DEBUG) <<
" d0, z0, phi0, theta, q/p: " << d0Low <<
", " << z0Low <<
", "
1186 << phi0Low <<
", " << thetaLow <<
", " << qOverPtLow <<
endmsg;
1215 return StatusCode::SUCCESS;
def retrieve(aClass, aKey=None)
HepMC::ConstGenParticlePtr scptr() const
Dereference/smart pointer.
NTuple::Array< float > m_nt_Trk_phi0_Up
phi0 parameter (Up track)
@ numberOfPixelHoles
number of pixel layers on track with absence of hits [unit8_t].
NTuple::Array< float > m_nt_Trk_chi2_Low
number of chi2 (Low track)
bool m_doTruth
switch on/off the truth information
NTuple::Array< float > m_nt_mc_Trk_genParticleEta
generated eta
NTuple::Array< float > m_nt_Trk_d0_Up
d0 parameter (Up track)
NTuple::Array< float > m_nt_Trk_phi0
phi0 parameter
Const iterator class for DataVector/DataList.
NTuple::Array< float > m_nt_Trk_d0_Low
d0 parameter (Low track)
void bookMatchingNtuple()
NTuple::Array< float > m_nt_mc_Trk_theta0
MonteCarlo theta0 parameter.
NTuple::Array< float > m_nt_mc_Trk_qoverpt
MonteCarlo q/pt parameter.
@ numberOfSCTSharedHits
number of SCT hits shared by several tracks [unit8_t].
NTuple::Array< float > m_nt_mc_Trk_genParticlePt
generated pt
void dumpTrack(int, const Trk::Track *, const std::string &)
double charge() const
Returns the charge.
NTuple::Array< float > m_nt_Trk_theta0
theta0 parameter
Scalar eta() const
pseudorapidity method
std::string m_ntupleName
ntuple name
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Dummy class used to allow special convertors to be called for surfaces owned by a detector element.
NTuple::Array< float > m_nt_Trk_pt
pt parameter
NTuple::Array< long > m_nt_Trk_nsharedPixels_Up
number of Pixel shared hits (Up track)
NTuple::Array< float > m_nt_mc_Trk_phi0
MonteCarlo phi0 parameter.
Scalar theta() const
theta method
NTuple::Array< long > m_nt_Trk_nshared_Low
number of shared hits (Low track)
NTuple::Array< long > m_nt_Trk_nholes_Low
number of holes (Low track)
NTuple::Array< float > m_nt_Trk_z0_Low
z0 parameter (Low track)
NTuple::Array< long > m_nt_Trk_nholes
number of holes
ParametersT< 5, Charged, PerigeeSurface > Perigee
bool msgLvl(const MSG::Level lvl) const
NTuple::Array< float > m_nt_Trk_chi2Prob_Low
number of chi2 probability (Low track)
NTuple::Array< float > m_nt_mc_trkistruth
Has the Track an associated truth track?
@ numberOfPixelHits
these are the pixel hits, including the b-layer [unit8_t].
NTuple::Array< float > m_nt_Trk_delta_z0
z0 parameter
@ numberOfTRTHits
number of TRT hits [unit8_t].
ToolHandle< Trk::ITrackParticleCreatorTool > m_particleCreator
Pointer to track particle creator tool.
NTuple::Array< float > m_nt_mc_Trk_charge
MonteCarlo charge parameter.
#define ATH_MSG_VERBOSE(x)
NTuple::Array< float > m_nt_Trk_pt_Low
pt parameter (Low track)
const Amg::Vector3D & center() const
Returns the center position of the Surface.
virtual ~InDetAlignFillTrack()
NTuple::Item< long > m_nt_nmctracks
number of mc tracks
const_iterator begin() const
NTuple::Array< long > m_nt_Trk_nhitstrt
number of TRT hits
int dumpTrackCol(const TrackCollection *)
NTuple::Array< float > m_nt_Trk_delta_eta
eta parameter
NTuple::Array< float > m_nt_Trk_delta_qoverpt
q/pt parameter
NTuple::Array< long > m_nt_Trk_nhitstrt_Up
number of TRT hits (Up track)
NTuple::Array< float > m_nt_Trk_qoverp_Low
q/p parameter (Low track)
virtual StatusCode initialize()
CONT::const_iterator const_iterator
NTuple::Array< long > m_nt_Trk_nholesPixels_Low
number of Pixel holes (Low track)
std::string m_TruthTrkCol
NTuple::Array< long > m_nt_Trk_nHits_Up
number of hits (Up track)
void dumpPerigee(const Trk::TrackParameters *, int)
NTuple::Array< float > m_nt_Trk_pt_Up
pt parameter (Up track)
NTuple::Array< float > m_nt_mc_Trk_prob
MonteCarlo prob parameter.
NTuple::Array< long > m_nt_Trk_nholes_Up
number of holes (Up track)
HepMC::ConstGenParticlePtr cptr() const
Dereference.
std::string m_inputLowCol
ServiceHandle< StoreGateSvc > & evtStore()
The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
a link optimized in size for a GenParticle in a McEventCollection
@ numberOfPixelSharedHits
number of Pixel all-layer hits shared by several tracks [unit8_t].
NTuple::Array< float > m_nt_Trk_z0_Up
z0 parameter (Up track)
virtual StatusCode FillTrack()
NTuple::Array< long > m_nt_Trk_nhitssct_Low
number of SCT hits (Low track)
NTuple::Array< long > m_nt_Trk_ndof_Low
number of ndof (Low track)
NTuple::Array< long > m_nt_Trk_nholesSCT_Up
number of SCT holes (Up track)
NTuple::Array< long > m_nt_Trk_nsharedSCT_Up
number of SCT shared hits (Up track)
NTuple::Item< long > m_nt_nLowtracks
number of Low tracks
NTuple::Item< long > m_nt_ntracks
number of tracks
::StatusCode StatusCode
StatusCode definition for legacy code.
NTuple::Array< long > m_nt_Trk_nholesSCT
number of SCT holes
NTuple::Array< float > m_nt_Trk_qoverp
q/p parameter
NTuple::Array< long > m_nt_Trk_nholesSCT_Low
number of SCT holes (Low track)
Eigen::Affine3d Transform3D
NTuple::Array< float > m_nt_Trk_phi0_Low
phi0 parameter (Low track)
NTuple::Array< float > m_nt_Trk_chi2_Up
number of chi2 (Up track)
NTuple::Array< float > m_nt_Trk_d0
d0 parameter
NTuple::Array< long > m_nt_Trk_nhitspix_Low
number of Pixel hits (Low track)
Class to represent and store fit qualities from track reconstruction in terms of and number of degre...
NTuple::Item< long > m_nt_nUptracks
number of Up tracks
NTuple::Array< long > m_nt_Trk_nhitssct_Up
number of SCT hits (Up track)
NTuple::Array< long > m_nt_Trk_nhitspix_Up
number of Pixel hits (Up track)
@ numberOfSCTHoles
number of SCT holes [unit8_t].
NTuple::Array< float > m_nt_Trk_delta_pt
pt parameter
bool isValid() const
Validity check.
NTuple::Array< long > m_nt_Trk_nsharedSCT_Low
number of SCT shared hits (Low track)
NTuple::Array< long > m_nt_Trk_nhitstrt_Low
number of TRT hits (Low track)
const HepMcParticleLink & particleLink() const
ToolHandle< Trk::IExtrapolator > m_extrapolator
Pointer to IExtrapolator.
NTuple::Array< float > m_nt_Trk_chi2
number of chi2
NTuple::Array< float > m_nt_Trk_chi2Prob_Up
number of chi2 probability (Up track)
bool setElement(ElementType element)
Set to point to an element.
const Perigee * perigeeParameters() const
return Perigee.
const GenParticle * ConstGenParticlePtr
NTuple::Array< long > m_nt_Trk_nHits_Low
number of hits (Low track)
bool setStorableObject(BaseConstReference data, bool replace=false, IProxyDict *sg=0)
Set link to point to a new container (storable).
NTuple::Array< float > m_nt_Trk_delta_theta0
theta parameter
double charge(const T &p)
NTuple::Array< float > m_nt_mc_Trk_qoverp
MonteCarlo q/p parameter.
virtual StatusCode finalize()
NTuple::Array< long > m_nt_Trk_nsharedPixels
number of Pixel shared hits
Eigen::Matrix< double, 3, 1 > Vector3D
NTuple::Array< long > m_nt_Trk_nholesPixels_Up
number of Pixel holes (Up track)
NTuple::Array< float > m_nt_Trk_z0
z0 parameter
NTuple::Array< float > m_nt_mc_Trk_eta
MonteCarlo eta parameter.
NTuple::Array< float > m_nt_Trk_chi2Prob
number of chi2 probability
NTuple::Array< float > m_nt_mc_Trk_vtxX
MonteCarlo Vertex.X parameter.
NTuple::Array< float > m_nt_mc_Trk_genParticlePhi
generated phi
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
NTuple::Array< float > m_nt_Trk_delta_charge
charge parameter
ParametersT< 5, Charged, PlaneSurface > AtaPlane
An STL vector of pointers that by default owns its pointed-to elements.
MC particle associated with a reco track + the quality of match.
NTuple::Array< long > m_nt_Trk_nsharedSCT
number of SCT shared hits
NTuple::Array< long > m_nt_Trk_nshared
number of shared hits
NTuple::Array< long > m_nt_Trk_nholesPixels
number of Pixel holes
NTuple::Array< float > m_nt_Trk_delta_phi0
phi0 parameter
#define ATH_MSG_WARNING(x)
NTuple::Array< long > m_nt_Trk_nHits
number of hits
NTuple::Array< long > m_nt_Trk_ndof
number of ndof
const_iterator end() const
float probability() const
NTuple::Array< float > m_nt_Trk_theta0_Low
theta0 parameter (Low track)
InDetAlignFillTrack(const std::string &type, const std::string &name, const IInterface *parent)
Eigen::Translation< double, 3 > Translation3D
ToolHandle< Trk::ITruthToTrack > m_truthToTrack
Pointer to TruthToTrack.
NTuple::Array< float > m_nt_Trk_delta_d0
d0 parameter
NTuple::Array< float > m_nt_Trk_qoverp_Up
q/p parameter (Up track)
NTuple::Array< long > m_nt_Trk_nhitssct
number of SCT hits
@ numberOfSCTHits
number of hits in SCT [unit8_t].
double chiSquared() const
returns the of the overall track fit
double eta() const
Access method for pseudorapidity - from momentum.
StatusCode dumpMatching(const TrackCollection *, const TrackCollection *)
int numberDoF() const
returns the number of degrees of freedom of the overall track or vertex fit as integer
NTuple::Array< long > m_nt_Trk_nhitspix
number of Pixel hits
NTuple::Array< float > m_nt_mc_Trk_pdg
MonteCarlo pdg parameter.
NTuple::Array< float > m_nt_mc_Trk_vtxY
MonteCarlo Vertex.Y parameter.
Class describing a TrackParticle.
NTuple::Array< float > m_nt_mc_Trk_pt
MonteCarlo pt parameter.
NTuple::Array< float > m_nt_mc_Trk_z0
MonteCarlo z0 parameter.
bool m_doMatching
switch on/off the matching information
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
size_type size() const noexcept
Returns the number of elements in the collection.
NTuple::Array< long > m_nt_Trk_nshared_Up
number of shared hits (Up track)
NTuple::Array< long > m_nt_Trk_nsharedPixels_Low
number of Pixel shared hits (Low track)
bool empty() const noexcept
Returns true if the collection is empty.
NTuple::Array< float > m_nt_mc_Trk_vtxZ
MonteCarlo Vertex.Z parameter.
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
NTuple::Item< long > m_nt_matchingTrk
matching tracks
NTuple::Array< float > m_nt_Trk_theta0_Up
theta0 parameter (Up track)
NTuple::Array< float > m_nt_mc_Trk_d0
MonteCarlo d0 parameter.
NTuple::Array< long > m_nt_Trk_ndof_Up
number of ndof (Up track)