ATLAS Offline Software
Calib.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 //***********************************************************************************************
6 // PixelCalibration - Program to run the pixel charge calibration (not IBL)
7 // ---------------------------------------
8 // begin : 01 11 2023
9 // email : sergi.rodriguez@cern.ch
10 //***********************************************************************************************
11 
12 
13 #include "Calib.h"
14 
15 bool Calib::totFitting (const pix::PixelMapping &pm, const std::string &inTotFile, std::map<unsigned int , std::vector<std::unique_ptr<CalibFrontEndInfo>> > &map_info){
16 
17  if (inTotFile.empty()) return false;
18 
19  TFile totFile(inTotFile.c_str(),"READ");
20  if (not totFile.IsOpen()) {
21  printf("Error - File %s could not be opened.\n",inTotFile.c_str());
22  totFile.Close();
23  return false;
24  } else {
25  printf("File %s opened.\n",inTotFile.c_str());
26  printf("Running TOT calibration...\n");
27  }
28 
29  // Creating here the histograms with this scope ownership
30  std::vector< std::vector< std::unique_ptr<TH1F> > > histogramsTOT;
31  std::vector< std::vector< std::unique_ptr<TH1F> > > histogramsTOTsig;
32  for(unsigned int FE = 0; FE < m_nFE; FE++){
33 
34  histogramsTOT.push_back( std::vector< std::unique_ptr<TH1F> >() );
35  histogramsTOTsig.push_back( std::vector< std::unique_ptr<TH1F> >() );
36 
37  //Here we combine long and ganged pixels
38  for(unsigned int pixel=0; pixel<2; pixel++){
39  TString title = "FE"+std::to_string(FE)+"_pixType"+std::to_string(pixel);
40  std::unique_ptr<TH1F> htot = std::make_unique<TH1F>(title+"_tot", title+"_tot", m_totnbins, m_totLo, m_totHi);
41  htot->SetDirectory(0);
42  std::unique_ptr<TH1F> htotsig = std::make_unique<TH1F>(title+"_totsig", title+"_totsig", m_totsigNBins, m_totsigLo, m_totsigHi);
43  htotsig->SetDirectory(0);
44  //cppcheck-suppress containerOutOfBounds
45  histogramsTOT.at(FE).push_back(std::move(htot));
46  //cppcheck-suppress containerOutOfBounds
47  histogramsTOTsig.at(FE).push_back(std::move(htotsig));
48  }
49  }
50 
51  // Start looping over the ROD, MOD and charges
52  TIter rodItr=getRodIterator(totFile);
53  TKey* rodKey;
54 
55  while ((rodKey=static_cast<TKey*>(rodItr()))) {
56  TDirectoryFile* rodDir = static_cast<TDirectoryFile*>(rodKey->ReadObj());
57  TIter modItr = getModuleIterator(rodDir);
58  const TString rodName(rodKey->GetName());
59  printf("%s\n",rodName.Data());
60  TKey* modKey;
61  while ((modKey=static_cast<TKey*>(modItr()))) {
62  const TString modName(modKey->GetName());
63 
64  if (not moduleInPart(modName)) continue;
65  if (not pm.contains(std::string(modName))) continue;
66 
67  if( m_runOneMOD and strcmp(modName, m_testMOD) != 0){
68  continue;
69  }
70 
71  printf(" -> %s\n",modName.Data());
72 
73  //creates arrays for the Tgraph
74  std::array<std::array<float, m_ncharge>, m_nFE> totArrI{};
75  std::array<std::array<float, m_ncharge>, m_nFE> totErrArrI{};
76  std::array<std::array<float, m_ncharge>, m_nFE> totSigArrI{};
77  std::array<std::array<float, m_ncharge>, m_nFE> totSigErrArrI{};
78  std::array<std::array<float, m_ncharge>, m_nFE> totLongArrI{};
79  std::array<std::array<float, m_ncharge>, m_nFE> totErrLongArrI{};
80 
81  // loop over charges
82  for (int c=0; c<m_ncharge; ++c) {
83 
84  // Get TH2 for a given charge
85  std::unique_ptr<TH2F> h2dTOTmean(get2DHistogramFromPath(rodDir,modName, "TOT_MEAN", c));
86  std::unique_ptr<TH2F> h2dTOTsig(get2DHistogramFromPath(rodDir,modName, "TOT_SIGMA", c));
87  if(!h2dTOTmean or !h2dTOTsig) {
88  return false;
89  }
90  h2dTOTmean->SetDirectory(0);
91  h2dTOTsig->SetDirectory(0);
92  // loop over pixels
93  for (unsigned int ieta = 0; ieta < m_etaBins; ieta++) {
94  for (unsigned int iphi = 0; iphi < m_phiBins; iphi++) {
95  float totmean = h2dTOTmean->GetBinContent(ieta + 1, iphi + 1);
96  float totsig = h2dTOTsig ->GetBinContent(ieta + 1, iphi + 1);
97 
98  if (totmean<0.1) {
99  continue;
100  }
101 
102  int FE = chipId(iphi, ieta);
103  int pixel= pixelType(iphi, ieta, true);
104 
105  if(FE<0){
106  return false;
107  }
108 
109  histogramsTOT.at(FE).at(pixel)->Fill(totmean);
110  histogramsTOTsig.at(FE).at(pixel)->Fill(totsig);
111 
112  }
113  }
114 
115  // free memory
116  h2dTOTmean.reset();
117  h2dTOTsig.reset();
118 
119  //filling arrays
120  for(unsigned int FE = 0; FE < m_nFE; FE++){
121  for(unsigned int pixel = 0; pixel <2; pixel++){
122 
123  if(pixel == 0){
124  totArrI.at(FE).at(c) = histogramsTOT.at(FE).at(pixel)->GetMean();
125  totErrArrI.at(FE).at(c) = histogramsTOT.at(FE).at(pixel)->GetMeanError();
126  totSigArrI.at(FE).at(c) = std::sqrt(std::pow(histogramsTOTsig.at(FE).at(pixel)->GetMean() ,2)+std::pow(histogramsTOT.at(FE).at(pixel)->GetRMS() ,2));
127  totSigErrArrI.at(FE).at(c) = std::sqrt(std::pow(histogramsTOTsig.at(FE).at(pixel)->GetMeanError(),2)+std::pow(histogramsTOT.at(FE).at(pixel)->GetRMSError(),2));
128 
129  if(totSigErrArrI.at(FE).at(c) > 1.0){
130  totArrI.at(FE).at(c) = 0.0;
131  }
132 
133  }
134  else{
135  totLongArrI.at(FE).at(c) = histogramsTOT.at(FE).at(pixel)->GetMean();
136  totErrLongArrI.at(FE).at(c) = histogramsTOT.at(FE).at(pixel)->GetMeanError();
137 
138  if(totErrLongArrI.at(FE).at(c) > 1.0){
139  totLongArrI.at(FE).at(c) = 0.0;
140  }
141 
142  }
143 
144  //reset histogram for next iteration
145  histogramsTOT.at(FE).at(pixel)->Reset("ICESM");
146  histogramsTOTsig.at(FE).at(pixel)->Reset("ICESM");
147  }
148  }
149  } // End of charge.
150 
151  // loop over FE and create a graph for fitting
152  for(unsigned int FE = 0; FE < m_nFE; FE++) {
153 
154  TString subdir(((FE < 10) ? "FE0" : "FE") +std::to_string(FE));
155 
156  std::vector<float> v_Q;
157  std::vector<float> v_Qerr;
158  std::vector<float> v_TOT;
159  std::vector<float> v_TOTerr;
160  std::vector<float> v_TOTsig;
161  std::vector<float> v_TOTsigerr;
162  std::vector<float> v_TOTlong;
163  std::vector<float> v_TOTlongerr;
164 
165  std::copy(std::begin(m_chargeArr) , std::end(m_chargeArr) , std::back_inserter(v_Q) );
166  std::copy(std::begin(m_chargeErrArr) , std::end(m_chargeErrArr) , std::back_inserter(v_Qerr) );
167  std::copy(std::begin(totArrI[FE]) , std::end(totArrI[FE]) , std::back_inserter(v_TOT) );
168  std::copy(std::begin(totErrArrI[FE]) , std::end(totErrArrI[FE]) , std::back_inserter(v_TOTerr) );
169  std::copy(std::begin(totSigArrI[FE]) , std::end(totSigArrI[FE]) , std::back_inserter(v_TOTsig) );
170  std::copy(std::begin(totSigErrArrI[FE]) , std::end(totSigErrArrI[FE]) , std::back_inserter(v_TOTsigerr) );
171  std::copy(std::begin(totLongArrI[FE]) , std::end(totLongArrI[FE]) , std::back_inserter(v_TOTlong) );
172  std::copy(std::begin(totErrLongArrI[FE]) , std::end(totErrLongArrI[FE]) , std::back_inserter(v_TOTlongerr) );
173 
174  std::vector<float> pixNormalParams;
175  std::vector<float> pixNormalParamsQuality;
176  std::vector<float> pixSigParams;
177  std::vector<float> pixSigParamsQuality;
178  std::vector<float> pixLongParams;
179  std::vector<float> pixLongParamsQuality;
180 
181  //For normal pixels and sig
182  uint8_t n_fit = 0;
183  do{
184  int vecsize = v_Q.size();
185 
186  std::unique_ptr<TGraphErrors> graphnormal = std::make_unique<TGraphErrors>(vecsize, &v_Q.at(0), &v_TOT.at(0) , &v_Qerr.at(0), &v_TOTerr.at(0) );
187  std::unique_ptr<TGraphErrors> graphsig = std::make_unique<TGraphErrors>(vecsize, &v_Q.at(0), &v_TOTsig.at(0) , &v_Qerr.at(0), &v_TOTsigerr.at(0) );
188 
189  std::unique_ptr<TF1> functnormal = std::make_unique<TF1>("normal" ,new funcTot , m_chargeArr[m_qthresh]-100, m_chargeArr[m_ncharge-1]+100, 3);
190  std::unique_ptr<TF1> functnormalsig = std::make_unique<TF1>("normal_sig",new funcDisp, m_chargeArr[m_qthresh]-100, m_chargeArr[m_ncharge-1]+100, 2);
191 
192  graphnormal->Fit(functnormal.get() ,"MRQ");
193  graphsig ->Fit(functnormalsig.get(),"MRQ");
194 
195  pixNormalParams = getParams(functnormal.get() ,3 );
196  pixSigParams = getParams(functnormalsig.get(),2 );
197 
198  pixNormalParamsQuality = getParams_quality(functnormal.get() );
199  pixSigParamsQuality = getParams_quality(functnormalsig.get());
200 
201  if(m_savefile){
202 
203  m_wFile->cd();
204  if( !m_wFile->Get(rodName+"/"+modName+"/TOTfits/"+subdir) ){
205  m_wFile->mkdir(rodName+"/"+modName+"/TOTfits/"+subdir,rodName);
206  }
207 
208  m_wFile->cd(rodName+"/"+modName+"/TOTfits/"+subdir);
209 
210  graphTitles(graphnormal, TString(modName+" - "+subdir+" - normal pixels: Fit: "+std::to_string(n_fit)).Data(), "TOT");
211  graphTitles(graphsig , TString(modName+" - "+subdir+" - normal pixels: Fit: "+std::to_string(n_fit)).Data(), "Charge smearing");
212 
213  graphnormal->Write(TString("normal_fit_"+std::to_string(n_fit)), TObject::kWriteDelete);
214  graphsig->Write(TString("smearing_fit_"+std::to_string(n_fit)), TObject::kWriteDelete);
215  n_fit++;
216  }
217 
218  functnormal.reset();
219  functnormalsig.reset();
220 
221  graphnormal.reset();
222  graphsig.reset();
223 
224  }while(reFit_normalPix(pixNormalParams, v_Q, v_Qerr, v_TOT, v_TOTerr, v_TOTsig, v_TOTsigerr, FE ) );
225 
226 
227  // Since we have modified the vector size we need to clear it and refill it for the long and gange pixels
228  v_Q.clear();
229  v_Qerr.clear();
230 
231  std::copy(std::begin(m_chargeArr) , std::end(m_chargeArr) , std::back_inserter(v_Q) );
232  std::copy(std::begin(m_chargeErrArr) , std::end(m_chargeErrArr) , std::back_inserter(v_Qerr) );
233 
234 
235  //For long and ganged pixels (combined)
236  n_fit = 0;
237  do{
238 
239  //size will be modified in the condition..
240  int vecsize = v_Q.size();
241 
242  std::unique_ptr<TGraphErrors> graflong = std::make_unique<TGraphErrors>(vecsize, &v_Q.at(0), &v_TOTlong.at(0), &v_Qerr.at(0), &v_TOTlongerr.at(0) );
243 
244  std::unique_ptr<TF1> functlong = std::make_unique<TF1>("long" ,new funcTot , m_chargeArr[m_qthresh]-100, m_chargeArr[m_ncharge-1]+100, 3);
245 
246  graflong ->Fit(functlong.get() ,"MRQ");
247 
248  pixLongParams = getParams(functlong.get() ,3 );
249  pixLongParamsQuality = getParams_quality(functlong.get() );
250 
251  if(m_savefile){
252 
253  m_wFile->cd();
254  if( !m_wFile->Get(rodName+"/"+modName+"/TOTfits/"+subdir) ){
255  m_wFile->mkdir(rodName+"/"+modName+"/TOTfits/"+subdir,rodName);
256  }
257 
258  m_wFile->cd(rodName+"/"+modName+"/TOTfits/"+subdir);
259 
260  graphTitles(graflong, TString(modName+" - "+subdir+" - long+ganged pixels: Fit: "+std::to_string(n_fit)).Data(), "TOT");
261 
262  graflong->Write(TString("long_ganged_fit_"+std::to_string(n_fit)), TObject::kWriteDelete);
263  n_fit++;
264  }
265 
266  //delete the TF1
267  functlong.reset();
268 
269  //delete the graphs
270  graflong.reset();
271 
272  // no need to loop here.. leaving to improve it in the future - might need refitting
273  }while( false);
274 
275 
276  // Find the module it belong to
277  int modID = pm.getID(std::string(modName));
278  auto itr = map_info.find( modID );
279 
280  if (itr != map_info.end()) {
281  (itr->second).at(FE)->set_NormalParams( pixNormalParams);
282  (itr->second).at(FE)->set_LongParams ( pixLongParams );
283  (itr->second).at(FE)->set_SigParams ( pixSigParams );
284 
285  (itr->second).at(FE)->set_times_fitted ( n_fit );
286 
287  (itr->second).at(FE)->set_NormalParamsQuality( pixNormalParamsQuality);
288  (itr->second).at(FE)->set_LongParamsQuality ( pixLongParamsQuality );
289  (itr->second).at(FE)->set_SigParamsQuality ( pixSigParamsQuality );
290  }
291  else{
292  printf("Error - Module not found in fitting step... Skipping.\n");
293  return false;
294  }
295 
296  } // End of FE loop
297  } // End of MOD
298  } // End of ROD
299 
300  // remove from memory
301  for(unsigned int FE = 0; FE < m_nFE; FE++) {
302  for(unsigned int pixel=0; pixel<2; pixel++) {
303  histogramsTOT.at(FE).at(pixel).reset();
304  histogramsTOTsig.at(FE).at(pixel).reset();
305  }
306  }
307 
308  totFile.Close();
309  return true;
310 }
311 
312 
313 
314 bool Calib::reFit_normalPix(std::vector<float> &params, std::vector<float> &q, std::vector<float> &qerr, std::vector<float> &tot, std::vector<float> &toterr, std::vector<float> &sig, std::vector<float> &sigerr, const unsigned int fe){
315 
316  float vecFit_size = q.size()+1;
317  // float vecFit_size = q.size() - m_qthresh;
318  float stopFit = m_ncharge/2.0;
319  // float stopFit = (m_ncharge - m_qthresh)/2.0;
320  if(vecFit_size < stopFit) {
321 
322  // Default values for the fit
323  params.at(0) = 0;
324  params.at(1) = -28284.3;
325  params.at(2) = 0;
326 
327  printf("reFit_normalPix: Refitting skipped. Not enough points to fit.\n");
328 
329  return false;
330  }
331 
332  float parAI0 = params.at(0);
333  float parEI0 = params.at(1);
334  float parCI0 = params.at(2);
335 
336  std::vector<float> v_discrepancy;
337 
338  for(unsigned int i = 0; i < q.size(); i++){
339  float discrepancy = std::abs( 1 - ( (parAI0 * parEI0 - parCI0 * tot.at(i)) / (tot.at(i) - parAI0) ) / q.at(i) );
340 
341  if( i < m_qthresh ){
342  discrepancy = 0.0;
343  }
344  v_discrepancy.push_back(discrepancy);
345  }
346 
347  auto itr_max = std::max_element(v_discrepancy.begin(),v_discrepancy.end());
348 
349  if(*itr_max > m_chi_error){
350 
351  size_t n_max = std::distance(v_discrepancy.begin(), itr_max);
352  printf("FE %02u Refitted, removing charge %5.0f with chi_error %7.5f\n", fe ,q.at(n_max),*itr_max);
353  q.erase(q.begin()+n_max);
354  qerr.erase(qerr.begin()+n_max);
355  tot.erase(tot.begin()+n_max);
356  toterr.erase(toterr.begin()+n_max);
357  sig.erase(sig.begin()+n_max);
358  sigerr.erase(sigerr.begin()+n_max);
359 
360  return true;
361  }
362 
363  return false;
364 
365 }
366 
367 
368 bool Calib::fillTiming(const pix::PixelMapping &pm, const std::string &inTimFile, std::map<unsigned int, std::vector<std::unique_ptr<CalibFrontEndInfo>> > &map_info) {
369 
370  if (inTimFile.empty()) return false;
371 
372  TFile timFile(inTimFile.c_str(),"READ");
373  if (not timFile.IsOpen()) {
374  printf("Error - File %s could not be opened.\n",inTimFile.c_str());
375  return false;
376  } else {
377  printf("File %s opened.\n",inTimFile.c_str());
378  printf("Running timming calibration...\n");
379  }
380 
381  // Creating here the histograms with this scope owner
382  std::vector< std::vector< std::unique_ptr<TH1F> > > histogramsTIM;
383  for(unsigned int FE = 0; FE < m_nFE; FE++){
384 
385  histogramsTIM.push_back( std::vector< std::unique_ptr<TH1F> >() );
386 
387  for(unsigned int pixel=0; pixel<3; pixel++){
388  TString title = "FE"+std::to_string(FE)+"_pixType"+std::to_string(pixel);
389  std::unique_ptr<TH1F> h = std::make_unique<TH1F>(title+"_thr", title+"_thr", m_timnbins, m_timLo, m_timHi);
390  h->SetDirectory(0);
391  //cppcheck-suppress containerOutOfBounds
392  histogramsTIM.at(FE).push_back(std::move(h));
393  }
394  }
395 
396  //Will strat looping over the RODs
397  TIter rodItr = getRodIterator(timFile);
398  TKey* rodKey;
399  while ((rodKey=static_cast<TKey*>(rodItr()))) {
400  const TString rodName(rodKey->GetName());
401  TDirectoryFile* rodDir = static_cast<TDirectoryFile*>(rodKey->ReadObj());
402  TKey* modKey;
403  TIter modItr=getModuleIterator(rodDir);
404 
405  // Looping over the MODs of each ROD
406  while ((modKey=static_cast<TKey*>(modItr()))) {
407  TString modName(modKey->GetName());
408 
409  if ( not moduleInPart(modName)){
410  continue;
411  }
412  if ( not pm.contains(std::string(modName))){
413  printf("Error - Module %s not found in the PixelMapping tool\n",modName.Data());
414  continue;
415  }
416 
417  if( m_runOneMOD and strcmp(modName, m_testMOD) != 0){
418  continue;
419  }
420 
421  std::unique_ptr<TH2F> h2dTim(get2DHistogramFromPath(rodDir,modName, "SCURVE_MEAN"));
422  h2dTim->SetDirectory(0);
423 
424  for (unsigned int ieta = 0; ieta < m_etaBins; ieta++) {
425  for (unsigned int iphi = 0; iphi < m_phiBins; iphi++) {
426 
427  float tim = h2dTim->GetBinContent(ieta + 1, iphi + 1);
428 
429  if (tim<0.5) {
430  continue;
431  }
432 
433  int FE = chipId(iphi, ieta);
434  int pixel= pixelType(iphi, ieta);
435 
436  if(FE<0){
437  return false;
438  }
439 
440  histogramsTIM.at(FE).at(pixel)->Fill(tim);
441  }
442  }
443 
444  // Freeing memory of TH2F
445  h2dTim.reset();
446 
447  int modID = pm.getID(std::string(modName));
448  auto itr = map_info.find( modID );
449  if (itr == map_info.end()) {
450  printf("Calib::fillTiming: ERROR - Mod ID= %16s not found. Creating it -----> Inform Pixel Offline Software Experts... \n",std::string(modName).c_str());
451 
452  map_info[modID] = std::vector<std::unique_ptr<CalibFrontEndInfo>> ();
453 
454  for(unsigned int FE = 0; FE < m_nFE; FE++){
455 
456  map_info[modID].push_back( std::unique_ptr<CalibFrontEndInfo>() );
457  std::unique_ptr<CalibFrontEndInfo> p = std::make_unique<CalibFrontEndInfo>(modID,FE,std::string(modName),std::string(rodKey->GetName()));
458  map_info[modID].at(FE) = std::move(p);
459 
460  for(unsigned int pixel=0; pixel<3; pixel++){
461 
462  // Saving information for the calibration
463  int tim_mean = histogramsTIM.at(FE).at(pixel)->GetMean();
464 
465  // Reset histograms for next front end
466  histogramsTIM.at(FE).at(pixel)->Reset("ICESM");
467 
468  if(pixel == 0){ // normal
469  map_info[modID].at(FE)->set_NormalIntime(tim_mean);
470  }
471  else if(pixel == 1){ // long
472  map_info[modID].at(FE)->set_LongIntime(tim_mean);
473  }
474  else if(pixel == 2){ // ganged
475  map_info[modID].at(FE)->set_GangedIntime(tim_mean);
476  }
477  else{
478  printf("Error - Bad pixel in Calib::fillTiming\n");
479  return false;
480  }
481  }
482  }
483 
484  }
485  else{
486  for(unsigned int FE = 0; FE < m_nFE; FE++){
487  TString subdir(((FE < 10) ? "FE0" : "FE") +std::to_string(FE));
488  for(unsigned int pixel=0; pixel<3; pixel++){
489 
490  // Saving information for the calibration
491  int tim_mean = histogramsTIM.at(FE).at(pixel)->GetMean();
492 
493  if(m_savefile){
494 
495  m_wFile->cd();
496  if( !m_wFile->Get(rodName+"/"+modName+"/Thresholds/"+subdir) ){
497  m_wFile->mkdir(rodName+"/"+modName+"/Thresholds/"+subdir,rodName);
498  }
499 
500  m_wFile->cd(rodName+"/"+modName+"/Thresholds/"+subdir);
501 
502  histogramsTIM.at(FE).at(pixel)->SetTitle("Intime;Pixel intime;Counts");
503 
504  std::string type = "";
505  if(pixel == 0) type = "normal";
506  else if(pixel == 1) type = "long";
507  else if(pixel == 2) type = "ganged";
508 
509  histogramsTIM.at(FE).at(pixel)->Write(TString("intime_"+type).Data(), TObject::kWriteDelete);
510  }
511 
512  // Reset histograms for next front end
513  histogramsTIM.at(FE).at(pixel)->Reset("ICESM");
514 
515  if(pixel == 0){ // normal
516  (itr->second).at(FE)->set_NormalIntime(tim_mean);
517  }
518  else if(pixel == 1){ // long
519  (itr->second).at(FE)->set_LongIntime(tim_mean);
520  }
521  else if(pixel == 2){ // ganged
522  (itr->second).at(FE)->set_GangedIntime(tim_mean);
523  }
524  else{
525  printf("Error - Bad pixel in Calib::fillTiming\n");
526  return false;
527  }
528  } // End of pixel type loop
529  } // End of FE loop
530  }
531  } // End of MOD loop
532  }// End of ROD loop
533 
534  for(unsigned int FE = 0; FE < m_nFE; FE++) {
535  for(unsigned int pixel=0; pixel<3; pixel++) {
536  histogramsTIM.at(FE).at(pixel).reset();
537  }
538  }
539 
540  //Closing file - not needed anymore
541  timFile.Close();
542  printf("DONE with threshold calibration.\n");
543  return true;
544 }
545 
546 
547 
548 bool Calib::fillThresholds(const pix::PixelMapping &pm, const std::string &inThrFile, std::map<unsigned int, std::vector<std::unique_ptr<CalibFrontEndInfo>> > &map_info) {
549 
550  if (inThrFile.empty()) return false;
551 
552  TFile riThrFile(inThrFile.c_str(),"READ");
553  if (not riThrFile.IsOpen()) {
554  printf("Error - File %s could not be opened.\n",inThrFile.c_str());
555  return false;
556  } else {
557  printf("File %s opened.\n",inThrFile.c_str());
558  printf("Running threshold calibration...\n");
559  }
560 
561  // Creating here the histograms with this scope owner
562  std::vector< std::vector< std::unique_ptr<TH1F> > > histogramsTHR;
563  std::vector< std::vector< std::unique_ptr<TH1F> > > histogramsSIG;
564  for(unsigned int FE = 0; FE < m_nFE; FE++){
565 
566  histogramsTHR.push_back( std::vector< std::unique_ptr<TH1F> >() );
567  histogramsSIG.push_back( std::vector< std::unique_ptr<TH1F> >() );
568 
569  for(unsigned int pixel=0; pixel<3; pixel++){
570  TString title = "FE"+std::to_string(FE)+"_pixType"+std::to_string(pixel);
571  std::unique_ptr<TH1F> hthr = std::make_unique<TH1F>(title+"_thr", title+"_thr", m_thrnbins, m_thrLo, m_thrHi);
572  hthr->SetDirectory(0);
573  std::unique_ptr<TH1F> hsig = std::make_unique<TH1F>(title+"_sig", title+"_sig", m_thrnbins, m_sigLo, m_sigHi);
574  hsig->SetDirectory(0);
575  //cppcheck-suppress containerOutOfBounds
576  histogramsTHR.at(FE).push_back(std::move(hthr));
577  //cppcheck-suppress containerOutOfBounds
578  histogramsSIG.at(FE).push_back(std::move(hsig));
579  }
580  }
581 
582  //Will start looping over the RODs
583  TIter rodItr = getRodIterator(riThrFile);
584  TKey* rodKey;
585  while ((rodKey=static_cast<TKey*>(rodItr()))) {
586  const TString rodName(rodKey->GetName());
587  TDirectoryFile* rodDir = (TDirectoryFile*)rodKey->ReadObj();
588  TKey* modKey;
589  TIter modItr=getModuleIterator(rodDir);
590 
591  // Looping over the MODs of each ROD
592  while ((modKey=static_cast<TKey*>(modItr()))) {
593  TString modName(modKey->GetName());
594 
595  if ( not moduleInPart(modName)){
596  continue;
597  }
598  if ( not pm.contains(std::string(modName))){
599  printf("Error - Module %s not found in the PixelMapping tool\n",modName.Data());
600  continue;
601  }
602 
603  if( m_runOneMOD and strcmp(modName, m_testMOD) != 0){
604  continue;
605  }
606 
607  // pixel discriminator threshold
608  std::unique_ptr<TH2F> h2dThr(get2DHistogramFromPath(rodDir,modName, "SCURVE_MEAN"));
609  h2dThr->SetDirectory(0);
610 
611  // Getting histogram for noise
612  std::unique_ptr<TH2F>h2dSig(get2DHistogramFromPath(rodDir,modName, "SCURVE_SIGMA"));
613  h2dSig->SetDirectory(0);
614 
615  for (unsigned int ieta = 0; ieta < m_etaBins; ieta++) {
616  for (unsigned int iphi = 0; iphi < m_phiBins; iphi++) {
617 
618  float thr = h2dThr->GetBinContent(ieta + 1, iphi + 1);
619  float sig = h2dSig->GetBinContent(ieta + 1, iphi + 1);
620 
621  if (thr == 0 || thr > 10000 || sig == 0 || sig > 1000) {
622  continue;
623  }
624 
625  int FE = chipId(iphi, ieta);
626  int pixel= pixelType(iphi, ieta);
627 
628  if(FE<0){
629  return false;
630  }
631 
632  histogramsTHR.at(FE).at(pixel)->Fill(thr);
633  histogramsSIG.at(FE).at(pixel)->Fill(sig);
634 
635  }
636  }
637 
638  // Freeing memory of TH2F
639  h2dThr.reset();
640  h2dSig.reset();
641 
642  int modID = pm.getID(std::string(modName));
643  auto itr = map_info.find( modID );
644 
645  // Map should be empty and therefore we need to create the key - if the key is repeated then it will throw an error
646  if (itr == map_info.end()) {
647 
648  map_info[modID] = std::vector<std::unique_ptr<CalibFrontEndInfo>> ();
649 
650  for(unsigned int FE = 0; FE < m_nFE; FE++){
651  TString subdir(((FE < 10) ? "FE0" : "FE") +std::to_string(FE));
652  map_info[modID].push_back( std::unique_ptr<CalibFrontEndInfo>() );
653  std::unique_ptr<CalibFrontEndInfo> p = std::make_unique<CalibFrontEndInfo>(modID,FE,std::string(modName),std::string(rodKey->GetName()));
654  map_info[modID].at(FE) = std::move(p);
655 
656  for(unsigned int pixel=0; pixel<3; pixel++){
657 
658  // Saving information for the calibration
659  int thr_mean = histogramsTHR.at(FE).at(pixel)->GetMean();
660  int thr_rms = histogramsTHR.at(FE).at(pixel)->GetRMS();
661  int sig_mean = histogramsSIG.at(FE).at(pixel)->GetMean();
662 
663  if(m_savefile){
664 
665  m_wFile->cd();
666  if( !m_wFile->Get(rodName+"/"+modName+"/Thresholds/"+subdir) ){
667  m_wFile->mkdir(rodName+"/"+modName+"/Thresholds/"+subdir,rodName);
668  }
669 
670  m_wFile->cd(rodName+"/"+modName+"/Thresholds/"+subdir);
671 
672  histogramsTHR.at(FE).at(pixel)->SetTitle("Threshold;Pixel threshold;Counts");
673  histogramsSIG.at(FE).at(pixel)->SetTitle("Sigma;Pixel sigma;Counts");
674 
675  std::string type = "";
676  if(pixel == 0) type = "normal";
677  else if(pixel == 1) type = "long";
678  else if(pixel == 2) type = "ganged";
679 
680  histogramsTHR.at(FE).at(pixel)->Write(TString("thres_"+type).Data(), TObject::kWriteDelete);
681  histogramsSIG.at(FE).at(pixel)->Write(TString("sigma_"+type).Data(), TObject::kWriteDelete);
682  }
683 
684  // Reset histograms for next front end
685  histogramsTHR.at(FE).at(pixel)->Reset("ICESM");
686  histogramsSIG.at(FE).at(pixel)->Reset("ICESM");
687 
688  if(pixel == 0){ // normal
689  map_info[modID].at(FE)->set_NormalThreshold(thr_mean);
690  map_info[modID].at(FE)->set_NormalRms(thr_rms);
691  map_info[modID].at(FE)->set_NormalNoise(sig_mean);
692  }
693  else if(pixel == 1){ // long
694  map_info[modID].at(FE)->set_LongThreshold(thr_mean);
695  map_info[modID].at(FE)->set_LongRms(thr_rms);
696  map_info[modID].at(FE)->set_LongNoise(sig_mean);
697  }
698  else if(pixel == 2){ // ganged
699  map_info[modID].at(FE)->set_GangedThreshold(thr_mean);
700  map_info[modID].at(FE)->set_GangedRms(thr_rms);
701  map_info[modID].at(FE)->set_GangedNoise(sig_mean);
702  }
703  else{
704  printf("Calib::fillThresholds: ERROR - Bad pixel in Calib::fillThresholds\n");
705  return false;
706  }
707  }
708  }
709 
710  }
711  else{
712  printf("Calib::fillThresholds: ERROR - REPEATED MOD ID: %s! Contact Offline team\n",std::string(modName).c_str());
713  return false;
714  }
715  } // End of MOD loop
716  } // End of ROD loop
717 
718  for(unsigned int FE = 0; FE < m_nFE; FE++) {
719  for(unsigned int pixel=0; pixel<3; pixel++) {
720  histogramsTHR.at(FE).at(pixel).reset();
721  histogramsSIG.at(FE).at(pixel).reset();
722  }
723  }
724 
725  //Closing file - not needed anymore
726  riThrFile.Close();
727  printf("DONE with threshold calibration.\n");
728  return true;
729 }
730 
731 
732 int Calib::chipId(int iphi, int ieta) {
733  int circ = -1;
734  if (iphi < 160) {
735  circ = (int)(ieta / 18);
736  } else {
737  circ = 15 - (int)(ieta / 18);
738  } // FE15, FE14, ... FE8
739 
740  if (circ>15){
741  printf("Error - FE id error: %d, setting it to -1",circ);
742  circ = -1;//error
743  }
744  return circ;
745 }
746 
747 
748 int Calib::pixelType(int iphi, int ieta, bool isForTOT) {
749 
750  // normal pixels ( by default )
751  int pixtype = 0;
752 
753  // define long pixels
754  if (ieta % 18 == 0 || ieta % 18 == 17) {
755  pixtype = 1;
756  }
757  // define ganged pixels
758  if (iphi > 152 && iphi < 160 && iphi % 2 == 1) {
759  pixtype = 2;
760  }
761  if (iphi > 159 && iphi < 167 && iphi % 2 == 0) {
762  pixtype = 2;
763  }
764 
765  if(isForTOT and pixtype == 2 ) pixtype=1;
766  return pixtype;
767 }
768 
769 TIter Calib::getRodIterator(const TFile & inputFile) {
770  TDirectoryFile* scanDir = static_cast<TDirectoryFile*>((static_cast<TKey*>(inputFile.GetListOfKeys()->First()))->ReadObj());
771  TList* rodKeyList = static_cast<TList*>(scanDir->GetListOfKeys());
772  return TIter(rodKeyList);
773 }
774 
775 TIter Calib::getModuleIterator( TDirectoryFile* rodDir) {
776  TList* modKeyList = static_cast<TList*>(rodDir->GetListOfKeys());
777  return TIter(modKeyList);
778 }
779 
780 TH2F* Calib::get2DHistogramFromPath( TDirectoryFile* rodDir, const TString & moduleName, const TString & histName, int charge) {
781  TString suffix = (charge<0) ? ("") : (TString("/C") + charge);
782  TString fullHistoPath = moduleName + "/" + histName + "/A0/B0" + suffix;
783  TDirectoryFile *histDir = static_cast<TDirectoryFile *>(rodDir->GetDirectory(fullHistoPath));
784 
785  if(!histDir){
786  printf("Error - Directory \"%s\" not found. Exiting..\n",fullHistoPath.Data());
787  return nullptr;
788  }
789  TH2F *pTH2 = static_cast<TH2F*>((static_cast<TKey*>(histDir->GetListOfKeys()->First()))->ReadObj());
790  pTH2->SetDirectory(0);
791 
792  return pTH2;
793 }
794 
795 bool Calib::moduleInPart(const TString & modName) {
796  if (modName == "DSP_ERRORS") {
797  return false;
798  }
799  return modName.BeginsWith(m_MODprefixes[m_whichPart]);
800 }
801 
802 
803 std::vector<float> Calib::getParams(const TF1 *f, unsigned int params){
804  std::vector<float> v;
805  for(unsigned int i = 0; i<params; i++){
806  v.push_back(f->GetParameter(i));
807  }
808 
809  return v;
810 }
811 
812 std::vector<float> Calib::getParams_quality(const TF1 *f){
813  std::vector<float> v;
814 
815  v.push_back(f->GetChisquare());
816  v.push_back(f->GetNDF());
817 
818  return v;
819 }
820 
821 void Calib::graphTitles(const std::unique_ptr<TGraphErrors> &graph, const std::string &name, const std::string &Yname){
822  graph->SetTitle(TString(name)+";Charge;"+TString(Yname));
823  graph->SetMarkerStyle(20);
824 }
825 
826 
Calib::funcDisp
Definition: Calib.h:126
ClusterMessageType::Data
@ Data
Calib::m_chi_error
static constexpr float m_chi_error
Definition: Calib.h:66
Calib::m_phiBins
static constexpr int m_phiBins
Definition: Calib.h:74
Calib::chipId
int chipId(int iphi, int ieta)
Definition: Calib.cxx:732
Calib::m_nFE
static constexpr int m_nFE
Definition: Calib.h:94
Calib::getRodIterator
TIter getRodIterator(const TFile &inputFile)
Definition: Calib.cxx:769
hotSpotInTAG.suffix
string suffix
Definition: hotSpotInTAG.py:185
AddEmptyComponent.histName
string histName
Definition: AddEmptyComponent.py:64
xAOD::uint8_t
uint8_t
Definition: Muon_v1.cxx:553
Calib::m_chargeArr
static constexpr std::array< float, m_ncharge > m_chargeArr
Definition: Calib.h:96
Calib::pixelType
int pixelType(int iphi, int ieta, bool isForTOT=false)
Definition: Calib.cxx:748
Calib::m_MODprefixes
const std::array< TString, 4 > m_MODprefixes
Definition: Calib.h:70
Calib::m_savefile
bool m_savefile
Definition: Calib.h:60
SH::scanDir
void scanDir(SampleHandler &sh, DiskList &list, const std::string &pattern, const std::string &samplePattern, const std::string &samplePostfix)
effects: scan the given directory and add all subdirectories as samples that contain root files.
Definition: ToolsDiscovery.cxx:45
ReadBchFromCrest.begin
begin
Definition: ReadBchFromCrest.py:80
Calib::m_sigLo
static constexpr float m_sigLo
Definition: Calib.h:79
Calib::m_thrLo
static constexpr float m_thrLo
Definition: Calib.h:77
Calib::m_totsigLo
static constexpr float m_totsigLo
Definition: Calib.h:90
pix::PixelMapping::getID
int getID(const std::string &geographicalID) const
Definition: PixelMapping.cxx:84
Calib::m_timnbins
static constexpr int m_timnbins
Definition: Calib.h:82
Calib::m_totHi
static constexpr float m_totHi
Definition: Calib.h:88
MuonR4::to_string
std::string to_string(const SectorProjector proj)
Definition: MsTrackSeeder.cxx:66
python.TrigEgammaMonitorHelper.TH2F
def TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:45
Calib::getParams_quality
std::vector< float > getParams_quality(const TF1 *f)
Definition: Calib.cxx:812
Calib::m_thrHi
static constexpr float m_thrHi
Definition: Calib.h:78
python.CaloAddPedShiftConfig.type
type
Definition: CaloAddPedShiftConfig.py:42
mergePhysValFiles.end
end
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:92
Calib::m_timHi
static constexpr float m_timHi
Definition: Calib.h:84
Calib::m_timLo
static constexpr float m_timLo
Definition: Calib.h:83
Calib::m_chargeErrArr
static constexpr std::array< float, m_ncharge > m_chargeErrArr
Definition: Calib.h:97
IBLCalibrationConfig.thr
thr
Definition: IBLCalibrationConfig.py:39
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:209
Calib::reFit_normalPix
bool reFit_normalPix(std::vector< float > &params, std::vector< float > &q, std::vector< float > &qerr, std::vector< float > &tot, std::vector< float > &toterr, std::vector< float > &sig, std::vector< float > &sigerr, const unsigned int fe)
Definition: Calib.cxx:314
CaloCondBlobAlgs_fillNoiseFromASCII.inputFile
string inputFile
Definition: CaloCondBlobAlgs_fillNoiseFromASCII.py:16
lumiFormat.i
int i
Definition: lumiFormat.py:85
Calib::m_wFile
std::unique_ptr< TFile > m_wFile
Definition: Calib.h:61
Calib::m_totsigHi
static constexpr float m_totsigHi
Definition: Calib.h:91
covarianceTool.title
title
Definition: covarianceTool.py:542
python.BuildSignatureFlags.sig
sig
Definition: BuildSignatureFlags.py:237
Calib::getModuleIterator
TIter getModuleIterator(TDirectoryFile *rodDir)
Definition: Calib.cxx:775
Calib::m_etaBins
static constexpr int m_etaBins
Definition: Calib.h:73
hist_file_dump.f
f
Definition: hist_file_dump.py:140
pix::PixelMapping
Definition: PixelMapping.h:18
ReadBchFromCrest.modName
modName
Definition: ReadBchFromCrest.py:380
Calib.h
Calib::m_qthresh
static constexpr int m_qthresh
Definition: Calib.h:100
PixelCalibrationConfig.tot
tot
Definition: PixelCalibrationConfig.py:28
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
charge
double charge(const T &p)
Definition: AtlasPID.h:997
Calib::m_thrnbins
static constexpr int m_thrnbins
Definition: Calib.h:76
Calib::m_totsigNBins
static constexpr int m_totsigNBins
Definition: Calib.h:89
Calib::m_sigHi
static constexpr float m_sigHi
Definition: Calib.h:80
python.PyAthena.v
v
Definition: PyAthena.py:154
Calib::graphTitles
void graphTitles(const std::unique_ptr< TGraphErrors > &graph, const std::string &name, const std::string &Yname)
Definition: Calib.cxx:821
Calib::m_totnbins
static constexpr int m_totnbins
Definition: Calib.h:86
h
python.CaloAddPedShiftConfig.int
int
Definition: CaloAddPedShiftConfig.py:45
Calib::m_runOneMOD
bool m_runOneMOD
Definition: Calib.h:63
Calib::fillThresholds
bool fillThresholds(const pix::PixelMapping &pm, const std::string &inThrFile, std::map< unsigned int, std::vector< std::unique_ptr< CalibFrontEndInfo >> > &map_info)
Definition: Calib.cxx:548
Calib::moduleInPart
bool moduleInPart(const TString &modName)
Definition: Calib.cxx:795
Calib::m_ncharge
static constexpr int m_ncharge
Definition: Calib.h:95
Calib::getParams
std::vector< float > getParams(const TF1 *f, unsigned int params)
Definition: Calib.cxx:803
extractSporadic.q
list q
Definition: extractSporadic.py:97
Calib::get2DHistogramFromPath
TH2F * get2DHistogramFromPath(TDirectoryFile *rodDir, const TString &moduleName, const TString &histName, int charge=-1)
Definition: Calib.cxx:780
Calib::totFitting
bool totFitting(const pix::PixelMapping &pm, const std::string &inTimFile, std::map< unsigned int, std::vector< std::unique_ptr< CalibFrontEndInfo >> > &map_info)
Definition: Calib.cxx:15
calibdata.copy
bool copy
Definition: calibdata.py:26
Calib::m_whichPart
int m_whichPart
Definition: Calib.h:69
PowhegControl_ttFCNC_NLO.params
params
Definition: PowhegControl_ttFCNC_NLO.py:226
python.LumiCalcRecover.subdir
subdir
Definition: LumiCalcRecover.py:25
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
Calib::funcTot
Definition: Calib.h:114
pix::PixelMapping::contains
bool contains(const std::string &geographicalID) const
Definition: PixelMapping.cxx:95
Amg::distance
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
Definition: GeoPrimitivesHelpers.h:54
Calib::m_totLo
static constexpr float m_totLo
Definition: Calib.h:87
python.compressB64.c
def c
Definition: compressB64.py:93
Calib::fillTiming
bool fillTiming(const pix::PixelMapping &pm, const std::string &inTimFile, std::map< unsigned int, std::vector< std::unique_ptr< CalibFrontEndInfo >> > &map_info)
Definition: Calib.cxx:368
SiliconTech::pixel
@ pixel
Calib::m_testMOD
TString m_testMOD
Definition: Calib.h:64