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