PWGLF/SPECTRA/Nuclei/B2/macros/DrawSec.C

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 // Draw corrections of secondaries
  17 // author: Eulogio Serradilla <eulogio.serradilla@cern.ch>
  18
  19 #if !defined(__CINT__) || defined(__MAKECINT__)
  20 #include <TROOT.h>
  21 #include <TStyle.h>
  22 #include <TFile.h>
  23 #include <TH1D.h>
  24 #include <TString.h>
  25 #include <TMath.h>
  26 #include <TCanvas.h>
  27 #endif
  28
  29 #include "B2.h"
  30
  31 Int_t GetNumberCols(Int_t lowbin, Int_t hibin);
  32
  33 TCanvas* DrawMC(TH1D** hData, TH1D** hPrim, TH1D** hMat, TH1D** hFdwn, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax);
  34
  35 TCanvas* DrawFit(TH1D** hData, TH1D** hDataFit, TH1D** hPrim, TH1D** hMat, TH1D** hFdwn, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax);
  36
  37 TCanvas* DrawGoF(TH1D** hData, TH1D** hDataFit, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax);
  38
  39 void DrawSec(const TString& inputFile="debug.root", const TString& tag="test", const TString& particle="Proton", Double_t ptmin=0.5, Double_t ptmax=3., Double_t dcaxyMin=-1.5, Double_t dcaxyMax=1.5)
  40 {
  41 //
  42 // draw corrections of secondaries
  43 //
  44         gStyle->SetPadTickX(1);
  45         gStyle->SetPadTickY(1);
  46         gStyle->SetOptStat(0);
  47         gStyle->SetOptTitle(1);
  48
  49         const Int_t kNBin = 50;
  50
  51         TFile* finput = new TFile(inputFile.Data());
  52         if (finput->IsZombie()) exit(1);
  53
  54         // data and templates
  55
  56         TH1D* hData[kNBin];
  57         TH1D* hDataMC[kNBin];
  58         TH1D* hPrim[kNBin];
  59         TH1D* hMat[kNBin];
  60         TH1D* hFdwn[kNBin];
  61
  62         // fitted data
  63
  64         TH1D* hDataFit[kNBin];
  65         TH1D* hPrimFit[kNBin];
  66         TH1D* hMatFit[kNBin];
  67         TH1D* hFdwnFit[kNBin];
  68
  69         TH1D* hPidPt = FindObj<TH1D>(finput, tag, particle + "_PID_Pt");
  70
  71         Int_t lowbin = hPidPt->GetXaxis()->FindFixBin(ptmin);
  72         Int_t hibin  = hPidPt->GetXaxis()->FindFixBin(ptmax);
  73
  74         for(Int_t i=lowbin; i<hibin && i<kNBin; ++i)
  75         {
  76                 // MC
  77                 hData[i]   = FindObj<TH1D>(finput, tag, particle + Form("_Data_DCAxy_%02d",i));
  78                 hDataMC[i] = FindObj<TH1D>(finput, tag, particle + Form("_SimData_DCAxy_%02d",i));
  79                 hPrim[i]   = FindObj<TH1D>(finput, tag, particle + Form("_Prim_DCAxy_%02d",i));
  80                 hMat[i]    = FindObj<TH1D>(finput, tag, particle + Form("_Mat_DCAxy_%02d",i));
  81
  82                 if(particle.Contains("Proton"))
  83                 {
  84                         hFdwn[i] = FindObj<TH1D>(finput, tag, particle + Form("_Fdwn_DCAxy_%02d",i));
  85                 }
  86                 else
  87                 {
  88                         hFdwn[i] = 0;
  89                 }
  90
  91                 // fitted data
  92                 hDataFit[i] = FindObj<TH1D>(finput, tag, particle + Form("_Fit_Data_DCAxy_%02d",i));
  93                 hPrimFit[i] = FindObj<TH1D>(finput, tag, particle + Form("_Fit_Prim_DCAxy_%02d",i));
  94                 hMatFit[i]  = FindObj<TH1D>(finput, tag, particle + Form("_Fit_Mat_DCAxy_%02d",i));
  95
  96                 if(particle.Contains("Proton"))
  97                 {
  98                         hFdwnFit[i] = FindObj<TH1D>(finput, tag, particle + Form("_Fit_Fdwn_DCAxy_%02d",i));
  99                 }
 100                 else
 101                 {
 102                         hFdwnFit[i] = 0;
 103                 }
 104         }
 105
 106         // draw
 107
 108         TCanvas* c0 = 0;
 109
 110         if(hDataMC[lowbin]!=0)
 111         {
 112                 c0 = DrawMC(hDataMC, hPrim, hMat, hFdwn, lowbin, hibin, particle + ".MC", Form("MC simulation for %s",particle.Data()), hPidPt, dcaxyMin, dcaxyMax);
 113                 c0->Update();
 114         }
 115
 116         TCanvas* c1;
 117
 118         if(hDataFit[lowbin]!=0)
 119         {
 120                 c1 = DrawFit(hData, hDataFit, hPrimFit, hMatFit, hFdwnFit, lowbin, hibin, particle + ".Fit", Form("Fitted DCA for %s",particle.Data()), hPidPt, dcaxyMin, dcaxyMax);
 121                 c1->Update();
 122         }
 123
 124         TCanvas* c2;
 125
 126         if(hDataFit[lowbin]!=0)
 127         {
 128                 c2 = DrawGoF(hData, hDataFit, lowbin, hibin, particle + ".GoF", Form("GoF for %s",particle.Data()), hPidPt, dcaxyMin, dcaxyMax);
 129                 c2->Update();
 130         }
 131 }
 132
 133 Int_t GetNumberCols(Int_t lowbin, Int_t hibin)
 134 {
 135 //
 136 // number of rows and cols for plotting the slices
 137 // square matrix (cols = rows)
 138 //
 139         return TMath::CeilNint(TMath::Sqrt(TMath::Abs(hibin - lowbin)));
 140 }
 141
 142 TCanvas* DrawMC(TH1D** hData, TH1D** hPrim, TH1D** hMat, TH1D** hFdwn, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax)
 143 {
 144 //
 145 // Draw MC DCA distributions for each pt bin
 146 //
 147         Int_t ncol = GetNumberCols(lowbin, hibin);
 148
 149         TCanvas* c0 = new TCanvas(name.Data(), title.Data());
 150         c0->Divide(ncol, ncol);
 151
 152         for(Int_t i=lowbin; i < hibin; ++i)
 153         {
 154                 c0->cd(i+1-lowbin);
 155
 156                 gPad->SetLogy();
 157
 158                 hData[i]->SetTitle(Form("%0.2f < p_{T} < %0.2f GeV/c", hPt->GetBinLowEdge(i), hPt->GetBinLowEdge(i)+hPt->GetBinWidth(i)));
 159                 hData[i]->SetMinimum(0.2);
 160
 161                 hData[i]->SetAxisRange(dcaxyMin, dcaxyMax, "X");
 162
 163                 hData[i]->SetMarkerStyle(kFullCircle);
 164                 hData[i]->Draw("E");
 165
 166                 if(hPrim[i] != 0)
 167                 {
 168                         hPrim[i]->SetLineColor(9);
 169                         hPrim[i]->Draw("sameHist");
 170                 }
 171
 172                 if(hMat[i] != 0)
 173                 {
 174                         hMat[i]->SetLineColor(46);
 175                         hMat[i]->Draw("sameHist");
 176                 }
 177
 178                 if(hFdwn != 0)
 179                 if(hFdwn[i] != 0)
 180                 {
 181                         hFdwn[i]->SetLineColor(8);
 182                         hFdwn[i]->Draw("sameHist");
 183                 }
 184         }
 185
 186         return c0;
 187 }
 188
 189 TCanvas* DrawFit(TH1D** hData, TH1D** hDataFit, TH1D** hPrim, TH1D** hMat, TH1D** hFdwn, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax)
 190 {
 191 //
 192 // Draw DCA fitted model for each pt bin in current canvas
 193 //
 194         Int_t ncol = GetNumberCols(lowbin, hibin);
 195
 196         TCanvas* c0 = new TCanvas(name.Data(), title.Data());
 197         c0->Divide(ncol, ncol);
 198
 199         for(Int_t i=lowbin; i < hibin && i < ncol*ncol; ++i)
 200         {
 201                 c0->cd(i+1-lowbin);
 202
 203                 gPad->SetLogy();
 204
 205                 hData[i]->SetTitle(Form("%0.2f < #it{p}_{T} < %0.2f GeV/c", hPt->GetBinLowEdge(i), hPt->GetBinLowEdge(i)+hPt->GetBinWidth(i)));
 206                 hData[i]->SetMinimum(0.2);
 207
 208                 hData[i]->SetAxisRange(dcaxyMin, dcaxyMax, "X");
 209
 210                 hData[i]->SetMarkerStyle(kFullCircle);
 211                 hData[i]->Draw("E");
 212
 213                 hDataFit[i]->SetLineColor(2);
 214                 hDataFit[i]->Draw("sameHist");
 215
 216                 if(hPrim[i] != 0)
 217                 {
 218                         hPrim[i]->SetLineColor(9);
 219                         hPrim[i]->Draw("sameHist");
 220                 }
 221
 222                 if(hMat[i] != 0)
 223                 {
 224                         hMat[i]->SetLineColor(46);
 225                         hMat[i]->Draw("sameHist");
 226                 }
 227
 228                 if(hFdwn != 0)
 229                 if(hFdwn[i] != 0)
 230                 {
 231                         hFdwn[i]->SetLineColor(8);
 232                         hFdwn[i]->Draw("sameHist");
 233                 }
 234         }
 235
 236         return c0;
 237 }
 238
 239 TCanvas* DrawGoF(TH1D** hData, TH1D** hDataFit, Int_t lowbin, Int_t hibin, const TString& name, const TString& title, const TH1D* hPt, Double_t dcaxyMin, Double_t dcaxyMax)
 240 {
 241 //
 242 // Draw a goodness of fit for each pt bin consisting of data/fit
 243 //
 244         Int_t ncol = GetNumberCols(lowbin, hibin);
 245
 246         TCanvas* c0 = new TCanvas(name.Data(), title.Data());
 247         c0->Divide(ncol, ncol);
 248
 249         for(Int_t i=lowbin; i < hibin; ++i)
 250         {
 251                 c0->cd(i+1-lowbin);
 252
 253                 TH1D* hDiv = (TH1D*)hData[i]->Clone("_Div_Data_Fit_");
 254                 hDiv->Sumw2();
 255                 hDiv->Divide(hDataFit[i]);
 256
 257                 hDiv->SetTitle(Form("%0.2f < p_{T} < %0.2f GeV/c", hPt->GetBinLowEdge(i), hPt->GetBinLowEdge(i)+hPt->GetBinWidth(i)));
 258                 hDiv->SetAxisRange(dcaxyMin, dcaxyMax, "X");
 259                 hDiv->SetAxisRange(0,3,"Y");
 260                 hDiv->SetMarkerStyle(kFullCircle);
 261                 hDiv->DrawCopy("E");
 262
 263                 delete hDiv;
 264         }
 265
 266         return c0;
 267 }