/************************************************************************** * 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. * **************************************************************************/ //------------------------------------------------------------------------------ // AlidNdPtAnalysisPbPb2011 class. // // a. functionality: // - fills analysis control histograms // - fills generic correction matrices // - generates correction matrices // // b. data members: // - generic correction matrices // - control histograms // // Author: J.Otwinowski 04/11/2008 // last change: 2011-04-04 by M.Knichel //------------------------------------------------------------------------------ #include "TH1.h" #include "TH2.h" #include "TCanvas.h" #include "THnSparse.h" #include "AliHeader.h" #include "AliGenEventHeader.h" #include "AliInputEventHandler.h" #include "AliAnalysisManager.h" #include "AliStack.h" #include "AliESDEvent.h" #include "AliMCEvent.h" #include "AliESDtrackCuts.h" #include "AliLog.h" #include "AliMultiplicity.h" #include "AliTracker.h" #include "AliCentrality.h" #include "AlidNdPtEventCuts.h" #include "AlidNdPtAcceptanceCuts.h" #include "AliPhysicsSelection.h" #include "AliTriggerAnalysis.h" #include "AliPWG0Helper.h" #include "AlidNdPtHelper.h" #include "AlidNdPtAnalysisPbPb2011.h" #include "TSystem.h" #include "TROOT.h" using namespace std; ClassImp(AlidNdPtAnalysisPbPb2011) //_____________________________________________________________________________ AlidNdPtAnalysisPbPb2011::AlidNdPtAnalysisPbPb2011(): AlidNdPt(), fAnalysisFolder(0), fHistogramsOn(kFALSE), // rec. track pt vs true track pt correlation matrix fTrackPtCorrelationMatrix(0), // event level correction fGenEventMatrix(0), fTriggerEventMatrix(0), fRecEventMatrix(0), // // track-event level correction // fGenTrackEventMatrix(0), fTriggerTrackEventMatrix(0), fRecTrackEventMatrix(0), // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix) fGenTrackMatrix(0), fGenPrimTrackMatrix(0), fRecPrimTrackMatrix(0), // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix) fRecTrackMatrix(0), fRecSecTrackMatrix(0), // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix) fRecMultTrackMatrix(0), // event control histograms fMCEventHist1(0), fRecEventHist1(0), fRecEventHist2(0), fRecMCEventHist1(0), fRecMCEventHist2(0), // rec. pt and eta resolution w.r.t MC fRecMCTrackHist1(0), //multple reconstructed tracks fMCMultRecTrackHist1(0), // rec. track control histograms fRecTrackHist3(0), fTriggerAnalysis(0), fCentralityEstimator(0), fMultNbins(0), fPtNbins(0), fPtCorrNbins(0), fEtaNbins(0), fZvNbins(0), fCentralityNbins(0), fBinsMult(0), fBinsPt(0), fBinsPtCorr(0), fBinsEta(0), fBinsZv(0), fBinsCentrality(0), fIsInit(kFALSE) { // default constructor for(Int_t i=0; iSetBinEdges(0,fBinsPt); fTrackPtCorrelationMatrix->SetBinEdges(1,fBinsPt); fTrackPtCorrelationMatrix->SetBinEdges(2,fBinsEta); fTrackPtCorrelationMatrix->SetBinEdges(3,fBinsCentrality); fTrackPtCorrelationMatrix->GetAxis(0)->SetTitle("Pt (GeV/c)"); fTrackPtCorrelationMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fTrackPtCorrelationMatrix->GetAxis(2)->SetTitle("mcEta"); fTrackPtCorrelationMatrix->GetAxis(3)->SetTitle("Centrality"); fTrackPtCorrelationMatrix->Sumw2(); // // Efficiency and contamination correction matrices // Int_t binsEventMatrix[3]={fZvNbins-1,fMultNbins-1,fCentralityNbins-1}; // Double_t minEventMatrix[3]={-30.,-0.5,0.}; // Double_t maxEventMatrix[3]={30.,10000.5,100.}; fGenEventMatrix = new THnSparseF("fGenEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix); fGenEventMatrix->SetBinEdges(0,fBinsZv); fGenEventMatrix->SetBinEdges(1,fBinsMult); fGenEventMatrix->SetBinEdges(2,fBinsCentrality); fGenEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fGenEventMatrix->GetAxis(1)->SetTitle("multiplicity MB"); fGenEventMatrix->GetAxis(2)->SetTitle("Centrality"); fGenEventMatrix->Sumw2(); // fTriggerEventMatrix = new THnSparseF("fTriggerEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix); fTriggerEventMatrix->SetBinEdges(0,fBinsZv); fTriggerEventMatrix->SetBinEdges(1,fBinsMult); fTriggerEventMatrix->SetBinEdges(2,fBinsCentrality); fTriggerEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fTriggerEventMatrix->GetAxis(1)->SetTitle("multiplicity MB"); fTriggerEventMatrix->GetAxis(2)->SetTitle("Centrality"); fTriggerEventMatrix->Sumw2(); // fRecEventMatrix = new THnSparseF("fRecEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix); //,minEventMatrix,maxEventMatrix); fRecEventMatrix->SetBinEdges(0,fBinsZv); fRecEventMatrix->SetBinEdges(1,fBinsMult); fRecEventMatrix->SetBinEdges(2,fBinsCentrality); fRecEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecEventMatrix->GetAxis(1)->SetTitle("multiplicity MB"); fRecEventMatrix->GetAxis(2)->SetTitle("Centrality"); fRecEventMatrix->Sumw2(); // // track to event corrections // fGenTrackEventMatrix = new THnSparseF("fGenTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fGenTrackEventMatrix->SetBinEdges(0,fBinsZv); fGenTrackEventMatrix->SetBinEdges(1,fBinsPtCorr); fGenTrackEventMatrix->SetBinEdges(2,fBinsEta); fGenTrackEventMatrix->SetBinEdges(3,fBinsCentrality); fGenTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fGenTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fGenTrackEventMatrix->GetAxis(2)->SetTitle("mcEta"); fGenTrackEventMatrix->GetAxis(3)->SetTitle("Centrality"); fGenTrackEventMatrix->Sumw2(); // fTriggerTrackEventMatrix = new THnSparseF("fTriggerTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fTriggerTrackEventMatrix->SetBinEdges(0,fBinsZv); fTriggerTrackEventMatrix->SetBinEdges(1,fBinsPtCorr); fTriggerTrackEventMatrix->SetBinEdges(2,fBinsEta); fTriggerTrackEventMatrix->SetBinEdges(3,fBinsCentrality); fTriggerTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fTriggerTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fTriggerTrackEventMatrix->GetAxis(2)->SetTitle("mcEta"); fTriggerTrackEventMatrix->GetAxis(3)->SetTitle("Centrality"); fTriggerTrackEventMatrix->Sumw2(); // fRecTrackEventMatrix = new THnSparseF("fRecTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fRecTrackEventMatrix->SetBinEdges(0,fBinsZv); fRecTrackEventMatrix->SetBinEdges(1,fBinsPtCorr); fRecTrackEventMatrix->SetBinEdges(2,fBinsEta); fRecTrackEventMatrix->SetBinEdges(3,fBinsCentrality); fRecTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fRecTrackEventMatrix->GetAxis(2)->SetTitle("mcEta"); fRecTrackEventMatrix->GetAxis(3)->SetTitle("Centrality"); fRecTrackEventMatrix->Sumw2(); // // tracks correction matrices // fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fGenTrackMatrix->SetBinEdges(0,fBinsZv); fGenTrackMatrix->SetBinEdges(1,fBinsPtCorr); fGenTrackMatrix->SetBinEdges(2,fBinsEta); fGenTrackMatrix->SetBinEdges(3,fBinsCentrality); fGenTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fGenTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fGenTrackMatrix->GetAxis(2)->SetTitle("mcEta"); fGenTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fGenTrackMatrix->Sumw2(); fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fGenPrimTrackMatrix->SetBinEdges(0,fBinsZv); fGenPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr); fGenPrimTrackMatrix->SetBinEdges(2,fBinsEta); fGenPrimTrackMatrix->SetBinEdges(3,fBinsCentrality); fGenPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fGenPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fGenPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta"); fGenPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fGenPrimTrackMatrix->Sumw2(); fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fRecPrimTrackMatrix->SetBinEdges(0,fBinsZv); fRecPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr); fRecPrimTrackMatrix->SetBinEdges(2,fBinsEta); fRecPrimTrackMatrix->SetBinEdges(3,fBinsCentrality); fRecPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fRecPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta"); fRecPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fRecPrimTrackMatrix->Sumw2(); // fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fRecTrackMatrix->SetBinEdges(0,fBinsZv); fRecTrackMatrix->SetBinEdges(1,fBinsPtCorr); fRecTrackMatrix->SetBinEdges(2,fBinsEta); fRecTrackMatrix->SetBinEdges(3,fBinsCentrality); fRecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fRecTrackMatrix->GetAxis(2)->SetTitle("mcEta"); fRecTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fRecTrackMatrix->Sumw2(); fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fRecSecTrackMatrix->SetBinEdges(0,fBinsZv); fRecSecTrackMatrix->SetBinEdges(1,fBinsPtCorr); fRecSecTrackMatrix->SetBinEdges(2,fBinsEta); fRecSecTrackMatrix->SetBinEdges(3,fBinsCentrality); fRecSecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecSecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)"); fRecSecTrackMatrix->GetAxis(2)->SetTitle("Eta"); fRecSecTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fRecSecTrackMatrix->Sumw2(); // fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix); fRecMultTrackMatrix->SetBinEdges(0,fBinsZv); fRecMultTrackMatrix->SetBinEdges(1,fBinsPtCorr); fRecMultTrackMatrix->SetBinEdges(2,fBinsEta); fRecMultTrackMatrix->SetBinEdges(3,fBinsCentrality); fRecMultTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)"); fRecMultTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)"); fRecMultTrackMatrix->GetAxis(2)->SetTitle("mcEta"); fRecMultTrackMatrix->GetAxis(3)->SetTitle("Centrality"); fRecMultTrackMatrix->Sumw2(); // // Control analysis histograms // Int_t binsMCEventHist1[4]={100,100,fZvNbins-1,fCentralityNbins-1}; Double_t minMCEventHist1[4]={-0.1,-0.1,-30.,0.}; Double_t maxMCEventHist1[4]={0.1,0.1,30.,100.}; fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv:Centrality",4,binsMCEventHist1,minMCEventHist1,maxMCEventHist1); fMCEventHist1->SetBinEdges(2,fBinsZv); fMCEventHist1->SetBinEdges(3,fBinsCentrality); fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)"); fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)"); fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)"); fMCEventHist1->GetAxis(3)->SetTitle("Centrality"); fMCEventHist1->Sumw2(); // Int_t binsRecEventHist1[4]={100,100,fZvNbins-1,fCentralityNbins-1}; Double_t minRecEventHist1[4]={-3.,-3.,-30.,0.}; Double_t maxRecEventHist1[4]={3.,3.,30.,100.}; fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv:Centrality",4,binsRecEventHist1,minRecEventHist1,maxRecEventHist1); fRecEventHist1->SetBinEdges(2,fBinsZv); fRecEventHist1->SetBinEdges(3,fBinsCentrality); fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)"); fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)"); fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)"); fRecEventHist1->GetAxis(3)->SetTitle("Centrality"); fRecEventHist1->Sumw2(); // Int_t binsRecEventHist2[3]={fZvNbins-1, fMultNbins-1, fCentralityNbins-1}; Double_t minRecEventHist2[3]={-30., -0.5, 0.}; Double_t maxRecEventHist2[3]={30., 10000.5, 100.}; fRecEventHist2 = new THnSparseF("fRecEventHist2","Zv:multMB:Centrality",3,binsRecEventHist2,minRecEventHist2,maxRecEventHist2); fRecEventHist2->SetBinEdges(0,fBinsZv); fRecEventHist2->SetBinEdges(1,fBinsMult); fRecEventHist2->SetBinEdges(2,fBinsCentrality); fRecEventHist2->GetAxis(0)->SetTitle("Zv (cm)"); fRecEventHist2->GetAxis(1)->SetTitle("multiplicity MB"); fRecEventHist2->GetAxis(2)->SetTitle("Centrality"); fRecEventHist2->Sumw2(); // Double_t kFact = 0.1; Int_t binsRecMCEventHist1[4]={100,100,100, fCentralityNbins-1}; Double_t minRecMCEventHist1[4]={-10.0*kFact,-10.0*kFact,-10.0*kFact, 0.}; Double_t maxRecMCEventHist1[4]={10.0*kFact,10.0*kFact,10.0*kFact, 100.}; fRecMCEventHist1 = new THnSparseF("fRecMCEventHist1","Xv-mcXv:Yv-mcYv:Zv-mcZv:Centrality",4,binsRecMCEventHist1,minRecMCEventHist1,maxRecMCEventHist1); fRecMCEventHist1->SetBinEdges(3,fBinsCentrality); fRecMCEventHist1->GetAxis(0)->SetTitle("Xv-mcXv (cm)"); fRecMCEventHist1->GetAxis(1)->SetTitle("Yv-mcYv (cm)"); fRecMCEventHist1->GetAxis(2)->SetTitle("Zv-mcZv (cm)"); fRecMCEventHist1->GetAxis(3)->SetTitle("Centrality"); fRecMCEventHist1->Sumw2(); // Int_t binsRecMCEventHist2[4]={100,100,fMultNbins-1, fCentralityNbins-1}; Double_t minRecMCEventHist2[4]={-10.0*kFact,-10.0*kFact,-0.5, 0.}; Double_t maxRecMCEventHist2[4]={10.0*kFact,10.0*kFact,10000.5, 100.}; fRecMCEventHist2 = new THnSparseF("fRecMCEventHist2","Xv-mcXv:Zv-mcZv:mult:Centrality",4,binsRecMCEventHist2,minRecMCEventHist2,maxRecMCEventHist2); fRecMCEventHist2->SetBinEdges(2,fBinsMult); fRecMCEventHist2->SetBinEdges(3,fBinsCentrality); fRecMCEventHist2->GetAxis(0)->SetTitle("Xv-mcXv (cm)"); fRecMCEventHist2->GetAxis(1)->SetTitle("Zv-mcZv (cm)"); fRecMCEventHist2->GetAxis(2)->SetTitle("multiplicity"); fRecMCEventHist2->GetAxis(3)->SetTitle("Centrality"); fRecMCEventHist2->Sumw2(); // char name[256]; char title[256]; for(Int_t i=0; iSetBinEdges(0,fBinsPt); fMCTrackHist1[i]->SetBinEdges(1,fBinsEta); fMCTrackHist1[i]->SetBinEdges(3,fBinsCentrality); fMCTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)"); fMCTrackHist1[i]->GetAxis(1)->SetTitle("mcEta"); fMCTrackHist1[i]->GetAxis(2)->SetTitle("mcPhi (rad)"); fMCTrackHist1[i]->GetAxis(3)->SetTitle("Centrality"); fMCTrackHist1[i]->Sumw2(); Int_t binsMCPrimTrackHist1[6]= {fPtNbins-1,fEtaNbins-1,6,20,4000, fCentralityNbins-1}; Double_t minMCPrimTrackHist1[6]={0.,-1.5,0.,0.,0., 0.}; Double_t maxMCPrimTrackHist1[6]={50.,1.5,6.,20.,4000., 100.}; snprintf(name,256,"fMCPrimTrackHist1_%d",i); snprintf(title,256,"mcPt:mcEta:pid:mech:mother:Centrality"); fMCPrimTrackHist1[i] = new THnSparseF(name,title,6,binsMCPrimTrackHist1,minMCPrimTrackHist1,maxMCPrimTrackHist1); fMCPrimTrackHist1[i]->SetBinEdges(0,fBinsPt); fMCPrimTrackHist1[i]->SetBinEdges(1,fBinsEta); fMCPrimTrackHist1[i]->SetBinEdges(5,fBinsCentrality); fMCPrimTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)"); fMCPrimTrackHist1[i]->GetAxis(1)->SetTitle("mcEta"); fMCPrimTrackHist1[i]->GetAxis(2)->SetTitle("pid"); fMCPrimTrackHist1[i]->GetAxis(3)->SetTitle("mech"); fMCPrimTrackHist1[i]->GetAxis(4)->SetTitle("mother"); fMCPrimTrackHist1[i]->GetAxis(5)->SetTitle("Centrality"); fMCPrimTrackHist1[i]->Sumw2(); Int_t binsMCPrimTrackHist2[4]= {4000,20,4000,fCentralityNbins-1}; Double_t minMCPrimTrackHist2[4]={0.,0.,0., 0.}; Double_t maxMCPrimTrackHist2[4]={4000.,20.,4000., 100.}; snprintf(name,256,"fMCPrimTrackHist2_%d",i); snprintf(title,256,"pdg:mech:mother:Centrality"); fMCPrimTrackHist2[i] = new THnSparseF(name,title,4,binsMCPrimTrackHist2,minMCPrimTrackHist2,maxMCPrimTrackHist2); fMCPrimTrackHist2[i]->SetBinEdges(3,fBinsCentrality); fMCPrimTrackHist2[i]->GetAxis(0)->SetTitle("pdg"); fMCPrimTrackHist2[i]->GetAxis(1)->SetTitle("mech"); fMCPrimTrackHist2[i]->GetAxis(2)->SetTitle("mother"); fMCPrimTrackHist2[i]->GetAxis(3)->SetTitle("Centrality"); fMCPrimTrackHist2[i]->Sumw2(); Int_t binsMCSecTrackHist1[6]= {fPtNbins-1,fEtaNbins-1,6,20,4000, fCentralityNbins-1}; Double_t minMCSecTrackHist1[6]={0.,-1.5,0.,0.,0., 0.}; Double_t maxMCSecTrackHist1[6]={50.,1.5,6.,20.,4000., 100.}; snprintf(name,256,"fMCSecTrackHist1_%d",i); snprintf(title,256,"mcPt:mcEta:pid:mech:mother:Centrality"); fMCSecTrackHist1[i] = new THnSparseF(name,title,6,binsMCSecTrackHist1,minMCSecTrackHist1,maxMCSecTrackHist1); fMCSecTrackHist1[i]->SetBinEdges(0,fBinsPt); fMCSecTrackHist1[i]->SetBinEdges(1,fBinsEta); fMCSecTrackHist1[i]->SetBinEdges(5,fBinsCentrality); fMCSecTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)"); fMCSecTrackHist1[i]->GetAxis(1)->SetTitle("mcEta"); fMCSecTrackHist1[i]->GetAxis(2)->SetTitle("pid"); fMCSecTrackHist1[i]->GetAxis(3)->SetTitle("mech"); fMCSecTrackHist1[i]->GetAxis(4)->SetTitle("mother"); fMCSecTrackHist1[i]->GetAxis(5)->SetTitle("Centrality"); fMCSecTrackHist1[i]->Sumw2(); Int_t binsRecTrackHist1[4]={fPtNbins-1,fEtaNbins-1,90, fCentralityNbins-1}; Double_t minRecTrackHist1[4]={0.,-1.5,0., 0.}; Double_t maxRecTrackHist1[4]={50.,1.5,2.*TMath::Pi(), 100.}; snprintf(name,256,"fRecTrackHist1_%d",i); snprintf(title,256,"Pt:Eta:Phi:Centrality"); fRecTrackHist1[i] = new THnSparseF(name,title,4,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1); fRecTrackHist1[i]->SetBinEdges(0,fBinsPt); fRecTrackHist1[i]->SetBinEdges(1,fBinsEta); fRecTrackHist1[i]->SetBinEdges(3,fBinsCentrality); fRecTrackHist1[i]->GetAxis(0)->SetTitle("p_{T} (GeV/c)"); fRecTrackHist1[i]->GetAxis(1)->SetTitle("#eta"); fRecTrackHist1[i]->GetAxis(2)->SetTitle("#phi (rad)"); fRecTrackHist1[i]->GetAxis(3)->SetTitle("Centrality"); fRecTrackHist1[i]->Sumw2(); snprintf(name,256,"fRecTrackHist2_%d",i); snprintf(title,256,"Zv:Pt:Eta:Centrality"); fRecTrackHist2[i] = new THnSparseF(name,title,4,binsTrackEvent); fRecTrackHist2[i]->SetBinEdges(0,fBinsZv); fRecTrackHist2[i]->SetBinEdges(1,fBinsPt); fRecTrackHist2[i]->SetBinEdges(2,fBinsEta); fRecTrackHist2[i]->SetBinEdges(3,fBinsCentrality); fRecTrackHist2[i]->GetAxis(0)->SetTitle("Zv (cm)"); fRecTrackHist2[i]->GetAxis(1)->SetTitle("p_{T} (GeV/c)"); fRecTrackHist2[i]->GetAxis(2)->SetTitle("#eta"); fRecTrackHist2[i]->GetAxis(3)->SetTitle("Centrality"); fRecTrackHist2[i]->Sumw2(); // Int_t binsRecTrackMultHist1[3]={fPtNbins-1,fMultNbins-1, fCentralityNbins-1}; Double_t minRecTrackMultHist1[3]={0.,-0.5, -0.}; Double_t maxRecTrackMultHist1[3]={50.,10000.5, 100.}; snprintf(name,256,"fRecTrackMultHist_%d",i); snprintf(title,256,"Pt:Mult:Centrality"); fRecTrackMultHist1[i] = new THnSparseF(name,title,3,binsRecTrackMultHist1,minRecTrackMultHist1,maxRecTrackMultHist1); fRecTrackMultHist1[i]->SetBinEdges(0,fBinsPt); fRecTrackMultHist1[i]->SetBinEdges(1,fBinsMult); fRecTrackMultHist1[i]->SetBinEdges(2,fBinsCentrality); fRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)"); fRecTrackMultHist1[i]->GetAxis(1)->SetTitle("multiplicity"); fRecTrackMultHist1[i]->GetAxis(2)->SetTitle("Centrality"); fRecTrackMultHist1[i]->Sumw2(); } Int_t binsRecMCTrackHist1[5] = {fPtNbins-1,fEtaNbins-1,100,100, fCentralityNbins-1}; Double_t minRecMCTrackHist1[5]={0.,-1.5,-0.5,-0.5, 0.}; Double_t maxRecMCTrackHist1[5]={50.,1.5,0.5,0.5, 100.}; snprintf(name,256,"fRecMCTrackHist1"); snprintf(title,256,"mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta):Centrality"); fRecMCTrackHist1 = new THnSparseF(name,title,5,binsRecMCTrackHist1,minRecMCTrackHist1,maxRecMCTrackHist1); fRecMCTrackHist1->SetBinEdges(0,fBinsPt); fRecMCTrackHist1->SetBinEdges(1,fBinsEta); fRecMCTrackHist1->SetBinEdges(4,fBinsCentrality); fRecMCTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)"); fRecMCTrackHist1->GetAxis(1)->SetTitle("mcEta"); fRecMCTrackHist1->GetAxis(2)->SetTitle("(Pt-mcPt)/mcPt"); fRecMCTrackHist1->GetAxis(3)->SetTitle("Eta-mcEta"); fRecMCTrackHist1->GetAxis(4)->SetTitle("Centrality"); Int_t binsMCMultRecTrackHist1[4] = {fPtNbins-1,fEtaNbins-1,6, fCentralityNbins-1}; Double_t minMCMultRecTrackHist1[4]={0.,-1.5,0., 0.}; Double_t maxMCMultRecTrackHist1[4]={50.,1.5,6., 100.}; snprintf(name,256,"fMCMultRecTrackHist1"); snprintf(title,256,"mcPt:mcEta:pid:Centrality"); fMCMultRecTrackHist1 = new THnSparseF(name,title,4,binsMCMultRecTrackHist1,minMCMultRecTrackHist1,maxMCMultRecTrackHist1); fMCMultRecTrackHist1->SetBinEdges(0,fBinsPt); fMCMultRecTrackHist1->SetBinEdges(1,fBinsEta); fMCMultRecTrackHist1->SetBinEdges(3,fBinsCentrality); fMCMultRecTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)"); fMCMultRecTrackHist1->GetAxis(1)->SetTitle("mcEta"); fMCMultRecTrackHist1->GetAxis(2)->SetTitle("pid"); fMCMultRecTrackHist1->GetAxis(3)->SetTitle("Centrality"); //nClust:chi2PerClust:pt:eta:phi:Centrality Int_t binsRecTrackHist3[6]={160,100,fPtNbins-1,fEtaNbins-1,90, fCentralityNbins-1}; Double_t minRecTrackHist3[6]={0., 0., 0., -1.5, 0., 0.}; Double_t maxRecRecTrackHist3[6]={160.,10., 50., 1.5, 2.*TMath::Pi(), 100.}; fRecTrackHist3 = new THnSparseF("fRecTrackHist3","nClust:chi2PerClust:pt:eta:phi:Centrality",6,binsRecTrackHist3,minRecTrackHist3,maxRecRecTrackHist3); fRecTrackHist3->SetBinEdges(2,fBinsPt); fRecTrackHist3->SetBinEdges(3,fBinsEta); fRecTrackHist3->SetBinEdges(5,fBinsCentrality); fRecTrackHist3->GetAxis(0)->SetTitle("nClust"); fRecTrackHist3->GetAxis(1)->SetTitle("chi2PerClust"); fRecTrackHist3->GetAxis(2)->SetTitle("p_{T} (GeV/c)"); fRecTrackHist3->GetAxis(3)->SetTitle("#eta"); fRecTrackHist3->GetAxis(4)->SetTitle("#phi (rad)"); fRecTrackHist3->GetAxis(5)->SetTitle("Centrality"); fRecTrackHist3->Sumw2(); // init folder fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder"); // init trigger analysis (for zdc cut) fTriggerAnalysis = new AliTriggerAnalysis; // set init flag fIsInit = kTRUE; } //_____________________________________________________________________________ void AlidNdPtAnalysisPbPb2011::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent) { // init if not done already if (!fIsInit) { Init(); } // // Process real and/or simulated events // if(!esdEvent) { AliDebug(AliLog::kError, "esdEvent not available"); return; } // get selection cuts AlidNdPtEventCuts *evtCuts = GetEventCuts(); AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); AlidNdPtAcceptanceCuts *recCuts = GetRecAcceptanceCuts(); AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); if(!evtCuts || !accCuts || !esdTrackCuts) { AliDebug(AliLog::kError, "cuts not available"); return; } if (0 == recCuts) { recCuts = accCuts;} // trigger selection Bool_t isEventTriggered = kTRUE; AliPhysicsSelection *physicsSelection = NULL; AliTriggerAnalysis* triggerAnalysis = NULL; // AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler(); if (!inputHandler) { Printf("ERROR: Could not receive input handler"); return; } if(evtCuts->IsTriggerRequired()) { // always MB //isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB; isEventTriggered = inputHandler->IsEventSelected() & GetTriggerMask(); physicsSelection = static_cast (inputHandler->GetEventSelection()); if(!physicsSelection) return; //SetPhysicsTriggerSelection(physicsSelection); if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) { // set trigger (V0AND) triggerAnalysis = physicsSelection->GetTriggerAnalysis(); if(!triggerAnalysis) return; isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger()); } } // centrality determination Float_t centralityF = -2.; AliCentrality *esdCentrality = esdEvent->GetCentrality(); centralityF = esdCentrality->GetCentralityPercentile((const char*)GetCentralityEstimator()); // if centrality percentile is larger then the highest bin, return if(centralityF > fBinsCentrality[fCentralityNbins-1]) { AliDebug(AliLog::kError, Form("Event out of centrality binning (event: %.2f, maxallowed: %.2f)",centralityF, fBinsCentrality[fCentralityNbins-1])); return; } // use MC information AliHeader* header = 0; AliGenEventHeader* genHeader = 0; AliStack* stack = 0; TArrayF vtxMC(3); Int_t multMCTrueTracks = 0; if(IsUseMCInfo()) { // if(!mcEvent) { AliDebug(AliLog::kError, "mcEvent not available"); return; } // get MC event header header = mcEvent->Header(); if (!header) { AliDebug(AliLog::kError, "Header not available"); return; } // MC particle stack stack = mcEvent->Stack(); if (!stack) { AliDebug(AliLog::kError, "Stack not available"); return; } // get MC vertex genHeader = header->GenEventHeader(); if (!genHeader) { AliDebug(AliLog::kError, "Could not retrieve genHeader from Header"); return; } genHeader->PrimaryVertex(vtxMC); Double_t vMCEventHist1[4]={vtxMC[0],vtxMC[1],vtxMC[2],centralityF}; fMCEventHist1->Fill(vMCEventHist1); // multipliticy of all MC primary tracks // in Zv, pt and eta ranges) multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts); } // end bUseMC // get reconstructed vertex const AliESDVertex* vtxESD = 0; if(GetAnalysisMode() == AlidNdPtHelper::kTPC) { vtxESD = esdEvent->GetPrimaryVertexTPC(); } else if(GetAnalysisMode() == AlidNdPtHelper::kTPCITS) { vtxESD = esdEvent->GetPrimaryVertexTracks(); } else { return; } if(!vtxESD) return; Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD); //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered); //printf("GetAnalysisMode() %d \n",GetAnalysisMode()); // vertex contributors Int_t multMBTracks = 0; if(GetAnalysisMode() == AlidNdPtHelper::kTPC || GetAnalysisMode() == AlidNdPtHelper::kTPCITS) { if(vtxESD->GetStatus()) { multMBTracks = vtxESD->GetNContributors(); } } else { AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode())); return; } TObjArray *allChargedTracks=0; //Int_t multAll=0, multAcc=0, multRec=0; Int_t multAll=0, multRec=0; Int_t *labelsAll=0, *labelsAcc=0, *labelsRec=0; // check event cuts if(isEventOK && isEventTriggered) { // get all charged tracks allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode()); if(!allChargedTracks) return; Int_t entries = allChargedTracks->GetEntries(); labelsAll = new Int_t[entries]; labelsAcc = new Int_t[entries]; labelsRec = new Int_t[entries]; for(Int_t i=0; iAt(i); if(!track) continue; if(track->Charge()==0) continue; if(IsUseMCInfo() && ( mcEvent->IsFromBGEvent(track->GetLabel()) ) && (track->GetLabel() >= 0) ) { // Printf("i = %i, label = %i", i, track->GetLabel()); continue; } // only postive charged if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) continue; // only negative charged if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) continue; // Double_t values[4] = {vtxESD->GetZv(),track->Pt(),track->Eta(), centralityF}; fRecTrackHist2[AlidNdPtHelper::kAllTracks]->Fill(values); FillHistograms(track,stack,AlidNdPtHelper::kAllTracks,centralityF); labelsAll[multAll] = TMath::Abs(track->GetLabel()); multAll++; if(esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track) && recCuts->AcceptTrackLocalTPC(track)) { fRecTrackHist2[AlidNdPtHelper::kRecTracks]->Fill(values); FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,centralityF); labelsRec[multRec] = TMath::Abs(track->GetLabel()); multRec++; } }//loop over entries // fill track multiplicity histograms //FillHistograms(allChargedTracks,labelsAll,multAll,labelsAcc,multAcc,labelsRec,multRec); Double_t vRecEventHist1[4] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),centralityF}; fRecEventHist1->Fill(vRecEventHist1); Double_t vRecEventHist2[3] = {vtxESD->GetZv(),Double_t(multMBTracks),centralityF}; fRecEventHist2->Fill(vRecEventHist2); } // triggered and event vertex if(IsUseMCInfo()) { // // event level corrections (zv,N_MB) // // all inelastic Double_t vEventMatrix[3] = {vtxMC[2],Double_t(multMBTracks),centralityF}; fGenEventMatrix->Fill(vEventMatrix); if(isEventTriggered) fTriggerEventMatrix->Fill(vEventMatrix); if(isEventOK && isEventTriggered) fRecEventMatrix->Fill(vEventMatrix); // // track-event level corrections (zv,pt,eta) // for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc) { TParticle* particle = stack->Particle(iMc); if (!particle) continue; // only charged particles if(!particle->GetPDG()) continue; Double_t charge = particle->GetPDG()->Charge()/3.; if ( TMath::Abs(charge) < 0.001 ) continue; // only postive charged if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.) continue; // only negative charged if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) continue; // physical primary Bool_t prim = stack->IsPhysicalPrimary(iMc); if(!prim) continue; // injected signal Bool_t isInjected = mcEvent->IsFromBGEvent(iMc); if(isInjected) continue; // checked accepted if(accCuts->AcceptTrack(particle)) { Double_t vTrackEventMatrix[4] = {vtxMC[2], particle->Pt(), particle->Eta(), centralityF}; fGenTrackEventMatrix->Fill(vTrackEventMatrix); if(!isEventTriggered) continue; fTriggerTrackEventMatrix->Fill(vTrackEventMatrix); if(!isEventOK) continue; fRecTrackEventMatrix->Fill(vTrackEventMatrix); }// if(accCuts->AcceptTrack(particle)) }// for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc) // // track-level corrections (zv,pt,eta) // if(isEventOK && isEventTriggered) { // fill MC and rec event control histograms if(fHistogramsOn) { Double_t vRecMCEventHist1[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2], centralityF}; fRecMCEventHist1->Fill(vRecMCEventHist1);// Double_t vRecMCEventHist2[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetZv()-vtxMC[2],Double_t(multMBTracks), centralityF}; fRecMCEventHist2->Fill(vRecMCEventHist2); }// // // MC histograms for track efficiency studies // for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc) { TParticle* particle = stack->Particle(iMc); if (!particle) continue; Double_t vTrackMatrix[4] = {vtxMC[2],particle->Pt(),particle->Eta(),centralityF}; // only charged particles if(!particle->GetPDG()) continue; Double_t charge = particle->GetPDG()->Charge()/3.; if (TMath::Abs(charge) < 0.001) continue; // only postive charged if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.) continue; // only negative charged if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) continue; // physical primary Bool_t prim = stack->IsPhysicalPrimary(iMc); // injected signal Bool_t isInjected = mcEvent->IsFromBGEvent(iMc); if(isInjected) continue; // check accepted if(accCuts->AcceptTrack(particle)) { if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) fGenPrimTrackMatrix->Fill(vTrackMatrix); // fill control histograms if(fHistogramsOn) FillHistograms(stack,iMc,AlidNdPtHelper::kAccTracks, centralityF); // check multiple found tracks Int_t multCount = 0; for(Int_t iRec=0; iRec1) { fRecMultTrackMatrix->Fill(vTrackMatrix); // fill control histogram if(fHistogramsOn) { Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle); Double_t vMCMultRecTrackHist1[4] = {particle->Pt(), particle->Eta(), Double_t(pid), centralityF}; fMCMultRecTrackHist1->Fill(vMCMultRecTrackHist1); } } } } // check reconstructed for(Int_t iRec=0; iRecFill(vTrackMatrix); if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) { fRecPrimTrackMatrix->Fill(vTrackMatrix); } if(!prim) fRecSecTrackMatrix->Fill(vTrackMatrix); // fill control histograms if(fHistogramsOn) FillHistograms(stack,iMc,AlidNdPtHelper::kRecTracks, centralityF); break; }//if(iMc == labelsRec[iRec]) }//reco tracks }//accepted tracks }//stack loop }//is triggered } // end bUseMC if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0; if(labelsAll) delete [] labelsAll; labelsAll = 0; if(labelsAcc) delete [] labelsAcc; labelsAcc = 0; if(labelsRec) delete [] labelsRec; labelsRec = 0; } //_____________________________________________________________________________ void AlidNdPtAnalysisPbPb2011::FillHistograms(TObjArray *const allChargedTracks,Int_t *const labelsAll,Int_t multAll,Int_t *const labelsAcc,Int_t multAcc,Int_t *const labelsRec,Int_t multRec, Float_t centralityF) { // multiplicity histograms if(!allChargedTracks) return; if(!labelsAll) return; if(!labelsAcc) return; if(!labelsRec) return; Int_t entries = allChargedTracks->GetEntries(); for(Int_t i=0; iAt(i); if(!track) continue; if(track->Charge() == 0) continue; Int_t label = TMath::Abs(track->GetLabel()); for(Int_t iAll=0; iAllPt(), Double_t(multAll), centralityF}; fRecTrackMultHist1[AlidNdPtHelper::kAllTracks]->Fill(v1); } } for(Int_t iAcc=0; iAccPt(), Double_t(multAcc), centralityF}; fRecTrackMultHist1[AlidNdPtHelper::kAccTracks]->Fill(v2); } } for(Int_t iRec=0; iRecPt(), Double_t(multRec), centralityF}; fRecTrackMultHist1[AlidNdPtHelper::kRecTracks]->Fill(v3); }//out } } } //_____________________________________________________________________________ void AlidNdPtAnalysisPbPb2011::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF) { // // Fill ESD track and MC histograms // if(!esdTrack) return; Float_t q = esdTrack->Charge(); if(TMath::Abs(q) < 0.001) return; Float_t pt = esdTrack->Pt(); Float_t eta = esdTrack->Eta(); Float_t phi = esdTrack->Phi(); Float_t dca[2], bCov[3]; esdTrack->GetImpactParameters(dca,bCov); Int_t nClust = esdTrack->GetTPCclusters(0); Float_t chi2PerCluster = 0.; if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust); // fill histograms Double_t values[4] = {pt,eta,phi,centralityF}; fRecTrackHist1[trackObj]->Fill(values); // // Fill rec vs MC information // if(!stack) return; Int_t label = TMath::Abs(esdTrack->GetLabel()); //if(label == 0) return; if(label > stack->GetNtrack()) return; TParticle* particle = stack->Particle(label); if(!particle) return; Int_t motherPdg = -1; TParticle* mother = 0; //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label); Int_t motherLabel = particle->GetMother(0); if(motherLabel>0) mother = stack->Particle(motherLabel); if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only //Int_t mech = particle->GetUniqueID(); // production mechanism if(!particle->GetPDG()) return; Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 if(TMath::Abs(gq)<0.001) return; Float_t gpt = particle->Pt(); Float_t geta = particle->Eta(); Double_t dpt=0; //printf("pt %f, gpt %f \n",pt,gpt); if(gpt) dpt = (pt-gpt)/gpt; Double_t deta = (eta-geta); // fill histograms if(trackObj == AlidNdPtHelper::kRecTracks) //RecTracks??? { Double_t vTrackPtCorrelationMatrix[4]={pt,gpt,geta,centralityF}; fTrackPtCorrelationMatrix->Fill(vTrackPtCorrelationMatrix); Double_t vRecMCTrackHist1[5]={gpt,geta,dpt,deta,centralityF}; fRecMCTrackHist1->Fill(vRecMCTrackHist1); } } //_____________________________________________________________________________ void AlidNdPtAnalysisPbPb2011::FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF) { // Fill MC histograms if(!stack) return; if(label > stack->GetNtrack()) return; TParticle* particle = stack->Particle(label); if(!particle) return; Int_t motherPdg = -1; TParticle* mother = 0; //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label); Int_t motherLabel = particle->GetMother(0); if(motherLabel>0) mother = stack->Particle(motherLabel); if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only Int_t mech = particle->GetUniqueID(); // production mechanism if(!particle->GetPDG()) return; Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 if(TMath::Abs(gq) < 0.001) return; Float_t gpt = particle->Pt(); //Float_t qgpt = particle->Pt() * gq; Float_t geta = particle->Eta(); Float_t gphi = particle->Phi(); //Float_t gpz = particle->Pz(); Bool_t prim = stack->IsPhysicalPrimary(label); //Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz(); Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle); //if(prim&&pid==5) printf("pdgcode %d, production mech %d \n",particle->GetPdgCode(),mech); //if(!prim) printf("motherPdg %d, particle %d, production mech %d\n",motherPdg, particle->GetPdgCode(),mech); // // fill histogram // Double_t vMCTrackHist1[4] = {gpt,geta,gphi,centralityF}; fMCTrackHist1[trackObj]->Fill(vMCTrackHist1); Double_t vMCPrimTrackHist1[6] = {gpt,geta,Double_t(pid),Double_t(mech),Double_t(motherPdg),centralityF}; Double_t vMCPrimTrackHist2[4] = {Double_t(TMath::Abs(particle->GetPdgCode())),Double_t(mech),Double_t(motherPdg),centralityF}; //if(prim && AliPWG0Helper::IsPrimaryCharged(particle, label)) fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1); if(prim) { fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1); if(pid == 5) fMCPrimTrackHist2[trackObj]->Fill(vMCPrimTrackHist2); } else { fMCSecTrackHist1[trackObj]->Fill(vMCPrimTrackHist1); } } //_____________________________________________________________________________ Long64_t AlidNdPtAnalysisPbPb2011::Merge(TCollection* const list) { // Merge list of objects (needed by PROOF) // init if not done already if (!fIsInit) { Init(); } if (!list) return 0; if (list->IsEmpty()) return 1; TIterator* iter = list->MakeIterator(); TObject* obj = 0; // //TList *collPhysSelection = new TList; // collection of generated histograms Int_t count=0; while((obj = iter->Next()) != 0) { AlidNdPtAnalysisPbPb2011* entry = dynamic_cast(obj); if (entry == 0) continue; // physics selection //printf("entry->GetPhysicsTriggerSelection() %p \n", entry->GetPhysicsTriggerSelection()); //AliPhysicsSelection *physSel = entry->GetPhysicsTriggerSelection(); //if( physSel ) collPhysSelection->Add(physSel); // fTrackPtCorrelationMatrix->Add(entry->fTrackPtCorrelationMatrix); // fGenEventMatrix->Add(entry->fGenEventMatrix); fTriggerEventMatrix->Add(entry->fTriggerEventMatrix); fRecEventMatrix->Add(entry->fRecEventMatrix); // fGenTrackEventMatrix->Add(entry->fGenTrackEventMatrix); fTriggerTrackEventMatrix->Add(entry->fTriggerTrackEventMatrix); fRecTrackEventMatrix->Add(entry->fRecTrackEventMatrix); // fGenTrackMatrix->Add(entry->fGenTrackMatrix); fGenPrimTrackMatrix->Add(entry->fGenPrimTrackMatrix); fRecPrimTrackMatrix->Add(entry->fRecPrimTrackMatrix); // fRecTrackMatrix->Add(entry->fRecTrackMatrix); fRecSecTrackMatrix->Add(entry->fRecSecTrackMatrix); // fRecMultTrackMatrix->Add(entry->fRecMultTrackMatrix); // // control analysis histograms // fMCEventHist1->Add(entry->fMCEventHist1); fRecEventHist1->Add(entry->fRecEventHist1); fRecEventHist2->Add(entry->fRecEventHist2); fRecMCEventHist1->Add(entry->fRecMCEventHist1); fRecMCEventHist2->Add(entry->fRecMCEventHist2); for(Int_t i=0; iAdd(entry->fMCTrackHist1[i]); fMCPrimTrackHist1[i]->Add(entry->fMCPrimTrackHist1[i]); fMCPrimTrackHist2[i]->Add(entry->fMCPrimTrackHist2[i]); fMCSecTrackHist1[i]->Add(entry->fMCSecTrackHist1[i]); fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]); fRecTrackHist2[i]->Add(entry->fRecTrackHist2[i]); fRecTrackMultHist1[i]->Add(entry->fRecTrackMultHist1[i]); } fRecMCTrackHist1->Add(entry->fRecMCTrackHist1); fMCMultRecTrackHist1->Add(entry->fMCMultRecTrackHist1); fRecTrackHist3->Add(entry->fRecTrackHist3); count++; } //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection(); //if( trigSelection ) trigSelection->Merge(collPhysSelection); //if(collPhysSelection) delete collPhysSelection; return count; } //_____________________________________________________________________________ void AlidNdPtAnalysisPbPb2011::Analyse() { // init if not done already if (!fIsInit) { Init(); } // Analyse histograms // TH1::AddDirectory(kFALSE); TH1 *h=0, *h1=0, *h2=0, *h2c = 0; THnSparse *hs=0; TH2 *h2D=0; char name[256]; TObjArray *aFolderObj = new TObjArray; if(!aFolderObj) return; // // LHC backgraund in all and 0-bins // //AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection(); //trigSel->SaveHistograms("physics_selection"); // // Reconstructed event vertex // h = fRecEventHist1->Projection(0); h->SetName("Xv"); aFolderObj->Add(h); h = fRecEventHist1->Projection(1); h->SetName("Yv"); aFolderObj->Add(h); h = fRecEventHist1->Projection(2); h->SetName("Zv"); aFolderObj->Add(h); // // multiplicity // h = fRecEventHist2->Projection(1); h->SetName("multMB"); aFolderObj->Add(h); h2D = fRecEventHist2->Projection(0,1); h2D->SetName("Zv_vs_multiplicity_MB"); aFolderObj->Add(h2D); // // reconstructed pt histograms // h = fRecTrackHist1[0]->Projection(0); h->Scale(1.,"width"); h->SetName("pt_all_ch"); aFolderObj->Add(h); h = fRecTrackHist1[1]->Projection(0); h->Scale(1.,"width"); h->SetName("pt_acc"); aFolderObj->Add(h); h = fRecTrackHist1[2]->Projection(0); h->Scale(1.,"width"); h->SetName("pt_rec"); aFolderObj->Add(h); // // reconstructed eta histograms // h = fRecTrackHist1[0]->Projection(1); h->SetName("eta_all_ch"); aFolderObj->Add(h); h = fRecTrackHist1[1]->Projection(1); h->SetName("eta_acc"); aFolderObj->Add(h); h = fRecTrackHist1[2]->Projection(1); h->SetName("eta_rec"); aFolderObj->Add(h); // // reconstructed phi histograms // h = fRecTrackHist1[0]->Projection(2); h->SetName("phi_all_ch"); aFolderObj->Add(h); h = fRecTrackHist1[1]->Projection(2); h->SetName("phi_acc"); aFolderObj->Add(h); h = fRecTrackHist1[2]->Projection(2); h->SetName("phi_rec"); aFolderObj->Add(h); // // reconstructed eta:pt histograms // h2D = fRecTrackHist1[0]->Projection(1,0); h2D->SetName("pt_eta_all_ch"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[1]->Projection(1,0); h2D->SetName("pt_eta_acc"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[2]->Projection(1,0); h2D->SetName("pt_eta_rec"); aFolderObj->Add(h2D); // // reconstructed phi:pt histograms // h2D = fRecTrackHist1[0]->Projection(2,0); h2D->SetName("pt_phi_all_ch"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[1]->Projection(2,0); h2D->SetName("pt_phi_acc"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[2]->Projection(2,0); h2D->SetName("pt_phi_rec"); aFolderObj->Add(h2D); // // reconstructed phi:eta histograms // h2D = fRecTrackHist1[0]->Projection(2,1); h2D->SetName("eta_phi_all_ch"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[1]->Projection(2,1); h2D->SetName("eta_phi_acc"); aFolderObj->Add(h2D); h2D = fRecTrackHist1[2]->Projection(2,1); h2D->SetName("eta_phi_rec"); aFolderObj->Add(h2D); // // reconstructed nClust, chi2 vs pt, eta, phi // if(fHistogramsOn) { h2D = fRecTrackHist3->Projection(0,1); h2D->SetName("nClust_chi2_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(0,2); h2D->SetName("nClust_pt_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(0,3); h2D->SetName("nClust_eta_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(0,4); h2D->SetName("nClust_phi_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(1,2); h2D->SetName("chi2_pt_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(1,3); h2D->SetName("chi2_eta_rec"); aFolderObj->Add(h2D); h2D = fRecTrackHist3->Projection(1,4); h2D->SetName("chi2_phi_rec"); aFolderObj->Add(h2D); } // // calculate corrections for empty events // with multMB==0 // // // normalised zv to generate zv for triggered events // h = fRecEventHist2->Projection(0); if( h->Integral() ) h->Scale(1./h->Integral()); h->SetName("zv_distribution_norm"); aFolderObj->Add(h); // // MC available // if(IsUseMCInfo()) { // // Event vertex resolution // h2D = fRecMCEventHist2->Projection(0,2); h2D->SetName("DeltaXv_vs_mult"); aFolderObj->Add(h2D); h2D = fRecMCEventHist2->Projection(1,2); h2D->SetName("DeltaZv_vs_mult"); aFolderObj->Add(h2D); // // normalised zv to get trigger/trigger+vertex event differences // F(zv) = E_trig(zv,0)/Int(E_trig(zv,0) / Sum(E_trigvtx(zv,n))/Sum(Int(E_trigvtx(zv,n))dzv) // fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0.,0.); h = fTriggerEventMatrix->Projection(0); h2D = fTriggerEventMatrix->Projection(0,1); if(h2D->Integral()) h->Scale(1./h2D->Integral()); h1 = fRecEventMatrix->Projection(0); h2D = fRecEventMatrix->Projection(0,1); if(h2D->Integral()) h1->Scale(1./h2D->Integral()); h->Divide(h1); h->SetName("zv_empty_events_norm"); aFolderObj->Add(h); fTriggerEventMatrix->GetAxis(1)->SetRange(1,fTriggerEventMatrix->GetAxis(1)->GetNbins()); // // rec. vs true track pt correlation matrix // hs = (THnSparse*)fTrackPtCorrelationMatrix->Clone("track_pt_correlation_matrix"); aFolderObj->Add(hs); // // trigger efficiency for INEL // h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fGenEventMatrix->Projection(0),"zv_trig_INEL_eff_matrix"); aFolderObj->Add(h); // // trigger bias correction (MB to INEL) // hs = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoInel_corr_matrix"); aFolderObj->Add(hs); h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoInel_corr_matrix"); aFolderObj->Add(h); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoInel_corr_matrix_2D"); aFolderObj->Add(h2D); h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoInel_corr_matrix"); aFolderObj->Add(h); // // event vertex reconstruction correction (MB) // hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix,fRecEventMatrix,"zv_mult_event_corr_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0,1),fRecEventMatrix->Projection(0,1),"zv_mult_event_corr_matrix_2D"); aFolderObj->Add(h2D); h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(1),fRecEventMatrix->Projection(1),"mult_event_corr_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fRecEventMatrix->Projection(0),"zv_event_corr_matrix"); aFolderObj->Add(h); // // track-event trigger bias correction (MB to INEL) // hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoInel_corr_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoInel_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoInel_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoInel_corr_matrix"); aFolderObj->Add(h2D); // efficiency h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1),fGenTrackEventMatrix->Projection(1),"pt_track_trig_MBtoInel_eff_matrix"); aFolderObj->Add(h); // // track-event vertex reconstruction correction (MB) // hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix,fRecTrackEventMatrix,"zv_pt_eta_track_event_corr_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,2),fRecTrackEventMatrix->Projection(1,2),"eta_pt_track_event_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,0),fRecTrackEventMatrix->Projection(1,0),"pt_zv_track_event_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(2,0),fRecTrackEventMatrix->Projection(2,0),"zv_eta_track_event_corr_matrix"); aFolderObj->Add(h2D); // efficiency h = AlidNdPtHelper::GenerateCorrMatrix(fRecTrackEventMatrix->Projection(1),fTriggerTrackEventMatrix->Projection(1),"pt_track_event_eff_matrix"); aFolderObj->Add(h); // // track rec. efficiency correction // hs = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix,fRecPrimTrackMatrix,"zv_pt_eta_track_corr_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,2),fRecPrimTrackMatrix->Projection(1,2),"eta_pt_track_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,0),fRecPrimTrackMatrix->Projection(1,0),"pt_zv_track_corr_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2,0),fRecPrimTrackMatrix->Projection(2,0),"zv_eta_track_corr_matrix"); aFolderObj->Add(h2D); h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(0),fRecPrimTrackMatrix->Projection(0),"zv_track_corr_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1),fRecPrimTrackMatrix->Projection(1),"pt_track_corr_matrix"); aFolderObj->Add(h); // efficiency h = AlidNdPtHelper::GenerateCorrMatrix(fRecPrimTrackMatrix->Projection(1), fGenPrimTrackMatrix->Projection(1),"pt_track_eff_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2),fRecPrimTrackMatrix->Projection(2),"eta_track_corr_matrix"); aFolderObj->Add(h); // // secondary track contamination correction // //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix"); hs = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_track_cont_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_track_cont_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_track_cont_matrix"); aFolderObj->Add(h2D); h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_track_cont_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_track_cont_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_track_cont_matrix"); aFolderObj->Add(h); // // multiple track reconstruction correction // //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix"); hs = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix"); aFolderObj->Add(hs); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_mult_track_cont_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_mult_track_cont_matrix"); aFolderObj->Add(h2D); h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_mult_track_cont_matrix"); aFolderObj->Add(h2D); h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_mult_track_cont_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_mult_track_cont_matrix"); aFolderObj->Add(h); h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_mult_track_cont_matrix"); aFolderObj->Add(h); // // Control histograms // if(fHistogramsOn) { // Efficiency electrons, muons, pions, kaons, protons, all fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,1); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,1); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_electrons"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(2,2); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(2,2); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_muons"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(3,3); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(3,3); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_pions"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(4,4); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(4,4); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_kaons"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(5,5); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(5,5); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_protons"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,5); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,5); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_selected"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,6); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,6); h1 = fMCPrimTrackHist1[1]->Projection(0); h2 = fMCPrimTrackHist1[2]->Projection(0); h2c = (TH1D *)h2->Clone(); h2c->Divide(h1); h2c->SetName("eff_pt_all"); aFolderObj->Add(h2c); fMCPrimTrackHist1[1]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[1]->GetAxis(1)->GetNbins()); fMCPrimTrackHist1[2]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[2]->GetAxis(1)->GetNbins()); // pt spetra // - rec, primaries, secondaries // - primaries (pid) // - secondaries (pid) // - secondaries (mech) // - secondaries (mother) // TH1D *mcPtAccall = fMCTrackHist1[1]->Projection(0); mcPtAccall->SetName("mc_pt_acc_all"); aFolderObj->Add(mcPtAccall); TH1D *mcPtAccprim = fMCPrimTrackHist1[1]->Projection(0); mcPtAccprim->SetName("mc_pt_acc_prim"); aFolderObj->Add(mcPtAccprim); TH1D *mcPtRecall = fMCTrackHist1[2]->Projection(0); mcPtRecall->SetName("mc_pt_rec_all"); aFolderObj->Add(mcPtRecall); TH1D *mcPtRecprim = fMCPrimTrackHist1[2]->Projection(0); mcPtRecprim->SetName("mc_pt_rec_prim"); aFolderObj->Add(mcPtRecprim); TH1D *mcPtRecsec = fMCSecTrackHist1[2]->Projection(0); mcPtRecsec->SetName("mc_pt_rec_sec"); aFolderObj->Add(mcPtRecsec); for(Int_t i = 0; i<6; i++) { snprintf(name,256,"mc_pt_rec_prim_pid_%d",i); fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1); h = fMCPrimTrackHist1[2]->Projection(0); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_pt_rec_sec_pid_%d",i); fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1); h = fMCSecTrackHist1[2]->Projection(0); h->SetName(name); aFolderObj->Add(h); // production mechanisms for given pid fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1); for(Int_t j=0; j<20; j++) { if(j == 4) { // decay snprintf(name,256,"mc_pt_rec_sec_pid_%d_decay",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(0); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_eta_rec_sec_pid_%d_decay",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(1); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_mother_rec_sec_pid_%d_decay",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(4); h->SetName(name); aFolderObj->Add(h); } else if (j == 5) { // conversion snprintf(name,256,"mc_pt_rec_sec_pid_%d_conv",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(0); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_eta_rec_sec_pid_%d_conv",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(1); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_mother_rec_sec_pid_%d_conv",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(4); h->SetName(name); aFolderObj->Add(h); } else if (j == 13) { // mat snprintf(name,256,"mc_pt_rec_sec_pid_%d_mat",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(0); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_eta_rec_sec_pid_%d_mat",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(1); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_eta_mother_rec_sec_pid_%d_mat",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(4,1); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_mother_rec_sec_pid_%d_mat",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(4); h->SetName(name); aFolderObj->Add(h); snprintf(name,256,"mc_pt_mother_rec_sec_pid_%d_mat",i); fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1); h = fMCSecTrackHist1[2]->Projection(4,0); h->SetName(name); aFolderObj->Add(h); } else { continue; } } } } // end fHistogramOn // // resolution histograms // only for reconstructed tracks // TH2F *h2F=0; TCanvas * c = new TCanvas("resol","resol"); c->cd(); // fRecMCTrackHist1->GetAxis(1)->SetRangeUser(-0.8,0.79); h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0); h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10); h->SetXTitle("p_{tmc} (GeV/c)"); h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution"); h->Draw(); h->SetName("pt_resolution_vs_mcpt"); aFolderObj->Add(h); h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0); h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10); h->SetXTitle("p_{tmc} (GeV/c)"); h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean"); h->Draw(); h->SetName("dpt_mean_vs_mcpt"); aFolderObj->Add(h); // h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0); h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10); h->SetXTitle("p_{tmc} (GeV/c)"); h->SetYTitle("(#eta-#eta_{mc}) resolution"); h->Draw(); h->SetName("eta_resolution_vs_mcpt"); aFolderObj->Add(h); h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0); h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10); h->SetXTitle("p_{tmc} (GeV/c)"); h->SetYTitle("(#eta-mc#eta) mean"); h->Draw(); h->SetName("deta_mean_vs_mcpt"); aFolderObj->Add(h); // fRecMCTrackHist1->GetAxis(1)->SetRange(1,fRecMCTrackHist1->GetAxis(1)->GetNbins()); h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1); h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10); h->SetXTitle("#eta_{mc}"); h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution"); h->Draw(); h->SetName("pt_resolution_vs_mceta"); aFolderObj->Add(h); h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1); h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10); h->SetXTitle("#eta_{mc}"); h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean"); h->Draw(); h->SetName("dpt_mean_vs_mceta"); aFolderObj->Add(h); // h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1); h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10); h->SetXTitle("#eta_{mc}"); h->SetYTitle("(#eta-#eta_{mc}) resolution"); h->Draw(); h->SetName("eta_resolution_vs_mceta"); aFolderObj->Add(h); h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1); h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10); h->SetXTitle("#eta_{mc}"); h->SetYTitle("(#eta-mc#eta) mean"); h->Draw(); h->SetName("deta_mean_vs_mceta"); aFolderObj->Add(h); fRecMCTrackHist1->GetAxis(0)->SetRange(1,fRecMCTrackHist1->GetAxis(0)->GetNbins()); } // end use MC info // export objects to analysis folder fAnalysisFolder = ExportToFolder(aFolderObj); // delete only TObjArray if(aFolderObj) delete aFolderObj; } //_____________________________________________________________________________ TFolder* AlidNdPtAnalysisPbPb2011::ExportToFolder(TObjArray * const array) { // recreate folder avery time and export objects to new one // AlidNdPtAnalysisPbPb2011 * 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* AlidNdPtAnalysisPbPb2011::CreateFolder(TString name,TString title) { // create folder for analysed histograms // TFolder *folder = 0; folder = new TFolder(name.Data(),title.Data()); return folder; }