2 // Implementation file for implementation of data analysis aft 900 GeV
4 // Author: A. Pulvirenti
12 #include "TParticle.h"
14 #include "TLorentzVector.h"
17 #include "AliESDpid.h"
18 #include "AliESDEvent.h"
19 #include "AliESDVertex.h"
20 #include "AliESDtrack.h"
22 #include "AliMCEvent.h"
23 #include "AliTOFT0maker.h"
24 #include "AliTOFcalib.h"
25 #include "AliCDBManager.h"
26 #include "AliITSPIDResponse.h"
27 #include "AliRsnMonitorTrack.h"
28 #include "AliRsnDaughter.h"
29 #include "AliRsnEvent.h"
31 #include "AliRsnAnalysisMonitorPairTask.h"
33 ClassImp(AliRsnAnalysisMonitorPairTask)
35 //__________________________________________________________________________________________________
36 AliRsnAnalysisMonitorPairTask::AliRsnAnalysisMonitorPairTask(const char *name) :
37 AliAnalysisTaskSE(name),
45 fTOFcalibrateESD(kFALSE),
46 fTOFcorrectTExp(kFALSE),
50 fEventCuts("eventCuts", AliRsnCut::kEvent),
51 fTrackCuts("trackCuts", AliRsnCut::kDaughter)
58 DefineOutput(1, TTree::Class());
61 //__________________________________________________________________________________________________
62 AliRsnAnalysisMonitorPairTask::AliRsnAnalysisMonitorPairTask(const AliRsnAnalysisMonitorPairTask& copy) :
63 AliAnalysisTaskSE(copy),
66 fRangeMin(copy.fRangeMin),
67 fRangeMax(copy.fRangeMax),
71 fTOFcalibrateESD(kFALSE),
72 fTOFcorrectTExp(kFALSE),
76 fEventCuts(copy.fEventCuts),
77 fTrackCuts(copy.fTrackCuts)
84 //__________________________________________________________________________________________________
85 AliRsnAnalysisMonitorPairTask& AliRsnAnalysisMonitorPairTask::operator=(const AliRsnAnalysisMonitorPairTask& copy)
88 // Assignment operator
91 fMass[0] = copy.fMass[0];
92 fMass[1] = copy.fMass[1];
94 fRangeMin = copy.fRangeMin;
95 fRangeMax = copy.fRangeMax;
97 fTOFcalibrateESD = copy.fTOFcalibrateESD;
98 fTOFcorrectTExp = copy.fTOFcorrectTExp;
99 fTOFuseT0 = copy.fTOFuseT0;
100 fTOFtuneMC = copy.fTOFtuneMC;
101 fTOFresolution = copy.fTOFresolution;
106 //__________________________________________________________________________________________________
107 AliRsnAnalysisMonitorPairTask::~AliRsnAnalysisMonitorPairTask()
113 if (fOut) delete fOut;
114 if (fESDpid) delete fESDpid;
115 if (fTrack) delete fTrack;
118 //__________________________________________________________________________________________________
119 void AliRsnAnalysisMonitorPairTask::UserCreateOutputObjects()
122 // Create the output data container
125 // setup TPC response
126 fESDpid = new AliESDpid;
127 fESDpid->GetTPCResponse().SetBetheBlochParameters(fTPCpar[0], fTPCpar[1], fTPCpar[2], fTPCpar[3], fTPCpar[4]);
129 // setup TOF maker & calibration
130 fTOFcalib = new AliTOFcalib;
131 fTOFmaker = new AliTOFT0maker(fESDpid, fTOFcalib);
132 fTOFmaker->SetTimeResolution(fTOFresolution);
134 // create output branch object
135 fTrack[0] = new AliRsnMonitorTrack;
136 fTrack[1] = new AliRsnMonitorTrack;
138 // create output tree
140 fOut = new TTree("rsnPairMonitor", "Informations on pairs for cut checking");
141 fOut->Branch("track1", "AliRsnMonitorTrack", &fTrack[0]);
142 fOut->Branch("track2", "AliRsnMonitorTrack", &fTrack[1]);
143 fOut->Branch("minv" , &fInvMass, "minv/F");
146 //__________________________________________________________________________________________________
147 void AliRsnAnalysisMonitorPairTask::UserExec(Option_t *)
150 // Main execution function.
151 // Fills the fHEvents data member with the following legenda:
152 // 0 -- event OK, prim vertex with tracks
153 // 1 -- event OK, prim vertex with SPD
154 // 2 -- event OK but vz large
158 // retrieve ESD event and related stack (if available)
159 AliESDEvent *esd = dynamic_cast<AliESDEvent*>(fInputEvent);
160 AliStack *stack = 0x0;
164 if (fMCEvent) stack = fMCEvent->Stack();
166 // create interface objects to AliRsnEvent to check event cuts
169 event.SetRefMC(fMCEvent);
170 if (!fEventCuts.IsSelected(&event)) return;
173 Int_t type = EventEval(esd);
175 // if processable, then process it
176 if (type == 0) ProcessESD(esd, esd->GetPrimaryVertexTracks());
177 else if (type == 1) ProcessESD(esd, esd->GetPrimaryVertexSPD());
180 // update histogram container
184 //__________________________________________________________________________________________________
185 void AliRsnAnalysisMonitorPairTask::Terminate(Option_t *)
192 //__________________________________________________________________________________________________
193 Int_t AliRsnAnalysisMonitorPairTask::EventEval(AliESDEvent *esd)
196 // Checks if the event is good for analysis.
198 // ---> 0 if a good primary vertex with tracks was found,
199 // ---> 1 if a good SPD primary vertex was found
200 // ---> 2 otherwise (event to be rejected)
201 // In any case, adds an entry to the TTree, to keep trace of all events.
204 // get the best primary vertex:
205 // first try that with tracks, then the SPD one
206 const AliESDVertex *vTrk = esd->GetPrimaryVertexTracks();
207 const AliESDVertex *vSPD = esd->GetPrimaryVertexSPD();
208 if(vTrk->GetNContributors() > 0)
212 else if (vSPD->GetNContributors() > 0)
222 //__________________________________________________________________________________________________
223 void AliRsnAnalysisMonitorPairTask::ProcessESD(AliESDEvent *esd, const AliESDVertex *v)
226 // Process the ESD container, to read all tracks and copy their useful values.
227 // All info are stored into an AliRsnMonitorTrack object and saved into the
228 // TClonesArray which is one of the branches of the output TTree.
231 // TOF stuff #1: init OCDB
232 Int_t run = esd->GetRunNumber();
233 AliCDBManager *cdb = AliCDBManager::Instance();
234 cdb->SetDefaultStorage("raw://");
236 // TOF stuff #2: init calibration
237 fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
239 // TOF stuff #3: calibrate
240 if (fTOFcalibrateESD) fTOFcalib->CalibrateESD(esd);
241 if (fTOFtuneMC) fTOFmaker->TuneForMC(esd);
244 fTOFmaker->ComputeT0TOF(esd);
245 fTOFmaker->ApplyT0TOF(esd);
246 fESDpid->MakePID(esd, kFALSE, 0.);
249 // loop on all tracks
250 Int_t i0 , i1, ch0, ch1, ntracks = esd->GetNumberOfTracks();
251 TLorentzVector v0, v1, sum;
253 // two nested loops on tracks
254 for (i0 = 0; i0 < ntracks; i0++)
256 if (!ProcessTrack(0, i0, esd, v)) continue;
257 ch0 = fTrack[0]->Charge();
258 v0.SetXYZM(fTrack[0]->PrecX(), fTrack[0]->PrecY(), fTrack[0]->PrecZ(), fMass[0]);
260 for (i1 = i0 + 1; i1 < ntracks; i1++)
262 if (!ProcessTrack(1, i1, esd, v)) continue;
263 ch1 = fTrack[1]->Charge();
265 // skip like-sign pairs
266 if (ch1 == ch0) continue;
268 // check invmass range
269 v1.SetXYZM(fTrack[1]->PrecX(), fTrack[1]->PrecY(), fTrack[1]->PrecZ(), fMass[1]);
272 if (fInvMass < fRangeMin || fInvMass > fRangeMax) continue;
274 // if here, add an entry to the tree
280 //__________________________________________________________________________________________________
281 Bool_t AliRsnAnalysisMonitorPairTask::ProcessTrack(Int_t myIndex, Int_t esdIndex, AliESDEvent *esd, const AliESDVertex *v)
284 // Process a single track, and stores all its info into one of the two branch objects
285 // defined by the first index, while the second chooses the track in the owner ESD
288 if (myIndex < 0 || myIndex > 1) return kFALSE;
291 // create the response function and initialize it to MC or not
292 // depending if the AliStack object is there or not
293 AliStack *stack = 0x0;
294 if (fMCEvent) stack = fMCEvent->Stack();
295 Bool_t isMC = (stack != 0x0);
296 AliITSPIDResponse itsrsp(isMC);
298 // create interfacr objects
301 event.SetRefMC(fMCEvent);
302 AliRsnDaughter daughter;
303 AliESDtrack *track = esd->GetTrack(esdIndex);
304 event.SetDaughter(daughter, esdIndex, AliRsnDaughter::kTrack);
306 // skip NULL pointers, kink daughters and tracks which
307 // cannot be propagated to primary vertex
308 if (!track) return kFALSE;
309 if ((Int_t)track->GetKinkIndex(0) > 0) return kFALSE;
310 if (!track->RelateToVertex(v, esd->GetMagneticField(), kVeryBig)) return kFALSE;
314 Bool_t isTPC, isITSSA, isTOF;
315 Float_t b[2], bCov[3];
318 // reset the output object
319 // 'usable' flag will need to be set to 'ok'
320 fTrack[myIndex]->Reset();
323 fTrack[myIndex]->CutsPassed() = fTrackCuts.IsSelected(&daughter);
325 // get MC info if possible
326 if (stack) fTrack[myIndex]->AdoptMC(TMath::Abs(track->GetLabel()), stack);
329 fTrack[myIndex]->Status() = (UInt_t)track->GetStatus();
330 fTrack[myIndex]->Length() = (Double_t)track->GetIntegratedLength();
331 fTrack[myIndex]->Charge() = (Int_t)track->Charge();
332 fTrack[myIndex]->PrecX() = (Double_t)track->Px();
333 fTrack[myIndex]->PrecY() = (Double_t)track->Py();
334 fTrack[myIndex]->PrecZ() = (Double_t)track->Pz();
336 // evaluate some flags from the status to decide what to do next in some points
337 isTPC = ((fTrack[myIndex]->Status() & AliESDtrack::kTPCin) != 0);
338 isITSSA = ((fTrack[myIndex]->Status() & AliESDtrack::kTPCin) == 0 && (fTrack[myIndex]->Status() & AliESDtrack::kITSrefit) != 0 && (fTrack[myIndex]->Status() & AliESDtrack::kITSpureSA) == 0 && (fTrack[myIndex]->Status() & AliESDtrack::kITSpid) != 0);
339 isTOF = ((fTrack[myIndex]->Status() & AliESDtrack::kTOFout) != 0 && (fTrack[myIndex]->Status() & AliESDtrack::kTIME) != 0);
341 // accept only tracks which are TPC+ITS or ITS standalone
344 fTrack[myIndex]->ITSsa() = kTRUE;
345 fTrack[myIndex]->TOFok() = kFALSE;
349 fTrack[myIndex]->ITSsa() = kFALSE;
350 fTrack[myIndex]->TOFok() = isTOF;
355 // get DCA to primary vertex
356 track->GetImpactParameters(b, bCov);
357 fTrack[myIndex]->DCAr() = (Double_t)b[0];
358 fTrack[myIndex]->DCAz() = (Double_t)b[1];
361 for (k = 0; k < 6; k++)
363 fTrack[myIndex]->ITSmap(k) = track->HasPointOnITSLayer(k);
367 fTrack[myIndex]->ITSchi2() = track->GetITSchi2();
368 fTrack[myIndex]->ITSsignal() = track->GetITSsignal();
369 nITS = fTrack[myIndex]->SSDcount() + fTrack[myIndex]->SDDcount();
370 for (k = 0; k < AliPID::kSPECIES; k++)
372 fTrack[myIndex]->ITSnsigma(k) = itsrsp.GetNumberOfSigmas(fTrack[myIndex]->Prec(), fTrack[myIndex]->ITSsignal(), (AliPID::EParticleType)k, nITS, kTRUE);
379 fTrack[myIndex]->TPCcount() = (Int_t)track->GetTPCclusters(0);
380 fTrack[myIndex]->TPCchi2() = (Double_t)track->GetTPCchi2();
381 fTrack[myIndex]->TPCsignal() = (Double_t)track->GetTPCsignal();
382 fTrack[myIndex]->PtpcX() = fTrack[myIndex]->PtpcY() = fTrack[myIndex]->PtpcZ() = 1E10;
383 if (track->GetInnerParam())
385 fTrack[myIndex]->PtpcX() = track->GetInnerParam()->Px();
386 fTrack[myIndex]->PtpcY() = track->GetInnerParam()->Py();
387 fTrack[myIndex]->PtpcZ() = track->GetInnerParam()->Pz();
388 for (k = 0; k < AliPID::kSPECIES; k++)
390 fTrack[myIndex]->TPCnsigma(k) = fESDpid->NumberOfSigmasTPC(track, (AliPID::EParticleType)k);
398 track->GetIntegratedTimes(time);
399 fTrack[myIndex]->TOFsignal() = (Double_t)track->GetTOFsignal();
400 for (k = 0; k < AliPID::kSPECIES; k++)
402 fTrack[myIndex]->TOFref(k) = time[k];
403 fTrack[myIndex]->TOFsigma(k) = (Double_t)fTOFmaker->GetExpectedSigma(fTrack[myIndex]->Prec(), time[k], AliPID::ParticleMass(k));