Retrofeed from the Release developement
[u/mrichter/AliRoot.git] / HBTAN / AliHBTPositionRandomizer.cxx
1 #include "AliHBTPositionRandomizer.h"
2 //___________________________________________________
3 ////////////////////////////////////////////////////////////////////////////////
4 // 
5 // class AliHBTPositionRandomizer
6 //
7 // These class randomizes particle vertex positions
8 // Piotr.Skowronski@cern.ch
9 //
10 ////////////////////////////////////////////////////////////////////////////////
11
12 #include <TRandom.h>
13 #include "AliAOD.h"
14 #include "AliVAODParticle.h"
15
16
17 ClassImp(AliHBTPositionRandomizer)
18
19 const Int_t AliHBTPositionRandomizer::fgkNumberOfPids = 10;
20
21 /*********************************************************************/
22
23 AliHBTPositionRandomizer::AliHBTPositionRandomizer():
24  fReader(0x0),
25  fDefaultRandomizer(0x0),
26  fRandomizers(0x0),
27  fNPid(0),
28  fPids(0x0),
29  fAddToExistingPos(kFALSE),
30  fOnlyParticlesFromVertex(kFALSE),
31  fRandomizeTracks(kFALSE),
32  fVX(0.0),
33  fVY(0.0),
34  fVZ(0.0)
35 {
36 //constructor
37 }
38 /*********************************************************************/
39
40 AliHBTPositionRandomizer::AliHBTPositionRandomizer(AliReader* reader):
41  fReader(reader),
42  fRandomizers(new TObjArray(fgkNumberOfPids)),
43  fNPid(1),
44  fPids(new Int_t[fgkNumberOfPids]),
45  fAddToExistingPos(kFALSE),
46  fOnlyParticlesFromVertex(kFALSE),
47  fRandomizeTracks(kFALSE),
48  fVX(0.0),
49  fVY(0.0),
50  fVZ(0.0)
51 {
52 //constructor
53  fRandomizers->AddAt(new AliHBTRndmGaussBall(8.0,0.0,0.0),0);
54
55 /*********************************************************************/
56
57 AliHBTPositionRandomizer::AliHBTPositionRandomizer(const AliHBTPositionRandomizer& in):
58  AliReader(in),
59  fReader(),
60  fRandomizers(0x0),
61  fNPid(0),
62  fPids(0x0),
63  fAddToExistingPos(kFALSE),
64  fOnlyParticlesFromVertex(kFALSE),
65  fRandomizeTracks(kFALSE),
66  fVX(0.0),
67  fVY(0.0),
68  fVZ(0.0)
69 {
70   //cpy constructor
71   in.Copy(*this);
72 }
73 /*********************************************************************/
74 AliHBTPositionRandomizer::~AliHBTPositionRandomizer()
75 {
76   //dtor
77   delete fReader;
78   delete fRandomizers;
79   delete [] fPids;
80 }
81 /*********************************************************************/
82 AliHBTPositionRandomizer& AliHBTPositionRandomizer::operator=(const AliHBTPositionRandomizer& in)
83 {
84   //assigment operator
85   in.Copy(*this);
86   return *this;
87 }
88 /*********************************************************************/
89
90 AliAOD* AliHBTPositionRandomizer::GetEventSim() const
91 {
92  // gets from fReader and randomizes current particle event
93  if (fReader == 0x0) 
94   {
95     Error("GetEventSim","Reader is null");
96     return 0x0;
97   } 
98  AliAOD *e =  fReader->GetEventSim();
99  if (e) 
100    if (e->IsRandomized() == kFALSE) 
101      Randomize(e);
102  return e;
103 }
104 /*********************************************************************/
105
106 AliAOD* AliHBTPositionRandomizer::GetEventRec() const
107 {
108  // gets from fReader and randomizes current track event
109  if (fReader == 0x0) 
110   {
111     Error("GetEventRec","Reader is null");
112     return 0x0;
113   }  
114  AliAOD *e =  fReader->GetEventRec();
115  if (fRandomizeTracks && e) if (e->IsRandomized() == kFALSE) Randomize(e);
116  return e;
117 }
118 /*********************************************************************/
119
120 AliAOD* AliHBTPositionRandomizer::GetEventSim(Int_t n)
121 {
122 //returns event n
123  if (fReader == 0x0) return 0x0;
124  AliAOD *e =  fReader->GetEventSim(n);
125  if (e->IsRandomized() == kFALSE) Randomize(e);
126  return e;
127 }
128
129 /*********************************************************************/
130 void AliHBTPositionRandomizer::Randomize(AliAOD* event) const
131 {
132 // randomizes postions of all particles in the event
133   static const Double_t kfmtocm = 1.e-13;
134   if (AliVAODParticle::GetDebug() > 5) Info("Randomize(AliAOD*)","");
135   if (event == 0x0) return;
136
137   for (Int_t i = 0; i < event->GetNumberOfParticles(); i++)
138    {
139      AliVAODParticle* p = event->GetParticle(i);
140      Double_t x,y,z,t=0.0;
141      AliHBTRndm* r = GetRandomizer(p->GetPdgCode());
142      r->Randomize(x,y,z,t,p);
143      
144      Double_t nx = x*kfmtocm;
145      Double_t ny = y*kfmtocm;
146      Double_t nz = z*kfmtocm;
147      Double_t nt = t*kfmtocm;
148      
149      if (fAddToExistingPos)
150       {
151        nx += p->Vx();
152        ny += p->Vy();
153        nz += p->Vz();
154        nt += p->T();
155       }
156      p->SetProductionVertex(nx,ny,nz,nt); 
157    }
158   event->SetRandomized();
159 }
160 /*********************************************************************/
161
162 AliHBTRndm* AliHBTPositionRandomizer::GetRandomizer(Int_t pdg) const
163 {
164   //returns randomizer for a given pdg 
165   Int_t idx = GetRandomizerIndex(pdg);//in most of cases 
166   if (idx < 0) idx = 0;//if not found return a default one
167   return (AliHBTRndm*)fRandomizers->At(idx);
168 }
169 /*********************************************************************/
170 Int_t AliHBTPositionRandomizer::GetRandomizerIndex(Int_t pdg) const
171 {
172   //returns randomizer index for a given pdg 
173
174   if (pdg == 0) return 0;
175   
176   for (Int_t i=1; i < fNPid; i++)
177    {
178      if (fPids[i] == pdg) 
179       return i;
180    }
181    
182   return -1;
183 }
184 /*********************************************************************/
185
186 void AliHBTPositionRandomizer::SetRandomizer(Int_t pid, AliHBTRndm* rndm)
187 {
188  //sets the randomizer for a given particle type
189   if (rndm == 0x0)
190    {
191      Error("SetRandomizer","Randomizer is null");
192      return;
193    }
194    
195   Int_t idx = GetRandomizerIndex(pid);
196   if (idx >= 0) 
197    {
198      delete fRandomizers->At(idx);
199      fRandomizers->AddAt(rndm,idx);
200    }  
201   
202   if (fNPid == fgkNumberOfPids)
203    {
204      Error("SetRandomizer","There is no more space in the array");
205      return;
206    }
207
208   fPids[fNPid] = pid;
209   fRandomizers->AddAt(rndm,fNPid);
210   fNPid++;
211 }
212 /*********************************************************************/
213
214 void AliHBTPositionRandomizer::SetGaussianBall(Int_t pid, Double_t r, Double_t meantime, Double_t sigmatime)
215 {
216  //Sets Gaussian Ball Model
217   SetGaussianBall(pid,r,r,r,meantime,sigmatime);
218 }
219 /*********************************************************************/
220
221 void AliHBTPositionRandomizer::SetGaussianBall(Int_t pid, Double_t rx, Double_t ry, Double_t rz, Double_t meantime, Double_t sigmatime)
222 {
223  //Sets Gaussian Ball Model
224   AliHBTRndm* rndm = new AliHBTRndmGaussBall(rx,ry,rz,meantime,sigmatime);
225   SetRandomizer(pid,rndm);
226 }
227 /*********************************************************************/
228
229 void AliHBTPositionRandomizer::SetCyllinderSurface(Int_t pid, Double_t r, Double_t l)
230 {
231  //Sets Cylinder Surface Model
232   AliHBTRndm* rndm = new  AliHBTRndmCyllSurf(r,l);
233   SetRandomizer(pid,rndm);
234 }
235 /*********************************************************************/
236
237 void AliHBTPositionRandomizer::SetEventVertex(Double_t x, Double_t y,Double_t z)
238 {
239 //sets event vertex position
240   fVX = x;
241   fVY = y;
242   fVZ = z;
243 }
244
245
246 void AliHBTPositionRandomizer::SetEllipse(Int_t pid, Double_t rx, Double_t ryz)
247 {
248 //sets the ellipse randomization for the given pid
249   AliHBTRndm* rndm = new AliHBTRndmEllipse(rx,ryz);
250   SetRandomizer(pid,rndm);   
251 }
252
253 /*********************************************************************/
254 //_____________________________________________________________________
255 ///////////////////////////////////////////////////////////////////////
256 //                                                                   //
257 //  class AliHBTRndmGaussBall                                        //
258 //                                                                   //
259 ///////////////////////////////////////////////////////////////////////
260
261 AliHBTRndmGaussBall::AliHBTRndmGaussBall():
262  fRx(0.0),
263  fRy(0.0),
264  fRz(0.0),
265  fTmean(0.0),
266  fTsigma(0.0)
267 {
268   //constructor
269 }
270 /*********************************************************************/
271
272 AliHBTRndmGaussBall::AliHBTRndmGaussBall(Float_t r, Double_t meantime, Double_t sigmatime):
273  fRx(r),
274  fRy(r),
275  fRz(r),
276  fTmean(meantime),
277  fTsigma(sigmatime)
278 {
279   //constructor
280 }
281 /*********************************************************************/
282
283 AliHBTRndmGaussBall::AliHBTRndmGaussBall(Float_t rx, Float_t ry, Float_t rz, Double_t meantime, Double_t sigmatime):
284  fRx(rx),
285  fRy(ry),
286  fRz(rz),
287  fTmean(meantime),
288  fTsigma(sigmatime)
289 {
290   //constructor
291 }
292 /*********************************************************************/
293
294
295 AliHBTRndmEllipse::AliHBTRndmEllipse(Float_t rmin, Float_t rmax):
296  fRmin(rmin),
297  fRmax(rmax)
298 {
299      //constructor
300 }
301
302 /*********************************************************************/
303
304 void AliHBTRndmGaussBall::Randomize(Double_t& x,Double_t& y,Double_t&z,Double_t&t, AliVAODParticle*/*particle*/) const
305 {
306 //randomizez gauss for each coordinate separately
307   x = gRandom->Gaus(0.0,fRx);
308   y = gRandom->Gaus(0.0,fRy);
309   z = gRandom->Gaus(0.0,fRz);
310   
311   if (fTsigma == 0.0)
312    {
313      t = fTmean;
314      return;
315    }
316   
317   t = gRandom->Gaus(fTmean,fTsigma);
318     
319 }
320 /*********************************************************************/
321 //_____________________________________________________________________
322 ///////////////////////////////////////////////////////////////////////
323 //                                                                   //
324 //  class AliHBTRndmGaussBall                                        //
325 //                                                                   //
326 ///////////////////////////////////////////////////////////////////////
327
328 void AliHBTRndmCyllSurf::Randomize(Double_t& x,Double_t& y,Double_t&z,Double_t&/*t*/, AliVAODParticle* particle) const
329 {
330 //Randomizes x,y,z
331    Double_t r = fR + gRandom->Gaus(0.0, 1.0);
332    Double_t sf = r/particle->Pt();//scaling factor for position transformation ->
333                              //we move direction of string momentum but legth defined by r
334    x = sf*particle->Px();
335    y = sf*particle->Py();
336    z = gRandom->Uniform(-fL,fL);
337 }
338
339 /*********************************************************************/
340 /*********************************************************************/
341
342 void AliHBTRndmEllipse::Randomize(Double_t& x, Double_t& y, Double_t& z,Double_t&/*t*/, AliVAODParticle*p) const
343 {
344     // p=0; //workaround - fix this damn little thingy
345    double R;
346      double phi=p->Phi();
347      
348      R=fRmin+(fRmax-fRmin)*TMath::Sin(phi);
349      x=R*TMath::Sin(phi);
350      y=R*TMath::Cos(phi);
351      z=z;
352 }