Yvonne's update
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliHFEtaggedTrackAnalysis.cxx
1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 * appear 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 // Class AliHFEtaggedTrackAnalysis
17 // Analyses tracks with an apriori PID information (i.e. using the daugther
18 // tracks from well-identified decays of neutral charged particles). Tracks
19 // are processed in the Process function, where given tracks are filtered 
20 // via the track cuts and used for PID later. The plugin fills Correction
21 // Framework containers and additional PID QA containers
22 //
23 // Author:
24 //   Markus Fasel <M.Fasel@gsi.de>
25 //
26 #include "AliCFCutBase.h"
27 #include "AliCFContainer.h"
28 #include "AliCFManager.h"
29 #include "AliESDtrack.h"
30 #include "AliLog.h"
31 #include "AliPID.h"
32 #include "AliPIDResponse.h"
33 #include "AliExternalTrackParam.h"
34
35 #include "AliHFEcollection.h"
36 #include "AliHFEcontainer.h"
37 #include "AliHFEcuts.h"
38 #include "AliHFEpid.h"
39 #include "AliHFEpidTPC.h"
40 #include "AliHFEpidQAmanager.h"
41 #include "AliHFEtaggedTrackAnalysis.h"
42 #include "AliHFEvarManager.h"
43
44 ClassImp(AliHFEtaggedTrackAnalysis)
45 //____________________________________________________________
46 AliHFEtaggedTrackAnalysis::AliHFEtaggedTrackAnalysis():
47     TNamed()
48   , fVarManager(NULL)
49   , fContainer(NULL)
50   , fPID(NULL)
51   , fPIDqa(NULL)
52   , fCuts(NULL)
53   , fCFM(NULL)
54   , fQAhistos(NULL)
55   , fCentralityF(0)
56   , fClean(kFALSE)
57   , fMagneticField(0.0)
58   , fVariablesTRD(kFALSE)
59   , fIsPbPb(kFALSE)
60   , fIspPb(kFALSE)
61   , fIsAOD(kFALSE) 
62 {
63   //
64   // Dummy constructor
65   //
66 }
67
68 //____________________________________________________________
69 AliHFEtaggedTrackAnalysis::AliHFEtaggedTrackAnalysis(const char *name):
70     TNamed(name, "")
71   , fVarManager(NULL)
72   , fContainer(NULL)
73   , fPID(NULL)
74   , fPIDqa(NULL)
75   , fCuts(NULL)
76   , fCFM(NULL)
77   , fQAhistos(NULL)
78   , fCentralityF(0)
79   , fClean(kFALSE)
80   , fMagneticField(0.0)
81   , fVariablesTRD(kFALSE)
82   , fIsPbPb(kFALSE)
83   , fIspPb(kFALSE)
84   , fIsAOD(kFALSE)
85 {
86   //
87   // Default constructor
88   //
89   fVarManager = new AliHFEvarManager("taggedTrackVarManager");
90   //fVarManager->AddVariable("pt");
91   //fVarManager->AddVariable("eta");
92   //fVarManager->AddVariable("phi");
93   //fVarManager->AddVariable("charge");
94   //fVarManager->AddVariable("species");
95   fPIDqa = new AliHFEpidQAmanager;
96   fCFM = new AliCFManager;
97   SetBit(kIsOwner, kTRUE);
98 }
99
100 //____________________________________________________________
101 AliHFEtaggedTrackAnalysis::AliHFEtaggedTrackAnalysis(const AliHFEtaggedTrackAnalysis &ref):
102     TNamed(ref)
103   , fVarManager(ref.fVarManager)
104   , fContainer(NULL)
105   , fPID(ref.fPID)
106   , fPIDqa(ref.fPIDqa)
107   , fCuts(ref.fCuts)
108   , fCFM(ref.fCFM)
109   , fQAhistos(ref.fQAhistos)
110   , fCentralityF(ref.fCentralityF)
111   , fClean(ref.fClean)
112   , fMagneticField(ref.fMagneticField)
113   , fVariablesTRD(ref.fVariablesTRD)
114   , fIsPbPb(ref.fIsPbPb)
115   , fIspPb(ref.fIspPb)
116   , fIsAOD(ref.fIsAOD)
117 {
118   //
119   // Copy constructor
120   //
121   if(ref.fContainer){
122     InitContainer();
123   }
124   SetBit(kIsOwner, kFALSE);
125 }
126
127 //____________________________________________________________
128 AliHFEtaggedTrackAnalysis &AliHFEtaggedTrackAnalysis::operator=(const AliHFEtaggedTrackAnalysis &ref){
129   //
130   // Assignment operator
131   //
132   TNamed::operator=(ref);
133   if(&ref != this){
134     fVarManager = ref.fVarManager;
135     fPID = ref.fPID;
136     fPIDqa = ref.fPIDqa;
137     fCuts = ref.fCuts;
138     fCFM = ref.fCFM;
139     fQAhistos = ref.fQAhistos;
140     fCentralityF = ref.fCentralityF;
141     fClean = ref.fClean;
142     fMagneticField = ref.fMagneticField;
143     fVariablesTRD = ref.fVariablesTRD;
144     fIsPbPb = ref.fIsPbPb;
145     fIspPb = ref.fIspPb;
146     fIsAOD = ref.fIsAOD;
147
148     if(ref.fContainer) InitContainer();
149    
150     SetBit(kIsOwner, kFALSE); 
151     SetBit(kIsOwnerCuts, kFALSE);
152   }
153   return *this;
154 }
155
156 //____________________________________________________________
157 AliHFEtaggedTrackAnalysis::~AliHFEtaggedTrackAnalysis(){
158   //
159   // Destructor
160   //
161   if(TestBit(kIsOwner)){
162     if(fVarManager) delete fVarManager;
163     if(fPIDqa) delete fPIDqa;
164   }
165   if(TestBit(kIsOwnerCuts)) delete fCuts;
166   if(fContainer) delete fContainer;
167 }
168
169 //____________________________________________________________
170 void AliHFEtaggedTrackAnalysis::InitContainer(){
171   //
172   // Initialize output container
173   //
174   if(fContainer) return;
175   Int_t nStepPID = 0;
176   if(!fPID){
177     AliError("No PID set - defining container without PID steps");
178   } else {
179     nStepPID = fPID->GetNumberOfPIDdetectors();
180   }
181   fContainer = new AliHFEcontainer("containerV0");
182   fVarManager->DefineVariables(fContainer);
183   fContainer->CreateContainer("taggedTrackContainerReco", "Container for Tagged Tracks", AliHFEcuts::kNcutStepsRecTrack + nStepPID);
184
185   // Set the step titles
186   for(Int_t istep = 0; istep < AliHFEcuts::kNcutStepsRecTrack; istep++)
187     fContainer->SetStepTitle("taggedTrackContainerReco", AliHFEcuts::RecoCutName(istep), istep);
188   for(Int_t ipid = 0; ipid < nStepPID; ipid++){
189     fContainer->SetStepTitle("taggedTrackContainerReco", fPID->SortedDetectorName(ipid), ipid + AliHFEcuts::kNcutStepsRecTrack);
190   }
191   fCFM->SetParticleContainer(fContainer->GetCFContainer("taggedTrackContainerReco"));
192
193   // temporarily special QA
194   fQAhistos = new AliHFEcollection("taggedTrackQA", "Special QA for the TaggedTrackAnalysis");
195   fQAhistos->CreateTH2F("TPCclusters2_1", "TPCclusterInfo for findable clusters for 2 neighbors", 30, 0.1, 10., 162, 0., 161.);
196   fQAhistos->CreateTH2F("TPCclusters2_0", "TPCclusterInfo for the ratio for 2 neighbors", 30, 0.1, 10., 110, 0., 1.1);
197   fQAhistos->CreateTH2F("TPCncls", "TPC number of clusters", 30, 0.1, 10., 162, 0., 161.);
198   fQAhistos->CreateTH2F("TPCclr", "TPC cluster ratio", 30, 0.1, 10., 110, 0., 1.1);
199   fQAhistos->BinLogAxis("TPCclusters2_1", 0);   // pt axis in logarithmic binning
200   fQAhistos->BinLogAxis("TPCclusters2_0", 0);   // pt axis in logarithmic binning
201   fQAhistos->BinLogAxis("TPCncls", 0);   // pt axis in logarithmic binning
202   fQAhistos->BinLogAxis("TPCclr", 0);   // pt axis in logarithmic binning
203 }
204
205 //____________________________________________________________
206 void AliHFEtaggedTrackAnalysis::ProcessTrack(AliVTrack *track, Int_t abinitioPID){
207   //
208   // Filter tracks tagged by V0 PID class
209   //
210   //
211   fVarManager->NewTrack(track, NULL, fCentralityF, abinitioPID, kTRUE);
212
213
214
215   // Phi Angle
216   if(fVariablesTRD) {
217     AliESDtrack *esdtrackc = dynamic_cast<AliESDtrack *>(track);
218     if(esdtrackc) {
219       
220       const AliExternalTrackParam *trueparam = NULL;
221       if(esdtrackc->GetOuterParam()) {
222               trueparam = esdtrackc->GetOuterParam();
223               fVarManager->NewTrack((AliVParticle *)trueparam, NULL, fCentralityF, abinitioPID, kTRUE);
224       }
225       else return;
226     }
227   }
228   
229
230   // Try a loose cut to reject pion contamination
231   if(fClean) {
232     if(abinitioPID == AliPID::kElectron){
233       AliHFEpidTPC *pidTPC = (AliHFEpidTPC *) fPID->GetDetPID(AliHFEpid::kTPCpid);
234       if(pidTPC) {
235               Double_t numberOfSigmaTPC = pidTPC->GetPIDResponse()->NumberOfSigmasTPC(track,AliPID::kElectron);
236               if(numberOfSigmaTPC < -5) return;
237       }
238     }
239   }
240   // temporarily monitoring of the number of TPC clusters 
241   AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
242   if(esdtrack && abinitioPID == AliPID::kElectron){
243     if((esdtrack->GetITSClusterMap() & (BIT(0) | BIT(1))) && (TMath::Abs(esdtrack->Eta()) < 0.8)){  // Only select quasi-primary tracks
244       fQAhistos->Fill("TPCclusters2_1", track->Pt(), esdtrack->GetTPCClusterInfo(2,1));
245       fQAhistos->Fill("TPCclusters2_0", track->Pt(), esdtrack->GetTPCNclsF() > 0 ? esdtrack->GetTPCClusterInfo(2,1)/esdtrack->GetTPCNclsF() : 0.);
246       fQAhistos->Fill("TPCncls", track->Pt(), esdtrack->GetTPCNcls());
247       fQAhistos->Fill("TPCclr", track->Pt(), esdtrack->GetTPCNclsF() > 0 ? static_cast<Double_t>(esdtrack->GetTPCNcls())/static_cast<Double_t>(esdtrack->GetTPCNclsF()) : 0.);
248     }
249    }
250    
251   Int_t offset = AliHFEcuts::kStepRecKineITSTPC;
252   fVarManager->FillContainer(fCFM->GetParticleContainer(), 0); // Fill Container without filtering
253   
254   Bool_t survived = kTRUE;
255   for(Int_t icut = AliHFEcuts::kStepRecKineITSTPC; icut <= AliHFEcuts::kStepHFEcutsTRD; icut++){
256     AliDebug(2, Form("Checking cut %d for species %s", icut + AliHFEcuts::kNcutStepsMCTrack, AliPID::ParticleName(abinitioPID)));
257     /*
258       TObjArray *cutlist = fCFM->GetParticleCutsList(icut + AliHFEcuts::kNcutStepsMCTrack);
259       if(!cutlist){
260       AliDebug(2, Form("No cuts for step %d set", icut + AliHFEcuts::kNcutStepsMCTrack));
261       } else {
262       AliDebug(2, Form("Cut Collection %s", cutlist->GetName()));
263       TIter cutiter(cutlist);
264       AliCFCutBase *cut;
265       while((cut = dynamic_cast<AliCFCutBase *>(cutiter()))){
266       AliDebug(2, Form("Cut object %s, QA on? %s", cut->GetName(), cut->IsQAOn() ? "yes" : "no"));
267       }
268       }
269     */
270     //if(!fCFM->CheckParticleCuts(icut + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)){
271     if(!fCuts->CheckParticleCuts(icut + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)){
272       AliDebug(2, Form("Track didn' survive cut %d", icut + AliHFEcuts::kNcutStepsMCTrack));
273       survived = kFALSE;
274       break;
275     }
276     AliDebug(2, Form("Cut passed, filling container %d", icut - offset + 1));
277     fVarManager->FillContainer(fCFM->GetParticleContainer(), icut - offset + 1);
278   }
279   
280    if(survived){
281      AliDebug(2, "Use track in the PID");
282      // Apply PID
283      AliHFEpidObject hfetrack;
284      if(!fIsAOD) hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
285      else hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);
286      hfetrack.SetRecTrack(track);
287      hfetrack.SetAbInitioPID(abinitioPID);
288      hfetrack.SetCentrality(fCentralityF);
289      if(fIsPbPb) hfetrack.SetPbPb();
290      else{
291       if(fIspPb) hfetrack.SetpPb();
292             else hfetrack.SetPP();
293      }
294      fPID->SetVarManager(fVarManager);
295      fPID->IsSelected(&hfetrack, fContainer, "taggedTrackContainer", fPIDqa);
296    }
297 }
298
299 //____________________________________________________________
300 void AliHFEtaggedTrackAnalysis::SetCuts(AliHFEcuts *cuts){
301   //
302   // Set HFE cuts to be used to filter the tagged tracks
303   //
304   if(!cuts){
305     AliWarning("Nob cuts provided - Using standard cuts");
306     fCuts = new AliHFEcuts("cutsTagged", "HFE Cuts for the V0 tagged tracks");
307     fCuts->CreateStandardCuts();
308     fCuts->SetQAOn();
309     SetBit(kIsOwnerCuts);
310   } else {
311     AliDebug(1, "Setting single track cuts");
312     fCuts = cuts;
313   }
314   const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
315   printf("Setting Number of cut steps %d\n", kNcutSteps);
316   fCFM->SetNStepParticle(kNcutSteps);
317   for(Int_t istep = 0; istep < kNcutSteps; istep++)
318     fCFM->SetParticleCutsList(istep, NULL);
319
320   fCuts->Initialize(fCFM); 
321 }
322
323 //____________________________________________________________
324 void AliHFEtaggedTrackAnalysis::SetPID(AliHFEpid *pid){
325   //
326   // Set the PID and initialize the the QA manager
327   //
328   fPID = pid;
329   fPIDqa->Initialize(fPID);
330 }
331
332 //____________________________________________________________
333 TList *AliHFEtaggedTrackAnalysis::GetPIDQA() const {
334   //
335   // return PID QA 
336   //
337   return fPIDqa->MakeList("PIDqa_taggedTracks");
338 }
339
340 //____________________________________________________________
341 TList *AliHFEtaggedTrackAnalysis::GetCutQA() const {
342   //
343   // return Cut QA
344   //
345   return fCuts->GetQAhistograms();
346 }
347
348