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