PWGHF/vertexingHF/macros/VertexResolutionsFromCmpHF.C

   1 void VertexResolutionsFromCmpHF(Int_t pdgSel=421,Int_t nprongsSel=2)
   2 {
   3   //
   4   // Computes secondary vertex resolutions from the ntuple
   5   // written to CmpHF.root by AliAnalysisTaskSECompareHF
   6   // A.Dainese
   7   //
   8
   9   gStyle->SetOptStat(0);
  10   gStyle->SetTitleXSize(0.05);
  11   gStyle->SetTitleYSize(0.05);
  12
  13   // input
  14   TFile *in = new TFile("CmpHF.root");
  15   TNtuple *nt = (TNtuple*)in->Get("fNtupleCmp");
  16   Float_t pdg,nprongs,ptRec,vxRec,vxTrue,errVx,vyRec,vyTrue,errVy,vzRec,vzTrue,errVz,mRec;
  17   nt->SetBranchAddress("pdg",&pdg);
  18   nt->SetBranchAddress("nprongs",&nprongs);
  19   nt->SetBranchAddress("PtRec",&ptRec);
  20   nt->SetBranchAddress("VxRec",&vxRec);
  21   nt->SetBranchAddress("VxTrue",&vxTrue);
  22   nt->SetBranchAddress("ErrVx",&errVx);
  23   nt->SetBranchAddress("VyRec",&vyRec);
  24   nt->SetBranchAddress("VyTrue",&vyTrue);
  25   nt->SetBranchAddress("ErrVy",&errVy);
  26   nt->SetBranchAddress("VzRec",&vzRec);
  27   nt->SetBranchAddress("VzTrue",&vzTrue);
  28   nt->SetBranchAddress("ErrVz",&errVz);
  29   nt->SetBranchAddress("Mrec",&mRec);
  30
  31   Int_t entries = (Int_t)nt->GetEntries();
  32
  33   // histograms for 15 pt bins
  34   TH1F *hV2x_Reco_ = new TH1F[15];
  35   TH1F *hV2y_Reco_ = new TH1F[15];
  36   TH1F *hV2z_Reco_ = new TH1F[15];
  37
  38   TH1F *hV2pullx_Reco_ = new TH1F[15];
  39   TH1F *hV2pully_Reco_ = new TH1F[15];
  40   TH1F *hV2pullz_Reco_ = new TH1F[15];
  41
  42   TH1F *hChi2_ = new TH1F[15];
  43
  44   TH1F *hDum = new TH1F("hDum","",100,-1000,1000);
  45   TH1F *hDum2 = new TH1F("hDum2","",100,-5,5);
  46   for(Int_t j=0; j<15; j++) {
  47     hV2x_Reco_[j] = *hDum;
  48     hV2y_Reco_[j] = *hDum;
  49     hV2z_Reco_[j] = *hDum;
  50     hV2pullx_Reco_[j] = *hDum2;
  51     hV2pully_Reco_[j] = *hDum2;
  52     hV2pullz_Reco_[j] = *hDum2;
  53   }
  54   delete hDum; hDum=0;
  55   delete hDum2; hDum2=0;
  56
  57
  58
  59   Int_t bin;
  60   // loop on candidates
  61   for(Int_t i=0; i<entries; i++) {
  62     nt->GetEvent(i);
  63     if(i%10000==0) cout<<" Processing entry "<<i<<" of "<<entries<<endl;
  64     if(TMath::Abs(pdg)!=pdgSel || nprongs!=nprongsSel) continue;
  65     bin = GetPtBin(ptRec);
  66     if(bin==-1) continue;
  67     hV2x_Reco_[bin].Fill(10000.*(vxRec-vxTrue));
  68     hV2y_Reco_[bin].Fill(10000.*(vyRec-vyTrue));
  69     hV2z_Reco_[bin].Fill(10000.*(vzRec-vzTrue));
  70     hV2pullx_Reco_[bin].Fill((vxRec-vxTrue)/errVx);
  71     hV2pully_Reco_[bin].Fill((vyRec-vyTrue)/errVy);
  72     hV2pullz_Reco_[bin].Fill((vzRec-vzTrue)/errVz);
  73   }
  74
  75
  76   // output histograms
  77   Int_t nbins = 15;
  78   Float_t xbins[16]={0,0.5,1,1.5,2,2.5,3,3.5,4,5,6,7,8,10,12,14};
  79   TH1F *hResV2x_Reco = new TH1F("hResV2x_Reco","",nbins,xbins);
  80   hResV2x_Reco->SetMaximum(300);
  81   hResV2x_Reco->SetMinimum(0);
  82   hResV2x_Reco->SetXTitle("p_{t} D^{0} [GeV/c]");
  83   hResV2x_Reco->SetYTitle("resolution D^{0}#rightarrow K#pi vertex [#mu m]");
  84   hResV2x_Reco->SetLineWidth(1);
  85   hResV2x_Reco->SetMarkerStyle(20);
  86
  87   TH1F *hResV2y_Reco = (TH1F*)hResV2x_Reco->Clone("hResV2y_Reco");
  88   hResV2y_Reco->SetYTitle("resolution y  D^{0} vertex [#mu m]");
  89
  90   TH1F *hResV2z_Reco = (TH1F*)hResV2x_Reco->Clone("hResV2z_Reco");
  91   hResV2z_Reco->SetYTitle("resolution z  D^{0} vertex [#mu m]");
  92   hResV2z_Reco->SetMaximum(200);
  93
  94
  95   TH1F *hPullV2x_Reco = new TH1F("hPullV2x_Reco","",nbins,xbins);
  96   hPullV2x_Reco->SetMaximum(4);
  97   hPullV2x_Reco->SetMinimum(0);
  98   hPullV2x_Reco->SetXTitle("p_{t} D^{0} [GeV/c]");
  99   hPullV2x_Reco->SetYTitle("pull D^{0}#rightarrow K#pi vertex");
 100   hPullV2x_Reco->SetLineWidth(1);
 101   hPullV2x_Reco->SetMarkerStyle(20);
 102
 103   TH1F *hPullV2y_Reco = (TH1F*)hPullV2x_Reco->Clone("hPullV2y_Reco");
 104   hPullV2y_Reco->SetYTitle("pull y  D^{0} vertex");
 105   hPullV2y_Reco->SetMaximum(4);
 106
 107   TH1F *hPullV2z_Reco = (TH1F*)hPullV2x_Reco->Clone("hPullV2z_Reco");
 108   hPullV2z_Reco->SetYTitle("pull z  D^{0} vertex");
 109   hPullV2z_Reco->SetMaximum(4);
 110
 111
 112   // fit gaussian
 113   TF1 *g = new TF1("g","gaus",-10000.,10000.);
 114   TCanvas *cccx = new TCanvas("cccx","cccx",0,0,800,800);
 115   cccx->Divide(5,3);
 116   TCanvas *cccz = new TCanvas("cccz","cccz",0,0,800,800);
 117   cccz->Divide(5,3);
 118   // fit
 119   for(Int_t j=0;j<15;j++) {
 120     cccz->cd(j+1);
 121     //  Secondary vertex   RECONSTRUCTED
 122
 123     g->SetRange(-3.*hV2y_Reco_[j].GetRMS(),+3.*hV2y_Reco_[j].GetRMS());
 124     if(j>-9) hV2y_Reco_[j].Rebin(4);
 125     hV2y_Reco_[j].Fit("g","R,Q");
 126     hResV2y_Reco->SetBinContent(j+1,hV2y_Reco_[j].GetRMS());
 127     hResV2y_Reco->SetBinContent(j+1,g->GetParameter(2));
 128     hResV2y_Reco->SetBinError(j+1,g->GetParError(2));
 129
 130     g->SetRange(-3.*hV2z_Reco_[j].GetRMS(),+3.*hV2z_Reco_[j].GetRMS());
 131     if(j>-9) hV2z_Reco_[j].Rebin(4);
 132     hV2z_Reco_[j].Fit("g","R,Q");
 133     hResV2z_Reco->SetBinContent(j+1,hV2z_Reco_[j].GetRMS());
 134     hResV2z_Reco->SetBinContent(j+1,g->GetParameter(2));
 135     hResV2z_Reco->SetBinError(j+1,g->GetParError(2));
 136
 137     cccx->cd(j+1);
 138     g->SetRange(-3.*hV2x_Reco_[j].GetRMS(),+3.*hV2x_Reco_[j].GetRMS());
 139     if(j>-9) hV2x_Reco_[j].Rebin(4);
 140     hV2x_Reco_[j].Draw();
 141     hV2x_Reco_[j].Fit("g","R,Q");
 142     hResV2x_Reco->SetBinContent(j+1,hV2x_Reco_[j].GetRMS());
 143     hResV2x_Reco->SetBinContent(j+1,g->GetParameter(2));
 144     hResV2x_Reco->SetBinError(j+1,g->GetParError(2));
 145
 146   }
 147
 148
 149
 150   TCanvas *ccccx = new TCanvas("ccccx","ccccx",0,0,800,800);
 151   ccccx->Divide(5,3);
 152   TCanvas *ccccz = new TCanvas("ccccz","ccccz",0,0,800,800);
 153   ccccz->Divide(5,3);
 154   // fit
 155   for(Int_t j=0;j<15;j++) {
 156     ccccz->cd(j+1);
 157     g->SetRange(-3.*hV2pully_Reco_[j].GetRMS(),+3.*hV2pully_Reco_[j].GetRMS());
 158     if(j>9) hV2pully_Reco_[j].Rebin(4);
 159     hV2pully_Reco_[j].Fit("g","R,Q");
 160     hPullV2y_Reco->SetBinContent(j+1,g->GetParameter(2));
 161     hPullV2y_Reco->SetBinError(j+1,g->GetParError(2));
 162
 163     g->SetRange(-3.*hV2pullz_Reco_[j].GetRMS(),+3.*hV2pullz_Reco_[j].GetRMS());
 164     if(j>9) hV2pullz_Reco_[j].Rebin(4);
 165     hV2pullz_Reco_[j].Fit("g","R,Q");
 166     hPullV2z_Reco->SetBinContent(j+1,g->GetParameter(2));
 167     hPullV2z_Reco->SetBinError(j+1,g->GetParError(2));
 168
 169     ccccx->cd(j+1);
 170     g->SetRange(-3.*hV2pullx_Reco_[j].GetRMS(),+3.*hV2pullx_Reco_[j].GetRMS());
 171     if(j>9) hV2pullx_Reco_[j].Rebin(4);
 172     hV2pullx_Reco_[j].Draw();
 173     hV2pullx_Reco_[j].Fit("g","R,Q");
 174     hPullV2x_Reco->SetBinContent(j+1,g->GetParameter(2));
 175     hPullV2x_Reco->SetBinError(j+1,g->GetParError(2));
 176   }
 177
 178
 179   // Draw Secondary vertex resolution
 180   TCanvas *cV2 = new TCanvas("cV2","cV2",0,0,700,900);
 181   cV2->Divide(1,3);
 182   cV2->cd(1);
 183   hResV2x_Reco->Draw("p,e");
 184   cV2->cd(2);
 185   hResV2y_Reco->Draw("p,e");
 186   cV2->cd(3);
 187   hResV2z_Reco->Draw("p,e");
 188
 189   // Draw Secondary vertex pulls
 190   TCanvas *cV2pull = new TCanvas("cV2pull","cV2pull",0,0,700,900);
 191   cV2pull->Divide(1,3);
 192   cV2pull->cd(1);
 193   hPullV2x_Reco->Draw("p,e");
 194   cV2pull->cd(2);
 195   hPullV2y_Reco->Draw("p,e");
 196   cV2pull->cd(3);
 197   hPullV2z_Reco->Draw("p,e");
 198
 199
 200   // Draw Secondary vertex resolution
 201   TCanvas *cSummary = new TCanvas("cSummary","cSummary",0,0,900,500);
 202   cSummary->Divide(2,1);
 203   cSummary->cd(1);
 204   hResV2x_Reco->SetMarkerStyle(20);
 205   hResV2x_Reco->SetMarkerColor(1);
 206   hResV2x_Reco->Draw("p,e");
 207   hResV2y_Reco->SetMarkerStyle(21);
 208   hResV2y_Reco->SetMarkerColor(2);
 209   hResV2y_Reco->Draw("p,e,same");
 210   hResV2z_Reco->SetMarkerStyle(22);
 211   hResV2z_Reco->SetMarkerColor(4);
 212   hResV2z_Reco->Draw("p,e,same");
 213   TLegend *legg = new TLegend(0.5,0.5,0.9,0.9);
 214   legg->SetFillStyle(0);
 215   legg->SetBorderSize(0);
 216   legg->AddEntry(hResV2x_Reco,"x coordinate","p");
 217   legg->AddEntry(hResV2y_Reco,"y coordinate","p");
 218   legg->AddEntry(hResV2z_Reco,"z coordinate","p");
 219   legg->Draw();
 220   cSummary->cd(2);
 221   hPullV2x_Reco->SetMarkerStyle(20);
 222   hPullV2x_Reco->SetMarkerColor(1);
 223   hPullV2x_Reco->Draw("p,e");
 224   hPullV2y_Reco->SetMarkerStyle(21);
 225   hPullV2y_Reco->SetMarkerColor(2);
 226   hPullV2y_Reco->Draw("p,e,same");
 227   hPullV2z_Reco->SetMarkerStyle(22);
 228   hPullV2z_Reco->SetMarkerColor(4);
 229   hPullV2z_Reco->Draw("p,e,same");
 230   legg->Draw();
 231
 232   return;
 233 }
 234 //---------------------------------------------------------------------------
 235 Int_t GetPtBin(Double_t pt) {
 236
 237   if(pt<0.5) return 0;
 238   if(pt<1) return 1;
 239   if(pt<1.5) return 2;
 240   if(pt<2) return 3;
 241   if(pt<2.5) return 4;
 242   if(pt<3) return 5;
 243   if(pt<3.5) return 6;
 244   if(pt<4) return 7;
 245   if(pt<5) return 8;
 246   if(pt<6) return 9;
 247   if(pt<7) return 10;
 248   if(pt<8) return 11;
 249   if(pt<10) return 12;
 250   if(pt<12) return 13;
 251   if(pt<14) return 14;
 252   return -1;
 253
 254 }
 255
 256
 257
 258
 259
 260
 261
 262
 263
 264
 265
 266