]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGUD/UPC/AliAnalysisTaskUpcPsi2s.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding Test bit 0 in AOD
[u/mrichter/AliRoot.git] / PWGUD / UPC / AliAnalysisTaskUpcPsi2s.cxx
CommitLineData
3d16cd00 1/*************************************************************************
2* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
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**************************************************************************/
15
16// c++ headers
17#include <iostream>
18#include <string.h>
19
20// root headers
21#include "TH1I.h"
22#include "TTree.h"
23#include "TClonesArray.h"
24#include "TParticle.h"
25#include "TObjString.h"
26#include "TFile.h"
17c65770 27#include "TDatabasePDG.h"
28#include "TLorentzVector.h"
3d16cd00 29
30// aliroot headers
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"
e63dc644 50#include "AliCentrality.h"
489951cd 51#include "AliKFVertex.h"
52#include "AliExternalTrackParam.h"
3d16cd00 53
54// my headers
55#include "AliAnalysisTaskUpcPsi2s.h"
56
57ClassImp(AliAnalysisTaskUpcPsi2s);
58
b8f92f9d 59using std::cout;
60using std::endl;
61
3d16cd00 62//trees for UPC analysis,
63// michal.broz@cern.ch
64
65//_____________________________________________________________________________
66AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s()
f052ef6f 67 : AliAnalysisTaskSE(),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
d457ec9e 68 fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
03cff2a0 69 fTOFtrig1(0), fTOFtrig2(0),
489951cd 70 fVtxContrib(0),fVtxChi2(0),fVtxNDF(0),
71 fBCrossNum(0),fNtracklets(0),fNLooseTracks(0),
3d16cd00 72 fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
73 fDataFilnam(0),fRecoPass(0),fEvtNum(0),
17c65770 74 fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
f052ef6f 75 fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
e63dc644 76 fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
f052ef6f 77 fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
78 fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
17c65770 79
3d16cd00 80{
81
82//Dummy constructor
83
84}//AliAnalysisTaskUpcPsi2s
85
86
87//_____________________________________________________________________________
88AliAnalysisTaskUpcPsi2s::AliAnalysisTaskUpcPsi2s(const char *name)
f052ef6f 89 : AliAnalysisTaskSE(name),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),fJPsiTree(0),fPsi2sTree(0),
d457ec9e 90 fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),
03cff2a0 91 fTOFtrig1(0), fTOFtrig2(0),
489951cd 92 fVtxContrib(0),fVtxChi2(0),fVtxNDF(0),
93 fBCrossNum(0),fNtracklets(0),fNLooseTracks(0),
3d16cd00 94 fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
95 fDataFilnam(0),fRecoPass(0),fEvtNum(0),
17c65770 96 fJPsiAODTracks(0),fJPsiESDTracks(0),fPsi2sAODTracks(0),fPsi2sESDTracks(0),
f052ef6f 97 fListTrig(0),fHistUpcTriggersPerRun(0),fHistZedTriggersPerRun(0),fHistCvlnTriggersPerRun(0),
e63dc644 98 fHistMBTriggersPerRun(0),fHistCentralTriggersPerRun(0),fHistSemiCentralTriggersPerRun(0),
f052ef6f 99 fListHist(0),fHistNeventsJPsi(0),fHistTPCsignalJPsi(0),fHistDiLeptonPtJPsi(0),fHistDiElectronMass(0),fHistDiMuonMass(0),
100 fHistNeventsPsi2s(0),fHistPsi2sMassVsPt(0),fHistPsi2sMassCoherent(0)
17c65770 101
3d16cd00 102{
103
104 // Constructor
105 if( strstr(name,"ESD") ) fType = 0;
106 if( strstr(name,"AOD") ) fType = 1;
107
108 Init();
109
110 DefineOutput(1, TTree::Class());
111 DefineOutput(2, TTree::Class());
f052ef6f 112 DefineOutput(3, TList::Class());
17c65770 113 DefineOutput(4, TList::Class());
3d16cd00 114
115}//AliAnalysisTaskUpcPsi2s
116
117//_____________________________________________________________________________
118void AliAnalysisTaskUpcPsi2s::Init()
119{
120
121 for(Int_t i=0; i<ntrg; i++) fTrigger[i] = kFALSE;
8b17ae4a 122 for(Int_t i=0; i<4; i++) fTOFphi[i] = -666;
489951cd 123 for(Int_t i=0; i<3; i++){
124 fVtxPos[i] = -666;
125 fVtxErr[i] = -666;
126 fKfVtxPos[i] = -666;
127 }
3d16cd00 128
129}//Init
130
131//_____________________________________________________________________________
132AliAnalysisTaskUpcPsi2s::~AliAnalysisTaskUpcPsi2s()
133{
134 // Destructor
135 if(fJPsiTree){
136 delete fJPsiTree;
137 fJPsiTree = 0x0;
138 }
139 if(fPsi2sTree){
140 delete fPsi2sTree;
141 fPsi2sTree = 0x0;
142 }
f052ef6f 143 if(fListTrig){
144 delete fListTrig;
145 fListTrig = 0x0;
3d16cd00 146 }
46e1d1dc 147 if(fListHist){
148 delete fListHist;
149 fListHist = 0x0;
17c65770 150 }
3d16cd00 151
152}//~AliAnalysisTaskUpcPsi2s
153
154
155//_____________________________________________________________________________
156void AliAnalysisTaskUpcPsi2s::UserCreateOutputObjects()
157{
3d16cd00 158 //input file
159 fDataFilnam = new TObjString();
160 fDataFilnam->SetString("");
161
162 //tracks
163 fJPsiAODTracks = new TClonesArray("AliAODTrack", 1000);
164 fJPsiESDTracks = new TClonesArray("AliESDtrack", 1000);
165 fPsi2sAODTracks = new TClonesArray("AliAODTrack", 1000);
166 fPsi2sESDTracks = new TClonesArray("AliESDtrack", 1000);
167
168 //output tree with JPsi candidate events
169 fJPsiTree = new TTree("fJPsiTree", "fJPsiTree");
170 fJPsiTree ->Branch("fRunNum", &fRunNum, "fRunNum/I");
171 fJPsiTree ->Branch("fPerNum", &fPerNum, "fPerNum/i");
172 fJPsiTree ->Branch("fOrbNum", &fOrbNum, "fOrbNum/i");
173
174 fJPsiTree ->Branch("fBCrossNum", &fBCrossNum, "fBCrossNum/s");
175 fJPsiTree ->Branch("fTrigger", &fTrigger[0], Form("fTrigger[%i]/O", ntrg));
176 fJPsiTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
177 fJPsiTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
178 fJPsiTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
489951cd 179 fJPsiTree ->Branch("fNLooseTracks", &fNLooseTracks, "fNLooseTracks/s");
3d16cd00 180 fJPsiTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
d457ec9e 181
03cff2a0
MB		 182 fJPsiTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
183 fJPsiTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
8b17ae4a 184 fJPsiTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
03cff2a0 185
489951cd 186 fJPsiTree ->Branch("fVtxPos", &fVtxPos[0], "fVtxPos[3]/D");
187 fJPsiTree ->Branch("fVtxErr", &fVtxErr[0], "fVtxErr[3]/D");
d457ec9e
M		 188 fJPsiTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
189 fJPsiTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
190
489951cd 191 fJPsiTree ->Branch("fKfVtxPos", &fKfVtxPos[0], "fKfVtxPos[3]/D");
192
3d16cd00 193 fJPsiTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
194 fJPsiTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
195 fJPsiTree ->Branch("fV0Adecision", &fV0Adecision, "fV0Adecision/I");
196 fJPsiTree ->Branch("fV0Cdecision", &fV0Cdecision, "fV0Cdecision/I");
197 fJPsiTree ->Branch("fDataFilnam", &fDataFilnam);
198 fJPsiTree ->Branch("fRecoPass", &fRecoPass, "fRecoPass/S");
199 fJPsiTree ->Branch("fEvtNum", &fEvtNum, "fEvtNum/L");
200 if( fType == 0 ) {
201 fJPsiTree ->Branch("fJPsiESDTracks", &fJPsiESDTracks);
202 }
203 if( fType == 1 ) {
204 fJPsiTree ->Branch("fJPsiAODTracks", &fJPsiAODTracks);
205 }
206
207 //output tree with Psi2s candidate events
208 fPsi2sTree = new TTree("fPsi2sTree", "fPsi2sTree");
209 fPsi2sTree ->Branch("fRunNum", &fRunNum, "fRunNum/I");
210 fPsi2sTree ->Branch("fPerNum", &fPerNum, "fPerNum/i");
211 fPsi2sTree ->Branch("fOrbNum", &fOrbNum, "fOrbNum/i");
212
213 fPsi2sTree ->Branch("fBCrossNum", &fBCrossNum, "fBCrossNum/s");
214 fPsi2sTree ->Branch("fTrigger", &fTrigger[0], Form("fTrigger[%i]/O", ntrg));
215 fPsi2sTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
216 fPsi2sTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
217 fPsi2sTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
489951cd 218 fPsi2sTree ->Branch("fNLooseTracks", &fNLooseTracks, "fNLooseTracks/s");
3d16cd00 219 fPsi2sTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
d457ec9e 220
03cff2a0
MB		 221 fPsi2sTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
222 fPsi2sTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
8b17ae4a 223 fPsi2sTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
03cff2a0 224
489951cd 225 fPsi2sTree ->Branch("fVtxPos", &fVtxPos[0], "fVtxPos[3]/D");
226 fPsi2sTree ->Branch("fVtxErr", &fVtxErr[0], "fVtxErr[3]/D");
d457ec9e
M		 227 fPsi2sTree ->Branch("fVtxChi2", &fVtxChi2, "fVtxChi2/D");
228 fPsi2sTree ->Branch("fVtxNDF", &fVtxNDF, "fVtxNDF/D");
229
489951cd 230 fPsi2sTree ->Branch("fKfVtxPos", &fKfVtxPos[0], "fKfVtxPos[3]/D");
231
3d16cd00 232 fPsi2sTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
233 fPsi2sTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
234 fPsi2sTree ->Branch("fV0Adecision", &fV0Adecision, "fV0Adecision/I");
235 fPsi2sTree ->Branch("fV0Cdecision", &fV0Cdecision, "fV0Cdecision/I");
236 fPsi2sTree ->Branch("fDataFilnam", &fDataFilnam);
237 fPsi2sTree ->Branch("fRecoPass", &fRecoPass, "fRecoPass/S");
238 fPsi2sTree ->Branch("fEvtNum", &fEvtNum, "fEvtNum/L");
239 if( fType == 0 ) {
240 fPsi2sTree ->Branch("fPsi2sESDTracks", &fPsi2sESDTracks);
241 }
242 if( fType == 1 ) {
243 fPsi2sTree ->Branch("fPsi2sAODTracks", &fPsi2sAODTracks);
244 }
245
f052ef6f 246 fListTrig = new TList();
247 fListTrig ->SetOwner();
248
249 fHistUpcTriggersPerRun = new TH1D("fHistUpcTriggersPerRun", "fHistUpcTriggersPerRun", 3000, 167000.5, 170000.5);
250 fListTrig->Add(fHistUpcTriggersPerRun);
251
252 fHistZedTriggersPerRun = new TH1D("fHistZedTriggersPerRun", "fHistZedTriggersPerRun", 3000, 167000.5, 170000.5);
253 fListTrig->Add(fHistZedTriggersPerRun);
254
255 fHistCvlnTriggersPerRun = new TH1D("fHistCvlnTriggersPerRun", "fHistCvlnTriggersPerRun", 3000, 167000.5, 170000.5);
256 fListTrig->Add(fHistCvlnTriggersPerRun);
257
e63dc644
MB		 258 fHistMBTriggersPerRun = new TH1D("fHistMBTriggersPerRun", "fHistMBTriggersPerRun", 3000, 167000.5, 170000.5);
259 fListTrig->Add(fHistMBTriggersPerRun);
260
261 fHistCentralTriggersPerRun = new TH1D("fHistCentralTriggersPerRun", "fHistCentralTriggersPerRun", 3000, 167000.5, 170000.5);
262 fListTrig->Add(fHistCentralTriggersPerRun);
263
264 fHistSemiCentralTriggersPerRun = new TH1D("fHistSemiCentralTriggersPerRun", "fHistSemiCentralTriggersPerRun", 3000, 167000.5, 170000.5);
265 fListTrig->Add(fHistSemiCentralTriggersPerRun);
266
46e1d1dc 267 fListHist = new TList();
268 fListHist ->SetOwner();
269
17c65770 270 TString CutNameJPsi[12] = {"Analyzed","Triggered","Vertex cut","V0 decision","Two good tracks",
271 "Like sign","Oposite sign","One p_{T}>1", "Both p_{T}>1","PID","Dimuom","Dielectron"};
272 fHistNeventsJPsi = new TH1D("fHistNeventsJPsi","fHistNeventsPsi2s",12,0.5,12.5);
273 for (Int_t i = 0; i<12; i++) fHistNeventsJPsi->GetXaxis()->SetBinLabel(i+1,CutNameJPsi[i].Data());
46e1d1dc 274 fListHist->Add(fHistNeventsJPsi);
17c65770 275
276 fHistTPCsignalJPsi = new TH2D("fHistTPCsignalJPsi","fHistTPCsignalJPsi",240,0,120,240,0,120);
46e1d1dc 277 fListHist->Add(fHistTPCsignalJPsi);
17c65770 278
279 fHistDiLeptonPtJPsi = new TH2D("fHistDiLeptonPtJPsi","fHistDiLeptonPtJPsi",350,0,3.5,350,0,3.5);
46e1d1dc 280 fListHist->Add(fHistDiLeptonPtJPsi);
17c65770 281
282 fHistDiElectronMass = new TH1D("fHistDiElectronMass","Invariant mass of J/#psi candidates",100,2,5);
283 fHistDiElectronMass->GetXaxis()->SetTitle("Invariant mass(e^{+}e^{-}) (GeV/c)");
46e1d1dc 284 fListHist->Add(fHistDiElectronMass);
17c65770 285
286 fHistDiMuonMass = new TH1D("fHistDiMuonMass","Invariant mass of J/#psi candidates",100,2,5);
287 fHistDiMuonMass->GetXaxis()->SetTitle("Invariant mass(#mu^{+}#mu^{-}) (GeV/c)");
46e1d1dc 288 fListHist->Add(fHistDiMuonMass);
17c65770 289
290 TString CutNamePsi2s[13] = {"Analyzed","Triggered","Vertex cut","V0 decision","Four good tracks",
291 "DiLepton - DiPion","Like sign leptons","Like sign pions","Like sign both","Oposite sign","PID","Dimuom","Dielectron"};
292
293 fHistNeventsPsi2s = new TH1D("fHistNeventsPsi2s","fHistNeventsPsi2s",13,0.5,13.5);
294 for (Int_t i = 0; i<13; i++) fHistNeventsPsi2s->GetXaxis()->SetBinLabel(i+1,CutNamePsi2s[i].Data());
46e1d1dc 295 fListHist->Add(fHistNeventsPsi2s);
17c65770 296
297 fHistPsi2sMassVsPt = new TH2D("fHistPsi2sMassVsPt","Mass vs p_{T} of #psi(2s) candidates",100,3,6,50,0,5);
298 fHistPsi2sMassVsPt->GetXaxis()->SetTitle("Invariant mass(l^{+}l^{-}#pi^{+}#pi^{-}) (GeV/c)");
299 fHistPsi2sMassVsPt->GetYaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
46e1d1dc 300 fListHist->Add(fHistPsi2sMassVsPt);
17c65770 301
302 fHistPsi2sMassCoherent = new TH1D("fHistPsi2sMassAllCoherent","Invariant mass of coherent #psi(2s) candidates",100,3,6);
303 fHistPsi2sMassCoherent->GetXaxis()->SetTitle("Invariant mass(l^{+}l^{-}#pi^{+}#pi^{-}) (GeV/c)");
46e1d1dc 304 fListHist->Add(fHistPsi2sMassCoherent);
17c65770 305
3d16cd00 306 PostData(1, fJPsiTree);
307 PostData(2, fPsi2sTree);
f052ef6f 308 PostData(3, fListTrig);
46e1d1dc 309 PostData(4, fListHist);
3d16cd00 310
311}//UserCreateOutputObjects
312
313//_____________________________________________________________________________
314void AliAnalysisTaskUpcPsi2s::UserExec(Option_t *)
315{
316
317 //cout<<"#################### Next event ##################"<<endl;
318
f052ef6f 319 if( fType == 0 ){
320 RunESDtrig();
321 if(fRunHist) RunESDhist();
322 if(fRunTree) RunESDtree();
323 }
324
325 if( fType == 1 ){
326 RunAODtrig();
17c65770 327 if(fRunHist) RunAODhist();
328 if(fRunTree) RunAODtree();
329 }
3d16cd00 330
331}//UserExec
17c65770 332//_____________________________________________________________________________
f052ef6f 333void AliAnalysisTaskUpcPsi2s::RunAODtrig()
334{
335
336 //input event
337 AliAODEvent *aod = (AliAODEvent*) InputEvent();
338 if(!aod) return;
339
340 fRunNum = aod ->GetRunNumber();
341 //Trigger
342 TString trigger = aod->GetFiredTriggerClasses();
343
344 if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
345
e63dc644
MB		 346 if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
347 if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
f052ef6f 348
349 fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
350 if(fL1inputs & (1 << 18)) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
e63dc644
MB		 351
352 //MB, Central and SemiCentral triggers
353 AliCentrality *centrality = aod->GetCentrality();
354 UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
355
03cff2a0
MB		 356 Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
357 //Double_t percentile = centrality->GetCentralityPercentile("V0M");
e63dc644 358
3c79e118 359 if( ((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) ||
360 ((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) ||
361 ((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) ){
362
363 if(percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
364 if(percentile<6 && percentile>0 && (trigger.Contains("CVHN_R2-B"))) fHistCentralTriggersPerRun->Fill(fRunNum);
365 if(percentile<50 && percentile>15 && (trigger.Contains("CVLN_B2-B"))) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
366 if(percentile<50 && percentile>15 && (trigger.Contains("CVLN_R1-B"))) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
367
368 }
369
f052ef6f 370PostData(3, fListTrig);
371
372}
373//_____________________________________________________________________________
17c65770 374void AliAnalysisTaskUpcPsi2s::RunAODhist()
375{
376
377 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
378
379 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
380 Double_t muonMass = partMuon->Mass();
381
382 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
383 Double_t electronMass = partElectron->Mass();
384
385 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
386 Double_t pionMass = partPion->Mass();
387
388 //input event
389 AliAODEvent *aod = (AliAODEvent*) InputEvent();
390 if(!aod) return;
391
392 fHistNeventsJPsi->Fill(1);
393 fHistNeventsPsi2s->Fill(1);
394
395 //Trigger
396 TString trigger = aod->GetFiredTriggerClasses();
397
8b17ae4a 398 if( !trigger.Contains("CCUP") ) return;
17c65770 399
400 fHistNeventsJPsi->Fill(2);
401 fHistNeventsPsi2s->Fill(2);
402
403 //primary vertex
404 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
405 fVtxContrib = fAODVertex->GetNContributors();
406 if(fVtxContrib < 2) return;
407
408 fHistNeventsJPsi->Fill(3);
409 fHistNeventsPsi2s->Fill(3);
410
17c65770 411 //VZERO, ZDC
412 AliAODVZERO *fV0data = aod ->GetVZEROData();
f052ef6f 413 AliAODZDC *fZDCdata = aod->GetZDCData();
17c65770 414
415 fV0Adecision = fV0data->GetV0ADecision();
416 fV0Cdecision = fV0data->GetV0CDecision();
417 if(fV0Adecision != AliAODVZERO::kV0Empty || fV0Cdecision != AliAODVZERO::kV0Empty) return;
418
f052ef6f 419 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
420 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
421
e63dc644 422 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
f052ef6f 423
17c65770 424 fHistNeventsJPsi->Fill(4);
425 fHistNeventsPsi2s->Fill(4);
426
17c65770 427 //Two tracks loop
fc8ccb59 428 Int_t nGoodTracks = 0;
17c65770 429 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
430
431 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
432 Short_t qLepton[4], qPion[4];
433 UInt_t nLepton=0, nPion=0, nHighPt=0;
489951cd 434 Double_t fRecTPCsignal[5];
17c65770 435 Int_t mass[3]={-1,-1,-1};
73c6b208 436
fc8ccb59 437
73c6b208 438 //Four track loop
439 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
440 AliAODTrack *trk = aod->GetTrack(itr);
441 if( !trk ) continue;
442
443 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
444 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
445 if(trk->GetTPCNcls() < 50)continue;
446 if(trk->Chi2perNDF() > 4)continue;
447 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB		 448 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
449 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
450 delete trk_clone;
451
73c6b208 452 if(TMath::Abs(dca[1]) > 2) continue;
453
454 TrackIndex[nGoodTracks] = itr;
455 nGoodTracks++;
456
457 if(nGoodTracks > 4) break;
458 }//Track loop
459
17c65770 460 nLepton=0; nPion=0; nHighPt=0;
461 mass[0]= -1; mass[1]= -1, mass[2]= -1;
462
463 if(nGoodTracks == 4){
464 fHistNeventsPsi2s->Fill(5);
465 for(Int_t i=0; i<4; i++){
466 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
467
468 if(trk->Pt() > 1){
489951cd 469 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
17c65770 470 qLepton[nLepton] = trk->Charge();
489951cd 471 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
17c65770 472 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
473 mass[nLepton] = 0;
474 }
489951cd 475 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
17c65770 476 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
477 mass[nLepton] = 1;
478 }
479 nLepton++;
480 }
481 else{
482 qPion[nPion] = trk->Charge();
483 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
484 nPion++;
485 }
486
487 if(nLepton > 2 || nPion > 2) break;
488 }
489 if((nLepton == 2) && (nPion == 2)){
490 fHistNeventsPsi2s->Fill(6);
491 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
492 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
493 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
494 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
495 fHistNeventsPsi2s->Fill(10);
496 if(mass[0] == mass[1]) {
497 fHistNeventsPsi2s->Fill(11);
498 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
499 vDilepton = vLepton[0]+vLepton[1];
500 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
501 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
502 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
503 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
504 }
505 }
506 }
507 }
508
fc8ccb59
MB		 509 nGoodTracks = 0;
510 //Two track loop
511 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
512 AliAODTrack *trk = aod->GetTrack(itr);
513 if( !trk ) continue;
514
515 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
516 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
517 if(trk->GetTPCNcls() < 70)continue;
518 if(trk->Chi2perNDF() > 4)continue;
519 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
520 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
521 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
522 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
523 delete trk_clone;
524 if(TMath::Abs(dca[1]) > 2) continue;
525 if(TMath::Abs(dca[0]) > 0.2) continue;
526
527 TrackIndex[nGoodTracks] = itr;
528 nGoodTracks++;
529
530 if(nGoodTracks > 2) break;
531 }//Track loop
532
533 nLepton=0; nPion=0; nHighPt=0;
534 mass[0]= -1; mass[1]= -1, mass[2]= -1;
535
536 if(nGoodTracks == 2){
537 fHistNeventsJPsi->Fill(5);
538 for(Int_t i=0; i<2; i++){
539 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
540 if(trk->Pt() > 1) nHighPt++;
489951cd 541 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
fc8ccb59 542 qLepton[nLepton] = trk->Charge();
489951cd 543 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
fc8ccb59
MB		 544 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
545 mass[nLepton] = 0;
546 }
489951cd 547 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
fc8ccb59
MB		 548 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
549 mass[nLepton] = 1;
550 }
551 nLepton++;
552 }
553 if(nLepton == 2){
554 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
555 if(qLepton[0]*qLepton[1] < 0){
556 fHistNeventsJPsi->Fill(7);
557 if(nHighPt > 0){
558 fHistNeventsJPsi->Fill(8);
489951cd 559 fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
fc8ccb59
MB		 560 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
561 if(mass[0] == mass[1] && mass[0] != -1) {
562 fHistNeventsJPsi->Fill(10);
563 vCandidate = vLepton[0]+vLepton[1];
564 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
565 if(mass[0] == 0) {
566 fHistDiMuonMass->Fill(vCandidate.M());
567 fHistNeventsJPsi->Fill(11);
568 }
569 if(mass[0] == 1) {
570 fHistDiElectronMass->Fill(vCandidate.M());
571 fHistNeventsJPsi->Fill(12);
572 }
573 }
574 }
575 }
576 }
577 }
578
579
46e1d1dc 580 PostData(4, fListHist);
17c65770 581
582}
3d16cd00 583
584//_____________________________________________________________________________
17c65770 585void AliAnalysisTaskUpcPsi2s::RunAODtree()
3d16cd00 586{
587 //input event
588 AliAODEvent *aod = (AliAODEvent*) InputEvent();
589 if(!aod) return;
590
591 //input data
592 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
593 fDataFilnam->Clear();
594 fDataFilnam->SetString(filnam);
595 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
596 fRunNum = aod ->GetRunNumber();
597
3d16cd00 598 //Trigger
599 TString trigger = aod->GetFiredTriggerClasses();
600
8b17ae4a
MB		 601 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
602 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
603 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
604
605 Bool_t isTriggered = kFALSE;
606 for(Int_t i=0; i<ntrg; i++) {
607 if( fTrigger[i] ) isTriggered = kTRUE;
608 }
609 if( !isTriggered ) return;
3d16cd00 610
611 //trigger inputs
612 fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
8b17ae4a 613 fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
3d16cd00 614
615 //Event identification
616 fPerNum = aod ->GetPeriodNumber();
617 fOrbNum = aod ->GetOrbitNumber();
618 fBCrossNum = aod ->GetBunchCrossNumber();
619
620 //primary vertex
621 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
622 fVtxContrib = fAODVertex->GetNContributors();
489951cd 623 fVtxPos[0] = fAODVertex->GetX();
624 fVtxPos[1] = fAODVertex->GetY();
625 fVtxPos[2] = fAODVertex->GetZ();
d457ec9e
M		 626 Double_t CovMatx[6];
627 fAODVertex->GetCovarianceMatrix(CovMatx);
489951cd 628 fVtxErr[0] = CovMatx[0];
629 fVtxErr[1] = CovMatx[1];
630 fVtxErr[2] = CovMatx[2];
d457ec9e
M		 631 fVtxChi2 = fAODVertex->GetChi2();
632 fVtxNDF = fAODVertex->GetNDF();
3d16cd00 633
634 //Tracklets
635 fNtracklets = aod->GetTracklets()->GetNumberOfTracklets();
636
637 //VZERO, ZDC
638 AliAODVZERO *fV0data = aod ->GetVZEROData();
639 AliAODZDC *fZDCdata = aod->GetZDCData();
640
641 fV0Adecision = fV0data->GetV0ADecision();
642 fV0Cdecision = fV0data->GetV0CDecision();
643 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
644 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
3d16cd00 645
489951cd 646 fNLooseTracks = 0;
647
648 //Track loop - loose cuts
649 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
650 AliAODTrack *trk = aod->GetTrack(itr);
651 if( !trk ) continue;
652
653 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
654 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
655 if(trk->GetTPCNcls() < 20)continue;
656 fNLooseTracks++;
657 }//Track loop -loose cuts
658
17c65770 659 Int_t nGoodTracks=0;
3d16cd00 660 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
fc8ccb59 661
489951cd 662 //Two track loop
3d16cd00 663 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
664 AliAODTrack *trk = aod->GetTrack(itr);
665 if( !trk ) continue;
3c79e118 666 if(!(trk->TestFilterBit(1<<0)))continue;
3d16cd00 667
489951cd 668 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
669 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
8b17ae4a 670 if(trk->GetTPCNcls() < 70)continue;
3d16cd00 671 if(trk->Chi2perNDF() > 4)continue;
8b17ae4a 672 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
3d16cd00 673 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB		 674 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
675 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
676 delete trk_clone;
3d16cd00 677 if(TMath::Abs(dca[1]) > 2) continue;
8b17ae4a
MB		 678 if(TMath::Abs(dca[0]) > 0.2) continue;
679
3d16cd00 680 TrackIndex[nGoodTracks] = itr;
681 nGoodTracks++;
682
8b17ae4a 683 if(nGoodTracks > 2) break;
3d16cd00 684 }//Track loop
8b17ae4a
MB		 685
686 if(nGoodTracks == 2){
489951cd 687
688 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
689 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
690 Double_t muonMass = partMuon->Mass();
691 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
692 Double_t electronMass = partElectron->Mass();
693 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
694 Double_t pionMass = partPion->Mass();
695
696 Double_t KFcov[21];
697 Double_t KFpar[6];
698 Double_t KFmass = pionMass;
699 Double_t fRecTPCsignal;
700 AliKFParticle *KFpart[2];
701 AliKFVertex *KFvtx = new AliKFVertex();
702 KFvtx->SetField(aod->GetMagneticField());
703
8b17ae4a 704 for(Int_t i=0; i<2; i++){
3d16cd00 705 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
706
707 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB		 708 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
709 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
710 delete trk_clone;
8b17ae4a
MB		 711
712 new((*fJPsiAODTracks)[i]) AliAODTrack(*trk);
713 ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
3d16cd00 714
489951cd 715 trk->GetPosition(KFpar);
716 trk->PxPyPz(KFpar+3);
717 trk->GetCovMatrix(KFcov);
718
719 if(trk->Pt() > 1){
720 fRecTPCsignal = trk->GetTPCsignal();
721 if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
722 if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
723 }
724 else KFmass = pionMass;
725
726 KFpart[i] = new AliKFParticle();
727 KFpart[i]->SetField(aod->GetMagneticField());
728 KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
729 KFvtx->AddDaughter(*KFpart[i]);
730
731
8b17ae4a 732 Double_t pos[3]={0,0,0};
489951cd 733 AliExternalTrackParam *parTrk = new AliExternalTrackParam();
734 parTrk->CopyFromVTrack((AliVTrack*) trk);
735 if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
8b17ae4a
MB		 736 else {
737 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
738 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
489951cd 739 }
740 delete parTrk;
3d16cd00 741 }
489951cd 742 fKfVtxPos[0]= KFvtx->GetX();
743 fKfVtxPos[1]= KFvtx->GetY();
744 fKfVtxPos[2]= KFvtx->GetZ();
745 for(UInt_t i=0; i<2; i++)delete KFpart[i];
746 delete KFvtx;
747
8b17ae4a 748 fJPsiTree ->Fill();
3d16cd00 749 }
8b17ae4a
MB		 750
751 nGoodTracks = 0;
752 //Four track loop
73c6b208 753 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
754 AliAODTrack *trk = aod->GetTrack(itr);
755 if( !trk ) continue;
3c79e118 756 if(!(trk->TestFilterBit(1<<0)))continue;
73c6b208 757
489951cd 758 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
759 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
8b17ae4a 760 if(trk->GetTPCNcls() < 50)continue;
73c6b208 761 if(trk->Chi2perNDF() > 4)continue;
762 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB		 763 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
764 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
765 delete trk_clone;
8b17ae4a 766 if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
73c6b208 767 if(TMath::Abs(dca[1]) > 2) continue;
8b17ae4a 768
73c6b208 769 TrackIndex[nGoodTracks] = itr;
770 nGoodTracks++;
771
8b17ae4a 772 if(nGoodTracks > 4) break;
73c6b208 773 }//Track loop
489951cd 774
8b17ae4a
MB		 775
776 if(nGoodTracks == 4){
489951cd 777
778 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
779 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
780 Double_t muonMass = partMuon->Mass();
781 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
782 Double_t electronMass = partElectron->Mass();
783 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
784 Double_t pionMass = partPion->Mass();
785
786 Double_t KFcov[21];
787 Double_t KFpar[6];
788 Double_t KFmass = pionMass;
789 Double_t fRecTPCsignal;
790 AliKFParticle *KFpart[4];
791 AliKFVertex *KFvtx = new AliKFVertex();
792 KFvtx->SetField(aod->GetMagneticField());
793
8b17ae4a 794 for(Int_t i=0; i<4; i++){
3d16cd00 795 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
796
797 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB		 798 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
799 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
800 delete trk_clone;
3d16cd00 801
8b17ae4a
MB		 802 new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
803 ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
804
489951cd 805 trk->GetPosition(KFpar);
806 trk->PxPyPz(KFpar+3);
807 trk->GetCovMatrix(KFcov);
808
809 if(trk->Pt() > 1){
810 fRecTPCsignal = trk->GetTPCsignal();
811 if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
812 if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
813 }
814 else KFmass = pionMass;
815
816 KFpart[i] = new AliKFParticle();
817 KFpart[i]->SetField(aod->GetMagneticField());
818 KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
819 KFvtx->AddDaughter(*KFpart[i]);
820
8b17ae4a 821 Double_t pos[3]={0,0,0};
489951cd 822 AliExternalTrackParam *parTrk = new AliExternalTrackParam();
823 parTrk->CopyFromVTrack((AliVTrack*) trk);
824 if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
8b17ae4a
MB		 825 else {
826 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
827 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
489951cd 828 }
829 delete parTrk;
3d16cd00 830 }
489951cd 831 fKfVtxPos[0]= KFvtx->GetX();
832 fKfVtxPos[1]= KFvtx->GetY();
833 fKfVtxPos[2]= KFvtx->GetZ();
834 for(UInt_t i=0; i<4; i++)delete KFpart[i];
835 delete KFvtx;
8b17ae4a 836 fPsi2sTree ->Fill();
3d16cd00 837 }
03cff2a0
MB		 838
839 PostData(1, fJPsiTree);
840 PostData(2, fPsi2sTree);
3d16cd00 841
842}//RunAOD
843
844//_____________________________________________________________________________
f052ef6f 845void AliAnalysisTaskUpcPsi2s::RunESDtrig()
846{
847
848 //input event
849 AliESDEvent *esd = (AliESDEvent*) InputEvent();
850 if(!esd) return;
851
852 fRunNum = esd ->GetRunNumber();
853 //Trigger
854 TString trigger = esd->GetFiredTriggerClasses();
855
856 if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
857
e63dc644
MB		 858 if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
859 if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
f052ef6f 860
861 if(esd->GetHeader()->IsTriggerInputFired("1ZED")) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
862
03cff2a0
MB		 863 //MB, Central and SemiCentral triggers
864 AliCentrality *centrality = esd->GetCentrality();
865 UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
866
867 //Double_t percentile = centrality->GetCentralityPercentile("V0M");
868 Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
869
870 if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
871
872 if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
873
874 if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
875
876
f052ef6f 877PostData(3, fListTrig);
878
879}
880//_____________________________________________________________________________
881void AliAnalysisTaskUpcPsi2s::RunESDhist()
882{
883
884 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
885
886 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
887 Double_t muonMass = partMuon->Mass();
888
889 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
890 Double_t electronMass = partElectron->Mass();
891
892 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
893 Double_t pionMass = partPion->Mass();
894
895 //input event
896 AliESDEvent *esd = (AliESDEvent*) InputEvent();
897 if(!esd) return;
898
899 fHistNeventsJPsi->Fill(1);
900 fHistNeventsPsi2s->Fill(1);
901
902 //Trigger
903 TString trigger = esd->GetFiredTriggerClasses();
904
8b17ae4a 905 if( !trigger.Contains("CCUP") ) return;
f052ef6f 906
907 fHistNeventsJPsi->Fill(2);
908 fHistNeventsPsi2s->Fill(2);
909
910 //primary vertex
911 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
912 fVtxContrib = fESDVertex->GetNContributors();
913 if(fVtxContrib < 2) return;
914
915 fHistNeventsJPsi->Fill(3);
916 fHistNeventsPsi2s->Fill(3);
917
918 //VZERO, ZDC
919 AliESDVZERO *fV0data = esd->GetVZEROData();
920 AliESDZDC *fZDCdata = esd->GetESDZDC();
921
922 fV0Adecision = fV0data->GetV0ADecision();
923 fV0Cdecision = fV0data->GetV0CDecision();
924 if(fV0Adecision != AliESDVZERO::kV0Empty || fV0Cdecision != AliESDVZERO::kV0Empty) return;
925
926 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
927 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
03cff2a0 928 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
f052ef6f 929
930 fHistNeventsJPsi->Fill(4);
931 fHistNeventsPsi2s->Fill(4);
932
933 Int_t nGoodTracks=0;
934 //Two tracks loop
935 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
936
937 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
938 Short_t qLepton[4], qPion[4];
939 UInt_t nLepton=0, nPion=0, nHighPt=0;
489951cd 940 Double_t fRecTPCsignal[5];
f052ef6f 941 Int_t mass[3]={-1,-1,-1};
942
03cff2a0 943 //Two Track loop
f052ef6f 944 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
945 AliESDtrack *trk = esd->GetTrack(itr);
946 if( !trk ) continue;
947
948 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
949 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
03cff2a0 950 if(trk->GetTPCNcls() < 70)continue;
f052ef6f 951 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
03cff2a0 952 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
f052ef6f 953 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
954 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
955 trk->GetImpactParameters(dca[0],dca[1]);
956 if(TMath::Abs(dca[1]) > 2) continue;
03cff2a0 957 if(TMath::Abs(dca[1]) > 0.2) continue;
f052ef6f 958
959 TrackIndex[nGoodTracks] = itr;
960 nGoodTracks++;
03cff2a0 961 if(nGoodTracks > 2) break;
f052ef6f 962 }//Track loop
03cff2a0 963
f052ef6f 964
965 if(nGoodTracks == 2){
966 fHistNeventsJPsi->Fill(5);
967 for(Int_t i=0; i<2; i++){
968 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
969 if(trk->Pt() > 1) nHighPt++;
489951cd 970 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
f052ef6f 971 qLepton[nLepton] = trk->Charge();
489951cd 972 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
f052ef6f 973 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
974 mass[nLepton] = 0;
975 }
489951cd 976 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
f052ef6f 977 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
978 mass[nLepton] = 1;
979 }
980 nLepton++;
981 }
982 if(nLepton == 2){
983 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
984 if(qLepton[0]*qLepton[1] < 0){
985 fHistNeventsJPsi->Fill(7);
986 if(nHighPt > 0){
987 fHistNeventsJPsi->Fill(8);
489951cd 988 fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
f052ef6f 989 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
990 if(mass[0] == mass[1] && mass[0] != -1) {
991 fHistNeventsJPsi->Fill(10);
992 vCandidate = vLepton[0]+vLepton[1];
993 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
994 if(mass[0] == 0) {
995 fHistDiMuonMass->Fill(vCandidate.M());
996 fHistNeventsJPsi->Fill(11);
997 }
998 if(mass[0] == 1) {
999 fHistDiElectronMass->Fill(vCandidate.M());
1000 fHistNeventsJPsi->Fill(12);
1001 }
1002 }
1003 }
1004 }
1005 }
1006 }
03cff2a0 1007 nGoodTracks = 0; nLepton=0; nPion=0; nHighPt=0;
f052ef6f 1008 mass[0]= -1; mass[1]= -1, mass[2]= -1;
1009
03cff2a0
MB		 1010 //Four Track loop
1011 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1012 AliESDtrack *trk = esd->GetTrack(itr);
1013 if( !trk ) continue;
1014
1015 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1016 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1017 if(trk->GetTPCNcls() < 50)continue;
1018 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1019 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1020 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1021 trk->GetImpactParameters(dca[0],dca[1]);
1022 if(TMath::Abs(dca[1]) > 2) continue;
1023
1024 TrackIndex[nGoodTracks] = itr;
1025 nGoodTracks++;
1026 if(nGoodTracks > 4) break;
1027 }//Track loop
1028
f052ef6f 1029 if(nGoodTracks == 4){
1030 fHistNeventsPsi2s->Fill(5);
1031 for(Int_t i=0; i<4; i++){
1032 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1033
1034 if(trk->Pt() > 1){
489951cd 1035 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
f052ef6f 1036 qLepton[nLepton] = trk->Charge();
489951cd 1037 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
f052ef6f 1038 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
1039 mass[nLepton] = 0;
1040 }
489951cd 1041 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
f052ef6f 1042 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
1043 mass[nLepton] = 1;
1044 }
1045 nLepton++;
1046 }
1047 else{
1048 qPion[nPion] = trk->Charge();
1049 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
1050 nPion++;
1051 }
1052
1053 if(nLepton > 2 || nPion > 2) break;
1054 }
1055 if((nLepton == 2) && (nPion == 2)){
1056 fHistNeventsPsi2s->Fill(6);
1057 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
1058 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
1059 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
1060 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
1061 fHistNeventsPsi2s->Fill(10);
1062 if(mass[0] == mass[1]) {
1063 fHistNeventsPsi2s->Fill(11);
1064 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
1065 vDilepton = vLepton[0]+vLepton[1];
1066 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
1067 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
1068 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
1069 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
1070 }
1071 }
1072 }
1073 }
1074
1075 PostData(4, fListHist);
1076
1077}
1078
1079//_____________________________________________________________________________
1080void AliAnalysisTaskUpcPsi2s::RunESDtree()
3d16cd00 1081{
1082
1083 //input event
1084 AliESDEvent *esd = (AliESDEvent*) InputEvent();
1085 if(!esd) return;
1086
1087 //input data
1088 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
1089 fDataFilnam->Clear();
1090 fDataFilnam->SetString(filnam);
1091 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
1092 fRunNum = esd->GetRunNumber();
1093
8b17ae4a 1094 //Trigger
3d16cd00 1095 TString trigger = esd->GetFiredTriggerClasses();
1096
8b17ae4a
MB		 1097 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
1098 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
1099 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
1100
1101 Bool_t isTriggered = kFALSE;
1102 for(Int_t i=0; i<ntrg; i++) {
1103 if( fTrigger[i] ) isTriggered = kTRUE;
1104 }
1105 if( !isTriggered ) return;
1106
3d16cd00 1107 //trigger inputs
1108 fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
1109 fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
03cff2a0
MB		 1110
1111 //TOF trigger info (0OMU)
1112 fTOFtrig1 = esd->GetHeader()->IsTriggerInputFired("0OMU");
1113 fTOFtrig2 = esd->GetHeader()->GetActiveTriggerInputs().Contains("0OMU") ? ((esd->GetHeader()) ? esd->GetHeader()->IsTriggerInputFired("0OMU") : kFALSE) : TESTBIT(esd->GetHeader()->GetL0TriggerInputs(), (kFALSE) ? 21 : 9);
3d16cd00 1114
1115 //Event identification
1116 fPerNum = esd->GetPeriodNumber();
1117 fOrbNum = esd->GetOrbitNumber();
1118 fBCrossNum = esd->GetBunchCrossNumber();
1119
1120 //primary vertex
1121 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
1122 fVtxContrib = fESDVertex->GetNContributors();
489951cd 1123 fVtxPos[0] = fESDVertex->GetX();
1124 fVtxPos[1] = fESDVertex->GetY();
1125 fVtxPos[2] = fESDVertex->GetZ();
03cff2a0
MB		 1126 Double_t CovMatx[6];
1127 fESDVertex->GetCovarianceMatrix(CovMatx);
489951cd 1128 fVtxErr[0] = CovMatx[0];
1129 fVtxErr[1] = CovMatx[1];
1130 fVtxErr[2] = CovMatx[2];
03cff2a0
MB		 1131 fVtxChi2 = fESDVertex->GetChi2();
1132 fVtxNDF = fESDVertex->GetNDF();
3d16cd00 1133
1134 //Tracklets
1135 fNtracklets = esd->GetMultiplicity()->GetNumberOfTracklets();
1136
1137 //VZERO, ZDC
1138 AliESDVZERO *fV0data = esd->GetVZEROData();
1139 AliESDZDC *fZDCdata = esd->GetESDZDC();
1140
1141 fV0Adecision = fV0data->GetV0ADecision();
1142 fV0Cdecision = fV0data->GetV0CDecision();
1143 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
1144 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
489951cd 1145
1146 fNLooseTracks = 0;
1147
1148 //Track loop - loose cuts
1149 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1150 AliESDtrack *trk = esd->GetTrack(itr);
1151 if( !trk ) continue;
1152
1153 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1154 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1155 if(trk->GetTPCNcls() < 20)continue;
1156 fNLooseTracks++;
1157 }//Track loop -loose cuts
3d16cd00 1158
17c65770 1159 Int_t nGoodTracks=0;
3d16cd00 1160 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
17c65770 1161
03cff2a0 1162 //Two Track loop
3d16cd00 1163 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1164 AliESDtrack *trk = esd->GetTrack(itr);
1165 if( !trk ) continue;
1166
1167 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1168 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
03cff2a0 1169 if(trk->GetTPCNcls() < 70)continue;
3d16cd00 1170 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
03cff2a0 1171 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
3d16cd00 1172 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1173 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1174 trk->GetImpactParameters(dca[0],dca[1]);
1175 if(TMath::Abs(dca[1]) > 2) continue;
03cff2a0 1176 if(TMath::Abs(dca[1]) > 0.2) continue;
3d16cd00 1177
1178 TrackIndex[nGoodTracks] = itr;
1179 nGoodTracks++;
03cff2a0 1180 if(nGoodTracks > 2) break;
3d16cd00 1181 }//Track loop
1182
1183 if(nGoodTracks == 2){
1184 for(Int_t i=0; i<2; i++){
1185 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1186
1187 AliExternalTrackParam cParam;
1188 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1189
1190 new((*fJPsiESDTracks)[i]) AliESDtrack(*trk);
1191
8b17ae4a
MB		 1192 Double_t pos[3]={0,0,0};
1193 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1194 else {
1195 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1196 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1197 }
3d16cd00 1198 }
1199 fJPsiTree ->Fill();
3d16cd00 1200 }
1201
03cff2a0
MB		 1202 nGoodTracks = 0;
1203 //Four track loop
1204 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1205 AliESDtrack *trk = esd->GetTrack(itr);
1206 if( !trk ) continue;
1207
1208 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1209 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1210 if(trk->GetTPCNcls() < 50)continue;
1211 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1212 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1213 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1214 trk->GetImpactParameters(dca[0],dca[1]);
1215 if(TMath::Abs(dca[1]) > 2) continue;
1216
1217 TrackIndex[nGoodTracks] = itr;
1218 nGoodTracks++;
1219 if(nGoodTracks > 4) break;
1220 }//Track loop
1221
3d16cd00 1222 if(nGoodTracks == 4){
1223 for(Int_t i=0; i<4; i++){
1224 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1225
1226 AliExternalTrackParam cParam;
1227 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1228
8b17ae4a
MB		 1229 new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
1230
1231 Double_t pos[3]={0,0,0};
1232 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1233 else {
1234 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1235 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1236 }
3d16cd00 1237 }
1238 fPsi2sTree ->Fill();
3d16cd00 1239 }
03cff2a0
MB		 1240
1241 PostData(1, fJPsiTree);
1242 PostData(2, fPsi2sTree);
3d16cd00 1243
1244}//RunESD
1245
1246//_____________________________________________________________________________
1247void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
1248{
1249
1250 cout<<"Analysis complete."<<endl;
1251}//Terminate
1252