Support files added
[u/mrichter/AliRoot.git] / MFT / AliMuonForwardTrackPair.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 //====================================================================================================================================================
17 //
18 //      Description of an ALICE muon forward track pair, i.e. a pair of AliMuonForwardTrack objects
19 //
20 //      Contact author: antonio.uras@cern.ch
21 //
22 //====================================================================================================================================================
23
24 #include "AliLog.h"
25 #include "AliMUONTrackParam.h"
26 #include "TParticle.h"
27 #include "AliMuonForwardTrack.h"
28 #include "AliMUONTrackExtrap.h"
29 #include "TClonesArray.h"
30 #include "TDatabasePDG.h"
31 #include "TLorentzVector.h"
32 #include "TRandom.h"
33 #include "AliMuonForwardTrackPair.h"
34
35 ClassImp(AliMuonForwardTrackPair)
36
37 //====================================================================================================================================================
38
39 AliMuonForwardTrackPair::AliMuonForwardTrackPair():
40   TObject(),
41   fMuonForwardTracks(0),
42   fKinemMC(0,0,0,0)
43 {
44
45   // default constructor
46
47   fMuonForwardTracks = new TClonesArray("AliMuonForwardTrack", 2);
48
49 }
50
51 //====================================================================================================================================================
52
53 AliMuonForwardTrackPair::AliMuonForwardTrackPair(AliMuonForwardTrack *track0, AliMuonForwardTrack *track1):
54   TObject(),
55   fMuonForwardTracks(0),
56   fKinemMC(0,0,0,0)
57 {
58
59   fMuonForwardTracks = new TClonesArray("AliMuonForwardTrack", 2);
60
61   new ((*fMuonForwardTracks)[0]) AliMuonForwardTrack(*track0);
62   new ((*fMuonForwardTracks)[1]) AliMuonForwardTrack(*track1);
63
64   SetKinemMC();
65
66 }
67
68 //====================================================================================================================================================
69
70 AliMuonForwardTrackPair::AliMuonForwardTrackPair(const AliMuonForwardTrackPair& trackPair): 
71   TObject(trackPair),
72   fMuonForwardTracks(trackPair.fMuonForwardTracks),
73   fKinemMC(trackPair.fKinemMC)
74 {
75
76   // copy constructor
77   
78 }
79
80 //====================================================================================================================================================
81
82 AliMuonForwardTrackPair& AliMuonForwardTrackPair::operator=(const AliMuonForwardTrackPair& trackPair) {
83
84   // Asignment operator
85
86   // check assignement to self
87   if (this == &trackPair) return *this;
88
89   // base class assignement
90   AliMuonForwardTrackPair::operator=(trackPair);
91   
92   // clear memory
93   Clear();
94   
95   fMuonForwardTracks = trackPair.fMuonForwardTracks;
96
97   return *this;
98
99 }
100
101 //====================================================================================================================================================
102
103 void AliMuonForwardTrackPair::SetTrack(Int_t iTrack, AliMuonForwardTrack *track) {
104
105   if (iTrack==0 || iTrack==1) {
106     new ((*fMuonForwardTracks)[iTrack]) AliMuonForwardTrack(*track);
107   }
108
109 }
110
111 //====================================================================================================================================================
112
113 Double_t AliMuonForwardTrackPair::GetWeightedOffset(Double_t x, Double_t y, Double_t z) {
114
115   Double_t weightedOffset[2]={0};
116
117   weightedOffset[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetWeightedOffset(x, y, z);
118   weightedOffset[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetWeightedOffset(x, y, z);
119
120   Double_t weightedOffsetDimuon = TMath::Sqrt(0.5 * (weightedOffset[0]*weightedOffset[0] + weightedOffset[1]*weightedOffset[1]));
121
122   return weightedOffsetDimuon;
123
124 }
125
126 //====================================================================================================================================================
127
128 Double_t AliMuonForwardTrackPair::GetMass(Double_t z, Int_t nClusters) {
129
130   Int_t idCluster[2] = {0};
131   if (nClusters>0) {
132     idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMFTClusters() - nClusters;
133     idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMFTClusters() - nClusters;
134   }
135   if (idCluster[0]<0) idCluster[0] = 0;
136   if (idCluster[1]<0) idCluster[1] = 0;
137
138   Double_t momentum[2] = {0};
139   
140   AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMFTCluster(idCluster[0]);
141   AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMFTCluster(idCluster[1]);
142
143   AliDebug(2, Form("MFT before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
144                    param0->Px(), param0->Py(), param0->Pz(), 
145                    param1->Px(), param1->Py(), param1->Pz()));
146
147   if (TMath::Abs(z)<1e6) {
148     AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
149     AliMUONTrackExtrap::ExtrapToZCov(param0, z);
150     AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
151     AliMUONTrackExtrap::ExtrapToZCov(param1, z);
152   }
153
154   AliDebug(2, Form("MFT after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
155                    param0->Px(), param0->Py(), param0->Pz(), 
156                    param1->Px(), param1->Py(), param1->Pz()));
157
158   momentum[0] = (param0->P());
159   momentum[1] = (param1->P());
160
161   Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
162
163   TLorentzVector dimu;
164
165   dimu.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
166
167   dimu.SetPx(param0->Px() + param1->Px());
168   dimu.SetPy(param0->Py() + param1->Py());
169   dimu.SetPz(param0->Pz() + param1->Pz());
170
171   return dimu.M();
172
173 }
174
175 //====================================================================================================================================================
176
177 Double_t AliMuonForwardTrackPair::GetMassWithoutMFT(Double_t x, Double_t y, Double_t z, Int_t nClusters) {
178
179   Int_t idCluster[2] = {0};
180   if (nClusters>0) {
181     idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMUONClusters() - nClusters;
182     idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMUONClusters() - nClusters;
183   }
184   if (idCluster[0]<0) idCluster[0] = 0;
185   if (idCluster[1]<0) idCluster[1] = 0;
186
187   AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMUONCluster(idCluster[0]);
188   AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMUONCluster(idCluster[1]);
189
190   AliDebug(2, Form("MUON before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
191                    param0->Px(), param0->Py(), param0->Pz(), 
192                    param1->Px(), param1->Py(), param1->Pz()));
193
194   AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
195   AliMUONTrackExtrap::ExtrapToVertex(param0, x, y, z, 0., 0.);   // this should reproduce what is done in AliMUONESDInterface::MUONToESD(...) 
196   AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
197   AliMUONTrackExtrap::ExtrapToVertex(param1, x, y, z, 0., 0.);   // this should reproduce what is done in AliMUONESDInterface::MUONToESD(...) 
198
199   AliDebug(2, Form("MUON after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
200                    param0->Px(), param0->Py(), param0->Pz(), 
201                    param1->Px(), param1->Py(), param1->Pz()));
202
203   Double_t momentum[2] = {0}; 
204
205   momentum[0] = param0->P();
206   momentum[1] = param1->P();
207
208   Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
209
210   TLorentzVector dimu;
211
212   dimu.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
213
214   dimu.SetPx(param0->Px() + param1->Px());
215   dimu.SetPy(param0->Py() + param1->Py());
216   dimu.SetPz(param0->Pz() + param1->Pz());
217
218   return dimu.M();
219
220 }
221
222 //====================================================================================================================================================
223
224 void AliMuonForwardTrackPair::SetKinemMC() {
225
226   AliDebug(2, Form("MC: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
227                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Px(),
228                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Py(),
229                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Pz(),
230                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Px(),
231                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Py(),
232                    ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Pz()));
233
234   fKinemMC.SetE(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Energy() +
235                 ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Energy());
236   
237   fKinemMC.SetPx(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Px() +
238                  ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Px());
239   
240   fKinemMC.SetPy(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Py() +
241                  ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Py());
242   
243   fKinemMC.SetPz(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Pz() +
244                  ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()->Pz());
245
246 }
247
248 //====================================================================================================================================================
249
250 Bool_t AliMuonForwardTrackPair::IsResonance() {
251
252   Bool_t result = kFALSE;
253
254   Int_t labelMC[2] = {0};
255   Int_t codePDG[2] = {0};
256   
257   for (Int_t iTrack=0; iTrack<2; iTrack++) {
258     labelMC[iTrack] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(iTrack))->GetParentMCLabel(0);
259     codePDG[iTrack] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(iTrack))->GetParentPDGCode(0);
260   }
261
262   AliDebug(1, Form("Muons' mothers: (%d, %d)", labelMC[0], labelMC[1]));
263
264   if (labelMC[0]==labelMC[1] && codePDG[0]==codePDG[1] && (codePDG[0]==   113 ||
265                                                            codePDG[0]==   223 ||
266                                                            codePDG[0]==   333 ||
267                                                            codePDG[0]==   443 ||
268                                                            codePDG[0]==100443 ||
269                                                            codePDG[0]==   553 ||
270                                                            codePDG[0]==100553 ) ) result = kTRUE;
271
272   if (result) AliDebug(1, Form("Pair is a resonance %d", codePDG[0]));
273
274   return result; 
275
276 }
277
278 //====================================================================================================================================================