AliPhysics  b752f14 (b752f14)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAODConversionPhoton.cxx
Go to the documentation of this file.
1 #include "AliAODMCParticle.h"
2 #include "AliAODMCHeader.h"
5 
6 using namespace std;
7 
9 
13 fDCArPrimVtx(0),
14 fDCAzPrimVtx(0),
15 fInvMassPair(0),
16 fCaloPhoton(0),
17 fCaloClusterRef(-1),
18 fNCaloPhotonMCLabels(0),
19 fNCaloPhotonMotherMCLabels(0),
20 fCaloPhotonMCFlags(0)
21 {
22  // initialize calo photon MC labels
23  for (Int_t i =0; i<50; i++){
24  fCaloPhotonMCLabels[i]=-1;
25  }
26  for (Int_t i =0; i<20; i++){
27  fCaloPhotonMotherMCLabels[i]=-1;
28  }
29  //Standard constructor
30 }
31 
33 AliAODConversionParticle(kfphoton),
35 fDCArPrimVtx(0),
36 fDCAzPrimVtx(0),
37 fInvMassPair(0),
38 fCaloPhoton(0),
39 fCaloClusterRef(-1),
40 fNCaloPhotonMCLabels(0),
41 fNCaloPhotonMotherMCLabels(0),
42 fCaloPhotonMCFlags(0)
43 {
44  //Constructor from kfphoton
45 
46  // puts the mass to zero and store dilepton mass
47  SetMass(kfphoton->M());
48 
49  //SetE(P());
50 
51  // initialize calo photon MC labels
52  for (Int_t i =0; i<50; i++){
53  fCaloPhotonMCLabels[i]=-1;
54  }
55  for (Int_t i =0; i<20; i++){
57  }
58 
59 }
60 
64 fDCArPrimVtx(0),
65 fDCAzPrimVtx(0),
66 fInvMassPair(0),
67 fCaloPhoton(0),
68 fCaloClusterRef(-1),
69 fNCaloPhotonMCLabels(0),
70 fNCaloPhotonMotherMCLabels(0),
71 fCaloPhotonMCFlags(0)
72 {
73  //Constructor from TLorentzVector
74 
75  // initialize calo photon MC labels
76  for (Int_t i =0; i<50; i++){
77  fCaloPhotonMCLabels[i]=-1;
78  }
79  for (Int_t i =0; i<20; i++){
81  }
82 }
83 
84 
85 
87 AliAODConversionParticle(original),
88 AliConversionPhotonBase(original),
89 fDCArPrimVtx(original.fDCArPrimVtx),
90 fDCAzPrimVtx(original.fDCAzPrimVtx),
91 fInvMassPair(original.fInvMassPair),
92 fCaloPhoton(original.fCaloPhoton),
93 fCaloClusterRef(original.fCaloClusterRef),
94 fNCaloPhotonMCLabels(original.fNCaloPhotonMCLabels),
95 fNCaloPhotonMotherMCLabels(original.fNCaloPhotonMotherMCLabels),
96 fCaloPhotonMCFlags(original.fCaloPhotonMCFlags)
97 {
98  //Copy constructor
99 
100  // initialize calo photon MC labels
101  for (Int_t i =0; i<50; i++){
103  }
104  for (Int_t i =0; i<20; i++){
106  }
107 }
108 
110 {
111  // empty standard destructor
112 }
113 
115 {
116  // assignment operator
117  return *this;
118 }
119 
122 
123  Double_t primCo[3] = {primVertex->GetX(),primVertex->GetY(),primVertex->GetZ()};
124 
125  Double_t absoluteP = TMath::Sqrt(TMath::Power(GetPx(),2) + TMath::Power(GetPy(),2) + TMath::Power(GetPz(),2));
126  Double_t p[3] = {GetPx()/absoluteP,GetPy()/absoluteP,GetPz()/absoluteP};
127  Double_t CP[3];
128 
129  CP[0] = fConversionPoint[0] - primCo[0];
130  CP[1] = fConversionPoint[1] - primCo[1];
131  CP[2] = fConversionPoint[2] - primCo[2];
132 
133  Double_t Lambda = - (CP[0]*p[0]+CP[1]*p[1]+CP[2]*p[2])/(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
134 
135  Double_t S[3];
136  S[0] = fConversionPoint[0] + p[0]*Lambda;
137  S[1] = fConversionPoint[1] + p[1]*Lambda;
138  S[2] = fConversionPoint[2] + p[2]*Lambda;
139 
140  fDCArPrimVtx = TMath::Sqrt( TMath::Power(primCo[0]-S[0],2) + TMath::Power(primCo[1]-S[1],2));
141  fDCAzPrimVtx = primCo[2]-S[2];
142 
143  return;
144 }
145 
146 
147 void AliAODConversionPhoton::SetCaloPhotonMCFlags(AliStack *MCStack, Bool_t enableSort){
148 
149  Bool_t isPhoton = kFALSE; // largest contribution to cluster is photon
150  Bool_t isElectron = kFALSE; // largest contribution to cluster is electron
151  Bool_t isConversion = kFALSE; // largest contribution to cluster is converted electron
152  Bool_t isConversionFullyContained = kFALSE; // largest contribution to cluster is converted electron, second electron has been found in same cluster
153  Bool_t isMerged = kFALSE; // cluster contains more than one particle from the same decay
154  Bool_t isMergedPartConv = kFALSE; // cluster contains more than one particle from the same decay and at least one of the particles came from a conversion
155  Bool_t isDalitz = kFALSE; // this cluster was created by a particle stemming from a dalitz decay
156  Bool_t isDalitzMerged = kFALSE; // this cluster was created by a particle stemming from a dalitz decay and more than one particle of the dalitz decay is contained in the cluster
157  Bool_t isPhotonWithElecMother = kFALSE; // this cluster is from a photon with an electron as mother
158  Bool_t isShower = kFALSE; // this cluster contains as a largest contribution a particle from a shower or radiative process
159  Bool_t isSubLeadingEM = kFALSE; // cluster contains at least one electron or photon from a pi0, eta or eta_prime in subleading contribution
160 
161 
162  TParticle* Photon = 0x0;
163  if (fNCaloPhotonMCLabels==0) return;
164 
165  if (enableSort){
166  // sort the array according to the energy of contributing particles
167  if (fNCaloPhotonMCLabels>1){
168  // cout << "start sorting" << endl;
169  Int_t* sortIdx = new Int_t[fNCaloPhotonMCLabels];
170  Double_t* energyPerPart = new Double_t[fNCaloPhotonMCLabels];
171  Long_t* orginalContrib = new Long_t[fNCaloPhotonMCLabels];
172  for (Int_t i = 0; i < fNCaloPhotonMCLabels; i++){
173  orginalContrib[i] = fCaloPhotonMCLabels[i];
174  if (fCaloPhotonMCLabels[i]> -1){
175  TParticle* dummy = MCStack->Particle(fCaloPhotonMCLabels[i]);
176  energyPerPart[i] = dummy->Energy();
177  // suppress energy of hadrons !!! DIRTY hack !!!
178  if (!(TMath::Abs(dummy->GetPdgCode())== 11 || TMath::Abs(dummy->GetPdgCode())== 22)){
179  // cout << "suppressed hadron energy for:" << dummy->GetPdgCode() << endl;
180  energyPerPart[i]= 0.25; // put energy to mip
181  }
182  } else {
183  energyPerPart[i] = 0;
184  }
185  }
186 
187  TMath::Sort(fNCaloPhotonMCLabels,energyPerPart,sortIdx);
188  for(Int_t index = 0; index < fNCaloPhotonMCLabels; index++) {
189  fCaloPhotonMCLabels[index] = orginalContrib [sortIdx[index]] ;
190 
191  }
192  delete [] sortIdx;
193  delete [] energyPerPart;
194  delete [] orginalContrib;
195  }
196  }
197 
198  if(GetCaloPhotonMCLabel(0)>-1) Photon = MCStack->Particle(GetCaloPhotonMCLabel(0));
199 
200  if(Photon == NULL){
201  return;
202  }
203 
204  Int_t particleMotherLabel = Photon->GetMother(0);
205  Int_t particleGrandMotherLabel = -1;
206  Int_t particleMotherPDG = -1;
207  Int_t particleGrandMotherPDG = -1;
208  Int_t particleMotherNDaugthers = 0;
209  Int_t particleGrandMotherNDaugthers = 0;
210  if (particleMotherLabel > -1){
211  particleMotherNDaugthers = MCStack->Particle(Photon->GetMother(0))->GetNDaughters();
212  particleGrandMotherLabel = MCStack->Particle(Photon->GetMother(0))->GetMother(0);
213  particleMotherPDG = MCStack->Particle(Photon->GetMother(0))->GetPdgCode();
214  if (particleGrandMotherLabel > -1){
215  particleGrandMotherPDG = MCStack->Particle(MCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetPdgCode();
216  particleGrandMotherNDaugthers = MCStack->Particle(MCStack->Particle(Photon->GetMother(0))->GetMother(0))->GetNDaughters();
217  }
218  }
219 
220  //determine mother/grandmother of leading particle and if it is pion/eta/eta_prime: fill array fCaloPhotonMotherMCLabels at position 0
221  //determine mother/grandmother of leading particle and if it is pion/eta/eta_prime: fill array fCaloPhotonMotherMCLabels at position 0
222  if (particleMotherLabel > -1){
223  if( TMath::Abs(particleMotherPDG) == 111 || TMath::Abs(particleMotherPDG) == 221 || TMath::Abs(particleMotherPDG) == 331 ){
224  fCaloPhotonMotherMCLabels[0] = particleMotherLabel;
226  } else if (TMath::Abs(particleMotherPDG) == 22 && particleGrandMotherLabel > -1){
227  if ( TMath::Abs(particleMotherPDG) == 22 && (TMath::Abs(particleGrandMotherPDG) == 111 || TMath::Abs(particleGrandMotherPDG) == 221 || TMath::Abs(particleGrandMotherPDG) == 331) ){
228  fCaloPhotonMotherMCLabels[0] = particleGrandMotherLabel;
230  }
231  } else {
232  fCaloPhotonMotherMCLabels[0] = particleMotherLabel;
234  }
235  }
236 
237  // Check whether the first contribution was photon
238  if(TMath::Abs(MCStack->Particle(GetCaloPhotonMCLabel(0))->GetPdgCode()) == 22){
239  isPhoton = kTRUE;
240  // did it decay via the dalitz channel
241  if (particleMotherLabel > -1 && particleMotherNDaugthers == 3)
242  isDalitz = kTRUE;
243  // Test whether particle stems from a shower or radiation
244  if (TMath::Abs(particleMotherPDG) == 11){ // check whether photon stems from electron
245  isPhotonWithElecMother = kTRUE;
246  if (particleGrandMotherLabel > -1){ // test whether first particle has a grandmother
247  if (TMath::Abs(particleGrandMotherPDG) == 22 )
248  isShower = kTRUE; // check whether grandmother is a photon (meaning this is most likely a shower)
249  }
250  }
251  }
252  // Check whether the first contribution was electron
253  if( TMath::Abs(MCStack->Particle(GetCaloPhotonMCLabel(0))->GetPdgCode()) == 11 ){
254  isElectron = kTRUE;
255  if (particleMotherLabel > -1) {
256  // was it a conversion
257  if (TMath::Abs(particleMotherPDG) == 22)
258  isConversion = kTRUE;
259  // did it decay via the dalitz channel
260  if (particleGrandMotherLabel > -1 && particleGrandMotherNDaugthers == 3 )
261  isDalitz = kTRUE;
262  }
263  if (particleGrandMotherLabel > -1){ // check whether electron has a grandmother
264  if (TMath::Abs(particleGrandMotherPDG) == 11 || TMath::Abs(particleGrandMotherPDG) == 22){ // test whether electron has photon or electron as grandmother (meaning will most likely be a shower)
265  isShower = kTRUE;
266  }
267  }
268  }
269 
270  Bool_t enablePrintOuts = kFALSE;
271 // if (fNCaloPhotonMCLabels>2)
272 // enablePrintOuts = kTRUE;
273 
274  // check whether there were other contributions to the cluster
275  if (fNCaloPhotonMCLabels>1){
276  TParticle* dummyPart =NULL;
277 
278  for (Int_t i = 1; i< fNCaloPhotonMCLabels; i++){
279  if (i > 49) continue; // abort if more than 50 entries to the cluster have been checked (more are not stored in these objects)
280  if (enablePrintOuts) cout << "checking particle: " << i << endl;
281  if (GetCaloPhotonMCLabel(i) < 0) continue;
282  dummyPart = MCStack->Particle(GetCaloPhotonMCLabel(i));
283  Int_t dummyPartMotherLabel = dummyPart->GetMother(0);
284  Int_t dummyPartGrandMotherLabel = -1;
285  Int_t dummyPartMotherPDG = -1;
286  Int_t dummyPartGrandMotherPDG = -1;
287  // check whether this particle has a mother & obtain the pdg code
288  if (dummyPartMotherLabel > -1){
289  dummyPartGrandMotherLabel = MCStack->Particle(dummyPart->GetMother(0))->GetMother(0);
290  dummyPartMotherPDG = MCStack->Particle(dummyPart->GetMother(0))->GetPdgCode();
291  // check whether this particle has a grandmother & obtain its pdg code
292  if (dummyPartGrandMotherLabel > -1){
293  dummyPartGrandMotherPDG = MCStack->Particle(MCStack->Particle(dummyPart->GetMother(0))->GetMother(0))->GetPdgCode();
294  }
295  }
296  // largest contribution was from photon and is not from shower or electron mother
297  if (isPhoton && (!isShower || !isPhotonWithElecMother )){
298  if (enablePrintOuts) cout << "lead gamma" << endl;
299  if (dummyPartMotherLabel > -1 && particleMotherLabel > -1){ // test whether first particle has a mother
300  if (dummyPartMotherLabel == particleMotherLabel)
301  isMerged = kTRUE; // test whether current and first particle have the same mother => i.e. other gamma from decay or dalitz electron
302  if (dummyPartGrandMotherLabel > -1){ // test whether first particle has a grandmother
303  // check whether particle is an electron from a conversion of a photon from the original mother
304  if (TMath::Abs(dummyPart->GetPdgCode()) == 11 && dummyPartGrandMotherLabel == particleMotherLabel )
305  isMergedPartConv = kTRUE;
306  // check whether particle is an electron from a dalitz decay from the original mother
307  if (TMath::Abs(dummyPart->GetPdgCode()) == 11 && dummyPartMotherLabel == particleMotherLabel )
308  isDalitzMerged = kTRUE;
309  }
310  }
311  }
312 
313  // largest contribution was from electron & not a from a shower
314  if (isElectron && !isShower){
315  if (enablePrintOuts) cout << "lead electron" << endl;
316  if (dummyPartMotherLabel > -1 && particleMotherLabel > -1){ // test whether first particle has a mother
317  if ( TMath::Abs(dummyPart->GetPdgCode()) == 11 && isConversion && dummyPartMotherLabel == particleMotherLabel ) {
318  isConversionFullyContained = kTRUE; // test whether conversion is fully contained in cluster
319  if (enablePrintOuts) cout << "found full conversion" << endl;
320  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 11) { // current particle is an electron
321  if (enablePrintOuts) cout << "current is electron" << endl;
322  // check whether current particle is a shower
323  Bool_t isShowerCurrent = kFALSE;
324  if (TMath::Abs(dummyPartGrandMotherPDG) == 11)
325  isShowerCurrent = kTRUE;
326  if (TMath::Abs(dummyPartGrandMotherPDG) == 22)
327  isShowerCurrent = kTRUE;
328 
329  // check whether current particle is a conversion
330  Bool_t isConvCurrent = kFALSE;
331  if (TMath::Abs(dummyPartMotherPDG) == 22 && !isShowerCurrent)
332  isConvCurrent = kTRUE;
333 
334  if (enablePrintOuts) cout << "conv current: " << isConvCurrent << "\t shower current: " << isShowerCurrent << endl;
335  // check whether orginal electron and this electron stem from the same particle and electron originated in dalitz decay
336  if( dummyPartMotherLabel == particleMotherLabel && TMath::Abs(particleMotherPDG) != 22 && !isShowerCurrent ) {
337  isDalitzMerged = kTRUE;
338  isMerged = kTRUE;
339  // both particles are from conversions
340  } else if (isConversion && isConvCurrent){
341  if (particleGrandMotherLabel > -1 && dummyPartGrandMotherLabel > -1){ // test whether first particle has a grandmother and current particle has one as well
342  // check whether orginal electron and this electron stem from the same particle and electron stems from conversion
343  if( dummyPartGrandMotherLabel == particleGrandMotherLabel && !isShowerCurrent )
344  isMergedPartConv = kTRUE;
345  }
346  // check whether current electron is from conversion & orginal e are from same mother
347  } else if (dummyPartGrandMotherLabel > -1 && isConvCurrent){
348  if (dummyPartGrandMotherLabel == particleMotherLabel && (TMath::Abs(particleMotherPDG) == 111 || TMath::Abs(particleMotherPDG) == 221 || TMath::Abs(particleMotherPDG) == 331) ){
349  isDalitzMerged = kTRUE;
350  isMerged = kTRUE;
351  }
352  // check whether current electron & orginal e is from conversion are from same mother
353  } else if (isConversion && particleGrandMotherLabel > -1){
354  if (dummyPartMotherLabel == particleGrandMotherLabel && (TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) ){
355  isDalitzMerged = kTRUE;
356  isMerged = kTRUE;
357  }
358  }
359  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 22){ // test whether this particle is a photon
360  // check whether current particle is a shower
361  Bool_t isShowerCurrent = kFALSE;
362  if (TMath::Abs(dummyPartMotherPDG) == 11)
363  isShowerCurrent = kTRUE;
364  if (TMath::Abs(dummyPartMotherPDG) == 22)
365  isShowerCurrent = kTRUE;
366 
367  // check whether orginal electron and this photon stem from the same particle and electron originated in dalitz
368  if( dummyPartMotherLabel == particleMotherLabel && !isShowerCurrent ) {
369  isDalitzMerged = kTRUE;
370  isMerged = kTRUE;
371  } else if (particleGrandMotherLabel > -1){ // test whether first particle has a grandmother
372  // check whether orginal electron and this photon stem from the same particle and electron stems from conversion
373  if( dummyPartMotherLabel == particleGrandMotherLabel && !isShowerCurrent )
374  isMergedPartConv = kTRUE;
375  }
376  }
377  }
378  }
379 
380 
381  if (dummyPartMotherLabel > -1){ // test whether particle has a mother
382  if (TMath::Abs(dummyPart->GetPdgCode()) == 22){ // test whether particle is a photon
383  //check if photon directly comes from a pion/eta/eta_prime decay
384  Bool_t helpN = true;
385  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
386  if (fCaloPhotonMotherMCLabels[j] == dummyPartMotherLabel){ //check if mother is already contained in fCaloPhotonMotherMCLabels
387  helpN = false;
388  }
389  }
390  if (helpN){
392  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
393  }
394  if ((TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) && !isPhoton && !isElectron)
395  isSubLeadingEM = kTRUE;
396  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 11){ //test whether particle is an electron
397  //check if electron comes from a pion decay
398  if ( TMath::Abs(dummyPartMotherPDG) != 22 ){
399  Bool_t helpN = true;
400  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
401  if (fCaloPhotonMotherMCLabels[j] == dummyPartMotherLabel){ //check if mother is already contained in fCaloPhotonMotherMCLabels
402  helpN = false;
403  }
404  }
405  if (helpN){
407  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
408  }
409  if ((TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) && !isPhoton && !isElectron)
410  isSubLeadingEM = kTRUE;
411  } else if (dummyPartGrandMotherLabel > -1){ //if it is not a dalitz decay, test whether particle has a grandmother
412  //check if it is a conversion electron that has pion/eta/eta_prime as grandmother
413  if ( (TMath::Abs(dummyPartGrandMotherPDG) == 111 || TMath::Abs(dummyPartGrandMotherPDG) == 221 || TMath::Abs(dummyPartGrandMotherPDG) == 331)){
414 
415  Bool_t helpN = true;
416  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
417  if (fCaloPhotonMotherMCLabels[j] == dummyPartGrandMotherLabel){ //check if grandmother is already contained in fCaloPhotonMotherMCLabels
418  helpN = false;
419  }
420  }
421  if (helpN){
422  fCaloPhotonMotherMCLabels[fNCaloPhotonMotherMCLabels] = dummyPartGrandMotherLabel;
423  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
424  }
425  if (!isPhoton && !isElectron)
426  isSubLeadingEM = kTRUE;
427  }
428  }
429  }
430  }
431  }
432  }
433  fCaloPhotonMCFlags = isPhoton *1 + isElectron *2 + isConversion*4+ isConversionFullyContained *8 + isMerged *16 + isMergedPartConv*32 + isDalitz *64 + isDalitzMerged *128 + isPhotonWithElecMother *256 + isShower * 512 + isSubLeadingEM * 1024;
434 }
435 
436 void AliAODConversionPhoton::SetCaloPhotonMCFlagsAOD(AliVEvent* event, Bool_t enableSort){
437 
438  TClonesArray *AODMCTrackArray = dynamic_cast<TClonesArray*>(event->FindListObject(AliAODMCParticle::StdBranchName()));
439  if (!AODMCTrackArray) return;
440 
441  Bool_t isPhoton = kFALSE; // largest contribution to cluster is photon
442  Bool_t isElectron = kFALSE; // largest contribution to cluster is electron
443  Bool_t isConversion = kFALSE; // largest contribution to cluster is converted electron
444  Bool_t isConversionFullyContained = kFALSE; // largest contribution to cluster is converted electron, second electron has been found in same cluster
445  Bool_t isMerged = kFALSE; // largest contribution to cluster is photon, second photon or electron from dalitz has been found in same cluster
446  Bool_t isMergedPartConv = kFALSE; // cluster contains more than one particle from the same decay and at least one of the particles came from a conversion
447  Bool_t isDalitz = kFALSE; // this cluster was created by a particle stemming from a dality decay
448  Bool_t isDalitzMerged = kFALSE; // this cluster was created by a particle stemming from a dality decay and more than one particle of the dalitz decay is contained in the cluster
449  Bool_t isPhotonWithElecMother = kFALSE; // this cluster is from a photon with an electron as mother
450  Bool_t isShower = kFALSE; // this cluster contains as a largest contribution a particle from a shower or radiative process
451  Bool_t isSubLeadingEM = kFALSE; // cluster contains at least one electron or photon from a pi0 or eta in subleading contribution
452 
453  AliAODMCParticle* Photon;
454  AliAODMCParticle* PhotonMother;
455  AliAODMCParticle* PhotonGrandMother;
456 
457  if (fNCaloPhotonMCLabels==0) return;
458  Photon = (AliAODMCParticle*) AODMCTrackArray->At(GetCaloPhotonMCLabel(0));
459 
460  if (enableSort){
461  // sort the array according to the energy of contributing particles
462  if (fNCaloPhotonMCLabels>1){
463  // cout << "start sorting" << endl;
464  Int_t* sortIdx = new Int_t[fNCaloPhotonMCLabels];
465  Double_t* energyPerPart = new Double_t[fNCaloPhotonMCLabels];
466  Long_t* orginalContrib = new Long_t[fNCaloPhotonMCLabels];
467  for (Int_t i = 0; i < fNCaloPhotonMCLabels; i++){
468  orginalContrib[i] = fCaloPhotonMCLabels[i];
469  if (fCaloPhotonMCLabels[i]> -1){
470  AliAODMCParticle* dummy = (AliAODMCParticle*) AODMCTrackArray->At(GetCaloPhotonMCLabel(i));
471  energyPerPart[i] = dummy->E();
472  // suppress energy of hadrons !!! DIRTY hack !!!
473  if (!(TMath::Abs(dummy->GetPdgCode())== 11 || TMath::Abs(dummy->GetPdgCode())== 22)){
474  // cout << "suppressed hadron energy" << endl;
475  energyPerPart[i]= 0.25; // put energy to mip
476  }
477  } else {
478  energyPerPart[i] = 0;
479  }
480  }
481 
482  TMath::Sort(fNCaloPhotonMCLabels,energyPerPart,sortIdx);
483  for(Int_t index = 0; index < fNCaloPhotonMCLabels; index++) {
484  fCaloPhotonMCLabels[index] = orginalContrib [sortIdx[index]] ;
485 
486  }
487  delete [] sortIdx;
488  delete [] energyPerPart;
489  delete [] orginalContrib;
490  }
491  }
492 
493  if(Photon == NULL){
494  return;
495  }
496 
497  Int_t particleMotherLabel = Photon->GetMother();
498  Int_t particleGrandMotherLabel = -1;
499  Int_t particleMotherPDG = -1;
500  Int_t particleGrandMotherPDG = -1;
501  Int_t particleMotherNDaugthers = 0;
502  Int_t particleGrandMotherNDaugthers = 0;
503  if (particleMotherLabel > -1){
504  PhotonMother = (AliAODMCParticle*) AODMCTrackArray->At(Photon->GetMother());
505  particleMotherNDaugthers = PhotonMother->GetNDaughters();
506  particleGrandMotherLabel = PhotonMother->GetMother();
507  particleMotherPDG = PhotonMother->GetPdgCode();
508  if (particleGrandMotherLabel > -1){
509  PhotonGrandMother = (AliAODMCParticle*) AODMCTrackArray->At(PhotonMother->GetMother());
510  particleGrandMotherPDG = PhotonGrandMother->GetPdgCode();
511  particleGrandMotherNDaugthers = PhotonGrandMother->GetNDaughters();
512  }
513  }
514 
515  //determine mother/grandmother of leading particle and if it is pion/eta/eta_prime: fill array fCaloPhotonMotherMCLabels at position 0
516  if (particleMotherLabel > -1){
517  if( TMath::Abs(particleMotherPDG) == 111 || TMath::Abs(particleMotherPDG) == 221 || TMath::Abs(particleMotherPDG) == 331 ){
518  fCaloPhotonMotherMCLabels[0] = particleMotherLabel;
520  } else if (TMath::Abs(particleMotherPDG) == 22 && particleGrandMotherLabel > -1){
521  if ( TMath::Abs(particleMotherPDG) == 22 && (TMath::Abs(particleGrandMotherPDG) == 111 || TMath::Abs(particleGrandMotherPDG) == 221 || TMath::Abs(particleGrandMotherPDG) == 331) ){
522  fCaloPhotonMotherMCLabels[0] = particleGrandMotherLabel;
524  }
525  } else {
526  fCaloPhotonMotherMCLabels[0] = particleMotherLabel;
528  }
529  }
530 
531  // Check whether the first contribution was photon
532  if(TMath::Abs(Photon->GetPdgCode()) == 22){
533  isPhoton = kTRUE;
534  // did it decay via the dalitz channel
535  if (particleMotherLabel > -1 && particleMotherNDaugthers == 3)
536  isDalitz = kTRUE;
537  // Test whether particle stems from a shower or radiation
538  if (TMath::Abs(particleMotherPDG) == 11){ // check whether photon stems from electron
539  isPhotonWithElecMother = kTRUE;
540  if (particleGrandMotherLabel > -1){ // test whether first particle has a grandmother
541  if (TMath::Abs(particleGrandMotherPDG) == 22 )
542  isShower = kTRUE; // check whether grandmother is a photon (meaning this is most likely a shower)
543  }
544  }
545  }
546  // Check whether the first contribution was electron
547  if(TMath::Abs(Photon->GetPdgCode()) == 11 ){
548  isElectron = kTRUE;
549  if (particleMotherLabel > -1) {
550  // was it a conversion
551  if (TMath::Abs(particleMotherPDG) == 22)
552  isConversion = kTRUE;
553  // did it decay via the dalitz channel
554  if (particleGrandMotherLabel > -1 && particleGrandMotherNDaugthers == 3 )
555  isDalitz = kTRUE;
556  }
557  if (particleGrandMotherLabel > -1){ // check whether electron has a grandmother
558  if (TMath::Abs(particleGrandMotherPDG) == 11 || TMath::Abs(particleGrandMotherPDG) == 22){ // test whether electron has photon or electron as grandmother (meaning will most likely be a shower)
559  isShower = kTRUE;
560  }
561  }
562  }
563 
564 
565  // check whether there were other contributions to the cluster
566  if (fNCaloPhotonMCLabels>1){
567  AliAODMCParticle* dummyPart = NULL;
568  AliAODMCParticle* dummyPartMother = NULL;
569  AliAODMCParticle* dummyPartGrandMother = NULL;
570  for (Int_t i = 1; i< fNCaloPhotonMCLabels; i++){
571  if (i > 49) continue; // abort if more than 50 entries to the cluster have been checked (more are not stored in these objects)
572  dummyPart = (AliAODMCParticle*) AODMCTrackArray->At(GetCaloPhotonMCLabel(i));
573  Int_t dummyPartMotherLabel = dummyPart->GetMother();
574  Int_t dummyPartGrandMotherLabel = -1;
575  Int_t dummyPartMotherPDG = -1;
576  Int_t dummyPartGrandMotherPDG = -1;
577 
578  // check whether this particle has a mother & obtain the pdg code
579  if (dummyPartMotherLabel > -1){
580  dummyPartMother = (AliAODMCParticle*) AODMCTrackArray->At(dummyPart->GetMother());
581  dummyPartGrandMotherLabel = dummyPartMother->GetMother();
582  dummyPartMotherPDG = dummyPartMother->GetPdgCode();
583  // check whether this particle has a grandmother & obtain its pdg code
584  if (dummyPartGrandMotherLabel > -1){
585  dummyPartGrandMother = (AliAODMCParticle*) AODMCTrackArray->At(dummyPartMother->GetMother());
586  dummyPartGrandMotherPDG = dummyPartGrandMother->GetPdgCode();
587  }
588  }
589 
590 
591  // largest contribution was from photon and is not from shower or electron mother
592  if (isPhoton && (!isShower || !isPhotonWithElecMother )){
593  if (dummyPartMotherLabel > -1 && particleMotherLabel > -1){ // test whether first particle has a mother
594  if (dummyPartMotherLabel == particleMotherLabel)
595  isMerged = kTRUE; // test whether current and first particle have the same mother => i.e. other gamma from decay or dalitz electron
596  if (dummyPartGrandMotherLabel > -1){ // test whether first particle has a grandmother
597  // check whether particle is an electron from a conversion of a photon from the original mother
598  if (TMath::Abs(dummyPart->GetPdgCode()) == 11 && dummyPartGrandMotherLabel == particleMotherLabel )
599  isMergedPartConv = kTRUE;
600  // check whether particle is an electron from a dalitz decay from the original mother
601  if (TMath::Abs(dummyPart->GetPdgCode()) == 11 && dummyPartMotherLabel == particleMotherLabel )
602  isDalitzMerged = kTRUE;
603  }
604  }
605  }
606 
607  // largest contribution was from electron & not a from a shower
608  if (isElectron && !isShower){
609  if (dummyPartMotherLabel > -1 && particleMotherLabel > -1){ // test whether first particle has a mother
610  if ( TMath::Abs(dummyPart->GetPdgCode()) == 11 && isConversion && dummyPartMotherLabel == particleMotherLabel ) {
611  isConversionFullyContained = kTRUE; // test whether conversion is fully contained in cluster
612  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 11) { // current particle is an electron
613 
614  // check whether current particle is a shower
615  Bool_t isShowerCurrent = kFALSE;
616  if (TMath::Abs(dummyPartGrandMotherPDG) == 11)
617  isShowerCurrent = kTRUE;
618  if (TMath::Abs(dummyPartGrandMotherPDG) == 22)
619  isShowerCurrent = kTRUE;
620 
621  // check whether current particle is a conversion
622  Bool_t isConvCurrent = kFALSE;
623  if (TMath::Abs(dummyPartMotherPDG) == 22 && !isShowerCurrent)
624  isConvCurrent = kTRUE;
625 
626  // check whether orginal electron and this electron stem from the same particle and electron originated in dalitz decay
627  if( dummyPartMotherLabel == particleMotherLabel && TMath::Abs(particleMotherPDG) != 22 && !isShowerCurrent ) {
628  isDalitzMerged = kTRUE;
629  isMerged = kTRUE;
630  // both particles are from conversions
631  } else if (isConversion && isConvCurrent){
632  if (particleGrandMotherLabel > -1 && dummyPartGrandMotherLabel > -1){ // test whether first particle has a grandmother and current particle has one as well
633  // check whether orginal electron and this electron stem from the same particle and electron stems from conversion
634  if( dummyPartGrandMotherLabel == particleGrandMotherLabel && !isShowerCurrent )
635  isMergedPartConv = kTRUE;
636  }
637  // check whether current electron is from conversion & orginal e are from same mother
638  } else if (dummyPartGrandMotherLabel > -1 && isConvCurrent){
639  if (dummyPartGrandMotherLabel == particleMotherLabel && (TMath::Abs(particleMotherPDG) == 111 || TMath::Abs(particleMotherPDG) == 221 || TMath::Abs(particleMotherPDG) == 331) ){
640  isDalitzMerged = kTRUE;
641  isMerged = kTRUE;
642  }
643  // check whether current electron & orginal e is from conversion are from same mother
644  } else if (isConversion && particleGrandMotherLabel > -1){
645  if (dummyPartMotherLabel == particleGrandMotherLabel && (TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) ){
646  isDalitzMerged = kTRUE;
647  isMerged = kTRUE;
648  }
649  }
650  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 22){ // test whether this particle is a photon
651  // check whether current particle is a shower
652  Bool_t isShowerCurrent = kFALSE;
653  if (TMath::Abs(dummyPartMotherPDG) == 11)
654  isShowerCurrent = kTRUE;
655  if (TMath::Abs(dummyPartMotherPDG) == 22)
656  isShowerCurrent = kTRUE;
657 
658  // check whether orginal electron and this photon stem from the same particle and electron originated in dalitz
659  if( dummyPartMotherLabel == particleMotherLabel && !isShowerCurrent ) {
660  isDalitzMerged = kTRUE;
661  isMerged = kTRUE;
662  } else if (particleGrandMotherLabel > -1){ // test whether first particle has a grandmother
663  // check whether orginal electron and this photon stem from the same particle and electron stems from conversion
664  if( dummyPartMotherLabel == particleGrandMotherLabel && !isShowerCurrent )
665  isMergedPartConv = kTRUE;
666  }
667  }
668  }
669  }
670 
671  if (dummyPartMotherLabel > -1){ // test whether particle has a mother
672  if (TMath::Abs(dummyPart->GetPdgCode()) == 22){ // test whether particle is a photon
673  //check if photon directly comes from a pion/eta/eta_prime decay
674  Bool_t helpN = true;
675  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
676  if (fCaloPhotonMotherMCLabels[j] == dummyPartMotherLabel){ //check if mother is already contained in fCaloPhotonMotherMCLabels
677  helpN = false;
678  }
679  }
680  if (helpN){
682  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
683  }
684  if ((TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) && !isPhoton && !isElectron)
685  isSubLeadingEM = kTRUE;
686  } else if (TMath::Abs(dummyPart->GetPdgCode()) == 11){ //test whether particle is an electron
687  //check if electron comes from a pion decay
688  if ( TMath::Abs(dummyPartMotherPDG) != 22 ){
689  Bool_t helpN = true;
690  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
691  if (fCaloPhotonMotherMCLabels[j] == dummyPartMotherLabel){ //check if mother is already contained in fCaloPhotonMotherMCLabels
692  helpN = false;
693  }
694  }
695  if (helpN){
697  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
698  }
699  if ((TMath::Abs(dummyPartMotherPDG) == 111 || TMath::Abs(dummyPartMotherPDG) == 221 || TMath::Abs(dummyPartMotherPDG) == 331) && !isPhoton && !isElectron)
700  isSubLeadingEM = kTRUE;
701  } else if (dummyPartGrandMotherLabel > -1){ //if it is not a dalitz decay, test whether particle has a grandmother
702  //check if it is a conversion electron that has pion/eta/eta_prime as grandmother
703  if ( (TMath::Abs(dummyPartGrandMotherPDG) == 111 || TMath::Abs(dummyPartGrandMotherPDG) == 221 || TMath::Abs(dummyPartGrandMotherPDG) == 331)){
704 
705  Bool_t helpN = true;
706  for(Int_t j=0; j<fNCaloPhotonMotherMCLabels; j++){
707  if (fCaloPhotonMotherMCLabels[j] == dummyPartGrandMotherLabel){ //check if grandmother is already contained in fCaloPhotonMotherMCLabels
708  helpN = false;
709  }
710  }
711  if (helpN){
712  fCaloPhotonMotherMCLabels[fNCaloPhotonMotherMCLabels] = dummyPartGrandMotherLabel;
713  fNCaloPhotonMotherMCLabels++; //only if particle label is not yet contained in array, count up fNCaloPhotonMotherMCLabels
714  }
715  if (!isPhoton && !isElectron)
716  isSubLeadingEM = kTRUE;
717  }
718  }
719  }
720  }
721  }
722  }
723  fCaloPhotonMCFlags = isPhoton *1 + isElectron *2 + isConversion*4+ isConversionFullyContained *8 + isMerged *16 + isMergedPartConv*32 + isDalitz *64 + isDalitzMerged *128 + isPhotonWithElecMother *256 + isShower * 512 + isSubLeadingEM * 1024;
724 }
725 
726 //_____________________________________________________________________________
727 // prints information of MC particles contributing to cluster
728 //_____________________________________________________________________________
730  cout << endl << endl << "particles contributing: " << endl;
731  for (Int_t i =0 ; i < fNCaloPhotonMCLabels; i++ ){
732  TParticle *dummy = NULL;
733  if (fCaloPhotonMCLabels[i]>0){
734  dummy = (TParticle*)MCStack->Particle(fCaloPhotonMCLabels[i]);
735  cout << i << "\t"<< fCaloPhotonMCLabels[i] << "\t pdg code: " <<dummy->GetPdgCode() << "\t prod radius: "<< dummy->R() << "\t energy: " << dummy->Energy() ;
736  if (dummy->GetMother(0) > -1){
737  TParticle* dummyMother = (TParticle*)MCStack->Particle(dummy->GetMother(0));
738  cout << "\t mother part: " << dummy->GetMother(0) << "\t mother pdg code: " << dummyMother->GetPdgCode() << "\t energy: " << dummyMother->Energy() << endl;
739  if (dummyMother->GetMother(0) > -1){
740  TParticle* dummyGrandMother = (TParticle*)MCStack->Particle(dummyMother->GetMother(0));
741  cout << "\t grandmother part: " << dummyMother->GetMother(0) << "\t grandmother pdg code: " << dummyGrandMother->GetPdgCode() << "\t energy: " << dummyGrandMother->Energy() << endl;
742  } else {
743  cout << endl;
744  }
745  } else {
746  cout << endl;
747  }
748  }
749 // if (dummy) delete dummy;
750  }
751  cout << "mothers contributing: " << endl;
752  for (Int_t i =0 ; i < fNCaloPhotonMotherMCLabels; i++ ){
753  TParticle *dummy = NULL;
754  if (fCaloPhotonMotherMCLabels[i]>0){
755  dummy = (TParticle*)MCStack->Particle(fCaloPhotonMotherMCLabels[i]);
756  cout << i << "\t"<< fCaloPhotonMotherMCLabels[i] << "\t pdg code: " <<dummy->GetPdgCode() << "\t prod radius: "<< dummy->R() << "\t energy: " << dummy->Energy() << endl;
757 
758 
759  }
760 // if (dummy) delete dummy;
761  }
762 }
763 
764 //_____________________________________________________________________________
765 // prints information of cluster as it is flagged
766 //_____________________________________________________________________________
768  cout << fCaloPhotonMCFlags << "\t photon: " << IsLargestComponentPhoton() << "\t electron: " << IsLargestComponentElectron() << "\t conversion: " << IsConversion() << "\t conversion full: "
769  << IsConversionFullyContained() << "\t merged: " << IsMerged() << "\t merged p.conv: " << IsMergedPartConv() << "\t Dalitz: " << IsDalitz() << "\t Dalitz merged: " << IsDalitzMerged()
770  << "\t rad: " << IsPhotonWithElecMother() << "\t shower: " << IsShower() << "\t no EM lead: " << IsEMNonLeading() << "\t EM sub lead: " << IsSubLeadingEM() << endl;
771 }
virtual Double_t GetPy() const
double Double_t
Definition: External.C:58
virtual Double_t GetPz() const
AliAODConversionPhoton & operator=(const AliAODConversionPhoton &g)
int Int_t
Definition: External.C:63
void PrintCaloMCLabelsAndInfo(AliStack *MCStack)
virtual Double_t GetPx() const
void SetCaloPhotonMCFlags(AliStack *MCStack, Bool_t enableSort)
ClassImp(AliAODConversionPhoton) AliAODConversionPhoton
bool Bool_t
Definition: External.C:53
void CalculateDistanceOfClossetApproachToPrimVtx(const AliVVertex *primVertex)
void SetCaloPhotonMCFlagsAOD(AliVEvent *event, Bool_t enableSort)