@@ -64,9 +64,11 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
6464 std ::vector < o2 ::itsmft ::ROFRecord > * digiRofRecordsArr {nullptr };
6565 digiTree -> SetBranchAddress ("TF3DigitROF" , & digiRofRecordsArr );
6666 auto& digiRofArr = * digiRofRecordsArr ;
67- o2 ::dataformats ::IOMCTruthContainerView * digiPlabelsArr {nullptr };
68- digiTree -> SetBranchAddress ("TF3DigitMCTruth" , & digiPlabelsArr );
67+ o2 ::dataformats ::IOMCTruthContainerView * digiLabelsArr {nullptr };
68+ digiTree -> SetBranchAddress ("TF3DigitMCTruth" , & digiLabelsArr );
6969 digiTree -> GetEntry (0 );
70+ o2 ::dataformats ::ConstMCTruthContainer < o2 ::MCCompLabel > digiLabels ;
71+ digiLabelsArr -> copyandflatten (digiLabels );
7072
7173 // Clusters
7274 TFile * clsFile = TFile ::Open (clsFilePath .data ());
@@ -76,13 +78,21 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
7678 std ::vector < o2 ::itsmft ::ROFRecord > * clsRofRecordsArr {nullptr };
7779 clsTree -> SetBranchAddress ("TF3ClusterROF" , & clsRofRecordsArr );
7880 auto& clsRofArr = * clsRofRecordsArr ;
79- o2 ::dataformats ::IOMCTruthContainerView * clsPlabelsArr {nullptr };
80- clsTree -> SetBranchAddress ("TF3ClusterMCTruth" , & clsPlabelsArr );
81+ o2 ::dataformats ::MCTruthContainer < o2 :: MCCompLabel > * clsLabels {nullptr };
82+ clsTree -> SetBranchAddress ("TF3ClusterMCTruth" , & clsLabels );
8183 clsTree -> GetEntry (0 );
82-
83- // Start script
84- o2 ::dataformats ::ConstMCTruthContainer < o2 ::MCCompLabel > labels ;
85- clsPlabelsArr -> copyandflatten (labels );
84+
85+ // Summary of entries in all branches
86+ std ::cout << std ::endl ;
87+ std ::cout << "---> Number of digits: " << digitsArray -> size () << std ::endl ;
88+ std ::cout << "---> Number of digit ROFs: " << digiRofArr .size () << std ::endl ;
89+ std ::cout << "---> Number of clusters: " << clsArray -> size () << std ::endl ;
90+ std ::cout << "---> Number of cluster ROFs: " << clsRofArr .size () << std ::endl ;
91+ std ::cout << "---> Number of digits with MC label: " << digiLabels .getNElements () << std ::endl ;
92+ std ::cout << "---> Number of digits with MC label: " << digiLabels .getIndexedSize () << std ::endl ;
93+ std ::cout << "---> Number of clusters with MC label: " << clsLabels -> getNElements () << std ::endl ;
94+ std ::cout << "---> Number of clusters with MC label: " << clsLabels -> getIndexedSize () << std ::endl ;
95+ std ::cout << std ::endl ;
8696
8797 auto clsTuple = new TNtuple ("clsTuple" , "clsTuple" , "chip_id:x:y:z:subdet_id:row:col:time" );
8898 clsTuple -> SetDirectory (nullptr );
@@ -93,6 +103,21 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
93103 TH1F * histXCoordDigit = new TH1F ("histXCoordDigit" , "histXCoordDigit" , 8000 , -100 , 100 );
94104 TH1F * histYCoordDigit = new TH1F ("histYCoordDigit" , "histYCoordDigit" , 8000 , -100 , 100 );
95105 TH1F * histZCoordDigit = new TH1F ("histZCoordDigit" , "histZCoordDigit" , 28000 , -400 , 400 );
106+ TH1F * histXCoordRes = new TH1F ("histXCoordRes" , "histXCoordRes" , 100 , -0.05 , 0.05 );
107+ TH1F * histYCoordRes = new TH1F ("histYCoordRes" , "histYCoordRes" , 100 , -0.05 , 0.05 );
108+ TH1F * histZCoordRes = new TH1F ("histZCoordRes" , "histZCoordRes" , 100 , -0.05 , 0.05 );
109+ TH1F * histTimeRes = new TH1F ("histTimeRes" , "histTimeRes" , 100 , -0.05 , 0.05 );
110+
111+ // Load all digits upfront and build a lookup map
112+ int nDigits = digiTree -> GetEntries ();
113+ std ::unordered_map < o2 ::MCCompLabel , int > digitsLabels ;
114+ for (int iDigit = 0 ; iDigit < digitsArray -> size (); ++ iDigit ) {
115+ auto label = digiLabels .getLabels (iDigit )[0 ];
116+ if (!label .isValid ()) {
117+ continue ;
118+ }
119+ digitsLabels .emplace (label , iDigit );
120+ }
96121
97122 // LOOP on : ROFRecord array
98123 for (unsigned int iROF = 0 ; iROF < clsRofArr .size (); ++ iROF ) {
@@ -101,8 +126,9 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
101126 const unsigned int rofNEntries = clsRofArr [iROF ].getNEntries ();
102127
103128 // LOOP on : digits array
129+ std ::cout << "\n\n ----> Starting loop on digits for ROF " << iROF << " with index " << rofIndex << " and nEntries " << rofNEntries << std ::endl ;
104130 for (unsigned int iDigit = rofIndex ; iDigit < rofIndex + rofNEntries ; iDigit ++ ) {
105- if (iDigit % 1000 == 0 ) {
131+ if (iDigit % 10000 == 0 ) {
106132 std ::cout << "Reading digit " << iDigit << " / " << digitsArray -> size () << std ::endl ;
107133 }
108134
@@ -127,15 +153,17 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
127153
128154
129155 // LOOP on : clusters array
156+ std ::cout << "\n\n ----> Starting loop on clusters for ROF " << iROF << " with index " << rofIndex << " and nEntries " << rofNEntries << std ::endl ;
130157 for (unsigned int iCls = rofIndex ; iCls < rofIndex + rofNEntries ; iCls ++ ) {
131- if (iCls % 1000 == 0 ) {
158+ if (iCls % 10000 == 0 ) {
132159 std ::cout << "Reading cluster " << iCls << " / " << clsArray -> size () << std ::endl ;
133160 }
134161
135162 Int_t iRow = (* clsArray )[iCls ].row ;
136163 Int_t iCol = (* clsArray )[iCls ].col ;
137164 Int_t chipID = (* clsArray )[iCls ].chipID ;
138165 Int_t subDetID = tofGeo -> getIOTOFLayer (chipID );
166+ Float_t time = (* clsArray )[iCls ].time ;
139167
140168 Float_t x = 0.f , y = 0.f , z = 0.f ;
141169 if (subDetID >= 0 ) {
@@ -144,7 +172,6 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
144172
145173 o2 ::math_utils ::Point3D < float > localClsCoords (x , y , z ); // local Digit
146174 const auto globalClsCoords = tofGeo -> getMatrixL2G (chipID )(localClsCoords ); // convert to global
147- std ::cout << "Cluster " << iCls << ": chipID = " << chipID << ", X=" << globalClsCoords .x () << ", Y=" << globalClsCoords .y () << ", Z=" << globalClsCoords .z () << ", time=" << (* clsArray )[iCls ].time << std ::endl ;
148175 clsTuple -> Fill ((* clsArray )[iCls ].chipID ,
149176 globalClsCoords .x (),
150177 globalClsCoords .y (),
@@ -156,6 +183,34 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
156183 histXCoordCls -> Fill (globalClsCoords .x ());
157184 histYCoordCls -> Fill (globalClsCoords .y ());
158185 histZCoordCls -> Fill (globalClsCoords .z ());
186+
187+ // Match to digit
188+ auto digitLabelFromCls = (clsLabels -> getLabels (iCls ))[0 ];
189+ auto digitEntry = digitsLabels .find (digitLabelFromCls );
190+
191+ if (digitEntry == digitsLabels .end ()) {
192+ LOG (error ) << "No matching digit for cluster " << iCls << " with label " << digitLabelFromCls .getRawValue ();
193+ continue ;
194+ }
195+
196+ int iDigit = digitEntry -> second ;
197+ Int_t iRowFromDigit = (* digitsArray )[iDigit ].getRow ();
198+ Int_t iColFromDigit = (* digitsArray )[iDigit ].getColumn ();
199+ Int_t iChipIDFromDigit = (* digitsArray )[iDigit ].getChipIndex ();
200+ Int_t iSubDetIDFromDigit = tofGeo -> getIOTOFLayer (iChipIDFromDigit );
201+ Float_t timeFromDigit = (* digitsArray )[iDigit ].getTime ();
202+
203+ float xFromDigit = 0.f , yFromDigit = 0.f , zFromDigit = 0.f ;
204+ if (iSubDetIDFromDigit >= 0 ) {
205+ segGeom -> detectorToLocal (iRowFromDigit , iColFromDigit , xFromDigit , zFromDigit , iSubDetIDFromDigit );
206+ }
207+
208+ o2 ::math_utils ::Point3D < float > localDigitCoordFromDigit (xFromDigit , yFromDigit , zFromDigit ); // local Digit
209+ const auto globalDigitCoordFromDigit = tofGeo -> getMatrixL2G (iChipIDFromDigit )(localDigitCoordFromDigit ); // convert to global
210+ histXCoordRes -> Fill (globalClsCoords .x () - globalDigitCoordFromDigit .X ());
211+ histYCoordRes -> Fill (globalClsCoords .y () - globalDigitCoordFromDigit .Y ());
212+ histZCoordRes -> Fill (globalClsCoords .z () - globalDigitCoordFromDigit .Z ());
213+ histTimeRes -> Fill (time - timeFromDigit );
159214 } // end loop on clusters array
160215 } // end loop on ROFRecords
161216
@@ -185,6 +240,10 @@ void CheckClustersIOTOF(std::string digiFilePath = "tf3digits.root", std::string
185240 histXCoordDigit -> Write ();
186241 histYCoordDigit -> Write ();
187242 histZCoordDigit -> Write ();
243+ histXCoordRes -> Write ();
244+ histYCoordRes -> Write ();
245+ histZCoordRes -> Write ();
246+ histTimeRes -> Write ();
188247 outFile -> Write ();
189248 outFile -> Close ();
190249
0 commit comments