+++ /dev/null
-//------------------------------------------------------------------------------
-// Implementation of AliComparisonDCA class. It keeps information from
-// comparison of reconstructed and MC particle tracks. In addtion,
-// it keeps selection cuts used during comparison. The comparison
-// information is stored in the ROOT histograms. Analysis of these
-// histograms can be done by using Analyse() class function. The result of
-// the analysis (histograms/graphs) are stored in the folder
-// which is a data member of AliComparisonDCA.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-/*
-
- // after running comparison task, read the file, and get component
- gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
- LoadMyLibs();
- TFile f("Output.root");
- //AliComparisonDCA * compObj = (AliComparisonDCA*)f.Get("AliComparisonDCA");
- AliComparisonDCA * compObj = (AliComparisonDCA*)cOutput->FindObject("AliComparisonDCA");
-
- // Analyse comparison data
- compObj->Analyse();
-
- // the output histograms/graphs will be stored in the folder "folderDCA"
- compObj->GetAnalysisFolder()->ls("*");
-
- // user can save whole comparison object (or only folder with anlysed histograms)
- // in the seperate output file (e.g.)
- TFile fout("Analysed_DCA.root","recreate");
- compObj->Write(); // compObj->GetAnalysisFolder()->Write();
- fout.Close();
-
-*/
-
-#include <TAxis.h>
-#include <TCanvas.h>
-#include <TGraph.h>
-#include <TGraph2D.h>
-#include <TH1.h>
-
-#include "AliComparisonDCA.h"
-#include "AliESDEvent.h"
-#include "AliESDRecInfo.h"
-#include "AliESDVertex.h"
-#include "AliLog.h"
-#include "AliMCInfo.h"
-#include "AliMCInfoCuts.h"
-#include "AliMathBase.h"
-#include "AliRecInfoCuts.h"
-#include "AliTracker.h"
-
-using namespace std;
-
-ClassImp(AliComparisonDCA)
-
-//_____________________________________________________________________________
-AliComparisonDCA::AliComparisonDCA():
- AliComparisonObject("AliComparisonDCA"),
-
- // DCA histograms
- fDCAHisto(0),
- /*
- fD0TanSPtTPCITS(0),
- fD1TanSPtTPCITS(0),
- fD0TanSPt(0),
- fD1TanSPt(0),
- fD0TanSPtTPC(0),
- fD1TanSPtTPC(0),
- */
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // default constructor
-}
-
-//_____________________________________________________________________________
-AliComparisonDCA::AliComparisonDCA(Char_t* name="AliComparisonDCA", Char_t* title="AliComparisonDCA",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
- AliComparisonObject(name,title),
-
- // DCA histograms
- fDCAHisto(0),
- /*
- fD0TanSPtTPCITS(0),
- fD1TanSPtTPCITS(0),
- fD0TanSPt(0),
- fD1TanSPt(0),
- fD0TanSPtTPC(0),
- fD1TanSPtTPC(0),
- */
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // named constructor
-
- SetAnalysisMode(analysisMode);
- SetHptGenerator(hptGenerator);
- Init();
-}
-
-
-//_____________________________________________________________________________
-AliComparisonDCA::~AliComparisonDCA()
-{
- // destructor
- if(fDCAHisto) delete fDCAHisto; fDCAHisto=0;
- /*
- if(fD0TanSPtTPCITS) delete fD0TanSPtTPCITS; fD0TanSPtTPCITS=0;
- if(fD1TanSPtTPCITS) delete fD1TanSPtTPCITS; fD1TanSPtTPCITS=0;
- if(fD0TanSPt) delete fD0TanSPt; fD0TanSPt=0;
- if(fD1TanSPt) delete fD1TanSPt; fD1TanSPt=0;
- if(fD0TanSPtTPC) delete fD0TanSPtTPC; fD0TanSPtTPC=0;
- if(fD1TanSPtTPC) delete fD1TanSPtTPC; fD1TanSPtTPC=0;
- */
- if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
-
-}
-
-//_____________________________________________________________________________
-void AliComparisonDCA::Init()
-{
- // DCA histograms
-
- Int_t nPBins = 31;
- Double_t binsP[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
- Double_t pMin = 0., pMax = 10.;
-
- if(IsHptGenerator() == kTRUE) {
- nPBins = 100;
- pMin = 0.; pMax = 100.;
- }
-
- //dca_r, dca_z, eta, pt
- Int_t binsQA[4] = {100,100,20,nPBins};
- Double_t xminQA[4] = {-10.,-10.,-1., pMin};
- Double_t xmaxQA[4] = {10.,10.,1., pMax};
-
- fDCAHisto = new THnSparseF("fDCAHisto","dca_r:dca_z:eta:pt",4,binsQA,xminQA,xmaxQA);
- if(!IsHptGenerator()) fDCAHisto->SetBinEdges(3,binsP);
-
- fDCAHisto->GetAxis(0)->SetTitle("dca_r (cm)");
- fDCAHisto->GetAxis(1)->SetTitle("dca_z (cm)");
- fDCAHisto->GetAxis(2)->SetTitle("eta");
- fDCAHisto->GetAxis(3)->SetTitle("pt (GeV/c)");
- fDCAHisto->Sumw2();
-
- /*
- fD0TanSPtTPCITS = new TH3F("DCAyTanSPtTPCITS","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
- fD0TanSPtTPCITS->SetXTitle("tan(#theta)");
- fD0TanSPtTPCITS->SetYTitle("#sqrt{p_{t}(GeV)}");
- fD0TanSPtTPCITS->SetZTitle("DCA_{xy}");
-
- fD1TanSPtTPCITS = new TH3F("DCAzTanSPtTPCITS","DCAzTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
- fD1TanSPtTPCITS->SetXTitle("tan(#theta)");
- fD1TanSPtTPCITS->SetYTitle("#sqrt(p_{t}(GeV))");
- fD1TanSPtTPCITS->SetZTitle("DCA_{z}");
-
- fD0TanSPt = new TH3F("DCAyTanSPt","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
- fD0TanSPt->SetXTitle("tan(#theta)");
- fD0TanSPt->SetYTitle("#sqrt{p_{t}(GeV)}");
- fD0TanSPt->SetZTitle("DCA_{xy}");
-
- fD1TanSPt = new TH3F("DCAzTanSPt","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
- fD1TanSPt->SetXTitle("tan(#theta)");
- fD1TanSPt->SetYTitle("#sqrt{p_{t}(GeV)}");
- fD1TanSPt->SetZTitle("DCA_{z}");
-
- fD0TanSPtTPC = new TH3F("DCAyTanSPtTPC","DCAyTanSPt",40,-2,2, 10,0.3,3, 100,-1,1);
- fD0TanSPtTPC->SetXTitle("tan(#theta)");
- fD0TanSPtTPC->SetYTitle("#sqrt{p_{t}(GeV)}");
- fD0TanSPtTPC->SetZTitle("DCA_{xy}");
-
- fD1TanSPtTPC = new TH3F("DCAzTanSPtTPC","DCAzTanSPt",40,-2,2, 10,0.3,3, 100, -1,1);
- fD1TanSPtTPC->SetXTitle("tan(#theta)");
- fD1TanSPtTPC->SetYTitle("#sqrt{p_{t}(GeV)}");
- fD1TanSPtTPC->SetZTitle("DCA_{z}");
- */
-
- // init cuts
- if(!fCutsMC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
- if(!fCutsRC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
-
- // init folder
- fAnalysisFolder = CreateFolder("folderDCA","Analysis DCA Folder");
-}
-
-//_____________________________________________________________________________
-void AliComparisonDCA::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
-{
- // Fill DCA comparison information
- AliExternalTrackParam *track = 0;
- //Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]
- Double_t kMaxD = 123456.0; // max distance
-
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
-
- Float_t mcpt = infoMC->GetParticle().Pt();
- Float_t mceta = infoMC->GetParticle().Eta();
- //Float_t spt = TMath::Sqrt(mcpt);
-
- // distance to Prim. vertex
- const Double_t* dv = infoMC->GetVDist();
-
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // Check selection cuts
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
- if (!isPrim) return;
- if (infoRC->GetStatus(1)!=3) return;
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
- //if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;
-
- // calculate and set prim. vertex
- AliESDVertex vertexMC;
- vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );
- vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );
- vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );
-
- // calculate track parameters at vertex
- if (infoRC->GetESDtrack()->GetTPCInnerParam())
- {
- if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )
- {
-
- //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);
- Double_t field[3]; track->GetBxByBz(field);
- Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);
-
- if(bDCAStatus) {
- Double_t vDCAHisto[4]={dca[0],dca[1],mceta,mcpt};
- fDCAHisto->Fill(vDCAHisto);
- }
- delete track;
- }
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonDCA::ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
-{
- // Fill DCA comparison information
- Int_t clusterITS[200];
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
-
- Float_t mcpt = infoMC->GetParticle().Pt();
- Float_t mceta = infoMC->GetParticle().Eta();
- //Float_t spt = TMath::Sqrt(mcpt);
-
- // distance to Prim. vertex
- const Double_t* dv = infoMC->GetVDist();
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // Check selection cuts
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
- if (!isPrim) return;
- if (infoRC->GetStatus(1)!=3) return;
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
- //if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;
-
- infoRC->GetESDtrack()->GetImpactParameters(dca,cov);
-
- // ITS + TPC
- if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)>fCutsRC->GetMinNClustersITS())
- {
- Double_t vDCAHisto[4]={dca[0],dca[1],mceta,mcpt};
- fDCAHisto->Fill(vDCAHisto);
- }
-}
-
-void AliComparisonDCA::ProcessConstrained(AliMCInfo* const /*infoMC*/, AliESDRecInfo * const /*infoRC*/)
-{
- // Fill DCA comparison information
-
- AliDebug(AliLog::kWarning, "Warning: Not implemented");
-}
-
-//_____________________________________________________________________________
-Long64_t AliComparisonDCA::Merge(TCollection* const list)
-{
- // Merge list of objects (needed by PROOF)
-
- if (!list)
- return 0;
-
- if (list->IsEmpty())
- return 1;
-
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- // collection of generated histograms
- Int_t count=0;
- while((obj = iter->Next()) != 0)
- {
- AliComparisonDCA* entry = dynamic_cast<AliComparisonDCA*>(obj);
- if (entry == 0) continue;
-
- fDCAHisto->Add(entry->fDCAHisto);
- /*
- fD0TanSPtTPCITS->Add(entry->fD0TanSPtTPCITS);
- fD1TanSPtTPCITS->Add(entry->fD1TanSPtTPCITS);
- fD0TanSPt->Add(entry->fD0TanSPt);
- fD1TanSPt->Add(entry->fD1TanSPt);
- fD0TanSPtTPC->Add(entry->fD0TanSPtTPC);
- fD1TanSPtTPC->Add(entry->fD1TanSPtTPC);
- */
-
- count++;
- }
-
-return count;
-}
-
-//_____________________________________________________________________________
-void AliComparisonDCA::Exec(AliMCInfo* const infoMC, AliESDRecInfo * const infoRC)
-{
- // Process comparison information
- if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);
- else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);
- else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);
- else {
- printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
- return;
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonDCA::Analyse()
-{
- //
- // Analyse comparison information and store output histograms
- // in the analysis folder "folderDCA"
- //
-
- TH1::AddDirectory(kFALSE);
- TObjArray *aFolderObj = new TObjArray;
-
- /*
- TGraph * gr[8]= { 0,0,0,0,0,0,0,0 };
- TGraph2D *gr2[8]= { 0,0,0,0,0,0,0,0};
- AliComparisonDCA * comp=this;
-
- // write results in the folder
- // Canvas to draw analysed histograms
- TCanvas * c = new TCanvas("canDCA","DCA resolution");
- c->Divide(2,4);
- //
- // DCA resolution
- //
- c->cd(1);
- gr[0] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPC,2,5);
- gr[0]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[0]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
- gr[0]->SetName("DCAXYResolTanTPC");
- gr[0]->SetTitle("resol. DCA_xy (TPC only)");
- gr[0]->Draw("Al*");
-
- aFolderObj->Add(gr[0]);
-
- c->cd(2);
- gr[1] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPC,2,5);
- gr[1]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[1]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
- gr[1]->SetName("DCAZResolTanTPC");
- gr[1]->SetTitle("resol. DCA_z (TPC only)");
- gr[1]->Draw("Al*");
-
- aFolderObj->Add(gr[1]);
-
- c->cd(3);
- gr[2] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPCITS,2,5);
- gr[2]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[2]->GetYaxis()->SetTitle("#sigmaDCA_xy (cm)");
- gr[2]->SetName("DCAXYResolTanTPCITS");
- gr[2]->SetTitle("resol. DCA_xy (TPC+ITS)");
- gr[2]->Draw("Al*");
-
- aFolderObj->Add(gr[2]);
-
- c->cd(4);
- gr[3] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPCITS,2,5);
- gr[3]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[3]->GetYaxis()->SetTitle("#sigmaDCA_z (cm)");
- gr[3]->SetName("DCAZResolTanTPCITS");
- gr[3]->SetTitle("resol. DCA_z (TPC+ITS)");
- gr[3]->Draw("Al*");
-
- aFolderObj->Add(gr[3]);
-
- //
- // DCA mean value
- //
- c->cd(5);
- gr[4] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPC,2,4);
- gr[4]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[4]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
- gr[4]->SetName("DCAXYMeanTanTPC");
- gr[4]->SetTitle("mean DCA_xy (TPC only)");
- gr[4]->Draw("Al*");
-
- aFolderObj->Add(gr[4]);
-
- c->cd(6);
- gr[5] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPC,2,4);
- gr[5]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[5]->GetYaxis()->SetTitle("mean DCA_z (cm)");
- gr[5]->SetName("DCAZMeanTanTPC");
- gr[5]->SetTitle("mean DCA_z (TPC only)");
- gr[5]->Draw("Al*");
-
- aFolderObj->Add(gr[5]);
-
- c->cd(7);
- gr[6] = AliMathBase::MakeStat1D(comp->fD0TanSPtTPCITS,2,4);
- gr[6]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[6]->GetYaxis()->SetTitle("mean DCA_xy (cm)");
- gr[6]->SetName("DCAXYMeanTanTPCITS");
- gr[6]->SetTitle("mean DCA_xy (TPC+ITS)");
- gr[6]->Draw("Al*");
-
- aFolderObj->Add(gr[6]);
-
- c->cd(8);
- gr[7] = AliMathBase::MakeStat1D(comp->fD1TanSPtTPCITS,2,4);
- gr[7]->GetXaxis()->SetTitle("Tan(#theta)");
- gr[7]->GetYaxis()->SetTitle("mean DCA_z (cm)");
- gr[7]->SetName("DCAZMeanTanTPCITS");
- gr[7]->SetTitle("mean DCA_z (TPC+ITS)");
- gr[7]->Draw("Al*");
-
- aFolderObj->Add(gr[7]);
-
- // 2D DCA resolution
- TCanvas * c1 = new TCanvas("canDCA1","2D DCA resolution");
- c1->Divide(2,4);
-
- // TPC only
- c1->cd(1);
- gr2[0] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPC,4,2,5);
- gr2[0]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[0]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[0]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
- gr2[0]->SetName("DCAXYResolSPTTanTPC");
- gr2[0]->SetTitle("#sigma DCA_xy (TPC only)");
- gr2[0]->GetHistogram()->Draw("colz");
-
- gr2[0]->GetHistogram()->SetName("DCAXYResolSPTTanTPC");
- aFolderObj->Add(gr2[0]->GetHistogram());
-
- c1->cd(2);
- gr2[1] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPC,4,2,5);
- gr2[1]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[1]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[1]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
- gr2[1]->SetName("DCAZResolSPTTanTPC");
- gr2[1]->SetTitle("#sigma DCA_z (TPC only)");
- gr2[1]->GetHistogram()->Draw("colz");
-
- gr2[1]->GetHistogram()->SetName("DCAZResolSPTTanTPC");
- aFolderObj->Add(gr2[1]->GetHistogram());
-
- // TPC+ITS
- c1->cd(3);
- gr2[2] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPCITS,4,2,5);
- gr2[2]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[2]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[2]->GetZaxis()->SetTitle("#sigmaDCA_xy (cm)");
- gr2[2]->SetName("DCAXYResolSPTTanTPCITS");
- gr2[2]->SetTitle("#sigma DCA_xy (TPC+ITS)");
- gr2[2]->GetHistogram()->Draw("colz");
-
- gr2[2]->GetHistogram()->SetName("DCAXYResolSPTTanTPCITS");
- aFolderObj->Add(gr2[2]->GetHistogram());
-
- c1->cd(4);
- gr2[3] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPCITS,4,2,5);
- gr2[3]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[3]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[3]->GetZaxis()->SetTitle("#sigmaDCA_z (cm)");
- gr2[3]->SetName("DCAZResolSPTTanTPCITS");
- gr2[3]->SetTitle("#sigma DCA_z (TPC+ITS)");
- gr2[3]->GetHistogram()->Draw("colz");
-
- gr2[3]->GetHistogram()->SetName("DCAZResolSPTTanTPCITS");
- aFolderObj->Add(gr2[3]->GetHistogram());
-
- // 2D DCA mean value
- c1->cd(5);
- gr2[4] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPC,4,2,4);
- gr2[4]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[4]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[4]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
- gr2[4]->SetName("DCAXYMeanSPTTanTPC");
- gr2[4]->SetTitle("mean DCA_xy (TPC only)");
- gr2[4]->GetHistogram()->Draw("colz");
-
- gr2[4]->GetHistogram()->SetName("DCAXYMeanSPTTanTPC");
- aFolderObj->Add(gr2[4]->GetHistogram());
-
- c1->cd(6);
- gr2[5] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPC,4,2,4);
- gr2[5]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[5]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[5]->GetZaxis()->SetTitle("mean DCA_z (cm)");
- gr2[5]->SetName("DCAZMeanSPTTanTPC");
- gr2[5]->SetTitle("mean DCA_z (TPC only)");
- gr2[5]->GetHistogram()->Draw("colz");
-
- gr2[5]->GetHistogram()->SetName("DCAZMeanSPTTanTPC");
- aFolderObj->Add(gr2[5]->GetHistogram());
-
- c1->cd(7);
- gr2[6] = AliMathBase::MakeStat2D(comp->fD0TanSPtTPCITS,4,2,4);
- gr2[6]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[6]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[6]->GetZaxis()->SetTitle("mean DCA_xy (cm)");
- gr2[6]->SetName("DCAXYMeanSPTTanTPCITS");
- gr2[6]->SetTitle("mean DCA_xy (TPC+ITS)");
- gr2[6]->GetHistogram()->Draw("colz");
-
- gr2[6]->GetHistogram()->SetName("DCAXYMeanSPTTanTPCITS");
- aFolderObj->Add(gr2[6]->GetHistogram());
-
- c1->cd(8);
- gr2[7] = AliMathBase::MakeStat2D(comp->fD1TanSPtTPCITS,4,2,4);
- gr2[7]->GetXaxis()->SetTitle("Tan(#theta)");
- gr2[7]->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr2[7]->GetZaxis()->SetTitle("mean DCA_z (cm)");
- gr2[7]->SetName("DCAZMeanSPTTanTPCITS");
- gr2[7]->SetTitle("mean DCA_z (TPC+ITS)");
- gr2[7]->GetHistogram()->Draw("colz");
-
- gr2[7]->GetHistogram()->SetName("DCAZMeanSPTTanTPCITS");
- aFolderObj->Add(gr2[7]->GetHistogram());
-
- */
- // export objects to analysis folder
- fAnalysisFolder = ExportToFolder(aFolderObj);
-
- // delete only TObjArray
- if(aFolderObj) delete aFolderObj;
-}
-
-//_____________________________________________________________________________
-TFolder* AliComparisonDCA::ExportToFolder(TObjArray * array)
-{
- // recreate folder avery time and export objects to new one
- //
- AliComparisonDCA * comp=this;
- TFolder *folder = comp->GetAnalysisFolder();
-
- TString name, title;
- TFolder *newFolder = 0;
- Int_t i = 0;
- Int_t size = array->GetSize();
-
- if(folder) {
- // get name and title from old folder
- name = folder->GetName();
- title = folder->GetTitle();
-
- // delete old one
- delete folder;
-
- // create new one
- newFolder = CreateFolder(name.Data(),title.Data());
- newFolder->SetOwner();
-
- // add objects to folder
- while(i < size) {
- newFolder->Add(array->At(i));
- i++;
- }
- }
-
-return newFolder;
-}
-
-
-//_____________________________________________________________________________
-TFolder* AliComparisonDCA::CreateFolder(TString name,TString title) {
-// create folder for analysed histograms
-TFolder *folder = 0;
- folder = new TFolder(name.Data(),title.Data());
-
- return folder;
-}
+++ /dev/null
-#ifndef ALICOMPARISONDCA_H
-#define ALICOMPARISONDCA_H
-
-//------------------------------------------------------------------------------
-// Class to keep information from comparison of
-// reconstructed and MC particle tracks (DCA - Distance of Closest Approach
-// to the vertex).
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-class AliMCInfo;
-class AliESDRecInfo;
-class AliESDEvent;
-class AliRecInfoCuts;
-class AliMCInfoCuts;
-class AliESDVertex;
-class TH3F;
-class TH3;
-class TString;
-class TNamed;
-
-#include "THnSparse.h"
-#include "AliComparisonObject.h"
-
-class AliComparisonDCA : public AliComparisonObject {
-public :
- AliComparisonDCA();
- AliComparisonDCA(Char_t* name, Char_t* title, Int_t analysisMode, Bool_t hptGenerator);
- ~AliComparisonDCA();
-
- // Init data members
- virtual void Init();
-
- // Execute analysis
- virtual void Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC);
-
- // Merge output objects (needed by PROOF)
- virtual Long64_t Merge(TCollection* const list);
-
- // Analyse output histograms
- virtual void Analyse();
-
- // Get analysis folder
- virtual TFolder* GetAnalysisFolder() const {return fAnalysisFolder;}
-
- // Create folder for analysed histograms
- TFolder *CreateFolder(TString folder = "folderDCA",TString title = "Analysed DCA histograms");
-
- // Export objects to folder
- TFolder *ExportToFolder(TObjArray * array=0);
-
- // Process events
- void ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC); // not implemented
-
- // Selection cuts
- void SetAliRecInfoCuts(AliRecInfoCuts* const cuts=0) {fCutsRC = cuts;}
- void SetAliMCInfoCuts(AliMCInfoCuts* const cuts=0) {fCutsMC = cuts;}
-
- AliRecInfoCuts* GetAliRecInfoCuts() const {return fCutsRC;}
- AliMCInfoCuts* GetAliMCInfoCuts() const {return fCutsMC;}
-
- // getters
- /*
- TH3F *GetD0TanSPtTPCITS() const {return fD0TanSPtTPCITS;}
- TH3F *GetD1TanSPtTPCITS() const {return fD1TanSPtTPCITS;}
- TH3F *GetD0TanSPt() const {return fD0TanSPt;}
- TH3F *GetD1TanSPt() const {return fD1TanSPt;}
- TH3F *GetD0TanSPtTPC() const {return fD0TanSPtTPC;}
- TH3F *GetD1TanSPtTPC() const {return fD1TanSPtTPC;}
- */
-
- // DCA
- THnSparse* GetDCAHisto() const {return fDCAHisto;}
-
-private:
-
- // DCA histograms
- THnSparseF *fDCAHisto; //-> dca_r:dca_z:eta:pt
-
- /*
- TH3F *fD0TanSPtTPCITS; //-> distance to vertex y (TPC+ITS clusters)
- TH3F *fD1TanSPtTPCITS; //-> distance to vertex z (TPC+ITS clusters)
- TH3F *fD0TanSPt; //-> distance to vertex y
- TH3F *fD1TanSPt; //-> distance to vertex z
- TH3F *fD0TanSPtTPC; //-> distance to vertex y (only TPC track parameters)
- TH3F *fD1TanSPtTPC; //-> distance to vertex z (only TPC track parameters)
- */
-
- // Global cuts objects
- AliRecInfoCuts* fCutsRC; // selection cuts for reconstructed tracks
- AliMCInfoCuts* fCutsMC; // selection cuts for MC tracks
-
- // analysis folder
- TFolder *fAnalysisFolder; // folder for analysed histograms
-
- AliComparisonDCA(const AliComparisonDCA&); // not implemented
- AliComparisonDCA& operator=(const AliComparisonDCA&); // not implemented
-
- ClassDef(AliComparisonDCA,1);
-};
-
-#endif
+++ /dev/null
-//------------------------------------------------------------------------------
-// Implementation of AliComparisonDEdx class. It keeps information from
-// comparison of reconstructed and MC particle tracks. In addtion,
-// it keeps selection cuts used during comparison. The comparison
-// information is stored in the ROOT histograms. Analysis of these
-// histograms can be done by using Analyse() class function. The result of
-// the analysis (histograms/graphs) are stored in the folder which is
-// a data of AliComparisonDEdx.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-/*
-
- // after running comparison task, read the file, and get component
- gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
- LoadMyLibs();
- TFile f("Output.root");
- //AliComparisonDEdx * compObj = (AliComparisonDEdx*)f.Get("AliComparisonDEdx");
- AliComparisonDEdx * compObj = (AliComparisonDEdx*)cOutput->FindObject("AliComparisonDEdx");
-
- // Analyse comparison data
- compObj->Analyse();
-
- // the output histograms/graphs will be stored in the folder "folderDEdx"
- compObj->GetAnalysisFolder()->ls("*");
-
- // user can save whole comparison object (or only folder with anlysed histograms)
- // in the seperate output file (e.g.)
- TFile fout("Analysed_DEdx.root"."recreate");
- compObj->Write(); // compObj->GetAnalysisFolder()->Write();
- fout.Close();
-
-*/
-
-#include <TAxis.h>
-#include <TCanvas.h>
-#include <TH1.h>
-
-#include "AliComparisonDEdx.h"
-#include "AliESDEvent.h"
-#include "AliESDRecInfo.h"
-#include "AliLog.h"
-#include "AliMCInfo.h"
-#include "AliMCInfoCuts.h"
-#include "AliMathBase.h"
-#include "AliRecInfoCuts.h"
-#include "AliTreeDraw.h"
-
-using namespace std;
-
-ClassImp(AliComparisonDEdx)
-
-//_____________________________________________________________________________
-AliComparisonDEdx::AliComparisonDEdx():
-// TNamed("AliComparisonDEdx","AliComparisonDEdx"),
- AliComparisonObject("AliComparisonDEdx"),
-
- // dEdx
- fDeDxHisto(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // default constructor
-}
-
-//_____________________________________________________________________________
-AliComparisonDEdx::AliComparisonDEdx(Char_t* name="AliComparisonDEdx", Char_t* title="AliComparisonDEdx",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):
- AliComparisonObject(name,title),
-
- // dEdx
- fDeDxHisto(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // named constructor
-
- SetAnalysisMode(analysisMode);
- SetHptGenerator(hptGenerator);
- Init();
-}
-
-
-//_____________________________________________________________________________
-AliComparisonDEdx::~AliComparisonDEdx()
-{
- // destructor
- if(fDeDxHisto) delete fDeDxHisto; fDeDxHisto=0;
-
- if(fCutsRC) delete fCutsRC; fCutsRC=0;
- if(fCutsMC) delete fCutsMC; fCutsMC=0;
-
- if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::Init()
-{
- // Init histograms
-
- // TPC dEdx
- Int_t nPBins = 31;
- Double_t binsP[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
- Double_t pMin = 0., pMax = 10.;
-
- if(IsHptGenerator() == kTRUE) {
- nPBins = 100;
- pMin = 0.; pMax = 100.;
- }
-
- //signal:alpha:y:z:snp:tgl:ncls:pid:p
- Int_t binsQA[9] = {600,50, 50, 50, 50, 50, 80, 5, nPBins};
- Double_t xminQA[9] = {0, -4,-20,-250, -1, -2, 0, 0., pMin};
- Double_t xmaxQA[9] = {300, 4, 20, 250, 1, 2, 160, 5., pMax};
-
- fDeDxHisto = new THnSparseF("fDeDxHisto","signal:alpha:y:z:snp:tgl:ncls:pid:momentum",9,binsQA,xminQA,xmaxQA);
- if(!IsHptGenerator()) fDeDxHisto->SetBinEdges(8,binsP);
-
- fDeDxHisto->GetAxis(0)->SetTitle("signal");
- fDeDxHisto->GetAxis(1)->SetTitle("alpha (rad)");
- fDeDxHisto->GetAxis(2)->SetTitle("y (cm)");
- fDeDxHisto->GetAxis(3)->SetTitle("z (cm)");
- fDeDxHisto->GetAxis(4)->SetTitle("snp");
- fDeDxHisto->GetAxis(5)->SetTitle("tgl");
- fDeDxHisto->GetAxis(6)->SetTitle("ncls");
- fDeDxHisto->GetAxis(6)->SetTitle("pid");
- fDeDxHisto->GetAxis(7)->SetTitle("p (GeV/c)");
- fDeDxHisto->Sumw2();
-
- // Init cuts
- if(!fCutsMC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
- if(!fCutsRC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
-
- // init folder
- fAnalysisFolder = CreateFolder("folderDEdx","Analysis de/dx Folder");
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
-{
- // Fill dE/dx comparison information
-
- // Check selection cuts
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
-
- Double_t pid = -1;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3;
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4;
-
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
-
- Float_t dedx = infoRC->GetESDtrack()->GetTPCsignal();
- Int_t ncls = infoRC->GetESDtrack()->GetTPCNcls();
-
- const AliExternalTrackParam *innerParam = 0;
- if ((innerParam = infoRC->GetESDtrack()->GetInnerParam()) == 0) return;
-
- Double_t pt = innerParam->Pt();
- Double_t lam = TMath::ATan2(innerParam->Pz(),innerParam->Pt());
- Double_t p = pt/TMath::Cos(lam);
- Double_t alpha = innerParam->GetAlpha();
- Double_t y = innerParam->GetY();
- Double_t z = innerParam->GetZ();
- Double_t snp = innerParam->GetSnp();
- Double_t tgl = innerParam->GetTgl();
-
- Double_t vDeDxHisto[9] = {dedx,alpha,y,z,snp,tgl,ncls,pid,p};
- fDeDxHisto->Fill(vDeDxHisto);
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::ProcessTPCITS(AliMCInfo* const /*infoMC*/, AliESDRecInfo *const /*infoRC*/)
-{
- // Fill dE/dx comparison information
-
- AliDebug(AliLog::kWarning, "Warning: Not implemented");
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::ProcessConstrained(AliMCInfo* const /*infoMC*/, AliESDRecInfo *const /*infoRC*/)
-{
- // Fill dE/dx comparison information
-
- AliDebug(AliLog::kWarning, "Warning: Not implemented");
-}
-
-//_____________________________________________________________________________
-Long64_t AliComparisonDEdx::Merge(TCollection* const list)
-{
- // Merge list of objects (needed by PROOF)
-
- if (!list)
- return 0;
-
- if (list->IsEmpty())
- return 1;
-
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- // collection of generated histograms
- Int_t count=0;
- while((obj = iter->Next()) != 0)
- {
- AliComparisonDEdx* entry = dynamic_cast<AliComparisonDEdx*>(obj);
- if (entry == 0) continue;
-
- fDeDxHisto->Add(entry->fDeDxHisto);
- count++;
- }
-
-return count;
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
-{
- // Process comparison information
-
- if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);
- else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);
- else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);
- else {
- printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
- return;
- }
-}
-
-//_____________________________________________________________________________
-TH1F* AliComparisonDEdx::MakeResol(TH2F * his, Int_t integ, Bool_t type)
-{
- // Make resolution histograms
- TH1F *hisr, *hism;
- if (!gPad) new TCanvas;
- hisr = AliTreeDraw::CreateResHistoI(his,&hism,integ);
- if (type) return hism;
- else
- return hisr;
-}
-
-//_____________________________________________________________________________
-void AliComparisonDEdx::Analyse()
-{
- // Analyze comparison information and store output histograms
- // in the folder "folderDEdx"
- //
- TH1::AddDirectory(kFALSE);
- TObjArray *aFolderObj = new TObjArray;
-
- /*
- TH1F *hiss=0;
- TGraph2D * gr=0;
-
- // write results in the folder
- TCanvas * c = new TCanvas("can","TPC dedx");
- c->cd();
-
- hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,0);
- hiss->SetXTitle("Tan(#theta)");
- hiss->SetYTitle("#sigma_{dEdx}");
- hiss->Draw();
- hiss->SetName("TPCdEdxResolTan");
-
- aFolderObj->Add(hiss);
- //
- hiss = comp->MakeResol(comp->fTPCSignalNormTan,4,1);
- hiss->SetXTitle("Tan(#theta)");
- hiss->SetYTitle("<dEdx>");
- hiss->Draw();
- hiss->SetName("TPCdEdxMeanTan");
-
- aFolderObj->Add(hiss);
- //
- gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,4);
- gr->GetXaxis()->SetTitle("Tan(#theta)");
- gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr->GetZaxis()->SetTitle("<dEdx>");
- gr->SetName("TPCdEdxMeanTanPt_1");
- gr->GetHistogram()->Draw("colz");
-
- gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_1");
- aFolderObj->Add(gr->GetHistogram());
- //
- gr = AliMathBase::MakeStat2D(comp->fTPCSignalNormTanSPt,3,1,5);
- gr->GetXaxis()->SetTitle("Tan(#theta)");
- gr->GetYaxis()->SetTitle("#sqrt{p_{t}(GeV)}");
- gr->GetZaxis()->SetTitle("#sigma_{dEdx}");
- gr->SetName("TPCdEdxMeanTanPt_2");
- gr->GetHistogram()->Draw("colz");
-
- gr->GetHistogram()->SetName("TPCdEdxMeanTanPt_2");
- aFolderObj->Add(gr->GetHistogram());
- */
-
- // export objects to analysis folder
- fAnalysisFolder = ExportToFolder(aFolderObj);
-
- // delete only TObjrArray
- if(aFolderObj) delete aFolderObj;
-}
-
-//_____________________________________________________________________________
-TFolder* AliComparisonDEdx::ExportToFolder(TObjArray * array)
-{
- // recreate folder avery time and export objects to new one
- //
- AliComparisonDEdx * comp=this;
- TFolder *folder = comp->GetAnalysisFolder();
-
- TString name, title;
- TFolder *newFolder = 0;
- Int_t i = 0;
- Int_t size = array->GetSize();
-
- if(folder) {
- // get name and title from old folder
- name = folder->GetName();
- title = folder->GetTitle();
-
- // delete old one
- delete folder;
-
- // create new one
- newFolder = CreateFolder(name.Data(),title.Data());
- newFolder->SetOwner();
-
- // add objects to folder
- while(i < size) {
- newFolder->Add(array->At(i));
- i++;
- }
- }
-
-return newFolder;
-}
-
-
-//_____________________________________________________________________________
-TFolder* AliComparisonDEdx::CreateFolder(TString name,TString title) {
-// create folder for analysed histograms
-TFolder *folder = 0;
- folder = new TFolder(name.Data(),title.Data());
-
- return folder;
-}
-
+++ /dev/null
-#ifndef ALICOMPARISONDEdx_H
-#define ALICOMPARISONDEdx_H
-
-//------------------------------------------------------------------------------
-// Class to keep information from comparison of
-// reconstructed and MC particle tracks (TPC dE/dx).
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-class TCanvas;
-class TH1F;
-class TH2F;
-class TNamed;
-class TString;
-
-class AliMCInfo;
-class AliESDRecInfo;
-class AliESDEvent;
-class AliRecInfoCuts;
-class AliMCInfoCuts;
-
-#include "THnSparse.h"
-#include "AliComparisonObject.h"
-
-class AliComparisonDEdx : public AliComparisonObject {
-public :
- AliComparisonDEdx();
- AliComparisonDEdx(Char_t* name, Char_t* title, Int_t analysisMode, Bool_t hptGenerator);
- ~AliComparisonDEdx();
-
- // Init data members
- virtual void Init();
-
- // Execute analysis
- virtual void Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC);
-
- // Merge output objects (needed by PROOF)
- virtual Long64_t Merge(TCollection* const list);
-
- // Analyse output histograms
- virtual void Analyse();
-
- // Get analysis folder
- virtual TFolder* GetAnalysisFolder() const {return fAnalysisFolder;}
-
- // Create folder for analysed histograms
- TFolder *CreateFolder(TString folder = "folderDEdx",TString title = "Analysed DEdx histograms");
-
- // Export objects to folder
- TFolder *ExportToFolder(TObjArray * array=0);
-
- // Process events
- void ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC);
- void ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC); // not implemented
- void ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC); // not implemented
-
- // Selection cuts
- void SetAliRecInfoCuts(AliRecInfoCuts* const cuts=0) {fCutsRC = cuts;}
- void SetAliMCInfoCuts(AliMCInfoCuts* const cuts=0) {fCutsMC = cuts;}
-
- AliRecInfoCuts* GetAliRecInfoCuts() const {return fCutsRC;}
- AliMCInfoCuts* GetAliMCInfoCuts() const {return fCutsMC;}
-
- static TH1F* MakeResol(TH2F * his, Int_t integ, Bool_t type);
-
- //
- // TPC dE/dx
- //
- THnSparse* GetDeDxHisto() const {return fDeDxHisto;}
-
-private:
-
- // TPC dE/dx
- THnSparseF *fDeDxHisto; //-> signal:alpha:y:z:snp:tgl:ncls:pid:p
-
- // Selection cuts
- AliRecInfoCuts* fCutsRC; // selection cuts for reconstructed tracks
- AliMCInfoCuts* fCutsMC; // selection cuts for MC tracks
-
- // analysis folder
- TFolder *fAnalysisFolder; // folder for analysed histograms
-
- AliComparisonDEdx(const AliComparisonDEdx&); // not implemented
- AliComparisonDEdx& operator=(const AliComparisonDEdx&); // not implemented
-
- ClassDef(AliComparisonDEdx,1);
-};
-
-#endif
+++ /dev/null
-//------------------------------------------------------------------------------\r
-// Implementation of AliComparisonEff class. It keeps information from \r
-// comparison of reconstructed and MC particle tracks. In addtion, \r
-// it keeps selection cuts used during comparison. The comparison \r
-// information is stored in the ROOT histograms. Analysis of these \r
-// histograms can be done by using Analyse() class function. The result of \r
-// the analysis (histograms/graphs) are stored in the folder which is \r
-// a data member of AliComparisonEff.\r
-// \r
-// Author: J.Otwinowski 04/02/2008 \r
-//------------------------------------------------------------------------------\r
-\r
-/*\r
- \r
- // after running comparison task, read the file, and get component\r
- gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");\r
- LoadMyLibs();\r
- TFile f("Output.root");\r
- //AliComparisonEff * compObj = (AliComparisonEff*)f.Get("AliComparisonEff");\r
- AliComparisonEff * compObj = (AliComparisonEff*)cOutput->FindObject("AliComparisonEff");\r
-\r
- // Analyse comparison data\r
- compObj->Analyse();\r
-\r
- // the output histograms/graphs will be stored in the folder "folderEff" \r
- compObj->GetAnalysisFolder()->ls("*");\r
-\r
- // user can save whole comparison object (or only folder with anlysed histograms) \r
- // in the seperate output file (e.g.)\r
- TFile fout("Analysed_Eff.root","recreate");\r
- compObj->Write(); // compObj->GetAnalysisFolder()->Write();\r
- fout.Close();\r
-\r
-*/\r
-\r
-#include <TAxis.h>\r
-#include <TH1D.h>\r
-\r
-// \r
-#include "AliESDtrack.h"\r
-#include "AliRecInfoCuts.h" \r
-#include "AliMCInfoCuts.h" \r
-#include "AliLog.h" \r
-#include "AliESDVertex.h" \r
-#include "AliExternalTrackParam.h" \r
-#include "AliTracker.h" \r
-#include "AliMCInfo.h" \r
-#include "AliESDRecInfo.h" \r
-#include "AliComparisonEff.h" \r
-\r
-using namespace std;\r
-\r
-ClassImp(AliComparisonEff)\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonEff::AliComparisonEff():\r
- AliComparisonObject("AliComparisonEff"),\r
-\r
- // histograms\r
- \r
- fEffHisto(0),\r
-\r
- // Cuts \r
- fCutsRC(0), \r
- fCutsMC(0),\r
-\r
- // histogram folder \r
- fAnalysisFolder(0)\r
-{\r
- // default consttructor \r
-}\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonEff::AliComparisonEff(Char_t* name="AliComparisonEff",Char_t*title="AliComparisonEff",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):\r
- AliComparisonObject(name,title),\r
-\r
- // histograms\r
- fEffHisto(0),\r
-\r
- // Cuts \r
- fCutsRC(0), \r
- fCutsMC(0),\r
-\r
- // histogram folder \r
- fAnalysisFolder(0)\r
-{\r
- // named constructor\r
- //\r
- SetAnalysisMode(analysisMode);\r
- SetHptGenerator(hptGenerator);\r
-\r
- Init();\r
-}\r
-\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonEff::~AliComparisonEff()\r
-{\r
-// destructor\r
-\r
- if(fEffHisto) delete fEffHisto; fEffHisto=0;\r
-\r
- if(fCutsRC) delete fCutsRC; fCutsRC=0;\r
- if(fCutsMC) delete fCutsMC; fCutsMC=0;\r
-\r
- if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonEff::Init()\r
-{\r
- // Init histograms\r
- //\r
- Int_t nPtBins = 31;\r
- Double_t binsPt[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};\r
- Double_t ptMin = 0., ptMax = 10.;\r
-\r
- if(IsHptGenerator() == kTRUE) {\r
- nPtBins = 100;\r
- ptMin = 0.; ptMax = 100.;\r
- }\r
-\r
- //mceta:mcphi:mcpt:pid:isPrim:recStatus:findable\r
- Int_t binsEffHisto[7]={30,90,nPtBins,5,2,2,2};\r
- Double_t minEffHisto[7]={-1.5,0.,ptMin,0.,0.,0.,0.};\r
- Double_t maxEffHisto[7]={ 1.5,2.*TMath::Pi(), ptMax,5.,2.,2.,2.};\r
-\r
- fEffHisto = new THnSparseF("fEffHisto","mceta:mcphi:mcpt:pid:isPrim:recStatus:findable",7,binsEffHisto,minEffHisto,maxEffHisto);\r
- if(!IsHptGenerator()) fEffHisto->SetBinEdges(2,binsPt);\r
-\r
- fEffHisto->GetAxis(0)->SetTitle("eta");\r
- fEffHisto->GetAxis(1)->SetTitle("phi (rad)");\r
- fEffHisto->GetAxis(2)->SetTitle("pt (GeV/c)");\r
- fEffHisto->GetAxis(3)->SetTitle("pid");\r
- fEffHisto->GetAxis(4)->SetTitle("isPrim");\r
- fEffHisto->GetAxis(5)->SetTitle("recStatus");\r
- fEffHisto->GetAxis(6)->SetTitle("findable");\r
- fEffHisto->Sumw2();\r
-\r
- // init cuts\r
- if(!fCutsMC) \r
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");\r
- if(!fCutsRC) \r
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");\r
-\r
- // init folder\r
- fAnalysisFolder = CreateFolder("folderEff","Analysis Efficiency Folder");\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonEff::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC)\r
-{\r
- // Fill efficiency comparison information\r
-\r
- AliExternalTrackParam *track = 0;\r
- //Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]\r
- Double_t kMaxD = 123456.0; // max distance\r
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
- AliESDVertex vertexMC; // MC primary vertex\r
-\r
- // distance to Prim. vertex \r
- const Double_t* dv = infoMC->GetVDist(); \r
-\r
- Float_t mcpt = infoMC->GetParticle().Pt();\r
- Float_t mceta = infoMC->GetParticle().Eta();\r
- Float_t mcphi = infoMC->GetParticle().Phi();\r
- if(mcphi<0) mcphi += 2.*TMath::Pi();\r
-\r
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
- Bool_t findable = (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
- Bool_t recStatus = kFALSE;\r
- \r
- // calculate and set prim. vertex\r
- vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
- vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
- vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
- \r
- // Only 5 charged particle species (e,mu,pi,K,p)\r
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
-\r
- // transform Pdg to Pid\r
- // Pdg convension is different for hadrons and leptons \r
- // (e.g. K+/K- = 321/-321; e+/e- = -11/11 ) \r
- Double_t pid = -1;\r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4; \r
-\r
- if (infoRC->GetESDtrack() && infoRC->GetESDtrack()->GetTPCInnerParam()) \r
- {\r
- if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0)\r
- {\r
- //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);\r
- Double_t field[3]; track->GetBxByBz(field); \r
- Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);\r
- if(bDCAStatus) {\r
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())\r
- {\r
- recStatus = kTRUE;\r
- }\r
- }\r
- delete track;\r
- }\r
- }\r
-\r
- // Fill histograms\r
- Double_t vEffHisto[7] = { mceta, mcphi, mcpt, pid, isPrim, recStatus, findable }; \r
- fEffHisto->Fill(vEffHisto);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonEff::ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC)\r
-{\r
- // Fill efficiency comparison information\r
- \r
- Int_t clusterITS[200];\r
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
-\r
- Float_t mcpt = infoMC->GetParticle().Pt();\r
- Float_t mceta = infoMC->GetParticle().Eta();\r
- Float_t mcphi = infoMC->GetParticle().Phi();\r
- if(mcphi<0) mcphi += 2.*TMath::Pi();\r
-\r
- // distance to Prim. vertex \r
- const Double_t* dv = infoMC->GetVDist(); \r
-\r
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
- Bool_t findable = (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
- Bool_t recStatus =kFALSE;\r
- \r
- // Only 5 charged particle species (e,mu,pi,K,p)\r
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
-\r
- // transform Pdg to Pid\r
- // Pdg convension is different for hadrons and leptons \r
- // (e.g. K+/K- = 321/-321; e+/e- = -11/11 ) \r
- Double_t pid = -1;\r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4; \r
-\r
- if(!infoRC->GetESDtrack()) return;\r
- if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)>fCutsRC->GetMinNClustersITS())\r
- {\r
- infoRC->GetESDtrack()->GetImpactParameters(dca,cov);\r
- if(TMath::Abs(dca[0]) < fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1]) < fCutsRC->GetMaxDCAToVertexZ()) \r
- {\r
- recStatus =(infoRC->GetStatus(1)==3);\r
- } \r
- }\r
-\r
- // fill histograms\r
- Double_t vEffHisto[7] = { mceta, mcphi, mcpt, pid, isPrim, recStatus, findable }; \r
- fEffHisto->Fill(vEffHisto);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonEff::ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC)\r
-{\r
- // Fill efficiency comparison information\r
- AliExternalTrackParam *track = 0;\r
- //Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]\r
- Double_t kMaxD = 123456.0; // max distance\r
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
- AliESDVertex vertexMC; // MC primary vertex\r
-\r
- // distance to Prim. vertex \r
- const Double_t* dv = infoMC->GetVDist(); \r
-\r
- Float_t mcpt = infoMC->GetParticle().Pt();\r
- Float_t mceta = infoMC->GetParticle().Eta();\r
- Float_t mcphi = infoMC->GetParticle().Phi();\r
- if(mcphi<0) mcphi += 2.*TMath::Pi();\r
-\r
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
- Bool_t findable = (infoMC->GetRowsWithDigits()>fCutsMC->GetMinRowsWithDigits());\r
- Bool_t recStatus = kFALSE;\r
- \r
- // calculate and set prim. vertex\r
- vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
- vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
- vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
- \r
- // Only 5 charged particle species (e,mu,pi,K,p)\r
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
-\r
- // transform Pdg to Pid\r
- // Pdg convension is different for hadrons and leptons \r
- // (e.g. K+/K- = 321/-321; e+/e- = -11/11 ) \r
- Double_t pid = -1;\r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetEM() ) pid = 0; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetMuM() ) pid = 1; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetPiP() ) pid = 2; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetKP() ) pid = 3; \r
- if( TMath::Abs(infoMC->GetParticle().GetPdgCode())==fCutsMC->GetProt() ) pid = 4; \r
-\r
- // constrained parameters resolution\r
- if (!infoRC->GetESDtrack()) return;\r
- const AliExternalTrackParam * cparam = infoRC->GetESDtrack()->GetConstrainedParam();\r
- if(!cparam) return;\r
-\r
- if ((track = new AliExternalTrackParam(*cparam)) != 0)\r
- {\r
- //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);\r
- Double_t field[3]; track->GetBxByBz(field); \r
- Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);\r
- if(bDCAStatus) {\r
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())\r
- {\r
- recStatus = (infoRC->GetStatus(1)!=3);\r
- }\r
- }\r
- delete track;\r
- }\r
-\r
- // Fill histograms\r
- Double_t vEffHisto[7] = { mceta, mcphi, mcpt, pid, isPrim, recStatus, findable }; \r
- fEffHisto->Fill(vEffHisto);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonEff::Exec(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC)\r
-{\r
- // Process comparison information\r
-\r
- if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);\r
- else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);\r
- else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);\r
- else {\r
- printf("ERROR: AnalysisMode %d \n",fAnalysisMode);\r
- return;\r
- }\r
-}\r
-\r
-//_____________________________________________________________________________\r
-Long64_t AliComparisonEff::Merge(TCollection* const list) \r
-{\r
- // Merge list of objects (needed by PROOF)\r
-\r
- if (!list)\r
- return 0;\r
-\r
- if (list->IsEmpty())\r
- return 1;\r
-\r
- TIterator* iter = list->MakeIterator();\r
- TObject* obj = 0;\r
-\r
- // collection of generated histograms\r
-\r
- Int_t count=0;\r
- while((obj = iter->Next()) != 0) \r
- {\r
- AliComparisonEff* entry = dynamic_cast<AliComparisonEff*>(obj);\r
- if (entry == 0) continue; \r
- \r
- fEffHisto->Add(entry->fEffHisto);\r
- count++;\r
- }\r
-\r
-return count;\r
-}\r
- \r
-//_____________________________________________________________________________\r
-void AliComparisonEff::Analyse() \r
-{\r
- // Analyse comparison information and store output histograms\r
- // in the folder "folderEff" \r
- //\r
- TH1::AddDirectory(kFALSE);\r
- TObjArray *aFolderObj = new TObjArray;\r
-\r
- // efficiency vs pt\r
- fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.9); // eta range\r
- fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // primary tracks\r
-\r
- // rec efficiency vs pt\r
- TH1D *ptAll = fEffHisto->Projection(2);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed \r
- TH1D *ptRec = fEffHisto->Projection(2);\r
- ptRec->Divide(ptAll);\r
- ptRec->SetName("ptRecEff");\r
-\r
- aFolderObj->Add(ptRec);\r
-\r
-\r
- // rec efficiency vs pid vs pt\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(2.,2.); // pions\r
- TH1D *ptAllPi = fEffHisto->Projection(2);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *ptRecPi = fEffHisto->Projection(2);\r
- ptRecPi->Divide(ptAllPi);\r
- ptRecPi->SetName("ptRecEffPi");\r
-\r
- aFolderObj->Add(ptRecPi);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(3.,3.); // kaons\r
- TH1D *ptAllK = fEffHisto->Projection(2);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *ptRecK = fEffHisto->Projection(2);\r
- ptRecK->Divide(ptAllK);\r
- ptRecK->SetName("ptRecEffK");\r
-\r
- aFolderObj->Add(ptRecK);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(4.,4.); // protons\r
- TH1D *ptAllP = fEffHisto->Projection(2);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *ptRecP = fEffHisto->Projection(2);\r
- ptRecP->Divide(ptAllP);\r
- ptRecP->SetName("ptRecEffP");\r
-\r
- aFolderObj->Add(ptRecP);\r
- \r
- // findable efficiency \r
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); \r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(6)->SetRangeUser(1.,1.); // findable\r
- TH1D *ptAllF = fEffHisto->Projection(2);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.);\r
- fEffHisto->GetAxis(6)->SetRangeUser(1.,1.);\r
- TH1D *ptRecF = fEffHisto->Projection(2); // rec findable\r
- ptRecF->Divide(ptAllF);\r
- ptRecF->SetName("ptRecFindableEff");\r
-\r
- aFolderObj->Add(ptRecF);\r
-\r
- //\r
- // efficiency vs eta\r
- //\r
-\r
- fEffHisto->GetAxis(0)->SetRangeUser(-1.5,1.5); // eta range\r
- fEffHisto->GetAxis(2)->SetRangeUser(0.2,10.); // pt range\r
- fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // primary tracks\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); // all\r
- fEffHisto->GetAxis(6)->SetRangeUser(0.,1.); // all\r
-\r
- // rec efficiency vs eta\r
- TH1D *etaAll = fEffHisto->Projection(0);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed \r
- TH1D *etaRec = fEffHisto->Projection(0);\r
- etaRec->Divide(etaAll);\r
- etaRec->SetName("etaRecEff");\r
-\r
- aFolderObj->Add(etaRec);\r
-\r
- // rec efficiency vs pid vs eta\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(2.,2.); // pions\r
- TH1D *etaAllPi = fEffHisto->Projection(0);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *etaRecPi = fEffHisto->Projection(0);\r
- etaRecPi->Divide(etaAllPi);\r
- etaRecPi->SetName("etaRecEffPi");\r
-\r
- aFolderObj->Add(etaRecPi);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(3.,3.); // kaons\r
- TH1D *etaAllK = fEffHisto->Projection(0);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *etaRecK = fEffHisto->Projection(0);\r
- etaRecK->Divide(etaAllK);\r
- etaRecK->SetName("etaRecEffK");\r
-\r
- aFolderObj->Add(etaRecK);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(3)->SetRangeUser(4.,4.); // protons\r
- TH1D *etaAllP = fEffHisto->Projection(0);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // reconstructed\r
- TH1D *etaRecP = fEffHisto->Projection(0);\r
- etaRecP->Divide(etaAllP);\r
- etaRecP->SetName("etaRecEffP");\r
-\r
- aFolderObj->Add(etaRecP);\r
- \r
- //\r
- // findable efficiency \r
- //\r
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); \r
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); \r
- fEffHisto->GetAxis(6)->SetRangeUser(1.,1.); // findable\r
- TH1D *etaAllF = fEffHisto->Projection(0);\r
-\r
- fEffHisto->GetAxis(5)->SetRangeUser(1.,1.);\r
- fEffHisto->GetAxis(6)->SetRangeUser(1.,1.);\r
- TH1D *etaRecF = fEffHisto->Projection(0); // rec findable\r
- etaRecF->Divide(etaAllF);\r
- etaRecF->SetName("etaRecFindableEff");\r
-\r
- aFolderObj->Add(etaRecF);\r
-\r
- // export objects to analysis folder\r
- fAnalysisFolder = ExportToFolder(aFolderObj);\r
-\r
- // delete only TObjArray\r
- if(aFolderObj) delete aFolderObj;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-TFolder* AliComparisonEff::ExportToFolder(TObjArray * array) \r
-{\r
- // recreate folder avery time and export objects to new one\r
- //\r
- AliComparisonEff * comp=this;\r
- TFolder *folder = comp->GetAnalysisFolder();\r
-\r
- TString name, title;\r
- TFolder *newFolder = 0;\r
- Int_t i = 0;\r
- Int_t size = array->GetSize();\r
-\r
- if(folder) { \r
- // get name and title from old folder\r
- name = folder->GetName(); \r
- title = folder->GetTitle(); \r
-\r
- // delete old one\r
- delete folder;\r
-\r
- // create new one\r
- newFolder = CreateFolder(name.Data(),title.Data());\r
- newFolder->SetOwner();\r
-\r
- // add objects to folder\r
- while(i < size) {\r
- newFolder->Add(array->At(i));\r
- i++;\r
- }\r
- }\r
-\r
-return newFolder;\r
-}\r
-\r
-\r
-//_____________________________________________________________________________\r
-TFolder* AliComparisonEff::CreateFolder(TString name,TString title) { \r
-// create folder for analysed histograms\r
-//\r
-TFolder *folder = 0;\r
- folder = new TFolder(name.Data(),title.Data());\r
-\r
- return folder;\r
-}\r
+++ /dev/null
-#ifndef ALICOMPARISONEFF_H
-#define ALICOMPARISONEFF_H
-
-//------------------------------------------------------------------------------
-// Class to keep information from comparison of
-// reconstructed and MC particle tracks (TPC efficiency).
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-class TFile;
-class AliMCInfo;
-class AliESDRecInfo;
-class AliESDEvent;
-class AliESD;
-class AliRecInfoCuts;
-class AliMCInfoCuts;
-class TString;
-class AliESDVertex;
-class TNamed;
-
-#include "THnSparse.h"
-#include "AliComparisonObject.h"
-
-class AliComparisonEff : public AliComparisonObject {
-public :
- AliComparisonEff();
- AliComparisonEff(Char_t* name, Char_t* title, Int_t analysisMode, Bool_t hptGenerator);
- ~AliComparisonEff();
-
- // Init data members
- virtual void Init();
-
- // Execute analysis
- virtual void Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC);
-
- // Merge output objects (needed by PROOF)
- virtual Long64_t Merge(TCollection* const list);
-
- // Analyse output histograms
- virtual void Analyse();
-
- // Get analysis folder
- virtual TFolder* GetAnalysisFolder() const {return fAnalysisFolder;}
-
- // Create folder for analysed histograms
- TFolder *CreateFolder(TString folder = "folderEff",TString title = "Analysed Efficiency histograms");
-
- // Export objects to folder
- TFolder *ExportToFolder(TObjArray *array=0);
-
- // Process events
- void ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
-
- // Selection cuts
- void SetAliRecInfoCuts(AliRecInfoCuts* const cuts=0) {fCutsRC = cuts;}
- void SetAliMCInfoCuts(AliMCInfoCuts* const cuts=0) {fCutsMC = cuts;}
-
- // Getters
- AliRecInfoCuts* GetAliRecInfoCuts() const {return fCutsRC;}
- AliMCInfoCuts* GetAliMCInfoCuts() const {return fCutsMC;}
-
- THnSparseF* GetEffHisto() const {return fEffHisto;}
-
-private:
-
- // Control histograms
- THnSparseF *fEffHisto; //-> mceta:mcphi:mcpt:pid:isPrim:recStatus:findable
-
- // Global cuts objects
- AliRecInfoCuts* fCutsRC; // selection cuts for reconstructed tracks
- AliMCInfoCuts* fCutsMC; // selection cuts for MC tracks
-
- // analysis folder
- TFolder *fAnalysisFolder; // folder for analysed histograms
-
- AliComparisonEff(const AliComparisonEff&); // not implemented
- AliComparisonEff& operator=(const AliComparisonEff&); // not implemented
-
- ClassDef(AliComparisonEff,1);
-};
-
-#endif
+++ /dev/null
-/**************************************************************************
-* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
-* *
-* Author: The ALICE Off-line Project. *
-* Contributors are mentioned in the code where appropriate. *
-* *
-* Permission to use, copy, modify and distribute this software and its *
-* documentation strictly for non-commercial purposes is hereby granted *
-* without fee, provided that the above copyright notice appears in all *
-* copies and that both the copyright notice and this permission notice *
-* appear in the supporting documentation. The authors make no claims *
-* about the suitability of this software for any purpose. It is *
-* provided "as is" without express or implied warranty. *
-**************************************************************************/
-
-//------------------------------------------------------------------------------
-// Implementation of abstract AliComparisonObject class. It keeps information from
-// comparison of reconstructed and MC particle tracks.
-//
-// Author: J.Otwinowski 14/04/2008
-//------------------------------------------------------------------------------
-
-#include <iostream>
-
-#include "AliMCInfo.h"
-#include "AliESDRecInfo.h"
-#include "AliComparisonObject.h"
-
-using namespace std;
-
-ClassImp(AliComparisonObject)
-
-//_____________________________________________________________________________
-AliComparisonObject::AliComparisonObject():
- TNamed("AliComparisonObject","AliComparisonObject"),
- fAnalysisMode(-1),
- fHptGenerator(kFALSE)
-{
- // constructor
-}
-
-//_____________________________________________________________________________
-AliComparisonObject::AliComparisonObject(const char* name, const char* title):
- TNamed(name,title),
- fAnalysisMode(-1),
- fHptGenerator(kFALSE)
-{
- // constructor
-}
-
-//_____________________________________________________________________________
-AliComparisonObject::~AliComparisonObject(){
- // destructor
-}
-
+++ /dev/null
-#ifndef ALICOMPARISONOBJECT_H\r
-#define ALICOMPARISONOBJECT_H\r
-\r
-//------------------------------------------------------------------------------\r
-// Abstract class to keep information from comparison of \r
-// reconstructed and MC particle tracks. \r
-// \r
-// Author: J.Otwinowski 04/14/2008 \r
-//------------------------------------------------------------------------------\r
-\r
-#include "TNamed.h"\r
-#include "TFolder.h"\r
-\r
-class AliMCInfo;\r
-class AliESDRecInfo;\r
-\r
-class AliComparisonObject : public TNamed {\r
-public :\r
- AliComparisonObject(); \r
- AliComparisonObject(const char* name="AliComparisonObject", const char* title="AliComparisonObject"); \r
- virtual ~AliComparisonObject();\r
-\r
- // Init data members\r
- // call once before event loop\r
- virtual void Init() = 0;\r
-\r
- // Execute analysis\r
- // call in the event loop \r
- virtual void Exec(AliMCInfo* const infoMC=0, AliESDRecInfo* const infoRC=0) = 0;\r
-\r
- // Merge output objects (needed by PROOF) \r
- virtual Long64_t Merge(TCollection* const list=0) = 0;\r
-\r
- // Analyse output histograms\r
- virtual void Analyse() = 0;\r
-\r
- // Get output folder for analysed histograms\r
- virtual TFolder* GetAnalysisFolder() const = 0;\r
-\r
- // set and get analysisMode\r
- void SetAnalysisMode(Int_t analysisMode=0) {fAnalysisMode = analysisMode;} \r
- Int_t GetAnalysisMode() {return fAnalysisMode;}\r
-\r
- // set and get hpt generator \r
- void SetHptGenerator(Bool_t hptGenerator=kFALSE) {fHptGenerator = hptGenerator;}\r
- Bool_t IsHptGenerator() {return fHptGenerator;}\r
-\r
-protected: \r
-\r
- // analysis mode\r
- Int_t fAnalysisMode; // 0-TPC, 1-TPCITS, 2-Constrained\r
-\r
- // hpt generator\r
- Bool_t fHptGenerator; // hpt event generator\r
-\r
- ClassDef(AliComparisonObject,1);\r
-};\r
-\r
-#endif\r
+++ /dev/null
-//------------------------------------------------------------------------------
-// Implementation of AliComparisonRes class. It keeps information from
-// comparison of reconstructed and MC particle tracks. In addtion,
-// it keeps selection cuts used during comparison. The comparison
-// information is stored in the ROOT histograms. Analysis of these
-// histograms can be done by using Analyse() class function. The result of
-// the analysis (histograms/graphs) are stored in the folder which is
-// a data member of AliComparisonRes.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-/*
-
- // after running comparison task, read the file, and get component
- gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");
- LoadMyLibs();
-
- TFile f("Output.root");
- //AliComparisonRes * compObj = (AliComparisonRes*)f.Get("AliComparisonRes");
- AliComparisonRes * compObj = (AliComparisonRes*)cOutput->FindObject("AliComparisonRes");
-
- // analyse comparison data
- compObj->Analyse();
-
- // the output histograms/graphs will be stored in the folder "folderRes"
- compObj->GetAnalysisFolder()->ls("*");
-
- // user can save whole comparison object (or only folder with anlysed histograms)
- // in the seperate output file (e.g.)
- TFile fout("Analysed_Res.root","recreate");
- compObj->Write(); // compObj->GetAnalysisFolder()->Write();
- fout.Close();
-
-*/
-
-#include "TCanvas.h"
-#include "TH1.h"
-#include "TAxis.h"
-
-#include "AliComparisonRes.h"
-#include "AliESDRecInfo.h"
-#include "AliESDVertex.h"
-#include "AliESDtrack.h"
-#include "AliLog.h"
-#include "AliMCInfo.h"
-#include "AliMCInfoCuts.h"
-#include "AliRecInfoCuts.h"
-#include "AliTracker.h"
-#include "AliTreeDraw.h"
-
-using namespace std;
-
-ClassImp(AliComparisonRes)
-
-//_____________________________________________________________________________
-AliComparisonRes::AliComparisonRes():
- AliComparisonObject("AliComparisonRes"),
- fResolHisto(0),
- fPullHisto(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- //Init();
-}
-
-//_____________________________________________________________________________
-AliComparisonRes::AliComparisonRes(Char_t* name="AliComparisonRes", Char_t* title="AliComparisonRes",Int_t analysisMode=0,Bool_t hptGenerator=kFALSE):
-//AliComparisonRes::AliComparisonRes(Char_t* name, Char_t* title,Int_t analysisMode,Bool_t hptGenerator):
- AliComparisonObject(name,title),
- fResolHisto(0),
- fPullHisto(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // named constructor
- //
- SetAnalysisMode(analysisMode);
- SetHptGenerator(hptGenerator);
-
- Init();
-}
-
-//_____________________________________________________________________________
-AliComparisonRes::~AliComparisonRes()
-{
- // destructor
-
- if(fResolHisto) delete fResolHisto; fResolHisto=0;
- if(fPullHisto) delete fPullHisto; fPullHisto=0;
-
- if(fCutsRC) delete fCutsRC; fCutsRC=0;
- if(fCutsMC) delete fCutsMC; fCutsMC=0;
-
- if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::Init(){
-
- //
- // histogram bining
- //
-
- // set pt bins
- Int_t nPtBins = 31;
- Double_t binsPt[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
- Double_t ptMin = 0., ptMax = 10.;
-
-
- if(IsHptGenerator() == kTRUE) {
- nPtBins = 100;
- ptMin = 0.; ptMax = 100.;
- }
-
- // res_y:res_z:res_phi,res_lambda:res_1pt:y:z:eta:phi:pt
- Int_t binsResolHisto[10]={100,100,100,100,100,50,100,30,90,nPtBins};
- Double_t minResolHisto[10]={-1.,-1.,-0.03,-0.03,-0.2,-0.03,-20.,-1.5,0.,ptMin};
- Double_t maxResolHisto[10]={ 1., 1., 0.03, 0.03, 0.2, 0.03, 20., 1.5,2.*TMath::Pi(), ptMax};
-
-
- fResolHisto = new THnSparseF("fResolHisto","res_y:res_z:res_phi:res_lambda:res_1pt:y:z:eta:phi:pt",10,binsResolHisto,minResolHisto,maxResolHisto);
- if(!IsHptGenerator()) fResolHisto->SetBinEdges(9,binsPt);
-
- fResolHisto->GetAxis(0)->SetTitle("res_y (cm)");
- fResolHisto->GetAxis(1)->SetTitle("res_z (cm)");
- fResolHisto->GetAxis(2)->SetTitle("res_phi (rad)");
- fResolHisto->GetAxis(3)->SetTitle("res_lambda (rad)");
- fResolHisto->GetAxis(4)->SetTitle("res_1pt (GeV/c)^{-1}");
- fResolHisto->GetAxis(5)->SetTitle("y (cm)");
- fResolHisto->GetAxis(6)->SetTitle("z (cm)");
- fResolHisto->GetAxis(7)->SetTitle("eta");
- fResolHisto->GetAxis(8)->SetTitle("phi (rad)");
- fResolHisto->GetAxis(9)->SetTitle("pt (GeV/c)");
- fResolHisto->Sumw2();
-
- //pull_y:pull_z:pull_lambda:pull_phi:pull_1pt:y:z:eta:phi:pt
- Int_t binsPullHisto[10]={100,100,100,100,100,50,100,30,90,nPtBins};
- Double_t minPullHisto[10]={-5.,-5.,-5.,-5.,-5.,-0.03, -20.,-1.5, 0., ptMin};
- Double_t maxPullHisto[10]={ 5., 5., 5., 5., 5., 0.03, 20., 1.5, 2.*TMath::Pi(),ptMax};
-
- fPullHisto = new THnSparseF("fPullHisto","pull_y:pull_z:pull_phi:pull_lambda:pull_1pt:y:z:eta:phi:pt",10,binsPullHisto,minPullHisto,maxPullHisto);
- if(!IsHptGenerator()) fPullHisto->SetBinEdges(9,binsPt);
-
- fPullHisto->GetAxis(0)->SetTitle("res_y (cm)");
- fPullHisto->GetAxis(1)->SetTitle("res_z (cm)");
- fPullHisto->GetAxis(2)->SetTitle("res_phi (rad)");
- fPullHisto->GetAxis(3)->SetTitle("res_lambda (rad)");
- fPullHisto->GetAxis(4)->SetTitle("res_1pt (GeV/c)^{-1}");
- fPullHisto->GetAxis(5)->SetTitle("y (rad)");
- fPullHisto->GetAxis(6)->SetTitle("z (rad)");
- fPullHisto->GetAxis(7)->SetTitle("eta");
- fPullHisto->GetAxis(8)->SetTitle("phi (rad)");
- fPullHisto->GetAxis(9)->SetTitle("pt (GeV/c)");
- fPullHisto->Sumw2();
-
- // Init cuts
- if(!fCutsMC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
- if(!fCutsRC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
-
- // init folder
- fAnalysisFolder = CreateFolder("folderRes","Analysis Resolution Folder");
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
-{
- // Fill resolution comparison information
- AliExternalTrackParam *track = 0;
- //Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]
- Double_t kMaxD = 123456.0; // max distance
- AliESDVertex vertexMC; // MC primary vertex
-
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
-
- Float_t delta1PtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
- Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
- //Float_t delta1Pt2TPC, deltaY1PtTPC, deltaZ1PtTPC, deltaPhi1PtTPC, deltaTheta1PtTPC;
-
- Float_t mceta = infoMC->GetParticle().Eta();
- Float_t mcphi = infoMC->GetParticle().Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = infoMC->GetParticle().Pt();
-
- // distance to Prim. vertex
- const Double_t* dv = infoMC->GetVDist();
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
- if (!isPrim) return;
- //if (infoRC->GetStatus(1)!=3) return; // TPC refit
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
-
- // calculate and set prim. vertex
- vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );
- vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );
- vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );
-
- // calculate track parameters at vertex
- const AliExternalTrackParam *innerTPC = 0;
- if ((innerTPC = infoRC->GetESDtrack()->GetTPCInnerParam()) != 0)
- {
- if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )
- {
- //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);
- Double_t field[3]; track->GetBxByBz(field);
- Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);
-
- // Fill parametrisation histograms (only TPC track)
- if(bDCAStatus)
- {
- // select primaries
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
-
- deltaYTPC= track->GetY()-infoMC->GetParticle().Vy();
- deltaZTPC = track->GetZ()-infoMC->GetParticle().Vz();
- deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(infoMC->GetParticle().Pz(),infoMC->GetParticle().Pt());
- deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(infoMC->GetParticle().Py(),infoMC->GetParticle().Px());
- delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
-
- pullYTPC= (track->GetY()-infoMC->GetParticle().Vy()) / TMath::Sqrt(track->GetSigmaY2());
- pullZTPC = (track->GetZ()-infoMC->GetParticle().Vz()) / TMath::Sqrt(track->GetSigmaZ2());
- pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
- pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
- pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
-
- Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,delta1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
- fResolHisto->Fill(vResolHisto);
-
- Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
- fPullHisto->Fill(vPullHisto);
- }
-
- /*
- delta1Pt2TPC = (1/mcpt-innerTPC->OneOverPt())/TMath::Power(1+1/mcpt,2);
- deltaY1PtTPC= (infoMC->GetParticle().Vy()-innerTPC->GetY()) / (0.2+1/mcpt);
- deltaZ1PtTPC = (infoMC->GetParticle().Vz()-innerTPC->GetZ()) / (0.2+1/mcpt);
- deltaPhi1PtTPC = deltaPhiTPC / (0.1+1/mcpt);
- deltaTheta1PtTPC = deltaLambdaTPC / (0.1+1/mcpt);
-
- fPhiEtaPtTPC->Fill(TMath::ATan2(innerTPC->Py(),innerTPC->Px()), innerTPC->Eta(), innerTPC->Pt());
-
- f1Pt2ResolS1PtTPC->Fill(s1mcpt,delta1Pt2TPC);
- fYResolS1PtTPC->Fill(s1mcpt,deltaY1PtTPC);
- fZResolS1PtTPC->Fill(s1mcpt,deltaZ1PtTPC);
- fPhiResolS1PtTPC->Fill(s1mcpt,deltaPhi1PtTPC);
- fThetaResolS1PtTPC->Fill(s1mcpt,deltaTheta1PtTPC);
-
- fPtResolPtTPC->Fill(mcpt,deltaPtTPC);
- fPtPullPtTPC->Fill(mcpt,pullPtTPC);
- fPhiResolEtaTPC->Fill(eta,deltaPhiTPC);
- fPhiPullEtaTPC->Fill(eta,pullPhiTPC);
- fLambdaResolEtaTPC->Fill(eta,deltaLambdaTPC);
- fLambdaPullEtaTPC->Fill(eta,pullLambdaTPC);
- */
- }
- delete track;
- }
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
-{
- Int_t clusterITS[200];
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
-
- Float_t delta1PtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
- Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
-
- Float_t mceta = infoMC->GetParticle().Eta();
- Float_t mcphi = infoMC->GetParticle().Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = infoMC->GetParticle().Pt();
-
- // distance to Prim. vertex
- const Double_t* dv = infoMC->GetVDist();
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
- if (!isPrim) return;
- if (infoRC->GetStatus(1)!=3) return; // TPC refit
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
-
- infoRC->GetESDtrack()->GetImpactParameters(dca,cov);
-
- // select primaries
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
- deltaYTPC= infoRC->GetESDtrack()->GetY()-infoMC->GetParticle().Vy();
- deltaZTPC = infoRC->GetESDtrack()->GetZ()-infoMC->GetParticle().Vz();
- deltaLambdaTPC = TMath::ATan2(infoRC->GetESDtrack()->Pz(),infoRC->GetESDtrack()->Pt())-TMath::ATan2(infoMC->GetParticle().Pz(),infoMC->GetParticle().Pt());
- deltaPhiTPC = TMath::ATan2(infoRC->GetESDtrack()->Py(),infoRC->GetESDtrack()->Px())-TMath::ATan2(infoMC->GetParticle().Py(),infoMC->GetParticle().Px());
- delta1PtTPC = (infoRC->GetESDtrack()->OneOverPt()-1./mcpt)*mcpt;
-
- pullYTPC= (infoRC->GetESDtrack()->GetY()-infoMC->GetParticle().Vy()) / TMath::Sqrt(infoRC->GetESDtrack()->GetSigmaY2());
- pullZTPC = (infoRC->GetESDtrack()->GetZ()-infoMC->GetParticle().Vz()) / TMath::Sqrt(infoRC->GetESDtrack()->GetSigmaZ2());
- pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(infoRC->GetESDtrack()->GetSigmaTgl2());
- pullPhiTPC = deltaPhiTPC / TMath::Sqrt(infoRC->GetESDtrack()->GetSigmaSnp2());
- pull1PtTPC = (infoRC->GetESDtrack()->OneOverPt()-1./mcpt) / TMath::Sqrt(infoRC->GetESDtrack()->GetSigma1Pt2());
-
- Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,delta1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
- Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
-
- // TPC and ITS clusters in the system
- if(infoRC->GetESDtrack()->GetITSclusters(clusterITS)>fCutsRC->GetMinNClustersITS())
- {
- fResolHisto->Fill(vResolHisto);
- fPullHisto->Fill(vPullHisto);
- }
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC)
-{
- // Fill resolution comparison information (constarained parameters)
- //
- AliExternalTrackParam *track = 0;
- //Double_t field = AliTracker::GetBz(); // nominal Bz field [kG]
- Double_t kMaxD = 123456.0; // max distance
- AliESDVertex vertexMC; // MC primary vertex
-
- Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
-
- Float_t delta1PtTPC, deltaYTPC, deltaZTPC, deltaPhiTPC, deltaLambdaTPC;
- Float_t pull1PtTPC, pullYTPC, pullZTPC, pullPhiTPC, pullLambdaTPC;
-
- Float_t mceta = infoMC->GetParticle().Eta();
- Float_t mcphi = infoMC->GetParticle().Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = infoMC->GetParticle().Pt();
-
- // distance to Prim. vertex
- const Double_t* dv = infoMC->GetVDist();
- Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();
-
- // calculate and set prim. vertex
- vertexMC.SetXv( infoMC->GetParticle().Vx() - dv[0] );
- vertexMC.SetYv( infoMC->GetParticle().Vy() - dv[1] );
- vertexMC.SetZv( infoMC->GetParticle().Vz() - dv[2] );
-
- // Check selection cuts
- if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return;
- if (!isPrim) return;
- if (infoRC->GetStatus(1)!=3) return; // TPC refit
- if (!infoRC->GetESDtrack()) return;
- if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;
-
- // constrained parameters resolution
- const AliExternalTrackParam * cparam = infoRC->GetESDtrack()->GetConstrainedParam();
- if(!cparam) return;
-
- if ((track = new AliExternalTrackParam(*cparam)) != 0)
- {
- //Bool_t bDCAStatus = track->PropagateToDCA(&vertexMC,field,kMaxD,dca,cov);
- Double_t field[3]; track->GetBxByBz(field);
- Bool_t bDCAStatus = track->PropagateToDCABxByBz(&vertexMC,field,kMaxD,dca,cov);
- if(bDCAStatus) {
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
- deltaYTPC= track->GetY()-infoMC->GetParticle().Vy();
- deltaZTPC = track->GetZ()-infoMC->GetParticle().Vz();
- deltaLambdaTPC = TMath::ATan2(track->Pz(),track->Pt())-TMath::ATan2(infoMC->GetParticle().Pz(),infoMC->GetParticle().Pt());
- deltaPhiTPC = TMath::ATan2(track->Py(),track->Px())-TMath::ATan2(infoMC->GetParticle().Py(),infoMC->GetParticle().Px());
- delta1PtTPC = (track->OneOverPt()-1./mcpt)*mcpt;
-
- pullYTPC= (track->GetY()-infoMC->GetParticle().Vy()) / TMath::Sqrt(track->GetSigmaY2());
- pullZTPC = (track->GetZ()-infoMC->GetParticle().Vz()) / TMath::Sqrt(track->GetSigmaZ2());
- pullLambdaTPC = deltaLambdaTPC / TMath::Sqrt(track->GetSigmaTgl2());
- pullPhiTPC = deltaPhiTPC / TMath::Sqrt(track->GetSigmaSnp2());
- pull1PtTPC = (track->OneOverPt()-1./mcpt) / TMath::Sqrt(track->GetSigma1Pt2());
-
- Double_t vResolHisto[10] = {deltaYTPC,deltaZTPC,deltaPhiTPC,deltaLambdaTPC,delta1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
- fResolHisto->Fill(vResolHisto);
-
- Double_t vPullHisto[10] = {pullYTPC,pullZTPC,pullPhiTPC,pullLambdaTPC,pull1PtTPC,infoMC->GetParticle().Vy(),infoMC->GetParticle().Vz(),mceta,mcphi,mcpt};
- fPullHisto->Fill(vPullHisto);
- }
- }
- delete track;
- }
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::Exec(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC){
-
- // Process comparison information
- if(GetAnalysisMode() == 0) ProcessTPC(infoMC,infoRC);
- else if(GetAnalysisMode() == 1) ProcessTPCITS(infoMC,infoRC);
- else if(GetAnalysisMode() == 2) ProcessConstrained(infoMC,infoRC);
- else {
- printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
- return;
- }
-}
-
-//_____________________________________________________________________________
-TH1F* AliComparisonRes::MakeResol(TH2F * his, Int_t integ, Bool_t type){
- // Create resolution histograms
-
- TH1F *hisr, *hism;
- if (!gPad) new TCanvas;
- hisr = AliTreeDraw::CreateResHistoI(his,&hism,integ);
- if (type) return hism;
- else
- return hisr;
-}
-
-//_____________________________________________________________________________
-void AliComparisonRes::Analyse(){
- // Analyse comparison information and store output histograms
- // in the folder "folderRes"
- //
- TH1::AddDirectory(kFALSE);
- TH1F *h=0;
- TH2F *h2D=0;
- TObjArray *aFolderObj = new TObjArray;
-
- // write results in the folder
- TCanvas * c = new TCanvas("Phi resol Tan","Phi resol Tan");
- c->cd();
-
- char name[256];
- char res_title[256] = {"res_y:res_z:res_lambda:res_phi:res_1pt:y:z:eta:phi:pt"} ;
- char pull_title[256] = {"pull_y:pull_z:pull_lambda:pull_phi:pull_1pt:y:z:eta:phi:pt"};
-
- for(Int_t i=0; i<5; i++)
- {
- for(Int_t j=5; j<10; j++)
- {
- if(j==7) fResolHisto->GetAxis(9)->SetRangeUser(0.2,10.); // pt threshold
- if(j==9) fResolHisto->GetAxis(7)->SetRangeUser(-0.9,0.9); // eta window
-
- h2D = (TH2F*)fResolHisto->Projection(i,j);
- h = AliComparisonRes::MakeResol(h2D,1,0);
- sprintf(name,"h_res_%d_vs_%d",i,j);
- h->SetName(name);
- h->SetTitle(res_title);
-
- aFolderObj->Add(h);
-
- h = AliComparisonRes::MakeResol(h2D,1,1);
- sprintf(name,"h_mean_res_%d_vs_%d",i,j);
- h->SetName(name);
- h->SetTitle(res_title);
-
- aFolderObj->Add(h);
-
- if(j==7) fResolHisto->GetAxis(9)->SetRangeUser(0.0,10.);
- if(j==9) fResolHisto->GetAxis(7)->SetRangeUser(-1.5,1.5);
-
- //
- if(j==7) fPullHisto->GetAxis(9)->SetRangeUser(0.2,10.);
- if(j==9) fPullHisto->GetAxis(7)->SetRangeUser(-0.9,0.9);
-
- h2D = (TH2F*)fPullHisto->Projection(i,j);
- h = AliComparisonRes::MakeResol(h2D,1,0);
- sprintf(name,"h_pull_%d_vs_%d",i,j);
- h->SetName(name);
- h->SetTitle(pull_title);
-
- aFolderObj->Add(h);
-
- h = AliComparisonRes::MakeResol(h2D,1,1);
- sprintf(name,"h_mean_pull_%d_vs_%d",i,j);
- h->SetName(name);
- h->SetTitle(pull_title);
-
- aFolderObj->Add(h);
-
- if(j==7) fPullHisto->GetAxis(9)->SetRangeUser(0.0,10.);
- if(j==9) fPullHisto->GetAxis(7)->SetRangeUser(-1.5,1.5);
- }
- }
-
- // export objects to analysis folder
- fAnalysisFolder = ExportToFolder(aFolderObj);
-
- // delete only TObjArray
- if(aFolderObj) delete aFolderObj;
-}
-
-//_____________________________________________________________________________
-TFolder* AliComparisonRes::ExportToFolder(TObjArray * array)
-{
- // recreate folder avery time and export objects to new one
- //
- AliComparisonRes * comp=this;
- TFolder *folder = comp->GetAnalysisFolder();
-
- TString name, title;
- TFolder *newFolder = 0;
- Int_t i = 0;
- Int_t size = array->GetSize();
-
- if(folder) {
- // get name and title from old folder
- name = folder->GetName();
- title = folder->GetTitle();
-
- // delete old one
- delete folder;
-
- // create new one
- newFolder = CreateFolder(name.Data(),title.Data());
- newFolder->SetOwner();
-
- // add objects to folder
- while(i < size) {
- newFolder->Add(array->At(i));
- i++;
- }
- }
-
-return newFolder;
-}
-
-//_____________________________________________________________________________
-Long64_t AliComparisonRes::Merge(TCollection* const list)
-{
- // Merge list of objects (needed by PROOF)
-
- if (!list)
- return 0;
-
- if (list->IsEmpty())
- return 1;
-
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- // collection of generated histograms
- Int_t count=0;
- while((obj = iter->Next()) != 0)
- {
- AliComparisonRes* entry = dynamic_cast<AliComparisonRes*>(obj);
- if (entry == 0) continue;
-
- fResolHisto->Add(entry->fResolHisto);
- fPullHisto->Add(entry->fPullHisto);
-
- count++;
- }
-
-return count;
-}
-
-//_____________________________________________________________________________
-TFolder* AliComparisonRes::CreateFolder(TString name,TString title) {
-// create folder for analysed histograms
-//
-TFolder *folder = 0;
- folder = new TFolder(name.Data(),title.Data());
-
- return folder;
-}
+++ /dev/null
-#ifndef ALICOMPARISONRES_H
-#define ALICOMPARISONRES_H
-
-//------------------------------------------------------------------------------
-// Class to keep information from comparison of
-// reconstructed and MC particle tracks (TPC resolution).
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-class TString;
-class TNamed;
-class TCanvas;
-class TH1F;
-class TH2F;
-
-class AliESDVertex;
-class AliESDtrack;
-class AliMCInfo;
-class AliESDRecInfo;
-class AliESDEvent;
-class AliMCInfoCuts;
-class AliRecInfoCuts;
-
-#include "THnSparse.h"
-#include "AliComparisonObject.h"
-
-class AliComparisonRes : public AliComparisonObject {
-public :
- AliComparisonRes();
- AliComparisonRes(Char_t* name, Char_t* title, Int_t analysisMode, Bool_t hptGenerator);
- virtual ~AliComparisonRes();
-
- // Init data members
- virtual void Init();
-
- // Execute analysis
- virtual void Exec(AliMCInfo* const infoMC, AliESDRecInfo *const infoRC);
-
- // Merge output objects (needed by PROOF)
- virtual Long64_t Merge(TCollection* const list);
-
- // Analyse output histograms
- virtual void Analyse();
-
- // Get analysis folder
- virtual TFolder* GetAnalysisFolder() const {return fAnalysisFolder;}
-
- // Process events
- void ProcessConstrained(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessTPC(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
- void ProcessTPCITS(AliMCInfo* const infoMC, AliESDRecInfo* const infoRC);
-
- // Create folder for analysed histograms
- TFolder *CreateFolder(TString folder = "folderRes",TString title = "Analysed Resolution histograms");
-
- // Export objects to folder
- TFolder *ExportToFolder(TObjArray * array=0);
-
- // Selection cuts
- void SetAliRecInfoCuts(AliRecInfoCuts* const cuts=0) {fCutsRC = cuts;}
- void SetAliMCInfoCuts(AliMCInfoCuts* const cuts=0) {fCutsMC = cuts;}
-
- AliRecInfoCuts* GetAliRecInfoCuts() const {return fCutsRC;}
- AliMCInfoCuts* GetAliMCInfoCuts() const {return fCutsMC;}
-
- static TH1F* MakeResol(TH2F * his, Int_t integ, Bool_t type);
-
- // getters
- //
- THnSparse *GetResolHisto() const { return fResolHisto; }
- THnSparse *GetPullHisto() const { return fPullHisto; }
-
-private:
- //
- // Control histograms
- // 5 track parameters (details in STEER/AliExternalTrackParam.h)
- //
-
- // resolution histogram
- THnSparseF *fResolHisto; //-> res_y:res_z:res_phi:res_lambda:res_1pt:y:z:eta:phi:pt
-
- // pull histogram
- THnSparseF *fPullHisto; //-> pull_y:pull_z:pull_phi:pull_lambda:pull_1pt:y:z:eta:phi:pt
-
- // Global cuts objects
- AliRecInfoCuts* fCutsRC; // selection cuts for reconstructed tracks
- AliMCInfoCuts* fCutsMC; // selection cuts for MC tracks
-
- // analysis folder
- TFolder *fAnalysisFolder; // folder for analysed histograms
-
- AliComparisonRes(const AliComparisonRes&); // not implemented
- AliComparisonRes& operator=(const AliComparisonRes&); // not implemented
-
- ClassDef(AliComparisonRes,1);
-};
-
-#endif
+++ /dev/null
-/**************************************************************************\r
-* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
-* *\r
-* Author: The ALICE Off-line Project. *\r
-* Contributors are mentioned in the code where appropriate. *\r
-* *\r
-* Permission to use, copy, modify and distribute this software and its *\r
-* documentation strictly for non-commercial purposes is hereby granted *\r
-* without fee, provided that the above copyright notice appears in all *\r
-* copies and that both the copyright notice and this permission notice *\r
-* appear in the supporting documentation. The authors make no claims *\r
-* about the suitability of this software for any purpose. It is *\r
-* provided "as is" without express or implied warranty. *\r
-**************************************************************************/\r
-\r
-//------------------------------------------------------------------------------\r
-// Implementation of the AliComparisonTask class. It compares properties of the \r
-// reconstructed and MC particle tracks under several conditions. \r
-// As the input it requires the tree with AliRecInfo and AliMCInfo branches. Such\r
-// tree can be prepared in advance by runing AliGenInfoMaker and then AliRecInfoMaker\r
-// (details in description of these classes).\r
-//\r
-// The comparison output objects deriving from AliComparisonObject \r
-// (e.g. AliComparisonRes, AliComparisonEff, AliComparisonDEdxA, AliComparisonDCA ...) \r
-// are stored in the output file (details in description of these classes).\r
-// \r
-// Author: J.Otwinowski 04/02/2008 \r
-//------------------------------------------------------------------------------\r
-\r
-#include "iostream"\r
-\r
-#include "TChain.h"\r
-#include "TTree.h"\r
-#include "TH1F.h"\r
-#include "TCanvas.h"\r
-#include "TList.h"\r
-#include "TFile.h"\r
-\r
-#include "AliAnalysisTask.h"\r
-#include "AliAnalysisManager.h"\r
-#include "AliESDEvent.h"\r
-#include "AliESDInputHandler.h"\r
-#include "AliESDVertex.h"\r
-#include "AliMagF.h"\r
-#include "AliTracker.h"\r
-#include "AliGeomManager.h"\r
-\r
-#include "AliMCInfo.h"\r
-#include "AliESDRecInfo.h"\r
-#include "AliMCInfoCuts.h"\r
-#include "AliRecInfoCuts.h"\r
-#include "AliComparisonRes.h"\r
-#include "AliComparisonEff.h"\r
-#include "AliComparisonDEdx.h"\r
-#include "AliComparisonDCA.h"\r
-#include "AliComparisonObject.h"\r
-#include "AliComparisonTask.h"\r
-\r
-using namespace std;\r
-\r
-ClassImp(AliComparisonTask)\r
-\r
-Int_t AliComparisonTask::fEvtNumber = 0;\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonTask::AliComparisonTask(const char *name) \r
- : AliAnalysisTask(name, "")\r
- , fTree(0)\r
- , fInfoMC(0)\r
- , fInfoRC(0)\r
- , fOutput(0)\r
- , fPitList(0)\r
- , fCompList(0)\r
-{\r
- // Constructor\r
-\r
- // Define input and output slots here\r
- DefineInput(0, TChain::Class());\r
- DefineOutput(0, TList::Class());\r
-\r
- // create the list for comparison objects\r
- fCompList = new TList;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-AliComparisonTask::~AliComparisonTask()\r
-{\r
- if(fOutput) delete fOutput; fOutput =0; \r
- if(fCompList) delete fCompList; fCompList =0; \r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonTask::ConnectInputData(Option_t *) \r
-{\r
- // Connect input data \r
- // Called once\r
-\r
- fTree = dynamic_cast<TTree*> (GetInputData(0));\r
- if (!fTree) {\r
- Printf("ERROR: Could not read chain from input slot 0");\r
- } else {\r
- fTree->SetBranchStatus("*",1);\r
- }\r
-\r
- if(fTree->GetBranch("MC") && fTree->GetBranch("RC")) {\r
- fTree->GetBranch("MC")->SetAddress(&fInfoMC);\r
- fTree->GetBranch("RC")->SetAddress(&fInfoRC);\r
- } else {\r
- Printf("ERROR: Could not get MC and RC branches");\r
- }\r
-}\r
-\r
-//_____________________________________________________________________________\r
-Bool_t AliComparisonTask::AddComparisonObject(AliComparisonObject *pObj) \r
-{\r
- // add comparison object to the list\r
- if(pObj == 0) {\r
- Printf("ERROR: Could not add comparison object");\r
- return kFALSE;\r
- }\r
-\r
- // add object to the list\r
- fCompList->AddLast(pObj);\r
- \r
-return kTRUE;\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonTask::CreateOutputObjects()\r
-{\r
- // Create histograms\r
- // Called once\r
-\r
- // create output list\r
- fOutput = new TList;\r
- fOutput->SetOwner();\r
- fPitList = fOutput->MakeIterator();\r
-\r
- AliComparisonObject *pObj=0;\r
- Int_t count=0;\r
-\r
- // add comparison objects to the output\r
- TIterator *pitCompList = fCompList->MakeIterator();\r
- pitCompList->Reset();\r
- while(( pObj = (AliComparisonObject *)pitCompList->Next()) != NULL) {\r
- fOutput->Add(pObj);\r
- count++;\r
- }\r
- Printf("CreateOutputObjects(): Number of output comparison objects: %d \n", count);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-Bool_t AliComparisonTask::ReadEntry(Int_t evt) \r
-{\r
-// Read entry from the tree\r
- Long64_t centry = fTree->LoadTree(evt);\r
- if(centry < 0) return kFALSE;\r
-\r
- if(fTree->GetBranch("MC") && fTree->GetBranch("RC")) {\r
- fTree->GetBranch("MC")->SetAddress(&fInfoMC);\r
- fTree->GetBranch("RC")->SetAddress(&fInfoRC);\r
- } else {\r
- Printf("ERROR: Could not get MC and RC branches");\r
- return kFALSE;\r
- }\r
- fTree->GetEntry(evt);\r
-\r
-return kTRUE;\r
-}\r
-//_____________________________________________________________________________\r
-void AliComparisonTask::Exec(Option_t *) \r
-{\r
- // Main loop\r
- // Called for each event\r
-\r
- AliComparisonObject *pObj=0;\r
-\r
- if (!fInfoMC && !fInfoRC) {\r
- Printf("ERROR: fInfoMC && fInfoRC not available");\r
- return;\r
- }\r
-\r
- // Process comparison\r
- Bool_t status = ReadEntry(fEvtNumber);\r
- if(status == kTRUE) \r
- {\r
- fPitList->Reset();\r
- while(( pObj = (AliComparisonObject *)fPitList->Next()) != NULL) {\r
- pObj->Exec(fInfoMC,fInfoRC);\r
- }\r
- }\r
-\r
- if( !( fEvtNumber % 10000) ) { \r
- cout << fEvtNumber << endl;\r
- }\r
- fEvtNumber++;\r
-\r
- // Post output data.\r
- PostData(0, fOutput);\r
-}\r
-\r
-//_____________________________________________________________________________\r
-void AliComparisonTask::Terminate(Option_t *) \r
-{\r
- // Called one at the end \r
- \r
- // check output data\r
- fOutput = dynamic_cast<TList*> (GetOutputData(0));\r
- if (!fOutput) {\r
- Printf("ERROR: AliComparisonTask::Terminate(): Output data not avaiable GetOutputData(0)==0x0 ..." );\r
- return;\r
- }\r
-}\r
+++ /dev/null
-#ifndef ALICOMPARISONRESTASK_H
-#define ALICOMPARISONRESTASK_H
-
-//------------------------------------------------------------------------------
-// Class to compare properties of reconstructed and MC particle tracks.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-class AliComparisonObject;
-class AliMagF;
-class TList;
-
-#include "AliAnalysisTask.h"
-#include "AliMCInfo.h"
-#include "AliESDRecInfo.h"
-
-class AliComparisonTask : public AliAnalysisTask {
- public:
- AliComparisonTask(const char *name = "AliComparisonTask");
- virtual ~AliComparisonTask();
-
- virtual void ConnectInputData(Option_t *);
- virtual void CreateOutputObjects();
- virtual void Exec(Option_t *option);
- virtual void Terminate(Option_t *);
-
- // Read TTree entry (event by event)
- Bool_t ReadEntry(Int_t evt);
-
- // Set comparison objects
- Bool_t AddComparisonObject(AliComparisonObject* comp);
-
- private:
- TTree* fTree; //! input tree
- AliMCInfo *fInfoMC; //! AliMCInfo object
- AliESDRecInfo *fInfoRC; //! AliESDRecInfo object
-
- TList* fOutput; //! list send on output slot 0
- static Int_t fEvtNumber; //! event number
- TIterator *fPitList; //! iterator over the output objetcs
- TList *fCompList; // list of comparison objects
-
- AliComparisonTask(const AliComparisonTask&); // not implemented
- AliComparisonTask& operator=(const AliComparisonTask&); // not implemented
-
- ClassDef(AliComparisonTask, 1); // example of analysis
-};
-
-#endif
AliTPCseed seed;
seed.Set(param->GetX(),param->GetAlpha(),param->GetParameter(),param->GetCovariance());
Float_t alpha= TMath::ATan2(refIn->Y(),refIn->X());
- seed.Rotate(alpha-seed.GetAlpha());
+ if(seed.Rotate(alpha-seed.GetAlpha())==kFALSE) return;
seed.SetMass(mass);
for (Float_t xlayer= seed.GetX(); xlayer>refIn->R(); xlayer-=step){
seed.PropagateTo(xlayer);
AliTPCseed seed;
seed.Set(param->GetX(),param->GetAlpha(),param->GetParameter(),param->GetCovariance());
Float_t alpha= TMath::ATan2(refIn->Y(),refIn->X());
- seed.Rotate(alpha-seed.GetAlpha());
+ if(seed.Rotate(alpha-seed.GetAlpha()) == kFALSE) return;
seed.SetMass(mass);
for (Float_t xlayer= seed.GetX(); xlayer>refIn->R(); xlayer-=step){
seed.PropagateTo(xlayer);
// \r
// MC histograms for efficiency studies\r
//\r
- if(!mcEvent) return;\r
+ if(mcEvent) { \r
\r
AliStack *stack = mcEvent->Stack();\r
- if (!stack) {\r
- AliDebug(AliLog::kError, "Stack not available");\r
- return;\r
- }\r
+ if (stack) {\r
\r
Int_t nPart = stack->GetNtrack();\r
//Int_t nPart = stack->GetNprimary();\r
Double_t vEffSecHisto[10] = { mceta, mcphi, mcpt, pid, recStatus, findable, mcR, mother_phi, mother_eta, charge }; \r
fEffSecHisto->Fill(vEffSecHisto);\r
}\r
+ }\r
+ }\r
\r
if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;\r
if (!outfile) return -1;
TChain* chain = new TChain("tpcQA");
+ if(!chain) return -1;
+
ifstream in;
in.open(infilelist);
in.close();
//TTree *tree = chain;
TTree *tree = chain->CopyTree("1");
+ if(!tree) return -1;
if (chain) { delete chain; chain=0; }
//TGraph* graph = dynamic_cast<TGraph*>(tree->DrawClone("run:run"));
//TGraph *graph = (TGraph*)gPad->GetPrimitive("Graph");
TFile* out = new TFile(outfile,"RECREATE");
+ if(!out) return -1;
+
out->cd();
const Char_t* condition = "meanTPCncl>0";
SaveGraph(tree,"meanTPCnclF","run",condition);
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_3_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_3_5_7"));
+ if(!his3D) return 8;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
}
} else {
his2D = pTPC->GetTPCTrackHisto()->Projection(3,5);
}
-
+ if(!his2D) return 8;
his2D->FitSlicesY(0,0,-1,10,"QNR",&arrayFit);
if (pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_3_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_3_5_7"));
+ if(!his3D) return 16;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
}
his2D = pTPC->GetTPCTrackHisto()->Projection(3,5);
pTPC->GetTPCTrackHisto()->GetAxis(8)->SetRangeUser(-1.5,1.5);
}
+ if(!his2D) return 16;
his2D->FitSlicesY(0,0,-1,10,"QNR",&arrayFit);
delete his2D;
if (pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_3_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_3_5_7"));
+ if(!his3D) return 32;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
}
his2D = pTPC->GetTPCTrackHisto()->Projection(3,5);
pTPC->GetTPCTrackHisto()->GetAxis(8)->SetRangeUser(-1.5,1.5);
}
+ if(!his2D) return 32;
+
his2D->FitSlicesY(0,0,-1,10,"QNR",&arrayFit);
delete his2D;
his1D = (TH1*) arrayFit.At(1);
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_0_5_7")) {
his3D0 = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_0_5_7"));
+ if(!his3D0) return 1;
his3D0->GetYaxis()->SetRangeUser(-1,1);
his3D0->GetZaxis()->SetRangeUser(0.25,10);
}
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_1_5_7")) {
his3D1 = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_1_5_7"));
+ if(!his3D1) return 1;
his3D1->GetYaxis()->SetRangeUser(-1,1);
his3D1->GetZaxis()->SetRangeUser(0.25,10);
}
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_2_5_7")) {
his3D2 = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_2_5_7"));
+ if(!his3D2) return 1;
his3D2->GetYaxis()->SetRangeUser(-1,1);
his3D2->GetZaxis()->SetRangeUser(0.25,10);
his3D2->GetXaxis()->SetRangeUser(0.4,1.1);
} else {
hprof = pTPC->GetTPCTrackHisto()->Projection(0,5)->ProfileX();
}
+ if(!hprof) return 1;
hprof->Fit(fpol1,"QNR","QNR",0.1,0.8);
slopeATPCncl= fpol1->GetParameter(1);
pTPC->GetTPCTrackHisto()->GetAxis(2)->SetRangeUser(0.4,1.1);
his1D = pTPC->GetTPCTrackHisto()->Projection(2,5)->ProfileX();
}
+ if(!his1D) return 1;
his1D->Fit(fpol1,"QNR","QNR",0.1,0.8);
slopeATPCnclF= fpol1->GetParameter(1);
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_4_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_4_5_7"));
+ if(!his3D) return 2;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
- }
+ }
+
if (his3D && !fgForceTHnSparse) {
his2D = dynamic_cast<TH2*>(his3D->Project3D("xy"));
} else {
his2D = pTPC->GetTPCTrackHisto()->Projection(4,5);
}
+ if(!his2D) return 2;
static TF1 *fpol1 = new TF1("fpol1","pol1");
TObjArray arrayFit;
if (pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_4_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_4_5_7"));
+ if(!his3D) return 64;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
- }
+ }
+
if (his3D && !fgForceTHnSparse) {
his2D = dynamic_cast<TH2*>(his3D->Project3D("xy"));
} else {
his2D = pTPC->GetTPCTrackHisto()->Projection(4,5);
}
+ if(!his2D) return 64;
static TF1 *fpol1 = new TF1("fpol1","pol1");
TObjArray arrayFit;
if (pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_4_5_7")) {
his3D = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_4_5_7"));
+ if(!his3D) return 128;
his3D->GetYaxis()->SetRangeUser(-1,1);
his3D->GetZaxis()->SetRangeUser(0.25,10);
}
} else {
his2D = pTPC->GetTPCTrackHisto()->Projection(4,5);
}
+ if(!his2D) return 128;
static TF1 *fpol1 = new TF1("fpol1","pol1");
TObjArray arrayFit;
} else {
his1D = pTPCgain->GetDeDxHisto()->Projection(0);
}
+ if(!his1D) return 4;
his1D->Fit(&gausFit,"QN","QN");
meanMIP = gausFit.GetParameter(1);
} else {
his2D = pTPCgain->GetDeDxHisto()->Projection(0,5);
}
+ if(!his2D) return 4;
+
TF1 * fpol = new TF1("fpol","pol1");
TObjArray arrayFit;
his2D->FitSlicesY(0,0,-1,10,"QN",&arrayFit);
} else {
his2D = pTPCgain->GetDeDxHisto()->Projection(0,5);
}
+ if(!his2D) return 4;
+
TF1 * fpolA = new TF1("fpolA","pol1");
TObjArray arrayFitA;
his2D->FitSlicesY(0,0,-1,10,"QN",&arrayFit);
} else {
his2D = pTPCgain->GetDeDxHisto()->Projection(0,1);
}
+ if(!his2D) return 4;
+
for(Int_t i = 0; i < 18; i++) { // loop over sectors; correct mapping to be checked!
//TH1* his1D=0;
Float_t phiLow = -TMath::Pi() + i*(20./360.)*(2*TMath::Pi());
} else {
his2D = pTPCgain->GetDeDxHisto()->Projection(0,1);
}
+ if(!his2D) return 4;
+
for(Int_t i = 0; i < 18; i++) { // loop over sectors; correct mapping to be checked!
//TH1* his1D=0;
Float_t phiLow = -TMath::Pi() + i*(20./360.)*(2*TMath::Pi());
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(6);
}
+ if(!his1D) return 1;
+
vertAll = his1D->GetEntries();
vertOK = his1D->GetBinContent(2);
if (vertAll>=1) {
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(0);
}
+ if(!his1D) return 1;
+
meanVertX = his1D->GetMean();
rmsVertX = his1D->GetRMS();
delete his1D;
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(1);
}
+ if(!his1D) return 1;
+
meanVertY = his1D->GetMean();
rmsVertY = his1D->GetRMS();
delete his1D;
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(2);
} meanVertZ = his1D->GetMean();
+ if(!his1D) return 1;
+
rmsVertZ = his1D->GetRMS();
delete his1D;
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(3);
}
+ if(!his1D) return 1;
+
meanMult = his1D->GetMean();
rmsMult = his1D->GetRMS();
delete his1D;
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(4);
}
+ if(!his1D) return 1;
+
meanMultPos = his1D->GetMean();
rmsMultPos = his1D->GetRMS();
delete his1D;
} else {
his1D = pTPC->GetTPCEventHisto()->Projection(5);
}
+ if(!his1D) return 1;
+
meanMultNeg = his1D->GetMean();
rmsMultNeg = his1D->GetRMS();
delete his1D;
if (pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_3_5_7")) {
his3D1 = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_pos_recvertex_3_5_7"));
+ if(!his3D1) return 256;
his3D1->GetYaxis()->SetRangeUser(0.1,0.8);
if (pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_3_5_7")) {
his3D2 = dynamic_cast<TH3*>(pTPC->GetHistos()->FindObject("h_tpc_track_neg_recvertex_3_5_7"));
+ if(!his3D2) return 256;
his3D2->GetYaxis()->SetRangeUser(0.1,0.8);
if (pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_5_8")) {
his2D = dynamic_cast<TH2*>(pTPC->GetHistos()->FindObject("h_tpc_track_all_recvertex_5_8"));
+ if(!his2D) return 512;
his1D1 = his2D->ProjectionX();
his1D1->Fit(fp1,"R");
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff