AliPhysics  a0db429 (a0db429)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros
AliAnalysisTaskSELc2pK0sfromAODtracks.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appeuear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id$ */
17 
18 //
19 //
20 // Lc->pK0s analysis code
21 //
22 // Input: AOD
23 // Output: TTree or THnSparse (mass vs pT vs Centrality)
24 //
25 // Cuts:
26 // TTree: very loose cut
27 // THnSparse: One THnSparse is created per cut. One cut is specified by
28 // an array of bits, each bit corresponds to a cut in "Cut" function.
29 // Use "AddCutStream" function to add a cut.
30 //
31 //-------------------------------------------------------------------------
32 //
33 // Authors: Y.S Watanabe(a)
34 // (a) CNS, the University of Tokyo
35 // Contatcs: wyosuke@cns.s.u-tokyo.ac.jp
36 //-------------------------------------------------------------------------
37 
38 #include <TSystem.h>
39 #include <TParticle.h>
40 #include <TParticlePDG.h>
41 #include <TH1F.h>
42 #include <TH1F.h>
43 #include <TH2F.h>
44 #include <TH3F.h>
45 #include <TRandom.h>
46 #include <THnSparse.h>
47 #include <TLorentzVector.h>
48 #include <TTree.h>
49 #include "TROOT.h"
50 #include <TDatabasePDG.h>
51 #include <AliAnalysisDataSlot.h>
52 #include <AliAnalysisDataContainer.h>
53 #include "AliStack.h"
54 #include "AliMCEvent.h"
55 #include "AliAnalysisManager.h"
56 #include "AliAODMCHeader.h"
57 #include "AliAODHandler.h"
58 #include "AliLog.h"
59 #include "AliExternalTrackParam.h"
60 #include "AliAODVertex.h"
61 #include "AliESDVertex.h"
62 #include "AliAODRecoDecay.h"
63 #include "AliAODRecoDecayHF.h"
64 #include "AliAODRecoCascadeHF.h"
65 #include "AliESDtrack.h"
66 #include "AliAODTrack.h"
67 #include "AliAODv0.h"
68 #include "AliAODcascade.h"
69 #include "AliAODMCParticle.h"
70 #include "AliAnalysisTaskSE.h"
72 #include "AliPIDResponse.h"
73 #include "AliPIDCombined.h"
74 #include "AliTOFPIDResponse.h"
75 #include "AliAODPidHF.h"
76 #include "AliInputEventHandler.h"
77 #include "AliESDtrackCuts.h"
78 #include "AliNeutralTrackParam.h"
79 #include "AliKFParticle.h"
80 #include "AliKFVertex.h"
81 #include "AliExternalTrackParam.h"
82 #include "AliESDtrack.h"
83 #include "AliCentrality.h"
84 #include "AliVertexerTracks.h"
86 #include <set>
87 
88 using std::cout;
89 using std::endl;
90 
94 
95 //__________________________________________________________________________
97  AliAnalysisTaskSE(),
98  fUseMCInfo(kFALSE),
99  fOutput(0),
100  fOutputAll(0),
101  fListCuts(0),
102  fCEvents(0),
103  fHTrigger(0),
104  fHCentrality(0),
105  fHReactionPlane(0),
106  fAnalCuts(0),
107  fIsEventSelected(kFALSE),
108  fWriteVariableTree(kFALSE),
109  fWriteEachVariableTree(kFALSE),
110  fWriteMCVariableTree(kFALSE),
111  fVariablesTree(0),
112  fProtonVariablesTree(0),
113  fV0VariablesTree(0),
114  fMCVariablesTree(0),
115  fMCProtonVariablesTree(0),
116  fMCV0VariablesTree(0),
117  fReconstructPrimVert(kFALSE),
118  fIsMB(kFALSE),
119  fIsSemi(kFALSE),
120  fIsCent(kFALSE),
121  fIsINT7(kFALSE),
122  fIsEMC7(kFALSE),
123  fCandidateVariables(),
124  fCandidateProtonVariables(),
125  fCandidateV0Variables(),
126  fCandidateMCVariables(),
127  fCandidateMCProtonVariables(),
128  fCandidateMCV0Variables(),
129  fVtx1(0),
130  fV1(0),
131  fVtxZ(0),
132  fBzkG(0),
133  fCentrality(0),
134  fReactionPlane(0),
135  fRunNumber(0),
136  fTriggerCheck(0),
137  fUseCentralityV0M(kFALSE),
138  fEvNumberCounter(0),
139  fCounter(0),
140  fHistonEvtvsRunNumber(0),
141  fHistonProtonvsRunNumber(0),
142  fHistonK0svsRunNumber(0),
143  fHistoLcMCGen(0),
144  fHistoLcK0SpMass(0),
145  fHistoLcK0SpMassMix(0),
146  fHistoLcK0SpMassCoarse(0),
147  fHistoLcK0SpMassMixCoarse(0),
148  fHistoK0spCorrelation(0),
149  fHistoK0spCorrelationMix(0),
150  fHistoK0spCorrelationMCS(0),
151  fHistoLcK0SpMassMCS(0),
152  fHistoLcK0SpPi0MassMCS(0),
153  fHistoLcKPluspMass(0),
154  fHistoLcKMinuspMass(0),
155  fHistoLcKPluspMassMix(0),
156  fHistoLcKMinuspMassMix(0),
157  fHistoBachPt(0),
158  fHistoBachPtMCS(0),
159  fHistoBachPtMCGen(0),
160  fHistoKaonPt(0),
161  fHistoKaonPtMCS(0),
162  fHistoKaonPtMCGen(0),
163  fHistoK0sMassvsPt(0),
164  fHistoK0sMassvsPtMCS(0),
165  fHistoK0sMassvsPtMCGen(0),
166  fHistod0Bach(0),
167  fHistod0V0(0),
168  fHistod0d0(0),
169  fHistoV0CosPA(0),
170  fHistoProbProton(0),
171  fHistoDecayLength(0),
172  fHistoK0SMass(0),
173  fHistoMassTagV0Min(0),
174  fHistoMassTagV0SameSignMin(0),
175  fDoEventMixing(0),
176  fNumberOfEventsForMixing (5),
177  fNzVtxBins (0),
178  fNCentBins (0),
179  fNRPBins (0),
180  fNOfPools(1),
181  fEventBuffer(0x0),
182  fEventInfo(0x0),
183  fProtonTracks(0x0),
184  fV0Tracks(0x0),
185  fProtonCutVarsArray(0x0),
186  fV0CutVarsArray(0x0)
187 {
188  //
189  // Default Constructor.
190  //
191 }
192 
193 //___________________________________________________________________________
196  Bool_t writeVariableTree) :
197  AliAnalysisTaskSE(name),
198  fUseMCInfo(kFALSE),
199  fOutput(0),
200  fOutputAll(0),
201  fListCuts(0),
202  fCEvents(0),
203  fHTrigger(0),
204  fHCentrality(0),
205  fHReactionPlane(0),
206  fAnalCuts(analCuts),
207  fIsEventSelected(kFALSE),
208  fWriteVariableTree(writeVariableTree),
209  fWriteEachVariableTree(kFALSE),
210  fWriteMCVariableTree(kFALSE),
211  fVariablesTree(0),
212  fProtonVariablesTree(0),
213  fV0VariablesTree(0),
214  fMCVariablesTree(0),
215  fMCProtonVariablesTree(0),
216  fMCV0VariablesTree(0),
217  fReconstructPrimVert(kFALSE),
218  fIsMB(kFALSE),
219  fIsSemi(kFALSE),
220  fIsCent(kFALSE),
221  fIsINT7(kFALSE),
222  fIsEMC7(kFALSE),
223  fCandidateVariables(),
224  fCandidateProtonVariables(),
225  fCandidateV0Variables(),
226  fCandidateMCVariables(),
227  fCandidateMCProtonVariables(),
228  fCandidateMCV0Variables(),
229  fVtx1(0),
230  fV1(0),
231  fVtxZ(0),
232  fBzkG(0),
233  fCentrality(0),
234  fReactionPlane(0),
235  fRunNumber(0),
236  fTriggerCheck(0),
237  fUseCentralityV0M(kFALSE),
238  fEvNumberCounter(0),
239  fCounter(0),
240  fHistonEvtvsRunNumber(0),
241  fHistonProtonvsRunNumber(0),
242  fHistonK0svsRunNumber(0),
243  fHistoLcMCGen(0),
244  fHistoLcK0SpMass(0),
245  fHistoLcK0SpMassMix(0),
246  fHistoLcK0SpMassCoarse(0),
247  fHistoLcK0SpMassMixCoarse(0),
248  fHistoK0spCorrelation(0),
249  fHistoK0spCorrelationMix(0),
250  fHistoK0spCorrelationMCS(0),
251  fHistoLcK0SpMassMCS(0),
252  fHistoLcK0SpPi0MassMCS(0),
253  fHistoLcKPluspMass(0),
254  fHistoLcKMinuspMass(0),
255  fHistoLcKPluspMassMix(0),
256  fHistoLcKMinuspMassMix(0),
257  fHistoBachPt(0),
258  fHistoBachPtMCS(0),
259  fHistoBachPtMCGen(0),
260  fHistoKaonPt(0),
261  fHistoKaonPtMCS(0),
262  fHistoKaonPtMCGen(0),
263  fHistoK0sMassvsPt(0),
264  fHistoK0sMassvsPtMCS(0),
265  fHistoK0sMassvsPtMCGen(0),
266  fHistod0Bach(0),
267  fHistod0V0(0),
268  fHistod0d0(0),
269  fHistoV0CosPA(0),
270  fHistoProbProton(0),
271  fHistoDecayLength(0),
272  fHistoK0SMass(0),
273  fHistoMassTagV0Min(0),
274  fHistoMassTagV0SameSignMin(0),
275  fDoEventMixing(0),
276  fNumberOfEventsForMixing (5),
277  fNzVtxBins (0),
278  fNCentBins (0),
279  fNRPBins (0),
280  fNOfPools(1),
281  fEventBuffer(0x0),
282  fEventInfo(0x0),
283  fProtonTracks(0x0),
284  fV0Tracks(0x0),
285  fProtonCutVarsArray(0x0),
286  fV0CutVarsArray(0x0)
287 {
288  //
289  // Constructor. Initialization of Inputs and Outputs
290  //
291  Info("AliAnalysisTaskSELc2pK0sfromAODtracks","Calling Constructor");
292 
293  DefineOutput(1,TList::Class()); //conters
294  DefineOutput(2,TList::Class());
295  DefineOutput(3,TList::Class()); //conters
296  DefineOutput(4,TTree::Class()); //My private output
297  DefineOutput(5,TTree::Class()); //My private output
298  DefineOutput(6,TTree::Class()); //My private output
299  DefineOutput(7,TTree::Class()); //My private output
300  DefineOutput(8,AliNormalizationCounter::Class());
301  DefineOutput(9,TTree::Class()); //My private output
302  DefineOutput(10,TTree::Class()); //My private output
303 }
304 
305 //___________________________________________________________________________
307  //
308  // destructor
309  //
310  Info("~AliAnalysisTaskSELc2pK0sfromAODtracks","Calling Destructor");
311 
312  if (fOutput) {
313  delete fOutput;
314  fOutput = 0;
315  }
316 
317  if (fOutputAll) {
318  delete fOutputAll;
319  fOutputAll = 0;
320  }
321 
322  if (fListCuts) {
323  delete fListCuts;
324  fListCuts = 0;
325  }
326 
327 
328  if (fAnalCuts) {
329  delete fAnalCuts;
330  fAnalCuts = 0;
331  }
332 
333  if (fVariablesTree) {
334  delete fVariablesTree;
335  fVariablesTree = 0;
336  }
337 
338  if (fProtonVariablesTree) {
339  delete fProtonVariablesTree;
341  }
342  if (fV0VariablesTree) {
343  delete fV0VariablesTree;
344  fV0VariablesTree = 0;
345  }
346  if (fMCVariablesTree) {
347  delete fMCVariablesTree;
348  fMCVariablesTree = 0;
349  }
351  delete fMCProtonVariablesTree;
353  }
354  if (fMCV0VariablesTree) {
355  delete fMCV0VariablesTree;
356  fMCV0VariablesTree = 0;
357  }
358 
359  if(fCounter){
360  delete fCounter;
361  fCounter = 0;
362  }
363 
364  if(fProtonTracks) fProtonTracks->Delete();
365  delete fProtonTracks;
366  if(fV0Tracks) fV0Tracks->Delete();
367  delete fV0Tracks;
369  delete fProtonCutVarsArray;
370  if(fV0CutVarsArray) fV0CutVarsArray->Delete();
371  delete fV0CutVarsArray;
372  if(fEventBuffer){
373  for(Int_t i=0; i<fNOfPools; i++) delete fEventBuffer[i];
374  delete fEventBuffer;
375  }
376  delete fEventInfo;
377 
378 }
379 
380 //_________________________________________________
382  //
383  // Initialization
384  //
385  //
386 
387  fIsEventSelected=kFALSE;
388 
389  if (fDebug > 1) AliInfo("Init");
390 
391  fListCuts = new TList();
392  fListCuts->SetOwner();
393  fListCuts->SetName("ListCuts");
395  PostData(2,fListCuts);
396 
397  return;
398 }
399 
400 //_________________________________________________
402 {
403  //
404  // UserExec
405  //
406 
407  if (!fInputEvent) {
408  AliError("NO EVENT FOUND!");
409  return;
410  }
411  AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
412  fCEvents->Fill(1);
414 
415  //------------------------------------------------
416  // First check if the event has proper vertex and B
417  //------------------------------------------------
418  fBzkG = (Double_t)aodEvent->GetMagneticField();
419  AliKFParticle::SetField(fBzkG);
420  if (TMath::Abs(fBzkG)<0.001) {
421  delete fV1;
422  return;
423  }
424  fCEvents->Fill(2);
425 
428 
429  //------------------------------------------------
430  // MC analysis setting
431  //------------------------------------------------
432  TClonesArray *mcArray = 0;
433  AliAODMCHeader *mcHeader=0;
434  if (fUseMCInfo) {
435  // MC array need for maching
436  mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
437  if (!mcArray) {
438  AliError("Could not find Monte-Carlo in AOD");
439  return;
440  }
441  fCEvents->Fill(6); // in case of MC events
442 
443  // load MC header
444  mcHeader = (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName());
445  if (!mcHeader) {
446  AliError("AliAnalysisTaskSELc2pK0sfromAODtracks::UserExec: MC header branch not found!\n");
447  return;
448  }
449  fCEvents->Fill(7); // in case of MC events
450 
451  Double_t zMCVertex = mcHeader->GetVtxZ();
452  if (TMath::Abs(zMCVertex) > fAnalCuts->GetMaxVtxZ()) {
453  AliDebug(2,Form("Event rejected: abs(zVtxMC)=%f > fAnalCuts->GetMaxVtxZ()=%f",zMCVertex,fAnalCuts->GetMaxVtxZ()));
454  return;
455  } else {
456  fCEvents->Fill(17); // in case of MC events
457  }
458  if ((TMath::Abs(zMCVertex) < fAnalCuts->GetMaxVtxZ()) && (!fAnalCuts->IsEventRejectedDuePhysicsSelection()) && (!fAnalCuts->IsEventRejectedDueToTrigger())) {
459  Bool_t selevt = MakeMCAnalysis(mcArray);
460  if(!selevt) return;
461  }
462  }
463 
464  //------------------------------------------------
465  // Event selection
466  //------------------------------------------------
467  fVtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
468  if (!fVtx1) return;
469 
470  Double_t pos[3],cov[6];
471  fVtx1->GetXYZ(pos);
472  fVtx1->GetCovarianceMatrix(cov);
473  fV1 = new AliESDVertex(pos,cov,100.,100,fVtx1->GetName());
474  fVtxZ = pos[2];
475 
476  Bool_t fIsTriggerNotOK = fAnalCuts->IsEventRejectedDueToTrigger();
477  if(!fIsTriggerNotOK) fCEvents->Fill(3);
478  fIsEventSelected = fAnalCuts->IsEventSelected(aodEvent); // better to initialize before CheckEventSelection call
479  if(!fIsEventSelected) {
480  delete fV1;
481  return;
482  }
483  fCEvents->Fill(4);
484 
485  fIsMB=(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kMB)==(AliVEvent::kMB);
486  fIsSemi=(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kSemiCentral)==(AliVEvent::kSemiCentral);
487  fIsCent=(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kCentral)==(AliVEvent::kCentral);
488  fIsINT7=(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kINT7)==(AliVEvent::kINT7);
489  fIsEMC7=(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kEMC7)==(AliVEvent::kEMC7);
491  if(fIsMB) fHTrigger->Fill(1);
492  if(fIsSemi) fHTrigger->Fill(2);
493  if(fIsCent) fHTrigger->Fill(3);
494  if(fIsINT7) fHTrigger->Fill(4);
495  if(fIsEMC7) fHTrigger->Fill(5);
496  if(fIsMB|fIsSemi|fIsCent) fHTrigger->Fill(7);
497  if(fIsINT7|fIsEMC7) fHTrigger->Fill(8);
498  if(fIsMB&fIsSemi) fHTrigger->Fill(10);
499  if(fIsMB&fIsCent) fHTrigger->Fill(11);
500  if(fIsINT7&fIsEMC7) fHTrigger->Fill(12);
501 
502  if(fUseCentralityV0M){
503  AliCentrality *cent = aodEvent->GetCentrality();
504  fCentrality = cent->GetCentralityPercentile("V0M");
505  }else{
506  fCentrality = 1.;
507  }
508  if(fCentrality<0.||fCentrality>100.-0.0000001) {
509  delete fV1;
510  return;
511  }
512  fHCentrality->Fill(fCentrality);
513  fRunNumber = aodEvent->GetRunNumber();
514 
515  Int_t runnumber_offset = 0;
516  Int_t runnumber = aodEvent->GetRunNumber();
517  if(runnumber<=131000&&runnumber>=114000){
518  runnumber_offset = 114000;//lhc10bcde
519  }else if(runnumber<=196000&&runnumber>=195000){
520  runnumber_offset = 195000;//lhc13bc
521  }else if(runnumber<=170593&&runnumber>=167902){
522  runnumber_offset = 167902;//lhc11h
523  }
524  fHistonEvtvsRunNumber->Fill(runnumber-runnumber_offset,1.);
525 
526  //------------------------------------------------
527  // Check if the event has v0 candidate
528  //------------------------------------------------
529  //Int_t nv0 = aodEvent->GetNumberOfV0s();
530  fCEvents->Fill(5);
531 
532  //------------------------------------------------
533  // Main analysis done in this function
534  //------------------------------------------------
535  MakeAnalysis(aodEvent,mcArray);
536 
537 
538  PostData(1,fOutput);
539  PostData(3,fOutputAll);
540  PostData(4,fVariablesTree);
541  PostData(5,fProtonVariablesTree);
542  PostData(6,fV0VariablesTree);
543  PostData(7,fMCVariablesTree);
544  PostData(8,fCounter);
545  PostData(9,fMCProtonVariablesTree);
546  PostData(10,fMCV0VariablesTree);
547 
548  fIsEventSelected=kFALSE;
549 
550  delete fV1;
551  return;
552 }
553 
554 //________________________________________ terminate ___________________________
556 {
557  // The Terminate() function is the last function to be called during
558  // a query. It always runs on the client, it can be used to present
559  // the results graphically or save the results to file.
560 
561  //AliInfo("Terminate","");
562  AliAnalysisTaskSE::Terminate();
563 
564  fOutput = dynamic_cast<TList*> (GetOutputData(1));
565  if (!fOutput) {
566  AliError("fOutput not available");
567  return;
568  }
569 
570  fOutputAll = dynamic_cast<TList*> (GetOutputData(3));
571  if (!fOutputAll) {
572  AliError("fOutputAll not available");
573  return;
574  }
575 
576  return;
577 }
578 
579 //___________________________________________________________________________
581 {
582  //
583  // UserCreateOutputObject
584  //
585  //AliInfo(Form("CreateOutputObjects of task %s\n", GetName()));
586 
587  //------------------------------------------------
588  // output object setting
589  //------------------------------------------------
590  fOutput = new TList();
591  fOutput->SetOwner();
592  fOutput->SetName("chist0");
593  DefineGeneralHistograms(); // define general histograms
594  PostData(1,fOutput);
595 
596  fOutputAll = new TList();
597  fOutputAll->SetOwner();
598  fOutputAll->SetName("anahisto");
599  DefineAnalysisHistograms(); // define general histograms
600  PostData(3,fOutputAll);
601 
603  PostData(4,fVariablesTree);
604 
606  PostData(5,fProtonVariablesTree);
607 
609  PostData(6,fV0VariablesTree);
610 
612  PostData(7,fMCVariablesTree);
613 
615  PostData(9,fMCProtonVariablesTree);
616 
618  PostData(10,fMCV0VariablesTree);
619 
620  //Counter for Normalization
621  TString normName="NormalizationCounter";
622  AliAnalysisDataContainer *cont = GetOutputSlot(8)->GetContainer();
623  if(cont)normName=(TString)cont->GetName();
624  fCounter = new AliNormalizationCounter(normName.Data());
625  fCounter->Init();
626  PostData(8,fCounter);
627 
628  if(fDoEventMixing){
629  fProtonTracks = new TObjArray();
630  fProtonTracks->SetOwner();
631  fV0Tracks = new TObjArray();
632  fV0Tracks->SetOwner();
633  fProtonCutVarsArray = new TObjArray();
634  fProtonCutVarsArray->SetOwner();
635  fV0CutVarsArray = new TObjArray();
636  fV0CutVarsArray->SetOwner();
637 
639  fEventBuffer = new TTree*[fNOfPools];
640  for(Int_t i=0; i<fNOfPools; i++){
641  fEventBuffer[i]=new TTree(Form("EventBuffer_%d",i), "Temporary buffer for event mixing");
642  fEventBuffer[i]->Branch("zVertex", &fVtxZ);
643  fEventBuffer[i]->Branch("centrality", &fCentrality);
644  fEventBuffer[i]->Branch("reactionplane", &fReactionPlane);
645  fEventBuffer[i]->Branch("eventInfo", "TObjString",&fEventInfo);
646  fEventBuffer[i]->Branch("v1array", "TObjArray", &fV0Tracks);
647  fEventBuffer[i]->Branch("v1varsarray", "TObjArray", &fV0CutVarsArray);
648  }
649  }
650 
651 
652  return;
653 }
654 
655 //-------------------------------------------------------------------------------
657 (
658  AliAODEvent *aodEvent, TClonesArray *mcArray
659  )
660 {
661  //
662  // Main Analysis part
663  //
664  if(fDoEventMixing){
665  if(fProtonTracks) fProtonTracks->Delete();
666  if(fV0Tracks) fV0Tracks->Delete();
667  if(fProtonCutVarsArray) fProtonCutVarsArray->Delete();
668  if(fV0CutVarsArray) fV0CutVarsArray->Delete();
669  }
670 
671  Int_t nV0s= aodEvent->GetNumberOfV0s();
672  Int_t nTracks= aodEvent->GetNumberOfTracks();
673 
674  Int_t seleTrkFlags[nTracks];
675  Int_t nSeleTrks=0;
676  SelectTrack(aodEvent,nTracks,nSeleTrks,seleTrkFlags,mcArray);
677 
678  Bool_t seleV0Flags[nV0s];
679  Int_t nSeleV0=0;
680  SelectV0(aodEvent,nV0s,nSeleV0,seleV0Flags,mcArray);
681 
682  //------------------------------------------------
683  // V0 loop
684  //------------------------------------------------
685  for (Int_t iv0 = 0; iv0<nV0s; iv0++) {
686  if(!seleV0Flags[iv0]) continue;
687  AliAODv0 *v0 = aodEvent->GetV0(iv0);
688  if(!v0) continue;
689 
690  AliAODTrack *cptrack = (AliAODTrack*)(v0->GetDaughter(0));
691  AliAODTrack *cntrack = (AliAODTrack*)(v0->GetDaughter(1));
692 
693  //------------------------------------------------
694  // track loop
695  //------------------------------------------------
696  for (Int_t itrk = 0; itrk<nTracks; itrk++) {
697  if(seleTrkFlags[itrk]!=1) continue;
698  AliAODTrack *trk = (AliAODTrack*)aodEvent->GetTrack(itrk);
699  if(trk->GetID()<0) continue;
700 
701  Int_t cpid = cptrack->GetID();
702  Int_t cnid = cntrack->GetID();
703  Int_t lpid = trk->GetID();
704  if((cpid==lpid)||(cnid==lpid)) continue;
705 
706  if(!fAnalCuts->SelectWithRoughCuts(v0,trk)) continue;
707 
708  AliAODVertex *secVert = ReconstructSecondaryVertex(v0,trk,aodEvent);
709  if(!secVert) continue;
710 
711  AliAODRecoCascadeHF *lcobj = MakeCascadeHF(v0,trk,aodEvent,secVert);
712  if(!lcobj) {
713  continue;
714  }
715 
716  FillROOTObjects(lcobj,v0,trk,aodEvent,mcArray);
717 
718  lcobj->GetSecondaryVtx()->RemoveDaughters();
719  lcobj->UnsetOwnPrimaryVtx();
720  delete lcobj;lcobj=NULL;
721  delete secVert;
722  }
723  }
724 
725  if(fDoEventMixing){
726  fEventInfo->SetString(Form("Ev%d_esd%d_E%d_V%d",AliAnalysisManager::GetAnalysisManager()->GetNcalls(),((AliAODHeader*)aodEvent->GetHeader())->GetEventNumberESDFile(),fProtonTracks->GetEntries(),fV0Tracks->GetEntries()));
727  Int_t ind=GetPoolIndex(fVtxZ,fCentrality,fReactionPlane);
728  if(ind>=0 && ind<fNOfPools){
729  if(fEventBuffer[ind]->GetEntries() >= fNumberOfEventsForMixing){
730  DoEventMixingWithPools(ind);
731  if(fEventBuffer[ind]->GetEntries() >= 20*fNumberOfEventsForMixing){
732  ResetPool(ind);
733  }
734  }
735  fEventBuffer[ind]->Fill();
736  }
737  }
738 
739 }
740 
742 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillROOTObjects(AliAODRecoCascadeHF *lcobj, AliAODv0 *v0, AliAODTrack *trk, AliAODEvent *aodEvent, TClonesArray *mcArray)
743 {
744  //
745  // Fill histograms or tree depending on fWriteVariableTree
746  //
747  if(!trk) return;
748  if(!v0) return;
749 
750  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
751  Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
752  Double_t mlamPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
753 
754  for(Int_t i=0;i<54;i++){
755  fCandidateVariables[i] = -9999.;
756  }
757 
758  Double_t pxp = trk->Px();
759  Double_t pyp = trk->Py();
760  Double_t pzp = trk->Pz();
761  Double_t momp = sqrt(pxp*pxp+pyp*pyp+pzp*pzp);
762  Double_t Ep = sqrt(momp*momp+mprPDG*mprPDG);
763 
764  Double_t pxv = v0->Px();
765  Double_t pyv = v0->Py();
766  Double_t pzv = v0->Pz();
767  Double_t momv = sqrt(pxv*pxv+pyv*pyv+pzv*pzv);
768  Double_t mv = v0->MassK0Short();
769  Double_t Ev = sqrt(momv*momv+mv*mv);
770 
771  Double_t cosoa = (pxp*pxv+pyp*pyv+pzp*pzv)/momp/momv;
772 
773  fCandidateVariables[ 0] = lcobj->InvMassLctoK0sP();
774  fCandidateVariables[ 1] = lcobj->Px();
775  fCandidateVariables[ 2] = lcobj->Py();
776  fCandidateVariables[ 3] = lcobj->Pz();
777  fCandidateVariables[ 4] = v0->MassK0Short();
778  fCandidateVariables[ 5] = lcobj->PxProng(0);
779  fCandidateVariables[ 6] = lcobj->PyProng(0);
780  fCandidateVariables[ 7] = lcobj->PzProng(0);
781  fCandidateVariables[ 8] = lcobj->PxProng(1);
782  fCandidateVariables[ 9] = lcobj->PyProng(1);
783  fCandidateVariables[10] = lcobj->PzProng(1);
784  fCandidateVariables[11] = fVtx1->GetX();
785  fCandidateVariables[12] = fVtx1->GetY();
786  fCandidateVariables[13] = fVtx1->GetZ();
788  fCandidateVariables[15] = lcobj->DecayLengthXY();
789  fCandidateVariables[16] = (Float_t) fAnalCuts->CalculateLcCosPAXY(lcobj);
790 
791  Double_t nSigmaTPCpr=-9999.;
792  Double_t nSigmaTOFpr=-9999.;
793  Double_t probProton=-9999.;
794  if(fAnalCuts->GetIsUsePID())
795  {
796  nSigmaTPCpr = fAnalCuts->GetPidHF()->GetPidResponse()->NumberOfSigmasTPC(trk,AliPID::kProton);
797  nSigmaTOFpr = fAnalCuts->GetPidHF()->GetPidResponse()->NumberOfSigmasTOF(trk,AliPID::kProton);
799  probProton = fAnalCuts->GetProtonProbabilityTPCTOF(trk);
800  }
801  fCandidateVariables[17] = nSigmaTPCpr;
802  fCandidateVariables[18] = nSigmaTOFpr;
803  fCandidateVariables[19] = probProton;
804  }
805  fCandidateVariables[20] = 0;
806  fCandidateVariables[21] = lcobj->Getd0Prong(0);
807  fCandidateVariables[22] = lcobj->Getd0Prong(1);
808 
809  Double_t pospx = v0->MomPosX();
810  Double_t pospy = v0->MomPosY();
811  Double_t pospz = v0->MomPosZ();
812  Double_t pose = sqrt(pospx*pospx+pospy*pospy+pospz*pospz+0.000511*0.000511);
813  Double_t negpx = v0->MomNegX();
814  Double_t negpy = v0->MomNegY();
815  Double_t negpz = v0->MomNegZ();
816  Double_t nege = sqrt(negpx*negpx+negpy*negpy+negpz*negpz+0.000511*0.000511);
817  Double_t massPhoton = sqrt(pow(pose+nege,2)-pow(pospx+negpx,2)-pow(pospy+negpy,2)-pow(pospz+negpz,2));
818  fCandidateVariables[23] = v0->MassLambda();
819  fCandidateVariables[24] = v0->MassAntiLambda();
820  fCandidateVariables[25] = massPhoton;
821  fCandidateVariables[26] = v0->DcaPosToPrimVertex();
822  fCandidateVariables[27] = v0->DcaNegToPrimVertex();
823  fCandidateVariables[28] = v0->DcaV0ToPrimVertex();
824  fCandidateVariables[29] = v0->DcaV0Daughters();
825 
826  AliAODMCParticle *mclc = 0;
827  AliAODMCParticle *mcpr = 0;
828  AliAODMCParticle *mcv0 = 0;
829  Int_t mclablc = 0;
830  Int_t mcpdgpr_array[100];
831  Int_t mcpdgv0_array[100];
832  Int_t mclabelpr_array[100];
833  Int_t mclabelv0_array[100];
834  Int_t mcngen_pr=-9999;
835  Int_t mcngen_v0=-9999;
836  Bool_t islambdac = kFALSE;
837  Bool_t islambdac3body = kFALSE;
838  Double_t decayvertx_mc = -9999.;
839  Double_t decayverty_mc = -9999.;
840  Double_t decayvertz_mc = -9999.;
841  Double_t genvertx_mc = -9999.;
842  Double_t genverty_mc = -9999.;
843  Double_t genvertz_mc = -9999.;
844  if(fUseMCInfo && mcArray){
845  mclablc = MatchToMC(lcobj,mcArray,mcpdgpr_array, mcpdgv0_array,mclabelpr_array,mclabelv0_array,mcngen_pr,mcngen_v0);
846  if(mclablc>-1){
847  mclc = (AliAODMCParticle*) mcArray->At(mclablc);
848  if(mclabelpr_array[0]>=0)
849  mcpr = (AliAODMCParticle*) mcArray->At(mclabelpr_array[0]);
850  if(mclabelv0_array[0]>=0)
851  mcv0 = (AliAODMCParticle*) mcArray->At(mclabelv0_array[0]);
852  }
853  if(mclc){
854  Int_t pdgcode = mclc->GetPdgCode();
855  if(abs(pdgcode)==4122 && abs(mcpdgpr_array[1])==4122 && abs(mcpdgv0_array[1])==311 && abs(mcpdgv0_array[2])==4122 && mclc->GetNDaughters()==2){
856  islambdac = kTRUE;
857  }
858  if(abs(pdgcode)==4122 && abs(mcpdgpr_array[1])==2214 && abs(mcpdgpr_array[2])==4122 && abs(mcpdgv0_array[1])==311 && abs(mcpdgv0_array[2])==4122 && mclc->GetNDaughters()==2){
859  islambdac3body = kTRUE;
860  }
861  if(abs(pdgcode)==4122 && abs(mcpdgpr_array[1])==3222 && abs(mcpdgpr_array[2])==4122 && abs(mcpdgv0_array[1])==311 && abs(mcpdgv0_array[2])==4122 && mclc->GetNDaughters()==2){
862  islambdac3body = kTRUE;
863  }
864  if(abs(pdgcode)==4122 && abs(mcpdgpr_array[1])==3222 && abs(mcpdgpr_array[2])==3224 && abs(mcpdgpr_array[3])==4122 && abs(mcpdgv0_array[1])==311 && abs(mcpdgv0_array[2])==4122 && mclc->GetNDaughters()==2){
865  islambdac3body = kTRUE;
866  }
867  if(abs(pdgcode)==4122 && abs(mcpdgpr_array[1])==4122 && abs(mcpdgv0_array[1])==311 && abs(mcpdgv0_array[2])==313 && abs(mcpdgv0_array[3])==4122 && mclc->GetNDaughters()==2){
868  islambdac3body = kTRUE;
869  }
870  genvertx_mc = mclc->Xv();
871  genverty_mc = mclc->Yv();
872  genvertz_mc = mclc->Zv();
873  }
874  if(mcpr){
875  decayvertx_mc = mcpr->Xv();
876  decayverty_mc = mcpr->Yv();
877  decayvertz_mc = mcpr->Zv();
878  }
879  }
880  fCandidateVariables[30] = (Float_t) islambdac;
881 
882  Double_t LcPx = lcobj->Px();
883  Double_t LcPy = lcobj->Py();
884  Double_t LcPt = TMath::Sqrt(LcPx*LcPx+LcPy*LcPy);
885 
886  Double_t d0z0[2],covd0z0[3];
887  trk->PropagateToDCA(fVtx1,fBzkG,kVeryBig,d0z0,covd0z0);
888  Double_t x0 = fVtx1->GetX();
889  Double_t y0 = fVtx1->GetY();
890  Double_t px0 = LcPx/LcPt;
891  Double_t py0 = LcPy/LcPt;
892  Double_t tx[3];
893  trk->GetXYZ(tx);
894  Double_t x1 = tx[0];
895  Double_t y1 = tx[1];
896  trk->GetPxPyPz(tx);
897  Double_t px1 = tx[0];
898  Double_t py1 = tx[1];
899  Double_t pt1 = sqrt(px1*px1+py1*py1);
900  px1 = px1/pt1;
901  py1 = py1/pt1;
902 
903  Double_t dx = x0 - x1;
904  Double_t dy = y0 - y1;
905 
906  Double_t Delta = -px0*py1+py0*px1;
907  Double_t a0 = -9999.;
908  if(Delta!=0)
909  {
910  a0 = 1./Delta * (py1 * dx - px1 * dy);
911  }
912  Double_t neovertx = x0 + a0 * px0;
913  Double_t neoverty = y0 + a0 * py0;
914  Double_t z0 = fVtx1->GetZ();
915  Double_t neovertz = z0 + TMath::Abs(a0)*trk->Pz()/trk->Pt();
916 
917  fCandidateVariables[31] = neovertx;
918  fCandidateVariables[32] = neoverty;
919  fCandidateVariables[33] = neovertz;
920  fCandidateVariables[34] = (Float_t) islambdac3body;
921  fCandidateVariables[35] = decayvertx_mc;
922  fCandidateVariables[36] = decayverty_mc;
923  fCandidateVariables[37] = decayvertz_mc;
924  if(mclc){
925  fCandidateVariables[38] = mclc->Px();
926  fCandidateVariables[39] = mclc->Py();
927  fCandidateVariables[40] = mclc->Pz();
928  }
929  fCandidateVariables[41] = genvertx_mc;
930  fCandidateVariables[42] = genverty_mc;
931  fCandidateVariables[43] = genvertz_mc;
932 
933  Double_t x1_k0s = v0->DecayVertexV0X();
934  Double_t y1_k0s = v0->DecayVertexV0Y();
935  Double_t px1_k0s = v0->Px();
936  Double_t py1_k0s = v0->Py();
937  Double_t pt1_k0s = sqrt(px1_k0s*px1_k0s+py1_k0s*py1_k0s);
938  px1_k0s = px1_k0s/pt1_k0s;
939  py1_k0s = py1_k0s/pt1_k0s;
940 
941  Double_t dx_k0s = x0 - x1_k0s;
942  Double_t dy_k0s = y0 - y1_k0s;
943 
944  Double_t Delta_k0s = -px0*py1_k0s+py0*px1_k0s;
945  Double_t a0_k0s = -9999.;
946  if(Delta_k0s!=0)
947  {
948  a0_k0s = 1./Delta_k0s * (py1_k0s * dx_k0s - px1_k0s * dy_k0s);
949  }
950  Double_t neovertx_k0s = x0 + a0_k0s * px0;
951  Double_t neoverty_k0s = y0 + a0_k0s * py0;
952  fCandidateVariables[44] = neovertx_k0s;
953  fCandidateVariables[45] = neoverty_k0s;
954  fCandidateVariables[46] = lcobj->GetDCA();
955  if(mclc){
956  fCandidateVariables[47] = mclc->GetPdgCode();
957  fCandidateVariables[48] = mcpdgpr_array[0];
958  fCandidateVariables[49] = mcpdgpr_array[1];
959  fCandidateVariables[50] = mcpdgpr_array[2];
960  fCandidateVariables[51] = mcpdgv0_array[0];
961  fCandidateVariables[52] = mcpdgv0_array[1];
962  fCandidateVariables[53] = mcpdgv0_array[2];
963  }
964 
965 
967  fVariablesTree->Fill();
968 
969  Double_t cont_correlation[5];
970  Double_t phi_trig = trk->Phi();
971  Double_t phi_assoc = v0->Phi();
972  Double_t eta_trig = trk->Eta();
973  Double_t eta_assoc = v0->Eta();
974  Double_t deltaphi = phi_assoc-phi_trig;
975  if(deltaphi<-M_PI/2.) deltaphi += 2 * M_PI;
976  if(deltaphi>3*M_PI/2.) deltaphi -= 2 * M_PI;
977  cont_correlation[0] = deltaphi;
978  cont_correlation[1] = eta_assoc-eta_trig;
979  cont_correlation[2] = trk->Pt();
980  cont_correlation[3] = v0->Pt();
981  cont_correlation[4] = fCentrality;
982  fHistoK0spCorrelation->Fill(cont_correlation);
983 
984  Double_t cont[3];
985  cont[0] = lcobj->InvMassLctoK0sP();
986  cont[1] = lcobj->Pt();
987  cont[2] = fCentrality;
988  if(cosoa>0.) fHistoLcK0SpMassCoarse->Fill(cont);
989 
991  {
992  fHistoLcK0SpMass->Fill(cont);
993 
994  fHistod0Bach->Fill(lcobj->Getd0Prong(0));
995  fHistod0V0->Fill(lcobj->Getd0Prong(1));
996  fHistod0d0->Fill(lcobj->Getd0Prong(0)*lcobj->Getd0Prong(1));
997  fHistoV0CosPA->Fill(lcobj->CosV0PointingAngle());
998  fHistoProbProton->Fill(probProton);
1000  fHistoK0SMass->Fill(v0->MassK0Short());
1001 
1002  if(fUseMCInfo){
1003  if(islambdac){
1004  fHistoLcK0SpMassMCS->Fill(cont);
1005  fHistoK0spCorrelationMCS->Fill(cont_correlation);
1006  }
1007  }
1008  }
1009  if(fUseMCInfo){
1010  if(islambdac3body){
1011  fHistoLcK0SpPi0MassMCS->Fill(cont);
1012  }
1013  }
1014 
1015  return;
1016 }
1018 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillMixROOTObjects(TLorentzVector *trkp, TLorentzVector *v0, TVector *prvars, TVector *v0vars)
1019 {
1023  if(!trkp) return;
1024  if(!v0) return;
1025 
1026  for(Int_t i=0;i<54;i++){
1027  fCandidateVariables[i] = -9999.;
1028  }
1029 
1030  Double_t pxp = trkp->Px();
1031  Double_t pyp = trkp->Py();
1032  Double_t pzp = trkp->Pz();
1033  Double_t momp = sqrt(pxp*pxp+pyp*pyp+pzp*pzp);
1034  Double_t Ep = sqrt(momp*momp+0.938272*0.938272);
1035 
1036  Double_t pxv = v0->Px();
1037  Double_t pyv = v0->Py();
1038  Double_t pzv = v0->Pz();
1039  Double_t momv = sqrt(pxv*pxv+pyv*pyv+pzv*pzv);
1040  Double_t mv = v0->M();
1041  Double_t Ev = sqrt(momv*momv+mv*mv);
1042 
1043  Double_t cosoa = (pxp*pxv+pyp*pyv+pzp*pzv)/momp/momv;
1044 
1045  Double_t pxsum = pxp + pxv;
1046  Double_t pysum = pyp + pyv;
1047  Double_t pzsum = pzp + pzv;
1048  Double_t Esum = Ep + Ev;
1049  Double_t mpk0s = sqrt(Esum*Esum-pxsum*pxsum-pysum*pysum-pzsum*pzsum);
1050 
1051 
1052  Double_t cont_correlation[5];
1053  Double_t phi_trig = trkp->Phi();
1054  Double_t phi_assoc = v0->Phi();
1055  Double_t eta_trig = trkp->Eta();
1056  Double_t eta_assoc = v0->Eta();
1057  Double_t deltaphi = phi_assoc-phi_trig;
1058  if(deltaphi<-M_PI/2.) deltaphi += 2 * M_PI;
1059  if(deltaphi>3*M_PI/2.) deltaphi -= 2 * M_PI;
1060  cont_correlation[0] = deltaphi;
1061  cont_correlation[1] = eta_assoc-eta_trig;
1062  cont_correlation[2] = trkp->Pt();
1063  cont_correlation[3] = v0->Pt();
1064  cont_correlation[4] = fCentrality;
1065  fHistoK0spCorrelationMix->Fill(cont_correlation);
1066 
1067  fCandidateVariables[ 0] = mpk0s;
1068  fCandidateVariables[ 1] = pxsum;
1069  fCandidateVariables[ 2] = pysum;
1070  fCandidateVariables[ 3] = pzsum;
1071  fCandidateVariables[ 4] = v0->M();
1072  fCandidateVariables[ 5] = trkp->Px();
1073  fCandidateVariables[ 6] = trkp->Py();
1074  fCandidateVariables[ 7] = trkp->Pz();
1075  fCandidateVariables[ 8] = v0->Px();
1076  fCandidateVariables[ 9] = v0->Py();
1077  fCandidateVariables[10] = v0->Pz();
1078  fCandidateVariables[20] = 1;
1079 
1080  Double_t LcPx = pxsum;
1081  Double_t LcPy = pysum;
1082  Double_t LcPt = TMath::Sqrt(LcPx*LcPx+LcPy*LcPy);
1083 
1084  Double_t x0 = (*prvars)[5];
1085  Double_t y0 = (*prvars)[6];
1086  Double_t px0 = LcPx/LcPt;
1087  Double_t py0 = LcPy/LcPt;
1088  Double_t tx[3];
1089  Double_t x1 = (*prvars)[1];
1090  Double_t y1 = (*prvars)[2];
1091  Double_t px1 = (*prvars)[3];
1092  Double_t py1 = (*prvars)[4];
1093 
1094  Double_t dx = x0 - x1;
1095  Double_t dy = y0 - y1;
1096 
1097  Double_t Delta = -px0*py1+py0*px1;
1098  Double_t a0 = -9999.;
1099  if(Delta!=0)
1100  {
1101  a0 = 1./Delta * (py1 * dx - px1 * dy);
1102  }
1103  Double_t neovertx = x0 + a0 * px0;
1104  Double_t neoverty = y0 + a0 * py0;
1105 
1106  fCandidateVariables[31] = neovertx;
1107  fCandidateVariables[32] = neoverty;
1108 
1109  if(fWriteVariableTree)
1110  fVariablesTree->Fill();
1111 
1112  Double_t rdhfcutvars[7];
1113  rdhfcutvars[0] = mpk0s;
1114  rdhfcutvars[1] = sqrt(pxsum*pxsum+pysum*pysum);
1115  rdhfcutvars[2] = (*prvars)[7];
1116  rdhfcutvars[3] = (*prvars)[8];
1117  rdhfcutvars[4] = a0;
1118  rdhfcutvars[5] = (*prvars)[0];
1119  rdhfcutvars[6] = (*v0vars)[0];
1120 
1121  if(fAnalCuts->IsSelected(trkp,v0,rdhfcutvars,AliRDHFCuts::kAll))
1122  {
1123  Double_t cont[3];
1124  cont[0] = mpk0s;
1125  cont[1] = sqrt(pxsum*pxsum+pysum*pysum);
1126  cont[2] = fCentrality;
1127  fHistoLcK0SpMassMixCoarse->Fill(cont);
1128  fHistoLcK0SpMassMix->Fill(cont);
1129  }
1130 
1131  return;
1132 }
1133 
1136 {
1137  //
1138  // Define tree variables
1139  //
1140 
1141  const char* nameoutput = GetOutputSlot(4)->GetContainer()->GetName();
1142  fVariablesTree = new TTree(nameoutput,"Candidates variables tree");
1143  Int_t nVar = 54;
1144  fCandidateVariables = new Float_t [nVar];
1145  TString * fCandidateVariableNames = new TString[nVar];
1146 
1147  fCandidateVariableNames[ 0]="InvMassLc2pK0s";
1148  fCandidateVariableNames[ 1]="LcPx";
1149  fCandidateVariableNames[ 2]="LcPy";
1150  fCandidateVariableNames[ 3]="LcPz";
1151  fCandidateVariableNames[ 4]="massK0Short";
1152  fCandidateVariableNames[ 5]="V0Px";
1153  fCandidateVariableNames[ 6]="V0Py";
1154  fCandidateVariableNames[ 7]="V0Pz";
1155  fCandidateVariableNames[ 8]="BachPx";
1156  fCandidateVariableNames[ 9]="BachPy";
1157  fCandidateVariableNames[10]="BachPz";
1158  fCandidateVariableNames[11]="PrimVertx";
1159  fCandidateVariableNames[12]="PrimVerty";
1160  fCandidateVariableNames[13]="PrimVertz";
1161  fCandidateVariableNames[14]="Centrality";
1162  fCandidateVariableNames[15]="DecayLengthXY";
1163  fCandidateVariableNames[16]="LcCosPAXY";
1164  fCandidateVariableNames[17]="nSigmaTPCpr";
1165  fCandidateVariableNames[18]="nSigmaTOFpr";
1166  fCandidateVariableNames[19]="probProton";
1167  fCandidateVariableNames[20]="Mixing";
1168  fCandidateVariableNames[21]="Bachd0";
1169  fCandidateVariableNames[22]="V0d0";
1170  fCandidateVariableNames[23]="massLambda";
1171  fCandidateVariableNames[24]="massAntiLambda";
1172  fCandidateVariableNames[25]="massPhoton";
1173  fCandidateVariableNames[26]="DcaPosToPrimVtx";
1174  fCandidateVariableNames[27]="DcaNegToPrimVtx";
1175  fCandidateVariableNames[28]="DcaV0ToPrimVtx";
1176  fCandidateVariableNames[29]="DcaV0Daughters";
1177  fCandidateVariableNames[30]="IsLambdac";
1178  fCandidateVariableNames[31]="SecVertX";
1179  fCandidateVariableNames[32]="SecVertY";
1180  fCandidateVariableNames[33]="SecVertZ";
1181  fCandidateVariableNames[34]="islambdac3body";
1182  fCandidateVariableNames[35]="SecVertXMC";
1183  fCandidateVariableNames[36]="SecVertYMC";
1184  fCandidateVariableNames[37]="SecVertZMC";
1185  fCandidateVariableNames[38]="LcPxMC";
1186  fCandidateVariableNames[39]="LcPyMC";
1187  fCandidateVariableNames[40]="LcPzMC";
1188  fCandidateVariableNames[41]="PrimVertXMC";
1189  fCandidateVariableNames[42]="PrimVertYMC";
1190  fCandidateVariableNames[43]="PrimVertZMC";
1191  fCandidateVariableNames[44]="SecVertXK0s";
1192  fCandidateVariableNames[45]="SecVertYK0s";
1193  fCandidateVariableNames[46]="DcaBachV0";
1194  fCandidateVariableNames[47]="MatchedPDG";
1195  fCandidateVariableNames[48]="ProtonPDG";
1196  fCandidateVariableNames[49]="MotherProtonPDG";
1197  fCandidateVariableNames[50]="GrMotherProtonPDG";
1198  fCandidateVariableNames[51]="V0PDG";
1199  fCandidateVariableNames[52]="MotherV0PDG";
1200  fCandidateVariableNames[53]="GrMotherV0PDG";
1201 
1202  for (Int_t ivar=0; ivar<nVar; ivar++) {
1203  fVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
1204  }
1205 
1206  return;
1207 }
1208 
1211  //
1212  // This is to define general histograms
1213  //
1214 
1215  fCEvents = new TH1F("fCEvents","conter",18,-0.5,17.5);
1216  fCEvents->SetStats(kTRUE);
1217  fCEvents->GetXaxis()->SetBinLabel(1,"X1");
1218  fCEvents->GetXaxis()->SetBinLabel(2,"Analyzed events");
1219  fCEvents->GetXaxis()->SetBinLabel(3,"AliAODVertex exists");
1220  fCEvents->GetXaxis()->SetBinLabel(4,"TriggerOK");
1221  fCEvents->GetXaxis()->SetBinLabel(5,"IsEventSelected");
1222  fCEvents->GetXaxis()->SetBinLabel(6,"CascadesHF exists");
1223  fCEvents->GetXaxis()->SetBinLabel(7,"MCarray exists");
1224  fCEvents->GetXaxis()->SetBinLabel(8,"MCheader exists");
1225  fCEvents->GetXaxis()->SetBinLabel(9,"triggerClass!=CINT1");
1226  fCEvents->GetXaxis()->SetBinLabel(10,"triggerMask!=kAnyINT");
1227  fCEvents->GetXaxis()->SetBinLabel(11,"triggerMask!=kAny");
1228  fCEvents->GetXaxis()->SetBinLabel(12,"vtxTitle.Contains(Z)");
1229  fCEvents->GetXaxis()->SetBinLabel(13,"vtxTitle.Contains(3D)");
1230  fCEvents->GetXaxis()->SetBinLabel(14,"vtxTitle.Doesn'tContain(Z-3D)");
1231  fCEvents->GetXaxis()->SetBinLabel(15,Form("zVtx<=%2.0fcm",fAnalCuts->GetMaxVtxZ()));
1232  fCEvents->GetXaxis()->SetBinLabel(16,"!IsEventSelected");
1233  fCEvents->GetXaxis()->SetBinLabel(17,"triggerMask!=kAnyINT || triggerClass!=CINT1");
1234  fCEvents->GetXaxis()->SetBinLabel(18,Form("zVtxMC<=%2.0fcm",fAnalCuts->GetMaxVtxZ()));
1235  //fCEvents->GetXaxis()->SetTitle("");
1236  fCEvents->GetYaxis()->SetTitle("counts");
1237 
1238  fHTrigger = new TH1F("fHTrigger","counter",18,-0.5,17.5);
1239  fHTrigger->SetStats(kTRUE);
1240  fHTrigger->GetXaxis()->SetBinLabel(1,"X1");
1241  fHTrigger->GetXaxis()->SetBinLabel(2,"kMB");
1242  fHTrigger->GetXaxis()->SetBinLabel(3,"kSemiCentral");
1243  fHTrigger->GetXaxis()->SetBinLabel(4,"kCentral");
1244  fHTrigger->GetXaxis()->SetBinLabel(5,"kINT7");
1245  fHTrigger->GetXaxis()->SetBinLabel(6,"kEMC7");
1246  //fHTrigger->GetXaxis()->SetBinLabel(7,"Space");
1247  fHTrigger->GetXaxis()->SetBinLabel(8,"kMB|kSemiCentral|kCentral");
1248  fHTrigger->GetXaxis()->SetBinLabel(9,"kINT7|kEMC7");
1249  fHTrigger->GetXaxis()->SetBinLabel(11,"kMB&kSemiCentral");
1250  fHTrigger->GetXaxis()->SetBinLabel(12,"kMB&kCentral");
1251  fHTrigger->GetXaxis()->SetBinLabel(13,"kINT7&kEMC7");
1252 
1253  fHCentrality = new TH1F("fHCentrality","conter",100,0.,100.);
1254  fHReactionPlane = new TH1F("fHReactionPlane","conter",100,-1.6,3.2);
1255 
1256  fOutput->Add(fCEvents);
1257  fOutput->Add(fHTrigger);
1258  fOutput->Add(fHCentrality);
1259  fOutput->Add(fHReactionPlane);
1260 
1261  return;
1262 }
1263 //__________________________________________________________________________
1265 {
1266  //
1267  // Define analyis histograms
1268  //
1269 
1270  //------------------------------------------------
1271  // Basic histogram
1272  //------------------------------------------------
1273  Int_t bins_base[3]= {200 ,20 ,10};
1274  Double_t xmin_base[3]={2.286-0.5,0 ,0.00};
1275  Double_t xmax_base[3]={2.286+0.5,20. ,100};
1276  fHistoLcK0SpMass = new THnSparseF("fHistoLcK0SpMass","",3,bins_base,xmin_base,xmax_base);
1278  fHistoLcK0SpMassMix = new THnSparseF("fHistoLcK0SpMassMix","",3,bins_base,xmin_base,xmax_base);
1280  fHistoLcK0SpMassMCS = new THnSparseF("fHistoLcK0SpMassMCS","",3,bins_base,xmin_base,xmax_base);
1282  fHistoLcKPluspMass = new THnSparseF("fHistoLcKPluspMass","",3,bins_base,xmin_base,xmax_base);
1284  fHistoLcKMinuspMass = new THnSparseF("fHistoLcKMinuspMass","",3,bins_base,xmin_base,xmax_base);
1286  fHistoLcKPluspMassMix = new THnSparseF("fHistoLcKPluspMassMix","",3,bins_base,xmin_base,xmax_base);
1288  fHistoLcKMinuspMassMix = new THnSparseF("fHistoLcKMinuspMassMix","",3,bins_base,xmin_base,xmax_base);
1290 
1291  Int_t bins_lcmcgen[3]= {100 ,20 ,10};
1292  Double_t xmin_lcmcgen[3]={0.,-1.0 ,0.0};
1293  Double_t xmax_lcmcgen[3]={20.,1.0 ,100};
1294  fHistoLcMCGen = new THnSparseF("fHistoLcMCGen","",3,bins_lcmcgen,xmin_lcmcgen,xmax_lcmcgen);
1295  fOutputAll->Add(fHistoLcMCGen);
1296 
1297  Int_t bins_coarse[3]= {160 ,20 ,10};
1298  Double_t xmin_coarse[3]={1.,0 ,0.00};
1299  Double_t xmax_coarse[3]={5.8,20. ,100};
1300  fHistoLcK0SpMassCoarse = new THnSparseF("fHistoLcK0SpMassCoarse","",3,bins_coarse,xmin_coarse,xmax_coarse);
1302  fHistoLcK0SpMassMixCoarse = new THnSparseF("fHistoLcK0SpMassMixCoarse","",3,bins_coarse,xmin_coarse,xmax_coarse);
1304 
1305  Int_t bins_correlation[5]= {40 ,20,20,20,10};
1306  Double_t xmin_correlation[5]={-M_PI/2.,-2.4,0.,0.,0.0};
1307  Double_t xmax_correlation[5]={1.5*M_PI,2.4,10.,10.,100};
1308  fHistoK0spCorrelation = new THnSparseF("fHistoK0spCorrelation","",5,bins_correlation,xmin_correlation,xmax_correlation);
1310  fHistoK0spCorrelationMix = new THnSparseF("fHistoK0spCorrelationMix","",5,bins_correlation,xmin_correlation,xmax_correlation);
1312  fHistoK0spCorrelationMCS = new THnSparseF("fHistoK0spCorrelationMCS","",5,bins_correlation,xmin_correlation,xmax_correlation);
1314 
1315 
1316  //------------------------------------------------
1317  // checking histograms
1318  //------------------------------------------------
1319  fHistoBachPt = new TH2F("fHistoBachPt","Bachelor p_{T}",100,0.0,10.0,20,-1.,1.);
1320  fOutputAll->Add(fHistoBachPt);
1321  fHistoBachPtMCS = new TH2F("fHistoBachPtMCS","Bachelor p_{T}",100,0.0,10.0,20,-1.,1.);
1323  fHistoBachPtMCGen = new TH2F("fHistoBachPtMCGen","Bachelor p_{T}",100,0.0,10.0,20,-1.,1.);
1325 
1326  fHistoKaonPt = new TH2F("fHistoKaonPt","Kaon p_{T}",100,0.0,10.0,20,-1.,1.);
1327  fOutputAll->Add(fHistoKaonPt);
1328  fHistoKaonPtMCS = new TH2F("fHistoKaonPtMCS","Kaon p_{T}",100,0.0,10.0,20,-1.,1.);
1330  fHistoKaonPtMCGen = new TH2F("fHistoKaonPtMCGen","Kaon p_{T}",100,0.0,10.0,20,-1.,1.);
1332 
1333  fHistoK0sMassvsPt=new TH3F("fHistoK0sMassvsPt","K0s mass",100,0.497-0.05,0.497+0.05,20,0.,10.,20,-1.,1.);
1335  fHistoK0sMassvsPtMCS=new TH3F("fHistoK0sMassvsPtMCS","K0s mass",100,0.497-0.05,0.497+0.05,20,0.,10.,20,-1.,1.);
1337  fHistoK0sMassvsPtMCGen=new TH3F("fHistoK0sMassvsPtMCGen","K0s mass",100,0.497-0.05,0.497+0.05,20,0.,10.,20,-1.,1.);
1339 
1340  fHistod0Bach = new TH1F("fHistod0Bach","Bachelor d_{0}",100,-0.5,0.5);
1341  fOutputAll->Add(fHistod0Bach);
1342  fHistod0V0 = new TH1F("fHistod0V0","V_{0} d_{0}",100,-0.5,0.5);
1343  fOutputAll->Add(fHistod0V0);
1344  fHistod0d0 = new TH1F("fHistod0d0","d_{0} d_{0}",100,-0.5,0.5);
1345  fOutputAll->Add(fHistod0d0);
1346  fHistoV0CosPA=new TH1F("fHistoV0CosPA","V0->Second vertex cospa",100,-1.,1.0);
1347  fOutputAll->Add(fHistoV0CosPA);
1348  fHistoProbProton=new TH1F("fHistoProbProton","ProbProton",100,0.,1.0);
1350  fHistoDecayLength=new TH1F("fHistoDecayLength","Decay Length",100,-0.1,0.1);
1352  fHistoK0SMass=new TH1F("fHistoK0SMass","K0S mass",100,0.497-0.05,0.497+0.05);
1353  fOutputAll->Add(fHistoK0SMass);
1354 
1355  fHistonEvtvsRunNumber=new TH1F("fHistonEvtvsRunNumber","",20000,-0.5,19999.5);
1357  fHistonProtonvsRunNumber=new TH1F("fHistonProtonvsRunNumber","",20000,-0.5,19999.5);
1359  fHistonK0svsRunNumber=new TH1F("fHistonK0svsRunNumber","",20000,-0.5,19999.5);
1361 
1362  fHistoMassTagV0Min=new TH1F("fHistoMassTagV0Min","",1500,0,1.5);
1364  fHistoMassTagV0SameSignMin=new TH1F("fHistoMassTagV0SameSignMin","",1500,0,1.5);
1366 
1367  return;
1368 }
1369 
1370 //________________________________________________________________________
1371 AliAODRecoCascadeHF* AliAnalysisTaskSELc2pK0sfromAODtracks::MakeCascadeHF(AliAODv0 *v0, AliAODTrack *part, AliAODEvent * aod, AliAODVertex *secVert)
1372 {
1373  //
1374  // Create AliAODRecoCascadeHF object from the argument
1375  //
1376 
1377  if(!v0) return 0x0;
1378  if(!part) return 0x0;
1379  if(!aod) return 0x0;
1380 
1381  //------------------------------------------------
1382  // PrimaryVertex
1383  //------------------------------------------------
1384  AliAODVertex *primVertexAOD;
1385  Bool_t unsetvtx = kFALSE;
1387  primVertexAOD = CallPrimaryVertex(v0,part,aod);
1388  if(!primVertexAOD){
1389  primVertexAOD = fVtx1;
1390  }else{
1391  unsetvtx = kTRUE;
1392  }
1393  }else{
1394  primVertexAOD = fVtx1;
1395  }
1396  if(!primVertexAOD) return 0x0;
1397  Double_t posprim[3]; primVertexAOD->GetXYZ(posprim);
1398 
1399  //------------------------------------------------
1400  // DCA between tracks
1401  //------------------------------------------------
1402  AliESDtrack *esdtrack = new AliESDtrack((AliVTrack*)part);
1403 
1404  AliNeutralTrackParam *trackV0=NULL;
1405  const AliVTrack *trackVV0 = dynamic_cast<const AliVTrack*>(v0);
1406  if(trackVV0) trackV0 = new AliNeutralTrackParam(trackVV0);
1407 
1408  Double_t xdummy, ydummy;
1409  Double_t dca = esdtrack->GetDCA(trackV0,fBzkG,xdummy,ydummy);
1410 
1411 
1412  //------------------------------------------------
1413  // Propagate all tracks to the secondary vertex and calculate momentum there
1414  //------------------------------------------------
1415 
1416  Double_t d0z0bach[2],covd0z0bach[3];
1417 // if(sqrt(pow(secVert->GetX(),2)+pow(secVert->GetY(),2))<1.){
1418 // part->PropagateToDCA(secVert,fBzkG,kVeryBig,d0z0bach,covd0z0bach);
1419 // trackV0->PropagateToDCA(secVert,fBzkG,kVeryBig);
1420 // }else{
1421 // part->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0bach,covd0z0bach);
1422 // trackV0->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig);
1423 // }
1424  part->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0bach,covd0z0bach);
1425  trackV0->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig);
1426  Double_t momv0_new[3]={-9999,-9999,-9999.};
1427  trackV0->GetPxPyPz(momv0_new);
1428 
1429  Double_t px[2],py[2],pz[2];
1430  px[0] = part->Px(); py[0] = part->Py(); pz[0] = part->Pz();
1431  px[1] = momv0_new[0]; py[1] = momv0_new[1]; pz[1] = momv0_new[2];
1432 
1433  //------------------------------------------------
1434  // d0
1435  //------------------------------------------------
1436  Double_t d0[3],d0err[3];
1437 
1438  part->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0bach,covd0z0bach);
1439  d0[0]= d0z0bach[0];
1440  d0err[0] = TMath::Sqrt(covd0z0bach[0]);
1441 
1442  Double_t d0z0v0[2],covd0z0v0[3];
1443  trackV0->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0v0,covd0z0v0);
1444  d0[1]= d0z0v0[0];
1445  d0err[1] = TMath::Sqrt(covd0z0v0[0]);
1446 
1447  //------------------------------------------------
1448  // Create AliAODRecoCascadeHF
1449  //------------------------------------------------
1450  Short_t charge = part->Charge();
1451  AliAODRecoCascadeHF *theCascade = new AliAODRecoCascadeHF(secVert,charge,px,py,pz,d0,d0err,dca);
1452  if(!theCascade)
1453  {
1454  if(unsetvtx) delete primVertexAOD; primVertexAOD=NULL;
1455  if(esdtrack) delete esdtrack;
1456  if(trackV0) delete trackV0;
1457  return 0x0;
1458  }
1459  theCascade->SetOwnPrimaryVtx(primVertexAOD);
1460  UShort_t id[2]={(UShort_t)part->GetID(),(UShort_t)trackV0->GetID()};
1461  theCascade->SetProngIDs(2,id);
1462 
1463  theCascade->GetSecondaryVtx()->AddDaughter(part);
1464  theCascade->GetSecondaryVtx()->AddDaughter(v0);
1465 
1466  if(unsetvtx) delete primVertexAOD; primVertexAOD=NULL;
1467  if(esdtrack) delete esdtrack;
1468  if(trackV0) delete trackV0;
1469 
1470  return theCascade;
1471 }
1472 
1473 //________________________________________________________________________
1474 AliAODVertex* AliAnalysisTaskSELc2pK0sfromAODtracks::CallPrimaryVertex(AliAODv0 *v0, AliAODTrack *trk, AliAODEvent* aod)
1475 {
1476  //
1477  // Make an array of tracks which should not be used in primary vertex calculation and
1478  // Call PrimaryVertex function
1479  //
1480 
1481  TObjArray *TrackArray = new TObjArray(3);
1482 
1483  AliESDtrack *cptrk1 = new AliESDtrack((AliVTrack*)trk);
1484  TrackArray->AddAt(cptrk1,0);
1485 
1486  AliESDtrack *cascptrack = new AliESDtrack((AliVTrack*)v0->GetDaughter(0));
1487  TrackArray->AddAt(cascptrack,1);
1488  AliESDtrack *cascntrack = new AliESDtrack((AliVTrack*)v0->GetDaughter(1));
1489  TrackArray->AddAt(cascntrack,2);
1490 
1491  AliAODVertex *newvert = PrimaryVertex(TrackArray,aod);
1492 
1493  for(Int_t i=0;i<3;i++)
1494  {
1495  AliESDtrack *tesd = (AliESDtrack*)TrackArray->UncheckedAt(i);
1496  delete tesd;
1497  }
1498  TrackArray->Clear();
1499  delete TrackArray;
1500 
1501  return newvert;
1502 }
1503 
1504 //________________________________________________________________________
1505 AliAODVertex* AliAnalysisTaskSELc2pK0sfromAODtracks::PrimaryVertex(const TObjArray *trkArray,
1506  AliVEvent *event)
1507 {
1508  //
1509  //Used only for pp
1510  //copied from AliAnalysisVertexingHF (except for the following 3 lines)
1511  //
1512 
1513  Bool_t fRecoPrimVtxSkippingTrks = kTRUE;
1514  Bool_t fRmTrksFromPrimVtx = kFALSE;
1515 
1516  AliESDVertex *vertexESD = 0;
1517  AliAODVertex *vertexAOD = 0;
1518 
1519  //vertexESD = new AliESDVertex(*fV1);
1520 
1521 
1522  if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx) {
1523  // primary vertex from the input event
1524 
1525  vertexESD = new AliESDVertex(*fV1);
1526 
1527  } else {
1528  // primary vertex specific to this candidate
1529 
1530  Int_t nTrks = trkArray->GetEntriesFast();
1531  AliVertexerTracks *vertexer = new AliVertexerTracks(event->GetMagneticField());
1532 
1533  if(fRecoPrimVtxSkippingTrks) {
1534  // recalculating the vertex
1535 
1536  if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraint")) {
1537  Float_t diamondcovxy[3];
1538  event->GetDiamondCovXY(diamondcovxy);
1539  Double_t pos[3]={event->GetDiamondX(),event->GetDiamondY(),0.};
1540  Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
1541  AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
1542  vertexer->SetVtxStart(diamond);
1543  delete diamond; diamond=NULL;
1544  if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraintOnlyFitter"))
1545  vertexer->SetOnlyFitter();
1546  }
1547  Int_t skipped[1000];
1548  Int_t nTrksToSkip=0,id;
1549  AliExternalTrackParam *t = 0;
1550  for(Int_t i=0; i<nTrks; i++) {
1551  t = (AliExternalTrackParam*)trkArray->UncheckedAt(i);
1552  id = (Int_t)t->GetID();
1553  if(id<0) continue;
1554  skipped[nTrksToSkip++] = id;
1555  }
1556  // TEMPORARY FIX
1557  // For AOD, skip also tracks without covariance matrix
1558  Double_t covtest[21];
1559  for(Int_t j=0; j<event->GetNumberOfTracks(); j++) {
1560  AliVTrack *vtrack = (AliVTrack*)event->GetTrack(j);
1561  if(!vtrack->GetCovarianceXYZPxPyPz(covtest)) {
1562  id = (Int_t)vtrack->GetID();
1563  if(id<0) continue;
1564  skipped[nTrksToSkip++] = id;
1565  }
1566  }
1567  for(Int_t ijk=nTrksToSkip; ijk<1000; ijk++) skipped[ijk]=-1;
1568  //
1569  vertexer->SetSkipTracks(nTrksToSkip,skipped);
1570  vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(event);
1571 
1572  } else if(fRmTrksFromPrimVtx && nTrks>0) {
1573  // removing the prongs tracks
1574 
1575  TObjArray rmArray(nTrks);
1576  UShort_t *rmId = new UShort_t[nTrks];
1577  AliESDtrack *esdTrack = 0;
1578  AliESDtrack *t = 0;
1579  for(Int_t i=0; i<nTrks; i++) {
1580  t = (AliESDtrack*)trkArray->UncheckedAt(i);
1581  esdTrack = new AliESDtrack(*t);
1582  rmArray.AddLast(esdTrack);
1583  if(esdTrack->GetID()>=0) {
1584  rmId[i]=(UShort_t)esdTrack->GetID();
1585  } else {
1586  rmId[i]=9999;
1587  }
1588  }
1589  Float_t diamondxy[2]={static_cast<Float_t>(event->GetDiamondX()),static_cast<Float_t>(event->GetDiamondY())};
1590  vertexESD = vertexer->RemoveTracksFromVertex(fV1,&rmArray,rmId,diamondxy);
1591  delete [] rmId; rmId=NULL;
1592  rmArray.Delete();
1593 
1594  }
1595 
1596  delete vertexer; vertexer=NULL;
1597  if(!vertexESD) return vertexAOD;
1598  if(vertexESD->GetNContributors()<=0) {
1599  //AliDebug(2,"vertexing failed");
1600  delete vertexESD; vertexESD=NULL;
1601  return vertexAOD;
1602  }
1603 
1604 
1605  }
1606 
1607  // convert to AliAODVertex
1608  Double_t pos[3],cov[6],chi2perNDF;
1609  vertexESD->GetXYZ(pos); // position
1610  vertexESD->GetCovMatrix(cov); //covariance matrix
1611  chi2perNDF = vertexESD->GetChi2toNDF();
1612  delete vertexESD; vertexESD=NULL;
1613 
1614  vertexAOD = new AliAODVertex(pos,cov,chi2perNDF);
1615 
1616  return vertexAOD;
1617 }
1618 
1619 //________________________________________________________________________
1620 AliAODVertex* AliAnalysisTaskSELc2pK0sfromAODtracks::ReconstructSecondaryVertex(AliAODv0 *v0, AliAODTrack *part, AliAODEvent * aod)
1621 {
1622  //
1623  // Reconstruct secondary vertex from trkArray (Copied from AliAnalysisVertexingHF)
1624  //
1625  //------------------------------------------------
1626  // PrimaryVertex
1627  //------------------------------------------------
1628  AliAODVertex *primVertexAOD;
1629  Bool_t unsetvtx = kFALSE;
1631  primVertexAOD = CallPrimaryVertex(v0,part,aod);
1632  if(!primVertexAOD){
1633  primVertexAOD = fVtx1;
1634  }else{
1635  unsetvtx = kTRUE;
1636  }
1637  }else{
1638  primVertexAOD = fVtx1;
1639  }
1640  if(!primVertexAOD) return 0x0;
1641 
1642  //------------------------------------------------
1643  // Secondary vertex
1644  //------------------------------------------------
1645 
1646  Double_t LcPx = part->Px()+v0->Px();
1647  Double_t LcPy = part->Py()+v0->Py();
1648  Double_t LcPt = TMath::Sqrt(LcPx*LcPx+LcPy*LcPy);
1649 
1650  Double_t d0z0[2],covd0z0[3];
1651  part->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1652  Double_t x0 = primVertexAOD->GetX();
1653  Double_t y0 = primVertexAOD->GetY();
1654  Double_t px0 = LcPx/LcPt;
1655  Double_t py0 = LcPy/LcPt;
1656  Double_t tx[3];
1657  part->GetXYZ(tx);
1658  Double_t x1 = tx[0];
1659  Double_t y1 = tx[1];
1660  part->GetPxPyPz(tx);
1661  Double_t px1 = tx[0];
1662  Double_t py1 = tx[1];
1663  Double_t pt1 = sqrt(px1*px1+py1*py1);
1664  px1 = px1/pt1;
1665  py1 = py1/pt1;
1666 
1667  Double_t dx = x0 - x1;
1668  Double_t dy = y0 - y1;
1669 
1670  Double_t Delta = -px0*py1+py0*px1;
1671  Double_t a0 = -9999.;
1672  if(Delta!=0)
1673  {
1674  a0 = 1./Delta * (py1 * dx - px1 * dy);
1675  }
1676  Double_t neovertx = x0 + a0 * px0;
1677  Double_t neoverty = y0 + a0 * py0;
1678  Double_t z0 = primVertexAOD->GetZ();
1679  Double_t neovertz = z0 + TMath::Abs(a0)*part->Pz()/part->Pt();
1680 
1681  if(unsetvtx) delete primVertexAOD; primVertexAOD=NULL;
1682 
1683  Double_t pos[3],cov[6],chi2perNDF;
1684  pos[0]=neovertx;
1685  pos[1]=neoverty;
1686  pos[2]=neovertz;
1687  cov[0]=0.0;
1688  cov[1]=0.0;
1689  cov[2]=0.0;
1690  cov[3]=0.0;
1691  cov[4]=0.0;
1692  cov[5]=0.0;
1693  chi2perNDF=0.0;
1694  AliAODVertex *secVert = new AliAODVertex(pos,cov,chi2perNDF);
1695  if(!secVert){
1696  return 0x0;
1697  }
1698  return secVert;
1699 }
1700 //________________________________________________________________________
1701 //AliAODVertex* AliAnalysisTaskSELc2pK0sfromAODtracks::ReconstructSecondaryVertex(AliAODv0 *v0, AliAODTrack *part, AliAODEvent * aod)
1702 //{
1703 // //
1704 // // Reconstruct secondary vertex from trkArray (Copied from AliAnalysisVertexingHF)
1705 // // Currently only returns Primary vertex (can we reconstruct secondary vertex from p - v0)
1706 // //
1707 //
1708 // AliAODVertex *primVertexAOD;
1709 // Bool_t unsetvtx = kFALSE;
1710 // if(fReconstructPrimVert){
1711 // primVertexAOD = CallPrimaryVertex(v0,part,aod);
1712 // if(!primVertexAOD){
1713 // primVertexAOD = fVtx1;
1714 // }else{
1715 // unsetvtx = kTRUE;
1716 // }
1717 // }else{
1718 // primVertexAOD = fVtx1;
1719 // }
1720 // if(!primVertexAOD) return 0x0;
1721 //
1722 // AliESDVertex * vertexESD = new AliESDVertex(*fV1);
1723 //
1724 // Double_t pos[3],cov[6],chi2perNDF;
1725 // vertexESD->GetXYZ(pos); // position
1726 // vertexESD->GetCovMatrix(cov); //covariance matrix
1727 // chi2perNDF = vertexESD->GetChi2toNDF();
1728 // delete vertexESD; vertexESD=NULL;
1729 //
1730 // AliAODVertex *secVert = new AliAODVertex(pos,cov,chi2perNDF);
1731 //
1732 // return secVert;
1733 //}
1734 //________________________________________________________________________
1735 void AliAnalysisTaskSELc2pK0sfromAODtracks::SelectTrack( const AliVEvent *event, Int_t trkEntries, Int_t &nSeleTrks,Int_t *seleFlags, TClonesArray *mcArray)
1736 {
1737  //
1738  // Select good tracks using fAnalCuts (AliRDHFCuts object) and return the array of their ids
1739  //
1740 
1741  if(trkEntries==0) return;
1742 
1743  nSeleTrks=0;
1744  for(Int_t i=0; i<trkEntries; i++) {
1745  seleFlags[i] =0;
1746 
1747  AliVTrack *track;
1748  track = (AliVTrack*)event->GetTrack(i);
1749 
1750  if(track->GetID()<0) continue;
1751  Double_t covtest[21];
1752  if(!track->GetCovarianceXYZPxPyPz(covtest)) continue;
1753 
1754  AliAODTrack *aodt = (AliAODTrack*)track;
1755  if(!fAnalCuts) continue;
1756  if(fAnalCuts->SingleTrkCuts(aodt,fVtx1)){
1757  seleFlags[i]=1;
1758  nSeleTrks++;
1759  FillProtonROOTObjects(aodt,mcArray);
1760 
1761  Double_t minmass;
1762  Bool_t isv0 = fAnalCuts->TagV0(aodt,(AliAODEvent*)event,trkEntries,minmass);
1763  fHistoMassTagV0Min->Fill(minmass);
1764  if(isv0) seleFlags[i] = 0;
1765 
1766 // Double_t minmasslike;
1767 // fAnalCuts->TagV0SameSign(aodt,(AliAODEvent*)event,trkEntries,minmasslike);
1768 // fHistoMassTagV0SameSignMin->Fill(minmasslike);
1769  }
1770  } // end loop on tracks
1771 }
1772 //________________________________________________________________________
1773 void AliAnalysisTaskSELc2pK0sfromAODtracks::SelectV0( const AliVEvent *event,Int_t nV0s,Int_t &nSeleV0, Bool_t *seleV0Flags, TClonesArray *mcArray)
1774 {
1775  //
1776  // Select good V0 using fAnalCuts (AliRDHFCuts object) and return the array of their ids
1777  //
1778 
1779  nSeleV0 = 0;
1780  for(Int_t iv0=0;iv0<nV0s;iv0++)
1781  {
1782  seleV0Flags[iv0] = kFALSE;
1783  AliAODv0 *v0 = ((AliAODEvent*)event)->GetV0(iv0);
1784 
1785  if(!fAnalCuts) continue;
1786  if(fAnalCuts->SingleV0Cuts(v0,fVtx1)){
1787  seleV0Flags[iv0] = kTRUE;
1788  nSeleV0++;
1789 
1790  FillV0ROOTObjects(v0, mcArray);
1791  }
1792  }
1793 }
1794 
1797 {
1798  //
1800  //
1801 
1802  const char* nameoutput = GetOutputSlot(5)->GetContainer()->GetName();
1803  fProtonVariablesTree = new TTree(nameoutput,"proton variables tree");
1804  Int_t nVar = 3;
1805  fCandidateProtonVariables = new Float_t [nVar];
1806  TString * fCandidateVariableNames = new TString[nVar];
1807 
1808  fCandidateVariableNames[ 0]="PrPx";
1809  fCandidateVariableNames[ 1]="PrPy";
1810  fCandidateVariableNames[ 2]="PrPz";
1811 
1812  for (Int_t ivar=0; ivar<nVar; ivar++) {
1813  fProtonVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateProtonVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
1814  }
1815 
1816  return;
1817 }
1818 //________________________________________________________________________
1819 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillProtonROOTObjects(AliAODTrack *trk, TClonesArray *mcArray)
1820 {
1821  //
1823  //
1824 
1825  if(!trk) return;
1826 
1827  TLorentzVector vpr;
1828  vpr.SetXYZM(trk->Px(),trk->Py(),trk->Pz(),0.938272);
1829  fHistoBachPt->Fill(vpr.Pt(),vpr.Rapidity());
1830  if(fDoEventMixing){
1831  fProtonTracks->AddLast(new TLorentzVector(trk->Px(),trk->Py(),trk->Pz(),trk->Charge()));
1832 
1833  Double_t d0z0[2],covd0z0[3];
1834  trk->PropagateToDCA(fVtx1,fBzkG,kVeryBig,d0z0,covd0z0);
1835  Double_t tx[3];
1836  trk->GetXYZ(tx);
1837  Double_t x1 = tx[0];
1838  Double_t y1 = tx[1];
1839  trk->GetPxPyPz(tx);
1840  Double_t px1 = tx[0];
1841  Double_t py1 = tx[1];
1842  Double_t pt1 = sqrt(px1*px1+py1*py1);
1843  px1 = px1/pt1;
1844  py1 = py1/pt1;
1845 
1846  TVector *varvec = new TVector(9);
1847  (*varvec)[0] = d0z0[0];
1848  (*varvec)[1] = x1;
1849  (*varvec)[2] = y1;
1850  (*varvec)[3] = px1;
1851  (*varvec)[4] = py1;
1852  (*varvec)[5] = fVtx1->GetX();
1853  (*varvec)[6] = fVtx1->GetY();
1854  (*varvec)[7] = fAnalCuts->GetPidHF()->GetPidResponse()->NumberOfSigmasTPC(trk,AliPID::kProton);
1855  (*varvec)[8] = fAnalCuts->GetPidHF()->GetPidResponse()->NumberOfSigmasTOF(trk,AliPID::kProton);
1856  fProtonCutVarsArray->AddLast(varvec);
1857  }
1858 
1859  Int_t pdgPr = -9999;
1860  if(fUseMCInfo)
1861  {
1862  Int_t labPr = trk->GetLabel();
1863  if(labPr<0) return;
1864  AliAODMCParticle *mcetrk = (AliAODMCParticle*)mcArray->At(labPr);
1865  if(!mcetrk) return;
1866  pdgPr = mcetrk->GetPdgCode();
1867  if(abs(pdgPr)!=2212) return;
1868  fHistoBachPtMCS->Fill(vpr.Pt(),vpr.Rapidity());
1869  }
1870 }
1873 {
1874  //
1876  //
1877 
1878  const char* nameoutput = GetOutputSlot(6)->GetContainer()->GetName();
1879  fV0VariablesTree = new TTree(nameoutput,"v0 variables tree");
1880  Int_t nVar = 3;
1881  fCandidateV0Variables = new Float_t [nVar];
1882  TString * fCandidateVariableNames = new TString[nVar];
1883 
1884  fCandidateVariableNames[ 0]="V0Px";
1885  fCandidateVariableNames[ 1]="V0Py";
1886  fCandidateVariableNames[ 2]="V0Pz";
1887 
1888  for (Int_t ivar=0; ivar<nVar; ivar++) {
1889  fV0VariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateV0Variables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
1890  }
1891 
1892  return;
1893 }
1894 //________________________________________________________________________
1895 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillV0ROOTObjects(AliAODv0 *v0, TClonesArray *mcArray)
1896 {
1897  //
1899  //
1900  if(!v0) return;
1901 
1902  TLorentzVector vk0;
1903  vk0.SetXYZM(v0->Px(),v0->Py(),v0->Pz(),v0->MassK0Short());
1904  fHistoK0sMassvsPt->Fill(v0->MassK0Short(),v0->Pt(),vk0.Rapidity());
1905 
1906  if(fDoEventMixing){
1907  Double_t ev0 = sqrt(v0->P2()+pow(v0->MassK0Short(),2));
1908  fV0Tracks->AddLast(new TLorentzVector(v0->Px(),v0->Py(),v0->Pz(),ev0));
1909 
1910  TVector *varvec = new TVector(1);
1911  (*varvec)[0] = v0->DcaV0ToPrimVertex();
1912  fV0CutVarsArray->AddLast(varvec);
1913  }
1914 
1915  Int_t v0pdgcode = -9999;
1916  if(fUseMCInfo)
1917  {
1918  Int_t pdgdgv0[2]={211,211};
1919  Int_t labV0 = v0->MatchToMC(310,mcArray,2,pdgdgv0); // the V0
1920  if(labV0<0) return;
1921  AliAODMCParticle *mcv0trk = (AliAODMCParticle*)mcArray->At(labV0);
1922  if(!mcv0trk) return;
1923  fHistoK0sMassvsPtMCS->Fill(v0->MassK0Short(),v0->Pt(),vk0.Rapidity());
1924  }
1925 }
1926 //_________________________________________________________________
1927 Int_t AliAnalysisTaskSELc2pK0sfromAODtracks::GetPoolIndex(Double_t zvert, Double_t mult, Double_t rp){
1928  //
1929  // check in which of the pools the current event falls
1930  //
1931 
1932  Int_t theBinZ=TMath::BinarySearch(fNzVtxBins,fZvtxBins,zvert);
1933  if(theBinZ<0 || theBinZ>=fNzVtxBins) return -1;
1934  Int_t theBinM=TMath::BinarySearch(fNCentBins,fCentBins,mult);
1935  if(theBinM<0 || theBinM>=fNCentBins) return -1;
1936  Int_t theBinR=TMath::BinarySearch(fNRPBins,fRPBins,rp);
1937  if(theBinR<0 || theBinR>=fNRPBins) return -1;
1938  return fNRPBins*fNCentBins*theBinZ+fNRPBins*theBinM+theBinR;
1939 }
1940 //_________________________________________________________________
1942  //
1943  // delete the contets of the pool
1944  //
1945  if(poolIndex<0 || poolIndex>=fNOfPools) return;
1946  delete fEventBuffer[poolIndex];
1947  fEventBuffer[poolIndex]=new TTree(Form("EventBuffer_%d",poolIndex), "Temporary buffer for event mixing");
1948 
1949  fEventBuffer[poolIndex]->Branch("zVertex", &fVtxZ);
1950  fEventBuffer[poolIndex]->Branch("centrality", &fCentrality);
1951  fEventBuffer[poolIndex]->Branch("eventInfo", "TObjString",&fEventInfo);
1952  fEventBuffer[poolIndex]->Branch("v1array", "TObjArray", &fV0Tracks);
1953  fEventBuffer[poolIndex]->Branch("v1varsarray", "TObjArray", &fV0CutVarsArray);
1954 
1955  return;
1956 }
1957 //_________________________________________________________________
1959 {
1960  //
1961  // perform mixed event analysis
1962  //
1963 
1964  if(poolIndex<0 || poolIndex>fNzVtxBins*fNCentBins*fNRPBins) return;
1965  if(fEventBuffer[poolIndex]->GetEntries()<fNumberOfEventsForMixing) return;
1966 
1967  Int_t nPro = fProtonTracks->GetEntries();
1968  Int_t nEvents=fEventBuffer[poolIndex]->GetEntries();
1969  Int_t nPro_test = fProtonCutVarsArray->GetEntries();
1970  if(nPro!=nPro_test){
1971  cout<<"Something wrong in mixing machinery"<<endl;
1972  exit(1);
1973  }
1974 
1975  TObjArray* v1array=0x0;
1976  TObjArray* v1varsarray=0x0;
1977  Float_t zVertex,cent;
1978  TObjString* eventInfo=0x0;
1979  fEventBuffer[poolIndex]->SetBranchAddress("eventInfo",&eventInfo);
1980  fEventBuffer[poolIndex]->SetBranchAddress("zVertex", &zVertex);
1981  fEventBuffer[poolIndex]->SetBranchAddress("centrality", &cent);
1982  fEventBuffer[poolIndex]->SetBranchAddress("v1array", &v1array);
1983  fEventBuffer[poolIndex]->SetBranchAddress("v1varsarray", &v1varsarray);
1984  for (Int_t i=0; i<nPro; i++)
1985  {
1986  TLorentzVector* trke=(TLorentzVector*) fProtonTracks->At(i);
1987  if(!trke)continue;
1988  TVector *prvarsarray = (TVector*)fProtonCutVarsArray->At(i);
1989 
1990  for(Int_t iEv=0; iEv<fNumberOfEventsForMixing; iEv++){
1991  fEventBuffer[poolIndex]->GetEvent(iEv + nEvents - fNumberOfEventsForMixing);
1992  Int_t nV01=v1array->GetEntries();
1993  for(Int_t iTr1=0; iTr1<nV01; iTr1++){
1994  TLorentzVector* v01=(TLorentzVector*)v1array->At(iTr1);
1995  if(!v01 ) continue;
1996  if(!fAnalCuts->SelectWithRoughCuts(v01,trke)) continue;
1997 
1998  TVector *v0varsarray = (TVector*) v1varsarray->At(iTr1);
1999  FillMixROOTObjects(trke,v01,prvarsarray,v0varsarray);
2000  }//v0 loop
2001 
2002  }//event loop
2003  }//track loop
2004 }
2005 
2008 {
2012 
2013  const char* nameoutput = GetOutputSlot(7)->GetContainer()->GetName();
2014  fMCVariablesTree = new TTree(nameoutput,"MC variables tree");
2015  Int_t nVar = 1;
2016  fCandidateMCVariables = new Float_t [nVar];
2017  TString * fCandidateVariableNames = new TString[nVar];
2018 
2019  fCandidateVariableNames[ 0]="Centrality";
2020 
2021  for (Int_t ivar=0; ivar<nVar; ivar++) {
2022  fMCVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateMCVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
2023  }
2024  return;
2025 }
2027 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillMCROOTObjects(AliAODMCParticle *mcpart, AliAODMCParticle *mcepart, AliAODMCParticle *mcv0part, Int_t decaytype)
2028 {
2029  //
2031  //
2032  if(!mcpart) return;
2033  if(!mcepart) return;
2034  if(!mcv0part) return;
2035 
2036 
2037  if(decaytype==0){
2038  Double_t contmc[3];
2039  contmc[0] = mcpart->Pt();
2040  contmc[1] = mcpart->Y();
2041  contmc[2] = fCentrality;
2042  fHistoLcMCGen->Fill(contmc);
2043  }
2044 }
2047 {
2048  //
2049  // Define proton tree variables
2050  //
2051 
2052  const char* nameoutput = GetOutputSlot(9)->GetContainer()->GetName();
2053  fMCProtonVariablesTree = new TTree(nameoutput,"MC Proton variables tree");
2054  Int_t nVar = 1;
2055  fCandidateMCProtonVariables = new Float_t [nVar];
2056  TString * fCandidateVariableNames = new TString[nVar];
2057 
2058  fCandidateVariableNames[ 0]="Centrality";
2059 
2060  for (Int_t ivar=0; ivar<nVar; ivar++) {
2061  fMCProtonVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateMCProtonVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
2062  }
2063  return;
2064 }
2066 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillMCProtonROOTObjects(AliAODMCParticle *mcepart, TClonesArray *mcArray)
2067 {
2068  //
2069  // Fill tree depending on fWriteMCVariableTree
2070  //
2071  if(!mcepart) return;
2072 }
2075 {
2076  //
2077  // Define Mc v0 tree variables
2078  //
2079 
2080  const char* nameoutput = GetOutputSlot(10)->GetContainer()->GetName();
2081  fMCV0VariablesTree = new TTree(nameoutput,"MC v0 variables tree");
2082  Int_t nVar = 1;
2083  fCandidateMCV0Variables = new Float_t [nVar];
2084  TString * fCandidateVariableNames = new TString[nVar];
2085 
2086  fCandidateVariableNames[ 0]="Centrality";
2087 
2088  for (Int_t ivar=0; ivar<nVar; ivar++) {
2089  fMCV0VariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateMCV0Variables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
2090  }
2091  return;
2092 }
2094 void AliAnalysisTaskSELc2pK0sfromAODtracks::FillMCV0ROOTObjects(AliAODMCParticle *mcv0part, TClonesArray *mcArray)
2095 {
2096  //
2097  // Fill histograms or tree depending on fWriteMCVariableTree
2098  //
2099  if(!mcv0part) return;
2100 }
2101 //_________________________________________________________________
2103 {
2104  //
2105  // Analyze AliAODmcparticle
2106  //
2107  Int_t nmcpart = mcArray->GetEntriesFast();
2108 
2109  for(Int_t i=0;i<nmcpart;i++)
2110  {
2111  AliAODMCParticle *mcpart = (AliAODMCParticle*) mcArray->At(i);
2112  if(TMath::Abs(mcpart->GetPdgCode())==4122){
2113  //cout<<"Lambdac"<<endl;
2114  Bool_t p_flag = kFALSE;
2115  Bool_t k0s_flag = kFALSE;
2116  AliAODMCParticle *mcepart = 0;
2117  AliAODMCParticle *mcv0part = 0;
2118  Int_t ndau = mcpart->GetLastDaughter()-mcpart->GetFirstDaughter()+1;
2119  if(ndau==2){
2120  for(Int_t idau=mcpart->GetFirstDaughter();idau<mcpart->GetLastDaughter()+1;idau++)
2121  {
2122  if(idau<0) break;
2123  AliAODMCParticle *mcdau = (AliAODMCParticle*) mcArray->At(idau);
2124  if(!mcdau) continue;
2125  if(TMath::Abs(mcdau->GetPdgCode())==2212){
2126  p_flag = kTRUE;
2127  mcepart = mcdau;
2128  }
2129  if(TMath::Abs(mcdau->GetPdgCode())==311){
2130  k0s_flag = kTRUE;
2131  mcv0part = mcdau;
2132  }
2133  }
2134  }
2135 
2136  Int_t decaytype = -9999;
2137  if(p_flag && k0s_flag) decaytype = 0;
2138 
2139  FillMCROOTObjects(mcpart,mcepart,mcv0part,decaytype);
2140  }
2141 
2142  if(TMath::Abs(mcpart->GetPdgCode())==2212 && mcpart->GetStatus()==1){
2143  AliESDtrackCuts *esdcuts = fAnalCuts->GetTrackCuts();
2144  Float_t etamin, etamax;
2145  esdcuts->GetEtaRange(etamin,etamax);
2146  if(fabs(mcpart->Eta())<etamax){
2147  fHistoBachPtMCGen->Fill(mcpart->Pt(),mcpart->Y());
2148  }
2149  FillMCProtonROOTObjects(mcpart, mcArray);
2150  }
2151  if(TMath::Abs(mcpart->GetPdgCode())==321 && mcpart->GetStatus()==1){
2152  AliESDtrackCuts *esdcuts = fAnalCuts->GetTrackCuts();
2153  Float_t etamin, etamax;
2154  esdcuts->GetEtaRange(etamin,etamax);
2155  if(fabs(mcpart->Eta())<etamax){
2156  fHistoKaonPtMCGen->Fill(mcpart->Pt(),mcpart->Y());
2157  }
2158  }
2159  if(TMath::Abs(mcpart->GetPdgCode())==310){
2160  Double_t etamin, etamax, rapmin, rapmax;
2161  fAnalCuts->GetProdV0EtaRange(etamin,etamax);
2162  fAnalCuts->GetProdV0RapRange(rapmin,rapmax);
2163  if((fabs(mcpart->Y())<rapmax) && (fabs(mcpart->Eta())<etamax)){
2164  fHistoK0sMassvsPtMCGen->Fill(0.497, mcpart->Pt(),mcpart->Y());
2165  }
2166  FillMCV0ROOTObjects(mcpart, mcArray);
2167  }
2168  }
2169 
2170  return kTRUE;
2171 }
2172 //________________________________________________________________________
2173 Int_t AliAnalysisTaskSELc2pK0sfromAODtracks::MatchToMC(AliAODRecoCascadeHF *elobj, TClonesArray *mcArray, Int_t *pdgarray_pr, Int_t *pdgarray_v0, Int_t *labelarray_pr, Int_t *labelarray_v0, Int_t &ngen_pr, Int_t &ngen_v0)
2174 {
2175  //
2176  // Match to MC
2177  //
2178  for(Int_t i=0;i<100;i++){
2179  pdgarray_pr[i] = -9999;
2180  labelarray_pr[i] = -9999;
2181  pdgarray_v0[i] = -9999;
2182  labelarray_v0[i] = -9999;
2183  }
2184  ngen_pr = 0;
2185  ngen_v0 = 0;
2186 
2187  AliVTrack *trk = dynamic_cast<AliVTrack*>(elobj->GetBachelor());
2188  if(!trk) return -1;
2189  Int_t labPr = trk->GetLabel();
2190  if(labPr<0) return -1;
2191  AliAODMCParticle *mcetrk = (AliAODMCParticle*)mcArray->At(labPr);
2192  if(!mcetrk) return -1;
2193  labelarray_pr[0] = labPr;
2194  pdgarray_pr[0] = mcetrk->GetPdgCode();
2195  ngen_pr ++;
2196 
2197  AliAODMCParticle *mcprimpr=0;
2198  mcprimpr = mcetrk;
2199  while(mcprimpr->GetMother()>=0) {
2200  Int_t labprim_pr=mcprimpr->GetMother();
2201  AliAODMCParticle *tmcprimpr = (AliAODMCParticle*)mcArray->At(labprim_pr);
2202  if(!tmcprimpr) {
2203  break;
2204  }
2205 
2206  mcprimpr = tmcprimpr;
2207  pdgarray_pr[ngen_pr] = mcprimpr->GetPdgCode();
2208  labelarray_pr[ngen_pr] = labprim_pr;
2209  ngen_pr ++;
2210  if(ngen_pr==100) break;
2211  }
2212 
2213  AliAODv0 *theV0 = dynamic_cast<AliAODv0*>(elobj->Getv0());
2214  if(!theV0) return -1;
2215  Int_t pdgdgv0[2]={211,211};
2216  Int_t labV0 = theV0->MatchToMC(310,mcArray,2,pdgdgv0); // the V0
2217  if(labV0<0) return -1;
2218  AliAODMCParticle *mcv0 = (AliAODMCParticle*)mcArray->At(labV0);
2219  if(!mcv0) return -1;
2220  labelarray_v0[0] = labV0;
2221  pdgarray_v0[0] = mcv0->GetPdgCode();
2222  ngen_v0 ++;
2223 
2224  AliAODMCParticle *mcprimv0=0;
2225  mcprimv0 = mcv0;
2226  while(mcprimv0->GetMother()>=0) {
2227  Int_t labprim_v0=mcprimv0->GetMother();
2228  AliAODMCParticle *tmcprimv0 = (AliAODMCParticle*)mcArray->At(labprim_v0);
2229  if(!tmcprimv0) {
2230  break;
2231  }
2232 
2233  mcprimv0 = tmcprimv0;
2234  pdgarray_v0[ngen_v0] = mcprimv0->GetPdgCode();
2235  labelarray_v0[ngen_v0] = labprim_v0;
2236  ngen_v0 ++;
2237  if(ngen_v0==100) break;
2238  }
2239 
2240  Bool_t same_flag = kFALSE;
2241  Int_t matchedlabel=-9999;
2242  for(Int_t iemc=0;iemc<ngen_pr;iemc++){
2243  for(Int_t ivmc=0;ivmc<ngen_v0;ivmc++){
2244  if(labelarray_pr[iemc]==labelarray_v0[ivmc]){
2245  same_flag = kTRUE;
2246  matchedlabel = labelarray_pr[iemc];
2247  break;
2248  }
2249  }
2250  if(same_flag) break;
2251  }
2252 
2253  return matchedlabel;
2254 
2255 }
2256 
Int_t charge
TTree * fMCVariablesTree
! tree of the candidate variables after track selection on output slot 7
void FillProtonROOTObjects(AliAODTrack *trk, TClonesArray *mcArray)
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Int_t MatchToMC(AliAODRecoCascadeHF *elobj, TClonesArray *mcArray, Int_t *pdgarray_pr, Int_t *pdgarray_v0, Int_t *labelarray_pr, Int_t *labelarray_v0, Int_t &ngen_pr, Int_t &ngen_v0)
TH2F * fHistoBachPtMCGen
! Bachelor pT histogram (efficiency denominator)
TObjArray * fProtonTracks
unique event id for mixed event check
Bool_t GetUseCombined()
Definition: AliAODPidHF.h:165
void SelectTrack(const AliVEvent *event, Int_t trkEntries, Int_t &nSeleTrks, Int_t *seleFlags, TClonesArray *mcArray)
AliAODVertex * PrimaryVertex(const TObjArray *trkArray, AliVEvent *event)
TH1F * fHistoProbProton
! Probability to be proton histogram
AliAODv0 * Getv0() const
TH1F * fHReactionPlane
! Histogram to check Reaction plane
Float_t * fCandidateProtonVariables
! variables to be written to the tree
TTree * fVariablesTree
flag to decide whether to write the candidate variables on a tree variables
void GetProdV0EtaRange(Double_t &a, Double_t &b)
void FillROOTObjects(AliAODRecoCascadeHF *lcobj, AliAODv0 *v0, AliAODTrack *trk, AliAODEvent *aod, TClonesArray *mcarray)
TH2F * fHistoKaonPtMCGen
! Kaon pT histogram (efficiency denominator)
Double_t GetMaxVtxZ() const
Definition: AliRDHFCuts.h:241
TH2F * fHistoKaonPtMCS
! Kaon pT histogram (efficiency numerator)
Bool_t SingleV0Cuts(AliAODv0 *v0, AliAODVertex *vert)
Double_t InvMassLctoK0sP() const
AliAODPidHF * GetPidHF() const
Definition: AliRDHFCuts.h:231
Float_t * fCandidateV0Variables
! variables to be written to the tree
const Double_t etamin
TH1F * fHCentrality
! Histogram to check Centrality
TTree * fMCProtonVariablesTree
! tree of the candidate variables after track selection on output slot 9
TTree * fMCV0VariablesTree
! tree of the candidate variables after track selection on output slot 10
TTree * fProtonVariablesTree
! tree of the candidate variables after track selection on output slot 5
Double_t CalculateLcCosPAXY(AliAODRecoDecayHF *obj)
virtual void UserCreateOutputObjects()
Implementation of interface methods.
void MakeAnalysis(AliAODEvent *aod, TClonesArray *mcArray)
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
AliAODTrack * GetBachelor() const
TH2F * fHistoBachPtMCS
! Bachelor pT histogram (efficiency numerator)
Bool_t SingleTrkCuts(AliAODTrack *trk, AliAODVertex *vtx)
TObjArray * fV0CutVarsArray
array of RDHF cut information
AliESDtrackCuts * GetTrackCuts() const
Definition: AliRDHFCuts.h:246
void FillMCV0ROOTObjects(AliAODMCParticle *mcv0part, TClonesArray *mcArray)
void SetProngIDs(Int_t nIDs, UShort_t *id)
void GetProdV0RapRange(Double_t &a, Double_t &b)
AliPIDResponse * GetPidResponse() const
Definition: AliAODPidHF.h:160
Int_t GetPoolIndex(Double_t zvert, Double_t mult, Double_t rp)
void SelectV0(const AliVEvent *event, Int_t nV0, Int_t &nSeleV0, Bool_t *seleV0Flags, TClonesArray *mcArray)
Double_t CosV0PointingAngle() const
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Bool_t IsEventRejectedDuePhysicsSelection() const
Definition: AliRDHFCuts.h:320
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
Double_t DecayLengthXY() const
TTree * fV0VariablesTree
! tree of the candidate variables after track selection on output slot 6
Int_t fNzVtxBins
maximum number of events to be used in event mixing
void FillV0ROOTObjects(AliAODv0 *v0, TClonesArray *mcArray)
Float_t * fCandidateMCVariables
! variables to be written to the tree
void FillMixROOTObjects(TLorentzVector *pt, TLorentzVector *ev, TVector *tinfo, TVector *vinfo)
Bool_t IsEventSelected(AliVEvent *event)
void StoreEvent(AliVEvent *, AliRDHFCuts *, Bool_t mc=kFALSE, Int_t multiplicity=-9999, Double_t spherocity=-99.)
void FillMCROOTObjects(AliAODMCParticle *part, AliAODMCParticle *mcepart, AliAODMCParticle *mcv0part, Int_t decaytype)
Float_t nEvents[nProd]
Bool_t fIsMB
Reconstruct primary vertex excluding candidate tracks.
TList * fOutputAll
! User Output slot 3 //analysis histograms
Bool_t TagV0(AliAODTrack *etrk, AliAODEvent *evt, Int_t ntrk, Double_t &minmass)
const Double_t etamax
Float_t * fCandidateMCV0Variables
! variables to be written to the tree
Bool_t SelectWithRoughCuts(AliAODv0 *v0, AliAODTrack *trk1)
Bool_t GetIsUsePID() const
Definition: AliRDHFCuts.h:254
TH1F * fCEvents
! Histogram to check selected events
Bool_t IsEventRejectedDueToTrigger() const
Definition: AliRDHFCuts.h:296
Float_t * fCandidateMCProtonVariables
! variables to be written to the tree
AliAODRecoCascadeHF * MakeCascadeHF(AliAODv0 *casc, AliAODTrack *trk, AliAODEvent *aod, AliAODVertex *vert)
void FillMCProtonROOTObjects(AliAODMCParticle *mcepart, TClonesArray *mcArray)
Double_t fCentBins[100]
number of centrality bins
AliNormalizationCounter * fCounter
EvNumber counter.
AliAODVertex * CallPrimaryVertex(AliAODv0 *v0, AliAODTrack *trk, AliAODEvent *evt)
TObjArray * fV0Tracks
array of electron-compatible tracks
AliAODVertex * ReconstructSecondaryVertex(AliAODv0 *casc, AliAODTrack *trk, AliAODEvent *aod)
TH1F * fHistoV0CosPA
! V0 cosine pointing angle to primary vertex
TObjArray * fProtonCutVarsArray
array of lambda-compatible tracks