ATLAS Offline Software
Loading...
Searching...
No Matches
TimingClass.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3*/
4
5#include "LArSamplesMon/TimingClass.h"
6#include "LArSamplesMon/Interface.h"
7#include "LArSamplesMon/History.h"
8#include "LArSamplesMon/AbsShapeErrorGetter.h"
9#include "LArSamplesMon/TreeShapeErrorGetter.h"
10#include "LArSamplesMon/Data.h"
11#include "LArSamplesMon/OFC.h"
12#include "LArCafJobs/Geometry.h"
13
14#include "TFile.h"
15
16using namespace TMath;
17using namespace std;
18
19// *************************************** //
20LArSamples::TimingClass::TimingClass()
21// *************************************** //
22{
23
24}
25
26
27
28// **************************************** //
29LArSamples::TimingClass::~TimingClass()
30// **************************************** //
31{
32
33}
34
35
36// ***************************************************************************** //
37void LArSamples::TimingClass::timePerFebAllFebs(const std::string& nrun, const std::string& name)
38// ***************************************************************************** //
39{
40 ofstream file;
41 string fname = "TimingFile" + nrun + "_" + name + ".txt";
42
43 if(FileEmptyCheck(fname))
44 {
45 std::cout<<" +++++ The file:" << fname <<" is empty, time-fit will not be performed +++"<<std::endl;
46 return;
47 }
48 else
49 file.open(fname.c_str(), ios::out);
50
51 for( uint i = 0; i < m_interface->nChannels(); i++ ){
52 if( i%10000 == 0 ) cout << "Processing entry " << i << endl;
53 if( m_interface->historySize(i) == 0 ) continue;
54 const History* history = m_interface->cellHistory(i); // extract the history of the cells of each channel
55 if( !history )continue;
56
57 int calo = Abs(history->cellInfo()->calo());
58 int caloPart = history->cellInfo()->calo();;
59
60 int det = 9;
61 if( caloPart > 0 ) det = 1; //sideA
62 if( caloPart < 0 ) det = 0; //sideC
63
64 int layer = history->cellInfo()->layer();
65 int ft = history->cellInfo()->feedThrough() ;
66 int slot = history->cellInfo()->slot() ;
67
68 for( uint j = 0; j < history->nData(); j++ ){
69 //int gain = history->data(j)->gain();
70 //int NumberRun = history->data(j)->run();
71 //remove bad channels
72 if( history->data(j)->problems(true) != "None" ) continue;
73 //noise cut
74 if( (history->data(j)->energy()/history->data(j)->noise()) < 5 ) continue;
75 //Time cut
76 if( fabs( history->data(j)->ofcTime() ) > 20.00 )continue;
77 //Energy Threshold
78 double t = history->data(j)->ofcTime(); //time
79 double e = history->data(j)->energy(); //energy
80 if( !EnergyThreshold( calo, layer, history->data(j)->quality(), ft, slot, e, t ) )continue;
81 //Fill the timing distributions with energy weighting
82 //if( calo <= 3 && gain == 0 ) h[det][ft][slot]->Fill( t,e );
83 //if( ( calo == 4 && gain == 1 ) || ( calo == 5 && gain == 1 && NumberRun < 216867 ) || ( calo == 5 && gain == 0 && NumberRun >= 216867) )h[det][ft][slot]->Fill( t,e );
84 //Print to file
85 int detec = -1;
86 if( name == "EMB" ) detec = 0;
87 else if( name == "EMEC" ) detec = 1;
88 else if( name == "HEC" ) detec = 2;
89 else if( name == "FCAL" ) detec = 3;
90 file << detec << " " << det << " " << ft << " " << slot << " " << t << " "<< e << endl;
91 }
92 }
93
94
95 return;
96
97}
98
99// ********************************************************************************** //
100bool LArSamples::TimingClass::FileEmptyCheck(const std::string& fname)
101// ********************************************************************************** //
102{
103 std::ifstream file(fname);
104 return file.peek() == std::ifstream::traits_type::eof();
105}
106
107// ********************************************************************************** //
108void LArSamples::TimingClass::fitTimePerFebAllFebs(const std::string& nrun, const std::string& name)
109// ********************************************************************************** //
110{
111
112 for( int d = 0; d < 2; d++ ){ //side
113 for( int ft = 0; ft < 32; ft++ ){ //feedthrough
114 for( int sl = 0; sl < 15; sl++ ) //slot
115 Median[d][ft][sl] = -99.;
116 }
117 }
118
119 TH1F *timefeb[2][32][16];
120 TF1 *fit1 = new TF1( "fit1" , "gaus" , 3 );
121
122 double Hist_entries, Fit_mean, Fit_sigma, Fit_error;
123
124 std::vector< std::vector<double> > myvec;
125 string tfilename;
126 tfilename = "TimingFile" + nrun + "_" + name + ".txt";
127
128 if(FileEmptyCheck(tfilename))
129 {
130 std::cout<<" +++++ The file:" << tfilename <<" is empty, time-fit will not be performed +++"<<std::endl;
131 return ;
132 }
133
134 myvec = readTimingFiles(tfilename);
135
136 string Filename = "OFCTime_PerFEB_" + name + ".root";
137 TFile *FebTime = new TFile( Filename.c_str(), "READ");
138
139 double median2, TotalE = 0.0;
140 int dete = -1, side, feedT, slot;
141
142 std::vector<double> myvector;
143 std::vector<double> myvector2;
144 std::vector<double> energy;
145
146 ofstream ff;
147 string ffname = "mediantest_" + nrun + "_" + name + ".txt";
148 ff.open(ffname.c_str(), ios::out);
149
150 for( uint l = 0; l < (myvec.size() -1); l++ ){
151 float mean = myvec[l][4];
152 double E = myvec[l][5];
153 myvector.push_back(mean);
154 myvector2.push_back(mean);
155 energy.push_back(E);
156 TotalE = TotalE + E;
157 //Fill the vector when the number of slot changes or when it reaches the end of the file
158 if( ( myvec[l][3] != myvec[l+1][3] ) || ( l == myvec.size()-2 ) ){
159 //Get the median of the time for each feb (defined by its slot and ft)
160 median2 = getTimeWeightedMedian(myvector, myvector2, energy, TotalE);
161 myvector.clear();
162 myvector2.clear();
163 energy.clear();
164 TotalE = 0;
165 side = int(myvec[l][1]);
166 feedT = int(myvec[l][2]);
167 slot = int(myvec[l][3]);
168 if( name == "EMB" ) dete = 0;
169 else if( name == "EMEC" ) dete = 1;
170 else if( name == "HEC" ) dete = 2;
171 else if( name == "FCAL" ) dete = 3;
172 ff << dete << " " << side << " " << feedT << " " << slot << " " << median2 << endl;
173 Median[side][feedT][slot] = median2;
174 }
175 }
176
177 string savename = "FEB_time_fitMean_" + nrun + "_" + name + ".txt";
178 ofstream file;
179 file.open(savename.c_str(), ios::out);
180
181 double median_error;
182 //double mean_error;
183 std::vector<double> timevector;
184
185 int det = -1;
186 if( name == "EMB" ) det = 0;
187 else if( name == "EMEC" ) det = 1;
188 else if( name == "HEC" ) det = 2;
189 else if( name == "FCAL" ) det = 3;
190
191 for( int d = 0; d < 2; d++ ){ //sides
192 for( int ft = 0; ft < 32; ft++ ){ //feedthrough
193 for( int sl = 0; sl < 15; sl++ ){ //slot
194
195 ostringstream histname;
196 histname << "h_" << name << "_" << d << "_" << ft << "_" << sl+1;
197 timefeb[d][ft][sl+1] = (TH1F*)FebTime->Get(histname.str().c_str());
198 if( !timefeb[d][ft][sl+1] )
199 break;
200
201 Hist_entries = timefeb[d][ft][sl+1]->GetEntries();
202 //Fit only is histograms are full
203 if( Hist_entries == 0 ) continue;
204
205 timefeb[d][ft][sl+1]->Fit( "fit1","mE","N",-10, 10 );
206
207 //double M = fit1->GetParameter(1);
208 double rms, rmsDist;
209 rms = fit1->GetParameter(2);
210 double n1 = -99., n2 = -99.;
211 rmsDist=timefeb[d][ft][sl+1]->GetRMS();
212 double max;
213 max = timefeb[d][ft][sl+1]->GetXaxis()->GetBinCenter(timefeb[d][ft][sl+1]->GetMaximumBin());
214
215 if( name == "FCAL" ) n1 = n2 = 0.8;
216 else if( name == "EMEC" ){
217 if( sl+1 == 13 && rmsDist > 3.0 ){
218 n1 = 1.5; n2 = 1.;
219 }
220 else if( rmsDist > 3. && sl+1 > 13 ){
221 n1 = 1.5; n2 = 0.3;
222 }
223 else if( sl+1 == 9 && ( ft == 2 || ft == 9 || ft == 15 || ft == 21 ) )
224 n1 = n2 = 3;
225 else if( ( sl+1 == 8 || sl+1 == 9 ) && rmsDist > 3.0 )
226 n1 = n2 = 2;
227 else if( sl+1 == 8 && ft == 13 )
228 n1 = n2 = 3;
229 else
230 n1 = n2 = 1.8;
231 }
232 else if( name == "HEC" ) n1 = n2 = 1.5;
233 if( name == "EMB"){
234 if( sl+1 == 9 || sl+1 == 10 ){
235 n1 = 3; n2 = 2;
236 }
237 else if( sl+1 >= 2 && sl+1 < 6 ){
238 n1 = 1.5; n2 = 2.5;
239 }
240 else if( sl+1 == 14 ){
241 n1 = 2; n2 = 1;
242 }
243 else
244 n1 = n2 = 1.8;
245 }
246
247 timefeb[d][ft][sl+1]->Fit("fit1","mE","N", max-n1*rms, max+n2*rms);
248 Fit_mean = fit1->GetParameter(1);
249 Fit_sigma = fit1->GetParameter(2);
250 Fit_error = fit1->GetParError(1);
251
252 median_error = 1.253*timefeb[d][ft][sl+1]->GetRMS()/sqrt(Hist_entries);
253 //mean_error = timefeb[d][ft][sl+1]->GetRMS()/sqrt(Hist_entries);
254
255 if( Hist_entries <= 50 ){
256 file << det << " " << d << " " << ft << " " << sl+1 << " " << Median[d][ft][sl+1] << " " << median_error << endl;
257 }
258 else if( gMinuit->fCstatu == "SUCCESSFUL" && Fit_error >= Fit_mean && Fit_sigma > 1.5*rmsDist ){
259 file << det << " " << d << " " << ft << " " << sl+1 << " " << Median[d][ft][sl+1] << " " << median_error << endl;
260 }
261 else if( gMinuit->fCstatu != "SUCCESSFUL" ){
262 file << det << " " << d << " " << ft << " " << sl+1 << " " << Median[d][ft][sl+1] << " " << median_error << endl;
263 }
264 else{
265 file << det << " " << d << " " << ft << " " << sl+1 << " " << Fit_mean << " " << Fit_error << endl;
266 }
267
268 }
269 }
270 }
271
272
273 ff.close();
274 file.close();
275
276 return;
277
278}
279
280
281/************************************************************/
282void LArSamples::TimingClass::Time(int dete, const std::string& nrun)
283/************************************************************/
284{
285
286 TH1F *h = new TH1F( Form("h_%d", dete) , Form("h_%d", dete) , 160, -20, 20 );
287 h->Sumw2();
288
289 TH1F *h1 = new TH1F( Form("h1_%d", dete) , Form("h1_%d", dete) , 100, 0, 15000 );
290 TH2F *h2 = new TH2F( Form("h2_%d", dete) , Form("h2_%d", dete) , 100, 0, 15000, 160, -20, 20 );
291
292
293 for( uint i = 0; i < m_interface->nChannels(); i++ ){
294 if( m_interface->historySize(i) == 0 )continue;
295 //Extract the history of the cells of each channel
296 const History* history = m_interface->cellHistory(i);
297 if( !history )continue;
298
299 int calo = Abs( history->cellInfo()->calo() );
300 //int caloPart = history->cellInfo()->calo();
301 int layer = history->cellInfo()->layer();
302 int slot = history->cellInfo()->slot();
303 int ft = history->cellInfo()->feedThrough();
304 //int detIndex = 1;
305 //if( calo == 1 )detIndex = 0;
306 //int side = 1;
307 //if( caloPart < 0 )side = 0;
308 //int channel = history->cellInfo()->channel();
309
310
311 for( uint j = 0; j < history->nData(); j++ ){
312 int Gain = history->data(j)->gain();
313 int NumberRun = history->data(j)->run();
314 //Bad Channels
315 if( history->data(j)->problems(true) != "None" ) continue;
316 //Noise cut
317 if( ( history->data(j)->energy() / history->data(j)->noise() ) < 5 )continue;
318 //Time cut
319 if( fabs( history->data(j)->ofcTime() ) > 20 )continue;
320 //Energy Threshold
321 double t = history->data(j)->ofcTime(); //time
322 double e = history->data(j)->energy(); //energy
323 if( !EnergyThreshold( calo, layer, history->data(j)->quality(), ft, slot, e, t ) )continue;
324 //Fill distributions
325 if( ( calo <= 3 && Gain == 0 ) || ( calo == 4 && Gain == 1 ) || ( calo == 5 && Gain == 1 && NumberRun < 216867 ) || ( calo == 5 && Gain == 0 && NumberRun >= 216867 ) ) h->Fill( t, e );
326 if( ( calo <= 3 && Gain == 0 ) || ( calo == 4 && Gain == 1 ) || ( calo == 5 && Gain == 1 && NumberRun < 216867 ) || ( calo == 5 && Gain == 0 && NumberRun >=216867 ) ) h1->Fill( history->data(j)->quality() );
327 if( ( calo <= 3 && Gain == 0 ) || ( calo == 4 && Gain == 1 ) || ( calo == 5 && Gain == 1 && NumberRun < 216867 ) || ( calo == 5 && Gain == 0 && NumberRun >= 216867) ) h2->Fill( history->data(j)->quality(), t );
328 }
329 }
330
331 TCanvas *c = new TCanvas( Form("c_%d_%s", dete, nrun.c_str()), Form("c_%d_%s", dete, nrun.c_str()), 129,165,700,600 );
332 c->cd();
333 h->SetTitle(Form("hist_%d", dete) );
334 h->SetLineColor(1);
335 h->SetLineWidth(2);
336 h->GetXaxis()->SetTitle("ofcTime [ns]");
337 h->SetMarkerSize(0);
338
339 TF1 *fit = new TF1("fit", "gaus", 3);
340 fit->SetParameter(1, h->GetMean());
341 fit->SetParameter(2, h->GetRMS());
342 fit->SetParLimits(1, -10, 10);
343 fit->SetParLimits(2, 0, 10);
344
345 double mean, rms;
346 mean = h->GetMean();
347 rms = h->GetRMS();
348 h->Fit("fit", "", "", mean-rms, mean+0.5*rms);
349 h->Draw("histsame""E2");
350
351 TFile *fi = new TFile( Form("Time_%d_%s.root", dete, nrun.c_str()), "RECREATE" );
352 h->Write();
353 fit->Write();
354 h1->Write();
355 h2->Write();
356 fi->Close();
357
358 TString path = "Plots/";
359 c->Print( path+Form("Time_%d_%s.eps", dete, nrun.c_str()) );
360 c->Print( path+Form("Time_%d_%s.png", dete, nrun.c_str()) );
361
362 return;
363
364}
365
366// ******************************************************************************* //
367void LArSamples::TimingClass::PlotFebAverageTime24(const std::string& nrun, const std::string& name)
368// ******************************************************************************* //
369{
370 std::string Filename = "FEB_time_fitMean_" + nrun + "_" + name + ".txt";
371
372 std::ifstream f(Filename.c_str(), ios::in);
373 std::vector<double> mean;
374 std::vector<int> side;
375
376 while (!f.eof())
377 {
378 int n1, n2, n3, n4;
379 double n5, n6;
380
381 f >> n1 >> n2 >> n3 >> n4 >> n5 >> n6;
382
383 if (f.good() && !f.bad() && !f.fail())
384 {
385 mean.push_back(n5);
386 side.push_back(n2);
387 }
388 }
389
390 f.close();
391
392 TH1F *h = new TH1F(Form("FEB_Av_%s", name.c_str()), Form("FEB_Av_%s", name.c_str()), 240, -30., 30.);
393 TH1F *sAh = new TH1F(Form("FEB_Av_%s_sideA", name.c_str()), Form("FEB_Av_%s_sideA", name.c_str()), 240, -30., 30.);
394 TH1F *sCh = new TH1F(Form("FEB_Av_%s_sideC", name.c_str()), Form("FEB_Av_%s_sideC", name.c_str()), 240, -30., 30.);
395
396 for (unsigned int i = 0; i < mean.size(); i++)
397 {
398 h->Fill(mean[i]);
399 // Fill sides
400 if (side.at(i) == 1) // A
401 sAh->Fill(mean[i]);
402 else // C
403 sCh->Fill(mean[i]);
404 }
405
406 h->Sumw2();
407 sAh->Sumw2();
408 sCh->Sumw2();
409
410 TCanvas *c = new TCanvas(Form("c_FEB_Av_%s", name.c_str()), Form("c_FEB_Av_%s", name.c_str()), 129, 165, 700, 600);
411 c->cd();
412 c->SetLogy();
413 h->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
414 h->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
415 h->Draw("hist");
416
417 TString path = "Plots/";
418 c->Print(path + Form("t_FEB_%s.png", name.c_str()));
419 c->Print(path + Form("t_FEB_%s.eps", name.c_str()));
420 c->Print(path + Form("t_FEB_%s.pdf", name.c_str()));
421
422 TCanvas *sAc = new TCanvas(Form("c_FEB_Av_%s_sideA", name.c_str()), Form("c_FEB_Av_%s_sideA", name.c_str()), 129, 165, 700, 600);
423 sAc->cd();
424 sAc->SetLogy();
425 sAh->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
426 sAh->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
427 sAh->Draw("hist");
428 sAc->Print(path + Form("t_FEB_%s_sideA.png", name.c_str()));
429
430 TCanvas *sCc = new TCanvas(Form("c_FEB_Av_%s_sideC", name.c_str()), Form("c_FEB_Av_%s_sideC", name.c_str()), 129, 165, 700, 600);
431 sCc->cd();
432 sCc->SetLogy();
433 sCh->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
434 sCh->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
435 sCh->Draw("hist");
436 sCc->Print(path + Form("t_FEB_%s_sideC.png", name.c_str()));
437
438 std::string plotname = "FEB_average_" + nrun + "_" + name + ".root";
439 TFile *fi = new TFile(plotname.c_str(), "RECREATE");
440 h->Write();
441 fi->Close();
442
443 std::string plotnamesA = "FEB_average_" + nrun + "_" + name + "_sideA.root";
444 TFile *fisA = new TFile(plotnamesA.c_str(), "RECREATE");
445 sAh->Write();
446 fisA->Close();
447
448 std::string plotnamesC = "FEB_average_" + nrun + "_" + name + "_sideC.root";
449 TFile *fisC = new TFile(plotnamesC.c_str(), "RECREATE");
450 sCh->Write();
451 fisC->Close();
452
453 return;
454}
455
456// ******************************************************************************* //
457void LArSamples::TimingClass::PlotFebAverageTime(const std::string& nrun, const std::string& name)
458// ******************************************************************************* //
459{
460 string Filename = "FEB_time_fitMean_" + nrun + "_" + name + ".txt";
461
462 ifstream f( Filename.c_str(), ios::in );
463 vector<double> mean;
464 vector<int> side;
465
466 while( !f.eof() ){
467 int n1, n2, n3, n4;
468 double n5, n6;
469
470 f >> n1 >> n2 >> n3 >> n4 >> n5 >> n6;
471
472 if( f.good() && !f.bad() && !f.fail() ){
473 mean.push_back(n5);
474 side.push_back( n2 );
475 }
476 }
477
478 f.close();
479
480 TH1F *h = new TH1F( Form("FEB_Av_%s", name.c_str()), Form("FEB_Av_%s", name.c_str()), 120, -15., 15. );
481 TH1F *sAh = new TH1F( Form("FEB_Av_%s_sideA", name.c_str()), Form("FEB_Av_%s_sideA", name.c_str()), 120, -15., 15. );
482 TH1F *sCh = new TH1F( Form("FEB_Av_%s_sideC", name.c_str()), Form("FEB_Av_%s_sideC", name.c_str()), 120, -15., 15. );
483
484 for( unsigned int i=0; i<mean.size(); i++ ){
485 h->Fill( mean[i] );
486 //Fill sides
487 if( side.at(i) == 1 ) //A
488 sAh->Fill( mean[i] );
489 else //C
490 sCh->Fill( mean[i] );
491 }
492
493 h->Sumw2();
494 sAh->Sumw2();
495 sCh->Sumw2();
496
497 TCanvas *c = new TCanvas( Form("c_FEB_Av_%s", name.c_str()), Form("c_FEB_Av_%s", name.c_str()), 129,165,700,600 );
498 c->cd();
499 c->SetLogy();
500 h->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
501 h->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
502 h->Draw("hist");
503
504 TString path = "Plots/";
505 c->Print( path+ Form("t_FEB_%s.png", name.c_str()) );
506 c->Print( path+ Form("t_FEB_%s.eps", name.c_str()) );
507 c->Print( path+ Form("t_FEB_%s.pdf", name.c_str()) );
508
509 TCanvas *sAc = new TCanvas( Form("c_FEB_Av_%s_sideA", name.c_str()), Form("c_FEB_Av_%s_sideA", name.c_str()), 129,165,700,600 );
510 sAc->cd();
511 sAc->SetLogy();
512 sAh->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
513 sAh->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
514 sAh->Draw("hist");
515 sAc->Print( path+ Form("t_FEB_%s_sideA.png", name.c_str()) );
516
517 TCanvas *sCc = new TCanvas( Form("c_FEB_Av_%s_sideC", name.c_str()), Form("c_FEB_Av_%s_sideC", name.c_str()), 129,165,700,600 );
518 sCc->cd();
519 sCc->SetLogy();
520 sCh->GetXaxis()->SetTitle("<t_{FEB}> [ns]");
521 sCh->GetYaxis()->SetTitle("Number of FEBs / 0.25 ns");
522 sCh->Draw("hist");
523 sCc->Print( path+ Form("t_FEB_%s_sideC.png", name.c_str()) );
524
525
526 std::string plotname = "FEB_average_" + nrun + "_" + name + ".root";
527 TFile *fi = new TFile( plotname.c_str(), "RECREATE" );
528 h->Write();
529 fi->Close();
530
531 return;
532
533}
534
535
536// ******************************************************* //
537void LArSamples::TimingClass::MergeFebTime( const std::string& nrun )
538// ******************************************************* //
539{
540 ofstream mergedfile;
541 string name = "FEB_time_fitMean_" + nrun + ".txt";
542 mergedfile.open( name.c_str(), ios::out );
543
544 string detparts[4] = {"EMB", "EMEC", "HEC", "FCAL"};
545
546 for( int i = 0; i < 4; i++ ){
547 string tmpname = detparts[i];
548 string file = "FEB_time_fitMean_" + nrun + "_" + tmpname + ".txt";
549
550 if(FileEmptyCheck(file))
551 {
552 std::cout<<" +++++ no information for " << tmpname <<", therefore it's timing information is not merged +++"<<std::endl;
553 continue;
554 }
555
556 ifstream f( file.c_str(), ios::in );
557 while( !f.eof() ){
558
559 int n1, n2, n3, n4;
560 double n5, n6;
561
562 f >> n1 >> n2 >> n3 >> n4 >> n5 >> n6;
563
564 if( f.good() && !f.bad() && !f.fail() )
565 mergedfile << n1 << " " << n2 << " " << n3 << " " << n4 << " " << n5 << " " << n6 << endl;
566 }
567
568 f.close();
569 }
570
571 mergedfile.close();
572
573 return;
574
575}
576
577
578// ****************************************************************************** //
579void LArSamples::TimingClass::getFebCorrection( const std::string& nrun )
580// ****************************************************************************** //
581{
582 // The structure of the Feb correction file is different than what we need
583 // The corrections are needed for each slot and not for each feed through
584 // This function rearranges the previously obtained feb correction file to put it in the right format
585
586 for( int i = 0; i < 4; i++ ){
587 for( int j = 0; j < 2; j++ ){
588 for( int k = 0; k < 16; k++ ){
589 for( int m = 0; m < 32; m++ ){
590 param[i][j][m][k] = -99.;
591 error[i][j][m][k] = -99;
592 }
593 }
594 }
595 }
596
597 ofstream febCorrfile;
598 string name = "FEB_Corr_" + nrun + ".txt";
599 febCorrfile.open( name.c_str(), ios::out );
600
601 string name2 = "FEB_time_fitMean_"+nrun+".txt";
602 ifstream f( name2.c_str(), ios::in );
603
604 if( f ){
605 while( !f.eof() ){
606 int n1, n2, n3, n4;
607 double n5, n6;
608
609 f >> n1 >> n2 >> n3 >> n4 >> n5 >> n6;
610 if (n1 < 0 || n1 >= 4 ||
611 n2 < 0 || n2 >= 2 ||
612 n3 < 0 || n3 >= 32 ||
613 n4 < 0 || n4 >= 16)
614 {
615 cerr << "Bad line from timing file " << n1 << " " << n2 << " " << n3 << " " << n4 << endl;
616 continue;
617 }
618
619 if( f.good() && !f.bad() && !f.fail() ){
620 param[n1][n2][n3][n4] = n5;
621 error[n1][n2][n3][n4] = n6;
622 }
623 }
624 }
625
626 f.close();
627
628 string HECFEBID[48] = {"0x3a1a0000", "0x3a520000", "0x3a820000", "0x3ab20000", "0x3a1a8000", "0x3a528000", "0x3a828000", "0x3ab28000", "0x3a1b0000", "0x3a530000", "0x3a830000", "0x3ab30000", "0x3a1b8000", "0x3a538000", "0x3a838000", "0x3ab38000", "0x3a1c0000", "0x3a540000", "0x3a840000", "0x3ab40000", "0x3a1c8000", "0x3a548000", "0x3a848000", "0x3ab48000", "0x3b1a0000", "0x3b520000", "0x3b820000", "0x3bb20000", "0x3b1a8000", "0x3b528000", "0x3b828000", "0x3bb28000", "0x3b1b0000", "0x3b530000", "0x3b830000", "0x3bb30000", "0x3b1b8000", "0x3b538000", "0x3b838000", "0x3bb38000", "0x3b1c0000", "0x3b540000", "0x3b840000", "0x3bb40000", "0x3b1c8000", "0x3b548000", "0x3b848000", "0x3bb48000" };
629 string FCALFEBID[28] = {"0x3a300000", "0x3a308000", "0x3a310000", "0x3a318000", "0x3a320000", "0x3a328000", "0x3a330000", "0x3a340000", "0x3a348000", "0x3a350000", "0x3a358000", "0x3a360000", "0x3a368000", "0x3a370000", "0x3b300000", "0x3b308000", "0x3b310000", "0x3b318000", "0x3b320000", "0x3b328000", "0x3b330000", "0x3b340000", "0x3b348000", "0x3b350000", "0x3b358000", "0x3b360000", "0x3b368000", "0x3b370000" };
630
631 ofstream corrfile;
632 string nameid = "FEB_Corr_FEBID_" + nrun + ".txt";
633 corrfile.open( nameid.c_str(), ios::out );
634
635 int kk = 0;
636 for( int i = 0; i < 4; i++ ){
637 kk = 0;
638 for( int j = 0; j < 2; j++ ){
639 for( int k = 0; k < 16; k++ ){
640 for( int m = 0; m < 32; m++ ){
641 if( param[i][j][m][k] != -99. && error[i][j][m][k] != -99.){
642 febCorrfile << i << " " << j << " " << k << " " << m << " " << param[i][j][m][k] << " " << error[i][j][m][k] << endl;
643 if( i < 2 )corrfile << "0x" << std::hex << (0x38000000|i<<25|j<<24|m<<19|(k-1)<<15) << " " << std::dec << param[i][j][m][k] << endl;
644 else if( i == 2 )corrfile << HECFEBID[kk] << " " << param[i][j][m][k] << endl;
645 else corrfile << FCALFEBID[kk] << " " << param[i][j][m][k] << endl;
646 kk++;
647 }
648 }
649 }
650 }
651 }
652
653 febCorrfile.close();
654 corrfile.close();
655
656 return;
657
658}
659
660
661// *************************************** //
662void LArSamples::TimingClass::PlotFebtime()
663// *************************************** //
664{
665 //Plot the average feb time per feed through and per slot
666 TGraphErrors **g = new TGraphErrors*[15];
667 TGraphErrors **g1 = new TGraphErrors*[15];
668 TGraphErrors **g2 = new TGraphErrors*[15];
669 TGraphErrors **g3 = new TGraphErrors*[15];
670 TGraphErrors **g4 = new TGraphErrors*[15];
671 TGraphErrors **g5 = new TGraphErrors*[15];
672 TGraphErrors **g6 = new TGraphErrors*[15];
673
674 TGraph** gzero = new TGraph*[15];
675 TCanvas **c = new TCanvas*[4];
676 TLatex *l = new TLatex();
677 l->SetTextSize(0.065);
678
679 int count = 0, count1 = 0, count2 = 0, count3 = 0, count4 = 0, count5 = 0, count6 = 0;
680 TString path = "Plots/";
681
682 for( int det = 0; det < 4; det++ ){
683 //Plot EMEC special crates
684 int special = 1;
685 if( det == 1 )special = 2;
686 for( int crates = 0; crates < special; crates++ ){
687 string cratename = "Standard";
688 if( crates == 1 ) cratename = "Special";
689 //
690 for( int side = 0; side < 2; side++ ){
691 c[side] = new TCanvas( Form("c%d%d%s", det, side, cratename.c_str() ), Form("c%d%d%s", det,side, cratename.c_str()), 1100, 600 );
692 c[side]->Divide( 5, 3 );
693
694 for( int i = 0; i < 15; i++ ){
695 g[i] = new TGraphErrors();
696 g1[i] = new TGraphErrors();
697 g2[i] = new TGraphErrors();
698 g3[i] = new TGraphErrors();
699 g4[i] = new TGraphErrors();
700 g5[i] = new TGraphErrors();
701 g6[i] = new TGraphErrors();
702
703 gzero[i] = new TGraph();
704
705 for( int j = 0; j < 32; j++ ){
706 if( det == 1 && crates == 0 && ( j == 2 || j == 3 || j == 9 || j == 10 || j == 15 || j == 16 || j ==21 || j == 22 ) )continue;
707 if( det == 1 && crates == 1 && !( j == 2 || j == 3 || j == 9 || j == 10 || j == 15 || j == 16 || j ==21 || j == 22 ) )continue;
708
709 gzero[i]->SetPoint( j, j, 0 );
710 if( param[det][side][j][i+1] == -99 )continue;
711 //
712 count++;
713 g[i]->SetPoint( count-1, j, param[det][side][j][i+1] );
714 g[i]->SetPointError( count-1, 0, error[det][side][j][i+1] );
715 //
716 if( det == 0 ){
717 g[i]->SetMarkerColor(4);
718 g[i]->GetYaxis()->SetRangeUser(-5, 5);
719 }
720 else{
721 if( j == 6 ){
722 count1++;
723 g1[i]->SetPoint(count1-1, j, param[det][side][j][i+1]);
724 g1[i]->SetPointError(count1-1, 0, error[det][side][j][i+1]);
725 g1[i]->SetMarkerColor(4);
726 }
727 else if( j == 3 || j == 10 || j == 16 || j == 22 ){
728 if( i >= 4 && i <= 9 ){
729 count2++;
730 g2[i]->SetPoint(count2-1, j, param[det][side][j][i+1]);
731 g2[i]->SetPointError(count2-1, 0, error[det][side][j][i+1]);
732 g2[i]->SetMarkerColor(3);
733 }
734 else if( i == 0 || i == 1){
735 count3++;
736 g3[i]->SetPoint(count3-1, j, param[det][side][j][i+1]);
737 g3[i]->SetPointError(count3-1, 0, error[det][side][j][i+1]);
738 g3[i]->SetMarkerColor(6);
739 }
740 }
741 else if( j == 2 || j == 9 || j == 15 || j == 21 ){
742 count4++;
743 g4[i]->SetPoint(count4-1, j, param[det][side][j][i+1]);
744 g4[i]->SetPointError(count4-1, 0, error[det][side][j][i+1]);
745 g4[i]->SetMarkerColor(2);
746 }
747 else if( i <= 12 ){
748 if( j%2 ){
749 count5++;
750 g5[i]->SetPoint(count5-1, j, param[det][side][j][i+1]);
751 g5[i]->SetPointError(count5-1, 0, error[det][side][j][i+1]);
752 g5[i]->SetMarkerColor(1+12*(j%2));
753 }
754 else{
755 count6++;
756 g6[i]->SetPoint(count6-1, j, param[det][side][j][i+1]);
757 g6[i]->SetPointError(count6-1, 0, error[det][side][j][i+1]);
758 g6[i]->SetMarkerColor(1+12*(j%2));
759 }
760 }
761 g[i]->GetYaxis()->SetRangeUser(-5, 5);
762 }
763 }
764
765 count = 0; count1 = 0; count2 = 0; count3 = 0; count4 = 0; count5 = 0; count6 = 0;
766
767 c[side]->cd(i+1);
768 g[i]->GetXaxis()->SetTitle("FT");
769 g[i]->GetYaxis()->SetTitle("ofcTime[ns]");
770 g[i]->GetXaxis()->SetLabelSize(0.065);
771 g[i]->GetYaxis()->SetLabelSize(0.065);
772 g[i]->SetMarkerSize(0.4);
773 g1[i]->SetMarkerSize(0.4);
774 g2[i]->SetMarkerSize(0.4);
775 g3[i]->SetMarkerSize(0.4);
776 g4[i]->SetMarkerSize(0.4);
777 g5[i]->SetMarkerSize(0.4);
778 g6[i]->SetMarkerSize(0.4);
779 g[i]->SetTitle( Form("#color[2]{#bf{Slot %d}}", i+1 ) );
780 g[i]->Draw("AP");
781 g1[i]->Draw("P");
782 g2[i]->Draw("P");
783 g3[i]->Draw("P");
784 g4[i]->Draw("P");
785 g5[i]->Draw("P");
786 g6[i]->Draw("P");
787
788 gzero[i]->SetMarkerSize(0.2);
789 gzero[i]->Draw("P");
790 if( i == 14 ){
791 string detname = "EMB";
792 if( det == 1 )detname = "EMEC";
793 else if( det == 2 )detname = "HEC";
794 else if( det == 3 )detname = "FCAL";
795 c[side]->Print( path+ Form("t_FEB_side%d_%s_%s.png", side, detname.c_str(), cratename.c_str()));
796 }
797 }
798 }
799 }
800 }
801
802 return;
803
804}
805
806
807// **************************************************************************************************************************** //
808bool LArSamples::TimingClass::EnergyThreshold( int calo, int layer, int quality, int ft, int slot, double energy, double time )
809// **************************************************************************************************************************** //
810{
811 bool pass = true;
812 if( quality > 4000 ) pass = false;
813
814 //
815 if( layer == 0 )
816 if( fabs( time ) > 10.00 )pass = false;
817
818 //
819 if( calo == 1 ){
820 if( layer == 0 && energy < 1000. ) pass = false;
821 if( layer == 1 && energy < 1000. ) pass = false;
822 if( layer == 2 && energy < 3000. ) pass = false;
823 if( layer == 3 && energy < 1500. ) pass = false;
824 }
825
826 //
827 if( calo == 2 || calo == 3 ){
828 if( layer == 0 && energy < 1500. ) pass = false;
829 if( layer == 1 ){
830 if( calo == 3 && energy < 3000. ) pass = false;
831 if( calo == 2 && energy < 1000. ) pass = false;
832 }
833
834 if( layer == 2 ){
835 if( calo == 3 && energy < 2000. ) pass = false;
836 if( calo == 2 && energy < 3000. ) pass = false;
837 }
838
839 if( layer == 3 && energy < 2000. ) pass = false;
840 if( layer == 3 && quality > 10000 ) pass = false;
841 if( slot == 13 && ( ft == 2 || ft == 9 || ft == 15 || ft == 21 ) && quality > 4000 ) pass = false;
842 if( slot == 14 && ( ft == 2 || ft == 9 || ft == 15 || ft == 21 ) && quality > 4000 ) pass = false;
843 if( slot == 15 && ( ft == 2 || ft == 9 || ft == 15 || ft == 21 ) && quality > 4000 ) pass = false;
844 }
845
846 if( calo == 4 && energy < 3500. ) pass = false;
847 if( calo == 5 && energy < 10000. ) pass = false;
848
849 return pass;
850
851}
852
853// ******************************************************************************* //
854vector< vector<double> > LArSamples::TimingClass::readTimingFiles(const std::string& file)
855// ******************************************************************************* //
856{
857 std::vector< std::vector <double> > Data;
858
859 ifstream f( file.c_str() );
860
861 TH1F * h[2][32][16];
862
863 string name;
864 if( file.find( "EMB" ) < file.length() )name = "EMB";
865 else if( file.find( "EMEC" ) < file.length() )name = "EMEC";
866 else if( file.find( "HEC" ) < file.length() )name = "HEC";
867 else if( file.find( "FCAL" ) < file.length() )name = "FCAL";
868
869 for( int d = 0; d < 2; d++ ){ //sides
870 for( int ft = 0; ft < 32; ft++ ){ //feedthrough
871 for( int sl = 0; sl < 15; sl++ ){ //slots
872 h[d][ft][sl+1] = new TH1F( Form("h_%s_%d_%d_%d", name.c_str(), d, ft, sl+1) , Form("h_%s_%d_%d_%d", name.c_str(),d,ft, sl+1), 200, -20, 20 );
873 h[d][ft][sl+1]->Sumw2();
874 }
875 }
876 }
877
878
879
880 if( !f )cerr << "Timing File " << file << "cannot be read" << endl;
881 else {
882 while( !f.eof() ){
883 int n1, n2, n3, n4;
884 double n5, n6;
885
886 if( f.good() && !f.bad() && !f.fail() ){
887 f >> n1 >> n2 >> n3 >> n4 >> n5 >> n6;
888 if (n2 < 0 || n2 >= 2 ||
889 n3 < 0 || n3 >= 32 ||
890 n4 < 0 || n4 >= 16)
891 {
892 cerr << "Bad line from timing file " << n2 << " " << n3 << " " << n4 << endl;
893 continue;
894 }
895 //
896 Data.push_back(std::vector<double>{static_cast<double>(n1),
897 static_cast<double>(n2),
898 static_cast<double>(n3),
899 static_cast<double>(n4),
900 static_cast<double>(n5),
901 static_cast<double>(n6)});
902 //
903 h[n2][n3][n4]->Fill( n5, n6 );
904 }
905 }
906 }
907
908 f.close();
909
910
911 //Save the Histograms
912 string filename = "OFCTime_PerFEB_" + name + ".root";
913 TFile *ff = new TFile( filename.c_str(), "RECREATE" );
914
915 for( int d = 0; d < 2; d++ ){ //side
916 for( int ft = 0; ft < 32; ft++ ){ //feedthrough
917 for( int sl = 0; sl < 15; sl++ ) //slot
918 h[d][ft][sl+1]->Write();
919 }
920 }
921
922 ff->Close();
923
924
925
926 return Data;
927
928}
929
930
931// ********************************************************************************************************************************************** //
932double LArSamples::TimingClass::getTimeWeightedMedian(std::vector<double> time, const std::vector<double>& time2, const std::vector<double>& weight, double totalW)
933// ********************************************************************************************************************************************* //
934{
935 TGraph *g = new TGraph();
936 double wmedian = -99., weights, w0 = -1., cumulWeights = 0.0;
937 int Size = time.size();
938 bool wcondition = false;
939
940 sort( time.begin(), time.end() );
941
942 for( int i = 0; i < Size; i++ ){
943 int pos = std::find( time2.begin(),time2.end(), time[i] ) - time2.begin();
944 weights = weight[pos] / totalW;
945 cumulWeights = cumulWeights + weights;
946 g->SetPoint( g->GetN(), cumulWeights, time[i] );
947 if( i == 0 ) w0 = weight[pos] / totalW;
948 wcondition = w0 > 0.5;
949 }
950
951 if( wcondition ) wmedian = g->Eval(w0+(1-w0)/2);
952 else wmedian= g->Eval(0.5);
953
954 g->Delete();
955
956 return wmedian;
957
958}
959
uint
unsigned int uint
Definition LArOFPhaseFill.cxx:20
max
#define max(a, b)
Definition cfImp.cxx:41
LArSamples::CellInfo::slot
short slot() const
Definition CellInfo.h:68
LArSamples::CellInfo::layer
short layer() const
Definition CellInfo.h:53
LArSamples::CellInfo::feedThrough
short feedThrough() const
Definition CellInfo.h:65
LArSamples::CellInfo::calo
CaloId calo() const
Definition CellInfo.h:50
LArSamples::Data::gain
CaloGain::CaloGain gain() const
Definition Data.h:85
LArSamples::Data::noise
double noise() const
Definition Data.cxx:64
LArSamples::Data::ofcTime
double ofcTime() const
Definition Data.h:111
LArSamples::Data::energy
double energy() const
Definition Data.h:108
LArSamples::Data::problems
TString problems(bool sayNone=false) const
Definition Data.cxx:135
LArSamples::Data::run
int run() const
Definition Data.cxx:27
LArSamples::Data::quality
double quality() const
Definition Data.h:114
LArSamples::History::cellInfo
const CellInfo * cellInfo() const
Definition History.h:56
LArSamples::History::data
const Data * data(unsigned int i) const
Definition History.cxx:91
LArSamples::History::nData
unsigned int nData() const
Definition History.h:51
LArSamples::TimingClass::readTimingFiles
std::vector< std::vector< double > > readTimingFiles(const std::string &file)
Definition TimingClass.cxx:854
LArSamples::TimingClass::Time
void Time(int dete, const std::string &nrun)
Definition TimingClass.cxx:282
LArSamples::TimingClass::PlotFebAverageTime24
void PlotFebAverageTime24(const std::string &nrun, const std::string &name)
Definition TimingClass.cxx:367
LArSamples::TimingClass::FileEmptyCheck
bool FileEmptyCheck(const std::string &fname)
Definition TimingClass.cxx:100
LArSamples::TimingClass::Median
double Median[2][32][16]
Definition TimingClass.h:50
LArSamples::TimingClass::fitTimePerFebAllFebs
void fitTimePerFebAllFebs(const std::string &nrun, const std::string &name)
Definition TimingClass.cxx:108
LArSamples::TimingClass::MergeFebTime
void MergeFebTime(const std::string &nrun)
Definition TimingClass.cxx:537
LArSamples::TimingClass::getFebCorrection
void getFebCorrection(const std::string &nrun)
Definition TimingClass.cxx:579
LArSamples::TimingClass::param
double param[4][2][32][16]
Definition TimingClass.h:51
LArSamples::TimingClass::PlotFebAverageTime
void PlotFebAverageTime(const std::string &nrun, const std::string &name)
Definition TimingClass.cxx:457
LArSamples::TimingClass::EnergyThreshold
bool EnergyThreshold(int calo, int layer, int quality, int ft, int slot, double energy, double time)
Definition TimingClass.cxx:808
LArSamples::TimingClass::timePerFebAllFebs
void timePerFebAllFebs(const std::string &nrun, const std::string &name)
Definition TimingClass.cxx:37
LArSamples::TimingClass::m_interface
const Interface * m_interface
Definition TimingClass.h:59
LArSamples::TimingClass::getTimeWeightedMedian
double getTimeWeightedMedian(std::vector< double > time, const std::vector< double > &time2, const std::vector< double > &weight, double totalW)
Definition TimingClass.cxx:932
mean
void mean(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
Definition dependence.cxx:254
count
int count(std::string s, const std::string &regx)
count how many occurances of a regx are in a string
Definition hcg.cxx:146
std
STL namespace.
std::sort
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
Definition DVL_algorithms.h:554
file
TFile * file
Definition tile_monitor.h:29
Generated on for ATLAS Offline Software by doxygen 1.14.0