StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
plotEEmcL2Timing.C
1 TFile *file = 0;
2 TChain *chain = 0;
3 
4 #include <vector>
5 
6 #include "/afs/rhic.bnl.gov/star/packages/DEV/StRoot/StEEmcUtil/EEfeeRaw/EEdims.h"
7 
8 // summary TTree branches created by StEEmcTimingMaker
9 Int_t mRunNumber;
10 Float_t mTowerDelay;
11 Float_t mMapmtDelay;
12 Int_t mTotalYield;
13 Int_t mTowerCrateYield[MaxTwCrates];
14 Int_t mMapmtCrateYield[MaxMapmtCrates];
15 Int_t mTowerChanYield[MaxTwCrates][MaxTwCrateCh];
16 Int_t mMapmtChanYield[MaxMapmtCrates][MaxMapmtCrateCh];
17 Int_t mTowerMin;
18 Int_t mTowerMax;
19 Int_t mMapmtMin;
20 Int_t mMapmtMax;
21 
22 // vectors to hold TGraphs for tower and mapmt crates
23 Int_t npoints = 0;
24 
25 TGraphErrors *towerCrateCurves[MaxTwCrates];
26 TGraphErrors *mapmtCrateCurves[MaxMapmtCrates];
27 
28 TGraphErrors *towerChanCurves[MaxTwCrates][MaxTwCrateCh];
29 TGraphErrors *mapmtChanCurves[MaxMapmtCrates][MaxMapmtCrateCh];
30 
31 // enables printing of output files (gif) for documentation
32 // figures will appear in same subdirector as input files,
33 // specified below.
34 //Bool_t doprint = true;
35 Bool_t doprint = false;
36 
37 void plotEEmcL2Timing( const Char_t *input_dir="timing_files/")
38 {
39 
40  // chain files
41  chainFiles(input_dir);
42 
43  // setup branches
44  setBranches(input_dir);
45 
46  // get total number of points
47  Long64_t nruns = chain -> GetEntries();
48  npoints=(Int_t)nruns;
49 
50  // setup the graphs for each crate and each channel
51  initGraphs();
52 
53  // loop over all runs
54  for ( Long64_t i = 0; i < nruns; i++ )
55  {
56  chain->GetEntry(i);
57 
58  fillCrates((int)i);
59  fillChannels((int)i);
60 
61  }
62 
63  // draw timing scan curves for tower crates
64  drawCrates();
65  // for ( Int_t ii=0;ii<MaxTwCrates;ii++ ) towerChannels(ii);
66  // for ( Int_t ii=0;ii<MaxMapmtCrates;ii++ ) mapmtChannels(ii);
67 
68  std::cout << "--------------------------------------------------------" << std::endl;
69  std::cout << "to view timing curves for any crate" << std::endl;
70  std::cout << std::endl;
71  std::cout << "towerChannels(icrate) -- icrate = 0-5 for tower crates 1-6"<<std::endl;
72  // std::cout << "mapmtChannels(icrate) -- icrate = 0-47 for mapmt crates 1-48"<<std::endl;
73  std::cout << "print() -- make gif/pdf files for all crates and channels"<<std::endl;
74 
75 }
76 // ----------------------------------------------------------------------------
77 void print()
78 {
79  doprint=true;
80  drawCrates();
81  for ( Int_t ii=0;ii<MaxTwCrates;ii++ ) drawChannels(ii);
82  //for ( Int_t ii=0;ii<MaxMapmtCrates;ii++ ) drawMapmt(ii);
83 }
84 
85 // ----------------------------------------------------------------------------
86 void drawCrates()
87 {
88 
89  // tower crates first
90  TCanvas *towers=new TCanvas("towers","towers",500,400);
91  const Char_t *opt[]={"ALP","LP","LP","LP","LP","LP"};
92 
93  TLegend *legend=new TLegend(0.125,0.6,0.325,0.85);
94 
95  Float_t ymax=0.;
96  for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
97  {
98  towerCrateCurves[icr]->Sort();
99  towerCrateCurves[icr]->Draw(opt[icr]);
100  TString crname="tw crate ";crname+=icr+1;
101  legend->AddEntry( towerCrateCurves[icr], crname, "lp" );
102  if ( towerCrateCurves[icr]->GetYaxis()->GetXmax() > ymax )
103  ymax=towerCrateCurves[icr]->GetYaxis()->GetXmax();
104  }
105  towerCrateCurves[0]->SetTitle("EEMC Tower Crate Timing Curves");
106  towerCrateCurves[0]->GetXaxis()->SetTitle("TCD phase[ns]");
107  TString ytitle=Form("Integral [ped+%i,ped+%i] / N_{events}",mTowerMin,mTowerMax);
108  towerCrateCurves[0]->GetYaxis()->SetTitle(ytitle);
109  towerCrateCurves[0]->GetYaxis()->SetRangeUser(0.,ymax);
110  legend->Draw();
111 
112  if(doprint)towers->Print("tower_crates.gif");
113 
114 #if 0
115  // next mapmt crates, 4 crates per canvas
116  TCanvas *mapmt = 0;
117 
118  for ( Int_t icr = 0; icr<MaxMapmtCrates; icr+=4 )
119  {
120 
121  TString cname="mapmt";cname+=icr+MinMapmtCrateID;
122  mapmt=new TCanvas(cname,cname,500,400);
123  TString title="EEMC Mapmt Crate Timing Curves ";
124  title+=icr+MinMapmtCrateID;
125  title+="-";
126  title+=icr+MinMapmtCrateID+3;
127  legend=new TLegend(0.625,0.12,0.825,0.375);
128 
129  Float_t ymax = 0.0;
130  for ( Int_t jcr=0;jcr<4;jcr++ ) {
131  Int_t index=icr+jcr;
132  Int_t crateid=MinMapmtCrateID+index;
133 
134  mapmtCrateCurves[index]->Sort();
135  mapmtCrateCurves[index]->Draw(opt[jcr]);
136  if ( mapmtCrateCurves[index]->GetYaxis()->GetXmax() > ymax )
137  ymax=mapmtCrateCurves[index]->GetYaxis()->GetXmax();
138  TString crname="mapmt crate ";crname+=crateid;
139  legend->AddEntry(mapmtCrateCurves[index],crname,"lp");
140 
141  }
142  mapmtCrateCurves[icr]->SetTitle(title);
143  mapmtCrateCurves[icr]->GetXaxis()->SetTitle("TCD phase[ns]");
144  mapmtCrateCurves[icr]->GetYaxis()->SetTitle("Integral [ped+25,ped+125] / Nevents");
145  TString ytitle=Form("Integral [ped+%i,ped+%i] / N_{events}",mMapmtMin,mMapmtMax);
146  mapmtCrateCurves[icr]->GetYaxis()->SetTitle(ytitle);
147  mapmtCrateCurves[icr]->GetYaxis()->SetRangeUser(0.,ymax*1.05);
148  legend->Draw();
149 
150  TString oname="mapmt_crates_";
151  oname+=icr+MinMapmtCrateID;
152  oname+="_";
153  oname+=icr+MinMapmtCrateID+4;
154  if(doprint)mapmt->Print(oname+".gif");
155 
156  }
157 #endif
158 
159 }
160 // ----------------------------------------------------------------------------
161 void mapmtChannels( Int_t crate )
162 {
163 
164  static const Int_t stride=16;
165  Int_t crateid = MinMapmtCrateID+crate;
166 
167 
168  TString fname="mapmt-crate-";fname+=crate+MinMapmtCrateID;fname+=".ps";
169  TCanvas *canvas=new TCanvas("canvas","canvas",850/2,1100/2);
170  canvas->Divide(1,2);
171  Int_t icanvas=0;
172 
173 
174  for ( Int_t ich=0;ich<192;ich+=stride )
175  {
176 
177  canvas->cd(1+icanvas%2);
178  icanvas++;
179 
180 
181  TString pname="crate";
182  pname+=crateid;
183  pname+="_ch";
184  pname+=ich;
185  pname+="-";
186  pname+=ich+stride-1;
187 
188  const Char_t *opts[]={"ALP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP"};
189 
190  // normalize
191  Float_t ymax=0.0;
192  Double_t sum[stride];for ( Int_t jj=0;jj<stride;jj++ )sum[jj]=0.;
193  Double_t max=0.;
194  for ( Int_t jch=0;jch<stride;jch++ ) // loop over channels in this one graph
195  {
196  Int_t index=ich+jch;
197  Double_t *Y=mapmtChanCurves[crate][index]->GetY();
198  for ( Int_t ipoint=0;ipoint<npoints;ipoint++ ) {
199  if ( Y[ipoint]>ymax ) ymax=Y[ipoint];
200  sum[jch]+=Y[ipoint];
201  }
202  if ( sum[jch]>max ) max=sum[jch];
203  }
204  if ( max <= 0. ) continue; // meh?
205 
206 
207 
208  TLegend *legend=new TLegend(0.55,0.11,0.85,0.525);
209  for ( Int_t jch=0;jch<stride;jch++ )
210  {
211  Int_t index=ich+jch;
212  mapmtChanCurves[crate][index]->SetMarkerSize(0.75);
213  // offset X axis of each of these
214  Double_t *X=mapmtChanCurves[crate][index]->GetX();
215  Double_t *Y=mapmtChanCurves[crate][index]->GetY();
216  Double_t *EY=mapmtChanCurves[crate][index]->GetEY();
217  if ( sum[jch]<= 0. ) continue;
218  // std::cout<<"before"<<std::endl;
219  for ( Int_t ip=0;ip<npoints;ip++ ){
220  Float_t shift = 0.5+ ((float)jch) - ((float)stride)/2.0;
221  Double_t yy=Y[ip];
222  X[ip]-= 0.1*shift;
223  Y[ip]*=max/sum[jch];
224  EY[ip]*=max/sum[jch];
225  // std::cout << "ip="<<ip<<" y="<<yy<<" y'="<<Y[ip]<<std::endl;
226  }
227  mapmtChanCurves[crate][index]->Sort();
228  if ( !jch )
229  mapmtChanCurves[crate][index]->GetXaxis()->SetRangeUser(0.,ymax*1.05);
230  mapmtChanCurves[crate][index]->SetMarkerColor(38+jch);
231  mapmtChanCurves[crate][index]->SetLineColor(38+jch);
232  mapmtChanCurves[crate][index]->SetMinimum(0.);
233  mapmtChanCurves[crate][index]->Draw(opts[jch]);
234 
235  TString label="crate ";label+=crate+1;label+=" chan ";label+=index;
236  legend->AddEntry(mapmtChanCurves[crate][index],label,"lp");
237 
238  }
239  legend->Draw();
240  canvas->Update();
241 
242  // if(doprint)c->Print(pname+".gif");
243  if ( !(icanvas%2) ){
244  canvas->Print(fname+"(");
245  canvas->Clear();
246  canvas->Divide(1,2);
247  }
248  // if(doprint)c->Print(pname+".gif");
249 
250  }
251  canvas->Print(fname+")");
252  gSystem->Exec(TString("ps2pdf ")+fname);
253 
254 }
255 // ----------------------------------------------------------------------------
256 void towerChannels( Int_t crate )
257 {
258 
259  static const Int_t stride=12;
260 
261  TString fname="tower-crate-";fname+=crate+1;fname+=".ps";
262  TCanvas *canvas=0;
263  canvas = new TCanvas("canvas","canvas",850/2,1100/2);
264 
265  canvas->Divide(1,2);
266  Int_t icanvas=0;
267 
268  for ( Int_t ich=0;ich<120;ich+=stride )
269  {
270 
271  canvas->cd(1+icanvas%2);
272  icanvas++;
273 
274  // TString aname="crate";aname+=crate;aname+=" channels ";aname+=ich;aname+=" to ";aname+=ich+stride-1;
275  // TCanvas *c = new TCanvas(aname,aname,400,300);
276 
277  TString pname="crate";
278  pname+=crate+1;
279  pname+="_ch";
280  pname+=ich;
281  pname+="-";
282  pname+=ich+stride-1;
283 
284  const Char_t *opts[]={"ALP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP"};
285 
286  // normalize
287  Double_t sum[stride];for ( Int_t jj=0;jj<stride;jj++ )sum[jj]=0.;
288  Double_t max=0.;
289  for ( Int_t jch=0;jch<stride;jch++ ) // loop over channels in this one graph
290  {
291 
292  Int_t index=ich+jch;
293  Double_t *Y=towerChanCurves[crate][index]->GetY();
294  for ( Int_t ipoint=0;ipoint<npoints;ipoint++ ) sum[jch]+=Y[ipoint];
295  if ( sum[jch]>max ) max=sum[jch];
296  }
297  if ( max <= 0. ) continue; // meh?
298 
299  TLegend *legend=new TLegend(0.125,0.15,0.325,0.45);
300  for ( Int_t jch=0;jch<stride;jch++ )
301  {
302 
303  Int_t index=ich+jch;
304  towerChanCurves[crate][index]->SetMarkerSize(0.75);
305  // offset X axis of each of these
306  Double_t *X=towerChanCurves[crate][index]->GetX();
307  Double_t *Y=towerChanCurves[crate][index]->GetY();
308  Double_t *EY=towerChanCurves[crate][index]->GetEY();
309  if ( sum[jch]<= 0. ) continue;
310  // std::cout<<"before"<<std::endl;
311  for ( Int_t ip=0;ip<npoints;ip++ ){
312  Float_t shift = 0.5+ ((float)jch) - ((float)stride)/2.0;
313  Double_t yy=Y[ip];
314  X[ip]-= 0.1*shift;
315  Y[ip]*=max/sum[jch];
316  EY[ip]*=max/sum[jch];
317  // std::cout << "ip="<<ip<<" y="<<yy<<" y'="<<Y[ip]<<std::endl;
318  }
319 
320  towerChanCurves[crate][index]->Sort();
321  towerChanCurves[crate][index]->SetMinimum(0.);
322  towerChanCurves[crate][index]->SetMarkerColor(38+jch);
323  towerChanCurves[crate][index]->SetLineColor(38+jch);
324  towerChanCurves[crate][index]->Draw(opts[jch]);
325 
326  TString label="crate ";label+=crate+1;label+=" chan ";label+=index;
327  legend->AddEntry(towerChanCurves[crate][index],label,"lp");
328 
329  }
330  legend->Draw();
331  canvas->Update();
332 
333  // if(doprint)c->Print(pname+".gif");
334  if ( !(icanvas%2) ){
335  if ( doprint ) canvas->Print(fname+"(");
336  canvas->Clear();
337  canvas->Divide(1,2);
338  }
339 
340  }
341  if ( doprint ) {
342  canvas->Print(fname+")");
343  gSystem->Exec(TString("ps2pdf ")+fname);
344  }
345 
346 }
347 
348 
349 
350 // ----------------------------------------------------------------------------
351 void fillCrates(Int_t ipoint)
352 {
353 #if 1
354  // loop over tower crates
355  for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
356  {
357  Float_t yield = (Float_t)mTowerCrateYield[icr];
358  Float_t total = (Float_t)mTotalYield;
359  if ( total > 10.0 ) {
360  Float_t eyield = TMath::Sqrt(yield);
361  Float_t etotal = TMath::Sqrt(total);
362  Float_t r = yield / total;
363  Float_t e1 = (yield>0)? eyield/yield : 0.0;
364  Float_t e2 = etotal/total;
365  Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
366  towerCrateCurves[icr]->SetPoint(ipoint, mTowerDelay, r );
367  towerCrateCurves[icr]->SetPointError( ipoint, 0., er );
368  }
369  else {
370  towerCrateCurves[icr]->SetPoint(ipoint, mTowerDelay, -1.0 );
371  towerCrateCurves[icr]->SetPointError( ipoint, 0., 0. );
372  }
373  }
374  // loop over mapmt crates
375  for ( Int_t icr=0;icr<MaxMapmtCrates;icr++ )
376  {
377  Float_t yield = (Float_t)mMapmtCrateYield[icr];
378  Float_t total = (Float_t)mTotalYield;
379  if ( total > 10.0 ) {
380  Float_t eyield = TMath::Sqrt(yield);
381  Float_t etotal = TMath::Sqrt(total);
382  Float_t r = yield / total;
383  Float_t e1 = (yield>0)? eyield/yield : 0.0;
384  Float_t e2 = etotal/total;
385  Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
386  mapmtCrateCurves[icr]->SetPoint(ipoint, mMapmtDelay, r );
387  mapmtCrateCurves[icr]->SetPointError( ipoint, 0., er );
388  }
389  else {
390  mapmtCrateCurves[icr]->SetPoint(ipoint, mMapmtDelay, -1. );
391  mapmtCrateCurves[icr]->SetPointError( ipoint, 0., 0. );
392  }
393  }
394 #endif
395 }
396 // ----------------------------------------------------------------------------
397 void fillChannels(Int_t ipoint)
398 {
399 
400 #if 1
401  // loop over tower crates
402  for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
403  {
404  for ( Int_t ich=0;ich<MaxTwCrateCh;ich++ )
405  {
406 
407  Float_t yield = (Float_t)mTowerChanYield[icr][ich];
408  Float_t total = (Float_t)mTotalYield;
409  if ( total > 10.0 ) {
410  Float_t eyield = TMath::Sqrt(yield);
411  Float_t etotal = TMath::Sqrt(total);
412  Float_t r = yield / total;
413  Float_t e1 = (yield>0)? eyield/yield : 0.0;
414  Float_t e2 = etotal/total;
415  Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
416  towerChanCurves[icr][ich]->SetPoint(ipoint, mTowerDelay, r );
417  towerChanCurves[icr][ich]->SetPointError( ipoint, 0., er );
418  }
419  else {
420  towerChanCurves[icr][ich]->SetPoint(ipoint, mTowerDelay, -1.0 );
421  towerChanCurves[icr][ich]->SetPointError( ipoint, 0., 0. );
422  }
423  }
424  }
425 #endif
426 #if 1
427  // loop over mapmt crates
428  for ( Int_t icr=0;icr<MaxMapmtCrates;icr++ )
429  {
430  for ( Int_t ich=0;ich<MaxMapmtCrateCh;ich++ )
431  {
432 
433  Float_t yield = (Float_t)mMapmtChanYield[icr][ich];
434  Float_t total = (Float_t)mTotalYield;
435  if ( total > 10.0 ) {
436  Float_t eyield = TMath::Sqrt(yield);
437  Float_t etotal = TMath::Sqrt(total);
438  Float_t r = yield / total;
439  Float_t e1 = (yield>0)? eyield/yield : 0.0;
440  Float_t e2 = etotal/total;
441  Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
442  mapmtChanCurves[icr][ich]->SetPoint(ipoint, mMapmtDelay, r );
443  mapmtChanCurves[icr][ich]->SetPointError( ipoint, 0., er );
444  }
445  else {
446  mapmtChanCurves[icr][ich]->SetPoint(ipoint, mMapmtDelay, -1.0 );
447  mapmtChanCurves[icr][ich]->SetPointError( ipoint, 0., 0. );
448  }
449  }
450  }
451 #endif
452 }
453 // ----------------------------------------------------------------------------
454 void initGraphs()
455 {
456 
457  for ( Int_t i=0;i<MaxTwCrates;i++ ){
458  towerCrateCurves[i] = new TGraphErrors(npoints);
459  towerCrateCurves[i]->SetMarkerStyle(20+i);
460  towerCrateCurves[i]->SetMarkerColor(i+1);
461  towerCrateCurves[i]->SetLineColor(i+1);
462 
463  for ( Int_t j=0;j<MaxTwCrateCh;j++ )
464  towerChanCurves[i][j]=(TGraphErrors*)towerCrateCurves[i]->Clone();
465 
466  }
467  for ( Int_t i=0;i<MaxMapmtCrates;i++ ){
468  mapmtCrateCurves[i]= new TGraphErrors(npoints);
469  mapmtCrateCurves[i]->SetMarkerStyle(20+i%4);
470  mapmtCrateCurves[i]->SetMarkerColor(1+i%4);
471  mapmtCrateCurves[i]->SetLineColor(1+i%4);
472 
473  for ( Int_t j=0;j<MaxMapmtCrateCh;j++ )
474  mapmtChanCurves[i][j]=(TGraphErrors*)mapmtCrateCurves[i]->Clone();
475 
476  }
477 
478 
479 
480 
481 }
482 // ----------------------------------------------------------------------------
483 void chainFiles(const Char_t *path)
484 {
485  chain=new TChain("timing","Timing summary");
486 
487  TFile *tfile = 0;
488  std::cout << "chaining files in " << path << std::endl;
489  TSystemDirectory *dir = new TSystemDirectory("dir",path);
490 
491  TIter next( dir->GetListOfFiles() );
492  TObject *file = 0;
493  while ( file = (TObject*)next() )
494  {
495  TString name=file->GetName();
496 
497  // sum the event counter histogram
498  if ( name.Contains(".root") ) {
499  // open the TFile and
500  std::cout << " + " << name << std::endl;
501  // tfile = TFile::Open(name);
502  chain->Add(name);
503  }
504  }
505 
506 }
507 // ----------------------------------------------------------------------------
508 void setBranches(const Char_t *dir)
509 {
510 
511  chain->SetBranchAddress("mRunNumber", &mRunNumber );
512  chain->SetBranchAddress("mTowerDelay", &mTowerDelay );
513  chain->SetBranchAddress("mMapmtDelay", &mMapmtDelay );
514 
515  chain->SetBranchAddress("mTotalYield", &mTotalYield );
516  chain->SetBranchAddress("mTowerCrateYield", &mTowerCrateYield );
517  chain->SetBranchAddress("mMapmtCrateYield", &mMapmtCrateYield );
518  chain->SetBranchAddress("mTowerChanYield", &mTowerChanYield );
519  chain->SetBranchAddress("mMapmtChanYield", &mMapmtChanYield );
520 
521  chain->SetBranchAddress("mTowerMin",&mTowerMin);
522  chain->SetBranchAddress("mTowerMax",&mTowerMax);
523  chain->SetBranchAddress("mMapmtMin",&mMapmtMin);
524  chain->SetBranchAddress("mMapmtMax",&mMapmtMax);
525 
526 }
527