PWGPP/TOF/trending/MakeTrendingTOFQAv2.C

   1 /*
   2   fbellini@cern.ch - last update on 09/09/2014
   3   Macro to run the TOF QA trending by accessing the std QA output,
   4   to be mainly used with the automatic scripts to fill the QA repository.
   5   Launch with
   6   aliroot -l -b -q "MakeTrendingTOFQA.C(\"${fullpath}/QAresults.root\", ${run}, ...)
   7   The macro produces a file containing the tree of trending variables and the main plots.
   8   A feature that displays the plots in canvases must be enable when needed.
   9 */
  10
  11
  12 Int_t MakeTrendingTOFQA(char * runlist,
  13                         Int_t year=2010,
  14                         char *period="LHC10c",
  15                         char* pass="cpass1_pass4",
  16                         char* nameSuffix ="_barrel",
  17                         Bool_t isMC=kFALSE,
  18                         Int_t trainId=0,
  19                         Bool_t saveHisto=kTRUE)
  20 {
  21   Int_t filesCounter=0;
  22   if (!runlist) {
  23     printf("Invalid list of runs given as input: nothing done\n");
  24     return 1;
  25   }
  26   Int_t runNumber;
  27   char infile[300];
  28
  29   char trendFileName[100];
  30   //Define trending output
  31   if (trainId==0){
  32     sprintf(trendFileName,"treeTOFQA_%s_%s.root",period,pass);
  33   } else {
  34     sprintf(trendFileName,"treeTOFQA_QA%i_%s_%s.root",trainId,period,pass);
  35   }
  36   TFile * trendFile=new TFile(trendFileName,"recreate");
  37   FILE * files = fopen(runlist, "r") ;
  38
  39   //create chain of treePostQA
  40   Long64_t nentries=100000, firstentry=0;
  41   TChain *chainTree = 0;
  42   chainTree=new TChain("trendTree");
  43
  44   while (fscanf(files,"%d",&runNumber)==1 ){
  45
  46     //get QAtrain output
  47     if (trainId==0){
  48       if (!isMC) sprintf(infile,"alien:///alice/data/%i/%s/000%d/%s/QAresults%s.root",year,period,runNumber,pass,nameSuffix);
  49       else sprintf(infile,"alien:///alice/sim/%i/%s/%d/QAresults%s.root",year,period,runNumber,nameSuffix);
  50     } else{
  51       if (!isMC) sprintf(infile,"alien:///alice/data/%i/%s/000%d/ESDs/%s/QA%i/QAresults%s.root",year,period,runNumber,pass,trainId,nameSuffix);
  52       else sprintf(infile,"alien:///alice/sim/%i/%s/%d/QA%i/QAresults%s.root",year,period,runNumber,trainId,nameSuffix);
  53     }
  54     Printf("============== Opening QA file(s) for run %i =======================\n",runNumber);
  55
  56     //run post-analysis
  57     if (MakeTrendingTOFQAv2(infile,runNumber,isMC,kFALSE,kTRUE,"raw://", saveHisto)==0){
  58       filesCounter++;
  59     } else Printf("Post analysis not run on QA output %s", infile);
  60   }
  61   Printf(":::: Processed %i runs", filesCounter);
  62   return;
  63 }
  64
  65 //---------------------------------------------------------------------------------
  66 Int_t MakeTrendingTOFQAv2(TString qafilename,       //full path of the QA output; set IsOnGrid to prepend "alien://"
  67                           Int_t runNumber,          // run number
  68                           Bool_t isMC=kFALSE,       //MC flag, to disable meaningless checks
  69                           Bool_t canvasE = kFALSE,  //enable display plots on canvas and save png
  70                           Bool_t IsOnGrid = kFALSE, //set to kTRUE to access files on the grid
  71                           TString ocdbStorage = "raw://",
  72                           Bool_t saveHisto=kTRUE) //set the default ocdb storage
  73 {
  74   // macro to generate tree with TOF QA trending variables
  75   // access qa PWGPP output files
  76   if (!qafilename) {
  77     Printf("Error - Invalid input file");
  78     return 1;
  79   }
  80
  81   /*set graphic style*/
  82   gStyle->SetCanvasColor(kWhite);
  83   gStyle->SetFrameFillColor(kWhite);
  84   gStyle->SetFrameBorderMode(0);
  85   gStyle->SetCanvasBorderMode(0);
  86   gStyle->SetTitleFillColor(kWhite);
  87   gStyle->SetTitleBorderSize(0)  ;
  88   gStyle->SetTitleFont(42);
  89   gStyle->SetTextFont(42);
  90   gStyle->SetStatColor(kWhite);
  91   gStyle->SetStatBorderSize(1);
  92   TGaxis::SetMaxDigits(3);
  93   gStyle->SetOptStat(10);
  94
  95   char defaultQAoutput[30]="QAresults.root";
  96   TString treePostFileName=Form("trending_%i.root",runNumber);
  97
  98   if (IsOnGrid) TGrid::Connect("alien://");
  99   TFile * fin = TFile::Open(qafilename,"r");
 100   if (!fin) {
 101     Printf("ERROR: QA output not found. Exiting...\n");
 102     return -1;
 103   } else {
 104     Printf("INFO: QA output file %s open. \n",fin->GetName());
 105   }
 106
 107   //access histograms lists
 108   char tofQAdirName[15]="TOF";
 109   char genListName[15]="base_noPID";
 110   char t0ListName[15]="timeZero_noPID";
 111   char pidListName[15]="pid_noPID";
 112   char trdListName[15]="trd_noPID";
 113   char trgListName[15]="trigger_noPID";
 114
 115   TDirectoryFile * tofQAdir=(TDirectoryFile*)fin->Get(tofQAdirName);
 116   if (!tofQAdir) {
 117     Printf("ERROR: TOF QA directory not present in input file.\n");
 118     return -1;
 119   }
 120   TList * generalList=(TList*)tofQAdir->Get(genListName);
 121   TList  *timeZeroList=(TList*)tofQAdir->Get(t0ListName);
 122   TList  *pidList=(TList*)tofQAdir->Get(pidListName);
 123   TList  *trdList=(TList*)tofQAdir->Get(trdListName);
 124   TList  *trgList=(TList*)tofQAdir->Get(trgListName);
 125
 126   if (!generalList) Printf("WARNING: general QA histograms absent or not accessible\n");
 127   if (!timeZeroList) Printf("WARNING: timeZero QA histograms absent or not accessible\n");
 128   if (!pidList) Printf("WARNING: PID QA histograms absent or not accessible\n");
 129   if (!trdList) Printf("WARNING: QA histograms for TRD checks absent or not accessible\n");
 130   if (!trgList) Printf("WARNING: QA histograms for trigger absent or not accessible\n");
 131
 132   if ( (!generalList) && (!timeZeroList) && (!pidList) ){
 133     printf("ERROR: no QA available \n");
 134     return -1;
 135   }
 136
 137   Printf(":::: Getting post-analysis info for run %i",runNumber);
 138   TFile * trendFile = new TFile(treePostFileName.Data(),"recreate");
 139
 140   Double_t avTime=-9999., peakTime=-9999., spreadTime=-9999., peakTimeErr=-9999., spreadTimeErr=-9999., negTimeRatio=-9999.,
 141     avRawTime=-9999., peakRawTime=-9999., spreadRawTime=-9999., peakRawTimeErr=-9999., spreadRawTimeErr=-9999., avTot=-9999., peakTot=-9999.,spreadTot=-9999.,  peakTotErr=-9999.,spreadTotErr=-9999.,
 142     orphansRatio=-9999., avL=-9999., negLratio=-9999.,
 143     effPt1=-9999., effPt2=-9999., matchEffLinFit1Gev=-9999.,matchEffLinFit1GevErr=-9999.;
 144
 145   Double_t avPiDiffTime=-9999.,peakPiDiffTime=-9999., spreadPiDiffTime=-9999.,peakPiDiffTimeErr=-9999., spreadPiDiffTimeErr=-9999.;
 146
 147   Double_t avT0A=-9999.,peakT0A=-9999., spreadT0A=-9999.,peakT0AErr=-9999., spreadT0AErr=-9999.;
 148   Double_t avT0C=-9999.,peakT0C=-9999., spreadT0C=-9999.,peakT0CErr=-9999., spreadT0CErr=-9999.;
 149   Double_t avT0AC=-9999.,peakT0AC=-9999., spreadT0AC=-9999.,peakT0ACErr=-9999., spreadT0ACErr=-9999.;
 150   Double_t avT0res=-9999.,peakT0res=-9999., spreadT0res=-9999.,peakT0resErr=-9999., spreadT0resErr=-9999.;
 151   Double_t avT0fillRes=-9999.;
 152
 153   Float_t avMulti=0;
 154   Float_t fractionEventsWHits=-9999.;
 155   /*number of good (HW ok && efficient && !noisy) TOF channels from OCDB*/
 156   Double_t goodChannelRatio=0.0;
 157
 158   TTree * ttree=new TTree("trending","tree of trending variables");
 159   ttree->Branch("run",&runNumber,"run/I");
 160   ttree->Branch("avMulti",&avMulti,"avMulti/F");
 161   ttree->Branch("fractionEventsWHits",&fractionEventsWHits,"fractionEventsWHits/F");
 162   ttree->Branch("goodChannelsRatio",&goodChannelRatio,"goodChannelRatio/D");
 163   ttree->Branch("avTime",&avTime,"avTime/D"); //mean time
 164   ttree->Branch("peakTime",&peakTime,"peakTime/D"); //main peak time after fit
 165   ttree->Branch("spreadTime",&spreadTime,"spreadTime/D"); //spread of main peak of time after fit
 166   ttree->Branch("peakTimeErr",&peakTimeErr,"peakTimeErr/D"); //main peak time after fit error
 167   ttree->Branch("spreadTimeErr",&spreadTimeErr,"spreadTimeErr/D"); //spread of main peak of time after fit error
 168   ttree->Branch("negTimeRatio",&negTimeRatio,"negTimeRatio/D"); //negative time ratio
 169
 170   ttree->Branch("avRawTime",&avRawTime,"avRawTime/D"); //mean raw time
 171   ttree->Branch("peakRawTime",&peakRawTime,"peakRawTime/D"); //mean peak of RAW TIME after fit
 172   ttree->Branch("spreadRawTime",&spreadRawTime,"spreadRawTime/D"); //spread of main peak of raw time after fit
 173   ttree->Branch("peakRawTimeErr",&peakRawTimeErr,"peakRawTimeErr/D"); //main peak raw  time after fit error
 174   ttree->Branch("spreadRawTimeErr",&spreadRawTimeErr,"spreadRawTimeErr/D"); //spread of  raw main peak of time after fit error
 175
 176   ttree->Branch("avTot",&avTot,"avTot/D"); //main peak tot
 177   ttree->Branch("peakTot",&peakTot,"peakTot/D"); // main peak of tot after fit
 178   ttree->Branch("spreadTot",&spreadTot,"spreadTot/D"); //spread of main peak of tot after fit
 179   ttree->Branch("peakTotErr",&peakTotErr,"peakTotErr/D"); // main peak of tot after fit
 180   ttree->Branch("spreadTotErr",&spreadTotErr,"spreadTotErr/D"); //spread of main peak of tot after fit
 181
 182   ttree->Branch("orphansRatio",&orphansRatio,"orphansRatio/D"); //orphans ratio
 183
 184   ttree->Branch("avL",&avL,"avL/D"); //mean track length
 185   ttree->Branch("negLratio",&negLratio,"negLratio/D");//ratio of tracks with track length <350 cm
 186   ttree->Branch("effPt1",&effPt1,"effPt1/D");//matching eff at 1 GeV/c
 187   ttree->Branch("effPt2",&effPt2,"effPt2/D"); //matching eff at 2 GeV/c
 188   ttree->Branch("matchEffLinFit1Gev",&matchEffLinFit1Gev,"matchEffLinFit1Gev/D");//matching eff fit param
 189   ttree->Branch("matchEffLinFit1GevErr",&matchEffLinFit1GevErr,"matchEffLinFit1GevErr/D");////matching eff fit param error
 190
 191   ttree->Branch("avPiDiffTime",&avPiDiffTime,"avPiDiffTime/D"); //mean t-texp
 192   ttree->Branch("peakPiDiffTime",&peakPiDiffTime,"peakPiDiffTime/D"); //main peak t-texp after fit
 193   ttree->Branch("spreadPiDiffTime",&spreadPiDiffTime,"spreadPiDiffTime/D"); //spread of main peak t-texp after fit
 194   ttree->Branch("peakPiDiffTimeErr",&peakPiDiffTimeErr,"peakPiDiffTimeErr/D"); //main peak t-texp after fit error
 195   ttree->Branch("spreadPiDiffTimeErr",&spreadPiDiffTimeErr,"spreadPiDiffTimeErr/D"); //spread of main peak of t-texp after fit error
 196
 197   ttree->Branch("avT0A",&avT0A,"avT0A/D"); //main peak t0A
 198   ttree->Branch("peakT0A",&peakT0A,"peakT0A/D"); // main peak of t0A after fit
 199   ttree->Branch("spreadT0A",&spreadT0A,"spreadTot/D"); //spread of main peak of t0A after fit
 200   ttree->Branch("peakT0AErr",&peakT0AErr,"peakT0AErr/D"); // main peak of t0A after fit
 201   ttree->Branch("spreadT0AErr",&spreadT0AErr,"spreadT0AErr/D"); //spread of main peak of t0A after fit
 202
 203   ttree->Branch("avT0C",&avT0C,"avT0C/D"); //main peak t0C
 204   ttree->Branch("peakT0C",&peakT0C,"peakT0C/D"); // main peak of t0C after fit
 205   ttree->Branch("spreadT0C",&spreadT0C,"spreadT0C/D"); //spread of main peak of t0C after fit
 206   ttree->Branch("peakT0CErr",&peakT0CErr,"peakT0CErr/D"); // main peak of t0C after fit
 207   ttree->Branch("spreadT0CErr",&spreadT0CErr,"spreadT0CErr/D"); //spread of main peak of t0C after fit
 208
 209   ttree->Branch("avT0AC",&avT0AC,"avT0AC/D"); //main peak t0AC
 210   ttree->Branch("peakT0AC",&peakT0AC,"peakT0AC/D"); // main peak of t0AC after fit
 211   ttree->Branch("spreadT0AC",&spreadT0AC,"spreadT0AC/D"); //spread of main peak of t0AC after fit
 212   ttree->Branch("peakT0ACErr",&peakT0ACErr,"peakT0ACErr/D"); // main peak of t0AC after fit
 213   ttree->Branch("spreadT0ACErr",&spreadT0ACErr,"spreadT0ACErr/D"); //spread of main peak of t0AC after fit
 214
 215   ttree->Branch("avT0res",&avT0res,"avT0res/D"); //main peak t0AC
 216   ttree->Branch("peakT0res",&peakT0res,"peakT0res/D"); // main peak of t0AC after fit
 217   ttree->Branch("spreadT0res",&spreadT0res,"spreadT0res/D"); //spread of main peak of t0AC after fit
 218   ttree->Branch("peakT0resErr",&peakT0resErr,"peakT0resErr/D"); // main peak of t0AC after fit
 219   ttree->Branch("spreadT0resErr",&spreadT0resErr,"spreadT0resErr/D"); //spread of main peak of t0AC after fit
 220   ttree->Branch("avT0fillRes",&avT0fillRes,"avT0fillRes/D"); //t0 fill res
 221
 222   //save quantities for trending
 223   goodChannelRatio=(Double_t)GetGoodTOFChannelsRatio(runNumber,kFALSE,ocdbStorage);
 224
 225
 226   //--------------------------------- Multiplicity ----------------------------------//
 227
 228   TH1F * hMulti = (TH1F*) generalList->FindObject("hTOFmulti_all");
 229   TH1F* hFractionEventsWhits = new TH1F("hFractionEventsWhits","hFractionEventsWhits;fraction of events with hits (%)",200,0.,100.);
 230   Float_t fraction=0.0;
 231   if (hMulti->GetEntries()>0.0) {
 232     fraction = ((Float_t) hMulti->GetBinContent(1))/((Float_t) hMulti->GetEntries());
 233     avMulti = hMulti->GetMean();
 234   } else { fraction=0.0; }
 235   hFractionEventsWhits->Fill(fraction*100.);
 236
 237   //--------------------------------- T0F signal ----------------------------------//
 238   TH1F * hRawTime = (TH1F*)generalList->FindObject("hRawTime_all");
 239   if ((hRawTime)&&(hRawTime->GetEntries()>0)){
 240     avRawTime=hRawTime->GetMean();
 241     if (!isMC){
 242       hRawTime->Fit("landau","RQ0","",200.,250.);
 243       if (hRawTime->GetFunction("landau")) {
 244         peakRawTime=(hRawTime->GetFunction("landau"))->GetParameter(1);
 245         spreadRawTime=(hRawTime->GetFunction("landau"))->GetParameter(2);
 246         peakRawTimeErr=(hRawTime->GetFunction("landau"))->GetParError(1);
 247         spreadRawTimeErr=(hRawTime->GetFunction("landau"))->GetParError(2);
 248       }
 249     } else {
 250       printf("Reminder: Raw time not available in MC simulated data.");
 251     }
 252   }
 253   MakeUpHisto(hRawTime, "matched tracks", 21, kGreen+2);
 254
 255   TH1F * hTime = (TH1F*)generalList->FindObject("hTime_all");
 256   if ((hTime)&&(hTime->GetEntries()>0)) {
 257     avTime=hTime->GetMean();
 258     hTime->Fit("landau","RQ0","",0.,50.);
 259     if (hTime->GetFunction("landau")) {
 260       peakTime=(hTime->GetFunction("landau"))->GetParameter(1);
 261       spreadTime=(hTime->GetFunction("landau"))->GetParameter(2);
 262       peakTimeErr=(hTime->GetFunction("landau"))->GetParError(1);
 263       spreadTimeErr=(hTime->GetFunction("landau"))->GetParError(2);
 264       negTimeRatio=((Double_t)hTime->Integral(1,3)*100.)/((Double_t)hTime->Integral());
 265     }
 266     MakeUpHisto(hTime, "matched tracks", 20, kBlue+2);
 267
 268     TLegend *lTime = new TLegend(0.7125881,0.6052519,0.979435,0.7408306,NULL,"brNDC");
 269     lTime->SetTextSize(0.04281433);
 270     lTime->AddEntry(hRawTime, "raw","L");
 271     lTime->AddEntry(hTime, "ESD","L");
 272     lTime->SetFillColor(kWhite);
 273     lTime->SetShadowColor(0);
 274   }
 275
 276   TH1F * hTot = (TH1F*)generalList->FindObject("hTot_all");
 277   if ((hTot)&&(hTot->GetEntries()>0)){
 278     avTot=hTot->GetMean();
 279     hTot->Fit("gaus","","",0.,50.);
 280     if (hTot->GetFunction("gaus")) {
 281       peakTot=(hTot->GetFunction("gaus"))->GetParameter(1);
 282       spreadTot=(hTot->GetFunction("gaus"))->GetParameter(2);
 283       peakTotErr=(hTot->GetFunction("gaus"))->GetParError(1);
 284       spreadTotErr=(hTot->GetFunction("gaus"))->GetParError(2);
 285     }
 286   }
 287   MakeUpHisto(hTot, "matched tracks", 8, kViolet-3);
 288
 289   char orphansTxt[200];
 290   if (hTot->GetEntries()>1){
 291     orphansRatio=((Float_t) hTot->GetBinContent(1))/((Float_t) hTot->GetEntries()) ;
 292   }
 293   sprintf(orphansTxt,"orphans/matched = %4.2f%%",orphansRatio*100.);
 294   TPaveText *tOrphans = new TPaveText(0.38,0.63,0.88,0.7, "NDC");
 295   tOrphans->SetBorderSize(0);
 296   tOrphans->SetTextSize(0.045);
 297   tOrphans->SetFillColor(0); //white background
 298   tOrphans->SetTextAlign(12);
 299   tOrphans->SetTextColor(kViolet-3);
 300   tOrphans->AddText(orphansTxt);
 301
 302   TH1F * hL=(TH1F*)generalList->FindObject("hMatchedL_all");
 303   char negLengthTxt[200];
 304   if (hL->GetEntries()>0){
 305     avL=hL->GetMean();
 306     negLratio=(hL->Integral(1,750))/((Float_t) hL->GetEntries()) ;
 307   }
 308   MakeUpHisto(hL, "matched tracks", 1, kBlue+2);
 309   sprintf(negLengthTxt,"trk with L<350cm /matched = %4.2f%%", negLratio*100.);
 310   TPaveText *tLength = new TPaveText(0.15,0.83,0.65,0.87, "NDC");
 311   tLength->SetBorderSize(0);
 312   tLength->SetTextSize(0.04);
 313   tLength->SetFillColor(0); //white background
 314   tLength->SetTextAlign(11);
 315   tLength->SetTextColor(kOrange-3);
 316   tLength->AddText(negLengthTxt);
 317
 318   //--------------------------------- residuals -------------------------------------//
 319   TH2F* hDxPos4profile = (TH2F*) generalList->FindObject("hMatchedDxVsPt_all");
 320   TH2F* hTOFmatchedDzVsStrip = (TH2F*)generalList->FindObject("hMatchedDzVsStrip_all");
 321
 322   //--------------------------------- matching eff ----------------------------------//
 323   //matching as function of pT
 324   TH1F * hMatchingVsPt = new TH1F("hMatchingVsPt","Matching probability vs. Pt; Pt(GeV/c); matching probability", 50, 0., 5. );
 325   TH1F * hDenom = (TH1F*)generalList->FindObject("hPrimaryPt_all");
 326   if (hDenom) {
 327     hDenom->Sumw2();
 328     hMatchingVsPt=(TH1F*) generalList->FindObject("hMatchedPt_all")->Clone();
 329     hMatchingVsPt->Sumw2();
 330     // hMatchingVsPt->Rebin(5);
 331     // hDenom->Rebin(5);
 332     hMatchingVsPt->Divide(hDenom);
 333     hMatchingVsPt->GetYaxis()->SetTitle("matching efficiency");
 334     hMatchingVsPt->SetTitle("TOF matching efficiency as function of transverse momentum");
 335     hMatchingVsPt->GetYaxis()->SetRangeUser(0,1.2);
 336   }
 337
 338   if (hMatchingVsPt->GetEntries()>0){
 339     hMatchingVsPt->Fit("pol0","","",1.0,10.);
 340     hMatchingVsPt->Draw();
 341     if (hMatchingVsPt->GetFunction("pol0")){
 342       matchEffLinFit1Gev=(hMatchingVsPt->GetFunction("pol0"))->GetParameter(0);
 343       matchEffLinFit1GevErr=(hMatchingVsPt->GetFunction("pol0"))->GetParError(0);
 344       //printf("Matching efficiency fit param is %f +- %f\n",matchEffLinFit1Gev,matchEffLinFit1GevErr );
 345     }
 346   } else {
 347     printf("WARNING: matching efficiency plot has 0 entries. Skipped!\n");
 348   }
 349   MakeUpHisto(hMatchingVsPt, "efficiency", 1, kBlue+2);
 350
 351   //matching as function of eta
 352   TH1F * hMatchingVsEta = new TH1F("hMatchingVsEta","Matching probability vs. #\Eta; #\Eta; matching probability", 20, -1., 1.);
 353   hDenom->Clear();
 354   hDenom=(TH1F*)generalList->FindObject("hPrimaryEta_all");
 355   if (hDenom) {
 356     hDenom->Sumw2();
 357     hMatchingVsEta=(TH1F*) generalList->FindObject("hMatchedEta_all")->Clone();
 358     hMatchingVsEta->Sumw2();
 359       // hMatchingVsEta->Rebin(5);
 360     // hDenom->Rebin(5);
 361     hMatchingVsEta->Divide(hDenom);
 362     hMatchingVsEta->GetXaxis()->SetRangeUser(-1,1);
 363     hMatchingVsEta->GetYaxis()->SetTitle("matching efficiency");
 364     hMatchingVsEta->GetYaxis()->SetRangeUser(0,1.2);
 365     hMatchingVsEta->SetTitle("TOF matching efficiency as function of pseudorapidity");
 366   }
 367   MakeUpHisto(hMatchingVsEta, "efficiency", 1, kBlue+2);
 368
 369   //matching as function of phi
 370   TH1F * hMatchingVsPhi = new TH1F("hMatchingVsPhi","Matching probability vs. Phi; Phi(rad); matching probability", 628, 0., 6.28);
 371   hDenom->Clear();
 372   hDenom=(TH1F*)generalList->FindObject("hPrimaryPhi_all");
 373   if (hDenom) {
 374     hDenom->Sumw2();
 375     hMatchingVsPhi=(TH1F*) generalList->FindObject("hMatchedPhi_all")->Clone();
 376     // hMatchingVsPhi->Rebin(2);
 377     // hDenom->Rebin(2);
 378     hMatchingVsPhi->Sumw2();
 379     hMatchingVsPhi->Divide(hDenom);
 380     hMatchingVsPhi->GetYaxis()->SetTitle("matching efficiency");
 381     hMatchingVsPhi->SetTitle("TOF matching efficiency as function of phi");
 382     hMatchingVsPhi->GetYaxis()->SetRangeUser(0,1.2);
 383   }
 384   MakeUpHisto(hMatchingVsPhi, "efficiency", 1, kBlue+2);
 385
 386    if (saveHisto) {
 387     trendFile->cd();
 388     hMulti->Write();
 389     hTime->Write();
 390     hRawTime->Write();
 391     hTot->Write();
 392     hL->Write();
 393     hDxPos4profile->Write();
 394     hTOFmatchedDzVsStrip->Write();
 395     hMatchingVsPt->Write();
 396     hMatchingVsEta->Write();
 397     hMatchingVsPhi->Write();
 398   }
 399
 400
 401
 402   //--------------------------------- t-texp ----------------------------------//
 403   TH2F * hBetaP=(TH2F*)pidList->FindObject("hMatchedBetaVsP_all");
 404   if (hBetaP) hBetaP->GetYaxis()->SetRangeUser(0.,1.2);
 405
 406   TH1F * hMass=(TH1F*)pidList->FindObject("hMatchedMass_all");
 407   MakeUpHisto(hMass, "tracks", 1, kBlue+2);
 408   // hMass->SetFillColor(kAzure+10);
 409   // hMass->SetFillStyle(1001);
 410   hMass->Rebin(2);
 411
 412   //pions
 413   TH1F * hPionDiff=(TH1F*)pidList->FindObject("hExpTimePi_all");
 414   if ((hPionDiff)&&(hPionDiff->GetEntries()>0)) {
 415     avPiDiffTime=hPionDiff->GetMean();
 416     hPionDiff->Fit("gaus","","",-1000.,500.);
 417     if (hPionDiff->GetFunction("gaus")){
 418       peakPiDiffTime=(hPionDiff->GetFunction("gaus"))->GetParameter(1);
 419       spreadPiDiffTime=(hPionDiff->GetFunction("gaus"))->GetParameter(2);
 420       peakPiDiffTimeErr=(hPionDiff->GetFunction("gaus"))->GetParError(1);
 421       spreadPiDiffTimeErr=(hPionDiff->GetFunction("gaus"))->GetParError(2);
 422       // printf("Main peak t-t_exp (gaus): mean = %f +- %f\n",peakPiDiffTime,peakPiDiffTimeErr );
 423       // printf("Main peak t-t_exp (gaus): spread = %f +- %f\n",spreadPiDiffTime,spreadPiDiffTimeErr );
 424     }
 425   }
 426
 427   TH1F * hDiffTimeT0fillPion=(TH1F*)pidList->FindObject("hExpTimePiFillSub_all");
 428   TH1F * hDiffTimeT0TOFPion1GeV=(TH1F*)pidList->FindObject("hExpTimePiT0Sub1GeV_all");
 429   TH2F * hDiffTimePi=(TH2F*)pidList->FindObject("hExpTimePiVsP_all");
 430   hDiffTimePi->GetYaxis()->SetRangeUser(-5000.,5000.);
 431   hDiffTimePi->SetTitle("PIONS t-t_{exp,#pi} (from tracking) vs. P");
 432
 433   //Kaon
 434   TH2F * hDiffTimeKa=(TH2F*)pidList->FindObject("hExpTimeKaVsP_all");
 435   hDiffTimeKa->SetTitle("KAONS t-t_{exp,K} (from tracking) vs. P");
 436   hDiffTimeKa->GetYaxis()->SetRangeUser(-5000.,5000.);
 437
 438   //Protons
 439   TH2F * hDiffTimePro=(TH2F*)pidList->FindObject("hExpTimeProVsP_all");
 440   hDiffTimePro->SetTitle("PROTONS t-t_{exp,p} (from tracking) vs. P");
 441   hDiffTimePro->GetYaxis()->SetRangeUser(-5000.,5000.);
 442
 443   TH2F * hSigmaPi=(TH2F*)pidList->FindObject("hTOFpidSigmaPi_all");
 444   hSigmaPi->GetYaxis()->SetRangeUser(-5.,5.);
 445   TProfile * profSigmaPi = (TProfile*)hSigmaPi->ProfileX("profSigmaPi");
 446   profSigmaPi->SetLineWidth(2);
 447   profSigmaPi->SetLineColor(kRed+2);
 448
 449   TH2F * hSigmaKa=(TH2F*)pidList->FindObject("hTOFpidSigmaKa_all");
 450   hSigmaKa->GetYaxis()->SetRangeUser(-5.,5.);
 451   TProfile * profSigmaKa = (TProfile*)hSigmaKa->ProfileX("profSigmaKa");
 452   profSigmaKa->SetLineWidth(2);
 453   profSigmaKa->SetLineColor(kBlue);
 454
 455   TH2F * hSigmaPro=(TH2F*)pidList->FindObject("hTOFpidSigmaPro_all");
 456   hSigmaPro->GetYaxis()->SetRangeUser(-5.,5.);
 457   TProfile * profSigmaPro = (TProfile*)hSigmaPro->ProfileX("profSigmaPro");
 458   profSigmaPro->SetLineWidth(2);
 459   profSigmaPro->SetLineColor(kGreen+2);
 460
 461    if (saveHisto) {
 462     trendFile->cd();
 463     hBetaP->Write();
 464     hMass->Write();
 465     hPionDiff->Write();
 466     hDiffTimeT0fillPion->Write();
 467     hDiffTimeT0TOFPion1GeV->Write();
 468     hDiffTimePi->Write();
 469     hDiffTimeKa->Write();
 470     hDiffTimePro->Write();
 471     hSigmaPi->Write();
 472     hSigmaKa->Write();
 473     hSigmaPro->Write();
 474    }
 475    //--------------------------------- T0 detector ----------------------------------//
 476
 477   TH1F*hT0A=(TH1F*)timeZeroList->FindObject("hT0A");
 478   if ((hT0A)&&(hT0A->GetEntries()>0)) {
 479     avT0A = hT0A->GetMean();
 480     hT0A->Fit("gaus","RQ0", "", -1000., 1000.);
 481     if (hT0A->GetFunction("gaus")) {
 482       peakT0A=(hT0A->GetFunction("gaus"))->GetParameter(1);
 483       spreadT0A=(hT0A->GetFunction("gaus"))->GetParameter(2);
 484       peakT0AErr=(hT0A->GetFunction("gaus"))->GetParError(1);
 485       spreadT0AErr=(hT0A->GetFunction("gaus"))->GetParError(2);
 486       // printf("Main peak T0A(gaus): mean = %f +- %f\n",peakT0A,peakT0AErr );
 487       // printf("Main peak T0A (gaus): spread = %f +- %f\n",spreadT0A,spreadT0AErr );
 488       //add integral of main peak over total
 489     }
 490   }
 491   MakeUpHisto(hT0A, "events", 8, kBlue);
 492   hT0A->Rebin(2);
 493
 494   TH1F*hT0C=(TH1F*)timeZeroList->FindObject("hT0C");
 495   if ((hT0C)&&(hT0C->GetEntries()>0)) {
 496     avT0C=hT0C->GetMean();
 497     hT0C->Fit("gaus","RQ0","", -1000., 1000.);
 498     if (hT0C->GetFunction("gaus")) {
 499       peakT0C=(hT0C->GetFunction("gaus"))->GetParameter(1);
 500       spreadT0C=(hT0C->GetFunction("gaus"))->GetParameter(2);
 501       peakT0CErr=(hT0C->GetFunction("gaus"))->GetParError(1);
 502       spreadT0CErr=(hT0C->GetFunction("gaus"))->GetParError(2);
 503       // printf("Main peak T0C(gaus): mean = %f +- %f\n",peakT0C,peakT0CErr );
 504       // printf("Main peak T0C (gaus): spread = %f +- %f\n",spreadT0C,spreadT0CErr );
 505       //add integral of main peak over total
 506     }
 507   }
 508   MakeUpHisto(hT0C, "events", 8, kGreen+1);
 509   hT0C->Rebin(2);
 510
 511   TH1F*hT0AC=(TH1F*)timeZeroList->FindObject("hT0AC");
 512   if ((hT0AC)&&(hT0AC->GetEntries()>0)) {
 513     avT0AC=hT0AC->GetMean();
 514     hT0AC->Fit("gaus","RQ0", "",-1000., 1000.);
 515     if (hT0AC->GetFunction("gaus")) {
 516       peakT0AC=(hT0AC->GetFunction("gaus"))->GetParameter(1);
 517       spreadT0AC=(hT0AC->GetFunction("gaus"))->GetParameter(2);
 518       peakT0ACErr=(hT0AC->GetFunction("gaus"))->GetParError(1);
 519       spreadT0ACErr=(hT0AC->GetFunction("gaus"))->GetParError(2);
 520       // printf("Main peak T0AC(gaus): mean = %f +- %f\n",peakT0AC,peakT0ACErr );
 521       // printf("Main peak T0AC (gaus): spread = %f +- %f\n",spreadT0AC,spreadT0ACErr );
 522     }
 523   }
 524   MakeUpHisto(hT0AC, "events", 8, kRed+1);
 525   hT0AC->Rebin(2);
 526
 527   TLegend *lT0 = new TLegend(0.7125881,0.6052519,0.979435,0.7408306,NULL,"brNDC");
 528   lT0->SetTextSize(0.041);
 529   lT0->AddEntry(hT0AC, "T0 A&C","L");
 530   lT0->AddEntry(hT0A, "T0 A","L");
 531   lT0->AddEntry(hT0C, "T0 C","L");
 532   lT0->SetFillColor(kWhite);
 533   lT0->SetShadowColor(0);
 534
 535   TH1F*hT0res=(TH1F*)timeZeroList->FindObject("hT0DetRes");
 536   if ((hT0res)&&(hT0res->GetEntries()>0)) {
 537     avT0res=hT0res->GetMean();
 538     hT0res->Fit("gaus");
 539     if (hT0res->GetFunction("gaus")) {
 540       peakT0res=(hT0res->GetFunction("gaus"))->GetParameter(1);
 541       spreadT0res=(hT0res->GetFunction("gaus"))->GetParameter(2);
 542       peakT0resErr=(hT0res->GetFunction("gaus"))->GetParError(1);
 543       spreadT0resErr=(hT0res->GetFunction("gaus"))->GetParError(2);
 544       // printf("Main peak T0res(gaus): mean = %f +- %f\n",peakT0res,peakT0resErr );
 545       // printf("Main peak T0res (gaus): spread = %f +- %f\n",spreadT0res,spreadT0resErr );
 546     //add integral of main peak over total
 547     }
 548   }
 549   TH1F*hT0fillRes=(TH1F*)timeZeroList->FindObject("hT0fillRes");
 550   if ((hT0fillRes)&&(hT0fillRes->GetEntries()>0)) {
 551     avT0fillRes=hT0fillRes->GetMean();
 552   }
 553
 554   if (saveHisto) {
 555     trendFile->cd();
 556     hT0AC->Write();
 557     hT0A->Write();
 558     hT0C->Write();
 559     hT0res->Write();
 560     hT0fillRes->Write();
 561   }
 562   //Fill tree and save to file
 563   ttree->Fill();
 564   printf("==============  Saving trending quantities in tree for run %i ===============\n",runNumber);
 565   trendFile->cd();
 566   ttree->Write();
 567   trendFile->Close();
 568
 569   if (canvasE){
 570     // TString plotDir(Form("Plots_run%d",runNumber));
 571     // gSystem->Exec(Form("mkdir %s",plotDir.Data()));
 572     TString plotDir(".");
 573
 574     TCanvas *cTrackProperties= new TCanvas("cTrackProperties","summary of matched tracks properties", 1200, 500);
 575     cTrackProperties->Divide(3,1);
 576     cTrackProperties->cd(1);
 577     gPad->SetLogy();
 578     hTime->Draw("");
 579     hRawTime ->Draw("same");
 580     lTime->Draw();
 581     cTrackProperties->cd(2);
 582     gPad->SetLogy();
 583     hTot->Draw("");
 584     tOrphans->Draw();
 585     cTrackProperties->cd(3);
 586     gPad->SetLogy();
 587     hL->Draw("");
 588     tLength->Draw();
 589
 590     TCanvas *cResiduals= new TCanvas("residuals","residuals", 900,500);
 591     cResiduals->Divide(2,1);
 592     cResiduals->cd(1);
 593     gPad->SetLogz();
 594     hDxPos4profile->GetYaxis()->SetRangeUser(-5.,5.);
 595     hDxPos4profile->Draw("colz");
 596     profDxPos->SetLineColor(kRed);
 597     profDxPos ->Draw("same");
 598     cResiduals->cd(2);
 599     gPad->SetLogz();
 600     hDxNeg4profile->GetYaxis()->SetRangeUser(-5.,5.);
 601     hDxNeg4profile->Draw("colz");
 602     profDxNeg->SetLineColor(kBlue);
 603     profDxNeg->Draw("same");
 604
 605     TCanvas* cProfile = new TCanvas("cProfile","cProfile",50,50, 750,550);
 606     cProfile->cd();
 607     gPad->SetLogz();
 608     hTOFmatchedDzVsStrip->Draw("colz");
 609     Int_t binmin = hTOFmatchedDzVsStrip->GetYaxis()->FindBin(-3);
 610     Int_t binmax = hTOFmatchedDzVsStrip->GetYaxis()->FindBin(3);
 611     TProfile* hDzProfile = (TProfile*)hTOFmatchedDzVsStrip->ProfileX("hDzProfile",binmin, binmax);
 612     hDzProfile->SetLineWidth(3);
 613     hDzProfile->Draw("same");
 614
 615     TCanvas *cMatchingPerformance= new TCanvas("cMatchingPerformance","summary of matching performance",1200,500);
 616     cMatchingPerformance->Divide(3,1);
 617     cMatchingPerformance->cd(1);
 618     hMatchingVsPt->Draw();
 619     cMatchingPerformance->cd(2);
 620     hMatchingVsEta->Draw();
 621     cMatchingPerformance->cd(3);
 622     hMatchingVsPhi->Draw();
 623
 624     TCanvas *cPidPerformance= new TCanvas("cPidPerformance","summary of pid performance", 900,500);
 625     cPidPerformance->Divide(2,1);
 626     cPidPerformance->cd(1);
 627     gPad->SetLogz();
 628     hBetaP->Draw("colz");
 629     cPidPerformance->cd(2);
 630     gPad->SetLogy();
 631     hMass->Draw("HIST ");
 632
 633     TCanvas *cPidPerformance2= new TCanvas("cPidPerformance2","summary of pid performance - expected times", 1200, 500);
 634     cPidPerformance2->Divide(3,1);
 635     cPidPerformance2->cd(1);
 636     gPad->SetLogz();
 637     hDiffTimePi->Draw("colz");
 638     //profDiffTimePi->Draw("same");
 639     cPidPerformance2->cd(2);
 640     gPad->SetLogz();
 641     hDiffTimeKa->Draw("colz");
 642     //profDiffTimeKa->Draw("same");
 643     cPidPerformance2->cd(3);
 644     gPad->SetLogz();
 645     hDiffTimePro->Draw("colz");
 646     //profDiffTimePro->Draw("same");
 647
 648     TLegend * lSigmaPid=new TLegend(0.75,0.75,0.95,0.95,"#sigma_{PID}");
 649     lSigmaPid->AddEntry(profSigmaPi,"#pi^{#pm}","l");
 650     lSigmaPid->AddEntry(profSigmaKa,"K^{#pm}","l");
 651     lSigmaPid->AddEntry(profSigmaPro,"p^{#pm}","l");
 652     TCanvas *cPidPerformance3= new TCanvas("cPidPerformance3","summary of pid performance - sigmas",1200,500);
 653     cPidPerformance3->Divide(3,1);
 654     cPidPerformance3->cd(1);
 655     gPad->SetLogz();
 656     hSigmaPi->Draw("colz");
 657     profSigmaPi->Draw("same");
 658     cPidPerformance3->cd(2);
 659     gPad->SetLogz();
 660     hSigmaKa->Draw("colz");
 661     profSigmaKa->Draw("same");
 662     cPidPerformance3->cd(3);
 663     gPad->SetLogz();
 664     hSigmaPro->Draw("colz");
 665     profSigmaPro->Draw("same");
 666
 667     TCanvas *cT0detector= new TCanvas("cT0detector","T0 detector",800,600);
 668     cT0detector->Divide(2,1);
 669     cT0detector->cd(1);
 670     gPad->SetGridx();
 671     hT0AC->Draw("");
 672     hT0AC->SetTitle("timeZero measured by T0 detector");
 673     hT0A ->Draw("same");
 674     hT0C ->Draw("same");
 675     lT0->Draw();
 676     cT0detector->cd(2);
 677     hT0res->Draw();
 678
 679     cPidPerformance3->Print(Form("%s/%i_PID_sigmas.png", plotDir.Data(), runNumber));
 680     cPidPerformance->Print(Form("%s/%i_PID.png",plotDir.Data(), runNumber));
 681     //cPidPerformanceTh->Print(Form("%s/PID_theoreticalTimes.png",plotDir.Data()));
 682     cPidPerformance2->Print(Form("%s/%i_PID_ExpTimes.png",plotDir.Data(), runNumber));
 683     cMatchingPerformance->Print(Form("%s/%i_Matching.png",plotDir.Data(), runNumber));
 684     cTrackProperties->Print(Form("%s/%i_TrackProperties.png",plotDir.Data(), runNumber));
 685     cResiduals->Print(Form("%s/%i_Residuals.png",plotDir.Data(), runNumber));
 686     cProfile->Print(Form("%s/%i_ProfileDZvsStripNumber.png",plotDir.Data(), runNumber));
 687     cT0detector->Print(Form("%s/%i_T0Detector.png",plotDir.Data(), runNumber));
 688   }
 689
 690   return  0;
 691 }
 692
 693
 694 //----------------------------------------------------------
 695 Double_t GetGoodTOFChannelsRatio(Int_t run = -1, Bool_t saveMap = kFALSE, TString OCDBstorage = "raw://")
 696 {
 697   /*
 698     It retrieves from OCDB the number of good (= efficient && not noisy && HW ok) TOF channels.
 699     Optionally is saves the channel map
 700   */
 701   if (run<=0) {
 702     printf("MakeTrendingTOFqa.C - ERROR in CheckCalibStatus(): invalid run number. Please set a run number.\n");
 703     return 0.0;
 704   }
 705
 706   AliCDBManager *cdb = AliCDBManager::Instance();
 707   cdb->SetDefaultStorage(OCDBstorage.Data());
 708   cdb->SetRun(run);
 709
 710   AliCDBEntry *cdbe = cdb->Get("TOF/Calib/Status");
 711   if (!cdbe) {
 712     printf("MakeTrendingTOFqa.C - ERROR in CheckCalibStatus(): OCDB entry not available. Please, try again.\n");
 713     return 0.0;
 714   }
 715
 716   AliTOFChannelOnlineStatusArray *array = (AliTOFChannelOnlineStatusArray *)cdbe->GetObject();
 717   TH2F *hOkMap = new TH2F("hOkMap", "Ok map (!noisy & !problematic & efficient);sector;strip", 72, 0., 18., 91, 0., 91.);
 718
 719   AliTOFcalibHisto calibHisto;
 720   calibHisto.LoadCalibHisto();
 721   AliTOFcalib calib;
 722   calib.Init();
 723   Int_t sector, sectorStrip, padx, fea;
 724   Float_t hitmapx, hitmapy;
 725   for (Int_t i = 0; i <  array->GetSize(); i++) {
 726     sector = calibHisto.GetCalibMap(AliTOFcalibHisto::kSector, i);
 727     sectorStrip = calibHisto.GetCalibMap(AliTOFcalibHisto::kSectorStrip, i);
 728     padx = calibHisto.GetCalibMap(AliTOFcalibHisto::kPadX, i);
 729     fea = padx / 12;
 730     hitmapx = sector + ((Double_t)(3 - fea) + 0.5) / 4.;
 731     hitmapy = sectorStrip;
 732     if ( !(array->GetNoiseStatus(i) == AliTOFChannelOnlineStatusArray::kTOFNoiseBad)   &&
 733          (calib.IsChannelEnabled(i,kTRUE,kTRUE)))
 734       hOkMap->Fill(hitmapx,hitmapy);
 735   }
 736   Int_t nOk=(Int_t) hOkMap->GetEntries();
 737   Double_t ratioOk=nOk/152928.;
 738   if (saveMap) hOkMap->SaveAs(Form("run%i_OKChannelsMap.root",run));
 739   cout << "###    Run " << run << ": TOF channels ok = " << nOk << "/ total 152928 channels = " << ratioOk*100. << "% of whole TOF" << endl;
 740   return ratioOk;
 741 }
 742
 743 //----------------------------------------------------------
 744 void MakeUpHisto(TH1* histo, TString titleY, Int_t marker=20, Color_t color=kBlue+2)
 745 {
 746   if (!histo) return;
 747   histo->SetMarkerStyle(marker);
 748   histo->SetMarkerSize(0.7);
 749   histo->SetMarkerColor(color);
 750   histo->SetLineColor(color);
 751   histo->SetFillColor(kWhite);
 752   histo->SetFillStyle(0);
 753   histo->GetYaxis()->SetTitle(titleY.Data());
 754   histo->GetYaxis()->SetTitleOffset(1.35);
 755   histo->GetXaxis()->SetLabelSize(0.03);
 756
 757   return;
 758 }