StRoot: StEEmcPool/StEEmcTimingMaker/macros/plotEEmcL2Timing.C Source File

00001 TFile  *file  = 0;
00002 TChain *chain = 0;
00003 
00004 #include <vector>
00005 
00006 #include "/afs/rhic.bnl.gov/star/packages/DEV/StRoot/StEEmcUtil/EEfeeRaw/EEdims.h"
00007 
00008 // summary TTree branches created by StEEmcTimingMaker
00009 Int_t   mRunNumber;
00010 Float_t mTowerDelay;
00011 Float_t mMapmtDelay;
00012 Int_t   mTotalYield;
00013 Int_t   mTowerCrateYield[MaxTwCrates];
00014 Int_t   mMapmtCrateYield[MaxMapmtCrates];
00015 Int_t   mTowerChanYield[MaxTwCrates][MaxTwCrateCh];
00016 Int_t   mMapmtChanYield[MaxMapmtCrates][MaxMapmtCrateCh];
00017 Int_t   mTowerMin;
00018 Int_t   mTowerMax;
00019 Int_t   mMapmtMin;
00020 Int_t   mMapmtMax;
00021 
00022 // vectors to hold TGraphs for tower and mapmt crates
00023 Int_t npoints = 0;
00024 
00025 TGraphErrors *towerCrateCurves[MaxTwCrates];
00026 TGraphErrors *mapmtCrateCurves[MaxMapmtCrates];
00027 
00028 TGraphErrors *towerChanCurves[MaxTwCrates][MaxTwCrateCh];
00029 TGraphErrors *mapmtChanCurves[MaxMapmtCrates][MaxMapmtCrateCh];
00030 
00031 // enables printing of output files (gif) for documentation
00032 // figures will appear in same subdirector as input files, 
00033 // specified below.
00034 //Bool_t doprint = true;
00035 Bool_t doprint = false;
00036 
00037 void plotEEmcL2Timing( const Char_t *input_dir="timing_files/")
00038 {  
00039   
00040   // chain files
00041   chainFiles(input_dir);
00042 
00043   // setup branches
00044   setBranches(input_dir);
00045 
00046   // get total number of points
00047   Long64_t nruns = chain -> GetEntries();
00048   npoints=(Int_t)nruns;
00049 
00050   // setup the graphs for each crate and each channel
00051   initGraphs();
00052 
00053   // loop over all runs
00054   for ( Long64_t i = 0; i < nruns; i++ ) 
00055     {
00056       chain->GetEntry(i);
00057             
00058       fillCrates((int)i);
00059       fillChannels((int)i);
00060 
00061     }
00062 
00063   // draw timing scan curves for tower crates
00064   drawCrates();
00065   // for ( Int_t ii=0;ii<MaxTwCrates;ii++ ) towerChannels(ii);
00066   // for ( Int_t ii=0;ii<MaxMapmtCrates;ii++ ) mapmtChannels(ii);
00067 
00068   std::cout << "--------------------------------------------------------" << std::endl;
00069   std::cout << "to view timing curves for any crate" << std::endl;
00070   std::cout << std::endl;
00071   std::cout << "towerChannels(icrate) -- icrate = 0-5 for tower crates 1-6"<<std::endl;
00072   //  std::cout << "mapmtChannels(icrate) -- icrate = 0-47 for mapmt crates 1-48"<<std::endl;
00073   std::cout << "print()               -- make gif/pdf files for all crates and channels"<<std::endl;
00074 
00075 }
00076 // ----------------------------------------------------------------------------
00077 void print()
00078 {
00079   doprint=true;
00080   drawCrates();
00081   for ( Int_t ii=0;ii<MaxTwCrates;ii++ ) drawChannels(ii);
00082   //for ( Int_t ii=0;ii<MaxMapmtCrates;ii++ ) drawMapmt(ii);
00083 }
00084 
00085 // ----------------------------------------------------------------------------
00086 void drawCrates()
00087 {
00088 
00089   // tower crates first
00090   TCanvas *towers=new TCanvas("towers","towers",500,400);
00091   const Char_t *opt[]={"ALP","LP","LP","LP","LP","LP"};
00092 
00093   TLegend *legend=new TLegend(0.125,0.6,0.325,0.85);
00094 
00095   Float_t ymax=0.;
00096   for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
00097     {
00098       towerCrateCurves[icr]->Sort();
00099       towerCrateCurves[icr]->Draw(opt[icr]);
00100       TString crname="tw crate ";crname+=icr+1;
00101       legend->AddEntry( towerCrateCurves[icr], crname, "lp" );
00102       if ( towerCrateCurves[icr]->GetYaxis()->GetXmax() > ymax )
00103         ymax=towerCrateCurves[icr]->GetYaxis()->GetXmax();
00104     }
00105   towerCrateCurves[0]->SetTitle("EEMC Tower Crate Timing Curves");
00106   towerCrateCurves[0]->GetXaxis()->SetTitle("TCD phase[ns]");
00107   TString ytitle=Form("Integral [ped+%i,ped+%i] / N_{events}",mTowerMin,mTowerMax);
00108   towerCrateCurves[0]->GetYaxis()->SetTitle(ytitle);
00109   towerCrateCurves[0]->GetYaxis()->SetRangeUser(0.,ymax);
00110   legend->Draw();
00111 
00112   if(doprint)towers->Print("tower_crates.gif");
00113 
00114 #if 0
00115   // next mapmt crates, 4 crates per canvas
00116   TCanvas *mapmt = 0;
00117 
00118   for ( Int_t icr = 0; icr<MaxMapmtCrates; icr+=4 ) 
00119     {
00120       
00121       TString cname="mapmt";cname+=icr+MinMapmtCrateID;
00122       mapmt=new TCanvas(cname,cname,500,400);
00123       TString title="EEMC Mapmt Crate Timing Curves ";
00124       title+=icr+MinMapmtCrateID;
00125       title+="-";
00126       title+=icr+MinMapmtCrateID+3;
00127       legend=new TLegend(0.625,0.12,0.825,0.375);
00128 
00129       Float_t ymax = 0.0;
00130       for ( Int_t jcr=0;jcr<4;jcr++ ) {
00131         Int_t index=icr+jcr;
00132         Int_t crateid=MinMapmtCrateID+index;
00133 
00134         mapmtCrateCurves[index]->Sort();
00135         mapmtCrateCurves[index]->Draw(opt[jcr]);
00136         if ( mapmtCrateCurves[index]->GetYaxis()->GetXmax() > ymax )
00137           ymax=mapmtCrateCurves[index]->GetYaxis()->GetXmax();
00138         TString crname="mapmt crate ";crname+=crateid;
00139         legend->AddEntry(mapmtCrateCurves[index],crname,"lp");
00140 
00141       }
00142       mapmtCrateCurves[icr]->SetTitle(title);
00143       mapmtCrateCurves[icr]->GetXaxis()->SetTitle("TCD phase[ns]");
00144       mapmtCrateCurves[icr]->GetYaxis()->SetTitle("Integral [ped+25,ped+125] / Nevents");
00145       TString ytitle=Form("Integral [ped+%i,ped+%i] / N_{events}",mMapmtMin,mMapmtMax);
00146       mapmtCrateCurves[icr]->GetYaxis()->SetTitle(ytitle);
00147       mapmtCrateCurves[icr]->GetYaxis()->SetRangeUser(0.,ymax*1.05);
00148       legend->Draw();
00149 
00150       TString oname="mapmt_crates_";
00151       oname+=icr+MinMapmtCrateID;
00152       oname+="_";
00153       oname+=icr+MinMapmtCrateID+4;
00154       if(doprint)mapmt->Print(oname+".gif");
00155 
00156     }
00157 #endif
00158 
00159 }
00160 // ----------------------------------------------------------------------------
00161 void mapmtChannels( Int_t crate )
00162 {
00163   
00164   static const Int_t stride=16;
00165   Int_t crateid = MinMapmtCrateID+crate;
00166 
00167 
00168   TString fname="mapmt-crate-";fname+=crate+MinMapmtCrateID;fname+=".ps";
00169   TCanvas *canvas=new TCanvas("canvas","canvas",850/2,1100/2);
00170   canvas->Divide(1,2);
00171   Int_t icanvas=0;
00172 
00173 
00174   for ( Int_t ich=0;ich<192;ich+=stride ) 
00175     {
00176 
00177       canvas->cd(1+icanvas%2);
00178       icanvas++;
00179 
00180 
00181       TString pname="crate";
00182       pname+=crateid;
00183       pname+="_ch";
00184       pname+=ich;
00185       pname+="-";
00186       pname+=ich+stride-1;
00187 
00188       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"};
00189 
00190       // normalize
00191       Float_t ymax=0.0;
00192       Double_t sum[stride];for ( Int_t jj=0;jj<stride;jj++ )sum[jj]=0.;
00193       Double_t max=0.;
00194       for ( Int_t jch=0;jch<stride;jch++ ) // loop over channels in this one graph
00195         {
00196           Int_t index=ich+jch;
00197           Double_t *Y=mapmtChanCurves[crate][index]->GetY();
00198           for ( Int_t ipoint=0;ipoint<npoints;ipoint++ ) {
00199             if ( Y[ipoint]>ymax ) ymax=Y[ipoint];
00200             sum[jch]+=Y[ipoint];
00201           }
00202           if ( sum[jch]>max ) max=sum[jch];
00203         }
00204       if ( max <= 0. ) continue; // meh?
00205 
00206       
00207 
00208       TLegend *legend=new TLegend(0.55,0.11,0.85,0.525);
00209       for ( Int_t jch=0;jch<stride;jch++ )
00210         {
00211           Int_t index=ich+jch;
00212           mapmtChanCurves[crate][index]->SetMarkerSize(0.75);
00213           // offset X axis of each of these
00214           Double_t *X=mapmtChanCurves[crate][index]->GetX();
00215           Double_t *Y=mapmtChanCurves[crate][index]->GetY();
00216           Double_t *EY=mapmtChanCurves[crate][index]->GetEY();
00217           if ( sum[jch]<= 0. ) continue;
00218           //      std::cout<<"before"<<std::endl;
00219           for ( Int_t ip=0;ip<npoints;ip++ ){
00220             Float_t shift = 0.5+ ((float)jch) - ((float)stride)/2.0;
00221             Double_t yy=Y[ip];
00222             X[ip]-= 0.1*shift;
00223             Y[ip]*=max/sum[jch];
00224             EY[ip]*=max/sum[jch];
00225             //      std::cout << "ip="<<ip<<" y="<<yy<<" y'="<<Y[ip]<<std::endl;
00226           }
00227           mapmtChanCurves[crate][index]->Sort();
00228           if ( !jch ) 
00229             mapmtChanCurves[crate][index]->GetXaxis()->SetRangeUser(0.,ymax*1.05);
00230           mapmtChanCurves[crate][index]->SetMarkerColor(38+jch);
00231           mapmtChanCurves[crate][index]->SetLineColor(38+jch);    
00232           mapmtChanCurves[crate][index]->SetMinimum(0.);          
00233           mapmtChanCurves[crate][index]->Draw(opts[jch]);
00234           
00235           TString label="crate ";label+=crate+1;label+=" chan ";label+=index;
00236           legend->AddEntry(mapmtChanCurves[crate][index],label,"lp");
00237 
00238         }
00239       legend->Draw();
00240       canvas->Update();
00241 
00242       //      if(doprint)c->Print(pname+".gif");
00243       if ( !(icanvas%2) ){
00244         canvas->Print(fname+"(");
00245         canvas->Clear();
00246         canvas->Divide(1,2);
00247       }
00248       //      if(doprint)c->Print(pname+".gif");
00249 
00250     }
00251   canvas->Print(fname+")");
00252   gSystem->Exec(TString("ps2pdf ")+fname);
00253 
00254 }
00255 // ----------------------------------------------------------------------------
00256 void towerChannels( Int_t crate )
00257 {
00258 
00259   static  const  Int_t stride=12;
00260 
00261   TString fname="tower-crate-";fname+=crate+1;fname+=".ps";
00262   TCanvas *canvas=0;
00263   canvas = new TCanvas("canvas","canvas",850/2,1100/2);
00264 
00265   canvas->Divide(1,2);
00266   Int_t icanvas=0;
00267 
00268   for ( Int_t ich=0;ich<120;ich+=stride )
00269     {
00270 
00271       canvas->cd(1+icanvas%2);
00272       icanvas++;
00273 
00274       //    TString aname="crate";aname+=crate;aname+=" channels ";aname+=ich;aname+=" to ";aname+=ich+stride-1;
00275       //    TCanvas *c = new TCanvas(aname,aname,400,300);
00276 
00277       TString pname="crate";
00278       pname+=crate+1;
00279       pname+="_ch";
00280       pname+=ich;
00281       pname+="-";
00282       pname+=ich+stride-1;
00283 
00284       const Char_t *opts[]={"ALP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP","LP"};
00285       
00286       // normalize
00287       Double_t sum[stride];for ( Int_t jj=0;jj<stride;jj++ )sum[jj]=0.;
00288       Double_t max=0.;
00289       for ( Int_t jch=0;jch<stride;jch++ ) // loop over channels in this one graph
00290         {
00291 
00292           Int_t index=ich+jch;
00293           Double_t *Y=towerChanCurves[crate][index]->GetY();
00294           for ( Int_t ipoint=0;ipoint<npoints;ipoint++ ) sum[jch]+=Y[ipoint];
00295           if ( sum[jch]>max ) max=sum[jch];
00296         }
00297       if ( max <= 0. ) continue; // meh?
00298 
00299       TLegend *legend=new TLegend(0.125,0.15,0.325,0.45);
00300       for ( Int_t jch=0;jch<stride;jch++ )
00301         {
00302 
00303           Int_t index=ich+jch;
00304           towerChanCurves[crate][index]->SetMarkerSize(0.75);
00305           // offset X axis of each of these
00306           Double_t *X=towerChanCurves[crate][index]->GetX();
00307           Double_t *Y=towerChanCurves[crate][index]->GetY();
00308           Double_t *EY=towerChanCurves[crate][index]->GetEY();
00309           if ( sum[jch]<= 0. ) continue;
00310           //      std::cout<<"before"<<std::endl;
00311           for ( Int_t ip=0;ip<npoints;ip++ ){
00312             Float_t shift = 0.5+ ((float)jch) - ((float)stride)/2.0;
00313             Double_t yy=Y[ip];
00314             X[ip]-= 0.1*shift;
00315             Y[ip]*=max/sum[jch];
00316             EY[ip]*=max/sum[jch];
00317             //      std::cout << "ip="<<ip<<" y="<<yy<<" y'="<<Y[ip]<<std::endl;
00318           }
00319 
00320           towerChanCurves[crate][index]->Sort();
00321           towerChanCurves[crate][index]->SetMinimum(0.);          
00322           towerChanCurves[crate][index]->SetMarkerColor(38+jch);
00323           towerChanCurves[crate][index]->SetLineColor(38+jch);    
00324           towerChanCurves[crate][index]->Draw(opts[jch]);
00325           
00326           TString label="crate ";label+=crate+1;label+=" chan ";label+=index;
00327           legend->AddEntry(towerChanCurves[crate][index],label,"lp");
00328 
00329         }
00330       legend->Draw();
00331       canvas->Update();
00332 
00333       //      if(doprint)c->Print(pname+".gif");
00334       if ( !(icanvas%2) ){
00335         if ( doprint ) canvas->Print(fname+"(");
00336         canvas->Clear();
00337         canvas->Divide(1,2);
00338       }
00339 
00340     }
00341   if ( doprint ) {
00342     canvas->Print(fname+")");
00343     gSystem->Exec(TString("ps2pdf ")+fname);
00344   }
00345 
00346 }
00347 
00348 
00349 
00350 // ----------------------------------------------------------------------------
00351 void fillCrates(Int_t ipoint)
00352 {
00353 #if 1
00354   // loop over tower crates
00355   for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
00356     {
00357       Float_t yield = (Float_t)mTowerCrateYield[icr];
00358       Float_t total = (Float_t)mTotalYield;
00359       if ( total > 10.0 ) {
00360         Float_t eyield = TMath::Sqrt(yield);
00361         Float_t etotal = TMath::Sqrt(total);
00362         Float_t r = yield / total;
00363         Float_t e1 = (yield>0)? eyield/yield : 0.0;
00364         Float_t e2 = etotal/total;
00365         Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
00366         towerCrateCurves[icr]->SetPoint(ipoint, mTowerDelay, r );
00367         towerCrateCurves[icr]->SetPointError( ipoint, 0., er );
00368       }
00369       else {
00370         towerCrateCurves[icr]->SetPoint(ipoint, mTowerDelay, -1.0 );
00371         towerCrateCurves[icr]->SetPointError( ipoint, 0., 0. );
00372       }
00373     }
00374   // loop over mapmt crates
00375   for ( Int_t icr=0;icr<MaxMapmtCrates;icr++ )
00376     {
00377       Float_t yield = (Float_t)mMapmtCrateYield[icr];
00378       Float_t total = (Float_t)mTotalYield;
00379       if ( total > 10.0 ) {
00380         Float_t eyield = TMath::Sqrt(yield);
00381         Float_t etotal = TMath::Sqrt(total);
00382         Float_t r = yield / total;
00383         Float_t e1 = (yield>0)? eyield/yield : 0.0;
00384         Float_t e2 = etotal/total;
00385         Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
00386         mapmtCrateCurves[icr]->SetPoint(ipoint, mMapmtDelay, r );
00387         mapmtCrateCurves[icr]->SetPointError( ipoint, 0., er );
00388       }
00389       else {
00390         mapmtCrateCurves[icr]->SetPoint(ipoint, mMapmtDelay, -1. );
00391         mapmtCrateCurves[icr]->SetPointError( ipoint, 0., 0. ); 
00392       }
00393     }
00394 #endif
00395 }
00396 // ----------------------------------------------------------------------------
00397 void fillChannels(Int_t ipoint)
00398 {
00399 
00400 #if 1
00401   // loop over tower crates
00402   for ( Int_t icr=0;icr<MaxTwCrates;icr++ )
00403     {
00404       for ( Int_t ich=0;ich<MaxTwCrateCh;ich++ ) 
00405         {
00406 
00407           Float_t yield = (Float_t)mTowerChanYield[icr][ich];
00408           Float_t total = (Float_t)mTotalYield;
00409           if ( total > 10.0 ) {
00410             Float_t eyield = TMath::Sqrt(yield);
00411             Float_t etotal = TMath::Sqrt(total);
00412             Float_t r = yield / total;
00413             Float_t e1 = (yield>0)? eyield/yield : 0.0;
00414             Float_t e2 = etotal/total;
00415             Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
00416             towerChanCurves[icr][ich]->SetPoint(ipoint, mTowerDelay, r );
00417             towerChanCurves[icr][ich]->SetPointError( ipoint, 0., er );
00418           }
00419           else {
00420             towerChanCurves[icr][ich]->SetPoint(ipoint, mTowerDelay, -1.0 );
00421             towerChanCurves[icr][ich]->SetPointError( ipoint, 0., 0. );
00422           }
00423         }
00424     }
00425 #endif
00426 #if 1
00427   // loop over mapmt crates
00428   for ( Int_t icr=0;icr<MaxMapmtCrates;icr++ )
00429     {
00430       for ( Int_t ich=0;ich<MaxMapmtCrateCh;ich++ ) 
00431         {
00432 
00433           Float_t yield = (Float_t)mMapmtChanYield[icr][ich];
00434           Float_t total = (Float_t)mTotalYield;
00435           if ( total > 10.0 ) {
00436             Float_t eyield = TMath::Sqrt(yield);
00437             Float_t etotal = TMath::Sqrt(total);
00438             Float_t r = yield / total;
00439             Float_t e1 = (yield>0)? eyield/yield : 0.0;
00440             Float_t e2 = etotal/total;
00441             Float_t er = r * TMath::Sqrt( e1*e1 + e2*e2 );
00442             mapmtChanCurves[icr][ich]->SetPoint(ipoint, mMapmtDelay, r );
00443             mapmtChanCurves[icr][ich]->SetPointError( ipoint, 0., er );
00444           }
00445           else {
00446             mapmtChanCurves[icr][ich]->SetPoint(ipoint, mMapmtDelay, -1.0 );
00447             mapmtChanCurves[icr][ich]->SetPointError( ipoint, 0., 0. );
00448           }
00449         }
00450     }
00451 #endif
00452 }
00453 // ----------------------------------------------------------------------------
00454 void initGraphs()
00455 {
00456 
00457   for ( Int_t i=0;i<MaxTwCrates;i++ ){
00458     towerCrateCurves[i] =  new TGraphErrors(npoints);
00459     towerCrateCurves[i]->SetMarkerStyle(20+i);
00460     towerCrateCurves[i]->SetMarkerColor(i+1);
00461     towerCrateCurves[i]->SetLineColor(i+1);
00462 
00463     for ( Int_t j=0;j<MaxTwCrateCh;j++ )
00464       towerChanCurves[i][j]=(TGraphErrors*)towerCrateCurves[i]->Clone();
00465 
00466   }
00467   for ( Int_t i=0;i<MaxMapmtCrates;i++ ){
00468     mapmtCrateCurves[i]= new TGraphErrors(npoints);
00469     mapmtCrateCurves[i]->SetMarkerStyle(20+i%4);
00470     mapmtCrateCurves[i]->SetMarkerColor(1+i%4);
00471     mapmtCrateCurves[i]->SetLineColor(1+i%4);
00472 
00473     for ( Int_t j=0;j<MaxMapmtCrateCh;j++ )
00474       mapmtChanCurves[i][j]=(TGraphErrors*)mapmtCrateCurves[i]->Clone();
00475 
00476   }
00477 
00478   
00479 
00480 
00481 }
00482 // ----------------------------------------------------------------------------
00483 void chainFiles(const Char_t *path)
00484 {
00485   chain=new TChain("timing","Timing summary");
00486 
00487   TFile *tfile = 0;
00488   std::cout << "chaining files in " << path << std::endl;
00489   TSystemDirectory *dir = new TSystemDirectory("dir",path);
00490   
00491   TIter next( dir->GetListOfFiles() );
00492   TObject *file = 0;
00493   while ( file = (TObject*)next() )
00494     {
00495       TString name=file->GetName();
00496 
00497       // sum the event counter histogram
00498       if ( name.Contains(".root") ) {
00499         // open the TFile and
00500         std::cout << " + " << name << std::endl;
00501         // tfile = TFile::Open(name);
00502         chain->Add(name);
00503       }
00504     }
00505 
00506 }
00507 // ----------------------------------------------------------------------------
00508 void setBranches(const Char_t *dir)
00509 {
00510   
00511   chain->SetBranchAddress("mRunNumber",  &mRunNumber );
00512   chain->SetBranchAddress("mTowerDelay", &mTowerDelay );
00513   chain->SetBranchAddress("mMapmtDelay", &mMapmtDelay );
00514 
00515   chain->SetBranchAddress("mTotalYield",      &mTotalYield );
00516   chain->SetBranchAddress("mTowerCrateYield", &mTowerCrateYield );
00517   chain->SetBranchAddress("mMapmtCrateYield", &mMapmtCrateYield );
00518   chain->SetBranchAddress("mTowerChanYield",  &mTowerChanYield  );
00519   chain->SetBranchAddress("mMapmtChanYield",  &mMapmtChanYield  );
00520 
00521   chain->SetBranchAddress("mTowerMin",&mTowerMin);
00522   chain->SetBranchAddress("mTowerMax",&mTowerMax);
00523   chain->SetBranchAddress("mMapmtMin",&mMapmtMin);
00524   chain->SetBranchAddress("mMapmtMax",&mMapmtMax);
00525 
00526 }
00527