1 #include "AliAODConversionMother.h"
2 #include "AliKFConversionMother.h"
3 #include "AliKFParticle.h"
6 // Author D. Lohner (Daniel.Lohner@cern.ch)
10 ClassImp(AliAODConversionMother)
12 AliAODConversionMother::AliAODConversionMother() :
13 AliAODConversionParticle(),
19 fdcaBetweenPhotons(1),
35 AliAODConversionMother::AliAODConversionMother(AliKFConversionMother *kf):
36 AliAODConversionParticle(),
37 fMCLabel(kf->GetMCLabel()),
39 fOpeningAngle(kf->GetOpeningAngle()),
40 fAlpha(kf->GetAlpha()),
42 fdcaBetweenPhotons(100),
48 SetPxPyPzE(kf->GetPx(),kf->GetPy(),kf->GetPz(),kf->GetE());
54 //Set Decay Photon Labels
55 fLabel[0]=kf->GetGammaLabel(0);
56 fLabel[1]=kf->GetGammaLabel(1);
60 AliAODConversionMother::AliAODConversionMother(AliAODConversionPhoton *y1,AliAODConversionPhoton *y2):
61 AliAODConversionParticle(),
67 fdcaBetweenPhotons(1),
73 SetPxPyPzE(y1->Px()+y2->Px(),y1->Py()+y2->Py(),y1->Pz()+y2->Pz(),y1->E()+y2->E());
75 // Calculate Opening Angle
76 TVector3 v1(y1->Px(),y1->Py(),y1->Pz());
77 TVector3 v2(y2->Px(),y2->Py(),y2->Pz());
78 fOpeningAngle=v1.Angle(v2);
79 fdcaBetweenPhotons = CalculateDistanceBetweenPhotons(y1,y2,fProductionVtx);
80 DetermineMesonQuality(y1,y2);
82 if((y1->E()+y2->E()) != 0){
83 fAlpha=TMath::Abs((y1->E()-y2->E())/(y1->E()+y2->E()));
86 // Set Chi2 to the mean chi2 of gammas
87 // fChi2=0.5*(y1->GetChi2perNDF()+y2->GetChi2perNDF());
89 //Set Decay Photon Labels
95 AliAODConversionMother::~AliAODConversionMother() {
96 // empty standard destructor
99 TParticle *AliAODConversionMother::GetMCParticle(AliStack *fMCStack){
100 if(!fMCStack){AliError("MC Stack not defined");return 0x0;}
103 return fMCStack->Particle(fMCLabel);
108 Bool_t AliAODConversionMother::IsTrueMeson(AliStack *fMCStack,Int_t pdgcode){
109 TParticle *part=GetMCParticle(fMCStack);
112 // Check if it is a true photon
113 if(part->GetPdgCode()==pdgcode){
120 Float_t AliAODConversionMother::CalculateDistanceBetweenPhotons(AliAODConversionPhoton* y1, AliAODConversionPhoton* y2 , Double_t prodPoint[3]){
122 TVector3 a(y1->GetConversionX(),y1->GetConversionY(),y1->GetConversionZ());
123 TVector3 b(y1->GetPx(),y1->GetPy(),y1->GetPz());
124 TVector3 c(y2->GetConversionX(),y2->GetConversionY(),y2->GetConversionZ());
125 TVector3 d(y2->GetPx(),y2->GetPy(),y2->GetPz());
127 TVector3 n = b.Cross(d);
128 TVector3 nn = n.Unit();
134 dist = TMath::Abs((e.Cross(d)).Mag())/TMath::Abs(d.Mag());
140 dist = TMath::Abs(n.Dot(c-a))/TMath::Abs(n.Mag());
141 Double_t lambda = (b.Dot(d) * (a-c).Dot(d) - d.Dot(d) * (a-c).Dot(b))/(b.Dot(b) * d.Dot(d) - TMath::Power(b.Dot(d),2));
142 Double_t mu = ((a-c).Dot(d) * b.Dot(b) - (a-c).Dot(b) * b.Dot(d) )/(b.Dot(b) * d.Dot(d) - TMath::Power(b.Dot(d),2));
144 TVector3 S1 = a + lambda* b;
145 TVector3 S2 = c + mu* d;
146 TVector3 Prod = S1 + 0.5*dist*(S2-S1).Unit();
147 prodPoint[0] = Prod(0);
148 prodPoint[1] = Prod(1);
149 prodPoint[2] = Prod(2);
152 if (dist > 1000) dist = 999.;
156 ///________________________________________________________________________
157 void AliAODConversionMother::CalculateDistanceOfClossetApproachToPrimVtx(const AliVVertex* primVertex){
159 Double_t primCo[3] = {primVertex->GetX(),primVertex->GetY(),primVertex->GetZ()};
161 Double_t absoluteP = TMath::Sqrt(TMath::Power(Px(),2) + TMath::Power(Py(),2) + TMath::Power(Pz(),2));
162 Double_t p[3] = {Px()/absoluteP,Py()/absoluteP,Pz()/absoluteP};
165 CP[0] = fProductionVtx[0] - primCo[0];
166 CP[1] = fProductionVtx[1] - primCo[1];
167 CP[2] = fProductionVtx[2] - primCo[2];
169 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]);
172 S[0] = fProductionVtx[0] + p[0]*Lambda;
173 S[1] = fProductionVtx[1] + p[1]*Lambda;
174 S[2] = fProductionVtx[2] + p[2]*Lambda;
176 fdcaRPrimVtx = TMath::Sqrt( TMath::Power(primCo[0]-S[0],2) + TMath::Power(primCo[1]-S[1],2));
177 fdcaZPrimVtx = primCo[2]-S[2];
180 // cout << "DCA z: " << dca[1] << "\t DCA r: " << dca[0] << "\t DCA 3d: " << TMath::Sqrt(dca[1]*dca[1] + dca[0]*dca[0]) << endl;
186 ///________________________________________________________________________
187 void AliAODConversionMother::DetermineMesonQuality(AliAODConversionPhoton* y1, AliAODConversionPhoton* y2){
188 UChar_t photonQA1 = y1->GetPhotonQuality();
189 UChar_t photonQA2 = y2->GetPhotonQuality();
191 if (photonQA1 == 0 || photonQA2 == 0){
195 if (photonQA1 == 1 && photonQA2 == 1){
199 if (photonQA1 == 2 && photonQA2 == 2){
203 if (photonQA1 == 3 && photonQA2 == 3){
227 if ((photonQA1 == 2 && photonQA2 == 3)|| (photonQA1 == 3 && photonQA2 == 2)){