]> git.uio.no Git - u/mrichter/AliRoot.git/blame - RICH/Opticals.h
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Protection against empty hit events
[u/mrichter/AliRoot.git] / RICH / Opticals.h
CommitLineData
d48cca74 1#ifdef __CINT__
2void Opticals()
3{
4 gROOT->Reset();
5#endif
6 int i;
7 const Int_t kNbins=26;
8 Float_t aPckov[kNbins];
9 for(i=0;i<kNbins;i++){ //Photons energy intervals
c60862bf 10 aPckov[i]=(0.1*i+5.5)*1e-9;
d48cca74 11 }
12
13 Float_t aIndexFreon[kNbins];
14 Float_t aIndexQuartz[kNbins];
15 Float_t aIndexOpaqueQuartz[kNbins];
16 Float_t aIndexCH4[kNbins];
17 Float_t aIndexGrid[kNbins];
18
c60862bf 19 Float_t e1= 10.666;Float_t e2= 18.125; Float_t f1= 46.411; Float_t f2= 228.71;//RICH TDR page 35
d48cca74 20 for (i=0;i<kNbins;i++){
21 aIndexFreon[i] = aPckov[i] * .0172 * 1e9 + 1.177;
22 Float_t ene=aPckov[i]*1e9;
c60862bf 23 aIndexQuartz[i] = TMath::Sqrt(1. + f1/(e1*e1 - ene*ene) + f2/(e2*e2 - ene*ene) );
d48cca74 24 aIndexOpaqueQuartz[i] =1;
25 aIndexCH4[i] =1.000444;
26 aIndexGrid[i] =1;
27 }
28
29 Float_t aAbsFreon[kNbins]={179.0987, 179.0987, 179.0987, 179.0987, 179.0987,
30 179.0987, 179.0987, 179.0987, 179.0987, 142.9206,
31 56.64957, 25.58622, 13.95293, 12.03905, 10.42953,
32 8.804196, 7.069031, 4.461292, 2.028366, 1.293013,
33 0.577267, 0.40746, 0.334964, 0.0, 0.0,
34 0};
35
36
37 Float_t aAbsQuartz[kNbins]={105.8, 65.52, 48.58, 42.85, 35.79,
38 31.262, 28.598, 27.527, 25.007, 22.815,
39 21.004, 19.266, 17.525, 15.878, 14.177,
40 11.719, 9.282, 6.62, 4.0925, 2.601,
41 1.149, 0.667, 0.3627, 0.192, 0.1497,
42 0.10857};
43
44 Float_t aAbsOpaqueQuartz[kNbins];
45 Float_t aAbsCH4[kNbins];
46 Float_t aAbsGrid[kNbins];
47 Float_t aAbsCsI[kNbins];
48 for (i=0;i<kNbins;i++){
49 aAbsCH4[i] =AbsoCH4(aPckov[i]*1e9);
50 aAbsOpaqueQuartz[i]=1e-5;
51 aAbsCsI[i] =1e-4;
52 aAbsGrid[i] =1e-4;
53 }
54
303f57b3 55 Float_t aQeCsI1[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
56 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
57 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
58 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
59 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
60 0.215799987};
61 Float_t aQeCsI2[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
62 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
63 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
64 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
65 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
66 0.215799987};
67 Float_t aQeCsI3[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
68 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
69 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
70 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
71 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
72 0.215799987};
73 Float_t aQeCsI4[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
74 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
75 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
76 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
77 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
78 0.215799987};
79 Float_t aQeCsI5[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
80 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
81 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
82 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
83 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
84 0.215799987};
85 Float_t aQeCsI6[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
86 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
87 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
88 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
89 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
90 0.215799987};
91 Float_t aQeCsI7[kNbins] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983,
d48cca74 92 0.010125, 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994,
93 0.121500008, 0.141749993, 0.157949999, 0.162, 0.166050002,
94 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992,
95 0.18592, 0.187579989, 0.189239994, 0.190899998, 0.207499996,
96 0.215799987};
97 Float_t aQeAll[kNbins];
98 for(i=0;i<kNbins;i++){
303f57b3 99 aQeCsI1[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
100 aQeCsI2[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
101 aQeCsI3[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
102 aQeCsI4[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
103 aQeCsI5[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
104 aQeCsI6[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
105 aQeCsI7[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
d48cca74 106 aQeAll[i]=1; //QE for all other materials except for PC must be 1.
107 }
108
109
110#ifdef __CINT__
111
112
113//Freon, Quartz, Opaque Quartz,Methane,CsI,Grid
114 const Int_t kFreonMarker= 24; const Int_t kFreonColor=kRed;
115 const Int_t kCH4Marker= 25; const Int_t kCH4Color=kGreen;
116 const Int_t kSiO2Marker= 26; const Int_t kSiO2Color=kBlue;
117 const Int_t kCsIMarker = 2; const Int_t kCsIColor =kMagenta;
118
119 TCanvas *pC=new TCanvas("c1","RICH optics to check",800,900);
120
121 pC->Divide(2,2);
122
123 pC->cd(1);
124 TGraph *pAbsFreonGr=new TGraph(kNbins,aPckov,aAbsFreon);
125 pAbsFreonGr->SetMarkerStyle(kFreonMarker); pAbsFreonGr->SetMarkerColor(kFreonColor);
126 TGraph *pAbsSiO2Gr=new TGraph(kNbins,aPckov,aAbsQuartz);
127 pAbsSiO2Gr->SetMarkerStyle(kSiO2Marker); pAbsSiO2Gr->SetMarkerColor(kSiO2Color);
128 TMultiGraph *pAbsMG=new TMultiGraph();
129 TLegend *pAbsLegend=new TLegend(0.6,0.3,0.85,0.5);
130 pAbsMG->Add(pAbsFreonGr); pAbsLegend->AddEntry(pAbsFreonGr, "freon","p"); //1
131 pAbsMG->Add(pAbsSiO2Gr); pAbsLegend->AddEntry(pAbsSiO2Gr, "quartz","p"); //2
132 pAbsMG->Draw("APL");
133 pAbsMG->GetXaxis()->SetTitle("energy, GeV");
134 pAbsMG->GetYaxis()->SetTitle("absorption length, cm");
135 pAbsMG->Draw("APL");
136 pAbsLegend->Draw();
137
138
139 pC->cd(2);
140 TGraph *pIndexFreonGr=new TGraph(kNbins,aPckov,aIndexFreon);
141 pIndexFreonGr->SetMarkerStyle(kFreonMarker); pIndexFreonGr->SetMarkerColor(kFreonColor);
142 TGraph *pIndexSiO2Gr=new TGraph(kNbins,aPckov,aIndexQuartz);
143 pIndexSiO2Gr->SetMarkerStyle(kSiO2Marker); pIndexSiO2Gr->SetMarkerColor(kSiO2Color);
144 TGraph *pIndexCH4Gr=new TGraph(kNbins,aPckov,aIndexCH4);
145 pIndexCH4Gr->SetMarkerStyle(kCH4Marker); pIndexCH4Gr->SetMarkerColor(kCH4Color);
146 TMultiGraph *pIndexMG=new TMultiGraph();
147 TLegend *pIndexLegend=new TLegend(0.6,0.2,0.85,0.4);
148 pIndexMG->Add(pIndexFreonGr); pIndexLegend->AddEntry(pIndexFreonGr, "freon","p"); //1
149 pIndexMG->Add(pIndexSiO2Gr); pIndexLegend->AddEntry(pIndexSiO2Gr, "quartz","p");
150 pIndexMG->Add(pIndexCH4Gr); pIndexLegend->AddEntry(pIndexCH4Gr, "CH4","p");
151 pIndexMG->Draw("APL");
152 pIndexMG->GetXaxis()->SetTitle("energy, GeV");
153 pIndexMG->GetYaxis()->SetTitle("refraction index");
154 pIndexMG->Draw("APL");
155 pIndexLegend->Draw();
156
157 pC->cd(3);
158 TGraph *pAbsCH4Gr=new TGraph(kNbins,aPckov,aAbsCH4);
159 pAbsCH4Gr->SetMarkerStyle(kCH4Marker); pAbsCH4Gr->SetMarkerColor(kCH4Color);
160 pAbsCH4Gr->Draw("APL");
161 pAbsCH4Gr->GetXaxis()->SetTitle("energy, GeV");
162 pAbsCH4Gr->SetTitle("CH4 absorption length, cm");
163 pAbsCH4Gr->Draw("APL");
164
303f57b3 165 TCanvas *pQeC=new TCanvas("pQeC","CsI QE currently all the same",800,900);
166 pQeC->Divide(2,4);
167 for(int i=1;i<=7;i++){
168 pQeC->cd(i);
169 switch(i){
170 case 1: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI1);pQeCsIGr->SetTitle("Module 1");break;
171 case 2: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI2);pQeCsIGr->SetTitle("Module 2");break;
172 case 3: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI3);pQeCsIGr->SetTitle("Module 3");break;
173 case 4: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI4);pQeCsIGr->SetTitle("Module 4");break;
174 case 5: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI5);pQeCsIGr->SetTitle("Module 5");break;
175 case 6: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI6);pQeCsIGr->SetTitle("Module 6");break;
176 case 7: TGraph *pQeCsIGr=new TGraph(kNbins,aPckov,aQeCsI7);pQeCsIGr->SetTitle("Module 7");break;
177 }
178 pQeCsIGr->SetMarkerStyle(kCsIMarker); pQeCsIGr->SetMarkerColor(kCsIColor);
179 pQeCsIGr->Draw("APL");
180 pQeCsIGr->GetXaxis()->SetTitle("energy, GeV");
181 pQeCsIGr->Draw("APL");
182 }
d48cca74 183}//main
184
185Float_t AbsoCH4(Float_t x)
186{
187
188 //KLOSCH,SCH4(9),WL(9),EM(9),ALENGTH(31)
189 Float_t sch4[9] = {.12,.16,.23,.38,.86,2.8,7.9,28.,80.}; //MB X 10^22
190 //Float_t wl[9] = {153.,152.,151.,150.,149.,148.,147.,146.,145};
191 Float_t em[9] = {8.1,8.158,8.212,8.267,8.322,8.378,8.435,8.493,8.55};
192 const Float_t kLosch=2.686763E19; // LOSCHMIDT NUMBER IN CM-3
303f57b3 193 const Float_t kIgas1=100, kIgas2=0, kOxy=10., kWater=5., kPressure=750.,kTemperature=283.;
d48cca74 194 Float_t pn=kPressure/760.;
195 Float_t tn=kTemperature/273.16;
196
197
198// ------- METHANE CROSS SECTION -----------------
199// ASTROPH. J. 214, L47 (1978)
200
201 Float_t sm=0;
202 if (x<7.75)
203 sm=.06e-22;
204
205 if(x>=7.75 && x<=8.1)
206 {
207 Float_t c0=-1.655279e-1;
208 Float_t c1=6.307392e-2;
209 Float_t c2=-8.011441e-3;
210 Float_t c3=3.392126e-4;
211 sm=(c0+c1*x+c2*x*x+c3*x*x*x)*1.e-18;
212 }
213
214 if (x> 8.1)
215 {
216 Int_t j=0;
217 while (x<=em[j] && x>=em[j+1])
218 {
219 j++;
220 Float_t a=(sch4[j+1]-sch4[j])/(em[j+1]-em[j]);
221 sm=(sch4[j]+a*(x-em[j]))*1e-22;
222 }
223 }
224
225 Float_t dm=(kIgas1/100.)*(1.-((kOxy+kWater)/1.e6))*kLosch*pn/tn;
226 Float_t abslm=1./sm/dm;
227
228// ------- ISOBUTHANE CROSS SECTION --------------
229// i-C4H10 (ai) abs. length from curves in
230// Lu-McDonald paper for BARI RICH workshop .
231// -----------------------------------------------------------
232
233 Float_t ai;
234 Float_t absli;
235 if (kIgas2 != 0)
236 {
237 if (x<7.25)
238 ai=100000000.;
239
240 if(x>=7.25 && x<7.375)
241 ai=24.3;
242
243 if(x>=7.375)
244 ai=.0000000001;
245
246 Float_t si = 1./(ai*kLosch*273.16/293.); // ISOB. CRO.SEC.IN CM2
247 Float_t di=(kIgas2/100.)*(1.-((kOxy+kWater)/1.e6))*kLosch*pn/tn;
248 absli =1./si/di;
249 }
250 else
251 absli=1.e18;
252// ---------------------------------------------------------
253//
254// transmission of O2
255//
256// y= path in cm, x=energy in eV
257// so= cross section for UV absorption in cm2
258// do= O2 molecular density in cm-3
259// ---------------------------------------------------------
260
261 Float_t abslo;
262 Float_t so=0;
263 if(x>=6.0)
264 {
265 if(x>=6.0 && x<6.5)
266 {
267 so=3.392709e-13 * TMath::Exp(2.864104 *x);
268 so=so*1e-18;
269 }
270
271 if(x>=6.5 && x<7.0)
272 {
273 so=2.910039e-34 * TMath::Exp(10.3337*x);
274 so=so*1e-18;
275 }
276
277
278 if (x>=7.0)
279 {
280 Float_t a0=-73770.76;
281 Float_t a1=46190.69;
282 Float_t a2=-11475.44;
283 Float_t a3=1412.611;
284 Float_t a4=-86.07027;
285 Float_t a5=2.074234;
286 so= a0+(a1*x)+(a2*x*x)+(a3*x*x*x)+(a4*x*x*x*x)+(a5*x*x*x*x*x);
287 so=so*1e-18;
288 }
289
290 Float_t dox=(kOxy/1e6)*kLosch*pn/tn;
291 abslo=1./so/dox;
292 }
293 else
294 abslo=1.e18;
295// ---------------------------------------------------------
296//
297// transmission of H2O
298//
299// y= path in cm, x=energy in eV
300// sw= cross section for UV absorption in cm2
301// dw= H2O molecular density in cm-3
302// ---------------------------------------------------------
303
304 Float_t abslw;
305
306 Float_t b0=29231.65;
307 Float_t b1=-15807.74;
308 Float_t b2=3192.926;
309 Float_t b3=-285.4809;
310 Float_t b4=9.533944;
311
312 if(x>6.75)
313 {
314 Float_t sw= b0+(b1*x)+(b2*x*x)+(b3*x*x*x)+(b4*x*x*x*x);
315 sw=sw*1e-18;
316 Float_t dw=(kWater/1e6)*kLosch*pn/tn;
317 abslw=1./sw/dw;
318 }
319 else
320 abslw=1.e18;
321
322// ---------------------------------------------------------
323
324 Float_t alength=1./(1./abslm+1./absli+1./abslo+1./abslw);
325 return (alength);
326}//AbsoCH4
327
328
329Float_t Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)
330{
331
332 //ENE(EV), PDOTI=COS(INC.ANG.), PDOTR=COS(POL.PLANE ROT.ANG.)
333
334 Float_t en[36] = {5.0,5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6.0,6.1,6.2,
335 6.3,6.4,6.5,6.6,6.7,6.8,6.9,7.0,7.1,7.2,7.3,7.4,7.5,7.6,7.7,
336 7.8,7.9,8.0,8.1,8.2,8.3,8.4,8.5};
337
338
339 Float_t csin[36] = {2.14,2.21,2.33,2.48,2.76,2.97,2.99,2.59,2.81,3.05,
340 2.86,2.53,2.55,2.66,2.79,2.96,3.18,3.05,2.84,2.81,2.38,2.11,
341 2.01,2.13,2.39,2.73,3.08,3.15,2.95,2.73,2.56,2.41,2.12,1.95,
342 1.72,1.53};
343
344 Float_t csik[36] = {0.,0.,0.,0.,0.,0.196,0.408,0.208,0.118,0.49,0.784,0.543,
345 0.424,0.404,0.371,0.514,0.922,1.102,1.139,1.376,1.461,1.253,0.878,
346 0.69,0.612,0.649,0.824,1.347,1.571,1.678,1.763,1.857,1.824,1.824,
347 1.714,1.498};
348 Float_t xe=ene;
349 Int_t j=Int_t(xe*10)-49;
350 Float_t cn=csin[j]+((csin[j+1]-csin[j])/0.1)*(xe-en[j]);
351 Float_t ck=csik[j]+((csik[j+1]-csik[j])/0.1)*(xe-en[j]);
352
353 //FORMULAE FROM HANDBOOK OF OPTICS, 33.23 OR
354 //W.R. HUNTER, J.O.S.A. 54 (1964),15 , J.O.S.A. 55(1965),1197
355
356 Float_t sinin=TMath::Sqrt(1-pdoti*pdoti);
357 Float_t tanin=sinin/pdoti;
358
359 Float_t c1=cn*cn-ck*ck-sinin*sinin;
360 Float_t c2=4*cn*cn*ck*ck;
361 Float_t aO=TMath::Sqrt(0.5*(TMath::Sqrt(c1*c1+c2)+c1));
362 Float_t b2=0.5*(TMath::Sqrt(c1*c1+c2)-c1);
363
364 Float_t rs=((aO-pdoti)*(aO-pdoti)+b2)/((aO+pdoti)*(aO+pdoti)+b2);
365 Float_t rp=rs*((aO-sinin*tanin)*(aO-sinin*tanin)+b2)/((aO+sinin*tanin)*(aO+sinin*tanin)+b2);
366
367
368 //CORRECTION FACTOR FOR SURFACE ROUGHNESS
369 //B.J. STAGG APPLIED OPTICS, 30(1991),4113
370
371 Float_t sigraf=18.;
372 Float_t lamb=1240/ene;
373 Float_t fresn;
374
375 Float_t rO=TMath::Exp(-(4*TMath::Pi()*pdoti*sigraf/lamb)*(4*TMath::Pi()*pdoti*sigraf/lamb));
376
377 if(pola)
378 {
379 Float_t pdotr=0.8; //DEGREE OF POLARIZATION : 1->P , -1->S
380 fresn=0.5*(rp*(1+pdotr)+rs*(1-pdotr));
381 }
382 else
383 fresn=0.5*(rp+rs);
384
385 fresn = fresn*rO;
386 return(fresn);
387}//Fresnel(...)
388#endif