Skip to content

Commit d9193dd

Browse files
committed
ITS3: split topo dict into IB and OB
Signed-off-by: Felix Schlepper <felix.schlepper@cern.ch>
1 parent b7aa66d commit d9193dd

File tree

15 files changed

+498
-307
lines changed

15 files changed

+498
-307
lines changed

Detectors/Upgrades/ITS3/base/include/ITS3Base/SpecsV2.h

Lines changed: 0 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -135,7 +135,6 @@ constexpr std::array<double, nLayers> radiiInner{radii[0] - silicon::thicknessIn
135135
constexpr std::array<double, nLayers> radiiOuter{radii[0] + silicon::thicknessOut, radii[1] + silicon::thicknessOut, radii[2] + silicon::thicknessOut}; // outer silicon radius
136136
constexpr std::array<double, nLayers> radiiMiddle{(radiiInner[0] + radiiOuter[0]) / 2., (radiiInner[1] + radiiOuter[1]) / 2., (radiiInner[2] + radiiOuter[2]) / 2.}; // middle silicon radius
137137
constexpr double nominalYShift{-metalstack::thickness / 2.}; // shift to position in silicion volume to the chip volume (silicon+metalstack)
138-
/*constexpr double nominalYShift{0}; // shift to position in silicion volume to the chip volume (silicon+metalstack)*/
139138

140139
// extra information of pixels and their response functions
141140
namespace pixelarray::pixels

Detectors/Upgrades/ITS3/macros/test/CheckClusterSize.C

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -255,13 +255,14 @@ void CheckClusterSize(std::string clusFileName = "o2clus_its.root",
255255

256256
auto pattId = cluster.getPatternID();
257257
auto id = cluster.getSensorID();
258+
auto ib = o2::its3::constants::detID::isDetITS3(id);
258259
int clusterSize{-1};
259-
if (pattId == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattId)) {
260+
if (pattId == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattId, ib)) {
260261
o2::itsmft::ClusterPattern patt(pattIt);
261262
clusterSize = patt.getNPixels();
262263
continue;
263264
} else {
264-
clusterSize = dict.getNpixels(pattId);
265+
clusterSize = dict.getNpixels(pattId, ib);
265266
}
266267

267268
const auto& label = (clusLabArr->getLabels(clEntry))[0];

Detectors/Upgrades/ITS3/macros/test/CheckClustersITS3.C

Lines changed: 16 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -42,7 +42,7 @@
4242
void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
4343
const std::string& hitfile = "o2sim_HitsIT3.root",
4444
const std::string& inputGeom = "",
45-
std::string dictfile = "../ccdb/IT3/Calib/ClusterDictionary/snapshot.root",
45+
std::string dictfile = "./ccdb/IT3/Calib/ClusterDictionary/snapshot.root",
4646
bool batch = false)
4747
{
4848
gROOT->SetBatch(batch);
@@ -63,7 +63,8 @@ void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
6363
std::vector<HitVec*> hitVecPool;
6464
std::vector<MC2HITS_map> mc2hitVec;
6565

66-
ULong_t cPattValid{0}, cPattInvalid{0}, cLabelInvalid{0}, cNoMC{0};
66+
ULong_t cPattValidIB{0}, cPattInvalidIB{0}, cLabelInvalidIB{0}, cNoMCIB{0};
67+
ULong_t cPattValidOB{0}, cPattInvalidOB{0}, cLabelInvalidOB{0}, cNoMCOB{0};
6768

6869
TFile fout("CheckClusters.root", "recreate");
6970
TNtuple nt("ntc", "cluster ntuple", "ev:lab:hlx:hlz:hgx:hgz:tx:tz:cgx:cgy:cgz:clx:cly:clz:dx:dy:dz:ex:ez:patid:rof:npx:id:eta:row:col:lay");
@@ -103,6 +104,7 @@ void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
103104
} else {
104105
LOG(info) << "Running without dictionary !";
105106
}
107+
dict.print();
106108

107109
// ROFrecords
108110
std::vector<ROFRec> rofRecVec, *rofRecVecP = &rofRecVec;
@@ -175,28 +177,26 @@ void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
175177
auto isIB = o2::its3::constants::detID::isDetITS3(chipID);
176178
auto layer = o2::its3::constants::detID::getDetID2Layer(chipID);
177179
auto clusterSize{-1};
178-
if (pattID == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattID)) {
180+
if (pattID == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattID, isIB)) {
179181
o2::itsmft::ClusterPattern patt(pattIt);
180182
locC = dict.getClusterCoordinates(cluster, patt, false);
181183
LOGP(debug, "I am invalid and I am on chip {}", chipID);
182-
++cPattInvalid;
184+
(isIB) ? ++cPattInvalidIB : ++cPattInvalidOB;
183185
continue;
184186
} else {
185187
locC = dict.getClusterCoordinates(cluster);
186-
errX = dict.getErrX(pattID);
187-
errZ = dict.getErrZ(pattID);
188-
errX *= (isIB) ? MosaixSegmentation::PitchRow : Segmentation::PitchRow;
189-
errZ *= (isIB) ? MosaixSegmentation::PitchCol : Segmentation::PitchCol;
190-
npix = dict.getNpixels(pattID);
191-
++cPattValid;
188+
errX = dict.getErrX(pattID, isIB);
189+
errZ = dict.getErrZ(pattID, isIB);
190+
npix = dict.getNpixels(pattID, isIB);
191+
(isIB) ? ++cPattValidIB : ++cPattValidOB;
192192
}
193193

194194
// Transformation to the local --> global
195195
auto gloC = gman->getMatrixL2G(chipID)(locC);
196196
const auto& lab = (clusLabArr->getLabels(clEntry))[0];
197197

198198
if (!lab.isValid()) {
199-
++cLabelInvalid;
199+
(isIB) ? ++cLabelInvalidIB : ++cLabelInvalidOB;
200200
continue;
201201
}
202202

@@ -208,7 +208,7 @@ void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
208208
auto hitEntry = mc2hit.find(key);
209209
if (hitEntry == mc2hit.end()) {
210210
LOG(debug) << "Failed to find MC hit entry for Tr" << trID << " chipID" << chipID;
211-
++cNoMC;
211+
(isIB) ? ++cNoMCIB : ++cNoMCOB;
212212
continue;
213213
}
214214
const auto& hit = (*hitArray)[hitEntry->second];
@@ -263,8 +263,10 @@ void CheckClustersITS3(const std::string& clusfile = "o2clus_its.root",
263263
}
264264
}
265265

266-
LOGP(info, "There were {} valid PatternIDs and {} ({:.1f}%) invalid ones", cPattValid, cPattInvalid, ((float)cPattInvalid / (float)(cPattInvalid + cPattValid)) * 100);
267-
LOGP(info, "There were {} invalid Labels and {} with No MC Hit information ", cLabelInvalid, cNoMC);
266+
LOGP(info, "IB {} valid PatternIDs and {} ({:.1f}%) invalid ones", cPattValidIB, cPattInvalidIB, ((float)cPattInvalidIB / (float)(cPattInvalidIB + cPattValidIB)) * 100);
267+
LOGP(info, "IB {} invalid Labels and {} with No MC Hit information ", cLabelInvalidIB, cNoMCIB);
268+
LOGP(info, "OB {} valid PatternIDs and {} ({:.1f}%) invalid ones", cPattValidOB, cPattInvalidOB, ((float)cPattInvalidOB / (float)(cPattInvalidOB + cPattValidOB)) * 100);
269+
LOGP(info, "OB {} invalid Labels and {} with No MC Hit information ", cLabelInvalidOB, cNoMCOB);
268270

269271
auto canvCgXCgY = new TCanvas("canvCgXCgY", "", 1600, 1600);
270272
canvCgXCgY->Divide(2, 2);

Detectors/Upgrades/ITS3/macros/test/CompareClustersAndDigits.C

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -260,7 +260,7 @@ void CompareClustersAndDigits(std::string clusfile = "o2clus_it3.root",
260260
const auto pattID = cluster.getPatternID();
261261
const auto isIB = o2::its3::constants::detID::isDetITS3(chipID);
262262
const auto layer = gman->getLayer(chipID);
263-
if (pattID == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattID)) {
263+
if (pattID == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattID, isIB)) {
264264
continue;
265265
}
266266
const auto& lab = (clusLabArr->getLabels(clEntry))[0];
@@ -316,7 +316,7 @@ void CompareClustersAndDigits(std::string clusfile = "o2clus_it3.root",
316316
data[chipID].cog->AddPoint(colC, rowC);
317317

318318
constexpr float delta = 1e-2;
319-
const auto& patt = dict.getPattern(cluster.getPatternID());
319+
const auto& patt = dict.getPattern(cluster.getPatternID(), isIB);
320320
auto box = new TBox(
321321
cluster.getCol() - delta - 0.5,
322322
cluster.getRow() - delta - 0.5,

Detectors/Upgrades/ITS3/macros/test/CreateDictionariesITS3.C

Lines changed: 68 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -60,7 +60,7 @@ void CreateDictionariesITS3(bool saveDeltas = false,
6060
std::string collContextfile = "collisioncontext.root",
6161
std::string inputGeom = "",
6262
float checkOutliers = 2., // reject outliers (MC dX or dZ exceeds row/col span by a factor above the threshold)
63-
float minPtMC = 0.01) // account only MC hits with pT above threshold
63+
float minPtMC = 0.1) // account only MC hits with pT above threshold
6464
{
6565
const int QEDSourceID = 99; // Clusters from this MC source correspond to QED electrons
6666

@@ -84,10 +84,11 @@ void CreateDictionariesITS3(bool saveDeltas = false,
8484
std::array<o2::its3::SegmentationMosaix, 3> mMosaixSegmentations{0, 1, 2};
8585
if (!clusDictFile.empty()) {
8686
clusDictOld.readFromFile(clusDictFile);
87-
LOGP(info, "Loaded external cluster dictionary with {} entries from {}", clusDictOld.getSize(), clusDictFile);
87+
LOGP(info, "Loaded external cluster dictionary with {} IB/{} OBentries from {}", clusDictOld.getSize(true), clusDictOld.getSize(false), clusDictFile);
8888
}
8989

90-
ULong_t cOk{0}, cOutliers{0}, cFailedMC{0};
90+
ULong_t cOkIB{0}, cOutliersIB{0}, cFailedMCIB{0};
91+
ULong_t cOkOB{0}, cOutliersOB{0}, cFailedMCOB{0};
9192

9293
TFile* fout = nullptr;
9394
TNtuple* nt = nullptr;
@@ -233,17 +234,18 @@ void CreateDictionariesITS3(bool saveDeltas = false,
233234

234235
const auto& cluster = (*clusArr)[clEntry];
235236
o2::itsmft::ClusterPattern pattern;
237+
bool ib = o2::its3::constants::detID::isDetITS3(cluster.getChipID());
236238

237239
if (cluster.getPatternID() != CompCluster::InvalidPatternID) {
238-
if (clusDictOld.getSize() == 0) {
240+
if (clusDictOld.getSize(ib) == 0) {
239241
LOG(error) << "Encountered patternID = " << cluster.getPatternID() << " != " << CompCluster::InvalidPatternID;
240242
LOG(error) << "Clusters have already been generated with a dictionary which was not provided";
241243
return;
242244
}
243-
if (clusDictOld.isGroup(cluster.getPatternID())) {
245+
if (clusDictOld.isGroup(cluster.getPatternID(), ib)) {
244246
pattern.acquirePattern(pattIdx);
245247
} else {
246-
pattern = clusDictOld.getPattern(cluster.getPatternID());
248+
pattern = clusDictOld.getPattern(cluster.getPatternID(), ib);
247249
}
248250
} else {
249251
pattern.acquirePattern(pattIdx);
@@ -270,9 +272,8 @@ void CreateDictionariesITS3(bool saveDeltas = false,
270272
o2::math_utils::Vector3D<float> xyzLocM;
271273
xyzLocM.SetCoordinates(0.5f * (xyzLocE.X() + xyzLocS.X()), 0.5f * (xyzLocE.Y() + xyzLocS.Y()), 0.5f * (xyzLocE.Z() + xyzLocS.Z()));
272274
auto locC = o2::its3::TopologyDictionary::getClusterCoordinates(cluster, pattern, false);
273-
bool isIB = o2::its3::constants::detID::isDetITS3(chipID);
274275
int layer = gman->getLayer(chipID);
275-
if (isIB) {
276+
if (ib) {
276277
float xFlat{0.}, yFlat{0.};
277278
mMosaixSegmentations[layer].curvedToFlat(xyzLocM.X(), xyzLocM.Y(), xFlat, yFlat);
278279
xyzLocM.SetCoordinates(xFlat, yFlat, xyzLocM.Z());
@@ -281,33 +282,33 @@ void CreateDictionariesITS3(bool saveDeltas = false,
281282
}
282283
dX = xyzLocM.X() - locC.X();
283284
dZ = xyzLocM.Z() - locC.Z();
284-
dX /= (isIB) ? o2::its3::SegmentationMosaix::PitchRow : o2::itsmft::SegmentationAlpide::PitchRow;
285-
dZ /= (isIB) ? o2::its3::SegmentationMosaix::PitchCol : o2::itsmft::SegmentationAlpide::PitchCol;
285+
dX /= (ib) ? o2::its3::SegmentationMosaix::PitchRow : o2::itsmft::SegmentationAlpide::PitchRow;
286+
dZ /= (ib) ? o2::its3::SegmentationMosaix::PitchCol : o2::itsmft::SegmentationAlpide::PitchCol;
286287
if (saveDeltas) {
287288
nt->Fill(topology.getHash(), dX, dZ);
288289
}
289290
if (checkOutliers > 0.) {
290291
if (bool bX = std::abs(dX) > topology.getRowSpan() * checkOutliers, bZ = std::abs(dZ) > topology.getColumnSpan() * checkOutliers; bX || bZ) { // ignore outlier
291-
++cOutliers;
292+
(ib) ? ++cOutliersIB : ++cOutliersOB;
292293
LOGP(debug, "Ignored Value dX={} > {} * {} -> {}", dX, topology.getRowSpan(), checkOutliers, bX);
293294
LOGP(debug, "Ignored Value dZ={} > {} * {} -> {}", dZ, topology.getColumnSpan(), checkOutliers, bZ);
294295
dX = dZ = BuildTopologyDictionary::IgnoreVal;
295296
} else {
296-
++cOk;
297+
(ib) ? ++cOkIB : ++cOkOB;
297298
}
298299
}
299300
}
300301
} else {
301302
/* LOGP(info, " Failed to find MC hit entry for Tr: {} chipID: {}", trID, chipID); */
302303
/* lab.print(); */
303-
++cFailedMC;
304+
(ib) ? ++cFailedMCIB : ++cFailedMCOB;
304305
}
305-
signalDictionary.accountTopology(topology, dX, dZ);
306+
signalDictionary.accountTopology(topology, ib, dX, dZ);
306307
} else {
307-
noiseDictionary.accountTopology(topology, dX, dZ);
308+
noiseDictionary.accountTopology(topology, ib, dX, dZ);
308309
}
309310
}
310-
completeDictionary.accountTopology(topology, dX, dZ);
311+
completeDictionary.accountTopology(topology, ib, dX, dZ);
311312
}
312313

313314
// clean MC cache for events which are not needed anymore
@@ -323,64 +324,92 @@ void CreateDictionariesITS3(bool saveDeltas = false,
323324
}
324325
}
325326

326-
LOGP(info, "Clusters: {} okay (failed MCHit2Clus {}); outliers {}", cOk, cFailedMC, cOutliers);
327+
LOGP(info, "IB Clusters: {} okay (failed MCHit2Clus {}); outliers {}", cOkIB, cFailedMCIB, cOutliersIB);
328+
LOGP(info, "OB Clusters: {} okay (failed MCHit2Clus {}); outliers {}", cOkOB, cFailedMCOB, cOutliersOB);
327329

328330
auto dID = o2::detectors::DetID::IT3;
329331

330332
LOGP(info, "Complete Dictionary:");
331-
completeDictionary.setThreshold(probThreshold);
333+
completeDictionary.setThreshold(probThreshold, true);
334+
completeDictionary.setThreshold(probThreshold, false);
332335
completeDictionary.groupRareTopologies();
333336
completeDictionary.printDictionaryBinary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, ""));
334337
completeDictionary.printDictionary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "", "txt"));
335338
completeDictionary.saveDictionaryRoot(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "", "root"));
336339

337340
TFile histogramOutput("histograms.root", "recreate");
338341
TCanvas* cComplete = new TCanvas("cComplete", "Distribution of all the topologies");
339-
cComplete->cd();
340-
cComplete->SetLogy();
341-
TH1F* hComplete = completeDictionary.getDictionary().getTopologyDistribution("hComplete");
342-
hComplete->SetDirectory(nullptr);
343-
hComplete->Draw("hist");
344-
hComplete->Write();
342+
cComplete->Divide(2, 1);
343+
cComplete->cd(1);
344+
TH1F* hCompleteIB = completeDictionary.getDictionary().getTopologyDistribution("hCompleteInnerBarrel", true);
345+
hCompleteIB->SetDirectory(nullptr);
346+
hCompleteIB->Draw("hist");
347+
gPad->SetLogy();
348+
cComplete->cd(2);
349+
TH1F* hCompleteOB = completeDictionary.getDictionary().getTopologyDistribution("hCompleteOuterBarrel", false);
350+
hCompleteOB->SetDirectory(nullptr);
351+
hCompleteOB->Draw("hist");
352+
gPad->SetLogy();
353+
histogramOutput.cd();
354+
hCompleteIB->Write();
355+
hCompleteOB->Write();
345356
cComplete->Write();
346357

347358
if (clusLabArr) {
348359
LOGP(info, "Noise Dictionary:");
349-
noiseDictionary.setThreshold(0.0001);
360+
noiseDictionary.setThreshold(0.0001, true);
361+
noiseDictionary.setThreshold(0.0001, false);
350362
noiseDictionary.groupRareTopologies();
351363
noiseDictionary.printDictionaryBinary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "noiseClusTopo"));
352364
noiseDictionary.printDictionary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "noiseClusTopo", "txt"));
353365
noiseDictionary.saveDictionaryRoot(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "noiseClusTopo", "root"));
354366

355367
LOGP(info, "Signal Dictionary:");
356-
signalDictionary.setThreshold(0.0001);
368+
signalDictionary.setThreshold(0.0001, true);
369+
signalDictionary.setThreshold(0.0001, false);
357370
signalDictionary.groupRareTopologies();
358371
signalDictionary.printDictionaryBinary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "signal"));
359372
signalDictionary.printDictionary(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "signal", "txt"));
360373
signalDictionary.saveDictionaryRoot(o2::base::DetectorNameConf::getAlpideClusterDictionaryFileName(dID, "signal", "root"));
361374

362375
LOGP(info, "Plotting Channels");
363376
auto cNoise = new TCanvas("cNoise", "Distribution of noise topologies");
364-
cNoise->cd();
365-
cNoise->SetLogy();
366-
auto hNoise = noiseDictionary.getDictionary().getTopologyDistribution("hNoise");
367-
hNoise->SetDirectory(nullptr);
368-
hNoise->Draw("hist");
377+
cNoise->Divide(2, 1);
378+
cNoise->cd(1);
379+
auto hNoiseIB = noiseDictionary.getDictionary().getTopologyDistribution("hNoiseInnerBarrel", true);
380+
hNoiseIB->SetDirectory(nullptr);
381+
hNoiseIB->Draw("hist");
382+
gPad->SetLogy();
383+
cNoise->cd(2);
384+
auto hNoiseOB = noiseDictionary.getDictionary().getTopologyDistribution("hNoiseOuterBarrel", false);
385+
hNoiseOB->SetDirectory(nullptr);
386+
hNoiseOB->Draw("hist");
387+
gPad->SetLogy();
369388
histogramOutput.cd();
370-
hNoise->Write();
389+
hNoiseIB->Write();
390+
hNoiseOB->Write();
371391
cNoise->Write();
392+
372393
auto cSignal = new TCanvas("cSignal", "cSignal");
373-
cSignal->cd();
394+
cSignal->Divide(2, 1);
395+
cSignal->cd(1);
396+
auto hSignalIB = signalDictionary.getDictionary().getTopologyDistribution("hSignalInnerBarrel", true);
397+
hSignalIB->SetDirectory(nullptr);
398+
hSignalIB->Draw("hist");
399+
gPad->SetLogy();
400+
cSignal->cd(2);
374401
cSignal->SetLogy();
375-
auto hSignal = signalDictionary.getDictionary().getTopologyDistribution("hSignal");
376-
hSignal->SetDirectory(nullptr);
377-
hSignal->Draw("hist");
402+
auto hSignalOB = signalDictionary.getDictionary().getTopologyDistribution("hSignalOuterBarrel", false);
403+
hSignalOB->SetDirectory(nullptr);
404+
hSignalOB->Draw("hist");
405+
gPad->SetLogy();
378406
histogramOutput.cd();
379-
hSignal->Write();
407+
hSignalIB->Write();
408+
hSignalOB->Write();
380409
cSignal->Write();
381-
sw.Stop();
382-
sw.Print();
383410
}
411+
sw.Stop();
412+
sw.Print();
384413
if (saveDeltas) {
385414
fout->cd();
386415
nt->Write();

0 commit comments

Comments
 (0)