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
60 //printf("Creating dimuon from %p %p\n",mu0,mu1);
65 //______________________________________________________________________________
66 AliAODDimuon::~AliAODDimuon()
73 //______________________________________________________________________________
74 void AliAODDimuon::BookP(){
75 // Fills the dimuon momentum if not filled yet
76 static UInt_t unID[2]={0,0};
78 fP=new TLorentzVector(Px(),Py(),Pz(),E());
79 unID[0]=fMu[0].GetUniqueID();
80 unID[1]=fMu[1].GetUniqueID();
82 // For efficiency reasons
83 if((unID[0]!=fMu[0].GetUniqueID())||(unID[1]!=fMu[1].GetUniqueID())){
84 fP->SetPxPyPzE(Px(),Py(),Pz(),E());
85 unID[0]=fMu[0].GetUniqueID();
86 unID[1]=fMu[1].GetUniqueID();
90 //______________________________________________________________________________
91 Double_t AliAODDimuon::Px() const {
93 if(this->CheckPointers())return -999999999;
94 return ((AliAODTrack*)fMu[0].GetObject())->Px()+
95 ((AliAODTrack*)fMu[1].GetObject())->Px();
98 //______________________________________________________________________________
99 Double_t AliAODDimuon::Py() const {
101 if(this->CheckPointers())return -999999999;
102 return ((AliAODTrack*)fMu[0].GetObject())->Py()+
103 ((AliAODTrack*)fMu[1].GetObject())->Py();
106 //______________________________________________________________________________
107 Double_t AliAODDimuon::Pz() const {
109 if(this->CheckPointers())return -999999999;
110 return ((AliAODTrack*)fMu[0].GetObject())->Pz()+
111 ((AliAODTrack*)fMu[1].GetObject())->Pz();
114 //______________________________________________________________________________
115 Double_t AliAODDimuon::Pt() const {
117 if(this->CheckPointers())return -999999999;
120 return TMath::Sqrt(px*px+py*py);
124 //______________________________________________________________________________
125 Double_t AliAODDimuon::E() const {
127 if(this->CheckPointers())return -999999999;
128 return ((AliAODTrack*)fMu[0].GetObject())->E()+
129 ((AliAODTrack*)fMu[1].GetObject())->E();
132 //______________________________________________________________________________
133 Double_t AliAODDimuon::P() const {
134 // This is just to override the virtual function
135 printf("You should never call: Double_t AliAODDimuon::P() const\n");
139 //______________________________________________________________________________
140 Double_t AliAODDimuon::P() {
142 if(this->CheckPointers())return -999999999;
147 //______________________________________________________________________________
148 Double_t AliAODDimuon::M() const {
149 // This is just to override the virtual function
150 printf("You should never call: Double_t AliAODDimuon::M() const\n");
154 //______________________________________________________________________________
155 Double_t AliAODDimuon::M() {
156 // Dimuon invariant mass
157 if(this->CheckPointers())return -999999999;
162 //______________________________________________________________________________
163 Double_t AliAODDimuon::Mass() {
164 // Dimuon invariant mass
165 if(this->CheckPointers())return -999999999;
170 //______________________________________________________________________________
171 Double_t AliAODDimuon::Eta() const {
172 // This is just to override the virtual function
173 printf("You should never call: Double_t AliAODDimuon::Eta() const\n");
177 //______________________________________________________________________________
178 Double_t AliAODDimuon::Eta() {
179 // Dimuon pseudorapidity
180 if(this->CheckPointers())return -999999999;
185 //______________________________________________________________________________
186 Double_t AliAODDimuon::Phi() const {
187 // This is just to override the virtual function
188 printf("You should never call: Double_t AliAODDimuon::Phi() const\n");
192 //______________________________________________________________________________
193 Double_t AliAODDimuon::Phi() {
194 // Dimuon asimuthal angle
195 if(this->CheckPointers())return -999999999;
199 //______________________________________________________________________________
200 Double_t AliAODDimuon::Theta() const {
201 // This is just to override the virtual function
202 printf("You should never call: Double_t AliAODDimuon::Theta() const\n");
206 //______________________________________________________________________________
207 Double_t AliAODDimuon::Theta() {
208 // Dimuon polar angle
209 if(this->CheckPointers())return -999999999;
214 //______________________________________________________________________________
215 Double_t AliAODDimuon::Y() const {
216 // This is just to override the virtual function
217 printf("You should never call: Double_t AliAODDimuon::Y() const\n");
221 //______________________________________________________________________________
222 Double_t AliAODDimuon::Y() {
224 if(this->CheckPointers())return -999999999;
226 return fP->Rapidity();
229 //______________________________________________________________________________
230 Short_t AliAODDimuon::Charge() const {
232 if(this->CheckPointers())return -999;
233 return ((AliAODTrack*)fMu[0].GetObject())->Charge()+
234 ((AliAODTrack*)fMu[1].GetObject())->Charge();
237 //______________________________________________________________________________
238 Int_t AliAODDimuon::CheckPointers() const{
239 // Checks if the track pointers have been initialized
240 if(fMu[0]==0||fMu[1]==0){
241 printf("Dimuon not initialized\n");
244 if((fMu[0].GetObject())==0||(fMu[1].GetObject())==0){
245 printf("Can not get objects. Got: %p %p\n",fMu[0].GetObject(),fMu[1].GetObject());
251 //______________________________________________________________________________
252 void AliAODDimuon::SetMu(Int_t imu, AliAODTrack *mu){
253 // Assign a track pointer
259 //______________________________________________________________________________
260 void AliAODDimuon::SetMuons(AliAODTrack *mu0, AliAODTrack *mu1){
261 // Assign the track pointers
266 //______________________________________________________________________________
267 Double_t AliAODDimuon::XF() {
269 //Double_t ebeam=((AliAODEventInfo*)fEi.GetObject())->EBeam();
270 Double_t ebeam = 3500.; // temporary
272 printf("AliAODDimuon::xf: can not compute xf with EBeam=%f\n",ebeam);
275 if(this->CheckPointers())return -999999999;
278 Double_t pMax=TMath::Sqrt(ebeam*ebeam-mDimu*mDimu);
282 //______________________________________________________________________________
283 // Calculation the Collins-Soper angle (adapted from code by R. Arnaldi)
284 Double_t AliAODDimuon::CostCS(){
285 // Cosinus of the Collins-Soper polar decay angle
286 if(CheckPointers())return -999999999;
287 Double_t ebeam=3500.; //temporary
289 printf("Can not compute costCS with EBeam=%f\n",ebeam);
292 Double_t mp=fMProton;
293 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-mp*mp);
294 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
295 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
296 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
297 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
298 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
299 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
300 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
301 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
302 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
303 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
305 // Fill the Lorentz vector for projectile and target
306 // For the moment we do not consider the crossing angle
307 // Projectile runs towards the MUON arm
308 TLorentzVector pProjLab(0.,0.,-pbeam,ebeam); // projectile
309 TLorentzVector pTargLab(0.,0., pbeam,ebeam); // target
311 // --- Get the muons parameters in the LAB frame
313 TLorentzVector pMu1Lab(pla10,pla11,pla12,e1);
314 TLorentzVector pMu2Lab(pla20,pla21,pla22,e2);
316 // --- Obtain the dimuon parameters in the LAB frame
318 TLorentzVector pDimuLab=pMu1Lab+pMu2Lab;
320 // --- Translate the dimuon parameters in the dimuon rest frame
322 TVector3 beta=(-1./pDimuLab.E())*pDimuLab.Vect();
323 TLorentzVector pMu1Dimu=pMu1Lab;
324 TLorentzVector pMu2Dimu=pMu2Lab;
325 TLorentzVector pProjDimu=pProjLab;
326 TLorentzVector pTargDimu=pTargLab;
327 pMu1Dimu.Boost(beta);
328 pMu2Dimu.Boost(beta);
329 pProjDimu.Boost(beta);
330 pTargDimu.Boost(beta);
332 // --- Determine the z axis for the CS angle
334 TVector3 zaxisCS=(((pProjDimu.Vect()).Unit())-((pTargDimu.Vect()).Unit())).Unit();
336 // --- Determine the CS angle (angle between mu+ and the z axis defined above)
340 cost = zaxisCS.Dot((pMu1Dimu.Vect()).Unit());
341 // Theta CS is not properly defined for Like-Sign muons
342 if(mu2Charge > 0 && cost<0) cost=-cost;
344 // Theta CS is not properly defined for Like-Sign muons
345 cost = zaxisCS.Dot((pMu2Dimu.Vect()).Unit());
346 if(mu2Charge < 0 && cost<0) cost=-cost;
351 //______________________________________________________________________________
352 // Calculation the Helicity polarization angle (adapted from code by R. Arnaldi)
353 Double_t AliAODDimuon::CostHe(){
354 // Cosinus of the polar decay angle in the Helicity reference frame
355 if(CheckPointers())return -999999999;
356 Double_t ebeam=3500; //temporary
358 printf("Can not compute costCS with EBeam=%f\n",ebeam);
361 Double_t pbeam=TMath::Sqrt(ebeam*ebeam-fMProton*fMProton);
362 Double_t pla10=((AliAODTrack*)fMu[0].GetObject())->Px();
363 Double_t pla11=((AliAODTrack*)fMu[0].GetObject())->Py();
364 Double_t pla12=((AliAODTrack*)fMu[0].GetObject())->Pz();
365 Double_t e1=((AliAODTrack*)fMu[0].GetObject())->E();
366 Double_t mu1Charge=((AliAODTrack*)fMu[0].GetObject())->Charge();
367 Double_t pla20=((AliAODTrack*)fMu[1].GetObject())->Px();
368 Double_t pla21=((AliAODTrack*)fMu[1].GetObject())->Py();
369 Double_t pla22=((AliAODTrack*)fMu[1].GetObject())->Pz();
370 Double_t e2=((AliAODTrack*)fMu[1].GetObject())->E();
371 Double_t mu2Charge=((AliAODTrack*)fMu[1].GetObject())->Charge();
373 // Fill the Lorentz vector for projectile and target
374 // For the moment we consider no crossing angle
375 // Projectile runs towards the MUON arm
376 TLorentzVector pProjLab(0.,0.,-pbeam,ebeam); // projectile
377 TLorentzVector pTargLab(0.,0., pbeam,ebeam); // target
379 // --- Get the muons parameters in the LAB frame
381 TLorentzVector pMu1Lab(pla10,pla11,pla12,e1);
382 TLorentzVector pMu2Lab(pla20,pla21,pla22,e2);
384 // --- Obtain the dimuon parameters in the LAB frame
386 TLorentzVector pDimuLab=pMu1Lab+pMu2Lab;
388 // --- Translate the dimuon parameters in the dimuon rest frame
390 TVector3 beta=(-1./pDimuLab.E())*pDimuLab.Vect();
391 TLorentzVector pMu1Dimu=pMu1Lab;
392 TLorentzVector pMu2Dimu=pMu2Lab;
393 pMu1Dimu.Boost(beta);
394 pMu2Dimu.Boost(beta);
396 // --- Translate the dimuon parameters in the CM frame
398 TLorentzVector pDimuCM; //CM frame
400 beta2=(-1./(fMProton+pProjLab.E()))*pProjLab.Vect();
402 pDimuCM.Boost(beta2);
404 // --- Determine the z axis for the calculation of the polarization angle
405 // (i.e. the direction of the dimuon in the CM system)
408 zaxis=(pDimuCM.Vect()).Unit();
410 // --- Calculation of the polarization angle (Helicity)
411 // (angle between mu+ and the z axis defined above)
415 cost = zaxis.Dot((pMu1Dimu.Vect()).Unit());
416 // Theta Helicity is not properly defined for Like-Sign muons
417 if(mu2Charge > 0 && cost<0) cost=-cost;
419 cost = zaxis.Dot((pMu2Dimu.Vect()).Unit());
420 // Theta Helicity is not properly defined for Like-Sign muons
421 if(mu2Charge < 0 && cost<0) cost=-cost;
426 //______________________________________________________________________________
427 Int_t AliAODDimuon::AnyPt(){
428 // Test if the two muons match two trigger tracks
429 if(this->CheckPointers())return 0;
430 return (((AliAODTrack*)fMu[0].GetObject())->MatchTrigger())&&
431 (((AliAODTrack*)fMu[0].GetObject())->MatchTrigger());
434 //______________________________________________________________________________
435 Int_t AliAODDimuon::LowPt(){
436 // Test if the two muons match two trigger tracks with a "Low Pt" cut
437 if(this->CheckPointers())return 0;
438 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt())&&
439 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerLowPt());
442 //______________________________________________________________________________
443 Int_t AliAODDimuon::HighPt(){
444 // Test if the two muons match two trigger tracks with a "High Pt" cut
445 if(this->CheckPointers())return 0;
446 return (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt())&&
447 (((AliAODTrack*)fMu[0].GetObject())->MatchTriggerHighPt());
450 //______________________________________________________________________________
451 Double_t AliAODDimuon::MaxChi2Match(){
452 // Maximum matching Chi2 between track and trigger track
453 if(this->CheckPointers())return -999999999;
454 return TMath::Max((((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()),
455 (((AliAODTrack*)fMu[0].GetObject())->GetChi2MatchTrigger()));