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"
64 #include "AliAnalysisTaskQAMultistrange.h"
66 ClassImp(AliAnalysisTaskQAMultistrange)
70 //________________________________________________________________________
71 AliAnalysisTaskQAMultistrange::AliAnalysisTaskQAMultistrange()
72 : AliAnalysisTaskSE (),
74 fAnalysisType ("ESD"),
75 fCollidingSystem ("PbPb"),
77 fkSDDSelectionOn (kTRUE),
78 fkQualityCutZprimVtxPos (kTRUE),
79 fkQualityCutNoTPConlyPrimVtx(kTRUE),
80 fkQualityCutTPCrefit (kTRUE),
81 fkQualityCutnTPCcls (kTRUE),
82 fkQualityCutPileup (kTRUE),
90 fMinPtCutOnDaughterTracks (0),
91 fEtaCutOnDaughterTracks (0),
93 fCFContCascadeCuts(0),
94 fCFContCascadeMCgen(0)
101 //________________________________________________________________________
102 AliAnalysisTaskQAMultistrange::AliAnalysisTaskQAMultistrange(const char *name)
103 : AliAnalysisTaskSE (name),
105 fAnalysisType ("ESD"),
106 fCollidingSystem ("PbPb"),
108 fkSDDSelectionOn (kTRUE),
109 fkQualityCutZprimVtxPos (kTRUE),
110 fkQualityCutNoTPConlyPrimVtx(kTRUE),
111 fkQualityCutTPCrefit (kTRUE),
112 fkQualityCutnTPCcls (kTRUE),
113 fkQualityCutPileup (kTRUE),
121 fMinPtCutOnDaughterTracks (0),
122 fEtaCutOnDaughterTracks (0),
124 fCFContCascadeCuts(0),
125 fCFContCascadeMCgen(0)
129 // Output slot #0 writes into a TList container (Cascade)
130 DefineOutput(1, AliCFContainer::Class());
131 DefineOutput(2, AliCFContainer::Class());
133 AliLog::SetClassDebugLevel("AliAnalysisTaskQAMultistrange",1);
137 AliAnalysisTaskQAMultistrange::~AliAnalysisTaskQAMultistrange() {
141 // For all TH1, 2, 3 HnSparse and CFContainer are in the fListCascade TList.
142 // They will be deleted when fListCascade is deleted by the TSelector dtor
143 // Because of TList::SetOwner() ...
144 if (fCFContCascadeCuts && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fCFContCascadeCuts; fCFContCascadeCuts = 0x0; }
145 if (fCFContCascadeMCgen && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fCFContCascadeMCgen; fCFContCascadeMCgen = 0x0; }
150 //_____________________________________________________________
151 void AliAnalysisTaskQAMultistrange::UserCreateOutputObjects() {
155 //----------------------------------------------------------------------------
156 // Option for AliLog: to suppress the extensive info prompted by a run with MC
157 //----------------------------------------------------------------------------
158 if (fisMC) AliLog::SetGlobalLogLevel(AliLog::kError);
160 //-----------------------------------------------
161 // Particle Identification Setup (new PID object)
162 //-----------------------------------------------
163 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
164 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
165 fPIDResponse = inputHandler->GetPIDResponse();
167 //---------------------------------------------------
168 // Define the container for the topological variables
169 //---------------------------------------------------
170 if(! fCFContCascadeCuts) {
171 // Container meant to store all the relevant distributions corresponding to the cut variables.
172 // NB: overflow/underflow of variables on which we want to cut later should be 0!!!
173 const Int_t lNbSteps = 4 ;
174 const Int_t lNbVariables = 20;
175 //Array for the number of bins in each dimension :
176 Int_t lNbBinsPerVar[lNbVariables] = {0};
177 lNbBinsPerVar[0] = 25; //DcaCascDaughters : [0.0,2.4,3.0] -> Rec.Cut = 2.0;
178 lNbBinsPerVar[1] = 25; //DcaBachToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.01;
179 lNbBinsPerVar[2] = 30; //CascCosineOfPointingAngle : [0.97,1.0] -> Rec.Cut = 0.98;
180 lNbBinsPerVar[3] = 40; //CascRadius : [0.0,3.9,1000.0] -> Rec.Cut = 0.2;
181 lNbBinsPerVar[4] = 30; //InvMassLambdaAsCascDghter : [1.1,1.3] -> Rec.Cut = 0.008;
182 lNbBinsPerVar[5] = 20; //DcaV0Daughters : [0.0,2.0] -> Rec.Cut = 1.5;
183 lNbBinsPerVar[6] = 201; //V0CosineOfPointingAngleToPV : [0.89,1.0] -> Rec.Cut = 0.9;
184 lNbBinsPerVar[7] = 40; //V0Radius : [0.0,3.9,1000.0] -> Rec.Cut = 0.2;
185 lNbBinsPerVar[8] = 40; //DcaV0ToPrimVertex : [0.0,0.39,110.0] -> Rec.Cut = 0.01;
186 lNbBinsPerVar[9] = 25; //DcaPosToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.05;
187 lNbBinsPerVar[10] = 25; //DcaNegToPrimVertex : [0.0,0.24,100.0] -> Rec.Cut = 0.05
188 lNbBinsPerVar[11] = 150; //InvMassXi : 2-MeV/c2 bins
189 lNbBinsPerVar[12] = 120; //InvMassOmega : 2-MeV/c2 bins
190 lNbBinsPerVar[13] = 100; //XiTransvMom : [0.0,10.0]
191 lNbBinsPerVar[14] = 110; //Y(Xi) : 0.02 in rapidity units
192 lNbBinsPerVar[15] = 110; //Y(Omega) : 0.02 in rapidity units
193 lNbBinsPerVar[16] = 112; //Proper lenght of cascade
194 lNbBinsPerVar[17] = 112; //Proper lenght of V0
195 lNbBinsPerVar[18] = 201; //V0CosineOfPointingAngleToXiV
196 lNbBinsPerVar[19] = 11; //Centrality
197 //define the container
198 fCFContCascadeCuts = new AliCFContainer("fCFContCascadeCuts","Container for Cascade cuts", lNbSteps, lNbVariables, lNbBinsPerVar );
199 //Setting the bin limits
201 Double_t *lBinLim0 = new Double_t[ lNbBinsPerVar[0] + 1 ];
202 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;
203 lBinLim0[ lNbBinsPerVar[0] ] = 3.0;
204 fCFContCascadeCuts -> SetBinLimits(0, lBinLim0);
206 //1 - DcaToPrimVertexXi
207 Double_t *lBinLim1 = new Double_t[ lNbBinsPerVar[1] + 1 ];
208 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;
209 lBinLim1[ lNbBinsPerVar[1] ] = 100.0;
210 fCFContCascadeCuts -> SetBinLimits(1, lBinLim1);
212 //2 - CascCosineOfPointingAngle
213 fCFContCascadeCuts->SetBinLimits(2, 0.97, 1.);
215 Double_t *lBinLim3 = new Double_t[ lNbBinsPerVar[3]+1 ];
216 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 ;
217 lBinLim3[ lNbBinsPerVar[3] ] = 1000.0;
218 fCFContCascadeCuts -> SetBinLimits(3, lBinLim3 );
220 //4 - InvMassLambdaAsCascDghter
221 fCFContCascadeCuts->SetBinLimits(4, 1.1, 1.13);
223 fCFContCascadeCuts -> SetBinLimits(5, 0., 2.);
224 //6 - V0CosineOfPointingAngleToPV
225 fCFContCascadeCuts -> SetBinLimits(6, 0.8, 1.001);
227 Double_t *lBinLim7 = new Double_t[ lNbBinsPerVar[7] + 1];
228 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;
229 lBinLim7[ lNbBinsPerVar[7] ] = 1000.0;
230 fCFContCascadeCuts -> SetBinLimits(7, lBinLim7);
232 //8 - DcaV0ToPrimVertex
233 Double_t *lBinLim8 = new Double_t[ lNbBinsPerVar[8]+1 ];
234 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 ;
235 lBinLim8[ lNbBinsPerVar[8] ] = 100.0;
236 fCFContCascadeCuts -> SetBinLimits(8, lBinLim8 );
238 //9 - DcaPosToPrimVertex
239 Double_t *lBinLim9 = new Double_t[ lNbBinsPerVar[9]+1 ];
240 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 ;
241 lBinLim9[ lNbBinsPerVar[9] ] = 100.0;
242 fCFContCascadeCuts -> SetBinLimits(9, lBinLim9 );
244 //10 - DcaNegToPrimVertex
245 Double_t *lBinLim10 = new Double_t[ lNbBinsPerVar[10]+1 ];
246 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 ;
247 lBinLim10[ lNbBinsPerVar[10] ] = 100.0;
248 fCFContCascadeCuts -> SetBinLimits(10, lBinLim10 );
251 fCFContCascadeCuts->SetBinLimits(11, 1.25, 1.40);
253 fCFContCascadeCuts->SetBinLimits(12, 1.62, 1.74);
255 fCFContCascadeCuts->SetBinLimits(13, 0.0, 10.0);
257 fCFContCascadeCuts->SetBinLimits(14, -1.1, 1.1);
259 fCFContCascadeCuts->SetBinLimits(15, -1.1, 1.1);
260 //16 - Proper time of cascade
261 Double_t *lBinLim16 = new Double_t[ lNbBinsPerVar[16]+1 ];
262 for(Int_t i=0; i< lNbBinsPerVar[16];i++) lBinLim16[i] = (Double_t) -1. + (110. + 1.0 ) / (lNbBinsPerVar[16] - 1) * (Double_t) i;
263 lBinLim16[ lNbBinsPerVar[16] ] = 2000.0;
264 fCFContCascadeCuts->SetBinLimits(16, lBinLim16);
265 //17 - Proper time of V0
266 fCFContCascadeCuts->SetBinLimits(17, lBinLim16);
267 //18 - V0CosineOfPointingAngleToXiV
268 fCFContCascadeCuts -> SetBinLimits(18, 0.8, 1.001);
270 Double_t *lBinLim19 = new Double_t[ lNbBinsPerVar[19]+1 ];
271 for(Int_t i=3; i< lNbBinsPerVar[19]+1;i++) lBinLim19[i] = (Double_t)(i-1)*10.;
275 fCFContCascadeCuts->SetBinLimits(19, lBinLim19 );
277 // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
278 fCFContCascadeCuts->SetStepTitle(0, "#Xi^{-} candidates");
279 fCFContCascadeCuts->SetStepTitle(1, "#bar{#Xi}^{+} candidates");
280 fCFContCascadeCuts->SetStepTitle(2, "#Omega^{-} candidates");
281 fCFContCascadeCuts->SetStepTitle(3, "#bar{#Omega}^{+} candidates");
282 // Setting the variable title, per axis
283 fCFContCascadeCuts->SetVarTitle(0, "Dca(cascade daughters) (cm)");
284 fCFContCascadeCuts->SetVarTitle(1, "ImpactParamToPV(bachelor) (cm)");
285 fCFContCascadeCuts->SetVarTitle(2, "cos(cascade PA)");
286 fCFContCascadeCuts->SetVarTitle(3, "R_{2d}(cascade decay) (cm)");
287 fCFContCascadeCuts->SetVarTitle(4, "M_{#Lambda}(as casc dghter) (GeV/c^{2})");
288 fCFContCascadeCuts->SetVarTitle(5, "Dca(V0 daughters) in Xi (cm)");
289 fCFContCascadeCuts->SetVarTitle(6, "cos(V0 PA) in cascade to PV");
290 fCFContCascadeCuts->SetVarTitle(7, "R_{2d}(V0 decay) (cm)");
291 fCFContCascadeCuts->SetVarTitle(8, "ImpactParamToPV(V0) (cm)");
292 fCFContCascadeCuts->SetVarTitle(9, "ImpactParamToPV(Pos) (cm)");
293 fCFContCascadeCuts->SetVarTitle(10, "ImpactParamToPV(Neg) (cm)");
294 fCFContCascadeCuts->SetVarTitle(11, "Inv. Mass(Xi) (GeV/c^{2})");
295 fCFContCascadeCuts->SetVarTitle(12, "Inv. Mass(Omega) (GeV/c^{2})");
296 fCFContCascadeCuts->SetVarTitle(13, "pt(cascade) (GeV/c)");
297 fCFContCascadeCuts->SetVarTitle(14, "Y(Xi)");
298 fCFContCascadeCuts->SetVarTitle(15, "Y(Omega)");
299 fCFContCascadeCuts->SetVarTitle(16, "mL/p (cascade) (cm)");
300 fCFContCascadeCuts->SetVarTitle(17, "mL/p (V0) (cm)");
301 fCFContCascadeCuts->SetVarTitle(18, "cos(V0 PA) in cascade to XiV");
302 fCFContCascadeCuts->SetVarTitle(19, "Centrality");
306 //-----------------------------------------------
307 // Define the Container for the MC generated info
308 //-----------------------------------------------
309 if(! fCFContCascadeMCgen) {
310 // Container meant to store general distributions for MC generated particles.
311 // NB: overflow/underflow of variables on which we want to cut later should be 0!!!
312 const Int_t lNbStepsMC = 4 ;
313 const Int_t lNbVariablesMC = 7;
314 //Array for the number of bins in each dimension :
315 Int_t lNbBinsPerVarMC[lNbVariablesMC] = {0};
316 lNbBinsPerVarMC[0] = 100; //Total momentum : [0.0,10.0]
317 lNbBinsPerVarMC[1] = 100; //Transverse momentum : [0.0,10.0]
318 lNbBinsPerVarMC[2] = 110; //Y : [-1.1,1.1]
319 lNbBinsPerVarMC[3] = 200; //eta : [-10, 10]
320 lNbBinsPerVarMC[4] = 200; //theta : [-10, 190]
321 lNbBinsPerVarMC[5] = 360; //Phi : [0., 360.]
322 lNbBinsPerVarMC[6] = 11; //Centrality
323 //define the container
324 fCFContCascadeMCgen = new AliCFContainer("fCFContCascadeMCgen","Container for MC gen cascade ", lNbStepsMC, lNbVariablesMC, lNbBinsPerVarMC );
325 //Setting the bin limits
327 fCFContCascadeMCgen->SetBinLimits(0, 0.0, 10.0);
328 //1 - Transverse Momentum
329 fCFContCascadeMCgen->SetBinLimits(1, 0.0, 10.0);
331 fCFContCascadeMCgen->SetBinLimits(2, -1.1, 1.1);
333 fCFContCascadeMCgen->SetBinLimits(3, -10, 10);
335 fCFContCascadeMCgen->SetBinLimits(4, -10, 190);
337 fCFContCascadeMCgen->SetBinLimits(5, 0.0, 360.0);
339 Double_t *lBinLim6MC = new Double_t[ lNbBinsPerVarMC[6]+1 ];
340 for(Int_t i=3; i< lNbBinsPerVarMC[6]+1;i++) lBinLim6MC[i] = (Double_t)(i-1)*10.;
343 lBinLim6MC[2] = 10.0;
344 fCFContCascadeMCgen->SetBinLimits(6, lBinLim6MC);
345 delete [] lBinLim6MC;
346 // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
347 fCFContCascadeMCgen->SetStepTitle(0, "#Xi^{-} candidates");
348 fCFContCascadeMCgen->SetStepTitle(1, "#bar{#Xi}^{+} candidates");
349 fCFContCascadeMCgen->SetStepTitle(2, "#Omega^{-} candidates");
350 fCFContCascadeMCgen->SetStepTitle(3, "#bar{#Omega}^{+} candidates");
351 // Setting the variable title, per axis
352 fCFContCascadeMCgen->SetVarTitle(0, "MC gen p_tot (GeV/c)");
353 fCFContCascadeMCgen->SetVarTitle(1, "MC gen p_T (GeV/c)");
354 fCFContCascadeMCgen->SetVarTitle(2, "MC gen Rapidity");
355 fCFContCascadeMCgen->SetVarTitle(3, "MC gen Pseudo-rapidity");
356 fCFContCascadeMCgen->SetVarTitle(4, "MC gen Theta");
357 fCFContCascadeMCgen->SetVarTitle(5, "MC gen Phi");
358 fCFContCascadeMCgen->SetVarTitle(6, "MC gen Centrality");
362 PostData(1, fCFContCascadeCuts);
363 PostData(2, fCFContCascadeMCgen);
365 }// end UserCreateOutputObjects
368 //________________________________________________________________________
369 void AliAnalysisTaskQAMultistrange::UserExec(Option_t *) {
372 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
373 // Main loop (called for each event)
374 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
376 //------------------------
377 // Define ESD/AOD handlers
378 //------------------------
379 AliESDEvent *lESDevent = 0x0;
380 AliAODEvent *lAODevent = 0x0;
381 AliMCEvent *lMCevent = 0x0; //for MC
382 AliStack *lMCstack = 0x0; //for MC
383 TClonesArray *arrayMC = 0; //for MC
389 //----------------------
390 // Check the PIDresponse
391 //----------------------
393 AliError("Cannot get pid response");
397 //---------------------
398 // Check the InputEvent
399 //---------------------
400 if (fAnalysisType == "ESD") {
401 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
403 AliWarning("ERROR: lESDevent not available \n");
407 lMCevent = MCEvent();
409 Printf("ERROR: Could not retrieve MC event \n");
410 cout << "Name of the file with pb :" << CurrentFileName() << endl;
413 lMCstack = lMCevent->Stack();
415 Printf("ERROR: Could not retrieve MC stack \n");
416 cout << "Name of the file with pb :" << CurrentFileName() << endl;
420 } else if (fAnalysisType == "AOD") {
421 lAODevent = dynamic_cast<AliAODEvent*>( InputEvent() );
423 AliWarning("ERROR: lAODevent not available \n");
427 arrayMC = (TClonesArray*) lAODevent->GetList()->FindObject(AliAODMCParticle::StdBranchName());
428 if (!arrayMC) AliFatal("Error: MC particles branch not found!\n");
431 Printf("Analysis type (ESD or AOD) not specified \n");
440 //----------------------------------------------------------
441 // Centrality/Multiplicity selection for PbPb/pPb collisions
442 //----------------------------------------------------------
444 AliCentrality* centrality = 0;
445 if (fAnalysisType == "ESD" && (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb")) centrality = lESDevent->GetCentrality();
446 else if (fAnalysisType == "AOD" && (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb")) centrality = lAODevent->GetCentrality();
447 Float_t lcentrality = 0.;
448 if (fCollidingSystem == "PbPb" || fCollidingSystem == "pPb") {
449 if (fkUseCleaning) lcentrality = centrality->GetCentralityPercentile(fCentrEstimator.Data());
451 lcentrality = centrality->GetCentralityPercentileUnchecked(fCentrEstimator.Data());
452 if (centrality->GetQuality()>1) {
453 PostData(1, fCFContCascadeCuts);
457 //if (lcentrality<fCentrLowLim||lcentrality>=fCentrUpLim) {
458 // PostData(1, fCFContCascadeCuts);
461 } else if (fCollidingSystem == "pp") lcentrality = 0.;
463 //---------------------
464 // SDD status selection
465 //---------------------
466 if (fkSDDSelectionOn && fAnalysisType == "ESD") {
467 TString trcl = lESDevent->GetFiredTriggerClasses();
468 if (fwithSDD) { if(!(trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
469 else if (!fwithSDD){ if((trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
470 } else if (fkSDDSelectionOn && fAnalysisType == "AOD") {
471 TString trcl = lAODevent->GetFiredTriggerClasses();
472 if (fwithSDD) { if(!(trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
473 else if (!fwithSDD) { if((trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
479 UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
480 Bool_t isSelected = 0;
481 if (fCollidingSystem == "pp" || fCollidingSystem == "PbPb") isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
482 else if (fCollidingSystem == "pPb") isSelected = (maskIsSelected & AliVEvent::kINT7) == AliVEvent::kINT7;
483 if (!isSelected){ PostData(1, fCFContCascadeCuts); return; }
485 //------------------------------
486 // Well-established PV selection
487 //------------------------------
488 Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
489 Double_t lMagneticField = -10.;
490 if (fCollidingSystem == "PbPb" || fCollidingSystem == "pp") {
491 if (fAnalysisType == "ESD") {
492 const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
493 const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
494 if (fkQualityCutNoTPConlyPrimVtx) {
495 const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
496 if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
497 AliWarning(" No SPD prim. vertex nor prim. Tracking vertex ... return !");
498 PostData(1, fCFContCascadeCuts);
502 lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
503 } else if (fAnalysisType == "AOD") {
504 const AliAODVertex *lPrimaryBestAODVtx = lAODevent->GetPrimaryVertex();
505 if (!lPrimaryBestAODVtx){
506 AliWarning("No prim. vertex in AOD... return!");
507 PostData(1, fCFContCascadeCuts);
510 lPrimaryBestAODVtx->GetXYZ( lBestPrimaryVtxPos );
512 } else if (fCollidingSystem == "pPb") {
513 if (fAnalysisType == "ESD") {
514 Bool_t fHasVertex = kFALSE;
515 const AliESDVertex *vertex = lESDevent->GetPrimaryVertexTracks();
516 if (vertex->GetNContributors() < 1) {
517 vertex = lESDevent->GetPrimaryVertexSPD();
518 if (vertex->GetNContributors() < 1) fHasVertex = kFALSE;
519 else fHasVertex = kTRUE;
520 TString vtxTyp = vertex->GetTitle();
522 vertex->GetCovarianceMatrix(cov);
523 Double_t zRes = TMath::Sqrt(cov[5]);
524 if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fHasVertex = kFALSE;
526 else fHasVertex = kTRUE;
527 if (fHasVertex == kFALSE) { //Is First event in chunk rejection: Still present!
528 AliWarning("Pb / | PV does not satisfy selection criteria!");
529 PostData(1, fCFContCascadeCuts);
532 vertex->GetXYZ( lBestPrimaryVtxPos );
533 } else if (fAnalysisType == "AOD") {
534 Bool_t fHasVertex = kFALSE;
535 const AliAODVertex *vertex = lAODevent->GetPrimaryVertex();
536 if (vertex->GetNContributors() < 1) {
537 vertex = lAODevent->GetPrimaryVertexSPD();
538 if (vertex->GetNContributors() < 1) fHasVertex = kFALSE;
539 else fHasVertex = kTRUE;
540 TString vtxTyp = vertex->GetTitle();
542 vertex->GetCovarianceMatrix(cov);
543 Double_t zRes = TMath::Sqrt(cov[5]);
544 if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fHasVertex = kFALSE;
546 else fHasVertex = kTRUE;
547 if (fHasVertex == kFALSE) { //Is First event in chunk rejection: Still present!
548 AliWarning("Pb / | PV does not satisfy selection criteria!");
549 PostData(1, fCFContCascadeCuts);
552 vertex->GetXYZ( lBestPrimaryVtxPos );
556 if (fAnalysisType == "ESD") lMagneticField = lESDevent->GetMagneticField();
557 else if (fAnalysisType == "AOD") lMagneticField = lAODevent->GetMagneticField();
559 //----------------------------
560 // Vertex Z position selection
561 //----------------------------
562 if (fkQualityCutZprimVtxPos) {
563 if (TMath::Abs(lBestPrimaryVtxPos[2]) > fVtxRange ) {
564 PostData(1, fCFContCascadeCuts);
569 //-----------------------
570 // Pilup selection for pp
571 //-----------------------
572 if (fCollidingSystem == "pp") {
573 if (fAnalysisType == "ESD") {
574 if (fkQualityCutPileup) { if(lESDevent->IsPileupFromSPD()){ PostData(1, fCFContCascadeCuts); return; } }
575 } else if (fAnalysisType == "AOD") {
576 if (fkQualityCutPileup) { if(lAODevent->IsPileupFromSPD()){ PostData(1, fCFContCascadeCuts); return; } }
581 ////////////////////////////
582 // MC GENERATED CASCADE PART
583 ////////////////////////////
589 Int_t lNbMCPrimary = 0;
590 if (fAnalysisType == "ESD") lNbMCPrimary = lMCstack->GetNprimary(); //lMCstack->GetNprimary(); or lMCstack->GetNtrack();
591 else if (fAnalysisType == "AOD") lNbMCPrimary = arrayMC->GetEntries();
593 for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++) {
596 Double_t partPt = 0.;
597 Double_t partEta = 0.;
598 Double_t partTheta = 0.;
599 Double_t partPhi = 0.;
600 Double_t partRap = 0.;
601 Double_t partEnergy = 0.; //for Rapidity
602 Double_t partPz = 0.; //for Rapidity
605 if ( fAnalysisType == "ESD" ) {
606 TParticle* lCurrentParticlePrimary = 0x0;
607 lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
608 if (!lCurrentParticlePrimary) {
609 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
612 if (!lMCstack->IsPhysicalPrimary(iCurrentLabelStack)) continue;
613 TParticle* cascMC = 0x0;
614 cascMC = lCurrentParticlePrimary;
616 Printf("MC TParticle pointer to Cascade = 0x0 ! Skip ...");
620 partPt = cascMC->Pt();
621 partEta = cascMC->Eta();
622 partTheta = cascMC->Theta()*180.0/TMath::Pi();
623 partPhi = cascMC->Phi()*180.0/TMath::Pi();
624 partEnergy = cascMC->Energy();
625 partPz = cascMC->Pz();
626 PDGcode = lCurrentParticlePrimary->GetPdgCode();
627 } else if ( fAnalysisType == "AOD" ) {
628 AliAODMCParticle *lCurrentParticleaod = (AliAODMCParticle*) arrayMC->At(iCurrentLabelStack);
629 if (!lCurrentParticleaod) {
630 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
633 if (!lCurrentParticleaod->IsPhysicalPrimary()) continue;
634 partP = lCurrentParticleaod->P();
635 partPt = lCurrentParticleaod->Pt();
636 partEta = lCurrentParticleaod->Eta();
637 partTheta = lCurrentParticleaod->Theta()*180.0/TMath::Pi();
638 partPhi = lCurrentParticleaod->Phi()*180.0/TMath::Pi();
639 partEnergy = lCurrentParticleaod->E();
640 partPz = lCurrentParticleaod->Pz();
641 PDGcode = lCurrentParticleaod->GetPdgCode();
643 partRap = 0.5*TMath::Log((partEnergy + partPz) / (partEnergy - partPz + 1.e-13));
645 // Fill the AliCFContainer
646 Double_t lContainerCutVarsMC[7] = {0.0};
647 lContainerCutVarsMC[0] = partP;
648 lContainerCutVarsMC[1] = partPt;
649 lContainerCutVarsMC[2] = partRap;
650 lContainerCutVarsMC[3] = partEta;
651 lContainerCutVarsMC[4] = partTheta;
652 lContainerCutVarsMC[5] = partPhi;
653 lContainerCutVarsMC[6] = lcentrality;
654 if (PDGcode == 3312) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,0); // for Xi-
655 if (PDGcode == -3312) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,1); // for Xi+
656 if (PDGcode == 3334) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,2); // for Omega-
657 if (PDGcode == -3334) fCFContCascadeMCgen->Fill(lContainerCutVarsMC,3); // for Omega+
663 //////////////////////////////
664 // CASCADE RECONSTRUCTION PART
665 //////////////////////////////
670 if (fAnalysisType == "ESD") ncascades = lESDevent->GetNumberOfCascades();
671 else if (fAnalysisType == "AOD") ncascades = lAODevent->GetNumberOfCascades();
673 for (Int_t iXi = 0; iXi < ncascades; iXi++) {// This is the begining of the Cascade loop (ESD or AOD)
675 // -------------------------------------
676 // - Initialisation of the local variables that will be needed for ESD/AOD
677 // -- Container variables (1st round)
678 Double_t lDcaXiDaughters = -1. ; //[Container]
679 Double_t lXiCosineOfPointingAngle = -1. ; //[Container]
680 Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 }; //Useful to define other variables: radius fid. vol., ctau, etc. for cascade
681 Double_t lXiRadius = -1000. ; //[Container]
682 UShort_t lPosTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
683 UShort_t lNegTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
684 UShort_t lBachTPCClusters = -1; //To check the quality of the tracks. For ESD only ...
685 Double_t lInvMassLambdaAsCascDghter = 0.; //[Container]
686 Double_t lDcaV0DaughtersXi = -1.; //[Container]
687 Double_t lDcaBachToPrimVertexXi = -1.; //[Container]
688 Double_t lDcaV0ToPrimVertexXi = -1.; //[Container]
689 Double_t lDcaPosToPrimVertexXi = -1.; //[Container]
690 Double_t lDcaNegToPrimVertexXi = -1.; //[Container]
691 Double_t lV0CosineOfPointingAngle = -1.; //[Container]
692 Double_t lV0toXiCosineOfPointingAngle = -1.; //[Container]
693 Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; //Useful to define other variables: radius fid. vol., ctau, etc. for VO
694 Double_t lV0RadiusXi = -1000.0; //[Container]
695 Double_t lV0quality = 0.; // ??
696 Double_t lInvMassXiMinus = 0.; //[Container]
697 Double_t lInvMassXiPlus = 0.; //[Container]
698 Double_t lInvMassOmegaMinus = 0.; //[Container]
699 Double_t lInvMassOmegaPlus = 0.; //[Container]
701 Bool_t lIsBachelorKaonForTPC = kFALSE;
702 Bool_t lIsBachelorPionForTPC = kFALSE;
703 Bool_t lIsNegPionForTPC = kFALSE;
704 Bool_t lIsPosPionForTPC = kFALSE;
705 Bool_t lIsNegProtonForTPC = kFALSE;
706 Bool_t lIsPosProtonForTPC = kFALSE;
707 // -- More container variables and quality checks
708 Double_t lXiMomX = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
709 Double_t lXiMomY = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
710 Double_t lXiMomZ = 0.; //Useful to define other variables: lXiTransvMom, lXiTotMom
711 Double_t lXiTransvMom = 0.; //[Container]
712 Double_t lXiTotMom = 0.; //Useful to define other variables: cTau
713 Double_t lV0PMomX = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
714 Double_t lV0PMomY = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
715 Double_t lV0PMomZ = 0.; //Useful to define other variables: lV0TotMom, lpTrackTransvMom
716 Double_t lV0NMomX = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
717 Double_t lV0NMomY = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
718 Double_t lV0NMomZ = 0.; //Useful to define other variables: lV0TotMom, lnTrackTransvMom
719 Double_t lV0TotMom = 0.; //Useful to define other variables: lctauV0
720 Double_t lBachMomX = 0.; //Useful to define other variables: lBachTransvMom
721 Double_t lBachMomY = 0.; //Useful to define other variables: lBachTransvMom
722 Double_t lBachMomZ = 0.; //Useful to define other variables: lBachTransvMom
723 Double_t lBachTransvMom = 0.; //Selection on the min bachelor pT
724 Double_t lpTrackTransvMom = 0.; //Selection on the min bachelor pT
725 Double_t lnTrackTransvMom = 0.; //Selection on the min bachelor pT
726 Short_t lChargeXi = -2; //Useful to select the particles based on the charge
727 Double_t lRapXi = -20.0; //[Container]
728 Double_t lRapOmega = -20.0; //[Container]
729 Float_t etaBach = 0.; //Selection on the eta range
730 Float_t etaPos = 0.; //Selection on the eta range
731 Float_t etaNeg = 0.; //Selection on the eta range
732 // -- variables for the AliCFContainer dedicated to cascade cut optmisiation: ESD and AOD
733 if (fAnalysisType == "ESD") {
735 // -------------------------------------
736 // - Load the cascades from the handler
737 AliESDcascade *xi = lESDevent->GetCascade(iXi);
740 // ---------------------------------------------------------------------------
741 // - Assigning the necessary variables for specific AliESDcascade data members
743 xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
744 lDcaXiDaughters = xi->GetDcaXiDaughters();
745 lXiCosineOfPointingAngle = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
746 // Take care : the best available vertex should be used (like in AliCascadeVertexer)
747 xi->GetXYZcascade( lPosXi[0], lPosXi[1], lPosXi[2] );
748 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
750 // -------------------------------------------------------------------------------------------------------------------------------
751 // - Around the tracks : Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
752 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( xi->GetPindex() );
753 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( xi->GetNindex() );
754 UInt_t lBachIdx = (UInt_t) TMath::Abs( xi->GetBindex() );
755 // Care track label can be negative in MC production (linked with the track quality)
756 // However = normally, not the case for track index ...
757 // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
758 if (lBachIdx == lIdxNegXi) continue;
759 if (lBachIdx == lIdxPosXi) continue;
760 // - Get the track for the daughters
761 AliESDtrack *pTrackXi = lESDevent->GetTrack( lIdxPosXi );
762 AliESDtrack *nTrackXi = lESDevent->GetTrack( lIdxNegXi );
763 AliESDtrack *bachTrackXi = lESDevent->GetTrack( lBachIdx );
764 if (!pTrackXi || !nTrackXi || !bachTrackXi ) continue;
765 // - Get the TPCnumber of cluster for the daughters
766 lPosTPCClusters = pTrackXi->GetTPCNcls();
767 lNegTPCClusters = nTrackXi->GetTPCNcls();
768 lBachTPCClusters = bachTrackXi->GetTPCNcls();
770 // ------------------------------------
771 // - Rejection of a poor quality tracks
772 if (fkQualityCutTPCrefit) {
773 // 1 - Poor quality related to TPCrefit
774 ULong_t pStatus = pTrackXi->GetStatus();
775 ULong_t nStatus = nTrackXi->GetStatus();
776 ULong_t bachStatus = bachTrackXi->GetStatus();
777 if ((pStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
778 if ((nStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
779 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" Bach. track has no TPCrefit ... continue!"); continue; }
781 if (fkQualityCutnTPCcls) {
782 // 2 - Poor quality related to TPC clusters
783 if (lPosTPCClusters < fMinnTPCcls) { AliWarning(" V0 Pos. track has less than minn TPC clusters ... continue!"); continue; }
784 if (lNegTPCClusters < fMinnTPCcls) { AliWarning(" V0 Neg. track has less than minn TPC clusters ... continue!"); continue; }
785 if (lBachTPCClusters < fMinnTPCcls) { AliWarning(" Bach. track has less than minn TPC clusters ... continue!"); continue; }
788 // ------------------------------
789 etaPos = pTrackXi->Eta();
790 etaNeg = nTrackXi->Eta();
791 etaBach = bachTrackXi->Eta();
792 lInvMassLambdaAsCascDghter = xi->GetEffMass();
793 lDcaV0DaughtersXi = xi->GetDcaV0Daughters();
794 lV0CosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
795 lDcaV0ToPrimVertexXi = xi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
796 lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField) );
797 xi->GetXYZ( lPosV0Xi[0], lPosV0Xi[1], lPosV0Xi[2] );
798 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
799 lDcaPosToPrimVertexXi = TMath::Abs( pTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) );
800 lDcaNegToPrimVertexXi = TMath::Abs( nTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) );
802 //----------------------------------------------------------------------------------------------------
803 // - Around effective masses. Change mass hypotheses to cover all the possibilities: Xi-/+, Omega -/+
804 if (bachTrackXi->Charge() < 0) {
805 // Calculate the effective mass of the Xi- candidate. pdg code 3312 = Xi-
807 xi->ChangeMassHypothesis(lV0quality , 3312);
808 lInvMassXiMinus = xi->GetEffMassXi();
809 // Calculate the effective mass of the Xi- candidate. pdg code 3334 = Omega-
811 xi->ChangeMassHypothesis(lV0quality , 3334);
812 lInvMassOmegaMinus = xi->GetEffMassXi();
813 // Back to default hyp.
815 xi->ChangeMassHypothesis(lV0quality , 3312);
816 }// end if negative bachelor
817 if ( bachTrackXi->Charge() > 0 ) {
818 // Calculate the effective mass of the Xi+ candidate. pdg code -3312 = Xi+
820 xi->ChangeMassHypothesis(lV0quality , -3312);
821 lInvMassXiPlus = xi->GetEffMassXi();
822 // Calculate the effective mass of the Xi+ candidate. pdg code -3334 = Omega+
824 xi->ChangeMassHypothesis(lV0quality , -3334);
825 lInvMassOmegaPlus = xi->GetEffMassXi();
826 // Back to "default" hyp.
828 xi->ChangeMassHypothesis(lV0quality , -3312);
829 }// end if positive bachelor
831 // ----------------------------------------------
832 // - TPC PID : 3-sigma bands on Bethe-Bloch curve
834 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
835 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
836 // Negative V0 daughter
837 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion )) < 4) lIsNegPionForTPC = kTRUE;
838 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
839 // Positive V0 daughter
840 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion )) < 4) lIsPosPionForTPC = kTRUE;
841 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
843 // ------------------------------
844 // - Miscellaneous pieces of info that may help regarding data quality assessment.
845 xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
846 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
847 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
848 xi->GetNPxPyPz(lV0NMomX,lV0NMomY,lV0NMomZ);
849 xi->GetPPxPyPz(lV0PMomX,lV0PMomY,lV0PMomZ);
850 lV0TotMom = TMath::Sqrt(TMath::Power(lV0NMomX+lV0PMomX,2)+TMath::Power(lV0NMomY+lV0PMomY,2)+TMath::Power(lV0NMomZ+lV0PMomZ,2));
851 xi->GetBPxPyPz( lBachMomX, lBachMomY, lBachMomZ );
852 lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
853 lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX + lV0NMomY*lV0NMomY );
854 lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX + lV0PMomY*lV0PMomY );
855 lChargeXi = xi->Charge();
856 lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
857 lRapXi = xi->RapXi();
858 lRapOmega = xi->RapOmega();
860 } else if (fAnalysisType == "AOD") {
862 // -------------------------------------
863 // - Load the cascades from the handler
864 const AliAODcascade *xi = lAODevent->GetCascade(iXi);
867 //----------------------------------------------------------------------------
868 // - Assigning the necessary variables for specific AliESDcascade data members
869 lDcaXiDaughters = xi->DcaXiDaughters();
870 lXiCosineOfPointingAngle = xi->CosPointingAngleXi( lBestPrimaryVtxPos[0],
871 lBestPrimaryVtxPos[1],
872 lBestPrimaryVtxPos[2] );
873 lPosXi[0] = xi->DecayVertexXiX();
874 lPosXi[1] = xi->DecayVertexXiY();
875 lPosXi[2] = xi->DecayVertexXiZ();
876 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
878 //-------------------------------------------------------------------------------------------------------------------------------
879 // - Around the tracks: Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
880 AliAODTrack *pTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );
881 AliAODTrack *nTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );
882 AliAODTrack *bachTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );
883 if (!pTrackXi || !nTrackXi || !bachTrackXi ) continue;
884 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( pTrackXi->GetID() );
885 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( nTrackXi->GetID() );
886 UInt_t lBachIdx = (UInt_t) TMath::Abs( bachTrackXi->GetID() );
887 // Care track label can be negative in MC production (linked with the track quality)
888 // However = normally, not the case for track index ...
890 // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
891 if (lBachIdx == lIdxNegXi) continue;
892 if (lBachIdx == lIdxPosXi) continue;
893 // - Get the TPCnumber of cluster for the daughters
894 lPosTPCClusters = pTrackXi->GetTPCNcls();
895 lNegTPCClusters = nTrackXi->GetTPCNcls();
896 lBachTPCClusters = bachTrackXi->GetTPCNcls();
898 // ------------------------------------
899 // - Rejection of a poor quality tracks
900 if (fkQualityCutTPCrefit) {
901 // - Poor quality related to TPCrefit
902 if (!(pTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
903 if (!(nTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
904 if (!(bachTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" Bach. track has no TPCrefit ... continue!"); continue; }
906 if (fkQualityCutnTPCcls) {
907 // - Poor quality related to TPC clusters
908 if (lPosTPCClusters < fMinnTPCcls) continue;
909 if (lNegTPCClusters < fMinnTPCcls) continue;
910 if (lBachTPCClusters < fMinnTPCcls) continue;
913 // ------------------------------------------------------------------------------------------------------------------------------
914 // - Around the tracks: Bach + V0 (AOD). Necessary variables for AODcascade data members coming from the AODv0 part (inheritance)
915 etaPos = pTrackXi->Eta();
916 etaNeg = nTrackXi->Eta();
917 etaBach = bachTrackXi->Eta();
918 lChargeXi = xi->ChargeXi();
919 if ( lChargeXi < 0) lInvMassLambdaAsCascDghter = xi->MassLambda();
920 else lInvMassLambdaAsCascDghter = xi->MassAntiLambda();
921 lDcaV0DaughtersXi = xi->DcaV0Daughters();
922 lDcaV0ToPrimVertexXi = xi->DcaV0ToPrimVertex();
923 lDcaBachToPrimVertexXi = xi->DcaBachToPrimVertex();
924 lPosV0Xi[0] = xi->DecayVertexV0X();
925 lPosV0Xi[1] = xi->DecayVertexV0Y();
926 lPosV0Xi[2] = xi->DecayVertexV0Z();
927 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
928 lV0CosineOfPointingAngle = xi->CosPointingAngle( lBestPrimaryVtxPos );
929 lDcaPosToPrimVertexXi = xi->DcaPosToPrimVertex();
930 lDcaNegToPrimVertexXi = xi->DcaNegToPrimVertex();
932 // ---------------------------------------------------------------------------------------------------
933 // - Around effective masses. Change mass hypotheses to cover all the possibilities: Xi-/+, Omega -/+
934 if ( lChargeXi < 0 ) lInvMassXiMinus = xi->MassXi();
935 if ( lChargeXi > 0 ) lInvMassXiPlus = xi->MassXi();
936 if ( lChargeXi < 0 ) lInvMassOmegaMinus = xi->MassOmega();
937 if ( lChargeXi > 0 ) lInvMassOmegaPlus = xi->MassOmega();
939 // ----------------------------------------------
940 // - TPC PID : 3-sigma bands on Bethe-Bloch curve
942 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
943 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
944 // Negative V0 daughter
945 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion )) < 4) lIsNegPionForTPC = kTRUE;
946 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
947 // Positive V0 daughter
948 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion )) < 4) lIsPosPionForTPC = kTRUE;
949 if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
951 //---------------------------------
952 // - Extra info for QA (AOD)
953 // Miscellaneous pieces of info that may help regarding data quality assessment.
954 // Cascade transverse and total momentum
955 lXiMomX = xi->MomXiX();
956 lXiMomY = xi->MomXiY();
957 lXiMomZ = xi->MomXiZ();
958 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
959 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
960 Double_t lV0MomX = xi->MomV0X();
961 Double_t lV0MomY = xi->MomV0Y();
962 Double_t lV0MomZ = xi->MomV0Z();
963 lV0TotMom = TMath::Sqrt(TMath::Power(lV0MomX,2)+TMath::Power(lV0MomY,2)+TMath::Power(lV0MomZ,2));
964 lBachMomX = xi->MomBachX();
965 lBachMomY = xi->MomBachY();
966 lBachMomZ = xi->MomBachZ();
967 lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
968 lV0NMomX = xi->MomNegX();
969 lV0NMomY = xi->MomNegY();
970 lV0PMomX = xi->MomPosX();
971 lV0PMomY = xi->MomPosY();
972 lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX + lV0NMomY*lV0NMomY );
973 lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX + lV0PMomY*lV0PMomY );
974 lV0toXiCosineOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );
975 lRapXi = xi->RapXi();
976 lRapOmega = xi->RapOmega();
978 }// end of AOD treatment
980 //---------------------------------------
981 // Cut on pt of the three daughter tracks
982 if (lBachTransvMom<fMinPtCutOnDaughterTracks) continue;
983 if (lpTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
984 if (lnTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
986 //---------------------------------------------------
987 // Cut on pseudorapidity of the three daughter tracks
988 if (TMath::Abs(etaBach)>fEtaCutOnDaughterTracks) continue;
989 if (TMath::Abs(etaPos)>fEtaCutOnDaughterTracks) continue;
990 if (TMath::Abs(etaNeg)>fEtaCutOnDaughterTracks) continue;
992 //----------------------------------
993 // Calculate proper time for cascade
994 Double_t cascadeMass = 0.;
995 if ( ( (lChargeXi<0) && lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC ) ||
996 ( (lChargeXi>0) && lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC ) ) cascadeMass = 1.321;
997 if ( ( (lChargeXi<0) && lIsBachelorKaonForTPC && lIsPosProtonForTPC && lIsNegPionForTPC ) ||
998 ( (lChargeXi>0) && lIsBachelorKaonForTPC && lIsNegProtonForTPC && lIsPosPionForTPC ) ) cascadeMass = 1.672;
999 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));
1000 if (lXiTotMom!=0) lctau = lctau*cascadeMass/lXiTotMom;
1002 // Calculate proper time for Lambda (reconstructed)
1003 Float_t lambdaMass = 1.115683; // PDG mass
1004 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));
1005 Float_t lctauV0 = -1.;
1006 if (lV0TotMom!=0) lctauV0 = distV0Xi*lambdaMass/lV0TotMom;
1009 // Fill the AliCFContainer (optimisation of topological selections)
1010 Double_t lContainerCutVars[21] = {0.0};
1011 lContainerCutVars[0] = lDcaXiDaughters;
1012 lContainerCutVars[1] = lDcaBachToPrimVertexXi;
1013 lContainerCutVars[2] = lXiCosineOfPointingAngle;
1014 lContainerCutVars[3] = lXiRadius;
1015 lContainerCutVars[4] = lInvMassLambdaAsCascDghter;
1016 lContainerCutVars[5] = lDcaV0DaughtersXi;
1017 lContainerCutVars[6] = lV0CosineOfPointingAngle;
1018 lContainerCutVars[7] = lV0RadiusXi;
1019 lContainerCutVars[8] = lDcaV0ToPrimVertexXi;
1020 lContainerCutVars[9] = lDcaPosToPrimVertexXi;
1021 lContainerCutVars[10] = lDcaNegToPrimVertexXi;
1022 lContainerCutVars[13] = lXiTransvMom;
1023 lContainerCutVars[16] = lctau;
1024 lContainerCutVars[17] = lctauV0;
1025 lContainerCutVars[18] = lV0toXiCosineOfPointingAngle;
1026 lContainerCutVars[19] = lcentrality;
1027 if ( lChargeXi < 0 ) {
1028 lContainerCutVars[11] = lInvMassXiMinus;
1029 lContainerCutVars[12] = lInvMassOmegaMinus;
1030 lContainerCutVars[14] = lRapXi;
1031 lContainerCutVars[15] = -1.;
1032 if (lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,0); // for Xi-
1033 lContainerCutVars[11] = lInvMassXiMinus;
1034 lContainerCutVars[12] = lInvMassOmegaMinus;
1035 lContainerCutVars[14] = -1.;
1036 lContainerCutVars[15] = lRapOmega;
1037 if (lIsBachelorKaonForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,2); // for Omega-
1039 lContainerCutVars[11] = lInvMassXiPlus;
1040 lContainerCutVars[12] = lInvMassOmegaPlus;
1041 lContainerCutVars[14] = lRapXi;
1042 lContainerCutVars[15] = -1.;
1043 if (lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,1); // for Xi+
1044 lContainerCutVars[11] = lInvMassXiPlus;
1045 lContainerCutVars[12] = lInvMassOmegaPlus;
1046 lContainerCutVars[14] = -1.;
1047 lContainerCutVars[15] = lRapOmega;
1048 if (lIsBachelorKaonForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,3); // for Omega+
1052 }// end of the Cascade loop (ESD or AOD)
1055 // Post output data.
1056 PostData(1, fCFContCascadeCuts);
1057 PostData(2, fCFContCascadeMCgen);
1060 //________________________________________________________________________
1061 void AliAnalysisTaskQAMultistrange::Terminate(Option_t *) {