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
54 //______________________________________________________________________________
55 AliAODDimuon::AliAODDimuon(TObject *mu0, TObject *mu1):
56 fP(0),fMProton(0.93827231)
58 // Creates a dimuon pair from two tracks
64 //______________________________________________________________________________
65 AliAODDimuon::~AliAODDimuon()
72 //______________________________________________________________________________
73 void AliAODDimuon::BookP(){
75 // build the TLorentz vector
77 fP=new TLorentzVector(Px(),Py(),Pz(),E());
78 fP->SetPxPyPzE(Px(),Py(),Pz(),E());
81 //______________________________________________________________________________
82 Double_t AliAODDimuon::Px() const {
84 if(this->CheckPointers())return -999999999;
85 return ((AliAODTrack*)fMu[0].GetObject())->Px()+
86 ((AliAODTrack*)fMu[1].GetObject())->Px();
89 //______________________________________________________________________________
90 Double_t AliAODDimuon::Py() const {
92 if(this->CheckPointers())return -999999999;
93 return ((AliAODTrack*)fMu[0].GetObject())->Py()+
94 ((AliAODTrack*)fMu[1].GetObject())->Py();
97 //______________________________________________________________________________
98 Double_t AliAODDimuon::Pz() const {
100 if(this->CheckPointers())return -999999999;
101 return ((AliAODTrack*)fMu[0].GetObject())->Pz()+
102 ((AliAODTrack*)fMu[1].GetObject())->Pz();
105 //______________________________________________________________________________
106 Double_t AliAODDimuon::Pt() const {
108 if(this->CheckPointers())return -999999999;
111 return TMath::Sqrt(px*px+py*py);
115 //______________________________________________________________________________
116 Double_t AliAODDimuon::E() const {
118 if(this->CheckPointers())return -999999999;
119 return ((AliAODTrack*)fMu[0].GetObject())->E()+
120 ((AliAODTrack*)fMu[1].GetObject())->E();
123 //______________________________________________________________________________
124 Double_t AliAODDimuon::P() const {
125 // This is just to override the virtual function
126 printf("You should never call: Double_t AliAODDimuon::P() const\n");
130 //______________________________________________________________________________
131 Double_t AliAODDimuon::P() {
133 if(this->CheckPointers())return -999999999;
138 //______________________________________________________________________________
139 Double_t AliAODDimuon::M() const {
140 // This is just to override the virtual function
141 printf("You should never call: Double_t AliAODDimuon::M() const\n");
145 //______________________________________________________________________________
146 Double_t AliAODDimuon::M() {
147 // Dimuon invariant mass
148 if(this->CheckPointers())return -999999999;
153 //______________________________________________________________________________
154 Double_t AliAODDimuon::Mass() {
155 // Dimuon invariant mass
156 if(this->CheckPointers())return -999999999;
161 //______________________________________________________________________________
162 Double_t AliAODDimuon::Eta() const {
163 // This is just to override the virtual function
164 printf("You should never call: Double_t AliAODDimuon::Eta() const\n");
168 //______________________________________________________________________________
169 Double_t AliAODDimuon::Eta() {
170 // Dimuon pseudorapidity
171 if(this->CheckPointers())return -999999999;
176 //______________________________________________________________________________
177 Double_t AliAODDimuon::Phi() const {
178 // This is just to override the virtual function
179 printf("You should never call: Double_t AliAODDimuon::Phi() const\n");
183 //______________________________________________________________________________
184 Double_t AliAODDimuon::Phi() {
185 // Dimuon asimuthal angle
186 if(this->CheckPointers())return -999999999;
190 //______________________________________________________________________________
191 Double_t AliAODDimuon::Theta() const {
192 // This is just to override the virtual function
193 printf("You should never call: Double_t AliAODDimuon::Theta() const\n");
197 //______________________________________________________________________________
198 Double_t AliAODDimuon::Theta() {
199 // Dimuon polar angle
200 if(this->CheckPointers())return -999999999;
205 //______________________________________________________________________________
206 Double_t AliAODDimuon::Y() const {
207 // This is just to override the virtual function
208 printf("You should never call: Double_t AliAODDimuon::Y() const\n");
212 //______________________________________________________________________________
213 Double_t AliAODDimuon::Y() {
215 if(this->CheckPointers())return -999999999;
217 return fP->Rapidity();
220 //______________________________________________________________________________
221 Short_t AliAODDimuon::Charge() const {
223 if(this->CheckPointers())return -999;
224 return ((AliAODTrack*)fMu[0].GetObject())->Charge()+
225 ((AliAODTrack*)fMu[1].GetObject())->Charge();
228 //______________________________________________________________________________
229 Int_t AliAODDimuon::CheckPointers() const{
230 // Checks if the track pointers have been initialized
231 if(fMu[0]==0||fMu[1]==0){
232 printf("Dimuon not initialized\n");
235 if((fMu[0].GetObject())==0||(fMu[1].GetObject())==0){
236 printf("Can not get objects\n");
242 //______________________________________________________________________________
243 void AliAODDimuon::SetMu(Int_t imu, AliAODTrack *mu){
244 // Assign a track pointer
250 //______________________________________________________________________________
251 void AliAODDimuon::SetMuons(AliAODTrack *mu0, AliAODTrack *mu1){
252 // Assign the track pointers
257 //______________________________________________________________________________
258 Double_t AliAODDimuon::XF() {
260 //Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
261 Double_t ebeam = 3500.; // temporary
263 printf("AliAODDimuon::xf: can not compute xf with EBeam=%f\n",ebeam);
266 if(this->CheckPointers())return -999999999;
269 Double_t pMax=TMath::Sqrt(ebeam*ebeam-mDimu*mDimu);
273 //______________________________________________________________________________
274 Double_t AliAODDimuon::CostCS(){
275 // Calculation the Collins-Soper angle (adapted from code by R. Arnaldi)
276 if(CheckPointers())return -999999999;
277 Double_t ebeam=3500.; //temporary
279 printf("Can not compute costCS with EBeam=%f\n",ebeam);
282 Double_t mp=fMProton;
283 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
284 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
285 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
286 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
287 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
288 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
289 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
290 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
291 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
292 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
293 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
295 // Fill the Lorentz vector for projectile and target
296 // For the moment we do not consider the crossing angle
297 // Projectile runs towards the MUON arm
298 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
299 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
301 // --- Get the muons parameters in the CM frame
303 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
304 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
306 // --- Obtain the dimuon parameters in the CM frame
308 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
310 // --- Translate the dimuon parameters in the dimuon rest frame
312 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
313 TLorentzVector pMu1Dimu=pMu1CM;
314 TLorentzVector pMu2Dimu=pMu2CM;
315 TLorentzVector pProjDimu=pProjCM;
316 TLorentzVector pTargDimu=pTargCM;
317 pMu1Dimu.Boost(beta);
318 pMu2Dimu.Boost(beta);
319 pProjDimu.Boost(beta);
320 pTargDimu.Boost(beta);
322 // --- Determine the z axis for the CS angle
324 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
326 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
330 cost = zaxisCS.Dot((pMu1Dimu.Vect()).Unit());
331 // Theta CS is not properly defined for Like-Sign muons
332 if(mu2Charge > 0 && cost<0) cost=-cost;
334 // Theta CS is not properly defined for Like-Sign muons
335 cost = zaxisCS.Dot((pMu2Dimu.Vect()).Unit());
336 if(mu2Charge < 0 && cost<0) cost=-cost;
341 //______________________________________________________________________________
342 Double_t AliAODDimuon::CostHe(){
343 // Calculation the Helicity polarization angle (adapted from code by R. Arnaldi)
344 if(CheckPointers())return -999999999;
345 Double_t ebeam=3500; //temporary
347 printf("Can not compute costCS with EBeam=%f\n",ebeam);
350 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
351 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
352 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
353 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
354 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
355 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
356 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
357 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
358 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
359 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
361 // --- Get the muons parameters in the CM frame
363 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
364 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
366 // --- Obtain the dimuon parameters in the CM frame
368 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
370 // --- Translate the dimuon parameters in the dimuon rest frame
372 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
373 TLorentzVector pMu1Dimu=pMu1CM;
374 TLorentzVector pMu2Dimu=pMu2CM;
375 pMu1Dimu.Boost(beta);
376 pMu2Dimu.Boost(beta);
378 // --- Determine the z axis for the calculation of the polarization angle
379 // (i.e. the direction of the dimuon in the CM system)
382 zaxis=(pDimuCM.Vect()).Unit();
384 // --- Calculation of the polarization angle (Helicity)
385 // (angle between mu+ and the z axis defined above)
389 cost = zaxis.Dot((pMu1Dimu.Vect()).Unit());
390 // Theta Helicity is not properly defined for Like-Sign muons
391 if(mu2Charge > 0 && cost<0) cost=-cost;
393 cost = zaxis.Dot((pMu2Dimu.Vect()).Unit());
394 // Theta Helicity is not properly defined for Like-Sign muons
395 if(mu2Charge < 0 && cost<0) cost=-cost;
400 //________________________________________________________________________
401 Double_t AliAODDimuon::PhiCS(){
402 // Cosinus of the Collins-Soper polar decay angle
403 if(CheckPointers())return -999999999;
404 Double_t ebeam=3500.; //temporary
406 printf("Can not compute phiCS with EBeam=%f\n",ebeam);
409 Double_t mp=fMProton;
410 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
411 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
412 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
413 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
414 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
415 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
416 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
417 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
418 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
419 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
421 // Fill the Lorentz vector for projectile and target
422 // For the moment we do not consider the crossing angle
423 // Projectile runs towards the MUON arm
424 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
425 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
427 // --- Get the muons parameters in the CM frame
429 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
430 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
432 // --- Obtain the dimuon parameters in the CM frame
434 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
436 // --- Translate the dimuon parameters in the dimuon rest frame
438 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
439 TLorentzVector pMu1Dimu=pMu1CM;
440 TLorentzVector pMu2Dimu=pMu2CM;
441 TLorentzVector pProjDimu=pProjCM;
442 TLorentzVector pTargDimu=pTargCM;
443 pMu1Dimu.Boost(beta);
444 pMu2Dimu.Boost(beta);
445 pProjDimu.Boost(beta);
446 pTargDimu.Boost(beta);
448 // --- Determine the z axis for the CS angle
450 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
452 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
454 TVector3 yaxisCS=(((pProjDimu.Vect()).Unit()).Cross((pTargDimu.Vect()).Unit())).Unit();
455 TVector3 xaxisCS=(yaxisCS.Cross(zaxisCS)).Unit();
458 if(mu1Charge>0) phi = TMath::ATan2((pMu1Dimu.Vect()).Dot(yaxisCS),((pMu1Dimu.Vect()).Dot(xaxisCS)));
459 else phi = TMath::ATan2((pMu2Dimu.Vect()).Dot(yaxisCS),((pMu2Dimu.Vect()).Dot(xaxisCS)));
464 //______________________________________________________________________________
465 Double_t AliAODDimuon::PhiHe(){
466 // Calculation the Helicity aimuthal angle (adapted from code by R. Arnaldi)
467 if(CheckPointers())return -999999999;
468 Double_t ebeam=3500; //temporary
470 printf("Can not compute phiHE with EBeam=%f\n",ebeam);
473 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-fMProton*fMProton);
474 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
475 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
476 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
477 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
478 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
479 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
480 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
481 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
482 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
484 // Fill the Lorentz vector for projectile and target
485 // For the moment we consider no crossing angle
486 // Projectile runs towards the MUON arm
487 TLorentzVector pProjCM(0.,0.,-pbeam,ebeam); // projectile
488 TLorentzVector pTargCM(0.,0., pbeam,ebeam); // target
490 // --- Get the muons parameters in the CM frame
492 TLorentzVector pMu1CM(pla10,pla11,pla12,e1);
493 TLorentzVector pMu2CM(pla20,pla21,pla22,e2);
495 // --- Obtain the dimuon parameters in the CM frame
497 TLorentzVector pDimuCM=pMu1CM+pMu2CM;
499 // --- Translate the muon parameters in the dimuon rest frame
501 TVector3 zaxis=(pDimuCM.Vect()).Unit();
503 // --- Translate the dimuon parameters in the dimuon rest frame
505 TVector3 beta=(-1./pDimuCM.E())*pDimuCM.Vect();
506 TLorentzVector pMu1Dimu=pMu1CM;
507 TLorentzVector pMu2Dimu=pMu2CM;
508 pMu1Dimu.Boost(beta);
509 pMu2Dimu.Boost(beta);
511 TLorentzVector pProjDimu=pProjCM;
512 TLorentzVector pTargDimu=pTargCM;
513 pProjDimu.Boost(beta);
514 pTargDimu.Boost(beta);
516 TVector3 yaxis=((pProjDimu.Vect()).Cross(pTargDimu.Vect())).Unit();
517 TVector3 xaxis=(yaxis.Cross(zaxis)).Unit();
519 // --- Calculation of the azimuthal angle (Helicity)
522 if(mu1Charge>0) phi = TMath::ATan2((pMu1Dimu.Vect()).Dot(yaxis),(pMu1Dimu.Vect()).Dot(xaxis));
523 else phi = TMath::ATan2((pMu2Dimu.Vect()).Dot(yaxis),(pMu2Dimu.Vect()).Dot(xaxis));
528 //______________________________________________________________________________
529 Int_t AliAODDimuon::AnyPt(){
530 // Test if the two muons match two trigger tracks
531 if(this->CheckPointers())return 0;
532 return (((AliAODTrack*)fMu[0].GetObject())->MatchTrigger())&&
533 (((AliAODTrack*)fMu[0].GetObject())->MatchTrigger());
536 //______________________________________________________________________________
537 Int_t AliAODDimuon::LowPt(){
538 // Test if the two muons match two trigger tracks with a "Low Pt" cut
539 if(this->CheckPointers())return 0;
540 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt())&&
541 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt());
544 //______________________________________________________________________________
545 Int_t AliAODDimuon::HighPt(){
546 // Test if the two muons match two trigger tracks with a "High Pt" cut
547 if(this->CheckPointers())return 0;
548 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt())&&
549 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt());
552 //______________________________________________________________________________
553 Double_t AliAODDimuon::MaxChi2Match(){
554 // Maximum matching Chi2 between track and trigger track
555 if(this->CheckPointers())return -999999999;
556 return TMath::Max((((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()),
557 (((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()));