]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/RESONANCES/extra/AliAnalysisTaskSigma1385.cxx
some coverity warnings fixed
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / extra / AliAnalysisTaskSigma1385.cxx
CommitLineData
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
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)
7356f978 51
35e49ca5 52//________________________________________________________________________
7356f978 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
35e49ca5 83{
7356f978 84 // Dummy Constructor
35e49ca5 85}
86
87//________________________________________________________________________
7356f978 88AliAnalysisTaskSigma1385::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//________________________________________________________________________
125void 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 167void 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
596Bool_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 738void 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