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. An instance of this class has also to be added to the AliAODEvent
11 // structure to provide additional information that is specific to MUON and
12 // therefore has not been included into the AOD header.
13 // Two transient data members are not stored on file as they can be recomputed
17 #include "AliAODDimuon.h"
18 #include "TLorentzVector.h"
19 #define AliAODDimuon_CXX
21 ClassImp(AliAODDimuon)
23 //______________________________________________________________________________
24 AliAODDimuon::AliAODDimuon():AliVParticle(),fP(0),fMProton(0.93827231)
26 // default constructor
31 //______________________________________________________________________________
32 AliAODDimuon::AliAODDimuon(const AliAODDimuon& dimu):
34 fP(0),fMProton(0.93827231)
41 //______________________________________________________________________________
42 AliAODDimuon &AliAODDimuon::operator=(const AliAODDimuon& dimu)
44 // assignment operator
55 //______________________________________________________________________________
56 AliAODDimuon::AliAODDimuon(TObject *mu0, TObject *mu1):
57 fP(0),fMProton(0.93827231)
59 // Creates a dimuon pair from two tracks
65 //______________________________________________________________________________
66 AliAODDimuon::~AliAODDimuon()
72 //______________________________________________________________________________
73 Double_t AliAODDimuon::Px() const {
75 if(CheckPointers())return -999999999;
76 return ((AliAODTrack*)fMu[0].GetObject())->Px()+
77 ((AliAODTrack*)fMu[1].GetObject())->Px();
80 //______________________________________________________________________________
81 Double_t AliAODDimuon::Py() const {
83 if(CheckPointers())return -999999999;
84 return ((AliAODTrack*)fMu[0].GetObject())->Py()+
85 ((AliAODTrack*)fMu[1].GetObject())->Py();
88 //______________________________________________________________________________
89 Double_t AliAODDimuon::Pz() const {
91 if(CheckPointers())return -999999999;
92 return ((AliAODTrack*)fMu[0].GetObject())->Pz()+
93 ((AliAODTrack*)fMu[1].GetObject())->Pz();
96 //______________________________________________________________________________
97 Double_t AliAODDimuon::Pt() const {
99 if(CheckPointers())return -999999999;
102 return TMath::Sqrt(px*px+py*py);
105 //______________________________________________________________________________
106 Double_t AliAODDimuon::E() const
110 if(CheckPointers())return -999999999;
112 return ((AliAODTrack*)fMu[0].GetObject())->E()+ ((AliAODTrack*)fMu[1].GetObject())->E();
115 //______________________________________________________________________________
116 Double_t AliAODDimuon::P() const {
118 if(CheckPointers())return -999999999;
122 //______________________________________________________________________________
123 Double_t AliAODDimuon::M() const {
124 // Dimuon invariant mass
125 if(CheckPointers())return -999999999;
129 //______________________________________________________________________________
130 Double_t AliAODDimuon::Eta() const {
131 // Dimuon pseudorapidity
132 if(CheckPointers())return -999999999;
136 //______________________________________________________________________________
137 Double_t AliAODDimuon::Phi() const {
138 // Dimuon asimuthal angle
139 if(CheckPointers())return -999999999;
144 //______________________________________________________________________________
145 Double_t AliAODDimuon::Theta() const {
146 // Dimuon polar angle
147 if(CheckPointers())return -999999999;
148 return TLV()->Theta();
151 //______________________________________________________________________________
152 Double_t AliAODDimuon::Y() const {
154 if(CheckPointers())return -999999999;
155 return TLV()->Rapidity();
158 //______________________________________________________________________________
159 Short_t AliAODDimuon::Charge() const {
161 if(CheckPointers())return -999;
162 return ((AliAODTrack*)fMu[0].GetObject())->Charge()+
163 ((AliAODTrack*)fMu[1].GetObject())->Charge();
166 //______________________________________________________________________________
167 Int_t AliAODDimuon::CheckPointers() const{
168 // Checks if the track pointers have been initialized
169 if(fMu[0]==0||fMu[1]==0){
170 printf("Dimuon not initialized\n");
173 if((fMu[0].GetObject())==0||(fMu[1].GetObject())==0){
174 printf("Can not get objects\n");
180 //______________________________________________________________________________
181 void AliAODDimuon::SetMu(Int_t imu, AliAODTrack *mu){
182 // Assign a track pointer
190 //______________________________________________________________________________
191 void AliAODDimuon::SetMuons(AliAODTrack *mu0, AliAODTrack *mu1){
192 // Assign the track pointers
199 //______________________________________________________________________________
200 Double_t AliAODDimuon::XF() {
203 if(CheckPointers())return -999999999;
205 //Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
206 Double_t ebeam = 3500.; // temporary
208 printf("AliAODDimuon::xf: can not compute xf with EBeam=%f\n",ebeam);
212 Double_t pMax=TMath::Sqrt(ebeam*ebeam-mDimu*mDimu);
216 //______________________________________________________________________________
217 Double_t AliAODDimuon::CostCS(){
218 // Calculation the Collins-Soper angle (adapted from code by R. Arnaldi)
219 if(CheckPointers())return -999999999;
220 Double_t ebeam=3500.; //temporary
222 printf("Can not compute costCS with EBeam=%f\n",ebeam);
225 Double_t mp=fMProton;
226 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
227 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
228 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
229 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
230 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
231 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
232 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
233 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
234 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
235 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
236 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
238 // Fill the Lorentz vector for projectile and target
239 // For the moment we do not consider the crossing angle
240 // Projectile runs towards the MUON arm
241 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
242 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
244 // --- Get the muons parameters in the CM frame
246 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
247 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
249 // --- Obtain the dimuon parameters in the CM frame
251 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
253 // --- Translate the dimuon parameters in the dimuon rest frame
255 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
256 TLorentzVector pMu1Dimu=pMu1CM;
257 TLorentzVector pMu2Dimu=pMu2CM;
258 TLorentzVector pProjDimu=pProjCM;
259 TLorentzVector pTargDimu=pTargCM;
260 pMu1Dimu.Boost(beta);
261 pMu2Dimu.Boost(beta);
262 pProjDimu.Boost(beta);
263 pTargDimu.Boost(beta);
265 // --- Determine the z axis for the CS angle
267 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
269 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
273 cost = zaxisCS.Dot((pMu1Dimu.Vect()).Unit());
274 // Theta CS is not properly defined for Like-Sign muons
275 if(mu2Charge > 0 && cost<0) cost=-cost;
277 // Theta CS is not properly defined for Like-Sign muons
278 cost = zaxisCS.Dot((pMu2Dimu.Vect()).Unit());
279 if(mu2Charge < 0 && cost<0) cost=-cost;
284 //______________________________________________________________________________
285 Double_t AliAODDimuon::CostHe(){
286 // Calculation the Helicity polarization angle (adapted from code by R. Arnaldi)
287 if(CheckPointers())return -999999999;
288 Double_t ebeam=3500; //temporary
290 printf("Can not compute costCS with EBeam=%f\n",ebeam);
293 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
294 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
295 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
296 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
297 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
298 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
299 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
300 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
301 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
302 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
304 // --- Get the muons parameters in the CM frame
306 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
307 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
309 // --- Obtain the dimuon parameters in the CM frame
311 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
313 // --- Translate the dimuon parameters in the dimuon rest frame
315 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
316 TLorentzVector pMu1Dimu=pMu1CM;
317 TLorentzVector pMu2Dimu=pMu2CM;
318 pMu1Dimu.Boost(beta);
319 pMu2Dimu.Boost(beta);
321 // --- Determine the z axis for the calculation of the polarization angle
322 // (i.e. the direction of the dimuon in the CM system)
325 zaxis=(pDimuCM.Vect()).Unit();
327 // --- Calculation of the polarization angle (Helicity)
328 // (angle between mu+ and the z axis defined above)
332 cost = zaxis.Dot((pMu1Dimu.Vect()).Unit());
333 // Theta Helicity is not properly defined for Like-Sign muons
334 if(mu2Charge > 0 && cost<0) cost=-cost;
336 cost = zaxis.Dot((pMu2Dimu.Vect()).Unit());
337 // Theta Helicity is not properly defined for Like-Sign muons
338 if(mu2Charge < 0 && cost<0) cost=-cost;
343 //________________________________________________________________________
344 Double_t AliAODDimuon::PhiCS(){
345 // Cosinus of the Collins-Soper polar decay angle
346 if(CheckPointers())return -999999999;
347 Double_t ebeam=3500.; //temporary
349 printf("Can not compute phiCS with EBeam=%f\n",ebeam);
352 Double_t mp=fMProton;
353 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
354 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
355 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
356 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
357 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
358 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
359 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
360 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
361 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
362 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
364 // Fill the Lorentz vector for projectile and target
365 // For the moment we do not consider the crossing angle
366 // Projectile runs towards the MUON arm
367 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
368 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
370 // --- Get the muons parameters in the CM frame
372 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
373 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
375 // --- Obtain the dimuon parameters in the CM frame
377 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
379 // --- Translate the dimuon parameters in the dimuon rest frame
381 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
382 TLorentzVector pMu1Dimu=pMu1CM;
383 TLorentzVector pMu2Dimu=pMu2CM;
384 TLorentzVector pProjDimu=pProjCM;
385 TLorentzVector pTargDimu=pTargCM;
386 pMu1Dimu.Boost(beta);
387 pMu2Dimu.Boost(beta);
388 pProjDimu.Boost(beta);
389 pTargDimu.Boost(beta);
391 // --- Determine the z axis for the CS angle
393 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
395 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
397 TVector3 yaxisCS=(((pProjDimu.Vect()).Unit()).Cross((pTargDimu.Vect()).Unit())).Unit();
398 TVector3 xaxisCS=(yaxisCS.Cross(zaxisCS)).Unit();
401 if(mu1Charge>0) phi = TMath::ATan2((pMu1Dimu.Vect()).Dot(yaxisCS),((pMu1Dimu.Vect()).Dot(xaxisCS)));
402 else phi = TMath::ATan2((pMu2Dimu.Vect()).Dot(yaxisCS),((pMu2Dimu.Vect()).Dot(xaxisCS)));
407 //______________________________________________________________________________
408 Double_t AliAODDimuon::PhiHe(){
409 // Calculation the Helicity aimuthal angle (adapted from code by R. Arnaldi)
410 if(CheckPointers())return -999999999;
411 Double_t ebeam=3500; //temporary
413 printf("Can not compute phiHE with EBeam=%f\n",ebeam);
416 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-fMProton*fMProton);
417 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
418 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
419 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
420 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
421 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
422 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
423 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
424 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
425 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
427 // Fill the Lorentz vector for projectile and target
428 // For the moment we consider no crossing angle
429 // Projectile runs towards the MUON arm
430 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
431 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
433 // --- Get the muons parameters in the CM frame
435 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
436 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
438 // --- Obtain the dimuon parameters in the CM frame
440 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
442 // --- Translate the muon parameters in the dimuon rest frame
444 TVector3 zaxis=(pDimuCM.Vect()).Unit();
446 // --- Translate the dimuon parameters in the dimuon rest frame
448 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
449 TLorentzVector pMu1Dimu=pMu1CM;
450 TLorentzVector pMu2Dimu=pMu2CM;
451 pMu1Dimu.Boost(beta);
452 pMu2Dimu.Boost(beta);
454 TLorentzVector pProjDimu=pProjCM;
455 TLorentzVector pTargDimu=pTargCM;
456 pProjDimu.Boost(beta);
457 pTargDimu.Boost(beta);
459 TVector3 yaxis=((pProjDimu.Vect()).Cross(pTargDimu.Vect())).Unit();
460 TVector3 xaxis=(yaxis.Cross(zaxis)).Unit();
462 // --- Calculation of the azimuthal angle (Helicity)
465 if(mu1Charge>0) phi = TMath::ATan2((pMu1Dimu.Vect()).Dot(yaxis),(pMu1Dimu.Vect()).Dot(xaxis));
466 else phi = TMath::ATan2((pMu2Dimu.Vect()).Dot(yaxis),(pMu2Dimu.Vect()).Dot(xaxis));
471 //______________________________________________________________________________
472 Int_t AliAODDimuon::AnyPt(){
473 // Test if the two muons match two trigger tracks
474 if(CheckPointers())return 0;
475 return (((AliAODTrack*)fMu[0].GetObject())->MatchTrigger())&&
476 (((AliAODTrack*)fMu[1].GetObject())->MatchTrigger());
479 //______________________________________________________________________________
480 Int_t AliAODDimuon::LowPt(){
481 // Test if the two muons match two trigger tracks with a "Low Pt" cut
482 if(CheckPointers())return 0;
483 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt())&&
484 (((AliAODTrack*)fMu[1].GetObject())->MatchTriggerLowPt());
487 //______________________________________________________________________________
488 Int_t AliAODDimuon::HighPt(){
489 // Test if the two muons match two trigger tracks with a "High Pt" cut
490 if(CheckPointers())return 0;
491 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt())&&
492 (((AliAODTrack*)fMu[1].GetObject())->MatchTriggerHighPt());
495 //______________________________________________________________________________
496 Double_t AliAODDimuon::MaxChi2Match(){
497 // Maximum matching Chi2 between track and trigger track
498 if(CheckPointers())return -999999999;
499 return TMath::Max((((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()),
500 (((AliAODTrack*)fMu[1].GetObject())->GetChi2MatchTrigger()));
503 //______________________________________________________________________________
504 TLorentzVector* AliAODDimuon::TLV() const
507 fP = new TLorentzVector(Px(),Py(),Pz(),E());