AliESDCentrality replaced by AliCentrality

[u/mrichter/AliRoot.git] / PWG3 / dielectron / AliDielectronPair.cxx

diff --git a/PWG3/dielectron/AliDielectronPair.cxx b/PWG3/dielectron/AliDielectronPair.cxx
index ad6994c8e41a2e0840c85110892bb51597bb4a0a..eb1d7f07cf759a9dedcbf071b92f4fbb7775fb51 100644 (file)
--- a/PWG3/dielectron/AliDielectronPair.cxx
+++ b/PWG3/dielectron/AliDielectronPair.cxx
@@ -98,22 +98,90 @@ void AliDielectronPair::SetTracks(AliVTrack * const particle1, Int_t pid1,
   }
 }
 
+//______________________________________________
+void AliDielectronPair::GetThetaPhiCM(Double_t &thetaHE, Double_t &phiHE, Double_t &thetaCS, Double_t &phiCS) const
+{
+  //
+  // Calculate theta and phi in helicity and Collins-Soper coordinate frame
+  //
+  const Double_t kBeamEnergy   = 3500.;
+  Double_t pxyz1[3]={0,0,0};
+  Double_t pxyz2[3]={0,0,0};
+  Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+  Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+  
+  AliVParticle *d1 = static_cast<AliVParticle*>(fRefD1.GetObject());
+  AliVParticle *d2 = static_cast<AliVParticle*>(fRefD2.GetObject());
+
+  d1->PxPyPz(pxyz1);
+  d2->PxPyPz(pxyz2);
+  
+  TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+  TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+  
+  // first & second daughter 4-mom
+  TLorentzVector p1Mom(pxyz1[0],pxyz1[1],pxyz1[2],
+                       TMath::Sqrt(pxyz1[0]*pxyz1[0]+pxyz1[1]*pxyz1[1]+pxyz1[2]*pxyz1[2]+eleMass*eleMass));
+  TLorentzVector p2Mom(pxyz2[0],pxyz2[1],pxyz2[2],
+                       TMath::Sqrt(pxyz2[0]*pxyz2[0]+pxyz2[1]*pxyz2[1]+pxyz2[2]*pxyz2[2]+eleMass*eleMass));
+  // J/Psi 4-momentum vector
+  TLorentzVector motherMom=p1Mom+p2Mom;
+  
+  // boost all the 4-mom vectors to the mother rest frame
+  TVector3 beta = (-1.0/motherMom.E())*motherMom.Vect();
+  p1Mom.Boost(beta);
+  p2Mom.Boost(beta);
+  projMom.Boost(beta);
+  targMom.Boost(beta);
+
+    // x,y,z axes
+  TVector3 zAxisHE = (motherMom.Vect()).Unit();
+  TVector3 zAxisCS = ((projMom.Vect()).Unit()-(targMom.Vect()).Unit()).Unit();
+  TVector3 yAxis = ((projMom.Vect()).Cross(targMom.Vect())).Unit();
+  TVector3 xAxisHE = (yAxis.Cross(zAxisHE)).Unit();
+  TVector3 xAxisCS = (yAxis.Cross(zAxisCS)).Unit();
+  
+  // fill theta and phi
+  if(d1->Charge()>0){
+    thetaHE = zAxisHE.Dot((p1Mom.Vect()).Unit());
+    thetaCS = zAxisCS.Dot((p1Mom.Vect()).Unit());
+    phiHE   = TMath::ATan2((p1Mom.Vect()).Dot(yAxis), (p1Mom.Vect()).Dot(xAxisHE));
+    phiCS   = TMath::ATan2((p1Mom.Vect()).Dot(yAxis), (p1Mom.Vect()).Dot(xAxisCS));
+  } else {
+    thetaHE = zAxisHE.Dot((p2Mom.Vect()).Unit());
+    thetaCS = zAxisCS.Dot((p2Mom.Vect()).Unit());
+    phiHE   = TMath::ATan2((p2Mom.Vect()).Dot(yAxis), (p2Mom.Vect()).Dot(xAxisHE));
+    phiCS   = TMath::ATan2((p2Mom.Vect()).Dot(yAxis), (p2Mom.Vect()).Dot(xAxisCS));
+  }
+}
+
 //______________________________________________
 Double_t AliDielectronPair::ThetaPhiCM(const AliVParticle* d1, const AliVParticle* d2, 
-                                      const Bool_t isHE, const Bool_t isTheta) {
-  // The function calculates theta and phi in the mother rest frame with 
+                                       const Bool_t isHE, const Bool_t isTheta)
+{
+  // The function calculates theta and phi in the mother rest frame with
   // respect to the helicity coordinate system and Collins-Soper coordinate system
   // TO DO: generalize for different decays (only J/Psi->e+e- now)
 
   // Laboratory frame 4-vectors:
   // projectile beam & target beam 4-mom
-  const Double_t kBeamEnergy   = 3500.;      //TODO: need to retrieve the beam energy from somewhere
-  TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton))); 
-  TLorentzVector targMom(0.,0.,kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton))); 
+  // TODO: need to retrieve the beam energy from somewhere
+  const Double_t kBeamEnergy   = 3500.;
+  Double_t px1=d1->Px();
+  Double_t py1=d1->Py();
+  Double_t pz1=d1->Pz();
+  Double_t px2=d2->Px();
+  Double_t py2=d2->Py();
+  Double_t pz2=d2->Pz();
+  Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+  Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+  
+  TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+  TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
   
   // first & second daughter 4-mom
-  TLorentzVector p1Mom(d1->Px(),d1->Py(),d1->Pz(),TMath::Sqrt(d1->Px()*d1->Px()+d1->Py()*d1->Py()+d1->Pz()*d1->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
-  TLorentzVector p2Mom(d2->Px(),d2->Py(),d2->Pz(),TMath::Sqrt(d2->Px()*d2->Px()+d2->Py()*d2->Py()+d2->Pz()*d2->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
+  TLorentzVector p1Mom(px1,py1,pz1,TMath::Sqrt(px1*px1+py1*py1+pz1*pz1+eleMass*eleMass));
+  TLorentzVector p2Mom(px2,py2,pz2,TMath::Sqrt(px2*px2+py2*py2+pz2*pz2+eleMass*eleMass));
   // J/Psi 4-momentum vector
   TLorentzVector motherMom=p1Mom+p2Mom;
   
@@ -155,15 +223,26 @@ Double_t AliDielectronPair::ThetaPhiCM(const Bool_t isHE, const Bool_t isTheta)
 
   // Laboratory frame 4-vectors:
   // projectile beam & target beam 4-mom
-  const Double_t kBeamEnergy   = 3500.;      //TODO: need to retrieve the beam energy from somewhere
-  TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton))); 
-  TLorentzVector targMom(0.,0.,kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton))); 
-
+  AliVParticle *d1 = static_cast<AliVParticle*>(fRefD1.GetObject());
+  AliVParticle *d2 = static_cast<AliVParticle*>(fRefD2.GetObject());
+  
+  const Double_t kBeamEnergy   = 3500.;
+  Double_t px1=d1->Px();
+  Double_t py1=d1->Py();
+  Double_t pz1=d1->Pz();
+  Double_t px2=d2->Px();
+  Double_t py2=d2->Py();
+  Double_t pz2=d2->Pz();
+  Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+  Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+  
+  TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+  TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+  
+  // first & second daughter 4-mom
   // first & second daughter 4-mom
-  AliVParticle *d1 = dynamic_cast<AliVParticle*>(fRefD1.GetObject());
-  AliVParticle *d2 = dynamic_cast<AliVParticle*>(fRefD2.GetObject());
-  TLorentzVector p1Mom(d1->Px(),d1->Py(),d1->Pz(),TMath::Sqrt(d1->Px()*d1->Px()+d1->Py()*d1->Py()+d1->Pz()*d1->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
-  TLorentzVector p2Mom(d2->Px(),d2->Py(),d2->Pz(),TMath::Sqrt(d2->Px()*d2->Px()+d2->Py()*d2->Py()+d2->Pz()*d2->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
+  TLorentzVector p1Mom(px1,py1,pz1,TMath::Sqrt(px1*px1+py1*py1+pz1*pz1+eleMass*eleMass));
+  TLorentzVector p2Mom(px2,py2,pz2,TMath::Sqrt(px2*px2+py2*py2+pz2*pz2+eleMass*eleMass));
   // J/Psi 4-momentum vector
   TLorentzVector motherMom=p1Mom+p2Mom;