]> git.uio.no Git - u/mrichter/AliRoot.git/blame - RICH/AliRICHResponseV0.cxx
Repositioned the RICH modules.
[u/mrichter/AliRoot.git] / RICH / AliRICHResponseV0.cxx
CommitLineData
237c933d 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 $Log$
d65c4a20 18 Revision 1.5 2001/02/23 17:25:08 jbarbosa
19 Re-definition of IntPH() to accomodate for wire sag effect.
20
733b4fa4 21 Revision 1.4 2000/12/01 17:37:44 morsch
22 Replace in argument of SetTrack(..) string constant by kPFeedBackPhoton.
23
9e4670f0 24 Revision 1.3 2000/10/03 21:44:09 morsch
25 Use AliSegmentation and AliHit abstract base classes.
26
a2f7eaf6 27 Revision 1.2 2000/10/02 21:28:12 fca
28 Removal of useless dependecies via forward declarations
29
94de3818 30 Revision 1.1 2000/06/12 15:29:37 jbarbosa
31 Cleaned up version.
32
237c933d 33*/
34
35#include "AliRICHResponseV0.h"
a2f7eaf6 36#include "AliSegmentation.h"
237c933d 37#include "AliRun.h"
94de3818 38#include "AliMC.h"
237c933d 39
40#include <TMath.h>
41#include <TRandom.h>
42#include <TParticle.h>
43//___________________________________________
44ClassImp(AliRICHResponseV0)
45
733b4fa4 46Float_t AliRICHResponseV0::IntPH(Float_t eloss, Float_t yhit)
237c933d 47{
48 // Get number of electrons and return charge
49
50 Int_t nel;
51 nel= Int_t(eloss/fEIonisation);
52
53 Float_t charge=0;
733b4fa4 54 Double_t gain_var=1;
55
237c933d 56 if (nel == 0) nel=1;
733b4fa4 57
58 if (fWireSag)
59 {
60 //printf("Voltage:%d, Yhit:%f\n",fVoltage, yhit);
61
62 if (fVoltage==2150)
63 {
64 gain_var = 9e-6*TMath::Power(yhit,4) + 2e-7*TMath::Power(yhit,3) - 0.0316*TMath::Power(yhit,2) - 3e-4*yhit + 25.367;
65 //gain_var = 9e-5*TMath::Power(yhit,4) + 2e-6*TMath::Power(yhit,3) - 0.316*TMath::Power(yhit,2) - 3e-3*yhit + 253.67;
66 }
67 if (fVoltage==2100)
68 gain_var = 8e-6*TMath::Power(yhit,4) + 2e-7*TMath::Power(yhit,3) - 0.0283*TMath::Power(yhit,2) - 2e-4*yhit + 23.015;
69 if (fVoltage==2050)
70 gain_var = 7e-6*TMath::Power(yhit,4) + 1e-7*TMath::Power(yhit,3) - 0.0254*TMath::Power(yhit,2) - 2e-4*yhit + 20.888;
71 if (fVoltage==2000)
72 gain_var = 6e-6*TMath::Power(yhit,4) + 8e-8*TMath::Power(yhit,3) - 0.0227*TMath::Power(yhit,2) - 1e-4*yhit + 18.961;
73
74 gain_var = gain_var/100;
75 //printf("Yhit:%f, Gain variation:%f\n",yhit,gain_var);
76
77 Float_t gain = (fChargeSlope + fChargeSlope*gain_var)*.9;
d65c4a20 78 //printf(" Yhit:%f, Gain variation:%f\n",yhit, gain);
733b4fa4 79
80 for (Int_t i=1;i<=nel;i++) {
81 charge -= gain*TMath::Log(gRandom->Rndm());
82 }
83 }
84 else
85 {
86 for (Int_t i=1;i<=nel;i++) {
87 charge -= fChargeSlope*TMath::Log(gRandom->Rndm());
88 }
89 }
90
237c933d 91 return charge;
92}
93
733b4fa4 94Float_t AliRICHResponseV0::IntPH(Float_t yhit)
237c933d 95{
96
97// Get number of electrons and return charge, for a single photon
98
733b4fa4 99 Float_t charge=0;
100 Double_t gain_var=1;
101
102 if (fWireSag)
103 {
104 if (fVoltage==2150)
105 {
106 gain_var = 9e-6*TMath::Power(yhit,4) + 2e-7*TMath::Power(yhit,3) - 0.0316*TMath::Power(yhit,2) - 3e-4*yhit + 25.367;
107 //gain_var = 9e-5*TMath::Power(yhit,4) + 2e-6*TMath::Power(yhit,3) - 0.316*TMath::Power(yhit,2) - 3e-3*yhit + 253.67;
108 }
109 if (fVoltage==2100)
110 gain_var = 8e-6*TMath::Power(yhit,4) + 2e-7*TMath::Power(yhit,3) - 0.0283*TMath::Power(yhit,2) - 2e-4*yhit + 23.015;
111 if (fVoltage==2050)
112 gain_var = 7e-6*TMath::Power(yhit,4) + 1e-7*TMath::Power(yhit,3) - 0.0254*TMath::Power(yhit,2) - 2e-4*yhit + 20.888;
113 if (fVoltage==2000)
114 gain_var = 6e-6*TMath::Power(yhit,4) + 8e-8*TMath::Power(yhit,3) - 0.0227*TMath::Power(yhit,2) - 1e-4*yhit + 18.961;
115
116 gain_var = gain_var/100;
117 //printf(" Yhit:%f, Gain variation:%f\n",yhit, gain_var);
118
119 Float_t gain = (fChargeSlope + fChargeSlope*gain_var)*.9;
120
121 charge -= gain*TMath::Log(gRandom->Rndm());
d65c4a20 122 //printf(" Yhit:%f, Gain variation:%f\n",yhit, gain);
733b4fa4 123 }
124 else
125 {
126 charge -= fChargeSlope*TMath::Log(gRandom->Rndm());
127 }
237c933d 128 return charge;
129}
130
131
132
133// -------------------------------------------
a2f7eaf6 134Float_t AliRICHResponseV0::IntXY(AliSegmentation * segmentation)
237c933d 135{
136
137 const Float_t kInversePitch = 1/fPitch;
138 Float_t response;
139//
140// Integration limits defined by segmentation model
141//
142
143 Float_t xi1, xi2, yi1, yi2;
144 segmentation->IntegrationLimits(xi1,xi2,yi1,yi2);
145
146 xi1=xi1*kInversePitch;
147 xi2=xi2*kInversePitch;
148 yi1=yi1*kInversePitch;
149 yi2=yi2*kInversePitch;
150
151 //printf("Integration Limits: %f-%f, %f-%f\n",xi1,xi2,yi1,yi2);
152
153 //printf("KInversePitch:%f\n",kInversePitch);
154
155 //
156// The Mathieson function
157 Double_t ux1=fSqrtKx3*TMath::TanH(fKx2*xi1);
158 Double_t ux2=fSqrtKx3*TMath::TanH(fKx2*xi2);
159
160 Double_t uy1=fSqrtKy3*TMath::TanH(fKy2*yi1);
161 Double_t uy2=fSqrtKy3*TMath::TanH(fKy2*yi2);
162
163 //printf("Integration Data: %f-%f, %f-%f\n",ux1,ux2,uy1,uy2);
164
165 //printf("%f %f %f %f\n",fSqrtKx3,fKx2,fKy4,fKx4);
166
167 response=4.*fKx4*(TMath::ATan(ux2)-TMath::ATan(ux1))*fKy4*(TMath::ATan(uy2)-TMath::ATan(uy1));
168
169 //printf("Response:%f\n",response);
170
171 return response;
172
173}
174
175Int_t AliRICHResponseV0::FeedBackPhotons(Float_t *source, Float_t qtot)
176{
177 //
178 // Generate FeedBack photons
179 //
180 Int_t j, ipart, nt;
181
182 Int_t sNfeed=0;
183
184
185 // Local variables
186 Float_t cthf, ranf[2], phif, enfp = 0, sthf;
187 Int_t i, ifeed;
188 Float_t e1[3], e2[3], e3[3];
189 Float_t vmod, uswop;
190 Float_t fp, random;
191 Float_t dir[3], phi;
192 Int_t nfp;
193 Float_t pol[3], mom[3];
194 TLorentzVector position;
195 //
196 // Determine number of feedback photons
197
198 // Get weight of current particle
199 TParticle *current = (TParticle*)
200 (*gAlice->Particles())[gAlice->CurrentTrack()];
201
202 ifeed = Int_t(current->GetWeight()/100+0.5);
203 ipart = gMC->TrackPid();
204 fp = fAlphaFeedback * qtot;
205 nfp = gRandom->Poisson(fp);
206
207 // This call to fill the time of flight
208 gMC->TrackPosition(position);
209 //
210 // Generate photons
211 for (i = 0; i <nfp; i++) {
212
213 // Direction
214 gMC->Rndm(ranf, 2);
215 cthf = ranf[0] * 2 - 1.;
216 if (cthf < 0) continue;
217 sthf = TMath::Sqrt((1 - cthf) * (1 + cthf));
218 phif = ranf[1] * 2 * TMath::Pi();
219 //
220 gMC->Rndm(&random, 1);
221 if (random <= .57) {
222 enfp = 7.5e-9;
223 } else if (random <= .7) {
224 enfp = 6.4e-9;
225 } else {
226 enfp = 7.9e-9;
227 }
228
229 dir[0] = sthf * TMath::Sin(phif);
230 dir[1] = cthf;
231 dir[2] = sthf * TMath::Cos(phif);
232 gMC->Gdtom(dir, mom, 2);
233 mom[0]*=enfp;
234 mom[1]*=enfp;
235 mom[2]*=enfp;
236
237 // Polarisation
238 e1[0] = 0;
239 e1[1] = -dir[2];
240 e1[2] = dir[1];
241
242 e2[0] = -dir[1];
243 e2[1] = dir[0];
244 e2[2] = 0;
245
246 e3[0] = dir[1];
247 e3[1] = 0;
248 e3[2] = -dir[0];
249
250 vmod=0;
251 for(j=0;j<3;j++) vmod+=e1[j]*e1[j];
252 if (!vmod) for(j=0;j<3;j++) {
253 uswop=e1[j];
254 e1[j]=e3[j];
255 e3[j]=uswop;
256 }
257 vmod=0;
258 for(j=0;j<3;j++) vmod+=e2[j]*e2[j];
259 if (!vmod) for(j=0;j<3;j++) {
260 uswop=e2[j];
261 e2[j]=e3[j];
262 e3[j]=uswop;
263 }
264
265 vmod=0;
266 for(j=0;j<3;j++) vmod+=e1[j]*e1[j];
267 vmod=TMath::Sqrt(1/vmod);
268 for(j=0;j<3;j++) e1[j]*=vmod;
269
270 vmod=0;
271 for(j=0;j<3;j++) vmod+=e2[j]*e2[j];
272 vmod=TMath::Sqrt(1/vmod);
273 for(j=0;j<3;j++) e2[j]*=vmod;
274
275 gMC->Rndm(ranf, 1);
276 phi = ranf[0] * 2 * TMath::Pi();
277 for(j=0;j<3;j++) pol[j]=e1[j]*TMath::Sin(phi)+e2[j]*TMath::Cos(phi);
278 gMC->Gdtom(pol, pol, 2);
279
280 // Put photon on the stack and label it as feedback (51, 52)
281 ++sNfeed;
282
283 gAlice->SetTrack(Int_t(1), gAlice->CurrentTrack(), Int_t(50000051),
284 mom,source,pol,position[3],
9e4670f0 285 kPFeedBackPhoton, nt, 1.);
237c933d 286 }
287 return(sNfeed);
288}
289
290
291
292