Fixed bug in the computation of dip angle value
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnAnalysisPhi7TeVNoPID.cxx
CommitLineData
2dab9030 1//
2// Implementation file for implementation of data analysis aft 900 GeV
3//
4// Author: A. Pulvirenti
5//
6
7#include "Riostream.h"
8#include <iomanip>
9
10#include "TH1.h"
11#include "TTree.h"
12#include "TParticle.h"
13#include "TRandom.h"
14#include "TLorentzVector.h"
15
16#include "AliLog.h"
17#include "AliESDpid.h"
18#include "AliESDEvent.h"
19#include "AliESDVertex.h"
20#include "AliESDtrack.h"
21#include "AliStack.h"
22#include "AliMCEvent.h"
23#include "AliTOFT0maker.h"
24#include "AliTOFcalib.h"
25#include "AliCDBManager.h"
26#include "AliITSPIDResponse.h"
27
28#include "AliRsnAnalysisPhi7TeVNoPID.h"
29
30//__________________________________________________________________________________________________
31AliRsnAnalysisPhi7TeVNoPID::AliRsnAnalysisPhi7TeVNoPID(const char *name) :
32 AliAnalysisTaskSE(name),
33 fUseMC(kFALSE),
34 fPDG(0),
35 fCh(0),
36 fIM(0.0),
37 fPt(0.0),
38 fY(0.0),
39 fEta(0.0),
40 fMaxVz(1E6),
41 fMaxITSband(1E6),
42 fTPCpLimit(0.35),
43 fMinTPCband(-1E6),
44 fMaxTPCband( 1E6),
45 fRsnTreeComp(0x0),
46 fRsnTreeTrue(0x0),
47 fOutList(0x0),
48 fHEvents(0x0),
49 fESDtrackCutsTPC(),
50 fESDtrackCutsITS(),
51 fESDpid(0x0),
52 fTOFmaker(0x0),
53 fTOFcalib(0x0),
54 fTOFcalibrateESD(kFALSE),
55 fTOFcorrectTExp(kFALSE),
56 fTOFuseT0(kFALSE),
57 fTOFtuneMC(kFALSE),
58 fTOFresolution(0.0)
59
60{
61//
62// Constructor
63//
64
65 DefineOutput(1, TTree::Class());
66 DefineOutput(2, TTree::Class());
67 DefineOutput(3, TList::Class());
68}
69
70//__________________________________________________________________________________________________
71AliRsnAnalysisPhi7TeVNoPID::AliRsnAnalysisPhi7TeVNoPID(const AliRsnAnalysisPhi7TeVNoPID& copy) :
72 AliAnalysisTaskSE(copy),
73 fUseMC(copy.fUseMC),
74 fPDG(0),
75 fCh(0),
76 fIM(0.0),
77 fPt(0.0),
78 fY(0.0),
79 fEta(0.0),
80 fMaxVz(copy.fMaxVz),
81 fMaxITSband(copy.fMaxITSband),
82 fTPCpLimit(copy.fTPCpLimit),
83 fMinTPCband(copy.fMinTPCband),
84 fMaxTPCband(copy.fMaxTPCband),
85 fRsnTreeComp(0x0),
86 fRsnTreeTrue(0x0),
87 fOutList(0x0),
88 fHEvents(0x0),
89 fESDtrackCutsTPC(copy.fESDtrackCutsTPC),
90 fESDtrackCutsITS(copy.fESDtrackCutsITS),
91 fESDpid(0x0),
92 fTOFmaker(0x0),
93 fTOFcalib(0x0),
94 fTOFcalibrateESD(kFALSE),
95 fTOFcorrectTExp(kFALSE),
96 fTOFuseT0(kFALSE),
97 fTOFtuneMC(kFALSE),
98 fTOFresolution(0.0)
99{
100//
101// Copy constructor
102//
103}
104
105//__________________________________________________________________________________________________
106AliRsnAnalysisPhi7TeVNoPID& AliRsnAnalysisPhi7TeVNoPID::operator=(const AliRsnAnalysisPhi7TeVNoPID& copy)
107{
108//
109// Assignment operator
110//
111
112 fUseMC = copy.fUseMC;
113
114 fMaxVz = copy.fMaxVz;
115 fMaxITSband = copy.fMaxITSband;
116
117 fTPCpLimit = copy.fTPCpLimit;
118 fMinTPCband = copy.fMinTPCband;
119 fMaxTPCband = copy.fMaxTPCband;
120
121 fESDtrackCutsTPC = copy.fESDtrackCutsTPC;
122 fESDtrackCutsITS = copy.fESDtrackCutsITS;
123
124 fTOFcalibrateESD = copy.fTOFcalibrateESD;
125 fTOFcorrectTExp = copy.fTOFcorrectTExp;
126 fTOFuseT0 = copy.fTOFuseT0;
127 fTOFtuneMC = copy.fTOFtuneMC;
128 fTOFresolution = copy.fTOFresolution;
129
130 return (*this);
131}
132
133//__________________________________________________________________________________________________
134AliRsnAnalysisPhi7TeVNoPID::~AliRsnAnalysisPhi7TeVNoPID()
135{
136//
137// Destructor
138//
139
140 if (fRsnTreeComp) delete fRsnTreeComp;
141 if (fRsnTreeTrue) delete fRsnTreeTrue;
142 if (fHEvents) delete fHEvents;
143 if (fESDpid) delete fESDpid;
144}
145
146//__________________________________________________________________________________________________
147void AliRsnAnalysisPhi7TeVNoPID::UserCreateOutputObjects()
148{
149//
150// Create the output data container
151//
152
153 // setup TPC response
154 fESDpid = new AliESDpid;
155 fESDpid->GetTPCResponse().SetBetheBlochParameters(fTPCpar[0],fTPCpar[1],fTPCpar[2],fTPCpar[3],fTPCpar[4]);
156
157 // setup TOF maker & calibration
158 fTOFcalib = new AliTOFcalib;
159 fTOFmaker = new AliTOFT0maker(fESDpid, fTOFcalib);
160 fTOFmaker->SetTimeResolution(fTOFresolution);
161
162 // initialize random
163 gRandom->SetSeed(0);
164
165 // create output trees
166 OpenFile(1);
167 fRsnTreeComp = new TTree("rsnTree", "Pairs");
168
169 fRsnTreeComp->Branch("pdg", &fPDG, "pdg/S" );
170 fRsnTreeComp->Branch("ch" , &fCh , "ch/S" );
171 fRsnTreeComp->Branch("im" , &fIM , "im/F" );
172 fRsnTreeComp->Branch("y" , &fY , "y/F" );
173 fRsnTreeComp->Branch("pt" , &fPt , "pt/F" );
174 fRsnTreeComp->Branch("eta", &fEta, "eta/F" );
175 fRsnTreeComp->Branch("its", &fITS, "its[2]/S");
176
177 fRsnTreeComp->Branch("p" , &fP , "p[2]/F");
178 fRsnTreeComp->Branch("ptpc" , &fPTPC , "ptpc[2]/F");
179 fRsnTreeComp->Branch("tpcpid", &fTPCnsigma, "tpcpid[2]/F");
180 fRsnTreeComp->Branch("itspid", &fITSnsigma, "itspid[2]/F");
181 fRsnTreeComp->Branch("tofpid", &fTOFdiff , "tofpid[2]/F");
182
183 OpenFile(2);
184 fRsnTreeTrue = new TTree("rsnTrue", "True pairs");
185
186 fRsnTreeTrue->Branch("im" , &fIM , "im/F" );
187 fRsnTreeTrue->Branch("y" , &fY , "y/F" );
188 fRsnTreeTrue->Branch("pt" , &fPt , "pt/F" );
189 fRsnTreeTrue->Branch("eta", &fEta, "eta/F");
190
191 OpenFile(3);
192 fOutList = new TList;
193 fHEvents = new TH1I("hEvents", "Event details", 5, 0, 5);
194 fVertexX[0] = new TH1F("hVertexTracksX", "X position of primary vertex (tracks)", 200, -2, 2);
195 fVertexY[0] = new TH1F("hVertexTracksY", "Y position of primary vertex (tracks)", 200, -2, 2);
196 fVertexZ[0] = new TH1F("hVertexTracksZ", "Z position of primary vertex (tracks)", 400, -40, 40);
197 fVertexX[1] = new TH1F("hVertexSPDX", "X position of primary vertex (SPD)", 200, -2, 2);
198 fVertexY[1] = new TH1F("hVertexSPDY", "Y position of primary vertex (SPD)", 200, -2, 2);
199 fVertexZ[1] = new TH1F("hVertexSPDZ", "Z position of primary vertex (SPD)", 400, -40, 40);
200
201 fHEvents->GetXaxis()->SetBinLabel(1, "Good vertex with tracks");
202 fHEvents->GetXaxis()->SetBinLabel(2, "Good vertex with SPD");
203 fHEvents->GetXaxis()->SetBinLabel(3, "Far vertex with tracks");
204 fHEvents->GetXaxis()->SetBinLabel(4, "Far vertex with SPD");
205 fHEvents->GetXaxis()->SetBinLabel(5, "No good vertex");
206
207 fOutList->Add(fHEvents);
208 fOutList->Add(fVertexX[0]);
209 fOutList->Add(fVertexY[0]);
210 fOutList->Add(fVertexZ[0]);
211 fOutList->Add(fVertexX[1]);
212 fOutList->Add(fVertexY[1]);
213 fOutList->Add(fVertexZ[1]);
214}
215
216//__________________________________________________________________________________________________
217void AliRsnAnalysisPhi7TeVNoPID::UserExec(Option_t *)
218{
219//
220// Main execution function.
221// Fills the fHEvents data member with the following legenda:
222// 0 -- event OK, prim vertex with tracks
223// 1 -- event OK, prim vertex with SPD
224// 2 -- event OK but vz large
225// 3 -- event bad
226//
227
228 static Int_t evNum = 0;
229 evNum++;
230
231 // retrieve ESD event and related stack (if available)
232 AliESDEvent *esd = dynamic_cast<AliESDEvent*>(fInputEvent);
233 AliStack *stack = (fMCEvent ? fMCEvent->Stack() : 0x0);
234
235 // check the event
236 Int_t eval = EventEval(esd);
237 fHEvents->Fill(eval);
238
239 // if the event is good for analysis, process it
240 if (eval == kGoodTracksPrimaryVertex || eval == kGoodSPDPrimaryVertex)
241 {
242 ProcessESD(esd, stack);
243 ProcessMC(stack);
244 }
245
246 // update histogram container
247 PostData(3, fOutList);
248}
249
250//__________________________________________________________________________________________________
251void AliRsnAnalysisPhi7TeVNoPID::Terminate(Option_t *)
252{
253//
254// Terminate
255//
256}
257
258//__________________________________________________________________________________________________
259Int_t AliRsnAnalysisPhi7TeVNoPID::EventEval(AliESDEvent *esd)
260{
261//
262// Checks if the event is good for analysis.
263// Returns one of the flag values defined in the header
264//
265
266 static Int_t evNum = 0;
267 evNum++;
268
269 // debug message
270 AliDebug(AliLog::kDebug + 1, Form("Event %d -- number of tracks = %d", evNum, esd->GetNumberOfTracks()));
271
272 // get the best primary vertex:
273 // first try the one with tracks
274 const AliESDVertex *vTrk = esd->GetPrimaryVertexTracks();
275 const AliESDVertex *vSPD = esd->GetPrimaryVertexSPD();
276 Double_t vzTrk = 1000.0;
277 Double_t vzSPD = 1000.0;
278 if (vTrk) vzTrk = TMath::Abs(vTrk->GetZv());
279 if (vSPD) vzSPD = TMath::Abs(vSPD->GetZv());
280 AliDebug(AliLog::kDebug + 1, Form("Event %d -- vertex with tracks: contributors = %d, abs(vz) = %f", evNum, vTrk->GetNContributors(), vzTrk));
281 AliDebug(AliLog::kDebug + 1, Form("Event %d -- vertex with SPD, contributors = %d, abs(vz) = %f", evNum, vSPD->GetNContributors(), vzSPD));
282 if(vTrk->GetNContributors() > 0)
283 {
284 // fill the histograms
285 fVertexX[0]->Fill(vTrk->GetXv());
286 fVertexY[0]->Fill(vTrk->GetYv());
287 fVertexZ[0]->Fill(vTrk->GetZv());
288
289 // check VZ position
290 if (vzTrk <= fMaxVz)
291 return kGoodTracksPrimaryVertex;
292 else
293 return kFarTracksPrimaryVertex;
294 }
295 else if (vSPD->GetNContributors() > 0)
296 {
297 // fill the histograms
298 fVertexX[1]->Fill(vSPD->GetXv());
299 fVertexY[1]->Fill(vSPD->GetYv());
300 fVertexZ[1]->Fill(vSPD->GetZv());
301
302 // check VZ position
303 if (vzSPD <= fMaxVz)
304 return kGoodSPDPrimaryVertex;
305 else
306 return kFarSPDPrimaryVertex;
307 }
308 else
309 return kNoGoodPrimaryVertex;
310}
311
312//__________________________________________________________________________________________________
313void AliRsnAnalysisPhi7TeVNoPID::ProcessESD
314(AliESDEvent *esd, AliStack *stack)
315{
316//
317// This function works with the ESD object
318//
319
320 // ITS stuff #1 create the response function
321 Bool_t isMC = (stack != 0x0);
322 AliITSPIDResponse itsrsp(isMC);
323
324 // TOF stuff #1: init OCDB
325 Int_t run = esd->GetRunNumber();
326 AliCDBManager *cdb = AliCDBManager::Instance();
327 cdb->SetDefaultStorage("raw://");
328 cdb->SetRun(run);
329 // TOF stuff #2: init calibration
330 fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
331 fTOFcalib->Init();
332 // TOF stuff #3: calibrate
333 if (fTOFcalibrateESD) fTOFcalib->CalibrateESD(esd);
334 if (fTOFtuneMC) fTOFmaker->TuneForMC(esd);
335 if (fTOFuseT0)
336 {
337 fTOFmaker->ComputeT0TOF(esd);
338 fTOFmaker->ApplyT0TOF(esd);
339 fESDpid->MakePID(esd, kFALSE, 0.);
340 }
341
342 // prepare to look on all tracks to select the ones
343 // which pass all the cuts
344 Int_t ntracks = esd->GetNumberOfTracks();
345 TArrayI pos(ntracks);
346 TArrayI neg(ntracks);
347 TArrayI itspos(ntracks);
348 TArrayI itsneg(ntracks);
349
350 // loop on all tracks
351 ULong_t status;
352 Int_t i, k, charge, npos = 0, nneg = 0, nITS;
353 Double_t times[10], itsSignal, mom, tofTime, tofRef;
354 Bool_t isTPC, isITSSA;
355 UChar_t itsCluMap;
356 for (i = 0; i < ntracks; i++)
357 {
358 AliESDtrack *track = esd->GetTrack(i);
359 if (!track) continue;
360
361 // get commonly used variables
362 status = (ULong_t)track->GetStatus();
363 mom = track->P();
364 isTPC = ((status & AliESDtrack::kTPCin) != 0);
365 isITSSA = ((status & AliESDtrack::kTPCin) == 0 && (status & AliESDtrack::kITSrefit) != 0 && (status & AliESDtrack::kITSpureSA) == 0 && (status & AliESDtrack::kITSpid) != 0);
366
367 // accept only tracks which are TPC+ITS or ITS standalone
368 if (!isTPC && !isITSSA) continue;
369
370 // check specific cuts for TPC and ITS-SA tracks
371 if (isTPC)
372 {
373 if (!fESDtrackCutsTPC.IsSelected(track)) continue;
374 }
375 else
376 {
377 if (!fESDtrackCutsITS.IsSelected(track)) continue;
378 itsCluMap = track->GetITSClusterMap();
379 nITS = 0;
380 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
381 if (nITS < 3) continue;
382 }
383
384 // if all checks are passed, add the track index in one of the
385 // charged tracks arrays
386 charge = (Int_t)track->Charge();
387 if (charge > 0)
388 {
389 pos[npos] = i;
390 if (isITSSA) itspos[npos] = 1; else itspos[npos] = 0;
391 npos++;
392 }
393 else if (charge < 0)
394 {
395 neg[nneg] = i;
396 if (isITSSA) itsneg[nneg] = 1; else itsneg[nneg] = 0;
397 nneg++;
398 }
399 }
400
401 // resize arrays accordingly
402 pos.Set(npos);
403 neg.Set(nneg);
404 itspos.Set(npos);
405 itsneg.Set(nneg);
406
407 // loop on unlike-sign pairs to compute invariant mass signal
408 Int_t ip, in, lp, ln;
409 AliPID pid;
410 Double_t kmass = pid.ParticleMass(AliPID::kKaon);
411 Double_t phimass = 1.019455;
412 TParticle *partp = 0x0, *partn = 0x0;
413 AliESDtrack *tp = 0x0, *tn = 0x0;
414 TLorentzVector vp, vn, vsum, vref;
415 for (ip = 0; ip < npos; ip++)
416 {
417 tp = esd->GetTrack(pos[ip]);
418 lp = TMath::Abs(tp->GetLabel());
419 if (stack) partp = stack->Particle(lp);
420
421 for (in = 0; in < nneg; in++)
422 {
423 if (pos[ip] == neg[in])
424 {
425 AliError("POS = NEG");
426 continue;
427 }
428 tn = esd->GetTrack(neg[in]);
429 ln = TMath::Abs(tn->GetLabel());
430 if (stack) partn = stack->Particle(ln);
431
432 fPDG = 0;
433 if (partp && partn)
434 {
435 if (partp->GetFirstMother() == partn->GetFirstMother())
436 {
437 if (partp->GetFirstMother() > 0)
438 {
439 TParticle *mum = stack->Particle(partp->GetFirstMother());
440 fPDG = mum->GetPdgCode();
441 }
442 }
443 }
444 fPDG = TMath::Abs(fPDG);
445
446 vp.SetXYZM(tp->Px(), tp->Py(), tp->Pz(), kmass);
447 vn.SetXYZM(tn->Px(), tn->Py(), tn->Pz(), kmass);
448 vsum = vp + vn;
449 vref.SetXYZM(vsum.X(), vsum.Y(), vsum.Z(), phimass);
450
451 fCh = 0;
452 fIM = (Float_t)vsum.M();
453 fPt = (Float_t)vsum.Perp();
454 fEta = (Float_t)vsum.Eta();
455 fY = (Float_t)vref.Rapidity();
456 fITS[0] = itspos[ip];
457 fITS[1] = itsneg[in];
458
459 if (fIM < 0.9 || fIM > 5.0) continue;
460 if (fPt < 0.0 || fPt > 20.0) continue;
461
462 // PID signal for track #1
463 // here it is enough to check if it is a TPC track
464 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
465 if ((tp->GetStatus() & AliESDtrack::kTPCin) != 0)
466 {
467 fP [0] = tp->P();
468 fPTPC [0] = tp->GetInnerParam()->P();
469 fTPCnsigma[0] = TMath::Abs(fESDpid->NumberOfSigmasTPC(tp, AliPID::kKaon));
470 fITSnsigma[0] = 1E6;
471 fTOFdiff [0] = 1E6;
472 // check TOF (only if momentum is large than function asymptote and flags are OK)
473 if (((tp->GetStatus() & AliESDtrack::kTOFout) != 0) && ((tp->GetStatus() & AliESDtrack::kTIME) != 0))
474 {
475 tp->GetIntegratedTimes(times);
476 tofTime = (Double_t)tp->GetTOFsignal();
477 tofRef = times[AliPID::kKaon];
478 if (tofRef > 0.0) fTOFdiff[0] = (tofTime - tofRef) / tofRef;
479 }
480 }
481 else
482 {
483 fP [0] = tp->P();
484 fPTPC [0] = 1E6;
485 fTPCnsigma[0] = 1E6;
486 fTOFdiff [0] = 1E6;
487 // check dE/dx
488 itsSignal = tp->GetITSsignal();
489 itsCluMap = tp->GetITSClusterMap();
490 nITS = 0;
491 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
492 fITSnsigma[0] = itsrsp.GetNumberOfSigmas(fP[0], itsSignal, AliPID::kKaon, nITS, kTRUE);
493 }
494
495 // PID signal for track #2
496 // here it is enough to check if it is a TPC track
497 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
498 if ((tp->GetStatus() & AliESDtrack::kTPCin) != 0)
499 {
500 fP [1] = tn->P();
501 fPTPC [1] = tn->GetInnerParam()->P();
502 fTPCnsigma[1] = TMath::Abs(fESDpid->NumberOfSigmasTPC(tn, AliPID::kKaon));
503 fITSnsigma[1] = 1E6;
504 fTOFdiff [1] = 1E6;
505 // check TOF (only if momentum is large than function asymptote and flags are OK)
506 if (((tn->GetStatus() & AliESDtrack::kTOFout) != 0) && ((tn->GetStatus() & AliESDtrack::kTIME) != 0))
507 {
508 tn->GetIntegratedTimes(times);
509 tofTime = (Double_t)tn->GetTOFsignal();
510 tofRef = times[AliPID::kKaon];
511 if (tofRef > 0.0) fTOFdiff[1] = (tofTime - tofRef) / tofRef;
512 }
513 }
514 else
515 {
516 fP [1] = tn->P();
517 fPTPC [1] = 1E6;
518 fTPCnsigma[1] = 1E6;
519 fTOFdiff [1] = 1E6;
520 // check dE/dx
521 itsSignal = tn->GetITSsignal();
522 itsCluMap = tn->GetITSClusterMap();
523 nITS = 0;
524 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
525 fITSnsigma[1] = itsrsp.GetNumberOfSigmas(fP[1], itsSignal, AliPID::kKaon, nITS, kTRUE);
526 }
527
528 fRsnTreeComp->Fill();
529 }
530 }
531
532 // loop on like-sign pairs to compute invariant mass background
533 Int_t i1, i2;
534 AliESDtrack *t1 = 0x0, *t2 = 0x0;
535 TLorentzVector v1, v2;
536
537 // pos-pos
538 for (i1 = 0; i1 < npos; i1++)
539 {
540 t1 = esd->GetTrack(pos[i1]);
541
542 for (i2 = i1+1; i2 < npos; i2++)
543 {
544 t2 = esd->GetTrack(pos[i2]);
545
546 v1.SetXYZM(t1->Px(), t1->Py(), t1->Pz(), kmass);
547 v2.SetXYZM(t2->Px(), t2->Py(), t2->Pz(), kmass);
548 vsum = v1 + v2;
549 vref.SetXYZM(vsum.X(), vsum.Y(), vsum.Z(), phimass);
550
551 fPDG = 0;
552 fCh = 1;
553 fIM = (Float_t)vsum.M();
554 fPt = (Float_t)vsum.Perp();
555 fEta = (Float_t)vsum.Eta();
556 fY = (Float_t)vref.Rapidity();
557 fITS[0] = itspos[i1];
558 fITS[1] = itspos[i2];
559
560 if (fIM < 0.9 || fIM > 5.0) continue;
561 if (fPt < 0.0 || fPt > 20.0) continue;
562
563 // PID signal for track #1
564 // here it is enough to check if it is a TPC track
565 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
566 if ((t1->GetStatus() & AliESDtrack::kTPCin) != 0)
567 {
568 fP [0] = t1->P();
569 fPTPC [0] = t1->GetInnerParam()->P();
570 fTPCnsigma[0] = TMath::Abs(fESDpid->NumberOfSigmasTPC(t1, AliPID::kKaon));
571 fITSnsigma[0] = 1E6;
572 fTOFdiff [0] = 1E6;
573 // check TOF (only if momentum is large than function asymptote and flags are OK)
574 if (((t1->GetStatus() & AliESDtrack::kTOFout) != 0) && ((t1->GetStatus() & AliESDtrack::kTIME) != 0))
575 {
576 t1->GetIntegratedTimes(times);
577 tofTime = (Double_t)t1->GetTOFsignal();
578 tofRef = times[AliPID::kKaon];
579 if (tofRef > 0.0) fTOFdiff[0] = (tofTime - tofRef) / tofRef;
580 }
581 }
582 else
583 {
584 fP [0] = t1->P();
585 fPTPC [0] = 1E6;
586 fTPCnsigma[0] = 1E6;
587 fTOFdiff [0] = 1E6;
588 // check dE/dx
589 itsSignal = t1->GetITSsignal();
590 itsCluMap = t1->GetITSClusterMap();
591 nITS = 0;
592 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
593 fITSnsigma[0] = itsrsp.GetNumberOfSigmas(fP[0], itsSignal, AliPID::kKaon, nITS, kTRUE);
594 }
595
596 // PID signal for track #2
597 // here it is enough to check if it is a TPC track
598 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
599 if ((t1->GetStatus() & AliESDtrack::kTPCin) != 0)
600 {
601 fP [1] = t2->P();
602 fPTPC [1] = t2->GetInnerParam()->P();
603 fTPCnsigma[1] = TMath::Abs(fESDpid->NumberOfSigmasTPC(t2, AliPID::kKaon));
604 fITSnsigma[1] = 1E6;
605 fTOFdiff [1] = 1E6;
606 // check TOF (only if momentum is large than function asymptote and flags are OK)
607 if (((t2->GetStatus() & AliESDtrack::kTOFout) != 0) && ((t2->GetStatus() & AliESDtrack::kTIME) != 0))
608 {
609 t2->GetIntegratedTimes(times);
610 tofTime = (Double_t)t2->GetTOFsignal();
611 tofRef = times[AliPID::kKaon];
612 if (tofRef > 0.0) fTOFdiff[1] = (tofTime - tofRef) / tofRef;
613 }
614 }
615 else
616 {
617 fP [1] = t2->P();
618 fPTPC [1] = 1E6;
619 fTPCnsigma[1] = 1E6;
620 fTOFdiff [1] = 1E6;
621 // check dE/dx
622 itsSignal = t2->GetITSsignal();
623 itsCluMap = t2->GetITSClusterMap();
624 nITS = 0;
625 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
626 fITSnsigma[1] = itsrsp.GetNumberOfSigmas(fP[1], itsSignal, AliPID::kKaon, nITS, kTRUE);
627 }
628
629 //fRsnTreeComp->Fill();
630 }
631 }
632 // neg-neg
633 for (i1 = 0; i1 < nneg; i1++)
634 {
635 t1 = esd->GetTrack(neg[i1]);
636
637 for (i2 = i1+1; i2 < nneg; i2++)
638 {
639 t2 = esd->GetTrack(neg[i2]);
640
641 v1.SetXYZM(t1->Px(), t1->Py(), t1->Pz(), kmass);
642 v2.SetXYZM(t2->Px(), t2->Py(), t2->Pz(), kmass);
643 vsum = v1 + v2;
644 vref.SetXYZM(vsum.X(), vsum.Y(), vsum.Z(), phimass);
645
646 fPDG = 0;
647 fCh = -1;
648 fIM = (Float_t)vsum.M();
649 fPt = (Float_t)vsum.Perp();
650 fEta = (Float_t)vsum.Eta();
651 fY = (Float_t)vref.Rapidity();
652 fITS[0] = itsneg[i1];
653 fITS[1] = itsneg[i2];
654
655 if (fIM < 0.9 || fIM > 5.0) continue;
656 if (fPt < 0.0 || fPt > 20.0) continue;
657
658 // PID signal for track #1
659 // here it is enough to check if it is a TPC track
660 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
661 if ((t1->GetStatus() & AliESDtrack::kTPCin) != 0)
662 {
663 fP [0] = t1->P();
664 fPTPC [0] = t1->GetInnerParam()->P();
665 fTPCnsigma[0] = TMath::Abs(fESDpid->NumberOfSigmasTPC(t1, AliPID::kKaon));
666 fITSnsigma[0] = 1E6;
667 fTOFdiff [0] = 1E6;
668 // check TOF (only if momentum is large than function asymptote and flags are OK)
669 if (((t1->GetStatus() & AliESDtrack::kTOFout) != 0) && ((t1->GetStatus() & AliESDtrack::kTIME) != 0))
670 {
671 t1->GetIntegratedTimes(times);
672 tofTime = (Double_t)t1->GetTOFsignal();
673 tofRef = times[AliPID::kKaon];
674 if (tofRef > 0.0) fTOFdiff[0] = (tofTime - tofRef) / tofRef;
675 }
676 }
677 else
678 {
679 fP [0] = t1->P();
680 fPTPC [0] = 1E6;
681 fTPCnsigma[0] = 1E6;
682 fTOFdiff [0] = 1E6;
683 // check dE/dx
684 itsSignal = t1->GetITSsignal();
685 itsCluMap = t1->GetITSClusterMap();
686 nITS = 0;
687 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
688 fITSnsigma[0] = itsrsp.GetNumberOfSigmas(fP[0], itsSignal, AliPID::kKaon, nITS, kTRUE);
689 }
690
691 // PID signal for track #2
692 // here it is enough to check if it is a TPC track
693 // since we excluded the case that it is neither a TPC+ITS nor an ITS-SA
694 if ((t1->GetStatus() & AliESDtrack::kTPCin) != 0)
695 {
696 fP [1] = t2->P();
697 fPTPC [1] = t2->GetInnerParam()->P();
698 fTPCnsigma[1] = TMath::Abs(fESDpid->NumberOfSigmasTPC(t2, AliPID::kKaon));
699 fITSnsigma[1] = 1E6;
700 fTOFdiff [1] = 1E6;
701 // check TOF (only if momentum is large than function asymptote and flags are OK)
702 if (((t2->GetStatus() & AliESDtrack::kTOFout) != 0) && ((t2->GetStatus() & AliESDtrack::kTIME) != 0))
703 {
704 t2->GetIntegratedTimes(times);
705 tofTime = (Double_t)t2->GetTOFsignal();
706 tofRef = times[AliPID::kKaon];
707 if (tofRef > 0.0) fTOFdiff[1] = (tofTime - tofRef) / tofRef;
708 }
709 }
710 else
711 {
712 fP [1] = t2->P();
713 fPTPC [1] = 1E6;
714 fTPCnsigma[1] = 1E6;
715 fTOFdiff [1] = 1E6;
716 // check dE/dx
717 itsSignal = t2->GetITSsignal();
718 itsCluMap = t2->GetITSClusterMap();
719 nITS = 0;
720 for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
721 fITSnsigma[1] = itsrsp.GetNumberOfSigmas(fP[1], itsSignal, AliPID::kKaon, nITS, kTRUE);
722 }
723
724 //fRsnTreeComp->Fill();
725 }
726 }
727
728 PostData(1, fRsnTreeComp);
729}
730
731//__________________________________________________________________________________________________
732void AliRsnAnalysisPhi7TeVNoPID::ProcessMC(AliStack *stack)
733{
734//
735// Function to process stack only
736//
737
738 if (!stack) return;
739 Int_t nPart = stack->GetNtrack();
740
741 // loop to compute invariant mass
742 Int_t ip, in;
743 AliPID pid;
744 Double_t kmass = pid.ParticleMass(AliPID::kKaon);
745 Double_t phimass = 1.019455;
746 TParticle *partp = 0x0, *partn = 0x0;
747 TLorentzVector vp, vn, vsum, vref;
748
749 for (ip = 0; ip < nPart; ip++)
750 {
751 partp = stack->Particle(ip);
752 if (partp->GetPdgCode() != 321) continue;
753
754 for (in = 0; in < nPart; in++)
755 {
756 partn = stack->Particle(in);
757 if (partn->GetPdgCode() != -321) continue;
758
759 fPDG = 0;
760 if (partp->GetFirstMother() == partn->GetFirstMother())
761 {
762 if (partp->GetFirstMother() > 0)
763 {
764 TParticle *mum = stack->Particle(partp->GetFirstMother());
765 fPDG = mum->GetPdgCode();
766 }
767 }
768 fPDG = TMath::Abs(fPDG);
769 if (fPDG != 333) continue;
770
771 vp.SetXYZM(partp->Px(), partp->Py(), partp->Pz(), kmass);
772 vn.SetXYZM(partn->Px(), partn->Py(), partn->Pz(), kmass);
773 vsum = vp + vn;
774 vref.SetXYZM(vsum.X(), vsum.Y(), vsum.Z(), phimass);
775
776 fIM = (Float_t)vsum.M();
777 fPt = (Float_t)vsum.Perp();
778 fEta = (Float_t)vsum.Eta();
779 fY = (Float_t)vref.Rapidity();
780
781 fRsnTreeTrue->Fill();
782 }
783 }
784
785 PostData(2, fRsnTreeTrue);
786}