AliPhysics  b752f14 (b752f14)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliV0ReaderStrange.cxx
Go to the documentation of this file.
1 #include <vector>
2 #include <TGeoGlobalMagField.h>
3 
4 #include "AliV0ReaderStrange.h"
5 #include "AliKFParticle.h"
6 #include "AliV0ParticleStrange.h"
7 #include "AliAODv0.h"
8 #include "AliESDv0.h"
9 #include "AliV0.h"
10 #include "AliAODEvent.h"
11 #include "AliESDEvent.h"
12 #include "AliVEvent.h"
13 #include "AliPID.h"
14 #include "AliMCEvent.h"
15 #include "AliStack.h"
16 #include "AliMCEventHandler.h"
17 #include "AliESDpid.h"
18 #include "AliESDtrackCuts.h"
19 #include "TRandom3.h"
20 #include "AliGenCocktailEventHeader.h"
21 #include "TList.h"
22 #include "AliKFConversionPhoton.h"
23 #include "AliAODConversionPhoton.h"
25 #include "TVector.h"
26 #include "AliKFVertex.h"
27 #include "AliAODTrack.h"
28 #include "AliESDtrack.h"
29 #include "AliAnalysisManager.h"
30 #include "AliInputEventHandler.h"
31 #include "AliAODHandler.h"
32 #include "AliPIDResponse.h"
33 #include "TChain.h"
34 #include "TFile.h"
35 #include "TString.h"
36 #include "TObjArray.h"
37 #include "AliVTrack.h"
38 #include "AliKFParticle.h"
39 
40 class iostream;
41 
42 using namespace std;
43 
45 
46 //________________________________________________________________________
48  fEventCuts(NULL),
49  fV0Cuts(NULL),
50  fVectorFoundGammas(0),
51  fConversionGammas(NULL),
52  fEventIsSelected(kFALSE)
53 {
54  // Default constructor
55 
56  DefineInput(0, TChain::Class());
57 }
58 
59 //________________________________________________________________________
61 {
62  // default deconstructor
63 
64  if(fConversionGammas){
65  fConversionGammas->Delete();// Clear Objects
66  delete fConversionGammas;
67  fConversionGammas=0x0;
68  }
69 }
70 
71 
72 //________________________________________________________________________
74 {
75  // Initialize function to be called once before analysis
76  if(fV0Cuts==NULL){
77  if(fV0Cuts==NULL)AliError("No V0 Cut Selection initialized");
78  }
79  if(fEventCuts==NULL){
80  if(fEventCuts==NULL)AliError("No Event Cut Selection initialized");
81  }
82 
83  if(fConversionGammas != NULL){
84  delete fConversionGammas;
85  fConversionGammas=NULL;
86  }
87 
88  if(fConversionGammas == NULL){
89  fConversionGammas = new TClonesArray("AliKFParticle",1000);
90 // if(kUseAODConversionPhoton){
91 // fConversionGammas = new TClonesArray("AliAODConversionPhoton",100);}
92 // else{
93 // fConversionGammas = new TClonesArray("AliKFConversionPhoton",100);}
94  }
95  fConversionGammas->Delete();//Reset the TClonesArray
96 }
97 
98 //________________________________________________________________________
100 {
101  fVectorFoundGammas.clear();
102 }
103 
104 //________________________________________________________________________
106 
107  AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
108  if(esdEvent) {
109  if (!TGeoGlobalMagField::Instance()->GetField()) esdEvent->InitMagneticField();
110  }
111 
112  // Check if correctly initialized
113  if(!fConversionGammas)Init();
114 
115  // User Exec
116  fEventIsSelected=ProcessEvent(fInputEvent,fMCEvent);
117 }
118 
119 //________________________________________________________________________
120 Bool_t AliV0ReaderStrange::ProcessEvent(AliVEvent *inputEvent,AliMCEvent *mcEvent)
121 {
122  //Reset the TClonesArray
123  fConversionGammas->Delete();
124 
125  fInputEvent=inputEvent;
126  fMCEvent=mcEvent;
127 
128  if(!fInputEvent){
129  AliError("No Input event");
130  return kFALSE;
131  }
132  if(!fEventCuts){AliError("No EventCuts");return kFALSE;}
133  if(!fV0Cuts){AliError("No V0 Cuts");return kFALSE;}
134 
135  // Event Cuts
136  if(!fEventCuts->EventIsSelected(fInputEvent,fMCEvent))return kFALSE;
137 
138  // Set Magnetic Field
139  AliKFParticle::SetField(fInputEvent->GetMagneticField());
140 
141  if(fInputEvent->IsA()==AliESDEvent::Class()){
142  ProcessESDV0s();
143  }
144  if(fInputEvent->IsA()==AliAODEvent::Class()){
145  ProcessAODV0s();
146  }
147 
148  return kTRUE;
149 }
150 
153 {
154  // Process ESD V0s for V0 reconstruction
155 
156  AliV0ParticleStrange *fCurrentMotherLambdaCandidate=NULL;
157 
158  AliESDEvent *fESDEvent=dynamic_cast<AliESDEvent*>(fInputEvent);
159 
160  for(Int_t iV0=0; iV0<fESDEvent->GetNumberOfV0s(); ++iV0){
161  AliESDv0 *fCurrentV0=(AliESDv0*)(fESDEvent->GetV0(iV0));
162  if(!fCurrentV0){
163  printf("Requested V0 does not exist");
164  continue;
165  }
166 
167  fCurrentMotherLambdaCandidate = ReconstructV0(fESDEvent, fCurrentV0, iV0);
168  if(fCurrentMotherLambdaCandidate){
169  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliV0ParticleStrange(*fCurrentMotherLambdaCandidate);
170  delete fCurrentMotherLambdaCandidate;
171  fCurrentMotherLambdaCandidate=NULL;
172  }
173  }
174 
175  return kTRUE;
176 }
177 
179 AliV0ParticleStrange *AliV0ReaderStrange::ReconstructV0(AliESDEvent *fESDEvent, AliESDv0 *fCurrentV0,Int_t currentV0Index)
180 {
181  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kV0In);
182 
183  //checks if on the fly mode is set
184  if(!fV0Cuts->SelectV0Finder(fCurrentV0->GetOnFlyStatus())){
185  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kOnFly);
186  return 0x0;
187  }
188 
189 // Double_t lV0CosineOfPointingAngle = fCurrentV0->GetV0CosineOfPointingAngle();
190 // GetV0CosineOfPointingAngle(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
191 
192  if(fCurrentV0->GetV0CosineOfPointingAngle() < 0.98) return 0x0;
193 
194  AliVTrack * pos = (AliVTrack*)fESDEvent->GetTrack(fCurrentV0->GetPindex());
195  AliVTrack * neg = (AliVTrack*)fESDEvent->GetTrack(fCurrentV0->GetNindex());
196  const AliExternalTrackParam * paramPos = fCurrentV0->GetParamP();
197  const AliExternalTrackParam * paramNeg = fCurrentV0->GetParamN();
198 
199  if(pos->GetSign() <0){//change tracks
200  pos=neg;
201  neg=fESDEvent->GetTrack(fCurrentV0->GetPindex());
202  paramPos=paramNeg;
203  paramNeg=fCurrentV0->GetParamP();
204  }
205 
206  if(!pos || !neg ) {
207  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kNoTracks);
208  return 0x0;
209 // return;
210  }
211  //cuts ------------ //remove like sign pairs
212  if(pos->GetSign() == neg->GetSign()){
213  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kSameSign);
214  return 0x0;
215  }
216 
217 // if(pos->Pt() < .16) return 0x0;
218 // if(neg->Pt() < .16) return 0x0;
219 // if(TMath::Abs(pos->Eta()) > 0.8) return 0x0;
220 // if(TMath::Abs(neg->Eta()) > 0.8) return 0x0;
221 
222 // if( neg->GetKinkIndex(0) > 0 || pos->GetKinkIndex(0) > 0) return 0x0 ;
223 // if( !(pos->GetStatus() & AliESDtrack::kTPCrefit) || !(neg->GetStatus() & AliESDtrack::kTPCrefit) ) return 0x0;
224 
225  Bool_t isPiPlus = fV0Cuts->GetPIDpion(pos);
226  Bool_t isPiMinus = fV0Cuts->GetPIDpion(neg);
227  Bool_t isProton = fV0Cuts->GetPIDproton(pos);
228  Bool_t isAntiProton = fV0Cuts->GetPIDproton(neg);
229 
230 
231  if(!((isProton && isPiMinus) || (isAntiProton && isPiPlus))) return 0x0;
232 
233  // pi plus: 211, pi minus: -211, proton: 2212, anti-proton: -2212
234 
235  Int_t pdgCodePos = 2211; //proton
236  Int_t pdgCodeNeg = -211; // pi minus
237  if(isPiPlus){
238  pdgCodePos = -211;
239  }
240  if(isAntiProton){
241  pdgCodeNeg = -2212;
242  }
243 
244 // AliVVertex *vertex = fESDEvent->GetPrimaryVertex();
245 // fCurrentV0->DecayLength(vertex);
246 
247  const AliKFParticle fCurrentPositiveKFParticle(*(paramPos), pdgCodePos);
248  const AliKFParticle fCurrentNegativeKFParticle(*(paramNeg), pdgCodeNeg);
249  const Bool_t gamma = kFALSE;
250  AliKFParticle *fCurrentMother = new AliKFParticle(fCurrentPositiveKFParticle, fCurrentNegativeKFParticle, gamma);
251 
252  AliV0ParticleStrange *fCurrentMotherV0 = new AliV0ParticleStrange(fCurrentMother);
253 
254  if(fMCEvent){
255 
256  AliStack *fMCStack= fMCEvent->Stack();
257 
258  Int_t labelp=TMath::Abs(fV0Cuts->GetTrack(fESDEvent,fCurrentMotherV0->GetTrackLabelPositive())->GetLabel());
259  Int_t labeln=TMath::Abs(fV0Cuts->GetTrack(fESDEvent,fCurrentMotherV0->GetTrackLabelNegative())->GetLabel());
260 
261 // cout << "rec: " << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
262 // cout << "recProp: " << fCurrentMotherKF->GetTrackLabelPositive() << "\t" << fCurrentMotherKF->GetTrackLabelNegative() << endl;
263 // cout << "MC: " << labeln << "\t" << labelp << endl;
264 
265  TParticle *fNegativeMCParticle = 0x0;
266  if(labeln>-1) fNegativeMCParticle = fMCStack->Particle(labeln);
267  TParticle *fPositiveMCParticle = 0x0;
268  if(labelp>-1) fPositiveMCParticle = fMCStack->Particle(labelp);
269 
270  if(fPositiveMCParticle&&fNegativeMCParticle){
271  fCurrentMotherV0->SetMCLabelPositive(labelp);
272  fCurrentMotherV0->SetMCLabelNegative(labeln);
273  }
274  }
275 
276  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kV0Out);
277  return fCurrentMotherV0;
278 }
279 
280 
283 {
284  // Process ESD V0s for V0 reconstruction
285 
286  AliV0ParticleStrange *fCurrentMotherLambdaCandidate=NULL;
287 
288  AliAODEvent *fAODEvent=dynamic_cast<AliAODEvent*>(fInputEvent);
289 
290  for(Int_t iV0=0; iV0<fAODEvent->GetNumberOfV0s(); ++iV0){
291  AliAODv0 *fCurrentV0=(AliAODv0*)(fAODEvent->GetV0(iV0));
292  if(!fCurrentV0){
293  printf("Requested V0 does not exist");
294  continue;
295  }
296 
297  fCurrentMotherLambdaCandidate = ReconstructV0(fAODEvent, fCurrentV0, iV0);
298  if(fCurrentMotherLambdaCandidate){
299  new((*fConversionGammas)[fConversionGammas->GetEntriesFast()]) AliV0ParticleStrange(*fCurrentMotherLambdaCandidate);
300  delete fCurrentMotherLambdaCandidate;
301  fCurrentMotherLambdaCandidate=NULL;
302  }
303  }
304 
305  return kTRUE;
306 }
307 
309 AliV0ParticleStrange *AliV0ReaderStrange::ReconstructV0(AliAODEvent *fAODEvent, AliAODv0 *fCurrentV0,Int_t currentV0Index)
310 {
311  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kV0In);
312 
313  //checks if on the fly mode is set
314  if(!fV0Cuts->SelectV0Finder(fCurrentV0->GetOnFlyStatus())){
315  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kOnFly);
316  return 0x0;
317  }
318 
319 // Double_t lV0CosineOfPointingAngle = fCurrentV0->GetV0CosineOfPointingAngle();
320 // GetV0CosineOfPointingAngle(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
321 
322  AliAODVertex *vertex = fAODEvent->GetPrimaryVertex();
323  if(fCurrentV0->CosPointingAngle(vertex) < 0.98) return 0x0;
324 
325  AliVTrack* pos = (AliVTrack*)fCurrentV0->GetDaughter(0);
326  AliVTrack* neg = (AliVTrack*)fCurrentV0->GetDaughter(1);
327 
328  AliExternalTrackParam paramPos;
329  paramPos.CopyFromVTrack(pos);
330  AliExternalTrackParam paramNeg;
331  paramPos.CopyFromVTrack(neg);
332 
333  if(pos->GetSign() <0){//change tracks
334  pos=neg;
335  neg=(AliVTrack*)fCurrentV0->GetDaughter(0);
336  paramPos=paramNeg;
337  paramNeg.CopyFromVTrack(pos);
338  }
339 
340  if(!pos || !neg ) {
341  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kNoTracks);
342  return 0x0;
343  }
344 
345  //cuts ------------ //remove like sign pairs
346  if(pos->GetSign() == neg->GetSign()){
347  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kSameSign);
348  return 0x0;
349  }
350 
351 // if(pos->Pt() < .16) return 0x0;
352 // if(neg->Pt() < .16) return 0x0;
353 // if(TMath::Abs(pos->Eta()) > 0.8) return 0x0;
354 // if(TMath::Abs(neg->Eta()) > 0.8) return 0x0;
355 
356 // if( neg->GetKinkIndex(0) > 0 || pos->GetKinkIndex(0) > 0) return 0x0 ;
357 // if( !(pos->GetStatus() & AliESDtrack::kTPCrefit) || !(neg->GetStatus() & AliESDtrack::kTPCrefit) ) return 0x0;
358 
359  Bool_t isPiPlus = fV0Cuts->GetPIDpion(pos);
360  Bool_t isPiMinus = fV0Cuts->GetPIDpion(neg);
361  Bool_t isProton = fV0Cuts->GetPIDproton(pos);
362  Bool_t isAntiProton = fV0Cuts->GetPIDproton(neg);
363 
364 
365  if(!((isProton && isPiMinus) || (isAntiProton && isPiPlus))) return 0x0;
366 
367  // pi plus: 211, pi minus: -211, proton: 2212, anti-proton: -2212
368 
369  Int_t pdgCodePos = 2211; //proton
370  Int_t pdgCodeNeg = -211; // pi minus
371  if(isPiPlus){
372  pdgCodePos = -211;
373  }
374  if(isAntiProton){
375  pdgCodeNeg = -2212;
376  }
377 
378  const AliKFParticle fCurrentPositiveKFParticle(paramPos, pdgCodePos);
379  const AliKFParticle fCurrentNegativeKFParticle(paramNeg, pdgCodeNeg);
380  const Bool_t gamma = kFALSE;
381  AliKFParticle *fCurrentMother = new AliKFParticle(fCurrentPositiveKFParticle, fCurrentNegativeKFParticle, gamma);
382 
383  AliV0ParticleStrange *fCurrentMotherV0 = new AliV0ParticleStrange(fCurrentMother);
384 
385  if(fMCEvent){
386  AliStack *fMCStack= fMCEvent->Stack();
387 
388  Int_t labelp=TMath::Abs(fV0Cuts->GetTrack(fAODEvent,fCurrentMotherV0->GetTrackLabelPositive())->GetLabel());
389  Int_t labeln=TMath::Abs(fV0Cuts->GetTrack(fAODEvent,fCurrentMotherV0->GetTrackLabelNegative())->GetLabel());
390 
391 // cout << "rec: " << currentTrackLabels[0] << "\t" << currentTrackLabels[1] << endl;
392 // cout << "recProp: " << fCurrentMotherKF->GetTrackLabelPositive() << "\t" << fCurrentMotherKF->GetTrackLabelNegative() << endl;
393 // cout << "MC: " << labeln << "\t" << labelp << endl;
394 
395  TParticle *fNegativeMCParticle = 0x0;
396  if(labeln>-1) fNegativeMCParticle = fMCStack->Particle(labeln);
397  TParticle *fPositiveMCParticle = 0x0;
398  if(labelp>-1) fPositiveMCParticle = fMCStack->Particle(labelp);
399 
400  if(fPositiveMCParticle&&fNegativeMCParticle){
401  fCurrentMotherV0->SetMCLabelPositive(labelp);
402  fCurrentMotherV0->SetMCLabelNegative(labeln);
403  }
404  }
405 
406 
407  fV0Cuts->FillV0CutIndex(AliV0CutsStrange::kV0Out);
408  return fCurrentMotherV0;
409 }
void SetMCLabelNegative(Int_t label)
int Int_t
Definition: External.C:63
void SetMCLabelPositive(Int_t label)
ClassImp(AliV0ReaderStrange) AliV0ReaderStrange
virtual void UserExec(Option_t *option)
Bool_t ProcessEvent(AliVEvent *inputEvent, AliMCEvent *mcEvent=NULL)
AliV0ParticleStrange * ReconstructV0(AliESDEvent *fESDEvent, AliESDv0 *fCurrentV0, Int_t currentV0Index)
const char Option_t
Definition: External.C:48
bool Bool_t
Definition: External.C:53