#include "TChain.h" #include "TTree.h" #include "TH1F.h" #include "TH2F.h" #include "TH3F.h" #include "TCanvas.h" #include "TList.h" #include "TParticle.h" #include "TParticlePDG.h" #include "TProfile.h" #include "AliAnalysisTask.h" #include "AliAnalysisManager.h" #include "AliESDEvent.h" #include "AliLog.h" #include "AliESDVertex.h" #include "AliESDInputHandler.h" #include "AliESDtrackCuts.h" #include "AliMultiplicity.h" #include "AliAnalysisTaskQASym.h" #include "AliExternalTrackParam.h" #include "AliTrackReference.h" #include "AliHeader.h" #include "AliGenEventHeader.h" #include "AliGenDPMjetEventHeader.h" // Analysis Task for basic QA on the ESD // Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing, // Andreas Morsch, Eva Sicking ClassImp(AliAnalysisTaskQASym) //________________________________________________________________________ AliAnalysisTaskQASym::AliAnalysisTaskQASym(const char *name) : AliAnalysisTaskSE(name) ,fTrackType(0) ,fStandAlone(0) ,fLow(0) ,fHigh(100) ,fFieldOn(kTRUE) ,fHists(0) ,fHistRECpt(0) ,fEta(0) ,fEtaWidth(0) ,fPhiWidth(0) ,fDcaWidth(0) ,fPtWidth(0) ,fEtaPhi(0) ,fEtaPt(0) ,fQPt(0) ,fDca(0) ,fDcaZ(0) ,fqRec(0) ,fSigmaPtHist(0) ,fRecPtPos(0) ,fRecPtNeg(0) ,fRecPhiPos(0) ,fRecPhiNeg(0) ,fRecEtaPos(0) ,fRecEtaNeg(0) ,fRecEtaPtPos(0) ,fRecEtaPtNeg(0) ,fRecDcaPos(0) ,fRecDcaNeg(0) ,fRecDcaNegInv(0) ,fRecDPos(0) ,fRecDNeg(0) ,fRecQPtPosEta(0) ,fRecQPtNegEta(0) ,fRecPtPosEta(0) ,fRecPtNegEta(0) ,fRecPhiPosEta(0) ,fRecPhiNegEta(0) ,fRecDcaPosEta(0) ,fRecDcaNegEta(0) ,fRecDPosEta(0) ,fRecDNegEta(0) ,fRecPtPosVz(0) ,fRecPtNegVz(0) ,fRecEtaPosVz(0) ,fRecEtaNegVz(0) ,fRecPhiPosVz(0) ,fRecPhiNegVz(0) ,fSignedDcaPosVz(0) ,fSignedDcaNegVz(0) ,fRecQPtPosEtaVz(0) ,fRecQPtNegEtaVz(0) ,fRecEtaPtPosVz(0) ,fRecEtaPtNegVz(0) ,fDeltaPhiAll(0) ,fDeltaPhiLeading(0) ,fDiffDcaD(0) ,fPhiRec(0) ,fThetaRec(0) ,fNumber(0) ,fNumberAfterCut(0) ,fVx(0) ,fVy(0) ,fVz(0) ,fNVertexSPD(0) ,fNVertexTracks(0) ,fRecDcaPosPhi(0) ,fRecDcaNegPhi(0) ,fRecPtPosPhi(0) ,fRecPtNegPhi(0) ,fRecEtaPosPhi(0) ,fRecEtaNegPhi(0) ,fRecQPtPhi(0) ,fRecEtaPtPosPhi(0) ,fRecEtaPtNegPhi(0) ,fRecPtPosEtaPos(0) ,fRecPtNegEtaPos(0) ,fRecPtPosEtaNeg(0) ,fRecPtNegEtaNeg(0) ,fRec1PtPosEtaPos(0) ,fRec1PtNegEtaPos(0) ,fRec1PtPosEtaNeg(0) ,fRec1PtNegEtaNeg(0) ,fRecPhiPosEtaPos(0) ,fRecPhiNegEtaPos(0) ,fRecPhiPosEtaNeg(0) ,fRecPhiNegEtaNeg(0) ,fRecDcaPosPhiEtaPos(0) ,fRecDcaNegPhiEtaPos(0) ,fRecDcaPosPhiEtaNeg(0) ,fRecDcaNegPhiEtaNeg(0) ,fRecDcaPosPtEtaPos(0) ,fRecDcaNegPtEtaPos(0) ,fRecDcaPosPtEtaNeg(0) ,fRecDcaNegPtEtaNeg(0) ,fRecPtPosPhiEtaPos(0) ,fRecPtNegPhiEtaPos(0) ,fRecPtPosPhiEtaNeg(0) ,fRecPtNegPhiEtaNeg(0) // ,fRecDcaPhiPtPosEtaPos(0) // ,fRecDcaPhiPtNegEtaPos(0) // ,fRecDcaPhiPtPosEtaNeg(0) // ,fRecDcaPhiPtNegEtaNeg(0) ,fEtavPt(0) ,fPhivPt(0) ,fCompareTPCparam(0) ,fITSlayer(0) ,fITSlayerEta(0) ,fITSlayerPhi(0) ,fCuts(0) { // Constructor for(Int_t i = 0;i<4;++i){ fVertexX[i]=0; fVertexY[i]=0; fVertexZ[i]=0; } for(Int_t i = 0;i<18;++i){ fRecPtTpcSector[i] = 0; fRecEtaTpcSector[i] = 0; fSignedDcaTpcSector[i] = 0; fRecQPtTpcSector[i] = 0; fRecEtaPtTpcSector[i] = 0; } for(Int_t i = 0;i< 7;++i){ fRecPtPosLadder[i] = 0; fRecPtNegLadder[i] = 0; fRecPhiPosLadder[i] = 0; fRecPhiNegLadder[i] = 0; fRecEtaPosLadder[i] = 0; fRecEtaNegLadder[i] = 0; fSignDcaPos[i] = 0; fSignDcaNeg[i] = 0; fSignDcaNegInv[i] = 0; fPtSigmaPos[i] =0; fPtSigmaNeg[i] =0; fqPtRec[i] =0; fDcaSigmaPos[i] =0; fDcaSigmaNeg[i] =0; } for(Int_t i = 0;i< 3;i++){ for(Int_t j = 0;j< 2;j++){ fEtaBinPt[i][j]=0; fPhiBinPt[i][j]=0; fDcaBinPt[i][j]=0; fEtaPhiBinPt[i][j]=0; } } DefineOutput(1, TList::Class()); } //________________________________________________________________________ void AliAnalysisTaskQASym::UserCreateOutputObjects() { // Create histograms // Called once Bool_t oldStatus = TH1::AddDirectoryStatus(); TH1::AddDirectory(kFALSE); Double_t range = 0.3; Double_t pt = 20.; fHists = new TList(); fHistRECpt = new TH1F("fHistRECpt", " p_{T}", 200, 0., pt); fEta = new TH1F("fEta", " #eta", 200, -2., 2.); fEtaWidth = new TH1F("fEtaWidth", " #eta", 200, -2., 2.); fPhiWidth = new TH1F("fPhiWidth", " #phi", 200, 0., 2*TMath::Pi()); fDcaWidth = new TH1F("fDcaWidth", "dca", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fPtWidth = new TH1F("fPtWidth", "p_{T}", 200, 0., pt); fEtavPt = new TH2F("fEtavPt", " #eta -p_{T}", 200, -2., 2., 100, -3, 4); fPhivPt = new TH2F("fPhivPt", " #phi -p_{T}", 200, 0, 2*TMath::Pi(), 100, -3, 5); fCompareTPCparam = new TH2F("fCompareTPCparam", "fCompareTPCparam", 100, -1., 1.,100,-5, 5); fITSlayer = new TH1F("fITSlayer", "fITSlayer", 8, -1.5, 6.5); fITSlayerEta = new TH2F("fITSlayerEta", "fITSlayerEta", 8, -1.5, 6.5, 200, -2.,2.); fITSlayerPhi = new TH2F("fITSlayerPhi", "fITSlayerPhi", 8, -1.5, 6.5, 200, 0,2*TMath::Pi()); fEtaPhi = new TH2F("fEtaPhi", " #eta - #phi", 200, -2., 2., 128, 0., 2. * TMath::Pi()); fThetaRec = new TH1F("fThetaRec", " #theta", 180, 0., TMath::Pi()); fPhiRec = new TH1F("fPhiRec", " #phi", 180, 0., 2*TMath::Pi()); fNumber = new TH1F("fNumber", "number of tracks per event", 500, -5, 4995); fNumberAfterCut = new TH1F("fNumberAfterCut", "number of tracks per event after cuts", 500, -5, 4995); fVx = new TH1F("fVx", "X of first track point", 100, -1., 1.); fVy = new TH1F("fVy", "Y of first track point", 100, -1., 1.); fVz = new TH1F("fVz", "Z of first track point", 200, -50., 50.); fNVertexSPD = new TH1F("fNVertexSPD", "Number of SPD vertices", 10, -0.5, 9.5); fNVertexTracks = new TH1F("fNVertexTracks", "Number of track vertices", 10, -0.5, 9.5); fEtaPt = new TH1F("fEtaPt", " #eta/p_{T} ", 100, -1., 1.); fQPt = new TH1F("fQPt", " charge/p_{T} ", 100, -1., 1.); fDca = new TH1F("fDca", " dca ", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fDcaZ = new TH1F("fDcaZ", "fDcaZ ",200, -3, 3);// limits fitting those of track cuts fqRec = new TH1F("fqRec", " charge all reconstructed particle", 21, -9.5, 10.5); fSigmaPtHist = new TH1F("fSigmaPtHist", "Log_{10}(#sigma_{p_{T}})", 200, -4., 8.); TString lable[4]={"", "SPD", "Track", "TPC"}; for(Int_t i=0;i<4;i++){ fVertexX[i] = new TH1F(Form("fVertexX%s",lable[i].Data()), Form("fVertexX%s",lable[i].Data()), 250, -1., 1.); fVertexY[i] = new TH1F(Form("fVertexY%s",lable[i].Data()), Form("fVertexY%s",lable[i].Data()), 250, -1., 1.); if(i==1 || i==2){ fVertexZ[i] = new TH1F(Form("fVertexZ%s",lable[i].Data()), Form("fVertexZ%s",lable[i].Data()), 200, -5., 5.); } else{ fVertexZ[i] = new TH1F(Form("fVertexZ%s",lable[i].Data()), Form("fVertexZ%s",lable[i].Data()), 200, -50., 50.); } } //------------ for(Int_t ITSlayer_case=0;ITSlayer_case<7;ITSlayer_case++){ fSignDcaPos[ITSlayer_case] = new TH1F(Form("fSignDcaPos%d", ITSlayer_case), " Signed dca", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fSignDcaPos[ITSlayer_case]->GetXaxis()->SetTitle("dca"); fSignDcaPos[ITSlayer_case]->GetYaxis()->SetTitle(""); fSignDcaNeg[ITSlayer_case] = new TH1F(Form("fSignDcaNeg%d", ITSlayer_case), " Signed dcas", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fSignDcaNeg[ITSlayer_case]->GetXaxis()->SetTitle("dca"); fSignDcaNeg[ITSlayer_case]->GetYaxis()->SetTitle(""); fSignDcaNegInv[ITSlayer_case] = new TH1F(Form("fSignDcaNegInv%d", ITSlayer_case), " inverse Signed dca ", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fSignDcaNegInv[ITSlayer_case]->GetXaxis()->SetTitle("-dca"); fSignDcaNegInv[ITSlayer_case]->GetYaxis()->SetTitle(""); fPtSigmaPos[ITSlayer_case] = new TH1F(Form("fPtSigmaPos%d", ITSlayer_case), " #sigma_{pT} ", 208, -4., 8.); fPtSigmaPos[ITSlayer_case]->GetXaxis()->SetTitle("Log_{10}(#sigma_{pT})"); fPtSigmaPos[ITSlayer_case]->GetYaxis()->SetTitle(""); fPtSigmaNeg[ITSlayer_case] = new TH1F(Form("fPtSigmaNeg%d",ITSlayer_case), " #sigma_{pT}", 208, -4., 8.); fPtSigmaNeg[ITSlayer_case]->GetXaxis()->SetTitle("Log_{10}(#sigma_{pT})"); fPtSigmaNeg[ITSlayer_case]->GetYaxis()->SetTitle(""); fqPtRec[ITSlayer_case] = new TH1F(Form("fqPtRec%d",ITSlayer_case), "q/ p_{T}", 200, -100., 100.); fqPtRec[ITSlayer_case]->GetXaxis()->SetTitle("q_{tr}/p_{T, tr} (GeV/c)"); fqPtRec[ITSlayer_case]->GetYaxis()->SetTitle(""); fDcaSigmaPos[ITSlayer_case] = new TH2F(Form("fDcaSigmaPos%d", ITSlayer_case), " p_{T} shift vs #sigma_{pT} ", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9),200, -4., 4. ); fDcaSigmaPos[ITSlayer_case]->GetXaxis()->SetTitle("signed DCA)"); fDcaSigmaPos[ITSlayer_case]->GetYaxis()->SetTitle("log_{10}(#sigma_{pT})"); fDcaSigmaNeg[ITSlayer_case] = new TH2F(Form("fDcaSigmaNeg%d", ITSlayer_case), " p_{T} shift vs #sigma_{pT} ", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9),200, -4., 4. ); fDcaSigmaNeg[ITSlayer_case]->GetXaxis()->SetTitle("signed DCA"); fDcaSigmaNeg[ITSlayer_case]->GetYaxis()->SetTitle("log_{10}(#sigma_{pT})"); } // YIELDs---------- positive and negative particles fRecPtPos = new TH1F("fRecPtPos", " p_{T}", 100, 0.,pt); fRecPtPos->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtNeg = new TH1F("fRecPtNeg", " p_{T} ", 100, 0., pt); fRecPtNeg->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPhiPos = new TH1F("fRecPhiPos", "#phi", 361, 0., 360.); fRecPhiPos->GetXaxis()->SetTitle("#phi (deg)"); fRecPhiNeg = new TH1F("fRecPhiNeg", "#phi ", 361, 0., 360.); fRecPhiNeg->GetXaxis()->SetTitle("#phi (deg)"); fRecEtaPos = new TH1F("fRecEtaPos", "#eta", 200, -2., 2.); fRecEtaPos->GetXaxis()->SetTitle("#eta"); fRecEtaNeg = new TH1F("fRecEtaNeg", "#eta", 200, -2., 2.); fRecEtaNeg->GetXaxis()->SetTitle("#eta"); fRecEtaPtPos = new TH1F("fRecEtaPtPos", "#eta/p_{T}", 200, -0.1, .1); fRecEtaPtPos->GetXaxis()->SetTitle("#eta/p_{T}"); fRecEtaPtNeg = new TH1F("fRecEtaPtNeg", "#eta/p_{T}", 200, -.1, .1); fRecEtaPtNeg->GetXaxis()->SetTitle("#eta/p_{T}"); fRecDcaPos = new TH1F("fRecDcaPos", " dca", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDcaPos->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNeg = new TH1F("fRecDcaNeg", " dca", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDcaNeg->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegInv = new TH1F("fRecDcaNegInv", " dca", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDcaNegInv->GetXaxis()->SetTitle("dca (cm)"); fRecDPos = new TH1F("fRecDPos", " d", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDPos->GetXaxis()->SetTitle("d (cm)"); fRecDNeg = new TH1F("fRecDNeg", "d", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDNeg->GetXaxis()->SetTitle("d (cm)"); // YIELDs ---------------- positive and negative eta fRecQPtPosEta = new TH1F("fRecQPtPosEta", "q/p_{T}", 200, -0.5, 0.5); fRecQPtPosEta->GetXaxis()->SetTitle("q/p_{T} "); fRecQPtNegEta = new TH1F("fRecQPtNegEta", "q/p_{T}", 200, -0.5, 0.5); fRecQPtNegEta->GetXaxis()->SetTitle("q/p_{T}"); fRecPtPosEta = new TH1F("fRecPtPosEta", " p_{T} ", 100, 0., pt); fRecPtPosEta->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtNegEta = new TH1F("fRecPtNegEta", " p_{T}", 100, 0., pt); fRecPtNegEta->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPhiPosEta = new TH1F("fRecPhiPosEta", "#phi", 361, 0., 360); fRecPhiPosEta->GetXaxis()->SetTitle("#phi (deg)"); fRecPhiNegEta = new TH1F("fRecPhiNegEta", "#phi ", 361, 0, 360); fRecPhiNegEta->GetXaxis()->SetTitle("#phi (deg)"); fRecDcaPosEta = new TH1F("fRecDcaPosEta", " dca ", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDcaPosEta->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegEta = new TH1F("fRecDcaNegEta", " dca", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDcaNegEta->GetXaxis()->SetTitle("dca (cm)"); fRecDPosEta = new TH1F("fRecDPosEta", " d", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); fRecDPosEta->GetXaxis()->SetTitle("d (cm)"); fRecDNegEta = new TH1F("fRecDNegEta", "d", 100, -5., 5.); fRecDNegEta->GetXaxis()->SetTitle("d (cm)"); fRecDcaPosPhi = new TH2F("fRecDcaPosPhi", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaPosPhi->GetXaxis()->SetTitle("dca (cm)"); fRecDcaPosPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecDcaNegPhi = new TH2F("fRecDcaNegPhi", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaNegPhi->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtPosPhi = new TH2F("fRecPtPosPhi", " log(p_T) vs. phi", 100, -2.5, 2., 180, 0, TMath::Pi()*2); fRecPtPosPhi->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtPosPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtNegPhi = new TH2F("fRecPtNegPhi", " log(p_T) vs. phi", 100,-2.5 , 2., 180, 0, TMath::Pi()*2); fRecPtNegPhi->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtNegPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecEtaPosPhi = new TH2F("fRecEtaPosPhi", " eta vs. phi", 100, -1.5, 1.5, 180, 0, TMath::Pi()*2); fRecEtaPosPhi->GetXaxis()->SetTitle("#eta"); fRecEtaPosPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecEtaNegPhi = new TH2F("fRecEtaNegPhi", " eta vs. phi", 100, -1.5, 1.5, 180, 0, TMath::Pi()*2); fRecEtaNegPhi->GetXaxis()->SetTitle("#eta"); fRecEtaNegPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecQPtPhi = new TH2F("fRecQPtPhi", " charge/p_T vs. phi", 100,-1. , 1., 180, 0, TMath::Pi()*2); fRecQPtPhi->GetXaxis()->SetTitle("charge/p_{T}"); fRecQPtPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecEtaPtPosPhi = new TH2F("fRecEtaPtPosPhi", " eta/p_T vs. phi", 100, -5, 5., 180, 0, TMath::Pi()*2); fRecEtaPtPosPhi->GetXaxis()->SetTitle("#eta/p_{T}"); fRecEtaPtPosPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecEtaPtNegPhi = new TH2F("fRecEtaPtNegPhi", " eta/p_T vs. phi", 100,-5 , 5., 180, 0, TMath::Pi()*2); fRecEtaPtNegPhi->GetXaxis()->SetTitle("#eta/p_{T}"); fRecEtaPtNegPhi->GetYaxis()->SetTitle("#phi (rad.)"); fRecDcaPosPhiEtaPos = new TH2F("fRecDcaPosPhiEtaPos", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaPosPhiEtaPos->GetXaxis()->SetTitle("dca (cm)"); fRecDcaPosPhiEtaPos->GetYaxis()->SetTitle("#phi (rad.)"); fRecDcaNegPhiEtaPos = new TH2F("fRecDcaNegPhiEtaPos", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaNegPhiEtaPos->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegPhiEtaPos->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtPosPhiEtaPos = new TH2F("fRecPtPosPhiEtaPos", " log(p_T) vs. phi", 100, -2.5, 2., 180, 0, TMath::Pi()*2); fRecPtPosPhiEtaPos->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtPosPhiEtaPos->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtNegPhiEtaPos = new TH2F("fRecPtNegPhiEtaPos", " log(p_T) vs. phi", 100,-2.5 , 2., 180, 0, TMath::Pi()*2); fRecPtNegPhiEtaPos->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtNegPhiEtaPos->GetYaxis()->SetTitle("#phi (rad.)"); fRecDcaPosPhiEtaNeg = new TH2F("fRecDcaPosPhiEtaNeg", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaPosPhiEtaNeg->GetXaxis()->SetTitle("dca (cm)"); fRecDcaPosPhiEtaNeg->GetYaxis()->SetTitle("#phi (rad.)"); fRecDcaNegPhiEtaNeg = new TH2F("fRecDcaNegPhiEtaNeg", " dca vs. phi", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 180, 0, TMath::Pi()*2); fRecDcaNegPhiEtaNeg->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegPhiEtaNeg->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtPosPhiEtaNeg = new TH2F("fRecPtPosPhiEtaNeg", " log(p_T) vs. phi", 100, -2.5, 2., 180, 0, TMath::Pi()*2); fRecPtPosPhiEtaNeg->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtPosPhiEtaNeg->GetYaxis()->SetTitle("#phi (rad.)"); fRecPtNegPhiEtaNeg = new TH2F("fRecPtNegPhiEtaNeg", " log(p_T) vs. phi", 100,-2.5 , 2., 180, 0, TMath::Pi()*2); fRecPtNegPhiEtaNeg->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtNegPhiEtaNeg->GetYaxis()->SetTitle("#phi (rad.)"); //new fRecDcaPosPtEtaPos = new TH2F("fRecDcaPosPtEtaPos", " dca vs. pt", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 200, -2, 2); fRecDcaPosPtEtaPos->GetXaxis()->SetTitle("dca (cm)"); fRecDcaPosPtEtaPos->GetYaxis()->SetTitle("log_{10}(p_{T})"); fRecDcaPosPtEtaNeg = new TH2F("fRecDcaPosPtEtaNeg", " dca vs. pt", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 200, -2, 2); fRecDcaPosPtEtaNeg->GetXaxis()->SetTitle("dca (cm)"); fRecDcaPosPtEtaNeg->GetYaxis()->SetTitle("log_{10}(p_{T})"); fRecDcaNegPtEtaPos = new TH2F("fRecDcaNegPtEtaPos", " dca vs. pt", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 200, -2, 2); fRecDcaNegPtEtaPos->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegPtEtaPos->GetYaxis()->SetTitle("log_{10}(p_{T})"); fRecDcaNegPtEtaNeg = new TH2F("fRecDcaNegPtEtaNeg", " dca vs. pt", 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 200, -2, 2); fRecDcaNegPtEtaNeg->GetXaxis()->SetTitle("dca (cm)"); fRecDcaNegPtEtaNeg->GetYaxis()->SetTitle("log_{10}(p_{T})"); // YIELDs ---------------- for TPC sectors for(Int_t sector=0; sector<18;sector++){ fRecPtTpcSector[sector] = new TH1F(Form("fRecPtTpcSector%02d",sector), Form("p_{T} distribution: TPC sector %d", sector),100, 0., pt); fRecPtTpcSector[sector]->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecEtaTpcSector[sector] = new TH1F(Form("fRecEtaTpcSector%02d",sector), Form("#eta distribution: TPC sector %d", sector),200, -2., 2.); fRecEtaTpcSector[sector]->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fSignedDcaTpcSector[sector] = new TH1F(Form("fSignedDcaTpcSector%02d",sector), Form("dca distribution: TPC sector %d", sector),200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9) ); fSignedDcaTpcSector[sector]->GetXaxis()->SetTitle("dca"); fRecQPtTpcSector[sector] = new TH1F(Form("fRecQPtTpcSector%02d",sector), Form("Q/ p_{T} distribution: TPC sector %d", sector),100, -1., 1.); fRecQPtTpcSector[sector]->GetXaxis()->SetTitle("Q/p_{T} (GeV/c)"); fRecEtaPtTpcSector[sector] = new TH1F(Form("fRecEtaPtTpcSector%02d",sector), Form("#eta/ p_{T} distribution: TPC sector %d", sector),100, -1., 1.); fRecEtaPtTpcSector[sector]->GetXaxis()->SetTitle("#eta/p_{T} (GeV/c)"); } // YIELDS ITS ladder for(Int_t i=0;i<7;i++){ fRecPtPosLadder[i] = new TH1F(Form("fRecPtPosLadder%d", i), " p_{T} distribution", 100, 0., pt); fRecPtPosLadder[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtNegLadder[i] = new TH1F(Form("fRecPtNegLadder%d",i), " p_{T} distribution ", 100, 0., pt); fRecPtNegLadder[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPhiPosLadder[i] = new TH1F(Form("fRecPhiPosLadder%d",i), "#phi distribution: all pos eta", 361, 0., 360); fRecPhiPosLadder[i]->GetXaxis()->SetTitle("#phi (deg)"); fRecPhiNegLadder[i] = new TH1F(Form("fRecPhiNegLadder%d", i), "#phi distribution: all neg eta", 361, 0, 360); fRecPhiNegLadder[i]->GetXaxis()->SetTitle("#phi (deg)"); fRecEtaPosLadder[i] = new TH1F(Form("fRecEtaPosLadder%d",i), "#eta distribution", 200, -2., 2.); fRecEtaPosLadder[i]->GetXaxis()->SetTitle("#eta)"); fRecEtaNegLadder[i] = new TH1F(Form("fRecEtaNegLadder%d", i), "#eta distribution", 200, -2., 2.); fRecEtaNegLadder[i]->GetXaxis()->SetTitle("#eta"); } Double_t vzmax = 15.; fRecPtPosVz = new TH2F("fRecPtPosVz", "p_{T} distribution vs Vz()", 100, -1., 2., 200,-vzmax,vzmax); fRecPtPosVz->GetXaxis()->SetTitle("log_{10}(p_{T})"); fRecPtNegVz = new TH2F("fRecPtNegVz", "p_{T} distribution vs Vz()", 100, -1., 2.,200,-vzmax,vzmax); fRecPtNegVz->GetXaxis()->SetTitle("Log_{10}(p_{T})"); fRecEtaPosVz= new TH2F("fRecEtaPosVz", "#eta distribution vs Vz()", 100, -2., 2., 200,-vzmax,vzmax); fRecEtaPosVz->GetXaxis()->SetTitle("#eta"); fRecEtaNegVz = new TH2F("fRecEtaNegVz", "#eta distribution vs Vz()", 100, -2., 2.,200,-vzmax,vzmax); fRecEtaNegVz->GetXaxis()->SetTitle("#eta"); fRecPhiPosVz= new TH2F("fRecPhiPosVz", "#eta distribution vs Vz()", 361, 0., 360., 200,-vzmax,vzmax); fRecPhiPosVz->GetXaxis()->SetTitle("#phi (deg)"); fRecPhiNegVz = new TH2F("fRecPhiNegVz", "dca vs Vz()", 361, 0., 360.,200,-vzmax,vzmax); fRecPhiNegVz->GetXaxis()->SetTitle("#phi (deg)"); fSignedDcaPosVz= new TH2F("fSignedDcaPosVz", "#eta distribution vs Vz()", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9), 200,-vzmax,vzmax); fSignedDcaPosVz->GetXaxis()->SetTitle("dca (cm)"); fSignedDcaNegVz = new TH2F("fSignedDcaNegVz", "dca vs Vz()", 200, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9),200,-vzmax,vzmax); fSignedDcaNegVz->GetXaxis()->SetTitle("dca (cm)"); fRecQPtPosEtaVz= new TH2F("fRecQPtPosEtaVz", " Q/p_{T} distribution vs Vz()", 100, -1., 1., 200,-vzmax,vzmax); fRecQPtPosEtaVz->GetXaxis()->SetTitle("Q/p_{T}"); fRecQPtNegEtaVz = new TH2F("fRecQPtNegEtaVz", " Q/p_{T} distribution vs Vz()", 100, -1., 1.,200,-vzmax,vzmax); fRecQPtNegEtaVz->GetXaxis()->SetTitle("Q/p_{T}"); fRecEtaPtPosVz= new TH2F("fRecEtaPtPosVz", " #eta/p_{T} distribution vs Vz()", 100, -1., 1., 200,-vzmax,vzmax); fRecEtaPtPosVz->GetXaxis()->SetTitle("#eta/p_{T"); fRecEtaPtNegVz = new TH2F("fRecEtaPtNegVz", " #eta/p_{T} distribution vs Vz()", 100, -1., 1.,200,-vzmax,vzmax); fRecEtaPtNegVz->GetXaxis()->SetTitle("#eta/p_{T}"); //new fDeltaPhiAll = new TH1F("fDeltaPhiAll", " #Delta #phi",200,-360,360); fDeltaPhiAll->GetXaxis()->SetTitle("#Delta #phi"); fDeltaPhiLeading = new TH2F("fDeltaPhiLeading", " #Delta #phi",361,-360,360, 361,0, 360); fDeltaPhiLeading->GetXaxis()->SetTitle("#Delta #phi (deg.)"); fDeltaPhiLeading->GetYaxis()->SetTitle("#phi_{leading particle} (deg.)"); fDiffDcaD = new TH1F("fDiffDcaD", "dca-d", 200, -1., 1.); fRecPtPosEtaPos = new TH1F("fRecPtPosEtaPos", "p_{T} distribution",100,0,pt); fRecPtPosEtaPos->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtNegEtaPos = new TH1F("fRecPtNegEtaPos", "p_{T} distribution",100,0,pt); fRecPtNegEtaPos->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtPosEtaNeg = new TH1F("fRecPtPosEtaNeg", "p_{T} distribution",100,0,pt); fRecPtPosEtaNeg->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRecPtNegEtaNeg = new TH1F("fRecPtNegEtaNeg", "p_{T} distribution",100,0,pt); fRecPtNegEtaNeg->GetXaxis()->SetTitle("p_{T} (GeV/c)"); fRec1PtPosEtaPos = new TH1F("fRec1PtPosEtaPos", "1/p_{T} distribution",100,0,0.5); fRec1PtPosEtaPos->GetXaxis()->SetTitle("p_{T} (c/GeV)"); fRec1PtNegEtaPos = new TH1F("fRec1PtNegEtaPos", "1/p_{T} distribution",100,0,0.5); fRec1PtNegEtaPos->GetXaxis()->SetTitle("p_{T} (c/GeV)"); fRec1PtPosEtaNeg = new TH1F("fRec1PtPosEtaNeg", "1/p_{T} distribution",100,0,0.5); fRec1PtPosEtaNeg->GetXaxis()->SetTitle("p_{T} (c/GeV)"); fRec1PtNegEtaNeg = new TH1F("fRec1PtNegEtaNeg", "1/p_{T} distribution",100,0,0.5); fRec1PtNegEtaNeg->GetXaxis()->SetTitle("1/p_{T} (c/GeV)"); fRecPhiPosEtaPos = new TH1F("fRecPhiPosEtaPos", "#phi",180,0,2*TMath::Pi()); fRecPhiPosEtaPos->GetXaxis()->SetTitle("#phi (rad.)"); fRecPhiNegEtaPos = new TH1F("fRecPhiNegEtaPos", "#phi",180,0,2*TMath::Pi()); fRecPhiNegEtaPos->GetXaxis()->SetTitle("#phi (rad.)"); fRecPhiPosEtaNeg = new TH1F("fRecPhiPosEtaNeg", "#phi",180,0,2*TMath::Pi()); fRecPhiPosEtaNeg->GetXaxis()->SetTitle("#phi (rad.)"); fRecPhiNegEtaNeg = new TH1F("fRecPhiNegEtaNeg", "#phi",180,0,2*TMath::Pi()); fRecPhiNegEtaNeg->GetXaxis()->SetTitle("#phi (rad.)"); // fRecDcaPhiPtPosEtaPos = new TH3F("fRecDcaPhiPtPosEtaPos", // "#phi- p_{T} - DCA", // 180,0,2*TMath::Pi(), // 100,0,pt, // 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); // fRecDcaPhiPtPosEtaPos->GetXaxis()->SetTitle("#phi (rad.)"); // fRecDcaPhiPtPosEtaPos->GetYaxis()->SetTitle("p_{T} (GeV/c)"); // fRecDcaPhiPtPosEtaPos->GetZaxis()->SetTitle("dca (cm)"); // fRecDcaPhiPtPosEtaNeg = new TH3F("fRecDcaPhiPtPosEtaNeg", // "#phi- p_{T} - DCA", // 180,0,2*TMath::Pi(), // 100,0,pt, // 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); // fRecDcaPhiPtPosEtaNeg->GetZaxis()->SetTitle("dca (cm)"); // fRecDcaPhiPtPosEtaNeg->GetXaxis()->SetTitle("#phi (rad.)"); // fRecDcaPhiPtPosEtaNeg->GetYaxis()->SetTitle("p_{T} (GeV/c)"); // fRecDcaPhiPtNegEtaPos = new TH3F("fRecDcaPhiPtNegEtaPos", // "#phi- p_{T} - DCA", // 180,0,2*TMath::Pi(), // 100,0,pt, // 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); // fRecDcaPhiPtNegEtaPos->GetZaxis()->SetTitle("dca (cm)"); // fRecDcaPhiPtNegEtaPos->GetXaxis()->SetTitle("#phi (rad.)"); // fRecDcaPhiPtNegEtaPos->GetYaxis()->SetTitle("p_{T} (GeV/c)"); // fRecDcaPhiPtNegEtaNeg = new TH3F("fRecDcaPhiPtNegEtaNeg", // "#phi- p_{T} - DCA", // 180,0,2*TMath::Pi(), // 100,0,pt, // 100, -range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9)); // fRecDcaPhiPtNegEtaNeg->GetZaxis()->SetTitle("dca (cm)"); // fRecDcaPhiPtNegEtaNeg->GetYaxis()->SetTitle("#phi (rad.)"); // fRecDcaPhiPtNegEtaNeg->GetXaxis()->SetTitle("p_{T} (GeV/c)"); TString charge[2]; charge[0]="Pos"; charge[1]="Neg"; for(Int_t i=0;i<3;i++){ for(Int_t j=0;j<2;j++){ fEtaBinPt[i][j] = new TH1F(Form("fEtaBinPt%d%s", i, charge[j].Data()), "eta", 200, -2., 2.); fPhiBinPt[i][j] = new TH1F(Form("fPhiBinPt%d%s", i,charge[j].Data() ), "phi", 181, 0, 2*TMath::Pi()); fDcaBinPt[i][j] = new TH1F(Form("fDcaBinPt%d%s", i, charge[j].Data()), "DCA", 200,-range*(1+Int_t(fTrackType/2)*9), range*(1+Int_t(fTrackType/2)*9) ); fEtaPhiBinPt[i][j]= new TH2F(Form("fEtaPhiBinPt%d%s", i, charge[j].Data()), "eta-phi", 200, -2., 2., 200, 0.,2*TMath::Pi()); } } fHists->SetOwner(); fHists->Add(fHistRECpt); fHists->Add(fEta); fHists->Add(fEtaWidth); fHists->Add(fPhiWidth); fHists->Add(fDcaWidth); fHists->Add(fPtWidth); fHists->Add(fEtavPt); fHists->Add(fPhivPt); fHists->Add(fCompareTPCparam); fHists->Add(fITSlayer); fHists->Add(fITSlayerEta); fHists->Add(fITSlayerPhi); fHists->Add(fEtaPhi); fHists->Add(fThetaRec); fHists->Add(fPhiRec); fHists->Add(fNumber); fHists->Add(fNumberAfterCut); fHists->Add(fVx); fHists->Add(fVy); fHists->Add(fVz); fHists->Add(fNVertexSPD); fHists->Add(fNVertexTracks); fHists->Add(fEtaPt); fHists->Add(fQPt); fHists->Add(fDca); fHists->Add(fDcaZ); fHists->Add(fDeltaPhiAll); fHists->Add(fDeltaPhiLeading); fHists->Add(fDiffDcaD); fHists->Add(fqRec); fHists->Add(fSigmaPtHist); fHists->Add(fRecPtPos); fHists->Add(fRecPtNeg); fHists->Add(fRecPhiPos); fHists->Add(fRecPhiNeg); fHists->Add(fRecEtaPos); fHists->Add(fRecEtaNeg); fHists->Add(fRecEtaPtPos); fHists->Add(fRecEtaPtNeg); fHists->Add(fRecDcaPos); fHists->Add(fRecDcaNeg); fHists->Add(fRecDcaNegInv); fHists->Add(fRecDPos); fHists->Add(fRecDNeg); fHists->Add(fRecQPtPosEta); fHists->Add(fRecQPtNegEta); fHists->Add(fRecPtPosEta); fHists->Add(fRecPtNegEta); fHists->Add(fRecPhiPosEta); fHists->Add(fRecPhiNegEta); fHists->Add(fRecDcaPosEta); fHists->Add(fRecDcaNegEta); fHists->Add(fRecDPosEta); fHists->Add(fRecDNegEta); for(Int_t i=0;i<4;i++){ fHists->Add(fVertexX[i]); fHists->Add(fVertexY[i]); fHists->Add(fVertexZ[i]); } for(Int_t i=0;i<18;i++){ fHists->Add(fRecPtTpcSector[i]); fHists->Add(fRecEtaTpcSector[i]); fHists->Add(fSignedDcaTpcSector[i]); fHists->Add(fRecQPtTpcSector[i]); fHists->Add(fRecEtaPtTpcSector[i]); } for(Int_t i=0;i<7;i++){ fHists->Add(fRecPtPosLadder[i]); fHists->Add(fRecPtNegLadder[i]); fHists->Add(fRecPhiPosLadder[i]); fHists->Add(fRecPhiNegLadder[i]); fHists->Add(fRecEtaPosLadder[i]); fHists->Add(fRecEtaNegLadder[i]); } fHists->Add(fRecPtPosVz); fHists->Add(fRecPtNegVz); fHists->Add(fRecEtaPosVz); fHists->Add(fRecEtaNegVz); fHists->Add(fRecPhiPosVz); fHists->Add(fRecPhiNegVz); fHists->Add(fSignedDcaPosVz); fHists->Add(fSignedDcaNegVz); fHists->Add(fRecQPtPosEtaVz); fHists->Add(fRecQPtNegEtaVz); fHists->Add(fRecEtaPtPosVz); fHists->Add(fRecEtaPtNegVz); for(Int_t i=0;i<7;i++){ fHists->Add(fSignDcaPos[i]); fHists->Add(fSignDcaNeg[i]); fHists->Add(fSignDcaNegInv[i]); fHists->Add(fPtSigmaPos[i]); fHists->Add(fPtSigmaNeg[i]); fHists->Add(fqPtRec[i]); fHists->Add(fDcaSigmaPos[i]); fHists->Add(fDcaSigmaNeg[i]); } fHists->Add(fRecDcaPosPhi); fHists->Add(fRecDcaNegPhi); fHists->Add(fRecPtPosPhi); fHists->Add(fRecPtNegPhi); fHists->Add(fRecEtaPosPhi); fHists->Add(fRecEtaNegPhi); fHists->Add(fRecQPtPhi); fHists->Add(fRecEtaPtPosPhi); fHists->Add(fRecEtaPtNegPhi); fHists->Add(fRecPtPosEtaPos); fHists->Add(fRecPtNegEtaPos); fHists->Add(fRecPtPosEtaNeg); fHists->Add(fRecPtNegEtaNeg); fHists->Add(fRec1PtPosEtaPos); fHists->Add(fRec1PtNegEtaPos); fHists->Add(fRec1PtPosEtaNeg); fHists->Add(fRec1PtNegEtaNeg); fHists->Add(fRecPhiPosEtaPos); fHists->Add(fRecPhiNegEtaPos); fHists->Add(fRecPhiPosEtaNeg); fHists->Add(fRecPhiNegEtaNeg); fHists->Add(fRecDcaPosPhiEtaPos); fHists->Add(fRecDcaNegPhiEtaPos); fHists->Add(fRecPtPosPhiEtaPos); fHists->Add(fRecPtNegPhiEtaPos); fHists->Add(fRecDcaPosPhiEtaNeg); fHists->Add(fRecDcaNegPhiEtaNeg); fHists->Add(fRecPtPosPhiEtaNeg); fHists->Add(fRecPtNegPhiEtaNeg); fHists->Add(fRecDcaPosPtEtaPos); fHists->Add(fRecDcaNegPtEtaPos); fHists->Add(fRecDcaPosPtEtaNeg); fHists->Add(fRecDcaNegPtEtaNeg); // fHists->Add(fRecDcaPhiPtPosEtaPos); // fHists->Add(fRecDcaPhiPtPosEtaNeg); // fHists->Add(fRecDcaPhiPtNegEtaPos); // fHists->Add(fRecDcaPhiPtNegEtaNeg); for(Int_t i=0;i<3;i++){ for(Int_t j=0;j<2;j++){ fHists->Add(fEtaBinPt[i][j]); fHists->Add(fPhiBinPt[i][j]); fHists->Add(fDcaBinPt[i][j]); fHists->Add(fEtaPhiBinPt[i][j]); } } // for (Int_t i=0; iGetEntries(); ++i) { // TH1 *h1 = dynamic_cast(fHists->At(i)); // if (h1){ // // Printf("%s ",h1->GetName()); // h1->Sumw2(); // } // } TH1::AddDirectory(oldStatus); PostData(1, fHists); } //__________________________________________________________ void AliAnalysisTaskQASym::UserExec(Option_t *) { // QA of global, TPC, ITS and ITS stand alone tracks // exploiting basic symmetries AliVEvent *event = InputEvent(); if (!event) { Printf("ERROR: Could not retrieve event"); return; } if(Entry()==0){ AliESDEvent* esd = static_cast(event); if(esd){ Printf("We are reading from ESD"); } } if(fDebug>1)Printf("There are %d tracks in this event", event->GetNumberOfTracks()); Int_t leadingTrack = 0; Float_t leadingEnergy = -20.; Float_t leadingPhi = 0;//TMath::Pi(); //check vertices AliESDEvent* esd = static_cast(event); Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks(); fNVertexSPD->Fill(nPileSPDVertices); fNVertexTracks->Fill(nPileTrkVertices); //check primary vertex Float_t vx = 0; Float_t vy = 0; Float_t vz = 0; //primary vertex: contribution from different vertexers const AliVVertex* vertex = event->GetPrimaryVertex(); if(!vertex) return; vx = vertex->GetX(); vy = vertex->GetY(); vz = vertex->GetZ(); if(vertex->GetNContributors()>0){ fVertexX[0]->Fill(vx); fVertexY[0]->Fill(vy); fVertexZ[0]->Fill(vz); } const AliVVertex* vertexSPD = esd->GetPrimaryVertexSPD(); if(vertexSPD){ if(vertexSPD->GetNContributors()>0){ fVertexX[1]->Fill(vertexSPD->GetX()); fVertexY[1]->Fill(vertexSPD->GetY()); fVertexZ[1]->Fill(vertexSPD->GetZ()); } } const AliVVertex* vertexTrack = esd->GetPrimaryVertexTracks(); if(vertexTrack){ if(vertexTrack->GetNContributors()>0){ fVertexX[2]->Fill(vertexTrack->GetX()); fVertexY[2]->Fill(vertexTrack->GetY()); fVertexZ[2]->Fill(vertexTrack->GetZ()); } } const AliVVertex* vertexTPC = esd->GetPrimaryVertexTPC(); if(vertexTPC){ if(vertexTPC->GetNContributors()>0){ fVertexX[3]->Fill(vertexTPC->GetX()); fVertexY[3]->Fill(vertexTPC->GetY()); fVertexZ[3]->Fill(vertexTPC->GetZ()); } } //cuts on general vertex if(vertex->GetNContributors()<1) return; if (TMath::Abs(vz) > 10.) return; fNumber->Fill(event->GetNumberOfTracks()); AliESDtrack *tpcP = 0x0; Int_t fNTracksAccepted=0; const Int_t arrSize = event->GetNumberOfTracks(); Float_t * phiArray = new Float_t[arrSize]; Float_t * etaArray = new Float_t[arrSize]; Float_t * ptArray = new Float_t[arrSize]; Float_t * dcaArray = new Float_t[arrSize]; Int_t * chargeArray = new Int_t [arrSize]; Bool_t * acceptedArray = new Bool_t [arrSize]; for (Int_t i = 0; i < event->GetNumberOfTracks(); i++) { phiArray[i] = 0.; etaArray[i] = 0.; ptArray[i] = 0.; dcaArray[i] = 0.; chargeArray[i] = 0; acceptedArray[i]= kFALSE; } for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++) { //prevent mem leak for TPConly track if(fTrackType==2&&tpcP){ delete tpcP; tpcP = 0; } AliVParticle *track = event->GetTrack(iTrack); AliESDtrack *esdtrack = static_cast(track); esdtrack->PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 10000.); if (!track) { Printf("ERROR: Could not receive track %d", iTrack); continue; } //__________ // run Task for global tracks or ITS tracks or TPC tracks const AliExternalTrackParam *tpcPin = 0x0; Double_t phiIn=0.; if(fTrackType==0){ //Fill all histograms with global tracks tpcP = esdtrack; if (!tpcP) continue; if (!fCuts->AcceptTrack(tpcP)) continue; phiIn = tpcP->Phi(); } else if(fTrackType==1){ //Fill all histograms with ITS tracks tpcP = esdtrack; phiIn = tpcP->Phi(); if (!tpcP) continue; if (!fCuts->AcceptTrack(tpcP)) continue; if(fStandAlone==kTRUE) { if(!(tpcP->GetStatus()&AliESDtrack::kITSpureSA))continue; } else if(fStandAlone==kFALSE){ if(tpcP->GetStatus()&AliESDtrack::kITSpureSA)continue; } } else if(fTrackType==2){ //Fill all histograms with TPC track information tpcPin = esdtrack->GetInnerParam(); if (!tpcPin) continue; phiIn=tpcPin->Phi(); tpcP = AliESDtrackCuts::GetTPCOnlyTrack(static_cast(event),esdtrack->GetID()); if (!tpcP) continue; if (!fCuts->AcceptTrack(tpcP)) continue; if(tpcP->GetNcls(1)>160)continue;//jacek's track cut if(tpcP->GetConstrainedChi2TPC()<0)continue; // jacek's track cut } else{ Printf("ERROR: wrong track type \n"); continue; } //___________ // fNTracksAccepted++; phiArray[iTrack] = phiIn; etaArray[iTrack] = tpcP->Eta(); ptArray[iTrack] = tpcP->Pt(); chargeArray[iTrack] = tpcP->Charge(); acceptedArray[iTrack]= kTRUE; if(tpcP->E()>leadingEnergy){ leadingTrack=iTrack; leadingEnergy=tpcP->E(); leadingPhi=phiIn; } fqRec->Fill(tpcP->Charge()); Double_t fSigmaPt = tpcP->GetSigma1Pt2(); fSigmaPt= sqrt(fSigmaPt); fSigmaPt= fSigmaPt *(tpcP->Pt()*tpcP->Pt()); if(TMath::Abs(fSigmaPt) < 1.e-10) continue; fSigmaPtHist->Fill(TMath::Log10(fSigmaPt)); // hits in ITS layer Int_t cas=-1; if(tpcP->HasPointOnITSLayer(0)) cas=0; else if(!tpcP->HasPointOnITSLayer(0) && tpcP->HasPointOnITSLayer(1)) cas=1; else if(!tpcP->HasPointOnITSLayer(0) && !tpcP->HasPointOnITSLayer(1) && tpcP->HasPointOnITSLayer(2)) cas=2; else if(!tpcP->HasPointOnITSLayer(0) && !tpcP->HasPointOnITSLayer(1) && !tpcP->HasPointOnITSLayer(2) && tpcP->HasPointOnITSLayer(3)) cas=3; else if(!tpcP->HasPointOnITSLayer(0) && !tpcP->HasPointOnITSLayer(1) && !tpcP->HasPointOnITSLayer(2) && !tpcP->HasPointOnITSLayer(3) && tpcP->HasPointOnITSLayer(4)) cas=4; else if( !tpcP->HasPointOnITSLayer(0) && !tpcP->HasPointOnITSLayer(1) && !tpcP->HasPointOnITSLayer(2) && !tpcP->HasPointOnITSLayer(3) && !tpcP->HasPointOnITSLayer(4) && tpcP->HasPointOnITSLayer(5)) cas=5; else cas=6; //------------------- Float_t fXVertexCor = 0.; Float_t fYVertexCor = 0.; fXVertexCor = tpcP->Xv() - vertex->GetX(); // coordinate corrected for vertex position fYVertexCor = tpcP->Yv() - vertex->GetY(); // " Double_t fSignedDca = (tpcP->Py()*fXVertexCor - tpcP->Px()*fYVertexCor)/tpcP->Pt(); fqPtRec[cas]->Fill(tpcP->Charge()/tpcP->Pt()); fHistRECpt->Fill(tpcP->Pt()); fEta->Fill(tpcP->Eta()); fEtavPt->Fill(tpcP->Eta(), TMath::Log(tpcP->Pt())); fPhivPt->Fill(phiIn, TMath::Log(tpcP->Pt())); fEtaPhi->Fill(tpcP->Eta(), phiIn); fThetaRec->Fill(tpcP->Theta()); fPhiRec->Fill(phiIn); fVx->Fill(tpcP->Xv()); fVy->Fill(tpcP->Yv()); fVz->Fill(tpcP->Zv()); fEtaPt->Fill(tpcP->Eta()/tpcP->Pt()); fQPt->Fill(tpcP->Charge()/tpcP->Pt()); fDca->Fill(fSignedDca); dcaArray[iTrack]=fSignedDca; fRecQPtPhi->Fill(tpcP->Charge()/tpcP->Pt(), phiIn); Float_t fXY = 0.; Float_t fZ = 0.; tpcP->GetImpactParameters(fXY,fZ); fDiffDcaD->Fill(fSignedDca+fXY); fDcaZ->Fill(fZ); if(fTrackType==2) fCompareTPCparam->Fill(fZ,tpcPin->GetTgl()); if(fTrackType!=2){//for global and ITS tracks for(Int_t itsLayer=0;itsLayer<6;itsLayer++){ if(tpcP->HasPointOnITSLayer(itsLayer)){ fITSlayer->Fill(itsLayer); fITSlayerEta->Fill(itsLayer, tpcP->Eta()); fITSlayerPhi->Fill(itsLayer, tpcP->Phi()); } } } //for positive particles if(tpcP->Charge()>0){ fRecPtPos->Fill(tpcP->Pt()); fRecPtPosLadder[cas]->Fill(tpcP->Pt()); fRecPtPosVz->Fill(TMath::Log10(tpcP->Pt()),tpcP->Zv()); fRecPhiPos->Fill(TMath::RadToDeg()*phiIn); fRecPhiPosLadder[cas]->Fill(TMath::RadToDeg()*phiIn); fRecPhiPosVz->Fill(TMath::RadToDeg()*phiIn,tpcP->Zv()); fSignedDcaPosVz->Fill(fSignedDca,tpcP->Zv()); fRecEtaPos->Fill(tpcP->Eta()); fRecEtaPosLadder[cas]->Fill(tpcP->Eta()); fRecEtaPtPos->Fill(tpcP->Eta()/tpcP->Pt()); fRecEtaPosVz->Fill(tpcP->Eta(),tpcP->Zv()); fRecEtaPtPosVz->Fill(tpcP->Eta()/tpcP->Pt(),tpcP->Zv()); fRecDcaPos->Fill(fSignedDca); fRecDcaPosPhi->Fill(fSignedDca, phiIn); fRecPtPosPhi->Fill(TMath::Log10(tpcP->Pt()), phiIn); fRecEtaPtPosPhi->Fill(tpcP->Eta()/tpcP->Pt(), phiIn); fRecEtaPosPhi->Fill(tpcP->Eta(), phiIn); fRecDPos->Fill(fXY); fSignDcaPos[cas]->Fill(fSignedDca); fDcaSigmaPos[cas]->Fill(fSignedDca, TMath::Log10(fSigmaPt)); fPtSigmaPos[cas]->Fill(TMath::Log10(fSigmaPt)); //pos eta if(tpcP->Eta()>0){ fRecPtPosEtaPos->Fill(tpcP->Pt()); fRec1PtPosEtaPos->Fill(1/tpcP->Pt()); fRecPhiPosEtaPos->Fill(phiIn); fRecDcaPosPhiEtaPos->Fill(fSignedDca, phiIn); fRecDcaPosPtEtaPos->Fill(fSignedDca, TMath::Log10(tpcP->Pt())); fRecPtPosPhiEtaPos->Fill(TMath::Log10(tpcP->Pt()), phiIn); //fRecDcaPhiPtPosEtaPos->Fill(phiIn, tpcP->Pt(), fSignedDca); } //neg eta else{ fRecPtPosEtaNeg->Fill(tpcP->Pt()); fRec1PtPosEtaNeg->Fill(1/tpcP->Pt()); fRecPhiPosEtaNeg->Fill(phiIn); fRecDcaPosPhiEtaNeg->Fill(fSignedDca, phiIn); fRecDcaPosPtEtaNeg->Fill(fSignedDca, TMath::Log10(tpcP->Pt())); fRecPtPosPhiEtaNeg->Fill(TMath::Log10(tpcP->Pt()), phiIn); //fRecDcaPhiPtPosEtaNeg->Fill(phiIn, tpcP->Pt(), fSignedDca); } } //and negative particles else { fRecPtNeg->Fill(tpcP->Pt()); fRecPtNegLadder[cas]->Fill(tpcP->Pt()); fRecPtNegVz->Fill(TMath::Log10(tpcP->Pt()),tpcP->Zv()); fRecPhiNeg->Fill(TMath::RadToDeg()*phiIn); fRecPhiNegLadder[cas]->Fill(TMath::RadToDeg()*phiIn); fRecPhiNegVz->Fill(TMath::RadToDeg()*phiIn,tpcP->Zv()); fSignedDcaNegVz->Fill(fSignedDca,tpcP->Zv()); fRecEtaPtNegVz->Fill(tpcP->Eta()/tpcP->Pt(),tpcP->Zv()); fRecEtaNeg->Fill(tpcP->Eta()); fRecEtaNegLadder[cas]->Fill(tpcP->Eta()); fRecEtaPtNeg->Fill(tpcP->Eta()/tpcP->Pt()); fRecEtaNegVz->Fill(tpcP->Eta(),tpcP->Zv()); fRecDcaNeg->Fill(fSignedDca); fRecDcaNegInv->Fill(-fSignedDca); fRecDcaNegPhi->Fill(fSignedDca, phiIn); fRecPtNegPhi->Fill(TMath::Log10(tpcP->Pt()), phiIn); fRecEtaNegPhi->Fill(tpcP->Eta(), phiIn); fRecEtaPtNegPhi->Fill(tpcP->Eta()/tpcP->Pt(), phiIn); fRecDNeg->Fill(fXY); fSignDcaNeg[cas]->Fill(fSignedDca); fSignDcaNegInv[cas]->Fill(-fSignedDca); fDcaSigmaNeg[cas]->Fill(fSignedDca,TMath::Log10(fSigmaPt)); fPtSigmaNeg[cas]->Fill(TMath::Log10(fSigmaPt)); //pos eta if(tpcP->Eta()>0){ fRecPtNegEtaPos->Fill(tpcP->Pt()); fRec1PtNegEtaPos->Fill(1/tpcP->Pt()); fRecPhiNegEtaPos->Fill(phiIn); fRecDcaNegPhiEtaPos->Fill(fSignedDca, phiIn); fRecDcaNegPtEtaPos->Fill(fSignedDca, TMath::Log10(tpcP->Pt())); fRecPtNegPhiEtaPos->Fill(TMath::Log10(tpcP->Pt()), phiIn); //fRecDcaPhiPtNegEtaPos->Fill(phiIn, tpcP->Pt(), fSignedDca); } //neg eta else{ fRecPtNegEtaNeg->Fill(tpcP->Pt()); fRec1PtNegEtaNeg->Fill(1/tpcP->Pt()); fRecPhiNegEtaNeg->Fill(phiIn); fRecDcaNegPhiEtaNeg->Fill(fSignedDca, phiIn); fRecDcaNegPtEtaNeg->Fill(fSignedDca, TMath::Log10(tpcP->Pt())); fRecPtNegPhiEtaNeg->Fill(TMath::Log10(tpcP->Pt()), phiIn); //fRecDcaPhiPtNegEtaNeg->Fill(phiIn, tpcP->Pt(), fSignedDca); } } //all particles with positive eta if(tpcP->Eta()>0){ fRecQPtPosEta->Fill(tpcP->Charge()/tpcP->Pt()); fRecPtPosEta->Fill(tpcP->Pt()); fRecPhiPosEta->Fill(TMath::RadToDeg()*phiIn); fRecQPtPosEtaVz->Fill(tpcP->Charge()/tpcP->Pt(),tpcP->Zv()); fRecDcaPosEta->Fill(fSignedDca); fRecDPosEta->Fill(fXY); } //all particles with negative eta (and eta==0) else{ fRecQPtNegEta->Fill(tpcP->Charge()/tpcP->Pt()); fRecPtNegEta->Fill(tpcP->Pt()); fRecPhiNegEta->Fill(TMath::RadToDeg()*phiIn); fRecQPtNegEtaVz->Fill(tpcP->Charge()/tpcP->Pt(),tpcP->Zv()); fRecDcaNegEta->Fill(fSignedDca); fRecDNegEta->Fill(fXY); } fRecPtTpcSector[Int_t(phiIn* TMath::RadToDeg()/20)]->Fill(tpcP->Pt()); fRecEtaTpcSector[Int_t(phiIn* TMath::RadToDeg()/20)]->Fill(tpcP->Eta()); fSignedDcaTpcSector[Int_t(phiIn* TMath::RadToDeg()/20)]->Fill(fSignedDca); fRecQPtTpcSector[Int_t(phiIn* TMath::RadToDeg()/20)]->Fill(tpcP->Charge()/tpcP->Pt()); fRecEtaPtTpcSector[Int_t(phiIn* TMath::RadToDeg()/20)]->Fill(tpcP->Eta()/tpcP->Pt()); // // another track loop // for (Int_t iTrack2 = 0; iTrack2 < event->GetNumberOfTracks(); iTrack2++) { // if(LeadingTrack==iTrack2) continue; // AliVParticle *track2 = event->GetTrack(iTrack2); // AliESDtrack* esdtrack2 = static_cast(track2); // if (!track2) { // Printf("ERROR: Could not receive track %d", iTrack); // continue; // } // if (!fCuts->AcceptTrack(esdtrack2)) continue; // //propagate to dca // esdtrack2->PropagateToDCA(event->GetPrimaryVertex(), // event->GetMagneticField(), 10000.); // fDeltaPhiLeading->Fill((LeadingPhi-esdtrack2->Phi())*TMath::RadToDeg(), // LeadingPhi*TMath::RadToDeg() ); // }//second track loop // if(fTrackType==2) delete tpcP; // delete in case of TPCOnlyTrack }//first track loop fNumberAfterCut->Fill(fNTracksAccepted); //second track loop for (Int_t iT = 0; iT < event->GetNumberOfTracks(); iT++) { if(acceptedArray[iT]){ if(ptArray[iT]>0.2 && ptArray[iT]<1. ){ fEtaBinPt[0][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT]); fDcaBinPt[0][Bool_t(chargeArray[iT]>0)]->Fill(dcaArray[iT]); fPhiBinPt[0][Bool_t(chargeArray[iT]>0)]->Fill(phiArray[iT]); fEtaPhiBinPt[0][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT], phiArray[iT]); } else if(ptArray[iT]>1. && ptArray[iT]<5.){ fEtaBinPt[1][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT]); fDcaBinPt[1][Bool_t(chargeArray[iT]>0)]->Fill(dcaArray[iT]); fPhiBinPt[1][Bool_t(chargeArray[iT]>0)]->Fill(phiArray[iT]); fEtaPhiBinPt[1][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT], phiArray[iT]); } else if (ptArray[iT]>5.){ fEtaBinPt[2][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT]); fDcaBinPt[2][Bool_t(chargeArray[iT]>0)]->Fill(dcaArray[iT]); fPhiBinPt[2][Bool_t(chargeArray[iT]>0)]->Fill(phiArray[iT]); fEtaPhiBinPt[2][Bool_t(chargeArray[iT]>0)]->Fill(etaArray[iT], phiArray[iT]); } if(fNTracksAccepted>=fLow&&fNTracksAccepted<=fHigh){ fEtaWidth->Fill(etaArray[iT]); fPhiWidth->Fill(phiArray[iT]); fDcaWidth->Fill(dcaArray[iT]); fPtWidth->Fill(ptArray[iT]); } } } //prevent mem leak for TPConly track if(fTrackType==2&&tpcP){ delete tpcP; tpcP = 0; } if(phiArray){ delete[] phiArray; phiArray=0; } if(etaArray){ delete[] etaArray; etaArray=0; } if(dcaArray){ delete[] dcaArray; dcaArray=0; } if(ptArray){ delete[] ptArray; ptArray=0; } if(chargeArray){ delete[] chargeArray; chargeArray=0; } if(acceptedArray){ delete[] acceptedArray; acceptedArray=0; } // Post output data. // PostData(1, fHistPt); PostData(1, fHists); } //________________________________________________________________________ void AliAnalysisTaskQASym::Terminate(Option_t *) { }