]>
Commit | Line | Data |
---|---|---|
35e49ca5 | 1 | /************************************************************************** |
2 | * Authors : Massimo Venaruzzo (massimo.venaruzzo@ts.infn.it) * | |
7356f978 | 3 | * Enrico Fragiacomo (enrico.fragiacomo@ts.infn.it) * |
35e49ca5 | 4 | * Contributors are mentioned in the code where appropriate. * |
5 | * * | |
6 | * Permission to use, copy, modify and distribute this software and its * | |
7 | * documentation strictly for non-commercial purposes is hereby granted * | |
8 | * without fee, provided that the above copyright notice appears in all * | |
9 | * copies and that both the copyright notice and this permission notice * | |
10 | * appear in the supporting documentation. The authors make no claims * | |
11 | * about the suitability of this software for any purpose. It is * | |
12 | * provided "as is" without express or implied warranty. * | |
13 | **************************************************************************/ | |
14 | ||
15 | //----------------------------------------------------------------- | |
16 | // AliAnalysisTaskSigma1385 class | |
17 | //----------------------------------------------------------------- | |
18 | ||
19 | class TTree; | |
20 | class TParticle; | |
21 | class TVector3; | |
22 | ||
23 | #include "AliAnalysisManager.h" | |
24 | #include <AliMCEventHandler.h> | |
25 | #include <AliMCEvent.h> | |
26 | #include <AliStack.h> | |
27 | ||
28 | class AliESDVertex; | |
29 | class AliESDv0; | |
30 | class AliAODv0; | |
31 | ||
32 | #include <iostream> | |
33 | ||
34 | #include "TList.h" | |
35 | #include "TH1.h" | |
36 | #include "TNtuple.h" | |
37 | #include "TGraph.h" | |
38 | #include "TCanvas.h" | |
39 | #include "TMath.h" | |
40 | #include "TChain.h" | |
41 | #include "AliLog.h" | |
42 | #include "AliESDEvent.h" | |
43 | #include "AliAODEvent.h" | |
44 | #include "AliCascadeVertexer.h" | |
45 | #include "AliESDcascade.h" | |
46 | #include "AliAODcascade.h" | |
47 | #include "AliAnalysisTaskSigma1385.h" | |
48 | #include "AliESDtrackCuts.h" | |
49 | ||
50 | ClassImp(AliAnalysisTaskSigma1385) | |
7356f978 | 51 | |
35e49ca5 | 52 | //________________________________________________________________________ |
7356f978 | 53 | AliAnalysisTaskSigma1385::AliAnalysisTaskSigma1385() |
54 | : AliAnalysisTaskSE(), fAnalysisType("ESD"), fCollidingSystems(0), fDataType("REAL"), | |
55 | fListHistCascade(0), fHistEventMultiplicity(0), fHistEventMultiplicityRAVS(0), fNtuple1(0), fNtuple2(0), fNtuple3(0), fNtuple4(0), | |
56 | ||
57 | //-------------------------------- For PID | |
58 | isMC(0), | |
59 | fIsMC(isMC), | |
60 | fCheckITS(kTRUE), | |
61 | fCheckTPC(kTRUE), | |
62 | fCheckTOF(kTRUE), | |
63 | fUseGlobal(kTRUE), | |
64 | fUseITSSA(kTRUE), | |
65 | fMaxITSband(3.0), | |
66 | fTPCpLimit(0.35), | |
67 | fMinTPCband(3.0), | |
68 | fMaxTPCband(5.0), | |
69 | fESDpid(0x0), | |
70 | fTOFmaker(0x0), | |
71 | fTOFcalib(0x0), | |
72 | fTOFcalibrateESD(!isMC), | |
73 | fTOFcorrectTExp(kTRUE), | |
74 | fTOFuseT0(kTRUE), | |
75 | fTOFtuneMC(isMC), | |
76 | fTOFresolution(100.0), | |
77 | fMinTOF(-2.5), | |
78 | fMaxTOF(3.5), | |
79 | fLastRun(-1) | |
80 | ||
81 | //-------------------------------- | |
82 | ||
35e49ca5 | 83 | { |
7356f978 | 84 | // Dummy Constructor |
35e49ca5 | 85 | } |
86 | ||
87 | //________________________________________________________________________ | |
7356f978 | 88 | AliAnalysisTaskSigma1385::AliAnalysisTaskSigma1385(const char *name) |
89 | : AliAnalysisTaskSE(name), fAnalysisType("ESD"), fCollidingSystems(0), fDataType("REAL"), | |
90 | fListHistCascade(0), fHistEventMultiplicity(0), fHistEventMultiplicityRAVS(0), fNtuple1(0), fNtuple2(0), fNtuple3(0), fNtuple4(0), | |
91 | ||
92 | //-------------------------------- For PID | |
93 | ||
94 | isMC(0), | |
95 | fIsMC(isMC), | |
96 | fCheckITS(kTRUE), | |
97 | fCheckTPC(kTRUE), | |
98 | fCheckTOF(kTRUE), | |
99 | fUseGlobal(kTRUE), | |
100 | fUseITSSA(kTRUE), | |
101 | fMaxITSband(3.0), | |
102 | fTPCpLimit(0.35), | |
103 | fMinTPCband(3.0), | |
104 | fMaxTPCband(5.0), | |
105 | fESDpid(0x0), | |
106 | fTOFmaker(0x0), | |
107 | fTOFcalib(0x0), | |
108 | fTOFcalibrateESD(!isMC), | |
109 | fTOFcorrectTExp(kTRUE), | |
110 | fTOFuseT0(kTRUE), | |
111 | fTOFtuneMC(isMC), | |
112 | fTOFresolution(100.0), | |
113 | fMinTOF(-2.5), | |
114 | fMaxTOF(3.5), | |
115 | fLastRun(-1) | |
116 | ||
117 | //-------------------------------- | |
35e49ca5 | 118 | { |
119 | ||
7356f978 | 120 | // Output slot #0 writes into a TList container (Cascade) |
121 | DefineOutput(1, TList::Class()); | |
35e49ca5 | 122 | } |
123 | ||
124 | //________________________________________________________________________ | |
125 | void AliAnalysisTaskSigma1385::UserCreateOutputObjects() | |
126 | { | |
7356f978 | 127 | fListHistCascade = new TList(); |
128 | ||
129 | if (! fHistEventMultiplicity) { | |
130 | fHistEventMultiplicity = new TH1F("fHistEventMultiplicity" , "Nb of Events" , 4, -1.0, 3.0); | |
131 | fListHistCascade->Add(fHistEventMultiplicity); | |
132 | } | |
133 | ||
134 | if (! fHistEventMultiplicityRAVS) { | |
135 | fHistEventMultiplicityRAVS = new TH1F("fHistEventMultiplicityRAVS" , "Nb of Events Rejected After Vertex selection" , 4, -1.0, 3.0); | |
136 | fListHistCascade->Add(fHistEventMultiplicityRAVS); | |
137 | } | |
138 | ||
139 | if (! fNtuple1) { | |
140 | fNtuple1 = new TNtuple("fNtuple1", "Ntuple1", "TrkNmb"); | |
141 | fNtuple1->SetDirectory(0); | |
142 | fListHistCascade->Add(fNtuple1); | |
143 | } | |
144 | ||
145 | if (! fNtuple2) { | |
146 | fNtuple2 = new TNtuple("fNtuple2", "Ntuple2", "s:dcal:lCosPoinAn:lDaugDCA:lambdap:lambdapt:lambdamass"); | |
147 | fNtuple2->SetDirectory(0); | |
148 | fListHistCascade->Add(fNtuple2); | |
149 | } | |
150 | ||
151 | if (! fNtuple3) { | |
152 | fNtuple3 = new TNtuple("fNtuple3", "Ntuple3", "c:dcapi:ppi:ptpi:bachphi:bachtheta:okPiTPC:okPiTOF"); | |
153 | fNtuple3->SetDirectory(0); | |
154 | fListHistCascade->Add(fNtuple3); | |
155 | } | |
156 | ||
157 | if (! fNtuple4) { | |
158 | fNtuple4 = new TNtuple("fNtuple4", "Ntuple4", "dca:mc:phi:theta:eta:y:pt:p:opang:invmass"); | |
159 | fListHistCascade->Add(fNtuple4); | |
160 | } | |
161 | ||
162 | ||
35e49ca5 | 163 | }// end UserCreateOutputObjects |
164 | ||
165 | ||
166 | //________________________________________________________________________ | |
7356f978 | 167 | void AliAnalysisTaskSigma1385::UserExec(Option_t *) |
35e49ca5 | 168 | { |
169 | ||
7356f978 | 170 | // Main loop |
171 | // Called for each event | |
172 | ||
173 | ||
174 | Info("AliAnalysisTaskSigma1385", "Starting UserExec"); | |
175 | ||
176 | AliMCEventHandler* eventHandler; | |
177 | AliMCEvent* mcEvent = 0; | |
178 | ||
179 | if (fDataType == "SIM") { | |
180 | ||
181 | eventHandler = dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); | |
182 | if (!eventHandler) { | |
183 | Printf("ERROR: Could not retrieve MC event handler"); | |
184 | return; | |
185 | } | |
186 | ||
187 | mcEvent = eventHandler->MCEvent(); | |
188 | if (!mcEvent) { | |
189 | Printf("ERROR: Could not retrieve MC event"); | |
190 | return; | |
191 | } | |
192 | ||
193 | } | |
194 | ||
195 | AliStack* stack = 0; | |
196 | if (fDataType == "SIM") {stack = mcEvent->Stack(); fIsMC = 1; isMC = 1;} | |
197 | ||
198 | AliESDEvent *lESDevent = 0x0; | |
199 | AliAODEvent *lAODevent = 0x0; | |
200 | ||
201 | ||
202 | // Connect to the InputEvent | |
203 | Int_t ncascades = -1; | |
204 | if (fAnalysisType == "ESD") { | |
205 | lESDevent = dynamic_cast<AliESDEvent*>(InputEvent()); | |
206 | if (!lESDevent) { | |
207 | Printf("ERROR: lESDevent not available \n"); | |
208 | return; | |
209 | } | |
210 | ncascades = lESDevent->GetNumberOfCascades(); | |
211 | } | |
212 | ||
213 | if (fAnalysisType == "AOD") { | |
214 | lAODevent = dynamic_cast<AliAODEvent*>(InputEvent()); | |
215 | if (!lAODevent) { | |
216 | Printf("ERROR: lAODevent not available \n"); | |
217 | return; | |
218 | } | |
219 | ncascades = lAODevent->GetNumberOfCascades(); | |
220 | } | |
221 | ||
222 | ||
223 | //------------------------------for PID | |
224 | ||
225 | SetCheckITS(kTRUE); | |
226 | SetCheckTPC(kTRUE); | |
227 | SetCheckTOF(kTRUE); | |
228 | ||
229 | // ----> set TPC range for PID and calibration | |
230 | SetTPCrange(5.0, 3.0); | |
231 | SetTPCpLimit(0.35); | |
232 | ||
233 | // ----> set ITS range for PID | |
234 | SetITSband(4.0); | |
235 | ||
236 | // ----> set TPC calibration | |
237 | if (fDataType == "SIM") SetTPCpar(2.15898 / 50.0, 1.75295E1, 3.40030E-9, 1.96178, 3.91720); | |
238 | else SetTPCpar(1.41543 / 50.0, 2.63394E1, 5.0411E-11, 2.12543, 4.88663); | |
239 | ||
240 | // ----> set the TOF calibration depending on type of input (sim/data) | |
241 | SetTOFcorrectTExp(kTRUE); | |
242 | SetTOFuseT0(kTRUE); | |
243 | SetTOFresolution(100.0); | |
244 | if (fDataType == "SIM") { | |
245 | SetTOFcalibrateESD(kFALSE); | |
246 | SetTOFtuneMC(kTRUE); | |
247 | } else { | |
248 | ||
249 | SetTOFcalibrateESD(kTRUE); | |
250 | SetTOFtuneMC(kFALSE); | |
251 | } | |
252 | ||
253 | ||
254 | if (!fESDpid) { | |
255 | fESDpid = new AliESDpid; | |
256 | fESDpid->GetTPCResponse().SetBetheBlochParameters(fTPCpar[0], fTPCpar[1], fTPCpar[2], fTPCpar[3], fTPCpar[4]); | |
257 | } | |
258 | ||
259 | // initialize DB to current run | |
260 | Int_t run = lESDevent->GetRunNumber(); | |
261 | if (run != fLastRun) { | |
262 | //cout << "Run = " << run << " -- LAST = " << fLastRun << endl; | |
263 | fLastRun = run; | |
264 | ||
265 | // setup TOF maker & calibration | |
266 | if (!fTOFcalib) fTOFcalib = new AliTOFcalib; | |
267 | fTOFcalib->SetCorrectTExp(fTOFcorrectTExp); | |
268 | if (!fTOFmaker) fTOFmaker = new AliTOFT0maker(fESDpid, fTOFcalib); | |
269 | fTOFmaker->SetTimeResolution(fTOFresolution); | |
270 | ||
271 | AliCDBManager *cdb = AliCDBManager::Instance(); | |
272 | cdb->ClearCache(); // suggestion by Annalisa | |
273 | cdb->Clear(); // suggestion by Annalisa | |
274 | cdb->SetDefaultStorage("raw://"); | |
275 | cdb->SetRun(run); | |
276 | fTOFcalib->SetCorrectTExp(fTOFcorrectTExp); | |
277 | fTOFcalib->Init(); | |
278 | } | |
279 | ||
280 | // if required, calibrate the TOF t0 maker with current event | |
281 | if (fTOFcalibrateESD) fTOFcalib->CalibrateESD(lESDevent); | |
282 | if (fTOFtuneMC) fTOFmaker->TuneForMC(lESDevent); | |
283 | if (fTOFuseT0) { | |
284 | fTOFmaker->ComputeT0TOF(lESDevent); | |
285 | fTOFmaker->ApplyT0TOF(lESDevent); | |
286 | fESDpid->MakePID(lESDevent, kFALSE, 0.); | |
287 | } | |
288 | ||
289 | ||
290 | //-------------------------------------------------------- | |
291 | ||
292 | ||
293 | fHistEventMultiplicity->Fill(1); | |
294 | ||
295 | //Some Quantities to characterize the event | |
35e49ca5 | 296 | |
7356f978 | 297 | Double_t b = lESDevent->GetMagneticField(); |
298 | Int_t TrackNumber = lESDevent->GetNumberOfTracks(); | |
299 | ||
300 | ||
301 | //--------------------------- | |
302 | // new part from AliCascadeVertexer (use vertexer for reconstructing Sigma(1385) | |
303 | // | |
304 | //const AliESDVertex *vtxT3D=lESDevent->GetPrimaryVertex(); | |
305 | const AliESDVertex *vtxT3D = lESDevent->GetPrimaryVertexTracks(); | |
306 | if (vtxT3D->GetNContributors() < 1) { | |
307 | // SPD vertex | |
308 | vtxT3D = lESDevent->GetPrimaryVertexSPD(); | |
309 | if (vtxT3D->GetNContributors() < 1) { | |
310 | ||
311 | ||
312 | fHistEventMultiplicityRAVS->Fill(1); | |
313 | return; | |
35e49ca5 | 314 | } |
7356f978 | 315 | |
316 | ||
317 | } | |
318 | ||
319 | Double_t xPrimaryVertex = vtxT3D->GetXv(); | |
320 | Double_t yPrimaryVertex = vtxT3D->GetYv(); | |
321 | Double_t zPrimaryVertex = vtxT3D->GetZv(); | |
322 | ||
323 | if (zPrimaryVertex > 10 || zPrimaryVertex < -10) return; | |
324 | ||
325 | //------------------------------------------------------- | |
326 | // New Part about tracks global selection criteria | |
327 | //------------------------------------------------------- | |
328 | ||
329 | AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts"); | |
330 | ||
331 | esdTrackCuts->SetAcceptKinkDaughters(0); // 0 = kFalse | |
332 | esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny); | |
333 | esdTrackCuts->SetMaxChi2PerClusterTPC(4); | |
334 | esdTrackCuts->SetMinNClustersTPC(70); | |
335 | ||
336 | ||
337 | //------------------------------------------------------- | |
338 | // loops over V0s | |
339 | //stores relevant V0s in an array | |
340 | Int_t nV0 = (Int_t)lESDevent->GetNumberOfV0s(); | |
341 | ||
342 | TObjArray vtcs(nV0); | |
343 | for (Int_t i = 0; i < nV0; i++) { | |
344 | AliESDv0 *v = lESDevent->GetV0(i); | |
345 | if (v->GetOnFlyStatus()) continue; // if kTRUE, then this V0 is recontructed | |
346 | ||
347 | vtcs.AddLast(v); | |
348 | } | |
349 | nV0 = vtcs.GetEntriesFast(); | |
350 | ||
351 | //------------------------------------------------------- | |
352 | // loops over bachelor tracks | |
353 | // stores relevant tracks in another array | |
354 | ||
355 | Int_t nentr = (Int_t)lESDevent->GetNumberOfTracks(); | |
356 | TArrayI trk(nentr); Int_t ntr = 0; | |
357 | ||
358 | for (Int_t i = 0; i < nentr; i++) { | |
359 | AliESDtrack *esdtr = lESDevent->GetTrack(i); | |
360 | ||
361 | if (!(esdtr->GetStatus() & AliESDtrack::kITSrefit)) continue; | |
362 | if (!(esdtr->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
363 | ||
364 | ||
365 | if (!(esdTrackCuts->AcceptTrack(esdtr))) continue; | |
366 | ||
367 | trk[ntr++] = i; | |
368 | } | |
369 | ||
370 | ||
371 | ||
372 | //----------------------------------------------------- | |
373 | // nested loops over V0s and bachelors | |
374 | Float_t massLambda = 1.11568; | |
375 | Int_t ncasc = 0; | |
376 | AliCascadeVertexer *cascvert = new AliCascadeVertexer(); | |
377 | for (Int_t i = 0; i < nV0; i++) { //loop on V0s | |
378 | ||
379 | // Lambda | |
380 | AliESDv0 *v = (AliESDv0*)vtcs.UncheckedAt(i); | |
381 | Int_t strness = 0; | |
382 | Double_t lambdaMass = 0; | |
383 | Float_t lambdaP = 0; | |
384 | Float_t lambdaPt = 0; | |
385 | Float_t LambdaDCA = 0; // DCA between Lambda and Primary Vertex | |
386 | Double_t v0cospointangle = 0; | |
387 | Float_t v0daughtersDCA = 0; | |
388 | ||
389 | ||
390 | // Lambda quality cuts | |
391 | UInt_t lIdxPosXi = (UInt_t) TMath::Abs(v->GetPindex()); | |
392 | UInt_t lIdxNegXi = (UInt_t) TMath::Abs(v->GetNindex()); | |
393 | ||
394 | AliESDtrack *pTrackXi = lESDevent->GetTrack(lIdxPosXi); | |
395 | AliESDtrack *nTrackXi = lESDevent->GetTrack(lIdxNegXi); | |
396 | ||
397 | // Filter like-sign V0 | |
398 | if (pTrackXi->GetSign() == nTrackXi->GetSign()) continue; | |
399 | ||
400 | // WARNING: the following selections cannot be done for AOD yet... | |
401 | ||
402 | ||
403 | v->ChangeMassHypothesis(kLambda0); // the v0 must be Lambda | |
404 | ||
405 | ||
406 | if ((TMath::Abs(v->GetEffMass() - massLambda)) < 0.1) { | |
407 | ||
408 | if (!(pTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
409 | if (!(nTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
410 | ||
411 | ||
412 | if ((pTrackXi->GetTPCNcls()) < 70) continue; | |
413 | if ((nTrackXi->GetTPCNcls()) < 70) continue; | |
414 | ||
415 | strness = 1; lambdaMass = v->GetEffMass(); lambdaP = v->P(); lambdaPt = v->Pt(); LambdaDCA = v->GetD(xPrimaryVertex, yPrimaryVertex, zPrimaryVertex); v0cospointangle = v->GetV0CosineOfPointingAngle(); v0daughtersDCA = v->GetDcaV0Daughters(); | |
416 | ||
417 | ||
418 | ||
419 | } | |
420 | ||
421 | ||
422 | v->ChangeMassHypothesis(kLambda0Bar); // the v0 must be Anti Lambda | |
423 | ||
424 | ||
425 | ||
426 | if ((TMath::Abs(v->GetEffMass() - massLambda)) < 0.1) { | |
427 | ||
428 | if (!(pTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
429 | if (!(nTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
430 | ||
431 | ||
432 | if ((pTrackXi->GetTPCNcls()) < 70) continue; | |
433 | if ((nTrackXi->GetTPCNcls()) < 70) continue; | |
434 | ||
435 | Int_t temp = strness + 1; strness = -1 * temp; lambdaMass = v->GetEffMass(); lambdaP = v->P(); lambdaPt = v->Pt(); LambdaDCA = v->GetD(xPrimaryVertex, yPrimaryVertex, zPrimaryVertex); v0cospointangle = v->GetV0CosineOfPointingAngle(); v0daughtersDCA = v->GetDcaV0Daughters(); | |
436 | ||
437 | ||
438 | } | |
439 | ||
440 | if (strness == 0) continue; | |
441 | ||
442 | ||
443 | for (Int_t j = 0; j < ntr; j++) { //loop on tracks | |
444 | ||
445 | // Pion bachelor | |
446 | Int_t bidx = trk[j]; | |
447 | Int_t bachTPCcls = 0; | |
448 | ||
449 | ||
450 | if ((bidx == v->GetIndex(0)) || (bidx == v->GetIndex(1))) continue; // bach and V0's daughter's must be different! | |
451 | ||
452 | AliESDtrack *btrk = lESDevent->GetTrack(bidx); | |
453 | ||
454 | if (!(btrk->GetStatus() & AliESDtrack::kTPCrefit)) continue; | |
455 | ||
456 | if ((bachTPCcls = btrk->GetTPCNcls()) < 70) continue; | |
457 | ||
458 | Bool_t *IsOkTrack = IsSelected(btrk); //for Alberto's PID | |
459 | ||
460 | Bool_t *okTrack = new Bool_t[5]; | |
461 | ||
462 | for (Int_t k = 0; k < 5; k++) okTrack[k] = IsOkTrack[k]; | |
463 | ||
464 | ||
465 | Int_t BachCharge = btrk->Charge(); | |
466 | Float_t PionDCA = TMath::Abs(btrk->GetD(xPrimaryVertex, yPrimaryVertex, b)); // DCA between Bachelor and Primary Vertex | |
467 | ||
468 | Double_t bachphi = btrk->Phi(); | |
469 | Double_t bachtheta = btrk->Theta(); | |
470 | Double_t bachmass = 0.13957; | |
471 | Double_t bachp = btrk->GetP(); | |
472 | Double_t bachpt = btrk->Pt(); | |
473 | ||
474 | ||
475 | // Distance between Lambda and Pion bachelor | |
476 | AliESDv0 v0(*v), *pv0 = &v0; | |
477 | AliExternalTrackParam bt(*btrk), *pbt = &bt; | |
478 | Double_t dca = cascvert->PropagateToDCA(pv0, pbt, b); // distance bwn V0 and bach (in cm) | |
479 | if (dca > 10.0) continue; // Note: was 0.1! Enlarged->further filter in second pass analysis | |
480 | ||
481 | ||
482 | AliESDcascade cascade(*pv0, *pbt, bidx); //constucts a sigma1385 candidate | |
483 | AliESDcascade *xi = &cascade; | |
484 | ||
485 | ||
486 | UInt_t lBachIdx = (UInt_t) TMath::Abs(xi->GetBindex()); | |
487 | AliESDtrack *bachTrackXi = lESDevent->GetTrack(lBachIdx); | |
488 | ||
489 | ||
490 | ||
491 | Bool_t MCtrue = 0; | |
492 | if (fDataType == "SIM") { | |
493 | ||
494 | Int_t pLabel = TMath::Abs(pTrackXi->GetLabel()); | |
495 | Int_t nLabel = TMath::Abs(nTrackXi->GetLabel()); | |
496 | Int_t bachLabel = TMath::Abs(bachTrackXi->GetLabel()); | |
497 | ||
498 | Int_t motherpLabel = stack->Particle(pLabel)->GetFirstMother(); | |
499 | Int_t mothernLabel = stack->Particle(nLabel)->GetFirstMother(); | |
500 | Int_t motherbachLabel = stack->Particle(bachLabel)->GetFirstMother(); | |
501 | ||
502 | ||
503 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
504 | //cout<< "plabel " << pLabel << " nlabel " << nLabel << " mother p " << motherpLabel << " mother n " << mothernLabel << " bachlabel " << bachLabel << " mother bach " << motherbachLabel << endl; | |
505 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
506 | ||
507 | if (motherpLabel > -1 && mothernLabel > -1) { | |
508 | TParticle *plambda = stack->Particle(motherpLabel); | |
509 | Int_t grandmother = plambda->GetFirstMother(); | |
510 | ||
511 | motherpLabel = TMath::Abs(motherpLabel); | |
512 | mothernLabel = TMath::Abs(mothernLabel); | |
513 | //motherbachLabel = TMath::Abs(motherbachLabel); | |
514 | ||
515 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
516 | //cout<< "plabel " << pLabel << " nlabel " << nLabel << " mother p " << motherpLabel << " mother n " << mothernLabel << " mother lambda " << grandmother << " bachlabel " << bachLabel << " mother bach " << motherbachLabel << endl; | |
517 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
518 | ||
519 | if (motherbachLabel > -1) { | |
520 | ||
521 | ||
522 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
523 | //cout<< "mother lambda " << grandmother << " mother bach " << motherbachLabel << " PDG Code "<< stack->Particle(grandmother)->GetPdgCode() << endl; | |
524 | //cout<< "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<endl; | |
525 | ||
526 | if ((motherpLabel == mothernLabel) && (grandmother == motherbachLabel) && | |
527 | ((TMath::Abs(stack->Particle(grandmother)->GetPdgCode()) == 3114) || | |
528 | (TMath::Abs(stack->Particle(grandmother)->GetPdgCode()) == 3224))) MCtrue = 1; | |
529 | ||
530 | } | |
531 | ||
532 | } | |
533 | ||
534 | } | |
535 | ||
536 | Double_t lBachMomX = 0., lBachMomY = 0., lBachMomZ = 0.; | |
537 | Double_t lPMom[3] = {0., 0., 0.,}; | |
538 | Double_t lNMom[3] = {0., 0., 0.,}; | |
539 | Float_t lLambdaMomX = 0., lLambdaMomY = 0., lLambdaMomZ = 0.; | |
540 | xi->GetBPxPyPz(lBachMomX, lBachMomY, lBachMomZ); | |
541 | xi->GetPPxPyPz(lPMom[0], lPMom[1], lPMom[2]); | |
542 | xi->GetNPxPyPz(lNMom[0], lNMom[1], lNMom[2]); | |
543 | lLambdaMomX = lPMom[0] + lNMom[0]; | |
544 | lLambdaMomY = lPMom[1] + lNMom[1]; | |
545 | lLambdaMomZ = lPMom[2] + lNMom[2]; | |
546 | ||
547 | Float_t lRapXi = xi->RapXi(); | |
548 | Float_t lEta = xi->Eta(); | |
549 | Float_t lTheta = xi->Theta(); | |
550 | Float_t lPhi = xi->Phi(); | |
551 | Float_t lPt = xi->Pt(); | |
552 | Float_t lP = xi->P(); | |
553 | ||
554 | // Support variables for invariant mass calculation | |
555 | TLorentzVector Lambda, Pion, Sigma; | |
556 | Double_t Elambda = TMath::Sqrt(lLambdaMomX * lLambdaMomX + lLambdaMomY * lLambdaMomY + lLambdaMomZ * lLambdaMomZ + lambdaMass * lambdaMass); | |
557 | Double_t Epion = TMath::Sqrt(lBachMomX * lBachMomX + lBachMomY * lBachMomY + lBachMomZ * lBachMomZ + bachmass * bachmass) ; | |
558 | ||
559 | Lambda.SetPxPyPzE(lLambdaMomX, lLambdaMomY, lLambdaMomZ, Elambda); | |
560 | Pion.SetPxPyPzE(lBachMomX, lBachMomY, lBachMomZ, Epion); | |
561 | ||
562 | Sigma = Lambda + Pion; | |
563 | ||
564 | Double_t openingangle, invmass = 0; | |
565 | ||
566 | openingangle = TMath::ACos((lBachMomX * lLambdaMomX + lBachMomY * lLambdaMomY + lBachMomZ * lLambdaMomZ) / ((TMath::Sqrt(lBachMomX * lBachMomX + lBachMomY * lBachMomY + lBachMomZ * lBachMomZ)) * (TMath::Sqrt(lLambdaMomX * lLambdaMomX + lLambdaMomY * lLambdaMomY + | |
567 | lLambdaMomZ * lLambdaMomZ)))); | |
568 | ||
569 | invmass = Sigma.M(); | |
570 | ||
571 | ||
572 | fNtuple1->Fill(TrackNumber); | |
573 | fNtuple2->Fill((1.*strness), LambdaDCA, v0cospointangle, v0daughtersDCA , lambdaP, lambdaPt, lambdaMass); | |
574 | fNtuple3->Fill((1.*BachCharge), PionDCA, bachp, bachpt, bachphi, bachtheta, okTrack[1], okTrack[2]); | |
575 | fNtuple4->Fill(dca, (1.*MCtrue), lPhi, lTheta, lEta, lRapXi, lPt, lP, openingangle, invmass); | |
576 | ||
577 | ncasc++; | |
578 | ||
579 | delete IsOkTrack; | |
580 | delete okTrack; | |
581 | ||
582 | ||
583 | } // end loop tracks | |
584 | } // end loop V0s | |
585 | ||
586 | Info("AliAnalysisTaskSigma1385", "Number of reconstructed Sigma(1385): %d", ncasc); | |
587 | ||
588 | ||
589 | // Post output data. | |
590 | PostData(1, fListHistCascade); | |
591 | ||
35e49ca5 | 592 | } |
593 | ||
594 | //________________________________________________________________________ | |
595 | ||
596 | Bool_t *AliAnalysisTaskSigma1385::IsSelected(AliESDtrack *track) | |
597 | { | |
598 | // | |
599 | // | |
7356f978 | 600 | Bool_t *okTrack = new Bool_t[5]; |
601 | ||
602 | for (Int_t i = 0; i < 5; i++) okTrack[i] = kFALSE; | |
603 | ||
604 | ||
605 | AliITSPIDResponse itsrsp(fIsMC); | |
606 | ||
607 | ||
608 | ULong_t status; | |
609 | Int_t k, nITS; | |
610 | Double_t times[10], tpcNSigma, tpcMaxNSigma, itsSignal, itsNSigma, mom, tofTime, tofSigma, tofRef, tofRel; | |
611 | Bool_t okTOF = kFALSE; | |
612 | Bool_t okTPC = kFALSE; | |
613 | Bool_t okITS = kFALSE; | |
614 | Bool_t isTPC = kFALSE; | |
615 | Bool_t isITSSA = kFALSE; | |
616 | Bool_t isTOF = kFALSE; | |
617 | UChar_t itsCluMap; | |
618 | ||
619 | // get commonly used variables | |
620 | status = (ULong_t)track->GetStatus(); | |
621 | mom = track->P(); | |
622 | isTPC = ((status & AliESDtrack::kTPCin) != 0); | |
623 | isITSSA = ((status & AliESDtrack::kTPCin) == 0 && (status & AliESDtrack::kITSrefit) != 0 && (status & AliESDtrack::kITSpureSA) == 0 && (status & AliESDtrack::kITSpid) != 0); | |
624 | isTOF = (((status & AliESDtrack::kTOFout) != 0) && ((status & AliESDtrack::kTIME) != 0) /* && mom > TMath::Max(b1, b2)*/); | |
625 | ||
626 | ||
627 | // check if the track type matches what is required | |
628 | if (!isTPC && !isITSSA) { | |
629 | AliDebug(AliLog::kDebug + 2, "Track is not either a TPC track or a ITS standalone. Rejected"); | |
630 | return okTrack; | |
631 | ||
632 | } else if (isTPC && !fUseGlobal) { | |
633 | AliDebug(AliLog::kDebug + 2, "Global tracks not used. Rejected"); | |
634 | return okTrack; | |
635 | ||
636 | } else if (isITSSA && !fUseITSSA) { | |
637 | AliDebug(AliLog::kDebug + 2, "ITS standalone not used. Rejected"); | |
638 | return okTrack; | |
639 | ||
640 | } | |
641 | ||
642 | // does a preliminary check on TOF values, if necessary | |
643 | // then, if the reference time or TOF signal are meaningless | |
644 | // even if the 'isTOF' flag is true, switch it to false | |
645 | if (isTOF) { | |
646 | track->GetIntegratedTimes(times); | |
647 | tofTime = (Double_t)track->GetTOFsignal(); | |
648 | tofSigma = fTOFmaker->GetExpectedSigma(mom, times[AliPID::kPion], AliPID::ParticleMass(AliPID::kPion)); | |
649 | tofRef = times[AliPID::kPion]; | |
650 | if (tofRef <= 0.0 && tofSigma <= 0.0) isTOF = kFALSE; | |
651 | } | |
652 | ||
653 | ||
654 | ||
655 | // check TPC dE/dx | |
656 | if (isTPC) { // this branch is entered by all global tracks | |
657 | // check TPC dE/dx: | |
658 | if (fCheckTPC) { | |
659 | tpcNSigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kPion)); | |
660 | if (track->GetInnerParam()->P() > fTPCpLimit) tpcMaxNSigma = fMinTPCband; else tpcMaxNSigma = fMaxTPCband; | |
661 | okTPC = (tpcNSigma <= tpcMaxNSigma); | |
662 | AliDebug(AliLog::kDebug + 2, Form("TPC nsigma = %f -- max = %f -- cut %s", tpcNSigma, tpcMaxNSigma, (okTPC ? "passed" : "failed"))); | |
663 | } else { | |
664 | // if TPC is not checked, it is as if all tracks do pass the cut | |
665 | okTPC = kTRUE; | |
666 | AliDebug(AliLog::kDebug + 2, "TPC not checked, track accepted"); | |
667 | } | |
668 | ||
669 | // check TOF (only if flags are OK) | |
670 | if (fCheckTOF) { | |
671 | if (isTOF) { | |
672 | // TOF can be checked only when track is matched there | |
673 | track->GetIntegratedTimes(times); | |
674 | tofTime = (Double_t)track->GetTOFsignal(); | |
675 | tofSigma = fTOFmaker->GetExpectedSigma(mom, times[AliPID::kPion], AliPID::ParticleMass(AliPID::kPion)); | |
676 | tofRef = times[AliPID::kPion]; | |
677 | ||
678 | tofRel = (tofTime - tofRef) / tofSigma; | |
679 | okTOF = ((tofRel >= fMinTOF) && (tofRel <= fMaxTOF)); | |
680 | AliDebug(AliLog::kDebug + 2, Form("TOF nsigma = %f -- range = %f %f -- cut %s", tofRel, fMinTOF, fMaxTOF, (okTOF ? "passed" : "failed"))); | |
681 | } else { | |
682 | // if TOF is not matched, the answer depends on TPC: | |
683 | // - if TPC is required, track is checked only there and TOF simply ignored | |
684 | // - if TPC is not required, track is rejected when TOF does not match it, if TOF check is required | |
685 | if (fCheckTPC) okTOF = kTRUE; else okTOF = kFALSE; | |
686 | } | |
687 | } else { | |
688 | okTOF = kTRUE; | |
35e49ca5 | 689 | } |
7356f978 | 690 | |
691 | } else if (isITSSA) { // this branch is entered by all ITS standalone tracks | |
692 | // check dE/dx only if this is required, otherwise ITS standalone are just added but not checked for PID | |
693 | if (fCheckITS) { | |
694 | itsSignal = track->GetITSsignal(); | |
695 | itsCluMap = track->GetITSClusterMap(); | |
696 | nITS = 0; | |
697 | for (k = 2; k < 6; k++) if (itsCluMap & (1 << k)) ++nITS; | |
698 | if (nITS < 3) return kFALSE; | |
699 | itsNSigma = itsrsp.GetNumberOfSigmas(mom, itsSignal, AliPID::kPion, nITS, kTRUE); | |
700 | okITS = ((TMath::Abs(itsNSigma)) <= fMaxITSband); | |
701 | AliDebug(AliLog::kDebug + 2, Form("ITS nsigma = %f -- max = %f -- cut %s", itsNSigma, fMaxITSband, (okITS ? "passed" : "failed"))); | |
702 | } else { | |
703 | okITS = kTRUE; | |
35e49ca5 | 704 | } |
7356f978 | 705 | } else { |
706 | // if we are here, the track is surely bad | |
707 | okITS = kFALSE; | |
708 | okTPC = kFALSE; | |
709 | okTOF = kFALSE; | |
35e49ca5 | 710 | } |
7356f978 | 711 | |
712 | ||
713 | okTrack[0] = okITS; | |
714 | okTrack[1] = okTPC; | |
715 | okTrack[2] = okTOF; | |
716 | okTrack[3] = isTPC; | |
717 | okTrack[4] = isITSSA; | |
718 | ||
719 | //cout<<"##########################################"<<endl; | |
720 | //cout<<"isITSSA "<<isITSSA<< " isTPC "<<isTPC<<endl; | |
721 | //cout<<"##########################################"<<endl; | |
722 | ||
723 | ||
724 | return okTrack; | |
35e49ca5 | 725 | } |
726 | ||
727 | ||
728 | ||
729 | ||
730 | ||
731 | ||
732 | ||
733 | ||
734 | ||
735 | ||
736 | ||
737 | //________________________________________________________________________ | |
7356f978 | 738 | void AliAnalysisTaskSigma1385::Terminate(Option_t *) |
35e49ca5 | 739 | { |
7356f978 | 740 | // Draw result to the screen |
741 | // Called once at the end of the query | |
35e49ca5 | 742 | } |
743 | ||
744 | ||
745 |