1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
23 #include "TClonesArray.h"
24 #include "TParticle.h"
25 #include "TObjString.h"
27 #include "TDatabasePDG.h"
28 #include "TLorentzVector.h"
31 #include "AliAnalysisManager.h"
32 #include "AliInputEventHandler.h"
33 #include "AliESDEvent.h"
34 #include "AliAODEvent.h"
35 #include "AliMCEvent.h"
36 #include "AliAODVZERO.h"
37 #include "AliAODZDC.h"
38 #include "AliESDVZERO.h"
39 #include "AliESDZDC.h"
40 #include "AliPIDResponse.h"
41 #include "AliAODTrack.h"
42 #include "AliAODPid.h"
43 #include "AliAODVertex.h"
44 #include "AliESDVertex.h"
45 #include "AliMultiplicity.h"
46 #include "AliESDtrack.h"
47 #include "AliESDMuonTrack.h"
48 #include "AliAODMCParticle.h"
49 #include "AliMCParticle.h"
50 #include "AliCentrality.h"
53 #include "AliAnalysisTaskUpcPsi2s.h"
55 ClassImp(AliAnalysisTaskUpcPsi2s);
60 //trees for UPC analysis,
61 // michal.broz@cern.ch
63 //_____________________________________________________________________________
64 AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s()
65 : AliAnalysisTaskSE(),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
66 fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
67 fTOFtrig1(0), fTOFtrig2(0),
68 fVtxContrib(0),fVtxPosX(0),fVtxPosY(0),fVtxPosZ(0),fVtxErrX(0),fVtxErrY(0),fVtxErrZ(0),fVtxChi2(0),fVtxNDF(0),
69 fBCrossNum(0),fNtracklets(0),
70 fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
71 fDataFilnam(0),fRecoPass(0),fEvtNum(0),
72 fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
73 fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
74 fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
75 fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
76 fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
82 }//AliAnalysisTaskUpcPsi2s
85 //_____________________________________________________________________________
86 AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s(const char *name)
87 : AliAnalysisTaskSE(name),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
88 fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
89 fTOFtrig1(0), fTOFtrig2(0),
90 fVtxContrib(0),fVtxPosX(0),fVtxPosY(0),fVtxPosZ(0),fVtxErrX(0),fVtxErrY(0),fVtxErrZ(0),fVtxChi2(0),fVtxNDF(0),
91 fBCrossNum(0),fNtracklets(0),
92 fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
93 fDataFilnam(0),fRecoPass(0),fEvtNum(0),
94 fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
95 fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
96 fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
97 fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
98 fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
103 if( strstr(name,"ESD") ) fType = 0;
104 if( strstr(name,"AOD") ) fType = 1;
108 DefineOutput(1, TTree::Class());
109 DefineOutput(2, TTree::Class());
110 DefineOutput(3, TList::Class());
111 DefineOutput(4, TList::Class());
113 }//AliAnalysisTaskUpcPsi2s
115 //_____________________________________________________________________________
116 void AliAnalysisTaskUpcPsi2s::Init()
119 for(Int_t i=0; i<ntrg; i++) fTrigger[i] = kFALSE;
120 for(Int_t i=0; i<4; i++) fTOFphi[i] = -666;
124 //_____________________________________________________________________________
125 AliAnalysisTaskUpcPsi2s::~AliAnalysisTaskUpcPsi2s()
145 }//~AliAnalysisTaskUpcPsi2s
148 //_____________________________________________________________________________
149 void AliAnalysisTaskUpcPsi2s::UserCreateOutputObjects()
152 fDataFilnam = new TObjString();
153 fDataFilnam->SetString("");
156 fJPsiAODTracks = new TClonesArray("AliAODTrack", 1000);
157 fJPsiESDTracks = new TClonesArray("AliESDtrack", 1000);
158 fPsi2sAODTracks = new TClonesArray("AliAODTrack", 1000);
159 fPsi2sESDTracks = new TClonesArray("AliESDtrack", 1000);
161 //output tree with JPsi candidate events
162 fJPsiTree = new TTree("fJPsiTree", "fJPsiTree");
163 fJPsiTree ->Branch("fRunNum", &fRunNum, "fRunNum/I");
164 fJPsiTree ->Branch("fPerNum", &fPerNum, "fPerNum/i");
165 fJPsiTree ->Branch("fOrbNum", &fOrbNum, "fOrbNum/i");
167 fJPsiTree ->Branch("fBCrossNum", &fBCrossNum, "fBCrossNum/s");
168 fJPsiTree ->Branch("fTrigger", &fTrigger[0], Form("fTrigger[%i]/O", ntrg));
169 fJPsiTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
170 fJPsiTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
171 fJPsiTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
172 fJPsiTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
174 fJPsiTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
175 fJPsiTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
176 fJPsiTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
178 fJPsiTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
179 fJPsiTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
180 fJPsiTree ->Branch("fVtxPosZ", &fVtxPosZ, "fVtxPosZ/D");
181 fJPsiTree ->Branch("fVtxErrX", &fVtxErrX, "fVtxErrX/D");
182 fJPsiTree ->Branch("fVtxErrY", &fVtxErrY, "fVtxErrY/D");
183 fJPsiTree ->Branch("fVtxErrZ", &fVtxErrZ, "fVtxErrZ/D");
184 fJPsiTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
185 fJPsiTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
187 fJPsiTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
188 fJPsiTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
189 fJPsiTree ->Branch("fV0Adecision", &fV0Adecision, "fV0Adecision/I");
190 fJPsiTree ->Branch("fV0Cdecision", &fV0Cdecision, "fV0Cdecision/I");
191 fJPsiTree ->Branch("fDataFilnam", &fDataFilnam);
192 fJPsiTree ->Branch("fRecoPass", &fRecoPass, "fRecoPass/S");
193 fJPsiTree ->Branch("fEvtNum", &fEvtNum, "fEvtNum/L");
195 fJPsiTree ->Branch("fJPsiESDTracks", &fJPsiESDTracks);
198 fJPsiTree ->Branch("fJPsiAODTracks", &fJPsiAODTracks);
201 //output tree with Psi2s candidate events
202 fPsi2sTree = new TTree("fPsi2sTree", "fPsi2sTree");
203 fPsi2sTree ->Branch("fRunNum", &fRunNum, "fRunNum/I");
204 fPsi2sTree ->Branch("fPerNum", &fPerNum, "fPerNum/i");
205 fPsi2sTree ->Branch("fOrbNum", &fOrbNum, "fOrbNum/i");
207 fPsi2sTree ->Branch("fBCrossNum", &fBCrossNum, "fBCrossNum/s");
208 fPsi2sTree ->Branch("fTrigger", &fTrigger[0], Form("fTrigger[%i]/O", ntrg));
209 fPsi2sTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
210 fPsi2sTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
211 fPsi2sTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
212 fPsi2sTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
214 fPsi2sTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
215 fPsi2sTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
216 fPsi2sTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
218 fPsi2sTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
219 fPsi2sTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
220 fPsi2sTree ->Branch("fVtxPosZ", &fVtxPosZ, "fVtxPosZ/D");
221 fPsi2sTree ->Branch("fVtxErrX", &fVtxErrX, "fVtxErrX/D");
222 fPsi2sTree ->Branch("fVtxErrY", &fVtxErrY, "fVtxErrY/D");
223 fPsi2sTree ->Branch("fVtxErrZ", &fVtxErrZ, "fVtxErrZ/D");
224 fPsi2sTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
225 fPsi2sTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
227 fPsi2sTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
228 fPsi2sTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
229 fPsi2sTree ->Branch("fV0Adecision", &fV0Adecision, "fV0Adecision/I");
230 fPsi2sTree ->Branch("fV0Cdecision", &fV0Cdecision, "fV0Cdecision/I");
231 fPsi2sTree ->Branch("fDataFilnam", &fDataFilnam);
232 fPsi2sTree ->Branch("fRecoPass", &fRecoPass, "fRecoPass/S");
233 fPsi2sTree ->Branch("fEvtNum", &fEvtNum, "fEvtNum/L");
235 fPsi2sTree ->Branch("fPsi2sESDTracks", &fPsi2sESDTracks);
238 fPsi2sTree ->Branch("fPsi2sAODTracks", &fPsi2sAODTracks);
241 fListTrig = new TList();
242 fListTrig ->SetOwner();
244 fHistUpcTriggersPerRun = new TH1D("fHistUpcTriggersPerRun", "fHistUpcTriggersPerRun", 3000, 167000.5, 170000.5);
245 fListTrig->Add(fHistUpcTriggersPerRun);
247 fHistZedTriggersPerRun = new TH1D("fHistZedTriggersPerRun", "fHistZedTriggersPerRun", 3000, 167000.5, 170000.5);
248 fListTrig->Add(fHistZedTriggersPerRun);
250 fHistCvlnTriggersPerRun = new TH1D("fHistCvlnTriggersPerRun", "fHistCvlnTriggersPerRun", 3000, 167000.5, 170000.5);
251 fListTrig->Add(fHistCvlnTriggersPerRun);
253 fHistMBTriggersPerRun = new TH1D("fHistMBTriggersPerRun", "fHistMBTriggersPerRun", 3000, 167000.5, 170000.5);
254 fListTrig->Add(fHistMBTriggersPerRun);
256 fHistCentralTriggersPerRun = new TH1D("fHistCentralTriggersPerRun", "fHistCentralTriggersPerRun", 3000, 167000.5, 170000.5);
257 fListTrig->Add(fHistCentralTriggersPerRun);
259 fHistSemiCentralTriggersPerRun = new TH1D("fHistSemiCentralTriggersPerRun", "fHistSemiCentralTriggersPerRun", 3000, 167000.5, 170000.5);
260 fListTrig->Add(fHistSemiCentralTriggersPerRun);
262 fListHist = new TList();
263 fListHist ->SetOwner();
265 TString CutNameJPsi[12] = {"Analyzed","Triggered","Vertex cut","V0 decision","Two good tracks",
266 "Like sign","Oposite sign","One p_{T}>1", "Both p_{T}>1","PID","Dimuom","Dielectron"};
267 fHistNeventsJPsi = new TH1D("fHistNeventsJPsi","fHistNeventsPsi2s",12,0.5,12.5);
268 for (Int_t i = 0; i<12; i++) fHistNeventsJPsi->GetXaxis()->SetBinLabel(i+1,CutNameJPsi[i].Data());
269 fListHist->Add(fHistNeventsJPsi);
271 fHistTPCsignalJPsi = new TH2D("fHistTPCsignalJPsi","fHistTPCsignalJPsi",240,0,120,240,0,120);
272 fListHist->Add(fHistTPCsignalJPsi);
274 fHistDiLeptonPtJPsi = new TH2D("fHistDiLeptonPtJPsi","fHistDiLeptonPtJPsi",350,0,3.5,350,0,3.5);
275 fListHist->Add(fHistDiLeptonPtJPsi);
277 fHistDiElectronMass = new TH1D("fHistDiElectronMass","Invariant mass of J/#psi candidates",100,2,5);
278 fHistDiElectronMass->GetXaxis()->SetTitle("Invariant mass(e^{+}e^{-}) (GeV/c)");
279 fListHist->Add(fHistDiElectronMass);
281 fHistDiMuonMass = new TH1D("fHistDiMuonMass","Invariant mass of J/#psi candidates",100,2,5);
282 fHistDiMuonMass->GetXaxis()->SetTitle("Invariant mass(#mu^{+}#mu^{-}) (GeV/c)");
283 fListHist->Add(fHistDiMuonMass);
285 TString CutNamePsi2s[13] = {"Analyzed","Triggered","Vertex cut","V0 decision","Four good tracks",
286 "DiLepton - DiPion","Like sign leptons","Like sign pions","Like sign both","Oposite sign","PID","Dimuom","Dielectron"};
288 fHistNeventsPsi2s = new TH1D("fHistNeventsPsi2s","fHistNeventsPsi2s",13,0.5,13.5);
289 for (Int_t i = 0; i<13; i++) fHistNeventsPsi2s->GetXaxis()->SetBinLabel(i+1,CutNamePsi2s[i].Data());
290 fListHist->Add(fHistNeventsPsi2s);
292 fHistPsi2sMassVsPt = new TH2D("fHistPsi2sMassVsPt","Mass vs p_{T} of #psi(2s) candidates",100,3,6,50,0,5);
293 fHistPsi2sMassVsPt->GetXaxis()->SetTitle("Invariant mass(l^{+}l^{-}#pi^{+}#pi^{-}) (GeV/c)");
294 fHistPsi2sMassVsPt->GetYaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
295 fListHist->Add(fHistPsi2sMassVsPt);
297 fHistPsi2sMassCoherent = new TH1D("fHistPsi2sMassAllCoherent","Invariant mass of coherent #psi(2s) candidates",100,3,6);
298 fHistPsi2sMassCoherent->GetXaxis()->SetTitle("Invariant mass(l^{+}l^{-}#pi^{+}#pi^{-}) (GeV/c)");
299 fListHist->Add(fHistPsi2sMassCoherent);
301 PostData(1, fJPsiTree);
302 PostData(2, fPsi2sTree);
303 PostData(3, fListTrig);
304 PostData(4, fListHist);
306 }//UserCreateOutputObjects
308 //_____________________________________________________________________________
309 void AliAnalysisTaskUpcPsi2s::UserExec(Option_t *)
312 //cout<<"#################### Next event ##################"<<endl;
316 if(fRunHist) RunESDhist();
317 if(fRunTree) RunESDtree();
322 if(fRunHist) RunAODhist();
323 if(fRunTree) RunAODtree();
327 //_____________________________________________________________________________
328 void AliAnalysisTaskUpcPsi2s::RunAODtrig()
332 AliAODEvent *aod = (AliAODEvent*) InputEvent();
335 fRunNum = aod ->GetRunNumber();
337 TString trigger = aod->GetFiredTriggerClasses();
339 if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
341 if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
342 if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
344 fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
345 if(fL1inputs & (1 << 18)) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
347 //MB, Central and SemiCentral triggers
348 AliCentrality *centrality = aod->GetCentrality();
349 UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
351 Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
352 //Double_t percentile = centrality->GetCentralityPercentile("V0M");
354 if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
356 if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
358 if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
360 PostData(3, fListTrig);
363 //_____________________________________________________________________________
364 void AliAnalysisTaskUpcPsi2s::RunAODhist()
367 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
369 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
370 Double_t muonMass = partMuon->Mass();
372 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
373 Double_t electronMass = partElectron->Mass();
375 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
376 Double_t pionMass = partPion->Mass();
379 AliAODEvent *aod = (AliAODEvent*) InputEvent();
382 fHistNeventsJPsi->Fill(1);
383 fHistNeventsPsi2s->Fill(1);
386 TString trigger = aod->GetFiredTriggerClasses();
388 if( !trigger.Contains("CCUP") ) return;
390 fHistNeventsJPsi->Fill(2);
391 fHistNeventsPsi2s->Fill(2);
394 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
395 fVtxContrib = fAODVertex->GetNContributors();
396 if(fVtxContrib < 2) return;
398 fHistNeventsJPsi->Fill(3);
399 fHistNeventsPsi2s->Fill(3);
402 AliAODVZERO *fV0data = aod ->GetVZEROData();
403 AliAODZDC *fZDCdata = aod->GetZDCData();
405 fV0Adecision = fV0data->GetV0ADecision();
406 fV0Cdecision = fV0data->GetV0CDecision();
407 if(fV0Adecision != AliAODVZERO::kV0Empty || fV0Cdecision != AliAODVZERO::kV0Empty) return;
409 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
410 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
412 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
414 fHistNeventsJPsi->Fill(4);
415 fHistNeventsPsi2s->Fill(4);
418 Int_t nGoodTracks = 0;
419 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
421 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
422 Short_t qLepton[4], qPion[4];
423 UInt_t nLepton=0, nPion=0, nHighPt=0;
424 Double_t jRecTPCsignal[5];
425 Int_t mass[3]={-1,-1,-1};
429 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
430 AliAODTrack *trk = aod->GetTrack(itr);
433 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
434 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
435 if(trk->GetTPCNcls() < 50)continue;
436 if(trk->Chi2perNDF() > 4)continue;
437 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
438 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
439 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
442 if(TMath::Abs(dca[1]) > 2) continue;
444 TrackIndex[nGoodTracks] = itr;
447 if(nGoodTracks > 4) break;
450 nLepton=0; nPion=0; nHighPt=0;
451 mass[0]= -1; mass[1]= -1, mass[2]= -1;
453 if(nGoodTracks == 4){
454 fHistNeventsPsi2s->Fill(5);
455 for(Int_t i=0; i<4; i++){
456 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
459 jRecTPCsignal[nLepton] = trk->GetTPCsignal();
460 qLepton[nLepton] = trk->Charge();
461 if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
462 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
465 if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
466 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
472 qPion[nPion] = trk->Charge();
473 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
477 if(nLepton > 2 || nPion > 2) break;
479 if((nLepton == 2) && (nPion == 2)){
480 fHistNeventsPsi2s->Fill(6);
481 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
482 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
483 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
484 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
485 fHistNeventsPsi2s->Fill(10);
486 if(mass[0] == mass[1]) {
487 fHistNeventsPsi2s->Fill(11);
488 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
489 vDilepton = vLepton[0]+vLepton[1];
490 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
491 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
492 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
493 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
501 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
502 AliAODTrack *trk = aod->GetTrack(itr);
505 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
506 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
507 if(trk->GetTPCNcls() < 70)continue;
508 if(trk->Chi2perNDF() > 4)continue;
509 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
510 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
511 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
512 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
514 if(TMath::Abs(dca[1]) > 2) continue;
515 if(TMath::Abs(dca[0]) > 0.2) continue;
517 TrackIndex[nGoodTracks] = itr;
520 if(nGoodTracks > 2) break;
523 nLepton=0; nPion=0; nHighPt=0;
524 mass[0]= -1; mass[1]= -1, mass[2]= -1;
526 if(nGoodTracks == 2){
527 fHistNeventsJPsi->Fill(5);
528 for(Int_t i=0; i<2; i++){
529 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
530 if(trk->Pt() > 1) nHighPt++;
531 jRecTPCsignal[nLepton] = trk->GetTPCsignal();
532 qLepton[nLepton] = trk->Charge();
533 if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
534 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
537 if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
538 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
544 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
545 if(qLepton[0]*qLepton[1] < 0){
546 fHistNeventsJPsi->Fill(7);
548 fHistNeventsJPsi->Fill(8);
549 fHistTPCsignalJPsi->Fill(jRecTPCsignal[0],jRecTPCsignal[1]);
550 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
551 if(mass[0] == mass[1] && mass[0] != -1) {
552 fHistNeventsJPsi->Fill(10);
553 vCandidate = vLepton[0]+vLepton[1];
554 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
556 fHistDiMuonMass->Fill(vCandidate.M());
557 fHistNeventsJPsi->Fill(11);
560 fHistDiElectronMass->Fill(vCandidate.M());
561 fHistNeventsJPsi->Fill(12);
570 PostData(4, fListHist);
574 //_____________________________________________________________________________
575 void AliAnalysisTaskUpcPsi2s::RunAODtree()
578 AliAODEvent *aod = (AliAODEvent*) InputEvent();
582 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
583 fDataFilnam->Clear();
584 fDataFilnam->SetString(filnam);
585 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
586 fRunNum = aod ->GetRunNumber();
589 TString trigger = aod->GetFiredTriggerClasses();
591 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
592 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
593 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
595 Bool_t isTriggered = kFALSE;
596 for(Int_t i=0; i<ntrg; i++) {
597 if( fTrigger[i] ) isTriggered = kTRUE;
599 if( !isTriggered ) return;
602 fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
603 fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
605 //Event identification
606 fPerNum = aod ->GetPeriodNumber();
607 fOrbNum = aod ->GetOrbitNumber();
608 fBCrossNum = aod ->GetBunchCrossNumber();
611 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
612 fVtxContrib = fAODVertex->GetNContributors();
613 fVtxPosX = fAODVertex->GetX();
614 fVtxPosY = fAODVertex->GetY();
615 fVtxPosZ = fAODVertex->GetZ();
617 fAODVertex->GetCovarianceMatrix(CovMatx);
618 fVtxErrX = CovMatx[0];
619 fVtxErrY = CovMatx[1];
620 fVtxErrZ = CovMatx[2];
621 fVtxChi2 = fAODVertex->GetChi2();
622 fVtxNDF = fAODVertex->GetNDF();
625 fNtracklets = aod->GetTracklets()->GetNumberOfTracklets();
628 AliAODVZERO *fV0data = aod ->GetVZEROData();
629 AliAODZDC *fZDCdata = aod->GetZDCData();
631 fV0Adecision = fV0data->GetV0ADecision();
632 fV0Cdecision = fV0data->GetV0CDecision();
633 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
634 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
637 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
640 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
641 AliAODTrack *trk = aod->GetTrack(itr);
644 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
645 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
646 if(trk->GetTPCNcls() < 70)continue;
647 if(trk->Chi2perNDF() > 4)continue;
648 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
649 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
650 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
651 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
653 if(TMath::Abs(dca[1]) > 2) continue;
654 if(TMath::Abs(dca[0]) > 0.2) continue;
656 TrackIndex[nGoodTracks] = itr;
659 if(nGoodTracks > 2) break;
662 if(nGoodTracks == 2){
663 for(Int_t i=0; i<2; i++){
664 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
666 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
667 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
668 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
671 new((*fJPsiAODTracks)[i]) AliAODTrack(*trk);
672 ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
674 Double_t pos[3]={0,0,0};
675 if(!trk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
677 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
678 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
686 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
687 AliAODTrack *trk = aod->GetTrack(itr);
690 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
691 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
692 if(trk->GetTPCNcls() < 50)continue;
693 if(trk->Chi2perNDF() > 4)continue;
694 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
695 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
696 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
698 if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
699 if(TMath::Abs(dca[1]) > 2) continue;
701 TrackIndex[nGoodTracks] = itr;
704 if(nGoodTracks > 4) break;
707 if(nGoodTracks == 4){
708 for(Int_t i=0; i<4; i++){
709 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
711 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
712 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
713 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
716 new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
717 ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
719 Double_t pos[3]={0,0,0};
720 if(!trk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
722 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
723 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
729 PostData(1, fJPsiTree);
730 PostData(2, fPsi2sTree);
734 //_____________________________________________________________________________
735 void AliAnalysisTaskUpcPsi2s::RunESDtrig()
739 AliESDEvent *esd = (AliESDEvent*) InputEvent();
742 fRunNum = esd ->GetRunNumber();
744 TString trigger = esd->GetFiredTriggerClasses();
746 if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
748 if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
749 if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
751 if(esd->GetHeader()->IsTriggerInputFired("1ZED")) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
753 //MB, Central and SemiCentral triggers
754 AliCentrality *centrality = esd->GetCentrality();
755 UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
757 //Double_t percentile = centrality->GetCentralityPercentile("V0M");
758 Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
760 if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
762 if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
764 if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
767 PostData(3, fListTrig);
770 //_____________________________________________________________________________
771 void AliAnalysisTaskUpcPsi2s::RunESDhist()
774 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
776 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
777 Double_t muonMass = partMuon->Mass();
779 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
780 Double_t electronMass = partElectron->Mass();
782 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
783 Double_t pionMass = partPion->Mass();
786 AliESDEvent *esd = (AliESDEvent*) InputEvent();
789 fHistNeventsJPsi->Fill(1);
790 fHistNeventsPsi2s->Fill(1);
793 TString trigger = esd->GetFiredTriggerClasses();
795 if( !trigger.Contains("CCUP") ) return;
797 fHistNeventsJPsi->Fill(2);
798 fHistNeventsPsi2s->Fill(2);
801 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
802 fVtxContrib = fESDVertex->GetNContributors();
803 if(fVtxContrib < 2) return;
805 fHistNeventsJPsi->Fill(3);
806 fHistNeventsPsi2s->Fill(3);
809 AliESDVZERO *fV0data = esd->GetVZEROData();
810 AliESDZDC *fZDCdata = esd->GetESDZDC();
812 fV0Adecision = fV0data->GetV0ADecision();
813 fV0Cdecision = fV0data->GetV0CDecision();
814 if(fV0Adecision != AliESDVZERO::kV0Empty || fV0Cdecision != AliESDVZERO::kV0Empty) return;
816 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
817 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
818 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
820 fHistNeventsJPsi->Fill(4);
821 fHistNeventsPsi2s->Fill(4);
825 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
827 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
828 Short_t qLepton[4], qPion[4];
829 UInt_t nLepton=0, nPion=0, nHighPt=0;
830 Double_t jRecTPCsignal[5];
831 Int_t mass[3]={-1,-1,-1};
834 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
835 AliESDtrack *trk = esd->GetTrack(itr);
838 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
839 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
840 if(trk->GetTPCNcls() < 70)continue;
841 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
842 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
843 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
844 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
845 trk->GetImpactParameters(dca[0],dca[1]);
846 if(TMath::Abs(dca[1]) > 2) continue;
847 if(TMath::Abs(dca[1]) > 0.2) continue;
849 TrackIndex[nGoodTracks] = itr;
851 if(nGoodTracks > 2) break;
855 if(nGoodTracks == 2){
856 fHistNeventsJPsi->Fill(5);
857 for(Int_t i=0; i<2; i++){
858 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
859 if(trk->Pt() > 1) nHighPt++;
860 jRecTPCsignal[nLepton] = trk->GetTPCsignal();
861 qLepton[nLepton] = trk->Charge();
862 if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
863 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
866 if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
867 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
873 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
874 if(qLepton[0]*qLepton[1] < 0){
875 fHistNeventsJPsi->Fill(7);
877 fHistNeventsJPsi->Fill(8);
878 fHistTPCsignalJPsi->Fill(jRecTPCsignal[0],jRecTPCsignal[1]);
879 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
880 if(mass[0] == mass[1] && mass[0] != -1) {
881 fHistNeventsJPsi->Fill(10);
882 vCandidate = vLepton[0]+vLepton[1];
883 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
885 fHistDiMuonMass->Fill(vCandidate.M());
886 fHistNeventsJPsi->Fill(11);
889 fHistDiElectronMass->Fill(vCandidate.M());
890 fHistNeventsJPsi->Fill(12);
897 nGoodTracks = 0; nLepton=0; nPion=0; nHighPt=0;
898 mass[0]= -1; mass[1]= -1, mass[2]= -1;
901 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
902 AliESDtrack *trk = esd->GetTrack(itr);
905 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
906 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
907 if(trk->GetTPCNcls() < 50)continue;
908 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
909 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
910 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
911 trk->GetImpactParameters(dca[0],dca[1]);
912 if(TMath::Abs(dca[1]) > 2) continue;
914 TrackIndex[nGoodTracks] = itr;
916 if(nGoodTracks > 4) break;
919 if(nGoodTracks == 4){
920 fHistNeventsPsi2s->Fill(5);
921 for(Int_t i=0; i<4; i++){
922 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
925 jRecTPCsignal[nLepton] = trk->GetTPCsignal();
926 qLepton[nLepton] = trk->Charge();
927 if(jRecTPCsignal[nLepton] > 40 && jRecTPCsignal[nLepton] < 70){
928 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
931 if(jRecTPCsignal[nLepton] > 70 && jRecTPCsignal[nLepton] < 100){
932 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
938 qPion[nPion] = trk->Charge();
939 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
943 if(nLepton > 2 || nPion > 2) break;
945 if((nLepton == 2) && (nPion == 2)){
946 fHistNeventsPsi2s->Fill(6);
947 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
948 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
949 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
950 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
951 fHistNeventsPsi2s->Fill(10);
952 if(mass[0] == mass[1]) {
953 fHistNeventsPsi2s->Fill(11);
954 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
955 vDilepton = vLepton[0]+vLepton[1];
956 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
957 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
958 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
959 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
965 PostData(4, fListHist);
969 //_____________________________________________________________________________
970 void AliAnalysisTaskUpcPsi2s::RunESDtree()
974 AliESDEvent *esd = (AliESDEvent*) InputEvent();
978 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
979 fDataFilnam->Clear();
980 fDataFilnam->SetString(filnam);
981 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
982 fRunNum = esd->GetRunNumber();
985 TString trigger = esd->GetFiredTriggerClasses();
987 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
988 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
989 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
991 Bool_t isTriggered = kFALSE;
992 for(Int_t i=0; i<ntrg; i++) {
993 if( fTrigger[i] ) isTriggered = kTRUE;
995 if( !isTriggered ) return;
998 fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
999 fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
1001 //TOF trigger info (0OMU)
1002 fTOFtrig1 = esd->GetHeader()->IsTriggerInputFired("0OMU");
1003 fTOFtrig2 = esd->GetHeader()->GetActiveTriggerInputs().Contains("0OMU") ? ((esd->GetHeader()) ? esd->GetHeader()->IsTriggerInputFired("0OMU") : kFALSE) : TESTBIT(esd->GetHeader()->GetL0TriggerInputs(), (kFALSE) ? 21 : 9);
1005 //Event identification
1006 fPerNum = esd->GetPeriodNumber();
1007 fOrbNum = esd->GetOrbitNumber();
1008 fBCrossNum = esd->GetBunchCrossNumber();
1011 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
1012 fVtxContrib = fESDVertex->GetNContributors();
1013 fVtxPosX = fESDVertex->GetX();
1014 fVtxPosY = fESDVertex->GetY();
1015 fVtxPosZ = fESDVertex->GetZ();
1016 Double_t CovMatx[6];
1017 fESDVertex->GetCovarianceMatrix(CovMatx);
1018 fVtxErrX = CovMatx[0];
1019 fVtxErrY = CovMatx[1];
1020 fVtxErrZ = CovMatx[2];
1021 fVtxChi2 = fESDVertex->GetChi2();
1022 fVtxNDF = fESDVertex->GetNDF();
1025 fNtracklets = esd->GetMultiplicity()->GetNumberOfTracklets();
1028 AliESDVZERO *fV0data = esd->GetVZEROData();
1029 AliESDZDC *fZDCdata = esd->GetESDZDC();
1031 fV0Adecision = fV0data->GetV0ADecision();
1032 fV0Cdecision = fV0data->GetV0CDecision();
1033 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
1034 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
1036 Int_t nGoodTracks=0;
1037 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
1040 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1041 AliESDtrack *trk = esd->GetTrack(itr);
1042 if( !trk ) continue;
1044 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1045 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1046 if(trk->GetTPCNcls() < 70)continue;
1047 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1048 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
1049 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1050 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1051 trk->GetImpactParameters(dca[0],dca[1]);
1052 if(TMath::Abs(dca[1]) > 2) continue;
1053 if(TMath::Abs(dca[1]) > 0.2) continue;
1055 TrackIndex[nGoodTracks] = itr;
1057 if(nGoodTracks > 2) break;
1060 if(nGoodTracks == 2){
1061 for(Int_t i=0; i<2; i++){
1062 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1064 AliExternalTrackParam cParam;
1065 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1067 new((*fJPsiESDTracks)[i]) AliESDtrack(*trk);
1069 Double_t pos[3]={0,0,0};
1070 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1072 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1073 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1081 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1082 AliESDtrack *trk = esd->GetTrack(itr);
1083 if( !trk ) continue;
1085 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1086 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1087 if(trk->GetTPCNcls() < 50)continue;
1088 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1089 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1090 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1091 trk->GetImpactParameters(dca[0],dca[1]);
1092 if(TMath::Abs(dca[1]) > 2) continue;
1094 TrackIndex[nGoodTracks] = itr;
1096 if(nGoodTracks > 4) break;
1099 if(nGoodTracks == 4){
1100 for(Int_t i=0; i<4; i++){
1101 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1103 AliExternalTrackParam cParam;
1104 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1106 new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
1108 Double_t pos[3]={0,0,0};
1109 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1111 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1112 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1115 fPsi2sTree ->Fill();
1118 PostData(1, fJPsiTree);
1119 PostData(2, fPsi2sTree);
1123 //_____________________________________________________________________________
1124 void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
1127 cout<<"Analysis complete."<<endl;