ATLAS Offline Software
Loading...
Searching...
No Matches
CalibrateIBL.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2026 CERN for the benefit of the ATLAS collaboration
3*/
4
5//======================================================================
6// PixelCalibIBL : Read calibration file, redo the fit.
7//
8// Note : this program only runs on IBL.
9//
10//
11// Mapping convention for FE in the calibration.
12//
13// ^
14// phi|
15// |
16// 336| sFE0 sFE1
17// |
18// 0+------------------------------>
19// 0 80 160
20// eta
21//======================================================================
22
23#include "CxxUtils/checker_macros.h"
24ATLAS_NO_CHECK_FILE_THREAD_SAFETY;
25
26
27// ignore the "demonstrate Xcheck" lines unless necessary
28// #define DEMOXCHECK true
29
30#include "TFile.h"
31#include "TDirectory.h"
32#include "TDirectoryFile.h"
33#include "TKey.h"
34#include "TH2D.h"
35#include "TH2F.h"
36#include "TProfile.h"
37#include "TGraph.h"
38#include "TGraphErrors.h"
39#include "TMultiGraph.h"
40#include "TMath.h"
41#include "TF1.h"
42#include "TStyle.h"
43#include "TVectorT.h"
44#include "THStack.h"
45#include "TRandom3.h"
46#include "ChargeCalibration/common/PixelMapping.h"
47#include "PathResolver/PathResolver.h"
48#include <cmath>
49#include <fstream>
50#include <iostream>
51#include <iterator>
52#include <string>
53#include <sstream>
54#include <utility>
55#include <map>
56#include <array>
57
58using namespace std;
59using pix::PixelMapping;
60
61// Lianyou: simplely capacitor charging curve
62double funcDisp(double *x, double *par)
63{
64 double ret = 9.9e10;
65 ret = par[0] * 0.1 * (1. - exp(-(x[0] - 1400.) * (x[0] - 1400.) / (1000. * par[1])) + par[2] * 0.0001 * sqrt(x[0]));
66 return ret;
67}
68
69double ChrgfromToT(float tot, double par[10])
70{
71 int lr = (tot <= 6. ? 0 : 5);
72 double chrg = 1400.;
73
74 double num = par[0 + lr] + tot * (par[1 + lr] + tot * (par[2 + lr] + tot * par[3 + lr]));
75 double denom = 1. + tot * par[4 + lr];
76 if (denom != 0.)
77 chrg = num / denom;
78
79 return chrg;
80}
81
82double ToTfromChrg(float chrg, double par[5])
83{
84 double tot = 0.;
85 double num = par[0] + chrg * (par[1] + chrg * (par[2] + chrg * par[3]));
86 double denom = 1. + chrg * par[4];
87 if (denom != 0.0)
88 tot = num / denom;
89 return tot;
90}
91
92double funcRation5(double *x, double *par)
93{
94 double ret = 9.9e10;
95 double num = par[0] + x[0] * (par[1] + x[0] * (par[2] + x[0] * par[3]));
96 double denom = 1. + x[0] * par[4];
97 if (denom != 0.0)
98 ret = num / denom;
99 return ret;
100}
101
102double funcDispConvoluted(double *x, double *par)
103{
104 double ret = 9.9e09;
105 ret = -par[0] * (x[0] - par[1]) + par[2];
106 return ret;
107}
108
109int iblCalib(const std::string& InDir, const std::string& THRscan, const std::string& scan, const std::string& scanLowCharge)
110{
111 // creating the object for the pixel mapping
112 PixelMapping pixmap(PathResolver::find_file("PixelCalibAlgs/mapping.csv", "DATAPATH"));
113
114 constexpr bool run3 = true; // new series of injection charges (22) and new format
115 const float badRDfracCut = 0.2; // Cut on the fraction of bad read out per frontend
116 float HDCshift = 2;
117
118 // Shall we try this later ?
119 constexpr int finerToT = 1; // smaller bin size for ToT if doing fit
120 constexpr bool moreFE = false;
121 constexpr int npsFEs = moreFE ? 8 : 2;
122
123 gStyle->SetOptFit(1100);
124
125#if defined(DEMOXCHECK)
126
127 std::cout << " DEMOXCHECK will run ... " << std::endl;
128
129 // For debugging a small area of pixels will print their readout and injection
130 int XcheckCharge[2] = {13, 14}, XcheckToT = 10;
131 int XcheckPhi[2] = {57, 63}, XcheckEta[2] = {23, 29};
132 std::string XcheckModule = "LI_S11_A_M3_A6";
133
134 // fill some histograms for monitoring
135 int ToTfill = 3 * finerToT, CHRGfill = 3;
136 const int HDC = 2; // specifying output
137
138#endif
139
140 Double_t tmpRCf_H[4] = {0.868295, 0.67523, 0.732279, 0.678458}; // initial values from history calibration
141 Double_t tmpRCf_T[4] = {1.02844, 1.10788, 1.07204, 1.11015}; // initial values from history calibration
142 Double_t reverseCF_H[4], reverseCF_T[4];
143 for (int t = 0; t < 4; t++)
144 {
145 reverseCF_H[t] = tmpRCf_H[t];
146 reverseCF_T[t] = tmpRCf_T[t];
147 }
148
149 Double_t RDentries[33] = {0, 10494.8, 14605.6, 62952.6, 8663.57, 17320.2, 20417, 49508.9,
150 21646.8, 31503.9, 14955.6, 16290.7, 11617.5, 16461.2, 11980.5, 17210.5,
151 12773.7, 18574, 13930.9, 19777, 14629.4, 20619.5, 15191.6, 20883,
152 15411.1, 19198, 12097.4, 12109.7, 8881.96, 17235.1, 20802.3, 11265.8,
153 0.};
154
155 std::vector<int> reversedModules = {182, 211, 213, 229, 223, 231, 246, 280, 295, 350, 357, 360, 371, 403, 416, 423};
156
157 // outputs : something@somewhere as one like ...
158 std::string Outdir = "./";
159 string spec = "ToTbin" + to_string(finerToT) + "_FrtEnd" + to_string(npsFEs) + "_";
160
161#if defined(DEMOXCHECK)
162 spec += "DEMOXCHECK_";
163#endif
164
165 std::string StrFileName = spec + scan;
166 std::string rootFileName = Outdir + "/TotChargeCalib_" + StrFileName + ".root";
167 std::string logFileName = Outdir + "/ChargeCalib_" + StrFileName + ".log";
168 std::string dbFileName = Outdir + "/ChargeCalib_" + StrFileName + ".TXT";
169 if (!run3)
170 {
171 rootFileName = Outdir + "/TotChargeCalib_run2_" + StrFileName + ".root";
172 logFileName = Outdir + "/ChargeCalib_run2_" + StrFileName + ".log";
173 dbFileName = Outdir + "/ChargeCalib_run2_" + StrFileName + ".TXT";
174 }
175
176 // fetch pictures from directory : ToT/ROD_L1_Sxx/LI_Sy_AC_Mz_ACk
177 ofstream logout(logFileName);
178 ofstream txtDB(dbFileName);
179
180 std::string inThrFile = "";
181 std::string inTotFile = "";
182 std::string rodPath = "";
183 std::string inTotFileAux = "";
184 std::string rodPathAux = "";
185
186 std::cout << " Running IBL calibration analysis ... " << endl;
187
188 inThrFile = InDir + THRscan + ".root";
189 inTotFile = InDir + scan + ".root";
190 rodPath = inTotFile + ":/" + scan.substr(scan.find("_")+1);
191 inTotFileAux = InDir + scanLowCharge + ".root";
192 rodPathAux = inTotFileAux + ":/" + scanLowCharge.substr(scanLowCharge.find("_")+1);
193
194 // historical record for run2 fomat
195 if (!run3)
196 {
197 inThrFile = "calib2018/IBL/THR_SCAN_S000083686.root";
198 inTotFile = "calib2018/IBL/TOT_SCAN_S000083690.root";
199 rodPath = inTotFile + ":/S000083690";
200 }
201
202 // selecting Q threshold
203 int nrow = 336, ncol = 160; // y-axis, x-axis respectively
204
205 const Int_t numChrgs = 22;
206 Double_t IBLchrgs[numChrgs] = {1400., 1500., 1750., 2000., 2500., 3000., 3500., 4000., 5000., 6000., 8000., 10000.,
207 12000., 14000., 16000., 18000., 20000., 22000., 24000., 26000., 28000., 30000.};
208
209 constexpr Int_t nchargeIBL = (run3 ? 22 : 19); //run3 is constexpr
210
211 int skip = numChrgs - nchargeIBL;
212 Double_t chrgAbaciIBL[nchargeIBL], chargeErrArrIBL[nchargeIBL];
213 for (int c = 0; c < nchargeIBL; c++)
214 {
215 chrgAbaciIBL[c] = IBLchrgs[c + skip];
216 chargeErrArrIBL[c] = 0.;
217 }
218
219 Double_t chrgsbins[nchargeIBL + 1];
220 for (int c = 1; c < nchargeIBL; c++)
221 chrgsbins[c] = 0.5 * (chrgAbaciIBL[c - 1] + chrgAbaciIBL[c]);
222
223 chrgsbins[0] = chrgAbaciIBL[0] - 0.5 * (chrgAbaciIBL[1] - chrgAbaciIBL[0]);
224 chrgsbins[nchargeIBL] = chrgAbaciIBL[nchargeIBL - 1] + 0.5 * (chrgAbaciIBL[nchargeIBL - 1] - chrgAbaciIBL[nchargeIBL - 2]);
225
226 // please note extra +1 for ending bins.
227 constexpr Int_t nToTibl = 16 * finerToT + 1;
228
229 Double_t totAbaci[nToTibl], totbins[nToTibl + 1];
230
231 totbins[0] = -0.25;
232 totbins[1] = 0.5;
233 totAbaci[0] = 0.0;
234 totAbaci[1] = 1.;
235 for (int t = 2; t < nToTibl; t++)
236 {
237 totAbaci[t] = totAbaci[1] + (t - 1) / (1. * finerToT);
238 totbins[t] = 0.5 * (totAbaci[t] + totAbaci[t - 1]);
239 }
240 totbins[nToTibl] = totAbaci[nToTibl - 1] + 0.5 * (totbins[nToTibl - 1] - totbins[nToTibl - 2]);
241
242 Double_t totErrArr[nToTibl];
243 for (int t = 0; t < nToTibl; t++)
244 totErrArr[t] = 0.;
245
246#if defined(DEMOXCHECK)
247 for (int t = 0; t < nToTibl; t++)
248 logout << " totAbaci : " << t << " " << totbins[t] << " < " << totAbaci[t] << " > " << totbins[t + 1] << endl;
249
250 if (finerToT == 1)
251 {
252 for (int t = 1; t < nToTibl; t++)
253 {
254 RDentries[t] = RDentries[2 * t - 1] + RDentries[2 * t];
255 logout << " finer RDentries : " << t << " " << 2 * t - 1 << " " << 2 * t << " " << RDentries[t] << endl;
256 }
257 }
258#endif
259
260 map<string, map<string, map<string, float>>> pcdMap;
261
262 TFile roFile(rootFileName.c_str(), "RECREATE");
263 TDirectory *roThrDir = roFile.mkdir("Threshold");
264 TDirectory *roTotDir = roFile.mkdir("ToT");
265
266 // ****************************** Threshold IBL ******************************
267
268 Double_t THR_avg[2][npsFEs], ThrSig_avg[2][npsFEs];
269 array<std::unique_ptr<TH2D> , npsFEs> h2_Thr{};
270 array<std::unique_ptr<TH2D> , npsFEs> h2_ThrSig{};
271
272 std::multimap<float, std::string, std::greater<float>> badThr_Order;
273 static const std::string knowModule{"LI_S06_C_M1_C1"};
274 for (int sfe = 0; sfe < npsFEs; sfe++)
275 {
276 stringstream ss;
277 ss << sfe;
278
279 string idx = "I" + ss.str();
280
281 pcdMap[knowModule][idx]["ThrNorm"] = -42.;
282 pcdMap[knowModule][idx]["ThrRmsNorm"] = -42.;
283 pcdMap[knowModule][idx]["ThrSigNorm"] = -42.;
284 pcdMap[knowModule][idx]["ThrLong"] = -42.;
285 pcdMap[knowModule][idx]["ThrRmsLong"] = -42.;
286 pcdMap[knowModule][idx]["ThrSigLong"] = -42.;
287
288 for (int nl = 0; nl < 2; nl++)
289 THR_avg[nl][sfe] = ThrSig_avg[nl][sfe] = -99.;
290
291 idx = "Threshold" + ss.str();
292 h2_Thr[sfe] = std::make_unique<TH2D>(idx.c_str(), " ", 2, 0, 2, 200, 0, 5000);
293 idx = "ThresholdSig" + ss.str();
294 h2_ThrSig[sfe] = std::make_unique<TH2D>(idx.c_str(), " ", 2, 0, 2, 200, 0, 500);
295 }
296
297 if (inThrFile.size() > 0)
298 {
299 std::cout << endl
300 << "INFO =>> [IBL] threshold scan analysis..." << endl;
301
302 TFile riThrFile(inThrFile.c_str(), "READ");
303 std::string chi2HistName = "SCURVE_CHI2";
304 std::string thrHistName = "SCURVE_MEAN";
305 std::string sigHistName = "SCURVE_SIGMA";
306
307 std::unique_ptr<TH1F> h1dChi2 = std::make_unique<TH1F>("h1dChi2", "", 200, 0, 1);
308 std::unique_ptr<TH1F> h1dThr = std::make_unique<TH1F>("h1dThr" , "", 200, 0, 5000);
309 std::unique_ptr<TH1F> h1dSig = std::make_unique<TH1F>("h1dSig" , "", 200, 0, 500);
310
311 TDirectoryFile *scanDir = (TDirectoryFile *)((TKey *)riThrFile.GetListOfKeys()->First())->ReadObj();
312 TList *rodKeyList = (TList *)scanDir->GetListOfKeys();
313 TIter rodItr(rodKeyList);
314 TKey *rodKey;
315
316 while ((rodKey = (TKey *)rodItr()))
317 {
318 std::string rodName(rodKey->GetName());
319 TDirectoryFile *rodDir = (TDirectoryFile *)rodKey->ReadObj();
320 TList *modKeyList = (TList *)rodDir->GetListOfKeys();
321 TIter modItr(modKeyList);
322 TKey *modKey;
323
324 TDirectory *dirRod = roThrDir->mkdir(rodName.c_str());
325
326 while ((modKey = (TKey *)modItr()))
327 {
328 std::string modName(modKey->GetName());
329 string modStr(modKey->GetName());
330 TDirectory *dirMod = dirRod->mkdir(modName.c_str());
331
332 std::string chi2HistDirPath = modName + "/" + chi2HistName + "/A0/B0";
333
334 TDirectory *chi2HistDir = (TDirectory *)rodDir->Get(chi2HistDirPath.c_str());
335 if (chi2HistDir == NULL)
336 {
337 cout << " Warning : NULL dir " << endl;
338 continue;
339 }
340
341 std::unique_ptr<TH2D> h2dChi2 (static_cast<TH2D*>(static_cast<TKey*>(chi2HistDir->GetListOfKeys()->First())->ReadObj()));
342 if (h2dChi2 == NULL)
343 {
344 cout << " Warning : NULL dir " << endl;
345 continue;
346 }
347 std::string thrHistDirPath = modName + "/" + thrHistName + "/A0/B0";
348 TDirectoryFile *thrHistDir = (TDirectoryFile *)rodDir->Get(thrHistDirPath.c_str());
349 std::unique_ptr<TH2D> h2dThr(static_cast<TH2D*>(static_cast<TKey*>(thrHistDir->GetListOfKeys()->First())->ReadObj()));
350 std::string sigHistDirPath = modName + "/" + sigHistName + "/A0/B0";
351
352 TDirectoryFile *sigHistDir = (TDirectoryFile *)rodDir->Get(sigHistDirPath.c_str());
353 std::unique_ptr<TH2D> h2dSig (static_cast<TH2D*>(static_cast<TKey*>(sigHistDir->GetListOfKeys()->First())->ReadObj()));
354
355 array<std::unique_ptr<TH1F>, npsFEs> h1_ThrNorm{};
356 array<std::unique_ptr<TH1F>, npsFEs> h1_ThrSigNorm{};
357 array<std::unique_ptr<TH1F>, npsFEs> h1_ThrLong{};
358 array<std::unique_ptr<TH1F>, npsFEs> h1_ThrSigLong{};
359 array<float, npsFEs> IlledThr{};
360
361 for (int sfe = 0; sfe < npsFEs; sfe++)
362 {
363 stringstream ss;
364 ss << sfe;
365
366 string hname = modStr + "ThrNorm" + ss.str();
367 h1_ThrNorm[sfe] = std::make_unique<TH1F>(hname.c_str(), "", 200, 0, 5000);
368 hname = modStr + "ThrSigNorm" + ss.str();
369 h1_ThrSigNorm[sfe] = std::make_unique<TH1F>(hname.c_str(), "", 200, 0, 500);
370 hname = modStr + "ThrLong" + ss.str();
371 h1_ThrLong[sfe] = std::make_unique<TH1F>(hname.c_str(), "", 200, 0, 5000);
372 hname = modStr + "ThrSigLong" + ss.str();
373 h1_ThrSigLong[sfe] = std::make_unique<TH1F>(hname.c_str(), "", 200, 0, 500);
374 }
375
376 for (int col = 1; col <= ncol; col++)
377 {
378 for (int row = 1; row <= nrow; row++)
379 {
380 float chi2 = h2dChi2->GetBinContent(col, row);
381 float thr = h2dThr->GetBinContent(col, row);
382 float sig = h2dSig->GetBinContent(col, row);
383
384 h1dChi2->Fill(chi2);
385 h1dThr->Fill(thr);
386 h1dSig->Fill(sig);
387 // Shall we measure/count the goodness of THR scan here ?
388 bool filled = true;
389 if (thr == 0 || thr > 10000 || sig == 0 || sig > 1000 || chi2 > 0.5 || chi2 <= 0)
390 {
391 filled = false;
392 }
393 int circ = -1;
394 if (filled)
395 {
396 if (col == 1 || col == ncol / 2 || col == ncol / 2 + 1 || col == ncol)
397 {
398 if (moreFE)
399 {
400 circ = (int)(row / 84.);
401 if (col > 80)
402 {
403 circ = 7 - circ;
404 }
405 }
406 else
407 {
408 if (col <= ncol / 2)
409 circ = 0;
410 else
411 circ = 1;
412 }
413
414 h1_ThrLong[circ]->Fill(thr);
415 h1_ThrSigLong[circ]->Fill(sig);
416 }
417 else
418 {
419 if (moreFE)
420 {
421 circ = (int)(row / 84.);
422 if (col > 80)
423 {
424 circ = 7 - circ;
425 }
426 }
427 else
428 {
429 if (col <= ncol / 2)
430 circ = 0;
431 else
432 circ = 1;
433 }
434
435 h1_ThrNorm[circ]->Fill(thr);
436 h1_ThrSigNorm[circ]->Fill(sig);
437 }
438 }
439 else
440 {
441 // Note: what is the intention? circ is -1 here
442 continue;
443 }
444 }
445 }
446 h2dSig.reset();
447 h2dThr.reset();
448 h2dChi2.reset();
449
450 for (int sfe = 0; sfe < npsFEs; sfe++)
451 {
452 stringstream ss;
453 ss << sfe;
454 string feName = "I" + ss.str();
455
456 float THRnorm = h1_ThrNorm[sfe]->GetMean();
457 bool valid = THRnorm > 100.;
458 pcdMap[modStr][feName]["ThrNorm"] = THRnorm;
459 pcdMap[modStr][feName]["ThrRmsNorm"] = h1_ThrNorm[sfe]->GetRMS();
460 float THRnormSig = h1_ThrSigNorm[sfe]->GetMean();
461 pcdMap[modStr][feName]["ThrSigNorm"] = THRnormSig;
462 float THRlong = h1_ThrLong[sfe]->GetMean();
463 pcdMap[modStr][feName]["ThrLong"] = THRlong;
464 pcdMap[modStr][feName]["ThrRmsLong"] = h1_ThrLong[sfe]->GetRMS();
465 float THRlongSig = h1_ThrSigLong[sfe]->GetMean();
466 pcdMap[modStr][feName]["ThrSigLong"] = THRlongSig;
467
468 float blank = 100. * IlledThr[sfe] * npsFEs / (1. * ncol * nrow);
469
470 if (blank > 0.01)
471 badThr_Order.insert(std::pair<float, std::string>(blank, (std::string)(modName + " : " + blank)));
472
473 if (valid)
474 {
475 h2_Thr[sfe]->Fill(0.5, THRnorm, 1);
476 h2_Thr[sfe]->Fill(1.5, THRlong, 1);
477 h2_ThrSig[sfe]->Fill(0.5, THRnormSig, 1);
478 h2_ThrSig[sfe]->Fill(1.5, THRlongSig, 1);
479 }
480 dirMod->WriteTObject(h1_ThrNorm[sfe].get());
481 dirMod->WriteTObject(h1_ThrLong[sfe].get());
482
483 h1_ThrNorm[sfe].reset();
484 h1_ThrSigNorm[sfe].reset();
485 h1_ThrLong[sfe].reset();
486 h1_ThrSigLong[sfe].reset();
487 }
488 }
489 }
490 roThrDir->WriteTObject(h1dChi2.get());
491 h1dChi2.reset();
492 roThrDir->WriteTObject(h1dThr.get());
493 h1dThr.reset();
494 roThrDir->WriteTObject(h1dSig.get());
495 h1dSig.reset();
496 }
497
498 for (int sfe = 0; sfe < npsFEs; sfe++)
499 {
500 stringstream ss;
501 ss << sfe;
502 string hName = "thrNorm" + ss.str();
503 THR_avg[0][sfe] = h2_Thr[sfe]->ProjectionY(hName.c_str(), 1, 1)->GetMean();
504 hName = "thrLong" + ss.str();
505 THR_avg[1][sfe] = h2_Thr[sfe]->ProjectionY(hName.c_str(), 2, 2)->GetMean();
506 hName = "thrSigNorm" + ss.str();
507 ThrSig_avg[0][sfe] = h2_ThrSig[sfe]->ProjectionY(hName.c_str(), 1, 1)->GetMean();
508 hName = "thrSigLong" + ss.str();
509 ThrSig_avg[1][sfe] = h2_ThrSig[sfe]->ProjectionY(hName.c_str(), 2, 2)->GetMean();
510 h2_Thr[sfe].reset();
511 h2_ThrSig[sfe].reset();
512 }
513 // ****************************** End - Threshold IBL ******************************
514
516 /************************* IBL ToT ***************************/
518
519 if (inTotFile.size() == 0 || inTotFileAux.size() == 0)
520 {
521 logout << " Missing ToT file from calib scan as input " << endl;
522 return 0;
523 }
524
525 std::cout << endl
526 << "INFO =>> [IBL] tot calib analysis..." << endl;
527
528 TFile riTotFile(inTotFile.c_str(), "READ");
529 TFile riTotFileAux(inTotFileAux.c_str(), "READ");
530
531 std::string totHistName = "TOT_MEAN";
532 std::string totSigHistName = "TOT_SIGMA";
533
534 TDirectoryFile *scanDir = (TDirectoryFile *)((TKey *)riTotFile.GetListOfKeys()->First())->ReadObj();
535 TList *rodKeyList = (TList *)scanDir->GetListOfKeys();
536 TIter rodItr(rodKeyList);
537 TKey *rodKey;
538
539 std::map<float, std::pair<vector<std::string>, vector<Double_t>>> ModuDataToPrint;
540
541#if defined(DEMOXCHECK)
542 vector<TH1F *> h1_ChrgEntry(nToTibl, nullptr);
543
544 vector<TH1F *> h1d_totSprdAll(nToTibl-1, nullptr);
545
546 for (int t = 0; t < nToTibl; t++)
547 {
548 stringstream tt;
549 tt << t;
550
551 string hname = "ChrgEntries_ToT_" + tt.str();
552 h1_ChrgEntry[t] = new TH1F(hname.c_str(), "RD entries ", 100, 0., 100000);
553
554 hname = "totSprdsFll_ToT_" + tt.str();
555 if (t < nToTibl - 1)
556 h1d_totSprdAll[t] = new TH1F(hname.c_str(), "ToT spread ", 36, 0.2, 1.1);
557 }
558#endif
559
560 constexpr Int_t totFE = 14 * 16 * npsFEs; // 16 modules on each of 14 ReadOutDisk
561 using Row = Double_t[nchargeIBL];
562 using Row2 = Double_t[nToTibl];
563 //avoid stack overflow, create on heap
564 std::unique_ptr<Row[]> TotArray{ new Row[totFE]{} };
565 std::unique_ptr<Row[]> TotErrArray{ new Row[totFE]{} };
566 std::unique_ptr<Row[]> TotSigArray { new Row[totFE]{} };
567 std::unique_ptr<Row[]> TotSigErrArray{ new Row[totFE]{} };
568 std::unique_ptr<Row2[]> ChrgArray{ new Row2[totFE]{} };
569 std::unique_ptr<Row[]> ChrgErrArray{ new Row[totFE]{} };
570
571 gRandom = new TRandom3(2203);
572 using Column = Double_t[totFE];
573 std::unique_ptr<Column[]> occuPhiEta{new Column[nchargeIBL]{}};
574
575 Int_t cntRod = 0;
576 std::map<float, std::string> devChrg_Order;
577 std::map<float, std::string> devToT_Order;
578
579 std::multimap<float, std::string, std::greater<float>> badModules_Order;
580 std::multimap<float, std::string, std::greater<float>> badModules_Order_detailed;
581
582 while ((rodKey = (TKey *)rodItr()))
583 {
584 std::string rodName(rodKey->GetName());
585 string rodStr(rodKey->GetName());
586 TDirectoryFile *rodDir = (TDirectoryFile *)rodKey->ReadObj();
587
588 TString path = rodDir->GetPath();
589 std::string pathAux = path.ReplaceAll(rodPath, rodPathAux).Data();
590
591 TDirectoryFile *rodDirAux = (TDirectoryFile *)riTotFileAux.Get(pathAux.c_str());
592 if (rodDirAux == NULL)
593 {
594 std::cout<< " Fail to get the rodPath in Aux: "<< pathAux << endl;
595 logout << " Fail to get the rodPath in Aux " << endl;
596 continue;
597 }
598
599 TDirectory *dirRod = roTotDir->mkdir(rodName.c_str());
600
601 TList *modKeyList = (TList *)rodDir->GetListOfKeys();
602 TIter modItr(modKeyList);
603
604 TKey *modKey;
605
606 int cntMod = 0;
607 string modNames[16 * npsFEs];
608
609 float occuChrgs[nToTibl][16 * npsFEs];
610
611#if defined(DEMOXCHECK)
612 std::string feName_maxDevChrg = "", feName_maxDevToT = "";
613 float maxDevChrg = -9., maxDevToT = -9., avgDevChrg = 0., avgDevToT = 0.;
614
615 vector<TH1F *> h1d_totSprd(nToTibl - 1, nullptr);
616
617 for (int t = 0; t < nToTibl - 1; t++)
618 {
619 stringstream tt;
620 tt << t + 1;
621
622 string prfmodname = "ToT_Sprd_" + tt.str();
623 h1d_totSprd[t] = new TH1F(prfmodname.c_str(), "ToT spread ", 28, 0.3, 1.0);
624 prfmodname = rodStr + "ToT spread @ ToT_" + tt.str();
625 h1d_totSprd[t]->SetTitle(prfmodname.c_str());
626 }
627#endif
628
629 while ((modKey = (TKey *)modItr()))
630 {
631 std::string modName(modKey->GetName());
632 string modStr(modKey->GetName());
633
634 TDirectory *dirMod = dirRod->mkdir(modName.c_str());
635 bool ibl3Dfe0 = false, ibl3Dfe1 = false;
636
637 int hashID = -1, hashIDL = -1, hashIDR = -1;
638 int bec = -1;
639 int layer = -1;
640 int phi_module = -1;
641 int eta_module = -1;
642 int eta_moduleL = -1, eta_moduleR = -1;
643 pixmap.mapping(modStr, &hashID, &bec, &layer, &phi_module, &eta_module);
644 if (hashID == -1)
645 {
646 pixmap.mapping(modStr + "_1", &hashID, &bec, &layer, &phi_module, &eta_module);
647 if (hashID != -1)
648 {
649 ibl3Dfe0 = true;
650 hashIDL = hashID;
651 eta_moduleL = eta_module;
652 }
653
654 pixmap.mapping(modStr + "_2", &hashID, &bec, &layer, &phi_module, &eta_module);
655 if (hashID != -1)
656 {
657 ibl3Dfe1 = true;
658 hashIDR = hashID;
659 eta_moduleR = eta_module;
660 }
661 }
662 if (hashID == -1){
663 std::cout<< "negative hash ID\n";
664 continue;
665 }
666
667
668 if (pcdMap.find(modStr) == pcdMap.end()){
669 std::cout<< "MOD not found "<< modStr << "\n";
670 continue;
671 }
672
673
674 // a stack of Charge-ToT data for ToT-charge converting
675 array<std::unique_ptr<TH2D>, npsFEs> h2d_XchrgYtot{};
676 array<std::unique_ptr<TH2D>, npsFEs> h2d_XchrgYToTSig{};
677
678 for (int sfe = 0; sfe < npsFEs; sfe++)
679 {
680 stringstream ss;
681 ss << sfe;
682
683 modNames[cntMod * npsFEs + sfe] = modStr + "-" + ss.str();
684
685 string prfmodname = modStr + "ToTvsChrg_FE" + ss.str();
686 h2d_XchrgYtot[sfe] = std::make_unique<TH2D>(prfmodname.c_str(), "ToT vs Charge", nchargeIBL, chrgsbins, nToTibl, totbins);
687
688 prfmodname = modStr + "ToT_Sig_vsChrg_FE" + ss.str();
689 h2d_XchrgYToTSig[sfe] = std::make_unique<TH2D>(prfmodname.c_str(), "ToT Sig vs Charge", nchargeIBL, chrgsbins, 100, 0., 1.);
690 }
691
692 // Fill ToT-Chrg for pixels per module
693 for (int c = 0; c < nchargeIBL; c++)
694 {
695 std::string totHistDirPath = modName + "/" + totHistName + "/A0/B0/C";
696 totHistDirPath += std::to_string(c);
697
698 std::unique_ptr<TH2F> h2dTot;
699 TDirectoryFile* totHistDir(static_cast<TDirectoryFile*> (rodDir->Get(totHistDirPath.c_str())));
700 if(!totHistDir){
701 std::cout<<" Missing totHistDir in : " << totHistDirPath << endl;
702 abort();
703 }
704 else {
705 h2dTot.reset(static_cast<TH2F*> ((static_cast<TKey*>(totHistDir->GetListOfKeys()->First()))->ReadObj()));
706 if (h2dTot) h2dTot->SetDirectory(0);
707 }
708 if (not h2dTot){
709 std::cout<<" Unrecoverable error in : " <<__LINE__ <<" of CalibrateIBL.cxx\n";
710 abort();
711 }
712
713 unique_ptr<TH2F> h2dTotAux;
714 TDirectoryFile* totHistDirAux(static_cast<TDirectoryFile*>(rodDirAux->Get(totHistDirPath.c_str())));
715 if (!totHistDirAux)
716 {
717 std::cout<<" Missing totHistDir in : " << totHistDirPath << endl;
718 logout << " Missing totHistDir in : " << totHistDirPath << endl;
719 abort();
720 }
721 else {
722 h2dTotAux.reset(static_cast<TH2F*> ((static_cast<TKey*>(totHistDirAux->GetListOfKeys()->First())->ReadObj())));
723
724 if(h2dTotAux) h2dTotAux->SetDirectory(0);
725 }
726 if (not h2dTotAux){
727 std::cout<<" Unrecoverable error in : " <<__LINE__ <<" of CalibrateIBL.cxx\n";
728 abort();
729 }
730 std::string totSigHistDirPath = modName + "/" + totSigHistName + "/A0/B0/C";
731 totSigHistDirPath += std::to_string(c);
732 TDirectoryFile* totSigHistDir (static_cast<TDirectoryFile*>(rodDir->Get(totSigHistDirPath.c_str())));
733 unique_ptr<TH2F> h2dTotSig (static_cast<TH2F*> ((static_cast<TKey*>(totSigHistDir->GetListOfKeys()->First()))->ReadObj()));
734 h2dTotSig->SetDirectory(0);
735
736 // Now one fill the calib-data from all pixles as basics.
737 for (int ieta = 0; ieta < ncol; ieta++)
738 {
739 for (int iphi = 0; iphi < nrow; iphi++)
740 {
741 int circ = -1;
742 if (moreFE)
743 {
744 // virtual FE with 80*84 pixels, so that sFE0, sFE1, sFE7, sFE6 will make up the physical 1'st FE
745 // while sFE2, sFE3, sFE5, sFE4 make up the physical 2'nd.
746 circ = (int)(iphi / 84.); // sFE0, sFE1, sFE2, sFE3
747 if (ieta > 80)
748 {
749 circ = 7 - circ;
750 } // sFE7, sFE6, sFE5, sFE4
751 }
752 else
753 {
754 if (ieta < ncol / 2)
755 circ = 0;
756 else
757 circ = 1;
758 }
759
760 float tot = h2dTot->GetBinContent(ieta + 1, iphi + 1) + HDCshift;
761 float totAux = h2dTotAux->GetBinContent(ieta + 1, iphi + 1);
762
763 // use ToT from higher HDC as default
764 float filltot = tot;
765
766 // use ToT from lower HDC for SMALL charges, dirty in studying
767 // May 4, 2022, please keep below strategy, it give reasonable merging till now
768
769 if (tot <= 5 || c <= 5)
770 {
771 if ((c == 0 && tot - totAux > 0.8) || (c >= 1 && c <= 3 && tot - totAux > 0.5) || (c > 4 && tot - totAux > 0.3))
772 filltot = totAux;
773 else
774 filltot = 0.5 * (totAux + tot);
775 if (c == 0 && totAux > 2.)
776 filltot = 1.;
777 }
778
779 if (c > 19)
780 filltot = 0.5 * (totAux + tot + 2.);
781
782 if (filltot > 0.){
783 h2d_XchrgYtot[circ]->Fill(chrgAbaciIBL[c], filltot, 1.);
784 }
785
786
787 float err = h2dTotSig->GetBinContent(ieta + 1, iphi + 1);
788 h2d_XchrgYToTSig[circ]->Fill(chrgAbaciIBL[c], err, 1.);
789 if (err == 0.)
790 err = 0.05;
791
792#if defined(DEMOXCHECK)
793 if (modName == XcheckModule && (iphi >= XcheckPhi[0] && iphi <= XcheckPhi[1]) && (ieta >= XcheckEta[0] && ieta <= XcheckEta[1]) && (c == XcheckCharge[0] || c == XcheckCharge[1] || (int)(tot + 0.5) == XcheckToT || (int)(totAux + 0.5) == XcheckToT))
794 logout << " tot = " << tot << ", totAux = " << totAux << " @injection " << chrgAbaciIBL[c] << " " << err << endl;
795#endif
796 }
797 } // all pixels done over a module
798
799 h2dTotSig.reset();
800 h2dTotAux.reset();
801 h2dTot.reset();
802
803 } // end looping over charges
804 float modHash = -1;
805 if (!(ibl3Dfe0 || ibl3Dfe1))
806 {
807 modHash = hashID;
808 }
809
810 array<TDirectory *, npsFEs> dirFE;
811
812 Int_t idxMod = cntRod * 16 * npsFEs + cntMod * npsFEs;
813
814 for (int sfe = 0; sfe < npsFEs; sfe++)
815 {
816 stringstream ss;
817 ss << sfe;
818
819 if (ibl3Dfe0 && sfe < npsFEs / 2)
820 {
821 modHash = hashIDL;
822 }
823 else if (ibl3Dfe1 && sfe >= npsFEs / 2)
824 {
825 modHash = hashIDR;
826 }
827 else
828 modHash += (sfe < npsFEs / 2 ? 0 : 1) * 0.8;
829
830 string prfmodname = "DirFE_" + ss.str();
831 dirFE[sfe] = dirMod->mkdir(prfmodname.c_str());
832
833 // standard deviation will be adopted for ERROR
834 prfmodname = modStr + "_profile_Tot_FE" + ss.str();
835 std::unique_ptr<TProfile> prfl_TotsFE(h2d_XchrgYtot[sfe]->ProfileX(prfmodname.c_str(), 1, -1, "s"));
836 prfl_TotsFE->SetTitle(prfmodname.c_str());
837 prfmodname = modStr + "_profile_Chrg_FE" + ss.str();
838 // shanly : seeming free from ROOT-7770 issue. 2021 Sep 17
839 std::unique_ptr<TProfile> prfl_ChrgsFE(h2d_XchrgYtot[sfe]->ProfileY(prfmodname.c_str(), 1, -1, "s"));
840 prfl_ChrgsFE->SetTitle(prfmodname.c_str());
841
842 Double_t TotArr[nchargeIBL], TotErrArr[nchargeIBL];
843 Double_t ChrgArr[nToTibl], ChrgErrArr[nToTibl];
844
845 Int_t idxFE = idxMod + sfe;
846 // ToT distribution along nchargeIBL bins charges
847 int dumbFE = 0;
848 for (int c = 0; c < nchargeIBL; c++)
849 {
850 TotArr[c] = TotArray[idxFE][c] = prfl_TotsFE->GetBinContent(c + 1);
851 TotErrArr[c] = TotErrArray[idxFE][c] = prfl_TotsFE->GetBinError(c + 1);
852
853 occuPhiEta[c][idxFE] = prfl_TotsFE->GetBinEntries(c + 1) * npsFEs / (1. * ncol * nrow);
854
855 if ((1. - occuPhiEta[c][idxFE]) > badRDfracCut)
856 dumbFE++;
857
858 stringstream sc;
859 sc << c;
860 prfmodname = modStr + "_ToTsig_FE" + ss.str() + "Chrg" + sc.str();
861
862 std::unique_ptr<TH1D> h_ToTsig(h2d_XchrgYToTSig[sfe]->ProjectionY(prfmodname.c_str(), c + 1, c + 1));
863 TotSigArray[idxFE][c] = TMath::Sqrt(h_ToTsig->GetMean() * h_ToTsig->GetMean() + h_ToTsig->GetRMS() * h_ToTsig->GetRMS());
864 TotSigErrArray[idxFE][c] = TMath::Sqrt(h_ToTsig->GetMeanError() * h_ToTsig->GetMeanError() + h_ToTsig->GetRMSError() * h_ToTsig->GetRMSError());
865
866 h_ToTsig.reset();
867 }
868 prfl_TotsFE.reset();
869
870 if (dumbFE > 5)
871 {
872 for (int c = 0; c < nchargeIBL; c++)
873 {
874 TotArr[c] = TotArray[idxFE][c] = 0.5 * (TotArray[idxFE - 1][c] + TotArray[idxFE - 2][c]);
875 TotErrArr[c] = TotErrArray[idxFE][c] = 0.5 * (TotErrArray[idxFE - 1][c] + TotErrArray[idxFE - 2][c]);
876
877 TotSigArray[idxFE][c] = 0.5 * (TotSigArray[idxFE - 1][c] + TotSigArray[idxFE - 2][c]);
878 TotSigErrArray[idxFE][c] = 0.5 * (TotSigErrArray[idxFE - 1][c] + TotSigErrArray[idxFE - 2][c]);
879
880 }
881 }
882 // Charges distribution along 16 bins ToT
883 Double_t errToT_overChrg[nToTibl];
884 for (int t = 0; t < nToTibl; t++)
885 {
886 ChrgArr[t] = ChrgArray[idxFE][t] = prfl_ChrgsFE->GetBinContent(t + 1);
887 ChrgErrArr[t] = ChrgErrArray[idxFE][t] = prfl_ChrgsFE->GetBinError(t + 1);
888
889 occuChrgs[t][cntMod * npsFEs + sfe] = prfl_ChrgsFE->GetBinEntries(t + 1) / RDentries[t];
890
891 stringstream st;
892 st << t;
893 prfmodname = modStr + "_Chrg_FE" + ss.str() + "ToT" + st.str();
894 std::unique_ptr<TH1D> h_chrg(h2d_XchrgYtot[sfe]->ProjectionX(prfmodname.c_str(), t + 1, t + 1));
895 st.str("");
896 st << totAbaci[t];
897 prfmodname = modStr + " Charge_FE" + ss.str() + " @ ToT : " + st.str();
898 h_chrg->SetTitle(prfmodname.c_str());
899
900 Double_t scl = h_chrg->Integral();
901 errToT_overChrg[t] = 0.;
902 if (scl > 0.)
903 for (int b = 0; b < nchargeIBL; b++)
904 {
905 errToT_overChrg[t] += TotErrArr[b] * (h_chrg->GetBinContent(b + 1)) / scl;
906 }
907 if (errToT_overChrg[t] == 0.)
908 errToT_overChrg[t] = 1.1;
909
910 h_chrg.reset();
911 }
912 prfl_ChrgsFE.reset();
913#if defined(DEMOXCHECK)
914 for (int b = 1; b < nToTibl; b++)
915 {
916 double err = errToT_overChrg[b];
917 h1d_totSprd[b - 1]->Fill(err);
918 h1d_totSprdAll[b - 1]->Fill(err);
919 }
920#endif
921
922 // try a correction upon dumb FE
923 if (dumbFE > 5)
924 {
925 for (int t = 0; t < nToTibl; t++)
926 {
927 ChrgArr[t] = ChrgArray[idxFE][t] = 0.5 * (ChrgArray[idxFE - 1][t] + ChrgArray[idxFE - 2][t]);
928 ChrgErrArr[t] = ChrgErrArray[idxFE][t] = 0.5 * (ChrgErrArray[idxFE - 1][t] + ChrgErrArray[idxFE - 2][t]);
929 }
930 }
931 bool xoticMod = std::find(reversedModules.begin(), reversedModules.end(), int(modHash)) != reversedModules.end();
932
933 // when charges' order are distorted along ToTs
934 bool reverseH = false, reverseT = false;
935 for (int t = 4; t >= 1; t--)
936 if (ChrgArr[t] > ChrgArr[t + 1])
937 {
938 reverseH = true;
939 break;
940 }
941 for (int t = 4; t >= 1; t--)
942 {
943 int tl = nToTibl - t;
944 if (ChrgArr[tl] < ChrgArr[tl - 1])
945 {
946 reverseT = true;
947 break;
948 }
949 }
950
951 if (reverseH)
952 {
953#if defined(DEMOXCHECK)
954 logout << " Head reversed : " << int(modHash) << " " << modName << " " << sfe << " " << ChrgArr[1] << " "
955 << ChrgArr[2] << " " << ChrgArr[3] << " " << ChrgArr[4] << " " << ChrgArr[5] << std::endl;
956#endif
957 // try a correction with some averaged factors
958 for (int t = 4; t >= 1; t--)
959 if (ChrgArr[t] > ChrgArr[t + 1] && (occuChrgs[t][cntMod * npsFEs + sfe] < 0.5 || xoticMod))
960 ChrgArr[t] = ChrgArray[idxFE][t] = reverseCF_H[t - 1] * ChrgArr[t + 1];
961 }
962 else // prepare the average factors along ToTs
963 {
964 for (int t = 4; t >= 1; t--)
965 if (ChrgArr[t + 1] != 0.)
966 reverseCF_H[t - 1] = (reverseCF_H[t - 1] == tmpRCf_H[t - 1] ? ChrgArr[t] / ChrgArr[t + 1] : 0.5 * (reverseCF_H[t - 1] + ChrgArr[t] / ChrgArr[t + 1]));
967 }
968
969 if (reverseT)
970 {
971#if defined(DEMOXCHECK)
972 logout << " Tail reversed : " << int(modHash) << " " << modName << " " << sfe << " " << ChrgArr[nToTibl - 5] << " "
973 << ChrgArr[nToTibl - 4] << " " << ChrgArr[nToTibl - 3] << " " << ChrgArr[nToTibl - 2] << " " << ChrgArr[nToTibl - 1] << std::endl;
974#endif
975 // try a correction with some averaged factors
976 for (int t = 4; t >= 1; t--)
977 {
978 int tl = nToTibl - t;
979 if (ChrgArr[tl] < ChrgArr[tl - 1] && (occuChrgs[tl][cntMod * npsFEs + sfe] < 0.5 || xoticMod))
980 ChrgArr[tl] = ChrgArray[idxFE][tl] = reverseCF_T[t - 1] * ChrgArr[tl - 1];
981 }
982 }
983 else // prepare the average factors along ToTs
984 {
985 for (int t = 4; t >= 1; t--)
986 {
987 int tl = nToTibl - t;
988 if (ChrgArr[tl - 1] != 0.)
989 reverseCF_T[t - 1] = (reverseCF_T[t - 1] == tmpRCf_T[t - 1] ? ChrgArr[tl] / ChrgArr[tl - 1] : 0.5 * (reverseCF_T[t - 1] + ChrgArr[tl] / ChrgArr[tl - 1]));
990 }
991 }
992
993 string gername = modStr + "_grToTsig_FE_" + ss.str();
994 std::unique_ptr<TGraphErrors> grTotSig = std::make_unique<TGraphErrors>(nchargeIBL, chrgAbaciIBL, TotSigArray[sfe], chargeErrArrIBL, TotSigErrArray[sfe]);
995 grTotSig->SetTitle(gername.c_str());
996 grTotSig->SetName(gername.c_str());
997
998 std::unique_ptr<TGraph> grTotSprd = std::make_unique<TGraph>(nToTibl, totAbaci, errToT_overChrg);
999 gername = "ToTsprd";
1000 grTotSprd->SetName(gername.c_str());
1001 gername = modStr + " ToT spread over FrontEnd @ FE " + ss.str();
1002 grTotSprd->SetTitle(gername.c_str());
1003 grTotSprd->SetLineColor(4);
1004 grTotSprd->SetMarkerStyle(4);
1005 grTotSprd->SetMarkerSize(1.);
1006
1007 TF1 ToTres("ToTres", funcDispConvoluted, 0.8, 16.5, 3);
1008 ToTres.SetParameter(1,13.);
1009 ToTres.SetParameter(2,0.65);
1010 ToTres.SetParLimits(1, 10., 16.);
1011 ToTres.SetParLimits(2, 0.3, 1.0);
1012
1013 grTotSprd->Fit(&ToTres, "MRQ");
1014 Double_t parP0 = ToTres.GetParameter(0);
1015 Double_t parP1 = ToTres.GetParameter(1);
1016
1017 dirFE[sfe]->WriteTObject(grTotSprd.get());
1018 grTotSprd.reset();
1019 // Lets' keep the bad 2-para fit untill the dataBase will be updated someday
1041
1042 gername = modStr + "_grToT_FE_" + ss.str();
1043 std::unique_ptr<TGraphErrors> grTot = std::make_unique<TGraphErrors>(nchargeIBL, chrgAbaciIBL, TotArr, chargeErrArrIBL, TotErrArr);
1044 grTot->SetTitle(gername.c_str());
1045 grTot->SetName(gername.c_str());
1046 TF1 *f1ToTfromCharge = new TF1("ToTfromCharge", funcRation5, chrgAbaciIBL[0] - 100., chrgAbaciIBL[nchargeIBL - 1] + 300., 5);
1047
1048 grTot->Fit(f1ToTfromCharge, "MRQ");
1049
1050 std::unique_ptr<TGraphErrors> grChrg = std::make_unique<TGraphErrors>(nToTibl, totAbaci, ChrgArr, totErrArr, ChrgErrArr);
1051
1052 gername = modStr + "_grChrg_FE_" + ss.str();
1053 grChrg->SetName(gername.c_str());
1054 grChrg->SetTitle(gername.c_str());
1055 grChrg->SetLineColor(3);
1056 grChrg->SetLineWidth(3);
1057
1058 TF1 *f1ChargefromToTLeft = new TF1("ChargefromToTL", funcRation5, 1., 7.5, 5);
1059 f1ChargefromToTLeft->SetLineColor(2);
1060 grChrg->Fit(f1ChargefromToTLeft, "MRQ");
1061 if (f1ChargefromToTLeft->GetChisquare() / 4. > 2.){
1062 logout << "bad fit Left... " << std::endl;
1063 }
1064
1065 TF1 *f1ChargefromToTRight = new TF1("ChargefromToTR", funcRation5, 5.5, 16., 5);
1066 f1ChargefromToTRight->SetLineColor(4);
1067 grChrg->Fit(f1ChargefromToTRight, "MRQ+");
1068
1069 dirFE[sfe]->WriteTObject(grTot.get());
1070 dirFE[sfe]->WriteTObject(grTotSig.get());
1071 dirFE[sfe]->WriteTObject(h2d_XchrgYtot[sfe].get());
1072
1073 h2d_XchrgYtot[sfe].reset();
1074 h2d_XchrgYToTSig[sfe].reset();
1075 grTot.reset();
1076 grTotSig.reset();
1077 grChrg.reset();
1078
1079 // prepare a map sor sorting
1080 string Idx = "I" + ss.str();
1081 vector<std::string> modName2prt;
1082 modName2prt.push_back(modName);
1083 modName2prt.push_back((std::string)(Idx));
1084 // unfortunately some tuning is necessary to satisfy dataBase reauirements
1085 // 5 for module ID, 4 for threshold, nToTibl - 1 for charges at ToT, 3 for ResToT
1086 std::vector<Double_t> prtAux;
1087 prtAux.reserve(5 + 4 + nToTibl - 1 + 3); // skip ToT = 0
1088 // @ 0 ... 4
1089 prtAux.push_back(modHash);
1090 prtAux.push_back(bec);
1091 prtAux.push_back(layer);
1092 prtAux.push_back(phi_module);
1093 prtAux.push_back(eta_module);
1094
1095 modName2prt[0] = (modStr + "_" + ss.str()).c_str();
1096 if (ibl3Dfe0 && sfe < npsFEs / 2)
1097 {
1098 if (!moreFE)
1099 modName2prt[0] = modName + "_0";
1100 prtAux[0] = hashIDL;
1101 prtAux[4] = eta_moduleL;
1102 }
1103 if (ibl3Dfe1 && sfe >= npsFEs / 2)
1104 {
1105 if (!moreFE)
1106 modName2prt[0] = modName + "_1";
1107 prtAux[0] = hashIDR;
1108 prtAux[4] = eta_moduleR;
1109 }
1110
1111 // @ 5 ... 8
1112 float ThrNorm = pcdMap[modStr][Idx]["ThrNorm"];
1113 if (ThrNorm == 0. || ThrNorm == -42.)
1114 {
1115 // retrieve something from average when it's missing in calibration
1116 prtAux.push_back(THR_avg[0][sfe]);
1117 prtAux.push_back(ThrSig_avg[0][sfe]);
1118 prtAux.push_back(THR_avg[1][sfe]);
1119 prtAux.push_back(ThrSig_avg[1][sfe]);
1120 }
1121 else
1122 {
1123 prtAux.push_back(pcdMap[modStr][Idx]["ThrNorm"]);
1124 prtAux.push_back(pcdMap[modStr][Idx]["ThrSigNorm"]);
1125 prtAux.push_back(pcdMap[modStr][Idx]["ThrLong"]);
1126 prtAux.push_back(pcdMap[modStr][Idx]["ThrSigLong"]);
1127 }
1128 // since this version, the goodness will be skipped
1129 // @ 10 ... nToTibl + 8, with ToT = 0 skipped
1130 for (int t = 1; t < nToTibl; t++)
1131 prtAux.push_back(ChrgArr[t]);
1132
1133 // @ nToTibl + 9, nToTibl + 10
1134 prtAux.push_back(parP0);
1135 prtAux.push_back(parP1);
1136
1137 std::pair<vector<std::string>, vector<Double_t>> payloadDB = std::pair<vector<std::string>, vector<Double_t>>(modName2prt, prtAux);
1138 ModuDataToPrint.insert(std::pair<float, std::pair<vector<std::string>, vector<Double_t>>>(modHash, payloadDB));
1139 } // end 1'st loop over FrontEnds
1140
1141 cntMod++;
1142 } // end of loop over Mods
1143
1144#if defined(DEMOXCHECK)
1145
1146 TDirectory *dirToTSprd = dirRod->mkdir("ToT_Spreads");
1147 for (int b = 0; b < nToTibl - 1; b++)
1148 {
1149 Int_t upper = h1d_totSprd[b]->FindLastBinAbove(1, 1);
1150 h1d_totSprd[b]->SetAxisRange(h1d_totSprd[b]->GetBinLowEdge(
1151 h1d_totSprd[b]->FindFirstBinAbove(1, 1)),
1152 h1d_totSprd[b]->GetBinLowEdge(upper) + h1d_totSprd[b]->GetBinWidth(upper), "X");
1153
1154 dirToTSprd->WriteTObject(h1d_totSprd[b]);
1155 delete h1d_totSprd[b];
1156 }
1157 h1d_totSprd.clear();
1158
1159 std::string hn = "";
1160
1161 TH2F *h2_badChrgs = new TH2F(rodName + "ChargeOccupancy", "Occupancy ", nchargeIBL, 0, nchargeIBL, 16 * npsFEs, 0, 16 * npsFEs);
1162 hn = "Occupancy along charges @ " + rodName;
1163 h2_badChrgs->SetTitle(hn);
1164
1165 TH2F *h2_badToT = new TH2F(rodName + "ToTOccupancy", "Occupancy ", nToTibl - 1, 0.5, nToTibl - 0.5, 16 * npsFEs, 0, 16 * npsFEs);
1166 hn = "Occupancy along ToT @ " + rodName;
1167 h2_badToT->SetTitle(hn);
1168
1169 bool fillChrg = false, fillToT = false;
1170
1171 for (int sfe = 0; sfe < 16 * npsFEs; sfe++)
1172 {
1173 float PhiEta = 0.;
1174 for (int c = 0; c < nchargeIBL; c++)
1175 {
1176 float occu = occuPhiEta[c][cntRod * 16 * npsFEs + sfe];
1177 float nt = 10. - floor(10. * occu + 0.5);
1178 PhiEta += 0.1 * nt;
1179 std::string modchrg = (std::string)(modNames[sfe] + "_Chrg_" + c);
1180 if (nt > 0.01)
1181 badModules_Order_detailed.insert(std::pair<float, std::string>(1. - occu, modchrg));
1182 if (nt > 5.)
1183 continue;
1184 if (nt < 0.1)
1185 nt = 0.01;
1186
1187 h2_badChrgs->SetBinContent(c + 1, sfe + 1, nt);
1188 fillChrg = true;
1189 }
1190 PhiEta /= nchargeIBL;
1191 if (PhiEta > 0.01)
1192 badModules_Order.insert(std::pair<float, std::string>(PhiEta,
1193 (std::string)(modNames[sfe])));
1194
1195 for (int t = 1; t < nToTibl; t++)
1196 {
1197 float nt = floor(10. * occuChrgs[t][sfe] + 0.5);
1198 if (nt > 5.)
1199 continue;
1200 h2_badToT->SetBinContent(t + 1, sfe + 1, nt);
1201 fillToT = true;
1202 }
1203
1204 std::string shortName = modNames[sfe].substr(9, 13).c_str();
1205 h2_badChrgs->GetYaxis()->SetBinLabel(sfe + 1, shortName);
1206 h2_badToT->GetYaxis()->SetBinLabel(sfe + 1, shortName);
1207 h2_badChrgs->SetStats(kFALSE);
1208 h2_badToT->SetStats(kFALSE);
1209 h2_badChrgs->SetTickLength(0.01, "Y");
1210 h2_badToT->SetTickLength(0.01, "Y");
1211 }
1212
1213 h2_badChrgs->SetOption("TEXT");
1214 h2_badToT->SetOption("TEXT");
1215
1216 if (fillChrg)
1217 {
1218 dirRod->WriteTObject(h2_badChrgs, "", "Overwrite");
1219 }
1220
1221 if (fillToT)
1222 {
1223 dirRod->WriteTObject(h2_badToT, "", "WriteDelete");
1224 }
1225
1226 delete h2_badChrgs;
1227 delete h2_badToT;
1228
1229#endif
1230
1231 cntRod++;
1232 } // end of loop over Rods
1233
1234#if defined(DEMOXCHECK)
1235
1236 for (int t = 0; t < 4; t++)
1237 {
1238 if (t == 0)
1239 logout << " correction factors for reversed charges : ";
1240 logout << reverseCF_H[t] << ", ";
1241 if (t == 3)
1242 logout << endl;
1243 }
1244 for (int t = 0; t < 4; t++)
1245 {
1246 if (t == 0)
1247 logout << " correction factors for reversed charges : ";
1248 logout << reverseCF_T[t] << ", ";
1249 if (t == 3)
1250 logout << endl;
1251 }
1252
1253 for (int t = 0; t < nToTibl; t++)
1254 {
1255 if (t == 0)
1256 logout << "RD Entries per ToT bin ";
1257 logout << h1_ChrgEntry[t]->GetMean() << ", ";
1258 if (t == nToTibl - 1)
1259 logout << std::endl;
1260 roTotDir->WriteTObject(h1_ChrgEntry[t]);
1261 delete h1_ChrgEntry[t];
1262 }
1263 h1_ChrgEntry.clear();
1264
1265 for (std::map<float, std::string>::const_iterator itr = devChrg_Order.begin();
1266 itr != devChrg_Order.end(); ++itr)
1267 logout << " Charge dev order : " << itr->second << " : " << itr->first << endl;
1268
1269 for (std::map<float, std::string>::const_iterator itr = devToT_Order.begin();
1270 itr != devToT_Order.end(); ++itr)
1271 logout << " ToT dev order : " << itr->second << " : " << itr->first << endl;
1272
1273 logout << "ToTRes (spread) & its RMS @ ToT from 1 to 16 : " << endl;
1274 for (int t = 0; t < nToTibl - 1; t++)
1275 {
1276 if (t == 0)
1277 logout << "[ ";
1278 logout << "[ " << h1d_totSprdAll[t]->GetMean() << " , " << h1d_totSprdAll[t]->GetStdDev() << " ] ";
1279 if (t < nToTibl - 2)
1280 logout << " , " << endl;
1281 if (t == nToTibl - 2)
1282 logout << " ] " << endl;
1283 delete h1d_totSprdAll[t];
1284 }
1285 h1d_totSprdAll.clear();
1286#endif
1287
1288 logout << " goto print " << std::endl;
1289
1290 // Comment out for information only
1291 // txtDB << "####################################################################################################" << endl;
1292 // txtDB << "=========== physics meaning of each column ====================" << endl;
1293 // txtDB << "hashId " << "Threshold " << "ThresholdRMS " << "ThresholdLong " << "ThresholdRMSLong " << "Chrg_ToT1 "
1294 // << "Chrg_ToT2 " << "Chrg_ToT3 " << "Chrg_ToT4 " << "Chrg_ToT5 " << "Chrg_ToT6 " << "Chrg_ToT7 " << "Chrg_ToT8 " << "Chrg_ToT9 "
1295 // << "Chrg_ToT10 " << "Chrg_ToT11 " << "Chrg_ToT12 " << "Chrg_ToT13 " << "Chrg_ToT14 " << "Chrg_ToT15 " << "Chrg_ToT16 "
1296 // << "ToT_sprd " << " ToT_sprd_var" << endl;
1297
1298 // txtDB << "######################## TXT for IBL calibration ____Start : _> ###########################" << endl;
1299
1300 // now print the dataBase payloads in a TXT format
1301 for (std::map<float, std::pair<vector<std::string>, vector<Double_t>>>::const_iterator itr = ModuDataToPrint.begin();
1302 itr != ModuDataToPrint.end(); ++itr)
1303 {
1304 int hash = floor(itr->first);
1305 std::pair<vector<std::string>, vector<Double_t>> payload = itr->second;
1306 vector<Double_t> fe = payload.second;
1307
1308 txtDB << hash << " ";
1309
1310 for (unsigned int t = 5; t < fe.size(); t++)
1311 {
1312 if (t < fe.size() - 2)
1313 {
1314 txtDB << (int)(fe[t]) << " ";
1315#if defined(DEMOXCHECK)
1316 if (fe[t] > fe[t + 1] && t >= 9 && t <= 23)
1317 logout << " Reversed charge ! " << t - 5 << " " << fe[t] << " " << fe[t + 1] << std::endl;
1318#endif
1319 }
1320 else
1321 {
1322 float p01 = abs(fe[t]);
1323 if (t < fe.size() - 1)
1324 txtDB << 0.001 * (int)(p01 * 1000.) << " ";
1325 else
1326 txtDB << 0.000000001 * (int)(p01 * 1000000000.) << " ";
1327 }
1328 }
1329 txtDB << endl;
1330 }
1331 // txtDB << "####################################################################################################" << endl;
1332 // txtDB << "######################## TXT for IBL calibration <___ : End _____ ###########################" << endl;
1333 // txtDB << "####################################################################################################" << endl;
1334
1335 std::cout << " Please find the file : " << dbFileName << " for dataBase payload " << endl;
1336
1337 // statistics for BAD frontends, only for ELOG reports
1338 logout << " modules lacking in RD during Threshold scan : " << std::endl;
1339 for (std::multimap<float, std::string, std::greater<float>>::const_iterator itr = badThr_Order.begin(); itr != badThr_Order.end(); ++itr)
1340 logout << " " << itr->second << " " << (itr->first) * 100. << "%" << endl;
1341
1342 logout << " modules lacking in RD during ToT scan : " << std::endl;
1343 for (std::multimap<float, std::string, std::greater<float>>::const_iterator itr = badModules_Order.begin(); itr != badModules_Order.end(); ++itr)
1344 logout << " " << itr->second << " " << (itr->first) * 100. << "%" << endl;
1345
1346 logout << " modules lacking in RD at certain charges during ToT scan : " << std::endl;
1347 for (std::multimap<float, std::string, std::greater<float>>::const_iterator itr = badModules_Order_detailed.begin();
1348 itr != badModules_Order_detailed.end(); ++itr)
1349 logout << " " << itr->second << " " << (itr->first) * 100. << "%" << endl;
1350
1351 roFile.Close();
1352 txtDB.close();
1353 logout.close();
1354
1355 return 0;
1356}
1357
1358void printError(){
1359 printf("ERROR - Argument not expected or wrongly set:\n\n");
1360 printf("Valid format is: ./IBLCalibration THR=SCAN_Sxxxxxxxxx TOT_HISDIS=SCAN_Sxxxxxxxxx TOT_LOWQ=SCAN_Sxxxxxxxxx directory_path=path/to/file/\n");
1361 printf("\n\t i.e: ./IBLCalibration THR=SCAN_S000087719 TOT_HISDIS=SCAN_S000087717 TOT_LOWQ=SCAN_S000087710 directory_path=/eos/user/x/xxxx/\n");
1362}
1363
1364int main(int argc, char *argv[]) {
1365
1366 bool saveInfo = false;
1367 std::string THR = "THR";
1368 std::string TOT_HISDIS = "TOT_HISDIS";
1369 std::string TOT_LOWQ = "TOT_LOWQ";
1370 std::string dpath = "directory_path";
1371
1372 for(int i=1; i<argc; i++){
1373 std::string aux(argv[i]);
1374 if (THR.compare(aux.substr(0,aux.find("="))) == 0) THR = aux.substr(aux.find("=")+1);
1375 else if(TOT_HISDIS.compare(aux.substr(0,aux.find("="))) == 0) TOT_HISDIS = aux.substr(aux.find("=")+1);
1376 else if(TOT_LOWQ.compare(aux.substr(0,aux.find("="))) == 0) TOT_LOWQ = aux.substr(aux.find("=")+1);
1377 else if(dpath.compare(aux.substr(0,aux.find("="))) == 0) dpath = aux.substr(aux.find("=")+1);
1378 else if(aux.compare("saveInfo") == 0) saveInfo = true;
1379 else{
1380 printError();
1381 return 1;
1382 }
1383 }
1384
1385 printf("%-14s = %s\n","Directory path",dpath.c_str());
1386 printf("%-14s = %s.root\n","THR",THR.c_str());
1387 printf("%-14s = %s.root\n","TOT_HISDIS",TOT_HISDIS.c_str());
1388 printf("%-14s = %s.root\n","TOT_LOWQ",TOT_LOWQ.c_str());
1389 printf("%-14s = %s\n\n\n" ,"Save root file",saveInfo ? "True" : "False" );
1390
1391 bool correctArgc = (THR.compare("THR") == 0) or (TOT_HISDIS.compare("TOT_HISDIS") == 0) or (TOT_LOWQ.compare("TOT_LOWQ") == 0) or (dpath.compare("directory_path") == 0);
1392
1393 if(correctArgc){
1394 printf("Cannot continue, one arguments is incorrect or not filled correctly...\n");
1395 printf("Helper below:\n**********************\n\n");
1396 printError();
1397 return 1;
1398 }
1399
1400 iblCalib(dpath, THR, TOT_HISDIS, TOT_LOWQ);
1401 return 0;
1402
1403}
1404
ToTfromChrg
double ToTfromChrg(float chrg, double par[5])
Definition CalibrateIBL.cxx:82
funcDispConvoluted
double funcDispConvoluted(double *x, double *par)
Definition CalibrateIBL.cxx:102
funcDisp
double funcDisp(double *x, double *par)
Definition CalibrateIBL.cxx:62
funcRation5
double funcRation5(double *x, double *par)
Definition CalibrateIBL.cxx:92
ChrgfromToT
double ChrgfromToT(float tot, double par[10])
Definition CalibrateIBL.cxx:69
printError
void printError()
Definition CalibrateIBL.cxx:1358
iblCalib
int iblCalib(const std::string &InDir, const std::string &THRscan, const std::string &scan, const std::string &scanLowCharge)
Definition CalibrateIBL.cxx:109
to_string
static std::string to_string(const std::vector< T > &v)
Definition FPGATrackSimGenScanTool.cxx:44
upper
int upper(int c)
Definition LArBadChannelParser.cxx:49
ss
static Double_t ss
Definition LArPhysWaveHECTool.cxx:37
sc
static Double_t sc
Definition LArPhysWaveHECTool.cxx:37
x
#define x
checker_macros.h
Define macros for attributes used to control the static checker.
ATLAS_NO_CHECK_FILE_THREAD_SAFETY
#define ATLAS_NO_CHECK_FILE_THREAD_SAFETY
Definition checker_macros.h:209
PathResolver::find_file
static std::string find_file(const std::string &logical_file_name, const std::string &search_path)
Definition PathResolver.cxx:220
pix::PixelMapping::mapping
void mapping(const std::string &geographicalID, int *hashID, int *bec, int *layer, int *phimod, int *etamod) const
Definition PixelMapping.cxx:68
chi2
double chi2(TH1 *h0, TH1 *h1)
Definition comparitor.cxx:525
get
T * get(TKey *tobj)
get a TObject* from a TKey* (why can't a TObject be a TKey?)
Definition hcg.cxx:132
main
int main()
Definition hello.cxx:18
scan
void scan(TDirectory *td=0, int depth=0)
Definition listroot.cxx:440
std
STL namespace.
Generated on for ATLAS Offline Software by doxygen 1.16.1