1 // AliAODDimuon: a class for AODs for the MUON Arm of the ALICE Experiment
2 // Author: P. Cortese, Universita' del Piemonte Orientale in Alessandria and
3 // INFN of Torino - Italy
5 // The class defines a dimuon pair object from two AliAODTrack objects.
6 // AliAODDimuon objects are supposed to be added to the AliAODEvent structure
7 // during analysis. They would then allow to calculate the dimuon-related
8 // kinematic variables with a minimal disk occupancy.
9 // The payload of the class has been reduced to two pointers to the two
10 // tracks with the addition of a pointer to the AliAODEventInfo. An instance of
11 // this class has also to be added to the AliAODEvent structure to provide
12 // additional information that is specific to MUON and therefore has not been
13 // included into the AOD header.
14 // Two transient data members are not stored on file as they can be recomputed
18 #include "AliAODDimuon.h"
19 #include "TLorentzVector.h"
20 #define AliAODDimuon_CXX
22 ClassImp(AliAODDimuon)
24 //______________________________________________________________________________
25 AliAODDimuon::AliAODDimuon():AliVParticle(),fEi(0),fP(0),fMProton(0.93827231)
27 // default constructor
32 //______________________________________________________________________________
33 AliAODDimuon::AliAODDimuon(const AliAODDimuon& dimu):
35 fEi(dimu.fEi),fP(0),fMProton(0.93827231)
42 //______________________________________________________________________________
43 AliAODDimuon &AliAODDimuon::operator=(const AliAODDimuon& dimu)
45 // assignment operator
56 //______________________________________________________________________________
57 AliAODDimuon::AliAODDimuon(TObject *mu0, TObject *mu1, TObject *ei):
58 fEi(ei),fP(0),fMProton(0.93827231)
60 // Creates a dimuon pair from two tracks and the EventInfo
62 //printf("Creating dimuon from %p %p\n",mu0,mu1);
67 //______________________________________________________________________________
68 AliAODDimuon::~AliAODDimuon()
75 //______________________________________________________________________________
76 void AliAODDimuon::BookP(){
77 // Fills the dimuon momentum if not filled yet
78 static UInt_t unID[2]={0,0};
80 fP=new TLorentzVector(Px(),Py(),Pz(),E());
81 unID[0]=fMu[0].GetUniqueID();
82 unID[1]=fMu[1].GetUniqueID();
84 // For efficiency reasons
85 if((unID[0]!=fMu[0].GetUniqueID())||(unID[1]!=fMu[1].GetUniqueID())){
86 fP->SetPxPyPzE(Px(),Py(),Pz(),E());
87 unID[0]=fMu[0].GetUniqueID();
88 unID[1]=fMu[1].GetUniqueID();
92 //______________________________________________________________________________
93 Double_t AliAODDimuon::Px() const {
95 if(this->CheckPointers())return -999999999;
96 return ((AliAODTrack*)fMu[0].GetObject())->Px()+
97 ((AliAODTrack*)fMu[1].GetObject())->Px();
100 //______________________________________________________________________________
101 Double_t AliAODDimuon::Py() const {
103 if(this->CheckPointers())return -999999999;
104 return ((AliAODTrack*)fMu[0].GetObject())->Py()+
105 ((AliAODTrack*)fMu[1].GetObject())->Py();
108 //______________________________________________________________________________
109 Double_t AliAODDimuon::Pz() const {
111 if(this->CheckPointers())return -999999999;
112 return ((AliAODTrack*)fMu[0].GetObject())->Pz()+
113 ((AliAODTrack*)fMu[1].GetObject())->Pz();
116 //______________________________________________________________________________
117 Double_t AliAODDimuon::Pt() const {
119 if(this->CheckPointers())return -999999999;
122 return TMath::Sqrt(px*px+py*py);
126 //______________________________________________________________________________
127 Double_t AliAODDimuon::E() const {
129 if(this->CheckPointers())return -999999999;
130 return ((AliAODTrack*)fMu[0].GetObject())->E()+
131 ((AliAODTrack*)fMu[1].GetObject())->E();
134 //______________________________________________________________________________
135 Double_t AliAODDimuon::P() const {
136 // This is just to override the virtual function
137 printf("You should never call: Double_t AliAODDimuon::P() const\n");
141 //______________________________________________________________________________
142 Double_t AliAODDimuon::P() {
144 if(this->CheckPointers())return -999999999;
149 //______________________________________________________________________________
150 Double_t AliAODDimuon::M() const {
151 // This is just to override the virtual function
152 printf("You should never call: Double_t AliAODDimuon::M() const\n");
156 //______________________________________________________________________________
157 Double_t AliAODDimuon::M() {
158 // Dimuon invariant mass
159 if(this->CheckPointers())return -999999999;
164 //______________________________________________________________________________
165 Double_t AliAODDimuon::Mass() {
166 // Dimuon invariant mass
167 if(this->CheckPointers())return -999999999;
172 //______________________________________________________________________________
173 Double_t AliAODDimuon::Eta() const {
174 // This is just to override the virtual function
175 printf("You should never call: Double_t AliAODDimuon::Eta() const\n");
179 //______________________________________________________________________________
180 Double_t AliAODDimuon::Eta() {
181 // Dimuon pseudorapidity
182 if(this->CheckPointers())return -999999999;
187 //______________________________________________________________________________
188 Double_t AliAODDimuon::Phi() const {
189 // This is just to override the virtual function
190 printf("You should never call: Double_t AliAODDimuon::Phi() const\n");
194 //______________________________________________________________________________
195 Double_t AliAODDimuon::Phi() {
196 // Dimuon asimuthal angle
197 if(this->CheckPointers())return -999999999;
201 //______________________________________________________________________________
202 Double_t AliAODDimuon::Theta() const {
203 // This is just to override the virtual function
204 printf("You should never call: Double_t AliAODDimuon::Theta() const\n");
208 //______________________________________________________________________________
209 Double_t AliAODDimuon::Theta() {
210 // Dimuon polar angle
211 if(this->CheckPointers())return -999999999;
216 //______________________________________________________________________________
217 Double_t AliAODDimuon::Y() const {
218 // This is just to override the virtual function
219 printf("You should never call: Double_t AliAODDimuon::Y() const\n");
223 //______________________________________________________________________________
224 Double_t AliAODDimuon::Y() {
226 if(this->CheckPointers())return -999999999;
228 return fP->Rapidity();
231 //______________________________________________________________________________
232 Short_t AliAODDimuon::Charge() const {
234 if(this->CheckPointers())return -999;
235 return ((AliAODTrack*)fMu[0].GetObject())->Charge()+
236 ((AliAODTrack*)fMu[1].GetObject())->Charge();
239 //______________________________________________________________________________
240 Int_t AliAODDimuon::CheckPointers() const{
241 // Checks if the track pointers have been initialized
242 if(fMu[0]==0||fMu[1]==0){
243 printf("Dimuon not initialized\n");
246 if((fMu[0].GetObject())==0||(fMu[1].GetObject())==0){
247 printf("Can not get objects. Got: %p %p\n",fMu[0].GetObject(),fMu[1].GetObject());
253 //______________________________________________________________________________
254 void AliAODDimuon::SetMu(Int_t imu, AliAODTrack *mu){
255 // Assign a track pointer
261 //______________________________________________________________________________
262 void AliAODDimuon::SetMuons(AliAODTrack *mu0, AliAODTrack *mu1){
263 // Assign the track pointers
268 //______________________________________________________________________________
269 Double_t AliAODDimuon::XF() {
271 Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
273 printf("AliAODDimuon::xf: can not compute xf with EBeam=%f\n",ebeam);
276 if(this->CheckPointers())return -999999999;
279 Double_t pMax=TMath::Sqrt(ebeam*ebeam-mDimu*mDimu);
283 //______________________________________________________________________________
284 // Calculation the Collins-Soper angle (adapted from code by R. Arnaldi)
285 Double_t AliAODDimuon::CostCS(){
286 // Cosinus of the Collins-Soper polar decay angle
287 if(CheckPointers())return -999999999;
289 printf("Pointer to MuonHeader not initialized\n");
292 if(fEi.GetObject()==0){
293 printf("Can not get MuonHeader object\n");
296 Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
298 printf("Can not compute costCS with EBeam=%f\n",ebeam);
301 Double_t mp=fMProton;
302 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
303 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
304 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
305 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
306 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
307 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
308 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
309 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
310 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
311 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
312 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
314 // Fill the Lorentz vector for projectile and target
315 // For the moment we do not consider the crossing angle
316 // Projectile runs towards the MUON arm
317 TLorentzVector pProjLab(0.,0.,-pbeam,ebeam); // projectile
318 TLorentzVector pTargLab(0.,0., pbeam,ebeam); // target
320 // --- Get the muons parameters in the LAB frame
322 TLorentzVector pMu1Lab(pla10,pla11,pla12,e1);
323 TLorentzVector pMu2Lab(pla20,pla21,pla22,e2);
325 // --- Obtain the dimuon parameters in the LAB frame
327 TLorentzVector pDimuLab=pMu1Lab+pMu2Lab;
329 // --- Translate the dimuon parameters in the dimuon rest frame
331 TVector3 beta=(-1./pDimuLab.E())*pDimuLab.Vect();
332 TLorentzVector pMu1Dimu=pMu1Lab;
333 TLorentzVector pMu2Dimu=pMu2Lab;
334 TLorentzVector pProjDimu=pProjLab;
335 TLorentzVector pTargDimu=pTargLab;
336 pMu1Dimu.Boost(beta);
337 pMu2Dimu.Boost(beta);
338 pProjDimu.Boost(beta);
339 pTargDimu.Boost(beta);
341 // --- Determine the z axis for the CS angle
343 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
345 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
349 cost = zaxisCS.Dot((pMu1Dimu.Vect()).Unit());
350 // Theta CS is not properly defined for Like-Sign muons
351 if(mu2Charge > 0 && cost<0) cost=-cost;
353 // Theta CS is not properly defined for Like-Sign muons
354 cost = zaxisCS.Dot((pMu2Dimu.Vect()).Unit());
355 if(mu2Charge < 0 && cost<0) cost=-cost;
360 //______________________________________________________________________________
361 // Calculation the Helicity polarization angle (adapted from code by R. Arnaldi)
362 Double_t AliAODDimuon::CostHe(){
363 // Cosinus of the polar decay angle in the Helicity reference frame
364 if(CheckPointers())return -999999999;
366 printf("Pointer to MuonHeader not initialized\n");
369 if(fEi.GetObject()==0){
370 printf("Can not get MuonHeader object\n");
373 Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
375 printf("Can not compute costCS with EBeam=%f\n",ebeam);
378 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-fMProton*fMProton);
379 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
380 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
381 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
382 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
383 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
384 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
385 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
386 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
387 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
388 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
390 // Fill the Lorentz vector for projectile and target
391 // For the moment we consider no crossing angle
392 // Projectile runs towards the MUON arm
393 TLorentzVector pProjLab(0.,0.,-pbeam,ebeam); // projectile
394 TLorentzVector pTargLab(0.,0., pbeam,ebeam); // target
396 // --- Get the muons parameters in the LAB frame
398 TLorentzVector pMu1Lab(pla10,pla11,pla12,e1);
399 TLorentzVector pMu2Lab(pla20,pla21,pla22,e2);
401 // --- Obtain the dimuon parameters in the LAB frame
403 TLorentzVector pDimuLab=pMu1Lab+pMu2Lab;
405 // --- Translate the dimuon parameters in the dimuon rest frame
407 TVector3 beta=(-1./pDimuLab.E())*pDimuLab.Vect();
408 TLorentzVector pMu1Dimu=pMu1Lab;
409 TLorentzVector pMu2Dimu=pMu2Lab;
410 pMu1Dimu.Boost(beta);
411 pMu2Dimu.Boost(beta);
413 // --- Translate the dimuon parameters in the CM frame
415 TLorentzVector pDimuCM; //CM frame
417 beta2=(-1./(fMProton+pProjLab.E()))*pProjLab.Vect();
419 pDimuCM.Boost(beta2);
421 // --- Determine the z axis for the calculation of the polarization angle
422 // (i.e. the direction of the dimuon in the CM system)
425 zaxis=(pDimuCM.Vect()).Unit();
427 // --- Calculation of the polarization angle (Helicity)
428 // (angle between mu+ and the z axis defined above)
432 cost = zaxis.Dot((pMu1Dimu.Vect()).Unit());
433 // Theta Helicity is not properly defined for Like-Sign muons
434 if(mu2Charge > 0 && cost<0) cost=-cost;
436 cost = zaxis.Dot((pMu2Dimu.Vect()).Unit());
437 // Theta Helicity is not properly defined for Like-Sign muons
438 if(mu2Charge < 0 && cost<0) cost=-cost;
443 //______________________________________________________________________________
444 Int_t AliAODDimuon::AnyPt(){
445 // Test if the two muons match two trigger tracks
446 if(this->CheckPointers())return 0;
447 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerAnyPt())&&
448 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerAnyPt());
451 //______________________________________________________________________________
452 Int_t AliAODDimuon::LowPt(){
453 // Test if the two muons match two trigger tracks with a "Low Pt" cut
454 if(this->CheckPointers())return 0;
455 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt())&&
456 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt());
459 //______________________________________________________________________________
460 Int_t AliAODDimuon::HighPt(){
461 // Test if the two muons match two trigger tracks with a "High Pt" cut
462 if(this->CheckPointers())return 0;
463 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt())&&
464 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt());
467 //______________________________________________________________________________
468 Double_t AliAODDimuon::MaxChi2Match(){
469 // Maximum matching Chi2 between track and trigger track
470 if(this->CheckPointers())return -999999999;
471 return TMath::Max((((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()),
472 (((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()));