88// In applying this license CERN does not waive the privileges and immunities
99// granted to it by virtue of its status as an Intergovernmental Organization
1010// or submit itself to any jurisdiction.
11+ // /
12+ // / \file ambiguousTrackPropagation.cxx
13+ // / \brief This code loops over central and MFT tracks and among the compatible
14+ // / collisions to this track, picks the one with the smallest DCAxy/DCAz and puts it
15+ // / in a table
16+ // / \author Anton Alkin <anton.alkin@cern.ch>
17+ // / \author Sarah Herrmann <sarah.herrmann@cern.ch>
18+ // / \author Gyula Bencedi <gyula.bencedi@cern.ch>
19+ // / \author Tulika Tripathy <tulika.tripathy@cern.ch>
1120
12- // \file trackPropagation.cxx
13- // \author Anton Alkin <anton.alkin@cern.ch>
14- // \author Sarah Herrmann <sarah.herrmann@cern.ch>
15- //
16- // \brief This code loops over central and MFT tracks and among the compatible
17- // collisions to this track, picks the one with the smallest DCAxy and puts it
18- // in a table
21+ #include " bestCollisionTable.h"
1922
20- #include " CCDB/BasicCCDBManager.h"
2123#include " Common/Core/trackUtilities.h"
24+ #include " Common/DataModel/CollisionAssociationTables.h"
2225#include " Common/DataModel/TrackSelectionTables.h"
26+
27+ #include " CCDB/BasicCCDBManager.h"
2328#include " CommonConstants/GeomConstants.h"
2429#include " DataFormatsParameters/GRPMagField.h"
2530#include " DetectorsBase/Propagator.h"
31+ #include " Field/MagneticField.h"
2632#include " Framework/AnalysisDataModel.h"
2733#include " Framework/AnalysisTask.h"
2834#include " Framework/runDataProcessing.h"
2935#include " ReconstructionDataFormats/TrackFwd.h"
36+
3037#include " Math/MatrixFunctions.h"
3138#include " Math/SMatrix.h"
32-
33- #include " Field/MagneticField.h"
3439#include " TGeoGlobalMagField.h"
3540
36- #include " Common/DataModel/CollisionAssociationTables.h "
37- #include " bestCollisionTable.h "
41+ #include < string >
42+ #include < vector >
3843
3944using SMatrix55 = ROOT::Math::SMatrix<double , 5 , 5 , ROOT::Math::MatRepSym<double , 5 >>;
4045using SMatrix5 = ROOT::Math::SVector<Double_t, 5 >;
@@ -46,33 +51,35 @@ using namespace o2;
4651using namespace o2 ::framework;
4752using namespace o2 ::aod::track;
4853
49- AxisSpec DCAxyAxis = {500 , -1 , 50 };
50-
5154struct AmbiguousTrackPropagation {
52- // Produces<aod::BestCollisions> tracksBestCollisions;
5355 Produces<aod::BestCollisionsFwd> fwdtracksBestCollisions;
56+ Produces<aod::BestCollisionsFwd3d> fwdtracksBestCollisions3d;
57+ Produces<aod::BestCollisionsFwd3dExtra> fwdtracksBestCollisions3dExtra;
5458 Produces<aod::BestCollFwdExtra> fwdtracksBestCollExtra;
5559 Produces<aod::ReassignedTracksCore> tracksReassignedCore;
5660 Produces<aod::ReassignedTracksExtra> tracksReassignedExtra;
5761 Service<o2::ccdb::BasicCCDBManager> ccdb;
5862
5963 int runNumber = -1 ;
60- float Bz = 0 ; // Magnetic field for MFT
61- static constexpr double centerMFT [3 ] = {0 , 0 , -61.4 }; // Field at center of MFT
64+ float bZ = 0 ; // Magnetic field for MFT
65+ static constexpr double kCenterMFT [3 ] = {0 , 0 , -61.4 }; // Field at center of MFT
6266
6367 o2::base::Propagator::MatCorrType matCorr =
6468 o2::base::Propagator::MatCorrType::USEMatCorrNONE;
6569
6670 o2::parameters::GRPMagField* grpmag = nullptr ;
6771
68- Configurable<std::string> ccdburl{" ccdb-url " , " http://alice-ccdb.cern.ch" , " url of the ccdb repository" };
72+ Configurable<std::string> ccdburl{" ccdburl " , " http://alice-ccdb.cern.ch" , " url of the ccdb repository" };
6973 Configurable<std::string> geoPath{" geoPath" , " GLO/Config/GeometryAligned" , " Path of the geometry file" };
7074 Configurable<std::string> grpmagPath{" grpmagPath" , " GLO/Config/GRPMagField" , " CCDB path of the GRPMagField object" };
7175 Configurable<std::string> mVtxPath {" mVtxPath" , " GLO/Calib/MeanVertex" , " Path of the mean vertex file" };
7276
7377 Configurable<bool > produceExtra{" produceExtra" , false , " Produce table with refitted track parameters" };
7478 Configurable<bool > produceHistos{" produceHistos" , false , " Produce control histograms" };
7579
80+ ConfigurableAxis binsDCAxy{" binsDCAxy" , {200 , -1 ., 1 .}, " " };
81+ ConfigurableAxis binsDCAz{" binsDCAz" , {200 , -1 ., 1 .}, " " };
82+
7683 HistogramRegistry registry{
7784 " registry" ,
7885 {
@@ -84,18 +91,27 @@ struct AmbiguousTrackPropagation {
8491
8592 void init (o2::framework::InitContext& /* initContext*/ )
8693 {
94+
95+ AxisSpec dcaXYAxis = {binsDCAxy, " dcaXYAxis" , " dcaXYAxis" };
96+ AxisSpec dcaZAxis = {binsDCAz, " dcaZAxis" , " dcaZAxis" };
97+
8798 ccdb->setURL (ccdburl);
8899 ccdb->setCaching (true );
89100 ccdb->setLocalObjectValidityChecking ();
90101
91102 if (produceHistos) {
92- if (doprocessMFT || doprocessMFTReassoc) {
103+ if (doprocessMFT || doprocessMFTReassoc || doprocessMFTReassoc3D ) {
93104 registry.add ({" DeltaZ" , " ; #Delta#it{z}" , {HistType::kTH1F , {{201 , -10.1 , 10.1 }}}});
94- registry.add ({" TracksDCAXY" , " ; DCA_{XY} (cm)" , {HistType::kTH1F , {DCAxyAxis }}});
95- registry.add ({" ReassignedDCAXY" , " ; DCA_{XY} (cm)" , {HistType::kTH1F , {DCAxyAxis }}});
96- registry.add ({" TracksOrigDCAXY" , " ; DCA_{XY} (wrt orig coll) (cm)" , {HistType::kTH1F , {DCAxyAxis }}});
105+ registry.add ({" TracksDCAXY" , " ; DCA_{XY} (cm)" , {HistType::kTH1F , {dcaXYAxis }}});
106+ registry.add ({" ReassignedDCAXY" , " ; DCA_{XY} (cm)" , {HistType::kTH1F , {dcaXYAxis }}});
107+ registry.add ({" TracksOrigDCAXY" , " ; DCA_{XY} (wrt orig coll) (cm)" , {HistType::kTH1F , {dcaXYAxis }}});
97108 registry.add ({" TracksAmbDegree" , " ; N_{coll}^{comp}" , {HistType::kTH1D , {{41 , -0.5 , 40.5 }}}});
98109 registry.add ({" TrackIsAmb" , " ; isAmbiguous" , {HistType::kTH1D , {{2 , -0.5 , 1.5 }}}});
110+ if (doprocessMFTReassoc3D) {
111+ registry.add ({" TracksDCAZ" , " ; DCA_{Z} (cm)" , {HistType::kTH1F , {dcaZAxis}}});
112+ registry.add ({" ReassignedDCAZ" , " ; DCA_{Z} (cm)" , {HistType::kTH1F , {dcaZAxis}}});
113+ registry.add ({" TracksOrigDCAZ" , " ; DCA_{Z} (wrt orig coll) (cm)" , {HistType::kTH1F , {dcaZAxis}}});
114+ }
99115 }
100116 if (doprocessCentral) {
101117 registry.add ({" PropagationFailures" , " " , {HistType::kTH1F , {{5 , 0.5 , 5.5 }}}});
@@ -122,18 +138,18 @@ struct AmbiguousTrackPropagation {
122138 o2::base::Propagator::initFieldFromGRP (grpmag);
123139 runNumber = bc.runNumber ();
124140
125- if (doprocessMFT || doprocessMFTReassoc) {
141+ if (doprocessMFT || doprocessMFTReassoc || doprocessMFTReassoc3D ) {
126142 o2::field::MagneticField* field = static_cast <o2::field::MagneticField*>(TGeoGlobalMagField::Instance ()->GetField ());
127- Bz = field->getBz (centerMFT );
128- LOG (info) << " The field at the center of the MFT is Bz = " << Bz ;
143+ bZ = field->getBz (kCenterMFT );
144+ LOG (info) << " The field at the center of the MFT is bZ = " << bZ ;
129145 }
130146 }
131147
132- static constexpr TrackSelectionFlags::flagtype trackSelectionITS =
148+ static constexpr TrackSelectionFlags::flagtype kTrackSelectionITS =
133149 TrackSelectionFlags::kITSNCls | TrackSelectionFlags::kITSChi2NDF |
134150 TrackSelectionFlags::kITSHits ;
135151
136- static constexpr TrackSelectionFlags::flagtype trackSelectionTPC =
152+ static constexpr TrackSelectionFlags::flagtype kTrackSelectionTPC =
137153 TrackSelectionFlags::kTPCNCls |
138154 TrackSelectionFlags::kTPCCrossedRowsOverNCls |
139155 TrackSelectionFlags::kTPCChi2NDF ;
@@ -153,18 +169,18 @@ struct AmbiguousTrackPropagation {
153169 std::array<float , 2 > dcaInfo;
154170 float bestDCA[2 ];
155171 o2::track::TrackParametrization<float > bestTrackPar;
156- for (auto & track : tracks) {
172+ for (auto const & track : tracks) {
157173 dcaInfo[0 ] = track.dcaXY (); // DCAxy
158174 dcaInfo[1 ] = track.dcaZ (); // DCAz
159175 bestDCA[0 ] = dcaInfo[0 ];
160176 bestDCA[1 ] = dcaInfo[1 ];
161177
162178 auto bestCol = track.has_collision () ? track.collisionId () : -1 ;
163- if ((track.trackCutFlag () & trackSelectionITS ) != trackSelectionITS ) {
179+ if ((track.trackCutFlag () & kTrackSelectionITS ) != kTrackSelectionITS ) {
164180 continue ;
165181 }
166182 if ((track.detectorMap () & (uint8_t )o2::aod::track::TPC) == (uint8_t )o2::aod::track::TPC) {
167- if ((track.trackCutFlag () & trackSelectionTPC ) != trackSelectionTPC ) {
183+ if ((track.trackCutFlag () & kTrackSelectionTPC ) != kTrackSelectionTPC ) {
168184 continue ;
169185 }
170186 }
@@ -181,7 +197,7 @@ struct AmbiguousTrackPropagation {
181197 }
182198 auto compatibleColls = track.compatibleColl ();
183199 int failures = 0 ;
184- for (auto & collision : compatibleColls) {
200+ for (auto const & collision : compatibleColls) {
185201 auto propagated = o2::base::Propagator::Instance ()->propagateToDCABxByBz ({collision.posX (), collision.posY (), collision.posZ ()}, trackPar, 2 .f , matCorr, &dcaInfo);
186202 if (!propagated) {
187203 ++failures;
@@ -239,7 +255,7 @@ struct AmbiguousTrackPropagation {
239255 float bestDCA = 0 .f , bestDCAx = 0 .f , bestDCAy = 0 .f ;
240256 o2::track::TrackParCovFwd bestTrackPar;
241257
242- for (auto & atrack : atracks) {
258+ for (auto const & atrack : atracks) {
243259 dcaInfo = 999 ; // DCAxy
244260 bestDCA = 999 ;
245261
@@ -254,14 +270,14 @@ struct AmbiguousTrackPropagation {
254270 int degree = 0 ; // degree of ambiguity of the track
255271
256272 auto compatibleBCs = atrack.bc_as <ExtBCs>();
257- for (auto & bc : compatibleBCs) {
273+ for (auto const & bc : compatibleBCs) {
258274 if (!bc.has_collisions ()) {
259275 continue ;
260276 }
261277 auto collisions = bc.collisions ();
262278 for (auto const & collision : collisions) {
263279 degree++;
264- trackPar.propagateToZhelix (collision.posZ (), Bz ); // track parameters propagation to the position of the z vertex
280+ trackPar.propagateToZhelix (collision.posZ (), bZ ); // track parameters propagation to the position of the z vertex
265281
266282 const auto dcaX (trackPar.getX () - collision.posX ());
267283 const auto dcaY (trackPar.getY () - collision.posY ());
@@ -325,7 +341,7 @@ struct AmbiguousTrackPropagation {
325341 float bestDCA = 0 .f , bestDCAx = 0 .f , bestDCAy = 0 .f ;
326342 o2::track::TrackParCovFwd bestTrackPar;
327343
328- for (auto & track : tracks) {
344+ for (auto const & track : tracks) {
329345 dcaInfo = 999 ; // DCAxy
330346 bestDCA = 999 ;
331347
@@ -345,9 +361,9 @@ struct AmbiguousTrackPropagation {
345361 SMatrix5 tpars (track.x (), track.y (), track.phi (), track.tgl (), track.signed1Pt ());
346362 o2::track::TrackParCovFwd trackPar{track.z (), tpars, tcovs, track.chi2 ()};
347363
348- for (auto & collision : compatibleColls) {
364+ for (auto const & collision : compatibleColls) {
349365
350- trackPar.propagateToZhelix (collision.posZ (), Bz ); // track parameters propagation to the position of the z vertex
366+ trackPar.propagateToZhelix (collision.posZ (), bZ ); // track parameters propagation to the position of the z vertex
351367
352368 const auto dcaX (trackPar.getX () - collision.posX ());
353369 const auto dcaY (trackPar.getY () - collision.posY ());
@@ -394,6 +410,88 @@ struct AmbiguousTrackPropagation {
394410 }
395411 }
396412 PROCESS_SWITCH (AmbiguousTrackPropagation, processMFTReassoc, " Fill BestCollisionsFwd for MFT ambiguous tracks with the new data model" , false );
413+
414+ void processMFTReassoc3D (MFTTracksWColls const & tracks, aod::Collisions const &, ExtBCs const & bcs)
415+ {
416+ if (bcs.size () == 0 ) {
417+ return ;
418+ }
419+ if (tracks.size () == 0 ) {
420+ return ;
421+ }
422+ auto bc = bcs.begin ();
423+ initCCDB (bc);
424+
425+ std::array<double , 3 > dcaInfOrig;
426+ std::array<double , 2 > dcaInfo;
427+ double bestDCA[2 ];
428+ o2::track::TrackParCovFwd bestTrackPar;
429+
430+ for (auto const & track : tracks) {
431+ dcaInfOrig[0 ] = 999 .f ; // original DCAx from propagation
432+ dcaInfOrig[1 ] = 999 .f ; // original DCAy from propagation
433+ dcaInfOrig[2 ] = 999 .f ; // original DCAz from propagation
434+ dcaInfo[0 ] = 999 .f ; // calcualted DCAxy
435+ dcaInfo[1 ] = 999 .f ; // calculated DCAz - same as original
436+ bestDCA[0 ] = 999 .f ; // minimal DCAxy
437+ bestDCA[1 ] = 999 .f ; // minimal DCAz
438+
439+ auto bestCol = track.has_collision () ? track.collisionId () : -1 ;
440+ auto compatibleColls = track.compatibleColl ();
441+
442+ std::vector<double > v1; // Temporary null vector for the computation of the covariance matrix
443+ SMatrix55 tcovs (v1.begin (), v1.end ());
444+ SMatrix5 tpars (track.x (), track.y (), track.phi (), track.tgl (), track.signed1Pt ());
445+ o2::track::TrackParCovFwd trackPar{track.z (), tpars, tcovs, track.chi2 ()};
446+
447+ for (auto const & collision : compatibleColls) {
448+
449+ trackPar.propagateToDCAhelix (bZ, {collision.posX (), collision.posY (), collision.posZ ()}, dcaInfOrig);
450+ dcaInfo[0 ] = std::sqrt (dcaInfOrig[0 ] * dcaInfOrig[0 ] + dcaInfOrig[1 ] * dcaInfOrig[1 ]);
451+ dcaInfo[1 ] = dcaInfOrig[2 ];
452+
453+ if ((std::abs (dcaInfo[0 ]) < std::abs (bestDCA[0 ])) && (std::abs (dcaInfo[1 ]) < std::abs (bestDCA[1 ]))) {
454+ bestCol = collision.globalIndex ();
455+ bestDCA[0 ] = dcaInfo[0 ];
456+ bestDCA[1 ] = dcaInfo[1 ];
457+ bestTrackPar = trackPar;
458+ }
459+ if ((track.collisionId () != collision.globalIndex ()) && produceHistos) {
460+ registry.fill (HIST (" DeltaZ" ), track.collision ().posZ () - collision.posZ ()); // deltaZ between the 1st coll zvtx and the other compatible ones
461+ }
462+ if (produceHistos) {
463+ registry.fill (HIST (" TracksDCAXY" ), dcaInfo[0 ]);
464+ registry.fill (HIST (" TracksDCAZ" ), dcaInfo[1 ]);
465+ }
466+
467+ if ((collision.globalIndex () == track.collisionId ()) && produceHistos) {
468+ registry.fill (HIST (" TracksOrigDCAXY" ), dcaInfo[0 ]);
469+ registry.fill (HIST (" TracksOrigDCAZ" ), dcaInfo[1 ]);
470+ }
471+ }
472+ if ((bestCol != track.collisionId ()) && produceHistos) {
473+ // reassigned
474+ registry.fill (HIST (" ReassignedDCAXY" ), bestDCA[0 ]);
475+ registry.fill (HIST (" ReassignedDCAZ" ), bestDCA[1 ]);
476+ }
477+ if (produceHistos) {
478+ int isAmbiguous = 0 ;
479+ registry.fill (HIST (" TracksAmbDegree" ), compatibleColls.size ());
480+ if (compatibleColls.size () > 1 ) {
481+ isAmbiguous = 1 ;
482+ }
483+ registry.fill (HIST (" TrackIsAmb" ), isAmbiguous);
484+ }
485+
486+ fwdtracksBestCollisions3d (track.globalIndex (), compatibleColls.size (), bestCol, bestDCA[0 ], bestDCA[1 ]);
487+ if (produceExtra) {
488+ fwdtracksBestCollisions3dExtra (bestTrackPar.getX (), bestTrackPar.getY (), bestTrackPar.getZ (),
489+ bestTrackPar.getTgl (), bestTrackPar.getInvQPt (), bestTrackPar.getPt (),
490+ bestTrackPar.getP (), bestTrackPar.getEta (), bestTrackPar.getPhi ());
491+ }
492+ }
493+ }
494+ PROCESS_SWITCH (AmbiguousTrackPropagation, processMFTReassoc3D, " Fill ReassignedTracks for MFT ambiguous tracks" , false );
397495};
398496
399497// ****************************************************************************************
0 commit comments