Coverity fixes
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / extra / AliAnalysisTaskSigma1385.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Authors : Massimo Venaruzzo (massimo.venaruzzo@ts.infn.it) *
3 * Enrico Fragiacomo (enrico.fragiacomo@ts.infn.it) *
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
19class TTree;
20class TParticle;
21class TVector3;
22
23#include "AliAnalysisManager.h"
24#include <AliMCEventHandler.h>
25#include <AliMCEvent.h>
26#include <AliStack.h>
27
28class AliESDVertex;
29class AliESDv0;
30class 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
50ClassImp(AliAnalysisTaskSigma1385)
51
52//________________________________________________________________________
53AliAnalysisTaskSigma1385::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
83{
84 // Dummy Constructor
85 Int_t i;
86 for (i = 0; i < 5; i++) fOkTrack[i] = kFALSE;
87}
88
89//________________________________________________________________________
90AliAnalysisTaskSigma1385::AliAnalysisTaskSigma1385(const char *name)
91 : AliAnalysisTaskSE(name), fAnalysisType("ESD"), fCollidingSystems(0), fDataType("REAL"),
92 fListHistCascade(0), fHistEventMultiplicity(0), fHistEventMultiplicityRAVS(0), fNtuple1(0), fNtuple2(0), fNtuple3(0), fNtuple4(0),
93
94 //-------------------------------- For PID
95
96 isMC(0),
97 fIsMC(isMC),
98 fCheckITS(kTRUE),
99 fCheckTPC(kTRUE),
100 fCheckTOF(kTRUE),
101 fUseGlobal(kTRUE),
102 fUseITSSA(kTRUE),
103 fMaxITSband(3.0),
104 fTPCpLimit(0.35),
105 fMinTPCband(3.0),
106 fMaxTPCband(5.0),
107 fESDpid(0x0),
108 fTOFmaker(0x0),
109 fTOFcalib(0x0),
110 fTOFcalibrateESD(!isMC),
111 fTOFcorrectTExp(kTRUE),
112 fTOFuseT0(kTRUE),
113 fTOFtuneMC(isMC),
114 fTOFresolution(100.0),
115 fMinTOF(-2.5),
116 fMaxTOF(3.5),
117 fLastRun(-1)
118
119 //--------------------------------
120{
121
122 // Output slot #0 writes into a TList container (Cascade)
123 DefineOutput(1, TList::Class());
124
125 Int_t i;
126 for (i = 0; i < 5; i++) fOkTrack[i] = kFALSE;
127}
128
129//________________________________________________________________________
130void AliAnalysisTaskSigma1385::UserCreateOutputObjects()
131{
132 fListHistCascade = new TList();
133
134 if (! fHistEventMultiplicity) {
135 fHistEventMultiplicity = new TH1F("fHistEventMultiplicity" , "Nb of Events" , 4, -1.0, 3.0);
136 fListHistCascade->Add(fHistEventMultiplicity);
137 }
138
139 if (! fHistEventMultiplicityRAVS) {
140 fHistEventMultiplicityRAVS = new TH1F("fHistEventMultiplicityRAVS" , "Nb of Events Rejected After Vertex selection" , 4, -1.0, 3.0);
141 fListHistCascade->Add(fHistEventMultiplicityRAVS);
142 }
143
144 if (! fNtuple1) {
145 fNtuple1 = new TNtuple("fNtuple1", "Ntuple1", "TrkNmb");
146 fNtuple1->SetDirectory(0);
147 fListHistCascade->Add(fNtuple1);
148 }
149
150 if (! fNtuple2) {
151 fNtuple2 = new TNtuple("fNtuple2", "Ntuple2", "s:dcal:lCosPoinAn:lDaugDCA:lambdap:lambdapt:lambdamass");
152 fNtuple2->SetDirectory(0);
153 fListHistCascade->Add(fNtuple2);
154 }
155
156 if (! fNtuple3) {
157 fNtuple3 = new TNtuple("fNtuple3", "Ntuple3", "c:dcapi:ppi:ptpi:bachphi:bachtheta:okPiTPC:okPiTOF");
158 fNtuple3->SetDirectory(0);
159 fListHistCascade->Add(fNtuple3);
160 }
161
162 if (! fNtuple4) {
163 fNtuple4 = new TNtuple("fNtuple4", "Ntuple4", "dca:mc:phi:theta:eta:y:pt:p:opang:invmass");
164 fListHistCascade->Add(fNtuple4);
165 }
166
167
168}// end UserCreateOutputObjects
169
170
171//________________________________________________________________________
172void AliAnalysisTaskSigma1385::UserExec(Option_t *)
173{
174
175 // Main loop
176 // Called for each event
177
178
179 Info("AliAnalysisTaskSigma1385", "Starting UserExec");
180
181 AliMCEventHandler* eventHandler;
182 AliMCEvent* mcEvent = 0;
183
184 if (fDataType == "SIM") {
185
186 eventHandler = dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
187 if (!eventHandler) {
188 Printf("ERROR: Could not retrieve MC event handler");
189 return;
190 }
191
192 mcEvent = eventHandler->MCEvent();
193 if (!mcEvent) {
194 Printf("ERROR: Could not retrieve MC event");
195 return;
196 }
197
198 }
199
200 AliStack* stack = 0;
201 if (fDataType == "SIM") {stack = mcEvent->Stack(); fIsMC = 1; isMC = 1;}
202
203 AliESDEvent *lESDevent = 0x0;
204 AliAODEvent *lAODevent = 0x0;
205
206
207 // Connect to the InputEvent
208 Int_t ncascades = -1;
209 if (fAnalysisType == "ESD") {
210 lESDevent = dynamic_cast<AliESDEvent*>(InputEvent());
211 if (lESDevent) ncascades = lESDevent->GetNumberOfCascades();
212 } else if (fAnalysisType == "AOD") {
213 lAODevent = dynamic_cast<AliAODEvent*>(InputEvent());
214 if (lAODevent) ncascades = lAODevent->GetNumberOfCascades();
215 } else {
216 Printf("ERROR: neither lESDevent nor lAODevent are available \n");
217 return;
218 }
219
220
221 //------------------------------for PID
222
223 SetCheckITS(kTRUE);
224 SetCheckTPC(kTRUE);
225 SetCheckTOF(kTRUE);
226
227// ----> set TPC range for PID and calibration
228 SetTPCrange(5.0, 3.0);
229 SetTPCpLimit(0.35);
230
231 // ----> set ITS range for PID
232 SetITSband(4.0);
233
234 // ----> set TPC calibration
235 if (fDataType == "SIM") SetTPCpar(2.15898 / 50.0, 1.75295E1, 3.40030E-9, 1.96178, 3.91720);
236 else SetTPCpar(1.41543 / 50.0, 2.63394E1, 5.0411E-11, 2.12543, 4.88663);
237
238 // ----> set the TOF calibration depending on type of input (sim/data)
239 SetTOFcorrectTExp(kTRUE);
240 SetTOFuseT0(kTRUE);
241 SetTOFresolution(100.0);
242 if (fDataType == "SIM") {
243 SetTOFcalibrateESD(kFALSE);
244 SetTOFtuneMC(kTRUE);
245 } else {
246
247 SetTOFcalibrateESD(kTRUE);
248 SetTOFtuneMC(kFALSE);
249 }
250
251
252 if (!fESDpid) {
253 fESDpid = new AliESDpid;
254 fESDpid->GetTPCResponse().SetBetheBlochParameters(fTPCpar[0], fTPCpar[1], fTPCpar[2], fTPCpar[3], fTPCpar[4]);
255 }
256
257 // initialize DB to current run
258 Int_t run = lESDevent->GetRunNumber();
259 if (run != fLastRun) {
260 //cout << "Run = " << run << " -- LAST = " << fLastRun << endl;
261 fLastRun = run;
262
263 // setup TOF maker & calibration
264 if (!fTOFcalib) fTOFcalib = new AliTOFcalib;
265 fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
266 if (!fTOFmaker) fTOFmaker = new AliTOFT0maker(fESDpid, fTOFcalib);
267 fTOFmaker->SetTimeResolution(fTOFresolution);
268
269 AliCDBManager *cdb = AliCDBManager::Instance();
270 cdb->ClearCache(); // suggestion by Annalisa
271 cdb->Clear(); // suggestion by Annalisa
272 cdb->SetDefaultStorage("raw://");
273 cdb->SetRun(run);
274 fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
275 fTOFcalib->Init();
276 }
277
278 // if required, calibrate the TOF t0 maker with current event
279 if (fTOFcalibrateESD) fTOFcalib->CalibrateESD(lESDevent);
280 if (fTOFtuneMC) fTOFmaker->TuneForMC(lESDevent);
281 if (fTOFuseT0) {
282 fTOFmaker->ComputeT0TOF(lESDevent);
283 fTOFmaker->ApplyT0TOF(lESDevent);
284 fESDpid->MakePID(lESDevent, kFALSE, 0.);
285 }
286
287
288 //--------------------------------------------------------
289
290
291 fHistEventMultiplicity->Fill(1);
292
293 //Some Quantities to characterize the event
294
295 Double_t b = lESDevent->GetMagneticField();
296 Int_t TrackNumber = lESDevent->GetNumberOfTracks();
297
298
299 //---------------------------
300 // new part from AliCascadeVertexer (use vertexer for reconstructing Sigma(1385)
301 //
302 //const AliESDVertex *vtxT3D=lESDevent->GetPrimaryVertex();
303 const AliESDVertex *vtxT3D = lESDevent->GetPrimaryVertexTracks();
304 if (vtxT3D->GetNContributors() < 1) {
305 // SPD vertex
306 vtxT3D = lESDevent->GetPrimaryVertexSPD();
307 if (vtxT3D->GetNContributors() < 1) {
308
309
310 fHistEventMultiplicityRAVS->Fill(1);
311 return;
312 }
313
314
315 }
316
317 Double_t xPrimaryVertex = vtxT3D->GetXv();
318 Double_t yPrimaryVertex = vtxT3D->GetYv();
319 Double_t zPrimaryVertex = vtxT3D->GetZv();
320
321 if (zPrimaryVertex > 10 || zPrimaryVertex < -10) return;
322
323 //-------------------------------------------------------
324 // New Part about tracks global selection criteria
325 //-------------------------------------------------------
326
327 AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
328
329 esdTrackCuts->SetAcceptKinkDaughters(0); // 0 = kFalse
330 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
331 esdTrackCuts->SetMaxChi2PerClusterTPC(4);
332 esdTrackCuts->SetMinNClustersTPC(70);
333
334
335 //-------------------------------------------------------
336 // loops over V0s
337 //stores relevant V0s in an array
338 Int_t nV0 = (Int_t)lESDevent->GetNumberOfV0s();
339
340 TObjArray vtcs(nV0);
341 for (Int_t i = 0; i < nV0; i++) {
342 AliESDv0 *v = lESDevent->GetV0(i);
343 if (v->GetOnFlyStatus()) continue; // if kTRUE, then this V0 is recontructed
344
345 vtcs.AddLast(v);
346 }
347 nV0 = vtcs.GetEntriesFast();
348
349 //-------------------------------------------------------
350 // loops over bachelor tracks
351 // stores relevant tracks in another array
352
353 Int_t nentr = (Int_t)lESDevent->GetNumberOfTracks();
354 TArrayI trk(nentr); Int_t ntr = 0;
355
356 for (Int_t i = 0; i < nentr; i++) {
357 AliESDtrack *esdtr = lESDevent->GetTrack(i);
358
359 if (!(esdtr->GetStatus() & AliESDtrack::kITSrefit)) continue;
360 if (!(esdtr->GetStatus() & AliESDtrack::kTPCrefit)) continue;
361
362
363 if (!(esdTrackCuts->AcceptTrack(esdtr))) continue;
364
365 trk[ntr++] = i;
366 }
367
368
369
370 //-----------------------------------------------------
371 // nested loops over V0s and bachelors
372 Float_t massLambda = 1.11568;
373 Int_t ncasc = 0;
374 AliCascadeVertexer *cascvert = new AliCascadeVertexer();
375 for (Int_t i = 0; i < nV0; i++) { //loop on V0s
376
377 // Lambda
378 AliESDv0 *v = (AliESDv0*)vtcs.UncheckedAt(i);
379 Int_t strness = 0;
380 Double_t lambdaMass = 0;
381 Float_t lambdaP = 0;
382 Float_t lambdaPt = 0;
383 Float_t LambdaDCA = 0; // DCA between Lambda and Primary Vertex
384 Double_t v0cospointangle = 0;
385 Float_t v0daughtersDCA = 0;
386
387
388 // Lambda quality cuts
389 UInt_t lIdxPosXi = (UInt_t) TMath::Abs(v->GetPindex());
390 UInt_t lIdxNegXi = (UInt_t) TMath::Abs(v->GetNindex());
391
392 AliESDtrack *pTrackXi = lESDevent->GetTrack(lIdxPosXi);
393 AliESDtrack *nTrackXi = lESDevent->GetTrack(lIdxNegXi);
394
395 // Filter like-sign V0
396 if (pTrackXi->GetSign() == nTrackXi->GetSign()) continue;
397
398 // WARNING: the following selections cannot be done for AOD yet...
399
400
401 v->ChangeMassHypothesis(kLambda0); // the v0 must be Lambda
402
403
404 if ((TMath::Abs(v->GetEffMass() - massLambda)) < 0.1) {
405
406 if (!(pTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue;
407 if (!(nTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue;
408
409
410 if ((pTrackXi->GetTPCNcls()) < 70) continue;
411 if ((nTrackXi->GetTPCNcls()) < 70) continue;
412
413 strness = 1; lambdaMass = v->GetEffMass(); lambdaP = v->P(); lambdaPt = v->Pt(); LambdaDCA = v->GetD(xPrimaryVertex, yPrimaryVertex, zPrimaryVertex); v0cospointangle = v->GetV0CosineOfPointingAngle(); v0daughtersDCA = v->GetDcaV0Daughters();
414
415
416
417 }
418
419
420 v->ChangeMassHypothesis(kLambda0Bar); // the v0 must be Anti Lambda
421
422
423
424 if ((TMath::Abs(v->GetEffMass() - massLambda)) < 0.1) {
425
426 if (!(pTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue;
427 if (!(nTrackXi->GetStatus() & AliESDtrack::kTPCrefit)) continue;
428
429
430 if ((pTrackXi->GetTPCNcls()) < 70) continue;
431 if ((nTrackXi->GetTPCNcls()) < 70) continue;
432
433 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();
434
435
436 }
437
438 if (strness == 0) continue;
439
440
441 for (Int_t j = 0; j < ntr; j++) { //loop on tracks
442
443 // Pion bachelor
444 Int_t bidx = trk[j];
445 Int_t bachTPCcls = 0;
446
447
448 if ((bidx == v->GetIndex(0)) || (bidx == v->GetIndex(1))) continue; // bach and V0's daughter's must be different!
449
450 AliESDtrack *btrk = lESDevent->GetTrack(bidx);
451
452 if (!(btrk->GetStatus() & AliESDtrack::kTPCrefit)) continue;
453
454 if ((bachTPCcls = btrk->GetTPCNcls()) < 70) continue;
455
456 //Bool_t *IsOkTrack = IsSelected(btrk); //for Alberto's PID
457
458 //Bool_t *okTrack = new Bool_t[5];
459
460 //for (Int_t k = 0; k < 5; k++) okTrack[k] = IsOkTrack[k];
461
462 IsSelected(btrk);
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 fNtuple3->Fill((1.*BachCharge), PionDCA, bachp, bachpt, bachphi, bachtheta, fOkTrack[1], fOkTrack[2]);
576 fNtuple4->Fill(dca, (1.*MCtrue), lPhi, lTheta, lEta, lRapXi, lPt, lP, openingangle, invmass);
577
578 ncasc++;
579
580 //delete IsOkTrack;
581 //delete okTrack;
582
583
584 } // end loop tracks
585 } // end loop V0s
586
587 Info("AliAnalysisTaskSigma1385", "Number of reconstructed Sigma(1385): %d", ncasc);
588
589
590 // Post output data.
591 PostData(1, fListHistCascade);
592
593}
594
595//________________________________________________________________________
596
597//Bool_t *AliAnalysisTaskSigma1385::IsSelected(AliESDtrack *track)
598void AliAnalysisTaskSigma1385::IsSelected(AliESDtrack *track)
599{
600//
601//
602 //Bool_t okTrack[5];
603
604 //for (Int_t i = 0; i < 5; i++) okTrack[i] = kFALSE;
605 for (Int_t i = 0; i < 5; i++) fOkTrack[i] = kFALSE;
606
607
608 AliITSPIDResponse itsrsp(fIsMC);
609
610
611 ULong_t status;
612 Int_t k, nITS;
613 Double_t times[10], tpcNSigma, tpcMaxNSigma, itsSignal, itsNSigma, mom, tofTime, tofSigma, tofRef, tofRel;
614 Bool_t okTOF = kFALSE;
615 Bool_t okTPC = kFALSE;
616 Bool_t okITS = kFALSE;
617 Bool_t isTPC = kFALSE;
618 Bool_t isITSSA = kFALSE;
619 Bool_t isTOF = kFALSE;
620 UChar_t itsCluMap;
621
622 // get commonly used variables
623 status = (ULong_t)track->GetStatus();
624 mom = track->P();
625 isTPC = ((status & AliESDtrack::kTPCin) != 0);
626 isITSSA = (!isTPC && (status & AliESDtrack::kITSrefit) != 0 && (status & AliESDtrack::kITSpureSA) == 0 && (status & AliESDtrack::kITSpid) != 0);
627 isTOF = (((status & AliESDtrack::kTOFout) != 0) && ((status & AliESDtrack::kTIME) != 0) /* && mom > TMath::Max(b1, b2)*/);
628
629
630 // check if the track type matches what is required
631 if (!isTPC && !isITSSA) {
632 AliDebug(AliLog::kDebug + 2, "Track is not either a TPC track or a ITS standalone. Rejected");
633 //return okTrack;
634 return;
635
636 } else if (isTPC && !fUseGlobal) {
637 AliDebug(AliLog::kDebug + 2, "Global tracks not used. Rejected");
638 //return okTrack;
639 return;
640
641 } else if (isITSSA && !fUseITSSA) {
642 AliDebug(AliLog::kDebug + 2, "ITS standalone not used. Rejected");
643 //return okTrack;
644 return;
645 }
646
647 // does a preliminary check on TOF values, if necessary
648 // then, if the reference time or TOF signal are meaningless
649 // even if the 'isTOF' flag is true, switch it to false
650 if (isTOF) {
651 track->GetIntegratedTimes(times);
652 tofTime = (Double_t)track->GetTOFsignal();
653 tofSigma = fTOFmaker->GetExpectedSigma(mom, times[AliPID::kPion], AliPID::ParticleMass(AliPID::kPion));
654 tofRef = times[AliPID::kPion];
655 if (tofRef <= 0.0 && tofSigma <= 0.0) isTOF = kFALSE;
656 }
657
658
659
660 // check TPC dE/dx
661 if (isTPC) { // this branch is entered by all global tracks
662 // check TPC dE/dx:
663 if (fCheckTPC) {
664 tpcNSigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kPion));
665 if (track->GetInnerParam()->P() > fTPCpLimit) tpcMaxNSigma = fMinTPCband; else tpcMaxNSigma = fMaxTPCband;
666 okTPC = (tpcNSigma <= tpcMaxNSigma);
667 AliDebug(AliLog::kDebug + 2, Form("TPC nsigma = %f -- max = %f -- cut %s", tpcNSigma, tpcMaxNSigma, (okTPC ? "passed" : "failed")));
668 } else {
669 // if TPC is not checked, it is as if all tracks do pass the cut
670 okTPC = kTRUE;
671 AliDebug(AliLog::kDebug + 2, "TPC not checked, track accepted");
672 }
673
674 // check TOF (only if flags are OK)
675 if (fCheckTOF) {
676 if (isTOF) {
677 // TOF can be checked only when track is matched there
678 track->GetIntegratedTimes(times);
679 tofTime = (Double_t)track->GetTOFsignal();
680 tofSigma = fTOFmaker->GetExpectedSigma(mom, times[AliPID::kPion], AliPID::ParticleMass(AliPID::kPion));
681 tofRef = times[AliPID::kPion];
682
683 tofRel = (tofTime - tofRef) / tofSigma;
684 okTOF = ((tofRel >= fMinTOF) && (tofRel <= fMaxTOF));
685 AliDebug(AliLog::kDebug + 2, Form("TOF nsigma = %f -- range = %f %f -- cut %s", tofRel, fMinTOF, fMaxTOF, (okTOF ? "passed" : "failed")));
686 } else {
687 // if TOF is not matched, the answer depends on TPC:
688 // - if TPC is required, track is checked only there and TOF simply ignored
689 // - if TPC is not required, track is rejected when TOF does not match it, if TOF check is required
690 if (fCheckTPC) okTOF = kTRUE; else okTOF = kFALSE;
691 }
692 } else {
693 okTOF = kTRUE;
694 }
695
696 } else if (isITSSA) { // this branch is entered by all ITS standalone tracks
697 // check dE/dx only if this is required, otherwise ITS standalone are just added but not checked for PID
698 if (fCheckITS) {
699 itsSignal = track->GetITSsignal();
700 itsCluMap = track->GetITSClusterMap();
701 nITS = 0;
702 for (k = 2; k < 6; k++) if (itsCluMap & (1 << k)) ++nITS;
703 //if (nITS < 3) return kFALSE;
704 //itsNSigma = itsrsp.GetNumberOfSigmas(mom, itsSignal, AliPID::kPion, nITS, kTRUE);
705 //okITS = ((TMath::Abs(itsNSigma)) <= fMaxITSband);
706 if (nITS < 3)
707 okITS = kFALSE;
708 else {
709 itsNSigma = itsrsp.GetNumberOfSigmas(mom, itsSignal, AliPID::kPion, nITS, kTRUE);
710 okITS = ((TMath::Abs(itsNSigma)) <= fMaxITSband);
711 AliDebug(AliLog::kDebug + 2, Form("ITS nsigma = %f -- max = %f -- cut %s", itsNSigma, fMaxITSband, (okITS ? "passed" : "failed")));
712 }
713 } else {
714 okITS = kTRUE;
715 }
716 } else {
717 // if we are here, the track is surely bad
718 okITS = kFALSE;
719 okTPC = kFALSE;
720 okTOF = kFALSE;
721 }
722
723
724 //okTrack[0] = okITS;
725 //okTrack[1] = okTPC;
726 //okTrack[2] = okTOF;
727 //okTrack[3] = isTPC;
728 //okTrack[4] = isITSSA;
729
730 fOkTrack[0] = okITS;
731 fOkTrack[1] = okTPC;
732 fOkTrack[2] = okTOF;
733 fOkTrack[3] = isTPC;
734 fOkTrack[4] = isITSSA;
735
736 //cout<<"##########################################"<<endl;
737 //cout<<"isITSSA "<<isITSSA<< " isTPC "<<isTPC<<endl;
738 //cout<<"##########################################"<<endl;
739
740
741 //return okTrack;
742}
743
744
745
746
747
748
749
750
751
752
753
754//________________________________________________________________________
755void AliAnalysisTaskSigma1385::Terminate(Option_t *)
756{
757 // Draw result to the screen
758 // Called once at the end of the query
759}
760
761
762