1 /**************************************************************************
2 * Authors : Domenico Colella *
6 * - Original AliAnalysisTaskCheckCascade (A. Maire, G. Hippolyte) *
7 * - Adapted to PbPb analysis (M. Nicassio) *
8 * - Adapted to work on all collisidng systems (pp, PbPb and pPb), *
9 * ESD/AOD and experimental/MC data (D. Colella) *
11 * Permission to use, copy, modify and distribute this software and its *
12 * documentation strictly for non-commercial purposes is hereby granted *
13 * without fee, provided that the above copyright notice appears in all *
14 * copies and that both the copyright notice and this permission notice *
15 * appear in the supporting documentation. The authors make no claims *
16 * about the suitability of this software for any purpose. It is *
17 * provided "as is" without express or implied warranty. *
18 **************************************************************************/
30 #include <Riostream.h>
31 #include "THnSparse.h"
36 #include "AliCentrality.h"
37 #include "AliHeader.h" //for MC
38 #include "AliMCEvent.h" //for MC
39 #include "AliStack.h" //for MC
40 #include "AliESDEvent.h"
41 #include "AliAODEvent.h"
42 #include "AliESDtrackCuts.h"
43 #include "AliPIDResponse.h"
45 #include "AliInputEventHandler.h"
46 #include "AliAnalysisManager.h"
47 #include "AliESDInputHandler.h"
48 #include "AliAODInputHandler.h"
49 #include "AliCFContainer.h"
50 #include "AliMultiplicity.h"
52 #include "AliGenEventHeader.h" //for MC
53 #include "AliGenCocktailEventHeader.h" //for MC
54 #include "AliGenHijingEventHeader.h" //for MC
55 #include "AliAODMCParticle.h" //for MC
57 #include "AliESDcascade.h"
58 #include "AliAODcascade.h"
59 #include "AliESDtrack.h"
60 #include "AliAODTrack.h"
61 #include "AliESDpid.h"
66 #include "AliAnalysisTaskQAMultistrange.h"
68 ClassImp(AliAnalysisTaskQAMultistrange)
72 //________________________________________________________________________
73 AliAnalysisTaskQAMultistrange::AliAnalysisTaskQAMultistrange()
74 : AliAnalysisTaskSE (),
76 fAnalysisType ("ESD"),
77 fCollidingSystem ("PbPb"),
79 fkQualityCutZprimVtxPos (kTRUE),
80 fkQualityCutNoTPConlyPrimVtx(kTRUE),
81 fkQualityCutTPCrefit (kTRUE),
82 fkQualityCutnTPCcls (kTRUE),
83 fkQualityCutPileup (kTRUE),
90 fMinPtCutOnDaughterTracks (0),
91 fEtaCutOnDaughterTracks (0),
94 fListHistMultistrangeQA(0),
96 fHistMassXiMinus(0), fHistMassXiPlus(0), fHistMassOmegaMinus(0), fHistMassOmegaPlus(0),
97 fCFContCascadeCuts(0),
98 fCFContCascadeMCgen(0)
105 //________________________________________________________________________
106 AliAnalysisTaskQAMultistrange::AliAnalysisTaskQAMultistrange(const char *name)
107 : AliAnalysisTaskSE (name),
109 fAnalysisType ("ESD"),
110 fCollidingSystem ("PbPb"),
112 fkQualityCutZprimVtxPos (kTRUE),
113 fkQualityCutNoTPConlyPrimVtx(kTRUE),
114 fkQualityCutTPCrefit (kTRUE),
115 fkQualityCutnTPCcls (kTRUE),
116 fkQualityCutPileup (kTRUE),
123 fMinPtCutOnDaughterTracks (0),
124 fEtaCutOnDaughterTracks (0),
127 fListHistMultistrangeQA(0),
129 fHistMassXiMinus(0), fHistMassXiPlus(0), fHistMassOmegaMinus(0), fHistMassOmegaPlus(0),
130 fCFContCascadeCuts(0),
131 fCFContCascadeMCgen(0)
135 // Output slot #0 writes into a TList container (Cascade)
136 DefineOutput(1, TList::Class());
137 DefineOutput(2, AliCFContainer::Class());
138 DefineOutput(3, AliCFContainer::Class());
140 AliLog::SetClassDebugLevel("AliAnalysisTaskQAMultistrange",1);
144 AliAnalysisTaskQAMultistrange::~AliAnalysisTaskQAMultistrange() {
148 // For all TH1, 2, 3 HnSparse and CFContainer are in the fListCascade TList.
149 // They will be deleted when fListCascade is deleted by the TSelector dtor
150 // Because of TList::SetOwner() ...
151 if (fListHistMultistrangeQA && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fListHistMultistrangeQA; fListHistMultistrangeQA = 0x0; }
152 if (fCFContCascadeCuts && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fCFContCascadeCuts; fCFContCascadeCuts = 0x0; }
153 if (fCFContCascadeMCgen && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fCFContCascadeMCgen; fCFContCascadeMCgen = 0x0; }
158 //_____________________________________________________________
159 void AliAnalysisTaskQAMultistrange::UserCreateOutputObjects() {
163 //----------------------------------------------------------------------------
164 // Option for AliLog: to suppress the extensive info prompted by a run with MC
165 //----------------------------------------------------------------------------
166 if (fisMC) AliLog::SetGlobalLogLevel(AliLog::kError);
168 //-----------------------------------------------
169 // Particle Identification Setup (new PID object)
170 //-----------------------------------------------
171 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
172 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
173 fPIDResponse = inputHandler->GetPIDResponse();
178 fListHistMultistrangeQA = new TList();
179 fListHistMultistrangeQA->SetOwner();
184 if(! fHistEventSel) {
185 fHistEventSel = new TH1F("fHistEventSel", "Event selection;Evt. Sel. Step;Count",6, 0, 6);
186 fHistEventSel->GetXaxis()->SetBinLabel(1, "Processed");
187 fHistEventSel->GetXaxis()->SetBinLabel(2, "PhysEvSel");
188 fHistEventSel->GetXaxis()->SetBinLabel(3, "Multiplicity");
189 fHistEventSel->GetXaxis()->SetBinLabel(4, "GoodPrVtx");
190 fHistEventSel->GetXaxis()->SetBinLabel(5, "PrVtxPosition");
191 fHistEventSel->GetXaxis()->SetBinLabel(6, "PileUpSel");
192 fListHistMultistrangeQA->Add(fHistEventSel);
194 if(! fHistMassXiMinus) {
195 fHistMassXiMinus = new TH1F("fHistMassXiMinus", "#Xi^{-} candidates;M(#Lambda,#pi^{-}) (GeV/c^{2}); Counts", 150, 1.25, 1.40);
196 fListHistMultistrangeQA->Add(fHistMassXiMinus);
198 if(! fHistMassXiPlus) {
199 fHistMassXiPlus = new TH1F("fHistMassXiPlus", "#Xi^{+} candidates; M(#bar{#Lambda}^{0},#pi^{+}) (GeV/c^{2}); Counts", 150, 1.25, 1.40);
200 fListHistMultistrangeQA->Add(fHistMassXiPlus);
202 if(! fHistMassOmegaMinus) {
203 fHistMassOmegaMinus = new TH1F("fHistMassOmegaMinus", "#Omega^{-} candidates; M(#Lambda,K^{-}) (GeV/c^{2}); Counts", 120, 1.62, 1.74);
204 fListHistMultistrangeQA->Add(fHistMassOmegaMinus);
206 if(! fHistMassOmegaPlus) {
207 fHistMassOmegaPlus = new TH1F("fHistMassOmegaPlus", "#Omega^{+} candidates;M(#bar{#Lambda}^{0},K^{+}) (GeV/c^{2}); Counts", 120, 1.62, 1.74);
208 fListHistMultistrangeQA->Add(fHistMassOmegaPlus);
212 //---------------------------------------------------
213 // Define the container for the topological variables
214 //---------------------------------------------------
215 if(! fCFContCascadeCuts) {
216 // Container meant to store all the relevant distributions corresponding to the cut variables.
217 // NB: overflow/underflow of variables on which we want to cut later should be 0!!!
218 const Int_t lNbSteps = 4 ;
219 const Int_t lNbVariables = 20;
220 //Array for the number of bins in each dimension :
221 Int_t lNbBinsPerVar[lNbVariables] = {0};
222 lNbBinsPerVar[0] = 25; //DcaCascDaughters : [0.0,2.4,3.0] -> Rec.Cut = 2.0;
223 lNbBinsPerVar[1] = 25; //DcaBachToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.01;
224 lNbBinsPerVar[2] = 30; //CascCosineOfPointingAngle : [0.97,1.0] -> Rec.Cut = 0.98;
225 lNbBinsPerVar[3] = 40; //CascRadius : [0.0,3.9,1000.0] -> Rec.Cut = 0.2;
226 lNbBinsPerVar[4] = 30; //InvMassLambdaAsCascDghter : [1.1,1.3] -> Rec.Cut = 0.008;
227 lNbBinsPerVar[5] = 20; //DcaV0Daughters : [0.0,2.0] -> Rec.Cut = 1.5;
228 lNbBinsPerVar[6] = 201; //V0CosineOfPointingAngleToPV : [0.89,1.0] -> Rec.Cut = 0.9;
229 lNbBinsPerVar[7] = 40; //V0Radius : [0.0,3.9,1000.0] -> Rec.Cut = 0.2;
230 lNbBinsPerVar[8] = 40; //DcaV0ToPrimVertex : [0.0,0.39,110.0] -> Rec.Cut = 0.01;
231 lNbBinsPerVar[9] = 25; //DcaPosToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.05;
232 lNbBinsPerVar[10] = 25; //DcaNegToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.05
233 lNbBinsPerVar[11] = 150; //InvMassXi : 2-MeV/c2 bins
234 lNbBinsPerVar[12] = 120; //InvMassOmega : 2-MeV/c2 bins
235 lNbBinsPerVar[13] = 100; //XiTransvMom : [0.0,10.0]
236 lNbBinsPerVar[14] = 110; //Y(Xi) : 0.02 in rapidity units
237 lNbBinsPerVar[15] = 110; //Y(Omega) : 0.02 in rapidity units
238 lNbBinsPerVar[16] = 112; //Proper lenght of cascade
239 lNbBinsPerVar[17] = 112; //Proper lenght of V0
240 lNbBinsPerVar[18] = 201; //V0CosineOfPointingAngleToXiV
241 lNbBinsPerVar[19] = 11; //Centrality
242 //define the container
243 fCFContCascadeCuts = new AliCFContainer("fCFContCascadeCuts","Container for Cascade cuts", lNbSteps, lNbVariables, lNbBinsPerVar );
244 //Setting the bin limits
246 Double_t *lBinLim0 = new Double_t[ lNbBinsPerVar[0] + 1 ];
247 for(Int_t i=0; i< lNbBinsPerVar[0]; i++) lBinLim0[i] = (Double_t)0.0 + (2.4 - 0.0)/(lNbBinsPerVar[0] - 1) * (Double_t)i;
248 lBinLim0[ lNbBinsPerVar[0] ] = 3.0;
249 fCFContCascadeCuts -> SetBinLimits(0, lBinLim0);
251 //1 - DcaToPrimVertexXi
252 Double_t *lBinLim1 = new Double_t[ lNbBinsPerVar[1] + 1 ];
253 for(Int_t i=0; i<lNbBinsPerVar[1]; i++) lBinLim1[i] = (Double_t)0.0 + (0.24 - 0.0)/(lNbBinsPerVar[1] - 1) * (Double_t)i;
254 lBinLim1[ lNbBinsPerVar[1] ] = 100.0;
255 fCFContCascadeCuts -> SetBinLimits(1, lBinLim1);
257 //2 - CascCosineOfPointingAngle
258 fCFContCascadeCuts->SetBinLimits(2, 0.97, 1.);
260 Double_t *lBinLim3 = new Double_t[ lNbBinsPerVar[3]+1 ];
261 for(Int_t i=0; i< lNbBinsPerVar[3]; i++) lBinLim3[i] = (Double_t)0.0 + (3.9 - 0.0 )/(lNbBinsPerVar[3] - 1) * (Double_t)i ;
262 lBinLim3[ lNbBinsPerVar[3] ] = 1000.0;
263 fCFContCascadeCuts -> SetBinLimits(3, lBinLim3 );
265 //4 - InvMassLambdaAsCascDghter
266 fCFContCascadeCuts->SetBinLimits(4, 1.1, 1.13);
268 fCFContCascadeCuts -> SetBinLimits(5, 0., 2.);
269 //6 - V0CosineOfPointingAngleToPV
270 fCFContCascadeCuts -> SetBinLimits(6, 0.8, 1.001);
272 Double_t *lBinLim7 = new Double_t[ lNbBinsPerVar[7] + 1];
273 for(Int_t i=0; i< lNbBinsPerVar[7];i++) lBinLim7[i] = (Double_t)0.0 + (3.9 - 0.0)/(lNbBinsPerVar[7] - 1) * (Double_t)i;
274 lBinLim7[ lNbBinsPerVar[7] ] = 1000.0;
275 fCFContCascadeCuts -> SetBinLimits(7, lBinLim7);
277 //8 - DcaV0ToPrimVertex
278 Double_t *lBinLim8 = new Double_t[ lNbBinsPerVar[8]+1 ];
279 for(Int_t i=0; i< lNbBinsPerVar[8];i++) lBinLim8[i] = (Double_t)0.0 + (0.39 - 0.0 )/(lNbBinsPerVar[8]-1) * (Double_t)i ;
280 lBinLim8[ lNbBinsPerVar[8] ] = 100.0;
281 fCFContCascadeCuts -> SetBinLimits(8, lBinLim8 );
283 //9 - DcaPosToPrimVertex
284 Double_t *lBinLim9 = new Double_t[ lNbBinsPerVar[9]+1 ];
285 for(Int_t i=0; i< lNbBinsPerVar[9];i++) lBinLim9[i] = (Double_t)0.0 + (0.24 - 0.0 )/(lNbBinsPerVar[9]-1) * (Double_t)i ;
286 lBinLim9[ lNbBinsPerVar[9] ] = 100.0;
287 fCFContCascadeCuts -> SetBinLimits(9, lBinLim9 );
289 //10 - DcaNegToPrimVertex
290 Double_t *lBinLim10 = new Double_t[ lNbBinsPerVar[10]+1 ];
291 for(Int_t i=0; i< lNbBinsPerVar[10];i++) lBinLim10[i] = (Double_t)0.0 + (0.24 - 0.0 )/(lNbBinsPerVar[10]-1) * (Double_t)i ;
292 lBinLim10[ lNbBinsPerVar[10] ] = 100.0;
293 fCFContCascadeCuts -> SetBinLimits(10, lBinLim10 );
296 fCFContCascadeCuts->SetBinLimits(11, 1.25, 1.40);
298 fCFContCascadeCuts->SetBinLimits(12, 1.62, 1.74);
300 fCFContCascadeCuts->SetBinLimits(13, 0.0, 10.0);
302 fCFContCascadeCuts->SetBinLimits(14, -1.1, 1.1);
304 fCFContCascadeCuts->SetBinLimits(15, -1.1, 1.1);
305 //16 - Proper time of cascade
306 Double_t *lBinLim16 = new Double_t[ lNbBinsPerVar[16]+1 ];
307 for(Int_t i=0; i< lNbBinsPerVar[16];i++) lBinLim16[i] = (Double_t) -1. + (110. + 1.0 ) / (lNbBinsPerVar[16] - 1) * (Double_t) i;
308 lBinLim16[ lNbBinsPerVar[16] ] = 2000.0;
309 fCFContCascadeCuts->SetBinLimits(16, lBinLim16);
310 //17 - Proper time of V0
311 fCFContCascadeCuts->SetBinLimits(17, lBinLim16);
312 //18 - V0CosineOfPointingAngleToXiV
313 fCFContCascadeCuts -> SetBinLimits(18, 0.8, 1.001);
315 Double_t *lBinLim19 = new Double_t[ lNbBinsPerVar[19]+1 ];
316 for(Int_t i=3; i< lNbBinsPerVar[19]+1;i++) lBinLim19[i] = (Double_t)(i-1)*10.;
320 fCFContCascadeCuts->SetBinLimits(19, lBinLim19 );
322 // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
323 fCFContCascadeCuts->SetStepTitle(0, "#Xi^{-} candidates");
324 fCFContCascadeCuts->SetStepTitle(1, "#bar{#Xi}^{+} candidates");
325 fCFContCascadeCuts->SetStepTitle(2, "#Omega^{-} candidates");
326 fCFContCascadeCuts->SetStepTitle(3, "#bar{#Omega}^{+} candidates");
327 // Setting the variable title, per axis
328 fCFContCascadeCuts->SetVarTitle(0, "Dca(cascade daughters) (cm)");
329 fCFContCascadeCuts->SetVarTitle(1, "ImpactParamToPV(bachelor) (cm)");
330 fCFContCascadeCuts->SetVarTitle(2, "cos(cascade PA)");
331 fCFContCascadeCuts->SetVarTitle(3, "R_{2d}(cascade decay) (cm)");
332 fCFContCascadeCuts->SetVarTitle(4, "M_{#Lambda}(as casc dghter) (GeV/c^{2})");
333 fCFContCascadeCuts->SetVarTitle(5, "Dca(V0 daughters) in Xi (cm)");
334 fCFContCascadeCuts->SetVarTitle(6, "cos(V0 PA) in cascade to PV");
335 fCFContCascadeCuts->SetVarTitle(7, "R_{2d}(V0 decay) (cm)");
336 fCFContCascadeCuts->SetVarTitle(8, "ImpactParamToPV(V0) (cm)");
337 fCFContCascadeCuts->SetVarTitle(9, "ImpactParamToPV(Pos) (cm)");
338 fCFContCascadeCuts->SetVarTitle(10, "ImpactParamToPV(Neg) (cm)");
339 fCFContCascadeCuts->SetVarTitle(11, "Inv. Mass(Xi) (GeV/c^{2})");
340 fCFContCascadeCuts->SetVarTitle(12, "Inv. Mass(Omega) (GeV/c^{2})");
341 fCFContCascadeCuts->SetVarTitle(13, "pt(cascade) (GeV/c)");
342 fCFContCascadeCuts->SetVarTitle(14, "Y(Xi)");
343 fCFContCascadeCuts->SetVarTitle(15, "Y(Omega)");
344 fCFContCascadeCuts->SetVarTitle(16, "mL/p (cascade) (cm)");
345 fCFContCascadeCuts->SetVarTitle(17, "mL/p (V0) (cm)");
346 fCFContCascadeCuts->SetVarTitle(18, "cos(V0 PA) in cascade to XiV");
347 fCFContCascadeCuts->SetVarTitle(19, "Centrality");
351 //-----------------------------------------------
352 // Define the Container for the MC generated info
353 //-----------------------------------------------
354 if(! fCFContCascadeMCgen) {
355 // Container meant to store general distributions for MC generated particles.
356 // NB: overflow/underflow of variables on which we want to cut later should be 0!!!
357 const Int_t lNbStepsMC = 4 ;
358 const Int_t lNbVariablesMC = 7;
359 //Array for the number of bins in each dimension :
360 Int_t lNbBinsPerVarMC[lNbVariablesMC] = {0};
361 lNbBinsPerVarMC[0] = 100; //Total momentum : [0.0,10.0]
362 lNbBinsPerVarMC[1] = 100; //Transverse momentum : [0.0,10.0]
363 lNbBinsPerVarMC[2] = 110; //Y : [-1.1,1.1]
364 lNbBinsPerVarMC[3] = 200; //eta : [-10, 10]
365 lNbBinsPerVarMC[4] = 200; //theta : [-10, 190]
366 lNbBinsPerVarMC[5] = 360; //Phi : [0., 360.]
367 lNbBinsPerVarMC[6] = 11; //Centrality
368 //define the container
369 fCFContCascadeMCgen = new AliCFContainer("fCFContCascadeMCgen","Container for MC gen cascade ", lNbStepsMC, lNbVariablesMC, lNbBinsPerVarMC );
370 //Setting the bin limits
372 fCFContCascadeMCgen->SetBinLimits(0, 0.0, 10.0);
373 //1 - Transverse Momentum
374 fCFContCascadeMCgen->SetBinLimits(1, 0.0, 10.0);
376 fCFContCascadeMCgen->SetBinLimits(2, -1.1, 1.1);
378 fCFContCascadeMCgen->SetBinLimits(3, -10, 10);
380 fCFContCascadeMCgen->SetBinLimits(4, -10, 190);
382 fCFContCascadeMCgen->SetBinLimits(5, 0.0, 360.0);
384 Double_t *lBinLim6MC = new Double_t[ lNbBinsPerVarMC[6]+1 ];
385 for(Int_t i=3; i< lNbBinsPerVarMC[6]+1;i++) lBinLim6MC[i] = (Double_t)(i-1)*10.;
388 lBinLim6MC[2] = 10.0;
389 fCFContCascadeMCgen->SetBinLimits(6, lBinLim6MC);
390 delete [] lBinLim6MC;
391 // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
392 fCFContCascadeMCgen->SetStepTitle(0, "#Xi^{-} candidates");
393 fCFContCascadeMCgen->SetStepTitle(1, "#bar{#Xi}^{+} candidates");
394 fCFContCascadeMCgen->SetStepTitle(2, "#Omega^{-} candidates");
395 fCFContCascadeMCgen->SetStepTitle(3, "#bar{#Omega}^{+} candidates");
396 // Setting the variable title, per axis
397 fCFContCascadeMCgen->SetVarTitle(0, "MC gen p_tot (GeV/c)");
398 fCFContCascadeMCgen->SetVarTitle(1, "MC gen p_T (GeV/c)");
399 fCFContCascadeMCgen->SetVarTitle(2, "MC gen Rapidity");
400 fCFContCascadeMCgen->SetVarTitle(3, "MC gen Pseudo-rapidity");
401 fCFContCascadeMCgen->SetVarTitle(4, "MC gen Theta");
402 fCFContCascadeMCgen->SetVarTitle(5, "MC gen Phi");
403 fCFContCascadeMCgen->SetVarTitle(6, "MC gen Centrality");
406 PostData(1, fListHistMultistrangeQA);
407 PostData(2, fCFContCascadeCuts);
408 PostData(3, fCFContCascadeMCgen);
411 }// end UserCreateOutputObjects
414 //________________________________________________________________________
415 void AliAnalysisTaskQAMultistrange::UserExec(Option_t *) {
418 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
419 // Main loop (called for each event)
420 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
422 //------------------------
423 // Define ESD/AOD handlers
424 //------------------------
425 AliESDEvent *lESDevent = 0x0;
426 AliAODEvent *lAODevent = 0x0;
427 AliMCEvent *lMCevent = 0x0; //for MC
428 AliStack *lMCstack = 0x0; //for MC
429 TClonesArray *arrayMC = 0; //for MC
435 //----------------------
436 // Check the PIDresponse
437 //----------------------
439 AliError("Cannot get pid response");
443 //---------------------
444 // Check the InputEvent
445 //---------------------
446 if (fAnalysisType == "ESD") {
447 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
449 AliWarning("ERROR: lESDevent not available \n");
453 lMCevent = MCEvent();
455 Printf("ERROR: Could not retrieve MC event \n");
456 cout << "Name of the file with pb :" << CurrentFileName() << endl;
459 lMCstack = lMCevent->Stack();
461 Printf("ERROR: Could not retrieve MC stack \n");
462 cout << "Name of the file with pb :" << CurrentFileName() << endl;
466 } else if (fAnalysisType == "AOD") {
467 lAODevent = dynamic_cast<AliAODEvent*>( InputEvent() );
469 AliWarning("ERROR: lAODevent not available \n");
473 arrayMC = (TClonesArray*) lAODevent->GetList()->FindObject(AliAODMCParticle::StdBranchName());
474 if (!arrayMC) AliFatal("Error: MC particles branch not found!\n");
477 Printf("Analysis type (ESD or AOD) not specified \n");
481 fHistEventSel->Fill(0.5);
490 UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
491 Bool_t isSelected = 0;
492 if (fCollidingSystem == "pp" || fCollidingSystem == "PbPb") isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
493 else if (fCollidingSystem == "pPb") isSelected = (maskIsSelected & AliVEvent::kINT7) == AliVEvent::kINT7;
495 PostData(1, fListHistMultistrangeQA);
496 PostData(2, fCFContCascadeCuts);
497 PostData(3, fCFContCascadeMCgen);
501 fHistEventSel->Fill(1.5);
503 //----------------------------------------------------------
504 // Centrality/Multiplicity selection for PbPb/pPb collisions
505 //----------------------------------------------------------
507 AliCentrality* centrality = 0;
508 if (fAnalysisType == "ESD" && (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb")) centrality = lESDevent->GetCentrality();
509 else if (fAnalysisType == "AOD" && (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb")) centrality = lAODevent->GetCentrality();
510 Float_t lcentrality = 0.;
511 if (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb") {
512 if (fkUseCleaning) lcentrality = centrality->GetCentralityPercentile(fCentrEstimator.Data());
514 lcentrality = centrality->GetCentralityPercentileUnchecked(fCentrEstimator.Data());
515 if (centrality->GetQuality()>1) {
516 PostData(1, fListHistMultistrangeQA);
517 PostData(2, fCFContCascadeCuts);
518 PostData(3, fCFContCascadeMCgen);
522 //if (lcentrality<fCentrLowLim||lcentrality>=fCentrUpLim) {
523 // PostData(1, fCFContCascadeCuts);
526 } else if (fCollidingSystem == "pp") lcentrality = 0.;
528 fHistEventSel->Fill(2.5);
530 //------------------------------
531 // Well-established PV selection
532 //------------------------------
533 Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
534 Double_t lMagneticField = -10.;
535 if (fCollidingSystem == "PbPb" || fCollidingSystem == "pp") {
536 if (fAnalysisType == "ESD") {
537 const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
538 const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
539 if (fkQualityCutNoTPConlyPrimVtx) {
540 const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
541 if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
542 AliWarning(" No SPD prim. vertex nor prim. Tracking vertex ... return !");
543 PostData(1, fListHistMultistrangeQA);
544 PostData(2, fCFContCascadeCuts);
545 PostData(3, fCFContCascadeMCgen);
549 lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
550 } else if (fAnalysisType == "AOD") {
551 const AliAODVertex *lPrimaryBestAODVtx = lAODevent->GetPrimaryVertex();
552 if (!lPrimaryBestAODVtx){
553 AliWarning("No prim. vertex in AOD... return!");
554 PostData(1, fListHistMultistrangeQA);
555 PostData(2, fCFContCascadeCuts);
556 PostData(3, fCFContCascadeMCgen);
559 lPrimaryBestAODVtx->GetXYZ( lBestPrimaryVtxPos );
561 } else if (fCollidingSystem == "pPb") {
562 if (fAnalysisType == "ESD") {
563 Bool_t fHasVertex = kFALSE;
564 const AliESDVertex *vertex = lESDevent->GetPrimaryVertexTracks();
565 if (vertex->GetNContributors() < 1) {
566 vertex = lESDevent->GetPrimaryVertexSPD();
567 if (vertex->GetNContributors() < 1) fHasVertex = kFALSE;
568 else fHasVertex = kTRUE;
569 TString vtxTyp = vertex->GetTitle();
571 vertex->GetCovarianceMatrix(cov);
572 Double_t zRes = TMath::Sqrt(cov[5]);
573 if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fHasVertex = kFALSE;
575 else fHasVertex = kTRUE;
576 if (fHasVertex == kFALSE) { //Is First event in chunk rejection: Still present!
577 AliWarning("Pb / | PV does not satisfy selection criteria!");
578 PostData(1, fListHistMultistrangeQA);
579 PostData(2, fCFContCascadeCuts);
580 PostData(3, fCFContCascadeMCgen);
583 vertex->GetXYZ( lBestPrimaryVtxPos );
584 } else if (fAnalysisType == "AOD") {
585 Bool_t fHasVertex = kFALSE;
586 const AliAODVertex *vertex = lAODevent->GetPrimaryVertex();
587 if (vertex->GetNContributors() < 1) {
588 vertex = lAODevent->GetPrimaryVertexSPD();
589 if (vertex->GetNContributors() < 1) fHasVertex = kFALSE;
590 else fHasVertex = kTRUE;
591 TString vtxTyp = vertex->GetTitle();
593 vertex->GetCovarianceMatrix(cov);
594 Double_t zRes = TMath::Sqrt(cov[5]);
595 if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fHasVertex = kFALSE;
597 else fHasVertex = kTRUE;
598 if (fHasVertex == kFALSE) { //Is First event in chunk rejection: Still present!
599 AliWarning("Pb / | PV does not satisfy selection criteria!");
600 PostData(1, fListHistMultistrangeQA);
601 PostData(2, fCFContCascadeCuts);
602 PostData(3, fCFContCascadeMCgen);
605 vertex->GetXYZ( lBestPrimaryVtxPos );
609 if (fAnalysisType == "ESD") lMagneticField = lESDevent->GetMagneticField();
610 else if (fAnalysisType == "AOD") lMagneticField = lAODevent->GetMagneticField();
612 fHistEventSel->Fill(3.5);
614 //----------------------------
615 // Vertex Z position selection
616 //----------------------------
617 if (fkQualityCutZprimVtxPos) {
618 if (TMath::Abs(lBestPrimaryVtxPos[2]) > fVtxRange ) {
619 PostData(1, fListHistMultistrangeQA);
620 PostData(2, fCFContCascadeCuts);
621 PostData(3, fCFContCascadeMCgen);
626 fHistEventSel->Fill(4.5);
628 //-----------------------
629 // Pilup selection for pp
630 //-----------------------
631 if (fCollidingSystem == "pp") {
632 if (fAnalysisType == "ESD") {
633 if (fkQualityCutPileup) { if(lESDevent->IsPileupFromSPD()){
634 PostData(1, fListHistMultistrangeQA);
635 PostData(2, fCFContCascadeCuts);
636 PostData(3, fCFContCascadeMCgen);
640 } else if (fAnalysisType == "AOD") {
641 if (fkQualityCutPileup) { if(lAODevent->IsPileupFromSPD()){
642 PostData(1, fListHistMultistrangeQA);
643 PostData(2, fCFContCascadeCuts);
644 PostData(3, fCFContCascadeMCgen);
651 fHistEventSel->Fill(5.5);
653 ////////////////////////////
654 // MC GENERATED CASCADE PART
655 ////////////////////////////
661 Int_t lNbMCPrimary = 0;
662 if (fAnalysisType == "ESD") lNbMCPrimary = lMCstack->GetNprimary(); //lMCstack->GetNprimary(); or lMCstack->GetNtrack();
663 else if (fAnalysisType == "AOD") lNbMCPrimary = arrayMC->GetEntries();
665 for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++) {
668 Double_t partPt = 0.;
669 Double_t partEta = 0.;
670 Double_t partTheta = 0.;
671 Double_t partPhi = 0.;
672 Double_t partRap = 0.;
673 Double_t partEnergy = 0.; //for Rapidity
674 Double_t partPz = 0.; //for Rapidity
677 if ( fAnalysisType == "ESD" ) {
678 TParticle* lCurrentParticlePrimary = 0x0;
679 lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
680 if (!lCurrentParticlePrimary) {
681 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
684 if (!lMCstack->IsPhysicalPrimary(iCurrentLabelStack)) continue;
685 TParticle* cascMC = 0x0;
686 cascMC = lCurrentParticlePrimary;
688 Printf("MC TParticle pointer to Cascade = 0x0 ! Skip ...");
692 partPt = cascMC->Pt();
693 partEta = cascMC->Eta();
694 partTheta = cascMC->Theta()*180.0/TMath::Pi();
695 partPhi = cascMC->Phi()*180.0/TMath::Pi();
696 partEnergy = cascMC->Energy();
697 partPz = cascMC->Pz();
698 PDGcode = lCurrentParticlePrimary->GetPdgCode();
699 } else if ( fAnalysisType == "AOD" ) {
700 AliAODMCParticle *lCurrentParticleaod = (AliAODMCParticle*) arrayMC->At(iCurrentLabelStack);
701 if (!lCurrentParticleaod) {
702 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
705 if (!lCurrentParticleaod->IsPhysicalPrimary()) continue;
706 partP = lCurrentParticleaod->P();
707 partPt = lCurrentParticleaod->Pt();
708 partEta = lCurrentParticleaod->Eta();
709 partTheta = lCurrentParticleaod->Theta()*180.0/TMath::Pi();
710 partPhi = lCurrentParticleaod->Phi()*180.0/TMath::Pi();
711 partEnergy = lCurrentParticleaod->E();
712 partPz = lCurrentParticleaod->Pz();
713 PDGcode = lCurrentParticleaod->GetPdgCode();
715 partRap = 0.5*TMath::Log((partEnergy + partPz) / (partEnergy - partPz + 1.e-13));
717 // Fill the AliCFContainer
718 Double_t lContainerCutVarsMC[7] = {0.0};
719 lContainerCutVarsMC[0] = partP;
720 lContainerCutVarsMC[1] = partPt;
721 lContainerCutVarsMC[2] = partRap;
722 lContainerCutVarsMC[3] = partEta;
723 lContainerCutVarsMC[4] = partTheta;
724 lContainerCutVarsMC[5] = partPhi;
725 lContainerCutVarsMC[6] = lcentrality;
726 if (PDGcode == 3312) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,0); // for Xi-
727 if (PDGcode == -3312) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,1); // for Xi+
728 if (PDGcode == 3334) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,2); // for Omega-
729 if (PDGcode == -3334) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,3); // for Omega+
735 //////////////////////////////
736 // CASCADE RECONSTRUCTION PART
737 //////////////////////////////
742 if (fAnalysisType == "ESD") ncascades = lESDevent->GetNumberOfCascades();
743 else if (fAnalysisType == "AOD") ncascades = lAODevent->GetNumberOfCascades();
745 for (Int_t iXi = 0; iXi < ncascades; iXi++) {// This is the begining of the Cascade loop (ESD or AOD)
747 // -------------------------------------
748 // - Initialisation of the local variables that will be needed for ESD/AOD
749 // -- Container variables (1st round)
750 Double_t lDcaXiDaughters = -1. ; //[Container]
751 Double_t lXiCosineOfPointingAngle = -1. ; //[Container]
752 Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 }; //Useful to define other variables: radius fid. vol., ctau, etc. for cascade
753 Double_t lXiRadius = -1000. ; //[Container]
754 UShort_t lPosTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
755 UShort_t lNegTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
756 UShort_t lBachTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
757 Double_t lInvMassLambdaAsCascDghter = 0.; //[Container]
758 Double_t lDcaV0DaughtersXi = -1.; //[Container]
759 Double_t lDcaBachToPrimVertexXi = -1.; //[Container]
760 Double_t lDcaV0ToPrimVertexXi = -1.; //[Container]
761 Double_t lDcaPosToPrimVertexXi = -1.; //[Container]
762 Double_t lDcaNegToPrimVertexXi = -1.; //[Container]
763 Double_t lV0CosineOfPointingAngle = -1.; //[Container]
764 Double_t lV0toXiCosineOfPointingAngle = -1.; //[Container]
765 Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; //Useful to define other variables: radius fid. vol., ctau, etc. for VO
766 Double_t lV0RadiusXi = -1000.0; //[Container]
767 Double_t lV0quality = 0.; // ??
768 Double_t lInvMassXiMinus = 0.; //[Container]
769 Double_t lInvMassXiPlus = 0.; //[Container]
770 Double_t lInvMassOmegaMinus = 0.; //[Container]
771 Double_t lInvMassOmegaPlus = 0.; //[Container]
773 Bool_t lIsBachelorKaonForTPC = kFALSE;
774 Bool_t lIsBachelorPionForTPC = kFALSE;
775 Bool_t lIsNegPionForTPC = kFALSE;
776 Bool_t lIsPosPionForTPC = kFALSE;
777 Bool_t lIsNegProtonForTPC = kFALSE;
778 Bool_t lIsPosProtonForTPC = kFALSE;
779 // -- More container variables and quality checks
780 Double_t lXiMomX = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
781 Double_t lXiMomY = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
782 Double_t lXiMomZ = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
783 Double_t lXiTransvMom = 0.; //[Container]
784 Double_t lXiTotMom = 0.; //Useful to define other variables: cTau
785 Double_t lV0PMomX = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
786 Double_t lV0PMomY = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
787 Double_t lV0PMomZ = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
788 Double_t lV0NMomX = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
789 Double_t lV0NMomY = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
790 Double_t lV0NMomZ = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
791 Double_t lV0TotMom = 0.; //Useful to define other variables: lctauV0
792 Double_t lBachMomX = 0.; //Useful to define other variables: lBachTransvMom
793 Double_t lBachMomY = 0.; //Useful to define other variables: lBachTransvMom
794 Double_t lBachMomZ = 0.; //Useful to define other variables: lBachTransvMom
795 Double_t lBachTransvMom = 0.; //Selection on the min bachelor pT
796 Double_t lpTrackTransvMom = 0.; //Selection on the min bachelor pT
797 Double_t lnTrackTransvMom = 0.; //Selection on the min bachelor pT
798 Short_t lChargeXi = -2; //Useful to select the particles based on the charge
799 Double_t lRapXi = -20.0; //[Container]
800 Double_t lRapOmega = -20.0; //[Container]
801 Float_t etaBach = 0.; //Selection on the eta range
802 Float_t etaPos = 0.; //Selection on the eta range
803 Float_t etaNeg = 0.; //Selection on the eta range
804 // -- variables for the AliCFContainer dedicated to cascade cut optmisiation: ESD and AOD
805 if (fAnalysisType == "ESD") {
807 // -------------------------------------
808 // - Load the cascades from the handler
809 AliESDcascade *xi = lESDevent->GetCascade(iXi);
812 // ---------------------------------------------------------------------------
813 // - Assigning the necessary variables for specific AliESDcascade data members
815 xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
816 lDcaXiDaughters = xi->GetDcaXiDaughters();
817 lXiCosineOfPointingAngle = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
818 // Take care : the best available vertex should be used (like in AliCascadeVertexer)
819 xi->GetXYZcascade( lPosXi[0], lPosXi[1], lPosXi[2] );
820 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
822 // -------------------------------------------------------------------------------------------------------------------------------
823 // - Around the tracks : Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
824 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( xi->GetPindex() );
825 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( xi->GetNindex() );
826 UInt_t lBachIdx = (UInt_t) TMath::Abs( xi->GetBindex() );
827 // Care track label can be negative in MC production (linked with the track quality)
828 // However = normally, not the case for track index ...
829 // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
830 if (lBachIdx == lIdxNegXi) continue;
831 if (lBachIdx == lIdxPosXi) continue;
832 // - Get the track for the daughters
833 AliESDtrack *pTrackXi = lESDevent->GetTrack( lIdxPosXi );
834 AliESDtrack *nTrackXi = lESDevent->GetTrack( lIdxNegXi );
835 AliESDtrack *bachTrackXi = lESDevent->GetTrack( lBachIdx );
836 if (!pTrackXi || !nTrackXi || !bachTrackXi ) continue;
837 // - Get the TPCnumber of cluster for the daughters
838 lPosTPCClusters = pTrackXi->GetTPCNcls();
839 lNegTPCClusters = nTrackXi->GetTPCNcls();
840 lBachTPCClusters = bachTrackXi->GetTPCNcls();
842 // ------------------------------------
843 // - Rejection of a poor quality tracks
844 if (fkQualityCutTPCrefit) {
845 // 1 - Poor quality related to TPCrefit
846 ULong_t pStatus = pTrackXi->GetStatus();
847 ULong_t nStatus = nTrackXi->GetStatus();
848 ULong_t bachStatus = bachTrackXi->GetStatus();
849 if ((pStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
850 if ((nStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
851 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" Bach. track has no TPCrefit ... continue!"); continue; }
853 if (fkQualityCutnTPCcls) {
854 // 2 - Poor quality related to TPC clusters
855 if (lPosTPCClusters < fMinnTPCcls) { AliWarning(" V0 Pos. track has less than minn TPC clusters ... continue!"); continue; }
856 if (lNegTPCClusters < fMinnTPCcls) { AliWarning(" V0 Neg. track has less than minn TPC clusters ... continue!"); continue; }
857 if (lBachTPCClusters < fMinnTPCcls) { AliWarning(" Bach. track has less than minn TPC clusters ... continue!"); continue; }
860 // ------------------------------
861 etaPos = pTrackXi->Eta();
862 etaNeg = nTrackXi->Eta();
863 etaBach = bachTrackXi->Eta();
864 lInvMassLambdaAsCascDghter = xi->GetEffMass();
865 lDcaV0DaughtersXi = xi->GetDcaV0Daughters();
866 lV0CosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
867 lDcaV0ToPrimVertexXi = xi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
868 lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField) );
869 xi->GetXYZ( lPosV0Xi[0], lPosV0Xi[1], lPosV0Xi[2] );
870 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
871 lDcaPosToPrimVertexXi = TMath::Abs( pTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) );
872 lDcaNegToPrimVertexXi = TMath::Abs( nTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) );
874 //----------------------------------------------------------------------------------------------------
875 // - Around effective masses. Change mass hypotheses to cover all the possibilities: Xi-/+, Omega -/+
876 if (bachTrackXi->Charge() < 0) {
877 // Calculate the effective mass of the Xi- candidate. pdg code 3312 = Xi-
879 xi->ChangeMassHypothesis(lV0quality , 3312);
880 lInvMassXiMinus = xi->GetEffMassXi();
881 // Calculate the effective mass of the Xi- candidate. pdg code 3334 = Omega-
883 xi->ChangeMassHypothesis(lV0quality , 3334);
884 lInvMassOmegaMinus = xi->GetEffMassXi();
885 // Back to default hyp.
887 xi->ChangeMassHypothesis(lV0quality , 3312);
888 }// end if negative bachelor
889 if ( bachTrackXi->Charge() > 0 ) {
890 // Calculate the effective mass of the Xi+ candidate. pdg code -3312 = Xi+
892 xi->ChangeMassHypothesis(lV0quality , -3312);
893 lInvMassXiPlus = xi->GetEffMassXi();
894 // Calculate the effective mass of the Xi+ candidate. pdg code -3334 = Omega+
896 xi->ChangeMassHypothesis(lV0quality , -3334);
897 lInvMassOmegaPlus = xi->GetEffMassXi();
898 // Back to "default" hyp.
900 xi->ChangeMassHypothesis(lV0quality , -3312);
901 }// end if positive bachelor
903 // ----------------------------------------------
904 // - TPC PID : 3-sigma bands on Bethe-Bloch curve
906 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
907 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
908 // Negative V0 daughter
909 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion )) < 4) lIsNegPionForTPC = kTRUE;
910 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
911 // Positive V0 daughter
912 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion )) < 4) lIsPosPionForTPC = kTRUE;
913 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
915 // ------------------------------
916 // - Miscellaneous pieces of info that may help regarding data quality assessment.
917 xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
918 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
919 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
920 xi->GetNPxPyPz(lV0NMomX,lV0NMomY,lV0NMomZ);
921 xi->GetPPxPyPz(lV0PMomX,lV0PMomY,lV0PMomZ);
922 lV0TotMom = TMath::Sqrt(TMath::Power(lV0NMomX+lV0PMomX,2)+TMath::Power(lV0NMomY+lV0PMomY,2)+TMath::Power(lV0NMomZ+lV0PMomZ,2));
923 xi->GetBPxPyPz( lBachMomX, lBachMomY, lBachMomZ );
924 lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
925 lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX + lV0NMomY*lV0NMomY );
926 lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX + lV0PMomY*lV0PMomY );
927 lChargeXi = xi->Charge();
928 lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
929 lRapXi = xi->RapXi();
930 lRapOmega = xi->RapOmega();
932 } else if (fAnalysisType == "AOD") {
934 // -------------------------------------
935 // - Load the cascades from the handler
936 const AliAODcascade *xi = lAODevent->GetCascade(iXi);
939 //----------------------------------------------------------------------------
940 // - Assigning the necessary variables for specific AliESDcascade data members
941 lDcaXiDaughters = xi->DcaXiDaughters();
942 lXiCosineOfPointingAngle = xi->CosPointingAngleXi( lBestPrimaryVtxPos[0],
943 lBestPrimaryVtxPos[1],
944 lBestPrimaryVtxPos[2] );
945 lPosXi[0] = xi->DecayVertexXiX();
946 lPosXi[1] = xi->DecayVertexXiY();
947 lPosXi[2] = xi->DecayVertexXiZ();
948 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
950 //-------------------------------------------------------------------------------------------------------------------------------
951 // - Around the tracks: Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
952 AliAODTrack *pTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );
953 AliAODTrack *nTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );
954 AliAODTrack *bachTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );
955 if (!pTrackXi || !nTrackXi || !bachTrackXi ) continue;
956 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( pTrackXi->GetID() );
957 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( nTrackXi->GetID() );
958 UInt_t lBachIdx = (UInt_t) TMath::Abs( bachTrackXi->GetID() );
959 // Care track label can be negative in MC production (linked with the track quality)
960 // However = normally, not the case for track index ...
962 // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
963 if (lBachIdx == lIdxNegXi) continue;
964 if (lBachIdx == lIdxPosXi) continue;
965 // - Get the TPCnumber of cluster for the daughters
966 lPosTPCClusters = pTrackXi->GetTPCNcls();
967 lNegTPCClusters = nTrackXi->GetTPCNcls();
968 lBachTPCClusters = bachTrackXi->GetTPCNcls();
970 // ------------------------------------
971 // - Rejection of a poor quality tracks
972 if (fkQualityCutTPCrefit) {
973 // - Poor quality related to TPCrefit
974 if (!(pTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
975 if (!(nTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
976 if (!(bachTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" Bach. track has no TPCrefit ... continue!"); continue; }
978 if (fkQualityCutnTPCcls) {
979 // - Poor quality related to TPC clusters
980 if (lPosTPCClusters < fMinnTPCcls) continue;
981 if (lNegTPCClusters < fMinnTPCcls) continue;
982 if (lBachTPCClusters < fMinnTPCcls) continue;
985 // ------------------------------------------------------------------------------------------------------------------------------
986 // - Around the tracks: Bach + V0 (AOD). Necessary variables for AODcascade data members coming from the AODv0 part (inheritance)
987 etaPos = pTrackXi->Eta();
988 etaNeg = nTrackXi->Eta();
989 etaBach = bachTrackXi->Eta();
990 lChargeXi = xi->ChargeXi();
991 if ( lChargeXi < 0) lInvMassLambdaAsCascDghter = xi->MassLambda();
992 else lInvMassLambdaAsCascDghter = xi->MassAntiLambda();
993 lDcaV0DaughtersXi = xi->DcaV0Daughters();
994 lDcaV0ToPrimVertexXi = xi->DcaV0ToPrimVertex();
995 lDcaBachToPrimVertexXi = xi->DcaBachToPrimVertex();
996 lPosV0Xi[0] = xi->DecayVertexV0X();
997 lPosV0Xi[1] = xi->DecayVertexV0Y();
998 lPosV0Xi[2] = xi->DecayVertexV0Z();
999 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
1000 lV0CosineOfPointingAngle = xi->CosPointingAngle( lBestPrimaryVtxPos );
1001 lDcaPosToPrimVertexXi = xi->DcaPosToPrimVertex();
1002 lDcaNegToPrimVertexXi = xi->DcaNegToPrimVertex();
1004 // ---------------------------------------------------------------------------------------------------
1005 // - Around effective masses. Change mass hypotheses to cover all the possibilities: Xi-/+, Omega -/+
1006 if ( lChargeXi < 0 ) lInvMassXiMinus = xi->MassXi();
1007 if ( lChargeXi > 0 ) lInvMassXiPlus = xi->MassXi();
1008 if ( lChargeXi < 0 ) lInvMassOmegaMinus = xi->MassOmega();
1009 if ( lChargeXi > 0 ) lInvMassOmegaPlus = xi->MassOmega();
1011 // ----------------------------------------------
1012 // - TPC PID : 3-sigma bands on Bethe-Bloch curve
1014 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
1015 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
1016 // Negative V0 daughter
1017 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion )) < 4) lIsNegPionForTPC = kTRUE;
1018 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
1019 // Positive V0 daughter
1020 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion )) < 4) lIsPosPionForTPC = kTRUE;
1021 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
1023 //---------------------------------
1024 // - Extra info for QA (AOD)
1025 // Miscellaneous pieces of info that may help regarding data quality assessment.
1026 // Cascade transverse and total momentum
1027 lXiMomX = xi->MomXiX();
1028 lXiMomY = xi->MomXiY();
1029 lXiMomZ = xi->MomXiZ();
1030 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
1031 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
1032 Double_t lV0MomX = xi->MomV0X();
1033 Double_t lV0MomY = xi->MomV0Y();
1034 Double_t lV0MomZ = xi->MomV0Z();
1035 lV0TotMom = TMath::Sqrt(TMath::Power(lV0MomX,2)+TMath::Power(lV0MomY,2)+TMath::Power(lV0MomZ,2));
1036 lBachMomX = xi->MomBachX();
1037 lBachMomY = xi->MomBachY();
1038 lBachMomZ = xi->MomBachZ();
1039 lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
1040 lV0NMomX = xi->MomNegX();
1041 lV0NMomY = xi->MomNegY();
1042 lV0PMomX = xi->MomPosX();
1043 lV0PMomY = xi->MomPosY();
1044 lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX + lV0NMomY*lV0NMomY );
1045 lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX + lV0PMomY*lV0PMomY );
1046 lV0toXiCosineOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );
1047 lRapXi = xi->RapXi();
1048 lRapOmega = xi->RapOmega();
1050 }// end of AOD treatment
1052 //---------------------------------------
1053 // Cut on pt of the three daughter tracks
1054 if (lBachTransvMom<fMinPtCutOnDaughterTracks) continue;
1055 if (lpTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
1056 if (lnTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
1058 //---------------------------------------------------
1059 // Cut on pseudorapidity of the three daughter tracks
1060 if (TMath::Abs(etaBach)>fEtaCutOnDaughterTracks) continue;
1061 if (TMath::Abs(etaPos)>fEtaCutOnDaughterTracks) continue;
1062 if (TMath::Abs(etaNeg)>fEtaCutOnDaughterTracks) continue;
1064 //----------------------------------
1065 // Calculate proper time for cascade
1066 Double_t cascadeMass = 0.;
1067 if ( ( (lChargeXi<0) && lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC ) ||
1068 ( (lChargeXi>0) && lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC ) ) cascadeMass = 1.321;
1069 if ( ( (lChargeXi<0) && lIsBachelorKaonForTPC && lIsPosProtonForTPC && lIsNegPionForTPC ) ||
1070 ( (lChargeXi>0) && lIsBachelorKaonForTPC && lIsNegProtonForTPC && lIsPosPionForTPC ) ) cascadeMass = 1.672;
1071 Double_t lctau = TMath::Sqrt(TMath::Power((lPosXi[0]-lBestPrimaryVtxPos[0]),2)+TMath::Power((lPosXi[1]-lBestPrimaryVtxPos[1]),2)+TMath::Power(( lPosXi[2]-lBestPrimaryVtxPos[2]),2));
1072 if (lXiTotMom!=0) lctau = lctau*cascadeMass/lXiTotMom;
1074 // Calculate proper time for Lambda (reconstructed)
1075 Float_t lambdaMass = 1.115683; // PDG mass
1076 Float_t distV0Xi = TMath::Sqrt(TMath::Power((lPosV0Xi[0]-lPosXi[0]),2)+TMath::Power((lPosV0Xi[1]-lPosXi[1]),2)+TMath::Power((lPosV0Xi[2]-lPosXi[2]),2));
1077 Float_t lctauV0 = -1.;
1078 if (lV0TotMom!=0) lctauV0 = distV0Xi*lambdaMass/lV0TotMom;
1080 // Fill the TH1F without PID info
1081 if ( lChargeXi < 0 ) {
1082 fHistMassXiMinus->Fill( lInvMassXiMinus );
1083 fHistMassOmegaMinus->Fill( lInvMassOmegaMinus );
1084 } else if ( lChargeXi > 0 ) {
1085 fHistMassXiPlus->Fill( lInvMassXiPlus );
1086 fHistMassOmegaPlus->Fill( lInvMassOmegaPlus );
1090 // Fill the AliCFContainer (optimisation of topological selections)
1091 Double_t lContainerCutVars[21] = {0.0};
1092 lContainerCutVars[0] = lDcaXiDaughters;
1093 lContainerCutVars[1] = lDcaBachToPrimVertexXi;
1094 lContainerCutVars[2] = lXiCosineOfPointingAngle;
1095 lContainerCutVars[3] = lXiRadius;
1096 lContainerCutVars[4] = lInvMassLambdaAsCascDghter;
1097 lContainerCutVars[5] = lDcaV0DaughtersXi;
1098 lContainerCutVars[6] = lV0CosineOfPointingAngle;
1099 lContainerCutVars[7] = lV0RadiusXi;
1100 lContainerCutVars[8] = lDcaV0ToPrimVertexXi;
1101 lContainerCutVars[9] = lDcaPosToPrimVertexXi;
1102 lContainerCutVars[10] = lDcaNegToPrimVertexXi;
1103 lContainerCutVars[13] = lXiTransvMom;
1104 lContainerCutVars[16] = lctau;
1105 lContainerCutVars[17] = lctauV0;
1106 lContainerCutVars[18] = lV0toXiCosineOfPointingAngle;
1107 lContainerCutVars[19] = lcentrality;
1108 if ( lChargeXi < 0 ) {
1109 lContainerCutVars[11] = lInvMassXiMinus;
1110 lContainerCutVars[12] = lInvMassOmegaMinus;
1111 lContainerCutVars[14] = lRapXi;
1112 lContainerCutVars[15] = -1.;
1113 if (lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,0); // for Xi-
1114 lContainerCutVars[11] = lInvMassXiMinus;
1115 lContainerCutVars[12] = lInvMassOmegaMinus;
1116 lContainerCutVars[14] = -1.;
1117 lContainerCutVars[15] = lRapOmega;
1118 if (lIsBachelorKaonForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,2); // for Omega-
1120 lContainerCutVars[11] = lInvMassXiPlus;
1121 lContainerCutVars[12] = lInvMassOmegaPlus;
1122 lContainerCutVars[14] = lRapXi;
1123 lContainerCutVars[15] = -1.;
1124 if (lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,1); // for Xi+
1125 lContainerCutVars[11] = lInvMassXiPlus;
1126 lContainerCutVars[12] = lInvMassOmegaPlus;
1127 lContainerCutVars[14] = -1.;
1128 lContainerCutVars[15] = lRapOmega;
1129 if (lIsBachelorKaonForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,3); // for Omega+
1133 }// end of the Cascade loop (ESD or AOD)
1136 // Post output data.
1137 PostData(1, fListHistMultistrangeQA);
1138 PostData(2, fCFContCascadeCuts);
1139 PostData(3, fCFContCascadeMCgen);
1143 //________________________________________________________________________
1144 void AliAnalysisTaskQAMultistrange::Terminate(Option_t *) {