Fixes for cmake
[u/mrichter/AliRoot.git] / PWG1 / AliGenKinkInfo.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 /*
19
20 Origin: marian.ivanov@cern.ch
21 Container classes with MC infomation
22
23 */
24
25 #if !defined(__CINT__) || defined(__MAKECINT__)
26 #include <stdio.h>
27 #include <string.h>
28 //ROOT includes
29 #include "TROOT.h"
30 #include "Rtypes.h"
31 #include "TFile.h"
32 #include "TTree.h"
33 #include "TChain.h"
34 #include "TCut.h"
35 #include "TString.h"
36 #include "TStopwatch.h"
37 #include "TParticle.h"
38 #include "TDatabasePDG.h"
39 #include "TSystem.h"
40 #include "TCanvas.h"
41 #include "TPolyLine3D.h"
42
43 //ALIROOT includes
44 #include "AliRun.h"
45 #include "AliStack.h"
46 #include "AliSimDigits.h"
47 #include "AliTPCParam.h"
48 #include "AliTPC.h"
49 #include "AliTPCLoader.h"
50 #include "AliDetector.h"
51 #include "AliTrackReference.h"
52 #include "AliTPCParamSR.h"
53 #include "AliTracker.h"
54 #include "AliMagF.h"
55 #include "AliHelix.h"
56 #include "AliTrackPointArray.h"
57
58 #endif
59 #include "AliGenKinkInfo.h"
60
61 //
62 // 
63
64 ClassImp(AliGenKinkInfo)
65
66
67
68 /////////////////////////////////////////////////////////////////////////////////
69
70 AliGenKinkInfo::AliGenKinkInfo():
71   fMCd(),       //info about daughter particle - second particle for V0
72   fMCm(),       //info about mother particle   - first particle for V0
73   fMCDist1(0),  //info about closest distance according closest MC - linear DCA
74   fMCDist2(0),  //info about closest distance parabolic DCA
75   fMCRr(0),     // rec position of the vertex 
76   fMCR(0)       //exact r position of the vertex
77 {
78   //
79   // default constructor
80   //
81   for (Int_t i=0;i<3;i++){
82     fMCPdr[i]=0;
83     fMCPd[i]=0;
84     fMCX[i]=0;
85     fMCPm[i]=0;
86     fMCAngle[i]=0;
87   }
88   for (Int_t i=0; i<2; i++) {
89     fPdg[i]= 0;
90     fLab[i]=0;
91   }
92 }
93
94 void AliGenKinkInfo::Update()
95 {
96   //
97   // Update information
98   //  some redundancy - faster acces to the values in analysis code
99   //
100   fMCPd[0] = fMCd.GetParticle().Px();
101   fMCPd[1] = fMCd.GetParticle().Py();
102   fMCPd[2] = fMCd.GetParticle().Pz();
103   fMCPd[3] = fMCd.GetParticle().P();
104   //
105   fMCX[0]  = fMCd.GetParticle().Vx();
106   fMCX[1]  = fMCd.GetParticle().Vy();
107   fMCX[2]  = fMCd.GetParticle().Vz();
108   fMCR       = TMath::Sqrt( fMCX[0]*fMCX[0]+fMCX[1]*fMCX[1]);
109   //
110   fPdg[0]    = fMCd.GetParticle().GetPdgCode();
111   fPdg[1]    = fMCm.GetParticle().GetPdgCode();
112   //
113   fLab[0]    = fMCd.GetParticle().GetUniqueID();
114   fLab[1]    = fMCm.GetParticle().GetUniqueID();
115   //
116   //
117   //
118   Double_t x1[3],p1[3];
119   TParticle& pdaughter = fMCd.GetParticle();
120   x1[0] = pdaughter.Vx();      
121   x1[1] = pdaughter.Vy();
122   x1[2] = pdaughter.Vz();
123   p1[0] = pdaughter.Px();
124   p1[1] = pdaughter.Py();
125   p1[2] = pdaughter.Pz();
126   Double_t sign = (pdaughter.GetPDG()->Charge()>0)? -1:1;
127   AliHelix dhelix1(x1,p1,sign);
128   //
129   //
130   Double_t x2[3],p2[3];            
131   //
132   TParticle& pmother = fMCm.GetParticle();
133   x2[0] = pmother.Vx();      
134   x2[1] = pmother.Vy();      
135   x2[2] = pmother.Vz();      
136   p2[0] = pmother.Px();
137   p2[1] = pmother.Py();
138   p2[2] = pmother.Pz();
139   //
140   const AliTrackReference & pdecay = fMCm.GetTRdecay();
141   x2[0] = pdecay.X();      
142   x2[1] = pdecay.Y();      
143   x2[2] = pdecay.Z();      
144   p2[0] = pdecay.Px();
145   p2[1] = pdecay.Py();
146   p2[2] = pdecay.Pz();
147   //
148   sign = (pmother.GetPDG()->Charge()>0) ? -1:1;
149   AliHelix mhelix(x2,p2,sign);
150   
151   //
152   //
153   //
154   //find intersection linear
155   //
156   Double_t distance1, distance2;
157   Double_t phase[2][2],radius[2];
158   Int_t  points = dhelix1.GetRPHIintersections(mhelix, phase, radius);
159   Double_t delta1=10000,delta2=10000;    
160   if (points>0){
161     dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
162     dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
163     dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
164   }
165   if (points==2){    
166     dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
167     dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
168     dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
169   }
170   distance1 = TMath::Min(delta1,delta2);
171   //
172   //find intersection parabolic
173   //
174   points = dhelix1.GetRPHIintersections(mhelix, phase, radius);
175   delta1=10000,delta2=10000;  
176   
177   if (points>0){
178     dhelix1.ParabolicDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
179   }
180   if (points==2){    
181     dhelix1.ParabolicDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
182   }
183   
184   distance2 = TMath::Min(delta1,delta2);
185   //
186   if (delta1<delta2){
187     //get V0 info
188     dhelix1.Evaluate(phase[0][0],fMCXr);
189     dhelix1.GetMomentum(phase[0][0],fMCPdr);
190     mhelix.GetMomentum(phase[0][1],fMCPm);
191     dhelix1.GetAngle(phase[0][0],mhelix,phase[0][1],fMCAngle);
192     fMCRr = TMath::Sqrt(radius[0]);
193   }
194   else{
195     dhelix1.Evaluate(phase[1][0],fMCXr);
196     dhelix1.GetMomentum(phase[1][0],fMCPdr);
197     mhelix.GetMomentum(phase[1][1],fMCPm);
198     dhelix1.GetAngle(phase[1][0],mhelix,phase[1][1],fMCAngle);
199     fMCRr = TMath::Sqrt(radius[1]);
200   }     
201   //            
202   //
203   fMCDist1 = TMath::Sqrt(distance1);
204   fMCDist2 = TMath::Sqrt(distance2);      
205     
206 }
207
208
209 Float_t AliGenKinkInfo::GetQt(){
210   //
211   //
212   Float_t momentumd = TMath::Sqrt(fMCPd[0]*fMCPd[0]+fMCPd[1]*fMCPd[1]+fMCPd[2]*fMCPd[2]);
213   return TMath::Sin(fMCAngle[2])*momentumd;
214 }
215
216