#include "AliAODMCParticle.h"
#include "AliMCParticle.h"
#include "AliCentrality.h"
+#include "AliKFVertex.h"
+#include "AliExternalTrackParam.h"
+#include "AliTriggerAnalysis.h"
// my headers
#include "AliAnalysisTaskUpcPsi2s.h"
//_____________________________________________________________________________
AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s()
- : AliAnalysisTaskSE(),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
+ : AliAnalysisTaskSE(),fType(0),isMC(kFALSE),fRunTree(kTRUE),fRunHist(kTRUE),fPIDResponse(0),fJPsiTree(0),fPsi2sTree(0),
fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
- fVtxContrib(0),fVtxPosX(0),fVtxPosY(0),fVtxPosZ(0),fVtxErrX(0),fVtxErrY(0),fVtxErrZ(0),fVtxChi2(0),fVtxNDF(0),fVtxType(0),
- fBCrossNum(0),fNtracklets(0),
+ fTOFtrig1(0), fTOFtrig2(0),
+ fVtxContrib(0),fVtxChi2(0),fVtxNDF(0),
+ fBCrossNum(0),fNtracklets(0),fNLooseTracks(0),
fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
fDataFilnam(0),fRecoPass(0),fEvtNum(0),
- fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
- fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
- fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
+ fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),fGenPart(0),
+ fListTrig(0),fHistCcup4TriggersPerRun(0), fHistCcup7TriggersPerRun(0), fHistCcup2TriggersPerRun(0),
+ fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0), fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
//_____________________________________________________________________________
AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s(const char *name)
- : AliAnalysisTaskSE(name),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
+ : AliAnalysisTaskSE(name),fType(0),isMC(kFALSE),fRunTree(kTRUE),fRunHist(kTRUE),fPIDResponse(0),fJPsiTree(0),fPsi2sTree(0),
fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
- fVtxContrib(0),fVtxPosX(0),fVtxPosY(0),fVtxPosZ(0),fVtxErrX(0),fVtxErrY(0),fVtxErrZ(0),fVtxChi2(0),fVtxNDF(0),fVtxType(0),
- fBCrossNum(0),fNtracklets(0),
+ fTOFtrig1(0), fTOFtrig2(0),
+ fVtxContrib(0),fVtxChi2(0),fVtxNDF(0),
+ fBCrossNum(0),fNtracklets(0),fNLooseTracks(0),
fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
fDataFilnam(0),fRecoPass(0),fEvtNum(0),
- fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
- fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
- fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
+ fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),fGenPart(0),
+ fListTrig(0),fHistCcup4TriggersPerRun(0), fHistCcup7TriggersPerRun(0), fHistCcup2TriggersPerRun(0),
+ fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0), fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
{
for(Int_t i=0; i<ntrg; i++) fTrigger[i] = kFALSE;
+ for(Int_t i=0; i<4; i++) {
+ fTOFphi[i] = -666;
+ fPIDMuon[i] = -666;
+ fPIDElectron[i] = -666;
+ fPIDPion[i] = -666;
+ fTriggerInputsMC[i] = kFALSE;
+ }
+ for(Int_t i=0; i<3; i++){
+ fVtxPos[i] = -666;
+ fVtxErr[i] = -666;
+ fKfVtxPos[i] = -666;
+ }
}//Init
//_____________________________________________________________________________
void AliAnalysisTaskUpcPsi2s::UserCreateOutputObjects()
{
+ //PID response
+ AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
+ AliInputEventHandler *inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
+ fPIDResponse = inputHandler->GetPIDResponse();
+
//input file
fDataFilnam = new TObjString();
fDataFilnam->SetString("");
fJPsiESDTracks = new TClonesArray("AliESDtrack", 1000);
fPsi2sAODTracks = new TClonesArray("AliAODTrack", 1000);
fPsi2sESDTracks = new TClonesArray("AliESDtrack", 1000);
+ fGenPart = new TClonesArray("TParticle", 1000);
//output tree with JPsi candidate events
fJPsiTree = new TTree("fJPsiTree", "fJPsiTree");
fJPsiTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
fJPsiTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
fJPsiTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
+ fJPsiTree ->Branch("fNLooseTracks", &fNLooseTracks, "fNLooseTracks/s");
fJPsiTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
- fJPsiTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
- fJPsiTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
- fJPsiTree ->Branch("fVtxPosZ", &fVtxPosZ, "fVtxPosZ/D");
- fJPsiTree ->Branch("fVtxErrX", &fVtxErrX, "fVtxErrX/D");
- fJPsiTree ->Branch("fVtxErrY", &fVtxErrY, "fVtxErrY/D");
- fJPsiTree ->Branch("fVtxErrZ", &fVtxErrZ, "fVtxErrZ/D");
+ fJPsiTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
+ fJPsiTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
+ fJPsiTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
+
+ fJPsiTree ->Branch("fPIDMuon", &fPIDMuon[0], "fPIDMuon[4]/D");
+ fJPsiTree ->Branch("fPIDElectron", &fPIDElectron[0], "fPIDElectron[4]/D");
+ fJPsiTree ->Branch("fPIDPion", &fPIDPion[0], "fPIDPion[4]/D");
+
+ fJPsiTree ->Branch("fVtxPos", &fVtxPos[0], "fVtxPos[3]/D");
+ fJPsiTree ->Branch("fVtxErr", &fVtxErr[0], "fVtxErr[3]/D");
fJPsiTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
fJPsiTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
- fJPsiTree ->Branch("fVtxType", &fVtxType, "fVtxType/B");
+
+ fJPsiTree ->Branch("fKfVtxPos", &fKfVtxPos[0], "fKfVtxPos[3]/D");
fJPsiTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
fJPsiTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
if( fType == 1 ) {
fJPsiTree ->Branch("fJPsiAODTracks", &fJPsiAODTracks);
}
+ if(isMC) {
+ fJPsiTree ->Branch("fGenPart", &fGenPart);
+ fJPsiTree ->Branch("fTriggerInputsMC", &fTriggerInputsMC[0], "fTriggerInputsMC[4]/O");
+ }
+
//output tree with Psi2s candidate events
fPsi2sTree = new TTree("fPsi2sTree", "fPsi2sTree");
fPsi2sTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
fPsi2sTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
fPsi2sTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
+ fPsi2sTree ->Branch("fNLooseTracks", &fNLooseTracks, "fNLooseTracks/s");
fPsi2sTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
- fPsi2sTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
- fPsi2sTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
- fPsi2sTree ->Branch("fVtxPosZ", &fVtxPosZ, "fVtxPosZ/D");
- fPsi2sTree ->Branch("fVtxErrX", &fVtxErrX, "fVtxErrX/D");
- fPsi2sTree ->Branch("fVtxErrY", &fVtxErrY, "fVtxErrY/D");
- fPsi2sTree ->Branch("fVtxErrZ", &fVtxErrZ, "fVtxErrZ/D");
+ fPsi2sTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
+ fPsi2sTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
+ fPsi2sTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
+
+ fPsi2sTree ->Branch("fPIDMuon", &fPIDMuon[0], "fPIDMuon[4]/D");
+ fPsi2sTree ->Branch("fPIDElectron", &fPIDElectron[0], "fPIDElectron[4]/D");
+ fPsi2sTree ->Branch("fPIDPion", &fPIDPion[0], "fPIDPion[4]/D");
+
+ fPsi2sTree ->Branch("fVtxPos", &fVtxPos[0], "fVtxPos[3]/D");
+ fPsi2sTree ->Branch("fVtxErr", &fVtxErr[0], "fVtxErr[3]/D");
fPsi2sTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
fPsi2sTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
- fPsi2sTree ->Branch("fVtxType", &fVtxType, "fVtxType/B");
+
+ fPsi2sTree ->Branch("fKfVtxPos", &fKfVtxPos[0], "fKfVtxPos[3]/D");
fPsi2sTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
fPsi2sTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
if( fType == 1 ) {
fPsi2sTree ->Branch("fPsi2sAODTracks", &fPsi2sAODTracks);
}
+ if(isMC) {
+ fPsi2sTree ->Branch("fGenPart", &fGenPart);
+ fPsi2sTree ->Branch("fTriggerInputsMC", &fTriggerInputsMC[0], "fTriggerInputsMC[4]/O");
+ }
+
fListTrig = new TList();
fListTrig ->SetOwner();
- fHistUpcTriggersPerRun = new TH1D("fHistUpcTriggersPerRun", "fHistUpcTriggersPerRun", 3000, 167000.5, 170000.5);
- fListTrig->Add(fHistUpcTriggersPerRun);
+ fHistCcup4TriggersPerRun = new TH1D("fHistCcup4TriggersPerRun", "fHistCcup4TriggersPerRun", 33000, 167000.5, 200000.5);
+ fListTrig->Add(fHistCcup4TriggersPerRun);
+
+ fHistCcup7TriggersPerRun = new TH1D("fHistCcup7TriggersPerRun", "fHistCcup7TriggersPerRun", 33000, 167000.5, 200000.5);
+ fListTrig->Add(fHistCcup7TriggersPerRun);
+
+ fHistCcup2TriggersPerRun = new TH1D("fHistCcup2TriggersPerRun", "fHistCcup2TriggersPerRun", 33000, 167000.5, 200000.5);
+ fListTrig->Add(fHistCcup2TriggersPerRun);
- fHistZedTriggersPerRun = new TH1D("fHistZedTriggersPerRun", "fHistZedTriggersPerRun", 3000, 167000.5, 170000.5);
+ fHistZedTriggersPerRun = new TH1D("fHistZedTriggersPerRun", "fHistZedTriggersPerRun", 33000, 167000.5, 200000.5);
fListTrig->Add(fHistZedTriggersPerRun);
- fHistCvlnTriggersPerRun = new TH1D("fHistCvlnTriggersPerRun", "fHistCvlnTriggersPerRun", 3000, 167000.5, 170000.5);
+ fHistCvlnTriggersPerRun = new TH1D("fHistCvlnTriggersPerRun", "fHistCvlnTriggersPerRun", 33000, 167000.5, 200000.5);
fListTrig->Add(fHistCvlnTriggersPerRun);
- fHistMBTriggersPerRun = new TH1D("fHistMBTriggersPerRun", "fHistMBTriggersPerRun", 3000, 167000.5, 170000.5);
+ fHistMBTriggersPerRun = new TH1D("fHistMBTriggersPerRun", "fHistMBTriggersPerRun", 33000, 167000.5, 200000.5);
fListTrig->Add(fHistMBTriggersPerRun);
- fHistCentralTriggersPerRun = new TH1D("fHistCentralTriggersPerRun", "fHistCentralTriggersPerRun", 3000, 167000.5, 170000.5);
+ fHistCentralTriggersPerRun = new TH1D("fHistCentralTriggersPerRun", "fHistCentralTriggersPerRun", 33000, 167000.5, 200000.5);
fListTrig->Add(fHistCentralTriggersPerRun);
- fHistSemiCentralTriggersPerRun = new TH1D("fHistSemiCentralTriggersPerRun", "fHistSemiCentralTriggersPerRun", 3000, 167000.5, 170000.5);
+ fHistSemiCentralTriggersPerRun = new TH1D("fHistSemiCentralTriggersPerRun", "fHistSemiCentralTriggersPerRun", 33000, 167000.5, 200000.5);
fListTrig->Add(fHistSemiCentralTriggersPerRun);
fListHist = new TList();
fListHist ->SetOwner();
- TString CutNameJPsi[12] = {"Analyzed","Triggered","Vertex cut","V0 decision","Two good tracks",
+ TString CutNameJPsi[13] = {"Analyzed","Triggered","Vertex cut","V0 decision","Neutron ZDC cut","Two good tracks",
"Like sign","Oposite sign","One p_{T}>1", "Both p_{T}>1","PID","Dimuom","Dielectron"};
- fHistNeventsJPsi = new TH1D("fHistNeventsJPsi","fHistNeventsPsi2s",12,0.5,12.5);
- for (Int_t i = 0; i<12; i++) fHistNeventsJPsi->GetXaxis()->SetBinLabel(i+1,CutNameJPsi[i].Data());
+ fHistNeventsJPsi = new TH1D("fHistNeventsJPsi","fHistNeventsPsi2s",13,0.5,13.5);
+ for (Int_t i = 0; i<13; i++) fHistNeventsJPsi->GetXaxis()->SetBinLabel(i+1,CutNameJPsi[i].Data());
fListHist->Add(fHistNeventsJPsi);
fHistTPCsignalJPsi = new TH2D("fHistTPCsignalJPsi","fHistTPCsignalJPsi",240,0,120,240,0,120);
fHistDiMuonMass->GetXaxis()->SetTitle("Invariant mass(#mu^{+}#mu^{-}) (GeV/c)");
fListHist->Add(fHistDiMuonMass);
- TString CutNamePsi2s[13] = {"Analyzed","Triggered","Vertex cut","V0 decision","Four good tracks",
+ TString CutNamePsi2s[14] = {"Analyzed","Triggered","Vertex cut","V0 decision","Neutron ZDC cut","Four good tracks",
"DiLepton - DiPion","Like sign leptons","Like sign pions","Like sign both","Oposite sign","PID","Dimuom","Dielectron"};
- fHistNeventsPsi2s = new TH1D("fHistNeventsPsi2s","fHistNeventsPsi2s",13,0.5,13.5);
- for (Int_t i = 0; i<13; i++) fHistNeventsPsi2s->GetXaxis()->SetBinLabel(i+1,CutNamePsi2s[i].Data());
+ fHistNeventsPsi2s = new TH1D("fHistNeventsPsi2s","fHistNeventsPsi2s",14,0.5,14.5);
+ for (Int_t i = 0; i<14; i++) fHistNeventsPsi2s->GetXaxis()->SetBinLabel(i+1,CutNamePsi2s[i].Data());
fListHist->Add(fHistNeventsPsi2s);
fHistPsi2sMassVsPt = new TH2D("fHistPsi2sMassVsPt","Mass vs p_{T} of #psi(2s) candidates",100,3,6,50,0,5);
//Trigger
TString trigger = aod->GetFiredTriggerClasses();
- if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
+ if(trigger.Contains("CCUP4-B")) fHistCcup4TriggersPerRun->Fill(fRunNum); //CCUP4 triggers
+ if(trigger.Contains("CCUP7-B")) fHistCcup7TriggersPerRun->Fill(fRunNum); //CCUP7 triggers
+ if(trigger.Contains("CCUP2-B")) fHistCcup2TriggersPerRun->Fill(fRunNum); //CCUP2 triggers
if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
AliCentrality *centrality = aod->GetCentrality();
UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
- Double_t percentile = centrality->GetCentralityPercentile("V0M");
+ Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
+ //Double_t percentile = centrality->GetCentralityPercentile("V0M");
- if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
+ if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<=80 && percentile>=0) fHistMBTriggersPerRun->Fill(fRunNum);
- if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
-
- if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
+ if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<=6 && percentile>=0 && (trigger.Contains("CVHN_R2-B"))) fHistCentralTriggersPerRun->Fill(fRunNum);
+ if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<=50 && percentile>=15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
+
PostData(3, fListTrig);
}
//Trigger
TString trigger = aod->GetFiredTriggerClasses();
- if( !trigger.Contains("CCUP4-B") ) return;
+ if(!isMC && !trigger.Contains("CCUP") ) return;
fHistNeventsJPsi->Fill(2);
fHistNeventsPsi2s->Fill(2);
fV0Cdecision = fV0data->GetV0CDecision();
if(fV0Adecision != AliAODVZERO::kV0Empty || fV0Cdecision != AliAODVZERO::kV0Empty) return;
+ fHistNeventsJPsi->Fill(4);
+ fHistNeventsPsi2s->Fill(4);
+
fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
- fHistNeventsJPsi->Fill(4);
- fHistNeventsPsi2s->Fill(4);
+ fHistNeventsJPsi->Fill(5);
+ fHistNeventsPsi2s->Fill(5);
- Int_t nGoodTracks=0;
//Two tracks loop
+ Int_t nGoodTracks = 0;
Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
Short_t qLepton[4], qPion[4];
UInt_t nLepton=0, nPion=0, nHighPt=0;
- Double_t jRecTPCsignal[5];
+ Double_t fRecTPCsignal[5];
Int_t mass[3]={-1,-1,-1};
- //Two track loop
- for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
- AliAODTrack *trk = aod->GetTrack(itr);
- if( !trk ) continue;
-
- if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
- if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
- if(trk->GetTPCNcls() < 70)continue;
- if(trk->Chi2perNDF() > 4)continue;
- if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
- Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
- if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
- if(TMath::Abs(dca[1]) > 2) continue;
- if(TMath::Abs(dca[0]) > 0.2) continue;
-
- TrackIndex[nGoodTracks] = itr;
- nGoodTracks++;
-
- if(nGoodTracks > 2) break;
- }//Track loop
-
- if(nGoodTracks == 2){
- fHistNeventsJPsi->Fill(5);
- for(Int_t i=0; i<2; i++){
- AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
- if(trk->Pt() > 1) nHighPt++;
- jRecTPCsignal[nLepton] = trk->GetTPCsignal();
- qLepton[nLepton] = trk->Charge();
- if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
- vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
- mass[nLepton] = 0;
- }
- if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
- vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
- mass[nLepton] = 1;
- }
- nLepton++;
- }
- if(nLepton == 2){
- if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
- if(qLepton[0]*qLepton[1] < 0){
- fHistNeventsJPsi->Fill(7);
- if(nHighPt > 0){
- fHistNeventsJPsi->Fill(8);
- fHistTPCsignalJPsi->Fill(jRecTPCsignal[0],jRecTPCsignal[1]);
- if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
- if(mass[0] == mass[1] && mass[0] != -1) {
- fHistNeventsJPsi->Fill(10);
- vCandidate = vLepton[0]+vLepton[1];
- if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
- if(mass[0] == 0) {
- fHistDiMuonMass->Fill(vCandidate.M());
- fHistNeventsJPsi->Fill(11);
- }
- if(mass[0] == 1) {
- fHistDiElectronMass->Fill(vCandidate.M());
- fHistNeventsJPsi->Fill(12);
- }
- }
- }
- }
- }
- }
-
- nGoodTracks = 0;
+
//Four track loop
for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
AliAODTrack *trk = aod->GetTrack(itr);
if( !trk ) continue;
+ if(!(trk->TestFilterBit(1<<0))) continue;
if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
if(trk->GetTPCNcls() < 50)continue;
if(trk->Chi2perNDF() > 4)continue;
Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
- if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
+
if(TMath::Abs(dca[1]) > 2) continue;
TrackIndex[nGoodTracks] = itr;
mass[0]= -1; mass[1]= -1, mass[2]= -1;
if(nGoodTracks == 4){
- fHistNeventsPsi2s->Fill(5);
+ fHistNeventsPsi2s->Fill(6);
for(Int_t i=0; i<4; i++){
AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
if(trk->Pt() > 1){
- jRecTPCsignal[nLepton] = trk->GetTPCsignal();
+ fRecTPCsignal[nLepton] = trk->GetTPCsignal();
qLepton[nLepton] = trk->Charge();
- if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
+ if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
mass[nLepton] = 0;
}
- if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
+ if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
mass[nLepton] = 1;
}
if(nLepton > 2 || nPion > 2) break;
}
if((nLepton == 2) && (nPion == 2)){
- fHistNeventsPsi2s->Fill(6);
- if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
- if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
- if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
+ fHistNeventsPsi2s->Fill(7);
+ if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(8);
+ if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(9);
+ if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(10);
if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
- fHistNeventsPsi2s->Fill(10);
+ fHistNeventsPsi2s->Fill(11);
if(mass[0] == mass[1]) {
- fHistNeventsPsi2s->Fill(11);
+ fHistNeventsPsi2s->Fill(12);
vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
vDilepton = vLepton[0]+vLepton[1];
fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
- if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
- if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
+ if(mass[0] == 0) fHistNeventsPsi2s->Fill(13);
+ if(mass[0] == 1) fHistNeventsPsi2s->Fill(14);
+ }
+ }
+ }
+ }
+
+ nGoodTracks = 0;
+ //Two track loop
+ for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
+ AliAODTrack *trk = aod->GetTrack(itr);
+ if( !trk ) continue;
+ if(!(trk->TestFilterBit(1<<0))) continue;
+
+ if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(trk->GetTPCNcls() < 70)continue;
+ if(trk->Chi2perNDF() > 4)continue;
+ if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
+ Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
+ if(TMath::Abs(dca[1]) > 2) continue;
+ if(TMath::Abs(dca[0]) > 0.2) continue;
+
+ TrackIndex[nGoodTracks] = itr;
+ nGoodTracks++;
+
+ if(nGoodTracks > 2) break;
+ }//Track loop
+
+ nLepton=0; nPion=0; nHighPt=0;
+ mass[0]= -1; mass[1]= -1, mass[2]= -1;
+
+ if(nGoodTracks == 2){
+ fHistNeventsJPsi->Fill(6);
+ for(Int_t i=0; i<2; i++){
+ AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
+ if(trk->Pt() > 1) nHighPt++;
+ fRecTPCsignal[nLepton] = trk->GetTPCsignal();
+ qLepton[nLepton] = trk->Charge();
+ if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
+ vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
+ mass[nLepton] = 0;
+ }
+ if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
+ vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
+ mass[nLepton] = 1;
+ }
+ nLepton++;
+ }
+ if(nLepton == 2){
+ if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(7);
+ if(qLepton[0]*qLepton[1] < 0){
+ fHistNeventsJPsi->Fill(8);
+ if(nHighPt > 0){
+ fHistNeventsJPsi->Fill(9);
+ fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
+ if(nHighPt == 2) fHistNeventsJPsi->Fill(10);
+ if(mass[0] == mass[1] && mass[0] != -1) {
+ fHistNeventsJPsi->Fill(11);
+ vCandidate = vLepton[0]+vLepton[1];
+ if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
+ if(mass[0] == 0) {
+ fHistDiMuonMass->Fill(vCandidate.M());
+ fHistNeventsJPsi->Fill(12);
+ }
+ if(mass[0] == 1) {
+ fHistDiElectronMass->Fill(vCandidate.M());
+ fHistNeventsJPsi->Fill(13);
+ }
+ }
}
}
}
}
+
PostData(4, fListHist);
}
AliAODEvent *aod = (AliAODEvent*) InputEvent();
if(!aod) return;
+ if(isMC) RunAODMC(aod);
+
//input data
const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
fDataFilnam->Clear();
//Trigger
TString trigger = aod->GetFiredTriggerClasses();
- if( !trigger.Contains("CCUP4-B") ) return;
+ fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
+ fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
+ fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
+
+ Bool_t isTriggered = kFALSE;
+ for(Int_t i=0; i<ntrg; i++) {
+ if( fTrigger[i] ) isTriggered = kTRUE;
+ }
+ if(!isMC && !isTriggered ) return;
//trigger inputs
fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
- fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
+ fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
//Event identification
fPerNum = aod ->GetPeriodNumber();
//primary vertex
AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
fVtxContrib = fAODVertex->GetNContributors();
- fVtxPosX = fAODVertex->GetX();
- fVtxPosY = fAODVertex->GetY();
- fVtxPosZ = fAODVertex->GetZ();
+ fVtxPos[0] = fAODVertex->GetX();
+ fVtxPos[1] = fAODVertex->GetY();
+ fVtxPos[2] = fAODVertex->GetZ();
Double_t CovMatx[6];
fAODVertex->GetCovarianceMatrix(CovMatx);
- fVtxErrX = CovMatx[0];
- fVtxErrY = CovMatx[1];
- fVtxErrZ = CovMatx[2];
+ fVtxErr[0] = CovMatx[0];
+ fVtxErr[1] = CovMatx[1];
+ fVtxErr[2] = CovMatx[2];
fVtxChi2 = fAODVertex->GetChi2();
fVtxNDF = fAODVertex->GetNDF();
- fVtxType = fAODVertex->GetType();
//Tracklets
fNtracklets = aod->GetTracklets()->GetNumberOfTracklets();
fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
+ fNLooseTracks = 0;
+
+ //Track loop - loose cuts
+ for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
+ AliAODTrack *trk = aod->GetTrack(itr);
+ if( !trk ) continue;
+ if(!(trk->TestFilterBit(1<<0))) continue;
+
+ if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
+ if(trk->GetTPCNcls() < 20)continue;
+ fNLooseTracks++;
+ }//Track loop -loose cuts
+
Int_t nGoodTracks=0;
Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
-
+
//Two track loop
for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
AliAODTrack *trk = aod->GetTrack(itr);
if( !trk ) continue;
-
- if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
- if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(!(trk->TestFilterBit(1<<0))) continue;
+
+ if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
if(trk->GetTPCNcls() < 70)continue;
if(trk->Chi2perNDF() > 4)continue;
if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
- if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
if(TMath::Abs(dca[1]) > 2) continue;
if(TMath::Abs(dca[0]) > 0.2) continue;
if(nGoodTracks > 2) break;
}//Track loop
-
+
+ fJPsiAODTracks->Clear("C");
if(nGoodTracks == 2){
+
+ TDatabasePDG *pdgdat = TDatabasePDG::Instance();
+ TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
+ Double_t muonMass = partMuon->Mass();
+ TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
+ Double_t electronMass = partElectron->Mass();
+ TParticlePDG *partPion = pdgdat->GetParticle( 211 );
+ Double_t pionMass = partPion->Mass();
+
+ Double_t KFcov[21];
+ Double_t KFpar[6];
+ Double_t KFmass = pionMass;
+ Double_t fRecTPCsignal;
+ AliKFParticle *KFpart[2];
+ AliKFVertex *KFvtx = new AliKFVertex();
+ KFvtx->SetField(aod->GetMagneticField());
+
for(Int_t i=0; i<2; i++){
AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
- trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov);
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
- trk->SetDCA(dca[0],dca[1]); //to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
new((*fJPsiAODTracks)[i]) AliAODTrack(*trk);
+ ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
+ fPIDMuon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kMuon);
+ fPIDElectron[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kElectron);
+ fPIDPion[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kPion);
+
+ trk->GetPosition(KFpar);
+ trk->PxPyPz(KFpar+3);
+ trk->GetCovMatrix(KFcov);
+
+ if(trk->Pt() > 1){
+ fRecTPCsignal = trk->GetTPCsignal();
+ if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
+ if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
+ }
+ else KFmass = pionMass;
+
+ KFpart[i] = new AliKFParticle();
+ KFpart[i]->SetField(aod->GetMagneticField());
+ KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
+ KFvtx->AddDaughter(*KFpart[i]);
+
+
+ Double_t pos[3]={0,0,0};
+ AliExternalTrackParam *parTrk = new AliExternalTrackParam();
+ parTrk->CopyFromVTrack((AliVTrack*) trk);
+ if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
+ else {
+ fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+ if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+ }
+ delete parTrk;
}
- fJPsiTree ->Fill();
- PostData(1, fJPsiTree);
- }
+ fKfVtxPos[0]= KFvtx->GetX();
+ fKfVtxPos[1]= KFvtx->GetY();
+ fKfVtxPos[2]= KFvtx->GetZ();
+ for(UInt_t i=0; i<2; i++)delete KFpart[i];
+ delete KFvtx;
- nGoodTracks = 0;
- //Four track loop
+ if(!isMC) fJPsiTree ->Fill();
+ }
+
+ nGoodTracks = 0;
+ //Four track loop
for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
AliAODTrack *trk = aod->GetTrack(itr);
if( !trk ) continue;
+ if(!(trk->TestFilterBit(1<<0))) continue;
- if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
- if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
if(trk->GetTPCNcls() < 50)continue;
if(trk->Chi2perNDF() > 4)continue;
Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
if(TMath::Abs(dca[1]) > 2) continue;
-
+
TrackIndex[nGoodTracks] = itr;
nGoodTracks++;
if(nGoodTracks > 4) break;
}//Track loop
-
-
+
+ fPsi2sAODTracks->Clear("C");
if(nGoodTracks == 4){
+
+ TDatabasePDG *pdgdat = TDatabasePDG::Instance();
+ TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
+ Double_t muonMass = partMuon->Mass();
+ TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
+ Double_t electronMass = partElectron->Mass();
+ TParticlePDG *partPion = pdgdat->GetParticle( 211 );
+ Double_t pionMass = partPion->Mass();
+
+ Double_t KFcov[21];
+ Double_t KFpar[6];
+ Double_t KFmass = pionMass;
+ Double_t fRecTPCsignal;
+ AliKFParticle *KFpart[4];
+ AliKFVertex *KFvtx = new AliKFVertex();
+ KFvtx->SetField(aod->GetMagneticField());
+
for(Int_t i=0; i<4; i++){
AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
- trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov);
+ AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
+ if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+ delete trk_clone;
+
+ new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
+ ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
- trk->SetDCA(dca[0],dca[1]); //to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
- new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
+ fPIDMuon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kMuon);
+ fPIDElectron[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kElectron);
+ fPIDPion[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kPion);
+
+ trk->GetPosition(KFpar);
+ trk->PxPyPz(KFpar+3);
+ trk->GetCovMatrix(KFcov);
+
+ if(trk->Pt() > 1){
+ fRecTPCsignal = trk->GetTPCsignal();
+ if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
+ if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
+ }
+ else KFmass = pionMass;
+
+ KFpart[i] = new AliKFParticle();
+ KFpart[i]->SetField(aod->GetMagneticField());
+ KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
+ KFvtx->AddDaughter(*KFpart[i]);
+
+ Double_t pos[3]={0,0,0};
+ AliExternalTrackParam *parTrk = new AliExternalTrackParam();
+ parTrk->CopyFromVTrack((AliVTrack*) trk);
+ if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
+ else {
+ fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+ if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+ }
+ delete parTrk;
}
- fPsi2sTree ->Fill();
- PostData(2, fPsi2sTree);
+ fKfVtxPos[0]= KFvtx->GetX();
+ fKfVtxPos[1]= KFvtx->GetY();
+ fKfVtxPos[2]= KFvtx->GetZ();
+ for(UInt_t i=0; i<4; i++)delete KFpart[i];
+ delete KFvtx;
+ if(!isMC) fPsi2sTree ->Fill();
+ }
+
+ if(isMC){
+ fJPsiTree ->Fill();
+ fPsi2sTree ->Fill();
}
+
+ PostData(1, fJPsiTree);
+ PostData(2, fPsi2sTree);
}//RunAOD
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcPsi2s::RunAODMC(AliAODEvent *aod)
+{
+
+ fGenPart->Clear("C");
+
+ TClonesArray *arrayMC = (TClonesArray*) aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
+ if(!arrayMC) return;
+
+ Int_t nmc=0;
+ //loop over mc particles
+ for(Int_t imc=0; imc<arrayMC->GetEntriesFast(); imc++) {
+ AliAODMCParticle *mcPart = (AliAODMCParticle*) arrayMC->At(imc);
+ if(!mcPart) continue;
+
+ if(mcPart->GetMother() >= 0) continue;
+
+ TParticle *part = (TParticle*) fGenPart->ConstructedAt(nmc++);
+ part->SetMomentum(mcPart->Px(), mcPart->Py(), mcPart->Pz(), mcPart->E());
+ part->SetPdgCode(mcPart->GetPdgCode());
+ part->SetUniqueID(imc);
+ }//loop over mc particles
+
+}//RunAODMC
+
+
//_____________________________________________________________________________
void AliAnalysisTaskUpcPsi2s::RunESDtrig()
{
//Trigger
TString trigger = esd->GetFiredTriggerClasses();
- if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
+ if(trigger.Contains("CCUP4-B")) fHistCcup4TriggersPerRun->Fill(fRunNum); //CCUP4 triggers
+ if(trigger.Contains("CCUP7-B")) fHistCcup7TriggersPerRun->Fill(fRunNum); //CCUP7 triggers
+ if(trigger.Contains("CCUP2-B")) fHistCcup2TriggersPerRun->Fill(fRunNum); //CCUP2 triggers
if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
if(esd->GetHeader()->IsTriggerInputFired("1ZED")) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
+ //MB, Central and SemiCentral triggers
+ AliCentrality *centrality = esd->GetCentrality();
+ UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
+
+ //Double_t percentile = centrality->GetCentralityPercentile("V0M");
+ Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
+
+ if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<=80 && percentile>=0) fHistMBTriggersPerRun->Fill(fRunNum);
+
+ if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<=6 && percentile>=0 && (trigger.Contains("CVHN_R2-B"))) fHistCentralTriggersPerRun->Fill(fRunNum);
+
+ if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<=50 && percentile>=15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
+
+
PostData(3, fListTrig);
}
//Trigger
TString trigger = esd->GetFiredTriggerClasses();
- if( !trigger.Contains("CCUP4-B") ) return;
+ if(!isMC && !trigger.Contains("CCUP") ) return;
fHistNeventsJPsi->Fill(2);
fHistNeventsPsi2s->Fill(2);
fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
- if( fZDCAenergy > 12000 || fZDCCenergy > 12000) return;
+ if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
fHistNeventsJPsi->Fill(4);
fHistNeventsPsi2s->Fill(4);
TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
Short_t qLepton[4], qPion[4];
UInt_t nLepton=0, nPion=0, nHighPt=0;
- Double_t jRecTPCsignal[5];
+ Double_t fRecTPCsignal[5];
Int_t mass[3]={-1,-1,-1};
- //Track loop
+ //Two Track loop
for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
AliESDtrack *trk = esd->GetTrack(itr);
if( !trk ) continue;
if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
- if(trk->GetTPCNcls() < 50)continue;
+ if(trk->GetTPCNcls() < 70)continue;
if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
+ if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
trk->GetImpactParameters(dca[0],dca[1]);
if(TMath::Abs(dca[1]) > 2) continue;
+ if(TMath::Abs(dca[1]) > 0.2) continue;
TrackIndex[nGoodTracks] = itr;
nGoodTracks++;
- if(nGoodTracks > 4) break;
+ if(nGoodTracks > 2) break;
}//Track loop
+
if(nGoodTracks == 2){
fHistNeventsJPsi->Fill(5);
for(Int_t i=0; i<2; i++){
AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
if(trk->Pt() > 1) nHighPt++;
- jRecTPCsignal[nLepton] = trk->GetTPCsignal();
+ fRecTPCsignal[nLepton] = trk->GetTPCsignal();
qLepton[nLepton] = trk->Charge();
- if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
+ if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
mass[nLepton] = 0;
}
- if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
+ if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
mass[nLepton] = 1;
}
fHistNeventsJPsi->Fill(7);
if(nHighPt > 0){
fHistNeventsJPsi->Fill(8);
- fHistTPCsignalJPsi->Fill(jRecTPCsignal[0],jRecTPCsignal[1]);
+ fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
if(mass[0] == mass[1] && mass[0] != -1) {
fHistNeventsJPsi->Fill(10);
}
}
}
- nLepton=0; nPion=0; nHighPt=0;
+ nGoodTracks = 0; nLepton=0; nPion=0; nHighPt=0;
mass[0]= -1; mass[1]= -1, mass[2]= -1;
+ //Four Track loop
+ for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
+ AliESDtrack *trk = esd->GetTrack(itr);
+ if( !trk ) continue;
+
+ if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(trk->GetTPCNcls() < 50)continue;
+ if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
+ Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
+ if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
+ trk->GetImpactParameters(dca[0],dca[1]);
+ if(TMath::Abs(dca[1]) > 2) continue;
+
+ TrackIndex[nGoodTracks] = itr;
+ nGoodTracks++;
+ if(nGoodTracks > 4) break;
+ }//Track loop
+
if(nGoodTracks == 4){
fHistNeventsPsi2s->Fill(5);
for(Int_t i=0; i<4; i++){
AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
if(trk->Pt() > 1){
- jRecTPCsignal[nLepton] = trk->GetTPCsignal();
+ fRecTPCsignal[nLepton] = trk->GetTPCsignal();
qLepton[nLepton] = trk->Charge();
- if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
+ if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
mass[nLepton] = 0;
}
- if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
+ if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
mass[nLepton] = 1;
}
//input event
AliESDEvent *esd = (AliESDEvent*) InputEvent();
if(!esd) return;
+
+ if(isMC) RunESDMC(esd);
//input data
const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
fRunNum = esd->GetRunNumber();
- //Trigger
+ //Trigger
TString trigger = esd->GetFiredTriggerClasses();
- if( !trigger.Contains("CCUP4-B") ) return;
-
+ fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
+ fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
+ fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
+
+ Bool_t isTriggered = kFALSE;
+ for(Int_t i=0; i<ntrg; i++) {
+ if( fTrigger[i] ) isTriggered = kTRUE;
+ }
+ if(!isMC && !isTriggered ) return;
+
//trigger inputs
fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
+
+ //TOF trigger info (0OMU)
+ fTOFtrig1 = esd->GetHeader()->IsTriggerInputFired("0OMU");
+ fTOFtrig2 = esd->GetHeader()->GetActiveTriggerInputs().Contains("0OMU") ? ((esd->GetHeader()) ? esd->GetHeader()->IsTriggerInputFired("0OMU") : kFALSE) : TESTBIT(esd->GetHeader()->GetL0TriggerInputs(), (kFALSE) ? 21 : 9);
//Event identification
fPerNum = esd->GetPeriodNumber();
//primary vertex
AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
fVtxContrib = fESDVertex->GetNContributors();
+ fVtxPos[0] = fESDVertex->GetX();
+ fVtxPos[1] = fESDVertex->GetY();
+ fVtxPos[2] = fESDVertex->GetZ();
+ Double_t CovMatx[6];
+ fESDVertex->GetCovarianceMatrix(CovMatx);
+ fVtxErr[0] = CovMatx[0];
+ fVtxErr[1] = CovMatx[1];
+ fVtxErr[2] = CovMatx[2];
+ fVtxChi2 = fESDVertex->GetChi2();
+ fVtxNDF = fESDVertex->GetNDF();
//Tracklets
fNtracklets = esd->GetMultiplicity()->GetNumberOfTracklets();
fV0Cdecision = fV0data->GetV0CDecision();
fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
+
+ fNLooseTracks = 0;
+
+ //Track loop - loose cuts
+ for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
+ AliESDtrack *trk = esd->GetTrack(itr);
+ if( !trk ) continue;
+
+ if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(trk->GetTPCNcls() < 20)continue;
+ fNLooseTracks++;
+ }//Track loop -loose cuts
Int_t nGoodTracks=0;
Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
- //Track loop
+ //Two Track loop
for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
AliESDtrack *trk = esd->GetTrack(itr);
if( !trk ) continue;
if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
- if(trk->GetTPCNcls() < 50)continue;
+ if(trk->GetTPCNcls() < 70)continue;
if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
+ if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
trk->GetImpactParameters(dca[0],dca[1]);
if(TMath::Abs(dca[1]) > 2) continue;
+ if(TMath::Abs(dca[1]) > 0.2) continue;
TrackIndex[nGoodTracks] = itr;
nGoodTracks++;
- if(nGoodTracks > 4) break;
+ if(nGoodTracks > 2) break;
}//Track loop
+ fJPsiESDTracks->Clear("C");
if(nGoodTracks == 2){
for(Int_t i=0; i<2; i++){
AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
new((*fJPsiESDTracks)[i]) AliESDtrack(*trk);
+ fPIDMuon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kMuon);
+ fPIDElectron[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kElectron);
+ fPIDPion[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kPion);
+
+ Double_t pos[3]={0,0,0};
+ if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
+ else {
+ fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+ if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+ }
}
- fJPsiTree ->Fill();
- PostData(1, fJPsiTree);
+ if(!isMC) fJPsiTree ->Fill();
}
+ nGoodTracks = 0;
+ //Four track loop
+ for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
+ AliESDtrack *trk = esd->GetTrack(itr);
+ if( !trk ) continue;
+
+ if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+ if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+ if(trk->GetTPCNcls() < 50)continue;
+ if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
+ Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
+ if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
+ trk->GetImpactParameters(dca[0],dca[1]);
+ if(TMath::Abs(dca[1]) > 2) continue;
+
+ TrackIndex[nGoodTracks] = itr;
+ nGoodTracks++;
+ if(nGoodTracks > 4) break;
+ }//Track loop
+
+ fPsi2sESDTracks->Clear("C");
if(nGoodTracks == 4){
for(Int_t i=0; i<4; i++){
AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
AliExternalTrackParam cParam;
trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
- new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
+ new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
+ fPIDMuon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kMuon);
+ fPIDElectron[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kElectron);
+ fPIDPion[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kPion);
+
+ Double_t pos[3]={0,0,0};
+ if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
+ else {
+ fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+ if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+ }
}
- fPsi2sTree ->Fill();
- PostData(2, fPsi2sTree);
+ if(!isMC) fPsi2sTree ->Fill();
}
+
+ if(isMC){
+ fJPsiTree ->Fill();
+ fPsi2sTree ->Fill();
+ }
+
+ PostData(1, fJPsiTree);
+ PostData(2, fPsi2sTree);
}//RunESD
+
//_____________________________________________________________________________
-void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
+void AliAnalysisTaskUpcPsi2s::RunESDMC(AliESDEvent* esd)
{
- cout<<"Analysis complete."<<endl;
-}//Terminate
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
+ AliTriggerAnalysis *fTrigAna = new AliTriggerAnalysis();
+ fTrigAna->SetAnalyzeMC(isMC);
+
+ if(fTrigAna->SPDFiredChips(esd,1,kFALSE,2) > 1) fTriggerInputsMC[0] = kTRUE;
+ if(fTrigAna->SPDFiredChips(esd,1,kFALSE,2) < 2) fTriggerInputsMC[0] = kFALSE;
+ fTriggerInputsMC[1] = esd->GetHeader()->IsTriggerInputFired("0OMU");
+ fTriggerInputsMC[2] = esd->GetHeader()->IsTriggerInputFired("0VBA");
+ fTriggerInputsMC[3] = esd->GetHeader()->IsTriggerInputFired("0VBC");
+ fGenPart->Clear("C");
+ AliMCEvent *mc = MCEvent();
+ if(!mc) return;
+ Int_t nmc = 0;
+ //loop over mc particles
+ for(Int_t imc=0; imc<mc->GetNumberOfTracks(); imc++) {
+ AliMCParticle *mcPart = (AliMCParticle*) mc->GetTrack(imc);
+ if(!mcPart) continue;
+ if(mcPart->GetMother() >= 0) continue;
+ TParticle *part = (TParticle*) fGenPart->ConstructedAt(nmc++);
+ part->SetMomentum(mcPart->Px(), mcPart->Py(), mcPart->Pz(), mcPart->E());
+ part->SetPdgCode(mcPart->PdgCode());
+ part->SetUniqueID(imc);
+ }//loop over mc particles
+}//RunESDMC
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
+{
+ cout<<"Analysis complete."<<endl;
+}//Terminate
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff