]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG1/global/AliAnalysisTaskVertexESD.cxx
fixed seeting of eta cut
[u/mrichter/AliRoot.git] / PWG1 / global / AliAnalysisTaskVertexESD.cxx
index bde90e0f0f2af5d47073dd9f107da69597187638..2ac17e1e2b9ea7b048347a8b20421ffce99496e7 100644 (file)
@@ -34,6 +34,8 @@
 #include <TH1F.h>
 #include <TH2F.h>  
 #include <TCanvas.h>
+#include <TGraphAsymmErrors.h>
+#include <TFile.h>
 
 #include "AliAnalysisTask.h"
 #include "AliAnalysisManager.h"
@@ -62,29 +64,49 @@ ClassImp(AliAnalysisTaskVertexESD)
 //________________________________________________________________________
 AliAnalysisTaskVertexESD::AliAnalysisTaskVertexESD(const char *name) : 
 AliAnalysisTaskSE(name), 
-fCheckEventType(kTRUE),
-fReadMC(kFALSE),
-fRecoVtxTPC(kFALSE),
-fRecoVtxITSTPC(kFALSE),
-fRecoVtxITSTPCHalfEvent(kFALSE),
-fOnlyITSTPCTracks(kFALSE),
-fOnlyITSSATracks(kFALSE),
-fFillNtuple(kFALSE),
-fFillTreeBeamSpot(kFALSE),
-fESD(0), 
-fOutput(0), 
-fNtupleVertexESD(0),
-fhSPDVertexX(0),
-fhSPDVertexY(0),
-fhSPDVertexZ(0),
-fhTRKVertexX(0),
-fhTRKVertexY(0),
-fhTRKVertexZ(0),
-fhTPCVertexX(0),
-fhTPCVertexY(0),
-fhTPCVertexZ(0),
-fhTrackRefs(0),
-fTreeBeamSpot(0)
+  fCheckEventType(kTRUE),
+  fReadMC(kFALSE),
+  fRecoVtxTPC(kTRUE),
+  fRecoVtxITSTPC(kTRUE),
+  fRecoVtxITSTPCHalfEvent(kFALSE),
+  fOnlyITSTPCTracks(kFALSE),
+  fOnlyITSSATracks(kFALSE),
+  fFillNtuple(kFALSE),
+  fFillTreeBeamSpot(kFALSE),
+  fESD(0), 
+  fOutput(0), 
+  fNtupleVertexESD(0),
+  fhSPDVertexX(0),
+  fhSPDVertexY(0),
+  fhSPDVertexZ(0),
+  fhTRKVertexX(0),
+  fhTRKVertexY(0),
+  fhTRKVertexZ(0),
+  fhTPCVertexX(0),
+  fhTPCVertexY(0),
+  fhTPCVertexZ(0),
+  fhTrackRefs(0),
+  fTreeBeamSpot(0),
+  fhTriggeredTrklets(0),
+  fhSPD3DTrklets(0),
+  fhSPDZTrklets(0),
+  fhTRKTrklets(0),
+  fhTRKcTrklets(0),
+  fhTRKncTrklets(0),
+  fhTPCTrklets(0),
+  fhTPCcTrklets(0),
+  fhTPCncTrklets(0),
+  fhSPD3DZreco(0),
+  fhSPDZZreco(0),
+  fhSPDVertexXPile(0),
+  fhSPDVertexYPile(0),
+  fhSPDVertexZPile(0),
+  fhSPDVertexDiffZPileContr2(0),
+  fhSPDVertexDiffZPileContr3(0),
+  fhSPDVertexDiffZPileContr4(0),
+  fhSPDVertexDiffZPileContr5(0),
+  fhSPDContributorsPile(0),
+  fhSPDDispContributors(0)
 {
   // Constructor
 
@@ -139,6 +161,23 @@ void AliAnalysisTaskVertexESD::UserCreateOutputObjects()
   fOutput->Add(fhTPCVertexY);
   fhTPCVertexZ = new TH1F("fhTPCVertexZ","TPCVertex z; z vertex [cm]; events",200,-20,20);
   fOutput->Add(fhTPCVertexZ);
+  
+  
+  fhSPDVertexXPile = new TH1F("fhSPDVertexXPile","SPDVertexPile x; x vertex [cm]; events",200,-20,20);
+  fOutput->Add(fhSPDVertexXPile);
+  fhSPDVertexYPile = new TH1F("fhSPDVertexYPile","SPDVertexPile y; y vertex [cm]; events",200,-20,20);
+  fOutput->Add(fhSPDVertexYPile);
+  fhSPDVertexZPile = new TH1F("fhSPDVertexZPile","SPDVertexPile z; z vertex [cm]; events",200,-40,40);
+  fOutput->Add(fhSPDVertexZPile);
+  
+  fhSPDVertexDiffZPileContr2 = new TH1F("fhSPDVertexDiffZPileContr2","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
+  fOutput->Add(fhSPDVertexDiffZPileContr2);
+  fhSPDVertexDiffZPileContr3 = new TH1F("fhSPDVertexDiffZPileContr3","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
+  fOutput->Add(fhSPDVertexDiffZPileContr3);
+  fhSPDVertexDiffZPileContr4 = new TH1F("fhSPDVertexDiffZPileContr4","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
+  fOutput->Add(fhSPDVertexDiffZPileContr4);
+  fhSPDVertexDiffZPileContr5 = new TH1F("fhSPDVertexDiffZPileContr5","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
+  fOutput->Add(fhSPDVertexDiffZPileContr5);
 
   fhTrackRefs = new TH2F("fhTrackRefs","Track references; x; y",1000,-4,4,1000,-4,4);
   fOutput->Add(fhTrackRefs);
@@ -158,6 +197,40 @@ void AliAnalysisTaskVertexESD::UserCreateOutputObjects()
   fTreeBeamSpot->Branch("ntrksTRKnc", &ntrksTRKnc, "ntrksTRKnc/s");
   fOutput->Add(fTreeBeamSpot);
 
+  Int_t nbinTrklets=29;
+  Float_t lowTrklets[30]={-0.5,0.5,1.5,2.5,3.5,4.5,5.5,7.5,10.5,25.5,50.5,75.5,100.5,150.5,200.5,300.5,400.5,500.5,600.5,800.5,1000.5,1500.5,2000.5,2500.5,3000.5,4000.5,5000.5,6000.5,8000.5,10000.5};
+  fhTriggeredTrklets = new TH1F("fhTriggeredTrklets","trklets dist for triggered ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTriggeredTrklets);
+  fhSPD3DTrklets = new TH1F("fhSPD3DTrklets","trklets dist for SPD3D ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhSPD3DTrklets);
+  fhSPDZTrklets = new TH1F("fhSPDZTrklets","trklets dist for SPDZ ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhSPDZTrklets);
+  fhTRKTrklets = new TH1F("fhTRKTrklets","trklets dist for TRK ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTRKTrklets);
+  fhTRKcTrklets = new TH1F("fhTRKcTrklets","trklets dist for TRKc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTRKcTrklets);
+  fhTRKncTrklets = new TH1F("fhTRKncTrklets","trklets dist for TRKnc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTRKncTrklets);
+  fhTPCTrklets = new TH1F("fhTPCTrklets","trklets dist for TPC ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTPCTrklets);
+  fhTPCcTrklets = new TH1F("fhTPCcTrklets","trklets dist for TPCc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTPCcTrklets);
+  fhTPCncTrklets = new TH1F("fhTPCncTrklets","trklets dist for TPCnc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
+  fOutput->Add(fhTPCncTrklets);
+  
+  Int_t nbinZreco = 16;
+  Double_t lowZreco[17]={-15.0,-10.0,-7.0,-5,-4,-3,-2,-1,0,1,2,3,4,5,7,10,15};
+  fhSPD3DZreco = new TH1F("fhSPD3DZreco","Zreco dist for SPD3D ev.; Zreco [cm]; entries",nbinZreco,lowZreco);
+  fOutput->Add(fhSPD3DZreco);
+  fhSPDZZreco = new TH1F("fhSPDZZreco","Zreco dist for SPDZ ev.; Zreco [cm]; entries",nbinZreco,lowZreco);
+  fOutput->Add(fhSPDZZreco);
+  
+  fhSPDContributorsPile = new TH1F("fhSPDContributorsPile","ncontributors pile up vertex; ncontributors; entries",200,-0.5,199.5);
+  fOutput->Add(fhSPDContributorsPile);
+  
+  fhSPDDispContributors = new TH2F("fhSPDDispContributors","ncontributors main-pile; ncontributors main; ncontributors pile",200,-0.5,199.5,200,-0.5,199.5);
+  fOutput->Add(fhSPDDispContributors);
+  
   PostData(1, fOutput);
 
   return;
@@ -175,7 +248,11 @@ void AliAnalysisTaskVertexESD::UserExec(Option_t *)
   }
   AliESDEvent* esdE = (AliESDEvent*) InputEvent();
   
-  if(fCheckEventType && (esdE->GetEventType())!=7) return; 
+  // Select PHYSICS events (type=7, for data) and MC events (type=0)
+  // fCheckEventType is kFALSE if fReadMC is kTRUE, hence check is skipped
+  if(fCheckEventType) {
+    if(esdE->GetEventType()!=7 && esdE->GetEventType()!=0) return; 
+  }
 
 
   TArrayF mcVertex(3);
@@ -343,6 +420,54 @@ void AliAnalysisTaskVertexESD::UserExec(Option_t *)
     }
   } 
 
+
+  
+   /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+  // filling Vertex reco efficiency plots
+  
+  if(eventTriggered ? 1. : 0.){
+    fhTriggeredTrklets->Fill(ntrklets);
+    
+    if(spdv){
+      if(spdv->GetNContributors()>0.5){
+       TString spdtitle = spdv->GetTitle();
+       if(spdtitle.Contains("vertexer: 3D") ? 1. : 0.){
+         fhSPD3DTrklets->Fill(ntrklets);
+         fhSPD3DZreco->Fill(spdv->GetZv());
+       }
+       else{
+         fhSPDZTrklets->Fill(ntrklets);
+         fhSPDZZreco->Fill(spdv->GetZv());
+       }
+      }
+    }
+    
+    if(trkv){
+      if(trkv->GetNContributors()>0.5)fhTRKTrklets->Fill(ntrklets);
+    }
+    if(fRecoVtxITSTPC) {
+      AliESDVertex *trkvc = ReconstructPrimaryVertexITSTPC(kTRUE);
+      if(trkvc->GetNContributors()>0.5)fhTRKcTrklets->Fill(ntrklets);
+      delete trkvc; trkvc=0;
+      AliESDVertex *trkvnc = ReconstructPrimaryVertexITSTPC(kFALSE);
+      if(trkvnc->GetNContributors()>0.5)fhTRKncTrklets->Fill(ntrklets);
+      delete trkvnc; trkvnc=0;
+    }
+    
+    if(tpcv){
+      if(tpcv->GetNContributors()>0.5)fhTPCTrklets->Fill(ntrklets);
+    }
+    if(fRecoVtxTPC) {
+       AliESDVertex *tpcvc = ReconstructPrimaryVertexTPC(kTRUE);
+       if(tpcvc->GetNContributors()>0.5)fhTPCcTrklets->Fill(ntrklets);
+       delete tpcvc; tpcvc=0;
+       AliESDVertex *tpcvnc = ReconstructPrimaryVertexTPC(kFALSE);
+       if(tpcvnc->GetNContributors()>0.5)fhTPCncTrklets->Fill(ntrklets);
+       delete tpcvnc; tpcvnc=0;
+      }
+  }
+  /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+  
   Float_t xpile=-999.;
   Float_t ypile=-999.;
   Float_t zpile=-999.;
@@ -351,7 +476,22 @@ void AliAnalysisTaskVertexESD::UserExec(Option_t *)
   Float_t ezpile=-999.;
   Int_t ntrkspile=-1;
 
-  if(esdE->GetNumberOfPileupVerticesSPD()>0 && spdvp){
+  if(esdE->GetNumberOfPileupVerticesSPD()>0 && spdvp && spdv){
+  
+    if(spdvp->GetNContributors()>0) {
+     
+         fhSPDVertexXPile->Fill(spdvp->GetXv());
+         fhSPDVertexYPile->Fill(spdvp->GetYv());  
+      fhSPDContributorsPile->Fill(spdvp->GetNContributors());
+      fhSPDDispContributors->Fill(spdv->GetNContributors(),spdvp->GetNContributors());
+      fhSPDVertexZPile->Fill(spdvp->GetZv());
+      if(spdvp->GetNContributors()>=2) {fhSPDVertexDiffZPileContr2->Fill(spdv->GetZv()-spdvp->GetZv());}
+      if(spdvp->GetNContributors()>=3) {fhSPDVertexDiffZPileContr3->Fill(spdv->GetZv()-spdvp->GetZv());}
+      if(spdvp->GetNContributors()>=4) {fhSPDVertexDiffZPileContr4->Fill(spdv->GetZv()-spdvp->GetZv());}
+      if(spdvp->GetNContributors()>=5) {fhSPDVertexDiffZPileContr5->Fill(spdv->GetZv()-spdvp->GetZv());}
+      
+    }
+   
     xpile=spdvp->GetXv();
     expile=spdvp->GetXRes();
     ypile=spdvp->GetYv();
@@ -559,6 +699,99 @@ void AliAnalysisTaskVertexESD::Terminate(Option_t *)
     return;
   }
   
+  //////////////////////////////////////////////////////
+  /*  
+  TH1F *fhTriggeredTrklets=(TH1F*)fOutput->FindObject("fhTriggeredTrklets");
+  TH1F *fhSPDZTrklets=(TH1F*)fOutput->FindObject("fhSPDZTrklets");
+  TH1F *fhSPD3DTrklets=(TH1F*)fOutput->FindObject("fhSPD3DTrklets");
+  TH1F *fhTRKTrklets=(TH1F*)fOutput->FindObject("fhTRKTrklets");
+  TH1F *fhTRKcTrklets=(TH1F*)fOutput->FindObject("fhTRKcTrklets");
+  TH1F *fhTRKncTrklets=(TH1F*)fOutput->FindObject("fhTRKncTrklets");
+  TH1F *fhSPDZZreco=(TH1F*)fOutput->FindObject("fhSPDZZreco");
+  TH1F *fhSPD3DZreco=(TH1F*)fOutput->FindObject("fhSPD3DZreco");
+  
+  TGraphAsymmErrors *fhSPDZEffTrklets=new TGraphAsymmErrors(fhSPDZTrklets,fhTriggeredTrklets,"w");
+  fhSPDZEffTrklets->SetName("fhSPDZEffTrklets");
+  fhSPDZEffTrklets->SetDrawOption("AP");
+  TGraphAsymmErrors *fhSPD3DEffTrklets=new TGraphAsymmErrors(fhSPD3DTrklets,fhTriggeredTrklets,"w");
+  fhSPD3DEffTrklets->SetName("fhSPD3DEffTrklets");
+  TH1F * fhSPDOverallTrklets=(TH1F*)fhSPDZTrklets->Clone("fhSPDOverallTrklets");
+  fhSPDOverallTrklets->Add(fhSPD3DTrklets);
+  TGraphAsymmErrors *fhSPDOverallEffTrklets=new TGraphAsymmErrors(fhSPDOverallTrklets,fhTriggeredTrklets,"w");
+  fhSPDOverallEffTrklets->SetName("fhSPDOverallEffTrklets");
+  TGraphAsymmErrors *fhTRKEffTrklets=new TGraphAsymmErrors(fhTRKTrklets,fhTriggeredTrklets,"w");
+  fhTRKEffTrklets->SetName("fhTRKEffTrklets");
+  TGraphAsymmErrors *fhTRKcEffTrklets=new TGraphAsymmErrors(fhTRKcTrklets,fhTriggeredTrklets,"w");
+  fhTRKcEffTrklets->SetName("fhTRKcEffTrklets");
+  TGraphAsymmErrors *fhTRKncEffTrklets=new TGraphAsymmErrors(fhTRKncTrklets,fhTriggeredTrklets,"w");
+  fhTRKncEffTrklets->SetName("fhTRKncEffTrklets");
+  TH1F * fhSPDOverallZreco=(TH1F*)fhSPDZZreco->Clone("fhSPDOverallZreco");
+  fhSPDOverallZreco->Add(fhSPD3DZreco);
+  TGraphAsymmErrors *fhSPDEffZreco=new TGraphAsymmErrors(fhSPD3DZreco,fhSPDOverallZreco,"w");
+  fhSPDEffZreco->SetName("fhSPDEffZreco");
+
+  TH1F *fhEff = new TH1F("hEff","hEff",6,0.5,6.5);
+  Int_t count=1;
+  if(fhSPDZTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhSPDZTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPDZTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"SPDZ");
+  
+  count++;
+  if(fhSPD3DTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhSPD3DTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPD3DTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"SPD3D");
+  
+  count++;
+  if(fhSPDOverallTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhSPDOverallTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPDOverallTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"SPD Overall");
+  
+  count++;
+  if(fhTRKTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhTRKTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"TRK");
+  
+  count++;
+  if(fhTRKcTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhTRKcTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKcTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"TRKc");
+  
+  count++;
+  if(fhTRKncTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
+    fhEff->Fill(count,fhTRKncTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
+    fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKncTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
+  }
+  fhEff->GetXaxis()->SetBinLabel(count,"TRKnc");
+  
+  count++;
+  fhEff->Print("all");
+  
+  TFile* fileEff = new TFile("VtxEff.root","recreate");
+  fhSPDZEffTrklets->Write();
+  fhSPD3DEffTrklets->Write();
+  fhSPDOverallEffTrklets->Write();
+  fhTRKEffTrklets->Write();
+  fhTRKcEffTrklets->Write();
+  fhTRKncEffTrklets->Write();
+  fhSPDEffZreco->Write();
+  fhEff->Write();
+  fileEff->Close();
+  delete fileEff;
+  
+  /////////////////////////////////////////
+  */
+
+
   if (!fNtupleVertexESD){
     Printf("ERROR: fNtuple not available");
     return;
@@ -575,7 +808,11 @@ AliESDVertex* AliAnalysisTaskVertexESD::ReconstructPrimaryVertexTPC(Bool_t const
   // On the fly reco of TPC vertex from ESD
   AliESDEvent* evt = (AliESDEvent*) fInputEvent;
   AliVertexerTracks vertexer(evt->GetMagneticField());
-  vertexer.SetTPCMode(); // defaults
+  if(evt->GetNumberOfTracks()<500) {
+    vertexer.SetTPCMode(); // defaults
+  } else { 
+    vertexer.SetTPCMode(0.1,1.0,5.0,10,1,3.,0.1,1.5,3.,30.,1,1);// PbPb
+  } 
   Float_t diamondcovxy[3]; evt->GetDiamondCovXY(diamondcovxy);
   Double_t pos[3]={evt->GetDiamondX(),evt->GetDiamondY(),0}; 
   Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
@@ -596,8 +833,12 @@ AliESDVertex* AliAnalysisTaskVertexESD::ReconstructPrimaryVertexITSTPC(Bool_t co
 
   AliESDEvent* evt = (AliESDEvent*) fInputEvent;
   AliVertexerTracks vertexer(evt->GetMagneticField());
-  vertexer.SetITSMode(); // defaults
-  vertexer.SetMinClusters(4); // default is 5
+  if(evt->GetNumberOfTracks()<500) {
+    vertexer.SetITSMode(); // defaults
+    vertexer.SetMinClusters(4); // default is 5
+  } else { 
+    vertexer.SetITSMode(0.1,0.1,0.5,5,1,3.,100.,1000.,3.,30.,1,1);// PbPb
+  } 
   //vertexer.SetITSpureSA();
   Float_t diamondcovxy[3]; evt->GetDiamondCovXY(diamondcovxy);
   Double_t pos[3]={evt->GetDiamondX(),evt->GetDiamondY(),0};