]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGUD/UPC/AliAnalysisTaskUpcPsi2s.cxx
trigger selection - crosscheck as a V0 daughter
[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
MB
358
359 if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
360
361 if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
362
363 if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
f052ef6f 364
365PostData(3, fListTrig);
366
367}
368//_____________________________________________________________________________
17c65770 369void AliAnalysisTaskUpcPsi2s::RunAODhist()
370{
371
372 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
373
374 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
375 Double_t muonMass = partMuon->Mass();
376
377 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
378 Double_t electronMass = partElectron->Mass();
379
380 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
381 Double_t pionMass = partPion->Mass();
382
383 //input event
384 AliAODEvent *aod = (AliAODEvent*) InputEvent();
385 if(!aod) return;
386
387 fHistNeventsJPsi->Fill(1);
388 fHistNeventsPsi2s->Fill(1);
389
390 //Trigger
391 TString trigger = aod->GetFiredTriggerClasses();
392
8b17ae4a 393 if( !trigger.Contains("CCUP") ) return;
17c65770 394
395 fHistNeventsJPsi->Fill(2);
396 fHistNeventsPsi2s->Fill(2);
397
398 //primary vertex
399 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
400 fVtxContrib = fAODVertex->GetNContributors();
401 if(fVtxContrib < 2) return;
402
403 fHistNeventsJPsi->Fill(3);
404 fHistNeventsPsi2s->Fill(3);
405
17c65770 406 //VZERO, ZDC
407 AliAODVZERO *fV0data = aod ->GetVZEROData();
f052ef6f 408 AliAODZDC *fZDCdata = aod->GetZDCData();
17c65770 409
410 fV0Adecision = fV0data->GetV0ADecision();
411 fV0Cdecision = fV0data->GetV0CDecision();
412 if(fV0Adecision != AliAODVZERO::kV0Empty || fV0Cdecision != AliAODVZERO::kV0Empty) return;
413
f052ef6f 414 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
415 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
416
e63dc644 417 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
f052ef6f 418
17c65770 419 fHistNeventsJPsi->Fill(4);
420 fHistNeventsPsi2s->Fill(4);
421
17c65770 422 //Two tracks loop
fc8ccb59 423 Int_t nGoodTracks = 0;
17c65770 424 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
425
426 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
427 Short_t qLepton[4], qPion[4];
428 UInt_t nLepton=0, nPion=0, nHighPt=0;
489951cd 429 Double_t fRecTPCsignal[5];
17c65770 430 Int_t mass[3]={-1,-1,-1};
73c6b208 431
fc8ccb59 432
73c6b208 433 //Four track loop
434 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
435 AliAODTrack *trk = aod->GetTrack(itr);
436 if( !trk ) continue;
437
438 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
439 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
440 if(trk->GetTPCNcls() < 50)continue;
441 if(trk->Chi2perNDF() > 4)continue;
442 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB
443 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
444 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
445 delete trk_clone;
446
73c6b208 447 if(TMath::Abs(dca[1]) > 2) continue;
448
449 TrackIndex[nGoodTracks] = itr;
450 nGoodTracks++;
451
452 if(nGoodTracks > 4) break;
453 }//Track loop
454
17c65770 455 nLepton=0; nPion=0; nHighPt=0;
456 mass[0]= -1; mass[1]= -1, mass[2]= -1;
457
458 if(nGoodTracks == 4){
459 fHistNeventsPsi2s->Fill(5);
460 for(Int_t i=0; i<4; i++){
461 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
462
463 if(trk->Pt() > 1){
489951cd 464 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
17c65770 465 qLepton[nLepton] = trk->Charge();
489951cd 466 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
17c65770 467 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
468 mass[nLepton] = 0;
469 }
489951cd 470 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
17c65770 471 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
472 mass[nLepton] = 1;
473 }
474 nLepton++;
475 }
476 else{
477 qPion[nPion] = trk->Charge();
478 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
479 nPion++;
480 }
481
482 if(nLepton > 2 || nPion > 2) break;
483 }
484 if((nLepton == 2) && (nPion == 2)){
485 fHistNeventsPsi2s->Fill(6);
486 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
487 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
488 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
489 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
490 fHistNeventsPsi2s->Fill(10);
491 if(mass[0] == mass[1]) {
492 fHistNeventsPsi2s->Fill(11);
493 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
494 vDilepton = vLepton[0]+vLepton[1];
495 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
496 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
497 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
498 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
499 }
500 }
501 }
502 }
503
fc8ccb59
MB
504 nGoodTracks = 0;
505 //Two track loop
506 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
507 AliAODTrack *trk = aod->GetTrack(itr);
508 if( !trk ) continue;
509
510 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
511 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
512 if(trk->GetTPCNcls() < 70)continue;
513 if(trk->Chi2perNDF() > 4)continue;
514 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
515 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
516 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
517 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
518 delete trk_clone;
519 if(TMath::Abs(dca[1]) > 2) continue;
520 if(TMath::Abs(dca[0]) > 0.2) continue;
521
522 TrackIndex[nGoodTracks] = itr;
523 nGoodTracks++;
524
525 if(nGoodTracks > 2) break;
526 }//Track loop
527
528 nLepton=0; nPion=0; nHighPt=0;
529 mass[0]= -1; mass[1]= -1, mass[2]= -1;
530
531 if(nGoodTracks == 2){
532 fHistNeventsJPsi->Fill(5);
533 for(Int_t i=0; i<2; i++){
534 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
535 if(trk->Pt() > 1) nHighPt++;
489951cd 536 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
fc8ccb59 537 qLepton[nLepton] = trk->Charge();
489951cd 538 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
fc8ccb59
MB
539 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
540 mass[nLepton] = 0;
541 }
489951cd 542 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
fc8ccb59
MB
543 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
544 mass[nLepton] = 1;
545 }
546 nLepton++;
547 }
548 if(nLepton == 2){
549 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
550 if(qLepton[0]*qLepton[1] < 0){
551 fHistNeventsJPsi->Fill(7);
552 if(nHighPt > 0){
553 fHistNeventsJPsi->Fill(8);
489951cd 554 fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
fc8ccb59
MB
555 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
556 if(mass[0] == mass[1] && mass[0] != -1) {
557 fHistNeventsJPsi->Fill(10);
558 vCandidate = vLepton[0]+vLepton[1];
559 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
560 if(mass[0] == 0) {
561 fHistDiMuonMass->Fill(vCandidate.M());
562 fHistNeventsJPsi->Fill(11);
563 }
564 if(mass[0] == 1) {
565 fHistDiElectronMass->Fill(vCandidate.M());
566 fHistNeventsJPsi->Fill(12);
567 }
568 }
569 }
570 }
571 }
572 }
573
574
46e1d1dc 575 PostData(4, fListHist);
17c65770 576
577}
3d16cd00 578
579//_____________________________________________________________________________
17c65770 580void AliAnalysisTaskUpcPsi2s::RunAODtree()
3d16cd00 581{
582 //input event
583 AliAODEvent *aod = (AliAODEvent*) InputEvent();
584 if(!aod) return;
585
586 //input data
587 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
588 fDataFilnam->Clear();
589 fDataFilnam->SetString(filnam);
590 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
591 fRunNum = aod ->GetRunNumber();
592
3d16cd00 593 //Trigger
594 TString trigger = aod->GetFiredTriggerClasses();
595
8b17ae4a
MB
596 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
597 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
598 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
599
600 Bool_t isTriggered = kFALSE;
601 for(Int_t i=0; i<ntrg; i++) {
602 if( fTrigger[i] ) isTriggered = kTRUE;
603 }
604 if( !isTriggered ) return;
3d16cd00 605
606 //trigger inputs
607 fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
8b17ae4a 608 fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
3d16cd00 609
610 //Event identification
611 fPerNum = aod ->GetPeriodNumber();
612 fOrbNum = aod ->GetOrbitNumber();
613 fBCrossNum = aod ->GetBunchCrossNumber();
614
615 //primary vertex
616 AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
617 fVtxContrib = fAODVertex->GetNContributors();
489951cd 618 fVtxPos[0] = fAODVertex->GetX();
619 fVtxPos[1] = fAODVertex->GetY();
620 fVtxPos[2] = fAODVertex->GetZ();
d457ec9e
M
621 Double_t CovMatx[6];
622 fAODVertex->GetCovarianceMatrix(CovMatx);
489951cd 623 fVtxErr[0] = CovMatx[0];
624 fVtxErr[1] = CovMatx[1];
625 fVtxErr[2] = CovMatx[2];
d457ec9e
M
626 fVtxChi2 = fAODVertex->GetChi2();
627 fVtxNDF = fAODVertex->GetNDF();
3d16cd00 628
629 //Tracklets
630 fNtracklets = aod->GetTracklets()->GetNumberOfTracklets();
631
632 //VZERO, ZDC
633 AliAODVZERO *fV0data = aod ->GetVZEROData();
634 AliAODZDC *fZDCdata = aod->GetZDCData();
635
636 fV0Adecision = fV0data->GetV0ADecision();
637 fV0Cdecision = fV0data->GetV0CDecision();
638 fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
639 fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
3d16cd00 640
489951cd 641 fNLooseTracks = 0;
642
643 //Track loop - loose cuts
644 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
645 AliAODTrack *trk = aod->GetTrack(itr);
646 if( !trk ) continue;
647
648 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
649 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
650 if(trk->GetTPCNcls() < 20)continue;
651 fNLooseTracks++;
652 }//Track loop -loose cuts
653
17c65770 654 Int_t nGoodTracks=0;
3d16cd00 655 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
fc8ccb59 656
489951cd 657 //Two track loop
3d16cd00 658 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
659 AliAODTrack *trk = aod->GetTrack(itr);
660 if( !trk ) continue;
661
489951cd 662 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
663 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
8b17ae4a 664 if(trk->GetTPCNcls() < 70)continue;
3d16cd00 665 if(trk->Chi2perNDF() > 4)continue;
8b17ae4a 666 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
3d16cd00 667 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB
668 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
669 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
670 delete trk_clone;
3d16cd00 671 if(TMath::Abs(dca[1]) > 2) continue;
8b17ae4a
MB
672 if(TMath::Abs(dca[0]) > 0.2) continue;
673
3d16cd00 674 TrackIndex[nGoodTracks] = itr;
675 nGoodTracks++;
676
8b17ae4a 677 if(nGoodTracks > 2) break;
3d16cd00 678 }//Track loop
8b17ae4a
MB
679
680 if(nGoodTracks == 2){
489951cd 681
682 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
683 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
684 Double_t muonMass = partMuon->Mass();
685 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
686 Double_t electronMass = partElectron->Mass();
687 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
688 Double_t pionMass = partPion->Mass();
689
690 Double_t KFcov[21];
691 Double_t KFpar[6];
692 Double_t KFmass = pionMass;
693 Double_t fRecTPCsignal;
694 AliKFParticle *KFpart[2];
695 AliKFVertex *KFvtx = new AliKFVertex();
696 KFvtx->SetField(aod->GetMagneticField());
697
8b17ae4a 698 for(Int_t i=0; i<2; i++){
3d16cd00 699 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
700
701 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB
702 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
703 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
704 delete trk_clone;
8b17ae4a
MB
705
706 new((*fJPsiAODTracks)[i]) AliAODTrack(*trk);
707 ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
3d16cd00 708
489951cd 709 trk->GetPosition(KFpar);
710 trk->PxPyPz(KFpar+3);
711 trk->GetCovMatrix(KFcov);
712
713 if(trk->Pt() > 1){
714 fRecTPCsignal = trk->GetTPCsignal();
715 if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
716 if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
717 }
718 else KFmass = pionMass;
719
720 KFpart[i] = new AliKFParticle();
721 KFpart[i]->SetField(aod->GetMagneticField());
722 KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
723 KFvtx->AddDaughter(*KFpart[i]);
724
725
8b17ae4a 726 Double_t pos[3]={0,0,0};
489951cd 727 AliExternalTrackParam *parTrk = new AliExternalTrackParam();
728 parTrk->CopyFromVTrack((AliVTrack*) trk);
729 if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
8b17ae4a
MB
730 else {
731 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
732 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
489951cd 733 }
734 delete parTrk;
3d16cd00 735 }
489951cd 736 fKfVtxPos[0]= KFvtx->GetX();
737 fKfVtxPos[1]= KFvtx->GetY();
738 fKfVtxPos[2]= KFvtx->GetZ();
739 for(UInt_t i=0; i<2; i++)delete KFpart[i];
740 delete KFvtx;
741
8b17ae4a 742 fJPsiTree ->Fill();
3d16cd00 743 }
8b17ae4a
MB
744
745 nGoodTracks = 0;
746 //Four track loop
73c6b208 747 for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
748 AliAODTrack *trk = aod->GetTrack(itr);
749 if( !trk ) continue;
750
489951cd 751 if(!(trk->GetStatus() & AliAODTrack::kTPCrefit) ) continue;
752 if(!(trk->GetStatus() & AliAODTrack::kITSrefit) ) continue;
8b17ae4a 753 if(trk->GetTPCNcls() < 50)continue;
73c6b208 754 if(trk->Chi2perNDF() > 4)continue;
755 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB
756 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
757 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
758 delete trk_clone;
8b17ae4a 759 if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
73c6b208 760 if(TMath::Abs(dca[1]) > 2) continue;
8b17ae4a 761
73c6b208 762 TrackIndex[nGoodTracks] = itr;
763 nGoodTracks++;
764
8b17ae4a 765 if(nGoodTracks > 4) break;
73c6b208 766 }//Track loop
489951cd 767
8b17ae4a
MB
768
769 if(nGoodTracks == 4){
489951cd 770
771 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
772 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
773 Double_t muonMass = partMuon->Mass();
774 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
775 Double_t electronMass = partElectron->Mass();
776 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
777 Double_t pionMass = partPion->Mass();
778
779 Double_t KFcov[21];
780 Double_t KFpar[6];
781 Double_t KFmass = pionMass;
782 Double_t fRecTPCsignal;
783 AliKFParticle *KFpart[4];
784 AliKFVertex *KFvtx = new AliKFVertex();
785 KFvtx->SetField(aod->GetMagneticField());
786
8b17ae4a 787 for(Int_t i=0; i<4; i++){
3d16cd00 788 AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
789
790 Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
fc8ccb59
MB
791 AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
792 if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
793 delete trk_clone;
3d16cd00 794
8b17ae4a
MB
795 new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
796 ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
797
489951cd 798 trk->GetPosition(KFpar);
799 trk->PxPyPz(KFpar+3);
800 trk->GetCovMatrix(KFcov);
801
802 if(trk->Pt() > 1){
803 fRecTPCsignal = trk->GetTPCsignal();
804 if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
805 if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
806 }
807 else KFmass = pionMass;
808
809 KFpart[i] = new AliKFParticle();
810 KFpart[i]->SetField(aod->GetMagneticField());
811 KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
812 KFvtx->AddDaughter(*KFpart[i]);
813
8b17ae4a 814 Double_t pos[3]={0,0,0};
489951cd 815 AliExternalTrackParam *parTrk = new AliExternalTrackParam();
816 parTrk->CopyFromVTrack((AliVTrack*) trk);
817 if(!parTrk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
8b17ae4a
MB
818 else {
819 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
820 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
489951cd 821 }
822 delete parTrk;
3d16cd00 823 }
489951cd 824 fKfVtxPos[0]= KFvtx->GetX();
825 fKfVtxPos[1]= KFvtx->GetY();
826 fKfVtxPos[2]= KFvtx->GetZ();
827 for(UInt_t i=0; i<4; i++)delete KFpart[i];
828 delete KFvtx;
8b17ae4a 829 fPsi2sTree ->Fill();
3d16cd00 830 }
03cff2a0
MB
831
832 PostData(1, fJPsiTree);
833 PostData(2, fPsi2sTree);
3d16cd00 834
835}//RunAOD
836
837//_____________________________________________________________________________
f052ef6f 838void AliAnalysisTaskUpcPsi2s::RunESDtrig()
839{
840
841 //input event
842 AliESDEvent *esd = (AliESDEvent*) InputEvent();
843 if(!esd) return;
844
845 fRunNum = esd ->GetRunNumber();
846 //Trigger
847 TString trigger = esd->GetFiredTriggerClasses();
848
849 if(trigger.Contains("CCUP4-B")) fHistUpcTriggersPerRun->Fill(fRunNum); //Upc triggers
850
e63dc644
MB
851 if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
852 if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
f052ef6f 853
854 if(esd->GetHeader()->IsTriggerInputFired("1ZED")) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
855
03cff2a0
MB
856 //MB, Central and SemiCentral triggers
857 AliCentrality *centrality = esd->GetCentrality();
858 UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
859
860 //Double_t percentile = centrality->GetCentralityPercentile("V0M");
861 Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
862
863 if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<80 && percentile>0) fHistMBTriggersPerRun->Fill(fRunNum);
864
865 if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<6 && percentile>0) fHistCentralTriggersPerRun->Fill(fRunNum);
866
867 if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<50 && percentile>15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
868
869
f052ef6f 870PostData(3, fListTrig);
871
872}
873//_____________________________________________________________________________
874void AliAnalysisTaskUpcPsi2s::RunESDhist()
875{
876
877 TDatabasePDG *pdgdat = TDatabasePDG::Instance();
878
879 TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
880 Double_t muonMass = partMuon->Mass();
881
882 TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
883 Double_t electronMass = partElectron->Mass();
884
885 TParticlePDG *partPion = pdgdat->GetParticle( 211 );
886 Double_t pionMass = partPion->Mass();
887
888 //input event
889 AliESDEvent *esd = (AliESDEvent*) InputEvent();
890 if(!esd) return;
891
892 fHistNeventsJPsi->Fill(1);
893 fHistNeventsPsi2s->Fill(1);
894
895 //Trigger
896 TString trigger = esd->GetFiredTriggerClasses();
897
8b17ae4a 898 if( !trigger.Contains("CCUP") ) return;
f052ef6f 899
900 fHistNeventsJPsi->Fill(2);
901 fHistNeventsPsi2s->Fill(2);
902
903 //primary vertex
904 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
905 fVtxContrib = fESDVertex->GetNContributors();
906 if(fVtxContrib < 2) return;
907
908 fHistNeventsJPsi->Fill(3);
909 fHistNeventsPsi2s->Fill(3);
910
911 //VZERO, ZDC
912 AliESDVZERO *fV0data = esd->GetVZEROData();
913 AliESDZDC *fZDCdata = esd->GetESDZDC();
914
915 fV0Adecision = fV0data->GetV0ADecision();
916 fV0Cdecision = fV0data->GetV0CDecision();
917 if(fV0Adecision != AliESDVZERO::kV0Empty || fV0Cdecision != AliESDVZERO::kV0Empty) return;
918
919 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
920 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
03cff2a0 921 if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
f052ef6f 922
923 fHistNeventsJPsi->Fill(4);
924 fHistNeventsPsi2s->Fill(4);
925
926 Int_t nGoodTracks=0;
927 //Two tracks loop
928 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
929
930 TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
931 Short_t qLepton[4], qPion[4];
932 UInt_t nLepton=0, nPion=0, nHighPt=0;
489951cd 933 Double_t fRecTPCsignal[5];
f052ef6f 934 Int_t mass[3]={-1,-1,-1};
935
03cff2a0 936 //Two Track loop
f052ef6f 937 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
938 AliESDtrack *trk = esd->GetTrack(itr);
939 if( !trk ) continue;
940
941 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
942 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
03cff2a0 943 if(trk->GetTPCNcls() < 70)continue;
f052ef6f 944 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
03cff2a0 945 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
f052ef6f 946 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
947 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
948 trk->GetImpactParameters(dca[0],dca[1]);
949 if(TMath::Abs(dca[1]) > 2) continue;
03cff2a0 950 if(TMath::Abs(dca[1]) > 0.2) continue;
f052ef6f 951
952 TrackIndex[nGoodTracks] = itr;
953 nGoodTracks++;
03cff2a0 954 if(nGoodTracks > 2) break;
f052ef6f 955 }//Track loop
03cff2a0 956
f052ef6f 957
958 if(nGoodTracks == 2){
959 fHistNeventsJPsi->Fill(5);
960 for(Int_t i=0; i<2; i++){
961 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
962 if(trk->Pt() > 1) nHighPt++;
489951cd 963 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
f052ef6f 964 qLepton[nLepton] = trk->Charge();
489951cd 965 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
f052ef6f 966 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
967 mass[nLepton] = 0;
968 }
489951cd 969 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
f052ef6f 970 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
971 mass[nLepton] = 1;
972 }
973 nLepton++;
974 }
975 if(nLepton == 2){
976 if(qLepton[0]*qLepton[1] > 0) fHistNeventsJPsi->Fill(6);
977 if(qLepton[0]*qLepton[1] < 0){
978 fHistNeventsJPsi->Fill(7);
979 if(nHighPt > 0){
980 fHistNeventsJPsi->Fill(8);
489951cd 981 fHistTPCsignalJPsi->Fill(fRecTPCsignal[0],fRecTPCsignal[1]);
f052ef6f 982 if(nHighPt == 2) fHistNeventsJPsi->Fill(9);
983 if(mass[0] == mass[1] && mass[0] != -1) {
984 fHistNeventsJPsi->Fill(10);
985 vCandidate = vLepton[0]+vLepton[1];
986 if( vCandidate.M() > 2.8 && vCandidate.M() < 3.2) fHistDiLeptonPtJPsi->Fill(vLepton[0].Pt(),vLepton[1].Pt());
987 if(mass[0] == 0) {
988 fHistDiMuonMass->Fill(vCandidate.M());
989 fHistNeventsJPsi->Fill(11);
990 }
991 if(mass[0] == 1) {
992 fHistDiElectronMass->Fill(vCandidate.M());
993 fHistNeventsJPsi->Fill(12);
994 }
995 }
996 }
997 }
998 }
999 }
03cff2a0 1000 nGoodTracks = 0; nLepton=0; nPion=0; nHighPt=0;
f052ef6f 1001 mass[0]= -1; mass[1]= -1, mass[2]= -1;
1002
03cff2a0
MB
1003 //Four Track loop
1004 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1005 AliESDtrack *trk = esd->GetTrack(itr);
1006 if( !trk ) continue;
1007
1008 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1009 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1010 if(trk->GetTPCNcls() < 50)continue;
1011 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1012 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1013 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1014 trk->GetImpactParameters(dca[0],dca[1]);
1015 if(TMath::Abs(dca[1]) > 2) continue;
1016
1017 TrackIndex[nGoodTracks] = itr;
1018 nGoodTracks++;
1019 if(nGoodTracks > 4) break;
1020 }//Track loop
1021
f052ef6f 1022 if(nGoodTracks == 4){
1023 fHistNeventsPsi2s->Fill(5);
1024 for(Int_t i=0; i<4; i++){
1025 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1026
1027 if(trk->Pt() > 1){
489951cd 1028 fRecTPCsignal[nLepton] = trk->GetTPCsignal();
f052ef6f 1029 qLepton[nLepton] = trk->Charge();
489951cd 1030 if(fRecTPCsignal[nLepton] > 40 && fRecTPCsignal[nLepton] < 70){
f052ef6f 1031 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), muonMass);
1032 mass[nLepton] = 0;
1033 }
489951cd 1034 if(fRecTPCsignal[nLepton] > 70 && fRecTPCsignal[nLepton] < 100){
f052ef6f 1035 vLepton[nLepton].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), electronMass);
1036 mass[nLepton] = 1;
1037 }
1038 nLepton++;
1039 }
1040 else{
1041 qPion[nPion] = trk->Charge();
1042 vPion[nPion].SetPtEtaPhiM(trk->Pt(), trk->Eta(), trk->Phi(), pionMass);
1043 nPion++;
1044 }
1045
1046 if(nLepton > 2 || nPion > 2) break;
1047 }
1048 if((nLepton == 2) && (nPion == 2)){
1049 fHistNeventsPsi2s->Fill(6);
1050 if(qLepton[0]*qLepton[1] > 0) fHistNeventsPsi2s->Fill(7);
1051 if(qPion[0]*qPion[1] > 0) fHistNeventsPsi2s->Fill(8);
1052 if((qLepton[0]*qLepton[1] > 0) && (qPion[0]*qPion[1] > 0)) fHistNeventsPsi2s->Fill(9);
1053 if((qLepton[0]*qLepton[1] < 0) && (qPion[0]*qPion[1] < 0)){
1054 fHistNeventsPsi2s->Fill(10);
1055 if(mass[0] == mass[1]) {
1056 fHistNeventsPsi2s->Fill(11);
1057 vCandidate = vLepton[0]+vLepton[1]+vPion[0]+vPion[1];
1058 vDilepton = vLepton[0]+vLepton[1];
1059 fHistPsi2sMassVsPt->Fill(vCandidate.M(),vCandidate.Pt());
1060 if(vCandidate.Pt() < 0.15) fHistPsi2sMassCoherent->Fill(vCandidate.M());
1061 if(mass[0] == 0) fHistNeventsPsi2s->Fill(12);
1062 if(mass[0] == 1) fHistNeventsPsi2s->Fill(13);
1063 }
1064 }
1065 }
1066 }
1067
1068 PostData(4, fListHist);
1069
1070}
1071
1072//_____________________________________________________________________________
1073void AliAnalysisTaskUpcPsi2s::RunESDtree()
3d16cd00 1074{
1075
1076 //input event
1077 AliESDEvent *esd = (AliESDEvent*) InputEvent();
1078 if(!esd) return;
1079
1080 //input data
1081 const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
1082 fDataFilnam->Clear();
1083 fDataFilnam->SetString(filnam);
1084 fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
1085 fRunNum = esd->GetRunNumber();
1086
8b17ae4a 1087 //Trigger
3d16cd00 1088 TString trigger = esd->GetFiredTriggerClasses();
1089
8b17ae4a
MB
1090 fTrigger[0] = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
1091 fTrigger[1] = trigger.Contains("CCUP2-B"); // Double gap
1092 fTrigger[2] = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
1093
1094 Bool_t isTriggered = kFALSE;
1095 for(Int_t i=0; i<ntrg; i++) {
1096 if( fTrigger[i] ) isTriggered = kTRUE;
1097 }
1098 if( !isTriggered ) return;
1099
3d16cd00 1100 //trigger inputs
1101 fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
1102 fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
03cff2a0
MB
1103
1104 //TOF trigger info (0OMU)
1105 fTOFtrig1 = esd->GetHeader()->IsTriggerInputFired("0OMU");
1106 fTOFtrig2 = esd->GetHeader()->GetActiveTriggerInputs().Contains("0OMU") ? ((esd->GetHeader()) ? esd->GetHeader()->IsTriggerInputFired("0OMU") : kFALSE) : TESTBIT(esd->GetHeader()->GetL0TriggerInputs(), (kFALSE) ? 21 : 9);
3d16cd00 1107
1108 //Event identification
1109 fPerNum = esd->GetPeriodNumber();
1110 fOrbNum = esd->GetOrbitNumber();
1111 fBCrossNum = esd->GetBunchCrossNumber();
1112
1113 //primary vertex
1114 AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
1115 fVtxContrib = fESDVertex->GetNContributors();
489951cd 1116 fVtxPos[0] = fESDVertex->GetX();
1117 fVtxPos[1] = fESDVertex->GetY();
1118 fVtxPos[2] = fESDVertex->GetZ();
03cff2a0
MB
1119 Double_t CovMatx[6];
1120 fESDVertex->GetCovarianceMatrix(CovMatx);
489951cd 1121 fVtxErr[0] = CovMatx[0];
1122 fVtxErr[1] = CovMatx[1];
1123 fVtxErr[2] = CovMatx[2];
03cff2a0
MB
1124 fVtxChi2 = fESDVertex->GetChi2();
1125 fVtxNDF = fESDVertex->GetNDF();
3d16cd00 1126
1127 //Tracklets
1128 fNtracklets = esd->GetMultiplicity()->GetNumberOfTracklets();
1129
1130 //VZERO, ZDC
1131 AliESDVZERO *fV0data = esd->GetVZEROData();
1132 AliESDZDC *fZDCdata = esd->GetESDZDC();
1133
1134 fV0Adecision = fV0data->GetV0ADecision();
1135 fV0Cdecision = fV0data->GetV0CDecision();
1136 fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
1137 fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
489951cd 1138
1139 fNLooseTracks = 0;
1140
1141 //Track loop - loose cuts
1142 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1143 AliESDtrack *trk = esd->GetTrack(itr);
1144 if( !trk ) continue;
1145
1146 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1147 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1148 if(trk->GetTPCNcls() < 20)continue;
1149 fNLooseTracks++;
1150 }//Track loop -loose cuts
3d16cd00 1151
17c65770 1152 Int_t nGoodTracks=0;
3d16cd00 1153 Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
17c65770 1154
03cff2a0 1155 //Two Track loop
3d16cd00 1156 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1157 AliESDtrack *trk = esd->GetTrack(itr);
1158 if( !trk ) continue;
1159
1160 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1161 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
03cff2a0 1162 if(trk->GetTPCNcls() < 70)continue;
3d16cd00 1163 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
03cff2a0 1164 if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
3d16cd00 1165 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1166 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1167 trk->GetImpactParameters(dca[0],dca[1]);
1168 if(TMath::Abs(dca[1]) > 2) continue;
03cff2a0 1169 if(TMath::Abs(dca[1]) > 0.2) continue;
3d16cd00 1170
1171 TrackIndex[nGoodTracks] = itr;
1172 nGoodTracks++;
03cff2a0 1173 if(nGoodTracks > 2) break;
3d16cd00 1174 }//Track loop
1175
1176 if(nGoodTracks == 2){
1177 for(Int_t i=0; i<2; i++){
1178 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1179
1180 AliExternalTrackParam cParam;
1181 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1182
1183 new((*fJPsiESDTracks)[i]) AliESDtrack(*trk);
1184
8b17ae4a
MB
1185 Double_t pos[3]={0,0,0};
1186 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1187 else {
1188 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1189 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1190 }
3d16cd00 1191 }
1192 fJPsiTree ->Fill();
3d16cd00 1193 }
1194
03cff2a0
MB
1195 nGoodTracks = 0;
1196 //Four track loop
1197 for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
1198 AliESDtrack *trk = esd->GetTrack(itr);
1199 if( !trk ) continue;
1200
1201 if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
1202 if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
1203 if(trk->GetTPCNcls() < 50)continue;
1204 if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
1205 Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
1206 if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
1207 trk->GetImpactParameters(dca[0],dca[1]);
1208 if(TMath::Abs(dca[1]) > 2) continue;
1209
1210 TrackIndex[nGoodTracks] = itr;
1211 nGoodTracks++;
1212 if(nGoodTracks > 4) break;
1213 }//Track loop
1214
3d16cd00 1215 if(nGoodTracks == 4){
1216 for(Int_t i=0; i<4; i++){
1217 AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
1218
1219 AliExternalTrackParam cParam;
1220 trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
1221
8b17ae4a
MB
1222 new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
1223
1224 Double_t pos[3]={0,0,0};
1225 if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
1226 else {
1227 fTOFphi[i] = TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
1228 if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
1229 }
3d16cd00 1230 }
1231 fPsi2sTree ->Fill();
3d16cd00 1232 }
03cff2a0
MB
1233
1234 PostData(1, fJPsiTree);
1235 PostData(2, fPsi2sTree);
3d16cd00 1236
1237}//RunESD
1238
1239//_____________________________________________________________________________
1240void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
1241{
1242
1243 cout<<"Analysis complete."<<endl;
1244}//Terminate
1245