Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Toulon.hh
Go to the documentation of this file.
1 
2 #ifndef __JPHYSICS__TOULON__
3 #define __JPHYSICS__TOULON__
4 
5 #include <map>
6 #include <cmath>
7 
8 #include "JTools/JConstants.hh"
9 #include "JTools/JFunction1D_t.hh"
10 
11 
12 /**
13  * \author mdejong
14  */
15 
16 namespace TOULON {
17 
18  namespace {
19  using JTOOLS::PI;
20  using JTOOLS::JPolint1Function1D_t;
21  using JTOOLS::JPolint2Function1D_t;
22  using JTOOLS::JPolint3Function1D_t;
23  using JTOOLS::JGridPolint1Function1D_t;
24  }
25 
26 
27  /**
28  * Photo-cathode area 10 inch PMT.
29  *
30  * \return photo-cathode area [m^2]
31  */
32  inline const double getPhotocathodeArea()
33  {
34  return 440e-4; // photo-cathode area [m^2]
35  }
36 
37 
38  /**
39  * Absoption length
40  *
41  * \param lambda wavelength of light [nm]
42  * \return absorption length [m]
43  */
44  inline double getAbsorptionLength(const double lambda)
45  {
46  typedef std::map<double, double> map_t;
47 
48  static map_t zmap;
49 
50  if (zmap.empty()) {
51 
52  zmap[290.0] = 0.0;
53  zmap[300.0] = 7.194;
54  zmap[310.0] = 9.174;
55  zmap[320.0] = 10.571;
56  zmap[330.0] = 12.626;
57  zmap[340.0] = 14.085;
58  zmap[350.0] = 15.898;
59  zmap[360.0] = 18.939;
60  zmap[370.0] = 21.142;
61  zmap[380.0] = 24.096;
62  zmap[390.0] = 27.548;
63  zmap[400.0] = 30.769;
64  zmap[410.0] = 34.364;
65  zmap[420.0] = 39.216;
66  zmap[430.0] = 42.194;
67  zmap[440.0] = 45.872;
68  zmap[450.0] = 50.000;
69  zmap[460.0] = 52.356;
70  zmap[470.0] = 54.945;
71  zmap[480.0] = 54.945;
72  zmap[490.0] = 51.020;
73  zmap[500.0] = 38.910;
74  zmap[510.0] = 28.011;
75  zmap[520.0] = 20.964;
76  zmap[530.0] = 19.724;
77  zmap[540.0] = 17.921;
78  zmap[550.0] = 15.674;
79  zmap[560.0] = 14.124;
80  zmap[570.0] = 12.516;
81  zmap[580.0] = 9.259;
82  zmap[590.0] = 6.369;
83  zmap[600.0] = 4.098;
84  zmap[610.0] = 3.460;
85  zmap[620.0] = 0.0;
86  }
87 
88  const double x = lambda;
89 
90  if (x > zmap.begin()->first && x < zmap.rbegin()->first) {
91 
92  map_t::const_iterator i = zmap.lower_bound(x);
93  map_t::const_iterator j = i;
94 
95  --j;
96 
97  if (i == zmap.begin()) {
98  ++i;
99  ++j;
100  }
101 
102  const double x1 = i->first;
103  const double x2 = j->first;
104 
105  const double y1 = i->second;
106  const double y2 = j->second;
107 
108  return y1 + (y2 - y1) * (x - x1) / (x2 - x1);
109 
110  } else
111 
112  return 0.0;
113  }
114 
115 
116  /**
117  * Mie scattering length
118  *
119  * \param lambda wavelength of light [nm]
120  * \return scattering length [m]
121  */
122  inline double getMieScatteringLength(const double lambda)
123  {
124  typedef std::map<double, double> map_t;
125 
126  static map_t zmap;
127 
128  if (zmap.empty()) {
129 
130  zmap[300.0] = 32.290;
131  zmap[310.0] = 34.000;
132  zmap[320.0] = 35.730;
133  zmap[330.0] = 37.480;
134  zmap[340.0] = 39.260;
135  zmap[350.0] = 41.070;
136  zmap[360.0] = 42.890;
137  zmap[370.0] = 44.740;
138  zmap[380.0] = 46.600;
139  zmap[390.0] = 48.490;
140  zmap[400.0] = 50.390;
141  zmap[410.0] = 52.310;
142  zmap[420.0] = 54.250;
143  zmap[430.0] = 56.210;
144  zmap[440.0] = 58.180;
145  zmap[450.0] = 60.160;
146  zmap[460.0] = 62.160;
147  zmap[470.0] = 64.180;
148  zmap[480.0] = 66.210;
149  zmap[490.0] = 68.240;
150  zmap[500.0] = 70.300;
151  zmap[510.0] = 72.360;
152  zmap[520.0] = 74.430;
153  zmap[530.0] = 76.510;
154  zmap[540.0] = 78.610;
155  zmap[550.0] = 80.710;
156  zmap[560.0] = 82.820;
157  zmap[570.0] = 84.940;
158  zmap[580.0] = 87.070;
159  zmap[590.0] = 89.200;
160  zmap[600.0] = 91.340;
161  zmap[610.0] = 93.490;
162  zmap[620.0] = 95.000;
163  }
164 
165  const double x = lambda;
166 
167  if (x > zmap.begin()->first && x < zmap.rbegin()->first) {
168 
169  map_t::const_iterator i = zmap.lower_bound(x);
170  map_t::const_iterator j = i;
171 
172  --j;
173 
174  if (i == zmap.begin()) {
175  ++i;
176  ++j;
177  }
178 
179  const double x1 = i->first;
180  const double x2 = j->first;
181 
182  const double y1 = i->second;
183  const double y2 = j->second;
184 
185  return y1 + (y2 - y1) * (x - x1) / (x2 - x1);
186 
187  } else
188 
189  return 0.0;
190  }
191 
192 
193  /**
194  * Rayleigh scattering length
195  *
196  * \param lambda wavelength of light [nm]
197  * \return scattering length [m]
198  */
199  inline double getRayleighScatteringLength(const double lambda)
200  {
201  typedef std::map<double, double> map_t;
202 
203  static map_t zmap;
204 
205  if (zmap.empty()) {
206 
207  zmap[300.0] = 43.410;
208  zmap[310.0] = 49.990;
209  zmap[320.0] = 57.300;
210  zmap[330.0] = 65.400;
211  zmap[340.0] = 74.360;
212  zmap[350.0] = 84.230;
213  zmap[360.0] = 95.080;
214  zmap[370.0] = 106.970;
215  zmap[380.0] = 119.970;
216  zmap[390.0] = 134.140;
217  zmap[400.0] = 149.570;
218  zmap[410.0] = 166.330;
219  zmap[420.0] = 184.490;
220  zmap[430.0] = 204.130;
221  zmap[440.0] = 225.340;
222  zmap[450.0] = 248.200;
223  zmap[460.0] = 272.800;
224  zmap[470.0] = 299.230;
225  zmap[480.0] = 327.590;
226  zmap[490.0] = 357.960;
227  zmap[500.0] = 390.450;
228  zmap[510.0] = 425.150;
229  zmap[520.0] = 462.170;
230  zmap[530.0] = 501.620;
231  zmap[540.0] = 543.610;
232  zmap[550.0] = 588.240;
233  zmap[560.0] = 635.620;
234  zmap[570.0] = 685.890;
235  zmap[580.0] = 739.150;
236  zmap[590.0] = 795.530;
237  zmap[600.0] = 855.150;
238  zmap[610.0] = 918.140;
239  zmap[620.0] = 979.000;
240  }
241 
242  const double x = lambda;
243 
244  if (x > zmap.begin()->first && x < zmap.rbegin()->first) {
245 
246  map_t::const_iterator i = zmap.lower_bound(x);
247  map_t::const_iterator j = i;
248 
249  --j;
250 
251  if (i == zmap.begin()) {
252  ++i;
253  ++j;
254  }
255 
256  const double x1 = i->first;
257  const double x2 = j->first;
258 
259  const double y1 = i->second;
260  const double y2 = j->second;
261 
262  return y1 + (y2 - y1) * (x - x1) / (x2 - x1);
263 
264  } else
265 
266  return 0.0;
267  }
268 
269 
270  /**
271  * Auxiliary method to describe light scattering in water (Rayleigh)
272  *
273  * \param a angular dependence parameter
274  * \param x cosine scattering angle
275  * \return probability
276  */
277  inline double rayleigh(const double a,
278  const double x)
279 
280  {
281  const double a0 = 1.0 / (1.0 + a/3.0) / (4*PI);
282 
283  return a0 * (1.0 + a*x*x);
284  }
285 
286 
287  /**
288  * Auxiliary method to describe light scattering in water (Rayleigh)
289  *
290  * \param x cosine scattering angle
291  * \return probability
292  */
293  inline double rayleigh(const double x)
294  {
295  return rayleigh(0.853, x);
296  }
297 
298 
299  /**
300  * Function to describe light scattering in water.
301  *
302  * \param ct cosine scattering angle
303  * \return probability
304  */
305  inline double getPetzholdScatteringProbability(const double ct)
306  {
307  static JPolint2Function1D_t f1;
308 
309  if (f1.empty()) {
310 
311  f1[+0.0017] = +19.377;
312  f1[+0.0022] = +17.907;
313  f1[+0.0028] = +16.464;
314  f1[+0.0035] = +15.333;
315  f1[+0.0044] = +14.148;
316  f1[+0.0055] = +13.043;
317  f1[+0.0069] = +12.059;
318  f1[+0.0087] = +12.021;
319  f1[+0.0110] = +9.992;
320  f1[+0.0139] = +8.898;
321  f1[+0.0175] = +7.853;
322  f1[+0.0220] = +6.845;
323  f1[+0.0277] = +5.900;
324  f1[+0.0348] = +4.989;
325  f1[+0.0438] = +4.175;
326  f1[+0.0552] = +3.473;
327  f1[+0.0695] = +2.870;
328  f1[+0.0875] = +2.358;
329  f1[+0.1101] = +1.938;
330  f1[+0.1386] = +1.579;
331  f1[+0.1745] = +1.258;
332  f1[+0.2618] = +0.796;
333  f1[+0.3491] = +0.525;
334  f1[+0.4363] = +0.391;
335  f1[+0.5236] = +0.270;
336  f1[+0.6109] = +0.214;
337  f1[+0.6981] = +0.170;
338  f1[+0.7854] = +0.136;
339  f1[+0.8727] = +0.110;
340  f1[+0.9599] = +0.087;
341  f1[+1.0472] = +0.071;
342  f1[+1.1345] = +0.060;
343  f1[+1.2217] = +0.050;
344  f1[+1.3090] = +0.042;
345  f1[+1.3963] = +0.036;
346  f1[+1.4835] = +0.031;
347  f1[+1.5708] = +0.027;
348  f1[+1.6581] = +0.024;
349  f1[+1.7453] = +0.021;
350  f1[+1.8326] = +0.019;
351  f1[+1.9199] = +0.017;
352  f1[+2.0071] = +0.016;
353  f1[+2.0944] = +0.015;
354  f1[+2.1817] = +0.013;
355  f1[+2.2689] = +0.012;
356  f1[+2.3562] = +0.011;
357  f1[+2.4435] = +0.009;
358  f1[+2.5307] = +0.008;
359  f1[+2.6180] = +0.007;
360  f1[+2.7053] = +0.007;
361  f1[+2.7925] = +0.006;
362  f1[+2.8798] = +0.005;
363  f1[+2.9671] = +0.003;
364  f1[+3.0543] = +0.002;
365  f1[+3.1416] = -0.000;
366 
367  for (JPolint1Function1D_t::iterator i = f1.begin(); i != f1.end(); ++i) {
368 
369  const double x = i->getX();
370  const double y = i->getY();
371 
372  i->setValue(x*y);
373  }
374 
375  f1.compile();
376 
377  f1.setExceptionHandler(new JPolint2Function1D_t::JDefaultResult(0.0));
378  }
379 
380  const double x = acos(ct);
381 
382  return f1(x) / x;
383  }
384 
385 
386  /**
387  * Angular acceptence of Antares PMT
388  *
389  * \param x cosine of angle of incidence
390  * \return probability
391  */
392  inline double getAngularAcceptance(const double x)
393  {
394  static JGridPolint1Function1D_t f1;
395 
396  if (f1.empty()) {
397 
398  f1[-1.000] = +0.997;
399  f1[-0.995] = +0.992;
400  f1[-0.990] = +0.986;
401  f1[-0.985] = +0.984;
402  f1[-0.980] = +0.983;
403  f1[-0.975] = +0.981;
404  f1[-0.970] = +0.976;
405  f1[-0.965] = +0.976;
406  f1[-0.960] = +0.974;
407  f1[-0.955] = +0.971;
408  f1[-0.950] = +0.970;
409  f1[-0.945] = +0.967;
410  f1[-0.940] = +0.964;
411  f1[-0.935] = +0.960;
412  f1[-0.930] = +0.956;
413  f1[-0.925] = +0.953;
414  f1[-0.920] = +0.949;
415  f1[-0.915] = +0.944;
416  f1[-0.910] = +0.940;
417  f1[-0.905] = +0.936;
418  f1[-0.900] = +0.934;
419  f1[-0.895] = +0.930;
420  f1[-0.890] = +0.926;
421  f1[-0.885] = +0.923;
422  f1[-0.880] = +0.918;
423  f1[-0.875] = +0.911;
424  f1[-0.870] = +0.905;
425  f1[-0.865] = +0.900;
426  f1[-0.860] = +0.895;
427  f1[-0.855] = +0.891;
428  f1[-0.850] = +0.886;
429  f1[-0.845] = +0.881;
430  f1[-0.840] = +0.876;
431  f1[-0.835] = +0.871;
432  f1[-0.830] = +0.866;
433  f1[-0.825] = +0.862;
434  f1[-0.820] = +0.859;
435  f1[-0.815] = +0.855;
436  f1[-0.810] = +0.850;
437  f1[-0.805] = +0.844;
438  f1[-0.800] = +0.838;
439  f1[-0.795] = +0.834;
440  f1[-0.790] = +0.831;
441  f1[-0.785] = +0.828;
442  f1[-0.780] = +0.823;
443  f1[-0.775] = +0.819;
444  f1[-0.770] = +0.815;
445  f1[-0.765] = +0.811;
446  f1[-0.760] = +0.808;
447  f1[-0.755] = +0.803;
448  f1[-0.750] = +0.798;
449  f1[-0.745] = +0.793;
450  f1[-0.740] = +0.789;
451  f1[-0.735] = +0.785;
452  f1[-0.730] = +0.782;
453  f1[-0.725] = +0.779;
454  f1[-0.720] = +0.776;
455  f1[-0.715] = +0.771;
456  f1[-0.710] = +0.767;
457  f1[-0.705] = +0.764;
458  f1[-0.700] = +0.763;
459  f1[-0.695] = +0.763;
460  f1[-0.690] = +0.761;
461  f1[-0.685] = +0.759;
462  f1[-0.680] = +0.755;
463  f1[-0.675] = +0.749;
464  f1[-0.670] = +0.743;
465  f1[-0.665] = +0.739;
466  f1[-0.660] = +0.736;
467  f1[-0.655] = +0.734;
468  f1[-0.650] = +0.733;
469  f1[-0.645] = +0.731;
470  f1[-0.640] = +0.728;
471  f1[-0.635] = +0.723;
472  f1[-0.630] = +0.719;
473  f1[-0.625] = +0.716;
474  f1[-0.620] = +0.714;
475  f1[-0.615] = +0.712;
476  f1[-0.610] = +0.710;
477  f1[-0.605] = +0.706;
478  f1[-0.600] = +0.703;
479  f1[-0.595] = +0.699;
480  f1[-0.590] = +0.696;
481  f1[-0.585] = +0.693;
482  f1[-0.580] = +0.690;
483  f1[-0.575] = +0.688;
484  f1[-0.570] = +0.685;
485  f1[-0.565] = +0.683;
486  f1[-0.560] = +0.681;
487  f1[-0.555] = +0.679;
488  f1[-0.550] = +0.677;
489  f1[-0.545] = +0.675;
490  f1[-0.540] = +0.672;
491  f1[-0.535] = +0.668;
492  f1[-0.530] = +0.665;
493  f1[-0.525] = +0.661;
494  f1[-0.520] = +0.658;
495  f1[-0.515] = +0.654;
496  f1[-0.510] = +0.651;
497  f1[-0.505] = +0.647;
498  f1[-0.500] = +0.643;
499  f1[-0.495] = +0.640;
500  f1[-0.490] = +0.637;
501  f1[-0.485] = +0.633;
502  f1[-0.480] = +0.629;
503  f1[-0.475] = +0.626;
504  f1[-0.470] = +0.622;
505  f1[-0.465] = +0.618;
506  f1[-0.460] = +0.615;
507  f1[-0.455] = +0.612;
508  f1[-0.450] = +0.608;
509  f1[-0.445] = +0.604;
510  f1[-0.440] = +0.600;
511  f1[-0.435] = +0.596;
512  f1[-0.430] = +0.592;
513  f1[-0.425] = +0.588;
514  f1[-0.420] = +0.584;
515  f1[-0.415] = +0.579;
516  f1[-0.410] = +0.574;
517  f1[-0.405] = +0.569;
518  f1[-0.400] = +0.564;
519  f1[-0.395] = +0.560;
520  f1[-0.390] = +0.557;
521  f1[-0.385] = +0.555;
522  f1[-0.380] = +0.553;
523  f1[-0.375] = +0.550;
524  f1[-0.370] = +0.546;
525  f1[-0.365] = +0.540;
526  f1[-0.360] = +0.534;
527  f1[-0.355] = +0.530;
528  f1[-0.350] = +0.527;
529  f1[-0.345] = +0.524;
530  f1[-0.340] = +0.521;
531  f1[-0.335] = +0.518;
532  f1[-0.330] = +0.514;
533  f1[-0.325] = +0.510;
534  f1[-0.320] = +0.506;
535  f1[-0.315] = +0.502;
536  f1[-0.310] = +0.498;
537  f1[-0.305] = +0.495;
538  f1[-0.300] = +0.491;
539  f1[-0.295] = +0.488;
540  f1[-0.290] = +0.485;
541  f1[-0.285] = +0.481;
542  f1[-0.280] = +0.478;
543  f1[-0.275] = +0.474;
544  f1[-0.270] = +0.471;
545  f1[-0.265] = +0.468;
546  f1[-0.260] = +0.465;
547  f1[-0.255] = +0.461;
548  f1[-0.250] = +0.458;
549  f1[-0.245] = +0.454;
550  f1[-0.240] = +0.451;
551  f1[-0.235] = +0.447;
552  f1[-0.230] = +0.445;
553  f1[-0.225] = +0.442;
554  f1[-0.220] = +0.439;
555  f1[-0.215] = +0.436;
556  f1[-0.210] = +0.433;
557  f1[-0.205] = +0.429;
558  f1[-0.200] = +0.425;
559  f1[-0.195] = +0.421;
560  f1[-0.190] = +0.418;
561  f1[-0.185] = +0.414;
562  f1[-0.180] = +0.411;
563  f1[-0.175] = +0.408;
564  f1[-0.170] = +0.405;
565  f1[-0.165] = +0.402;
566  f1[-0.160] = +0.399;
567  f1[-0.155] = +0.396;
568  f1[-0.150] = +0.393;
569  f1[-0.145] = +0.389;
570  f1[-0.140] = +0.386;
571  f1[-0.135] = +0.382;
572  f1[-0.130] = +0.379;
573  f1[-0.125] = +0.375;
574  f1[-0.120] = +0.372;
575  f1[-0.115] = +0.369;
576  f1[-0.110] = +0.366;
577  f1[-0.105] = +0.363;
578  f1[-0.100] = +0.359;
579  f1[-0.095] = +0.356;
580  f1[-0.090] = +0.353;
581  f1[-0.085] = +0.350;
582  f1[-0.080] = +0.347;
583  f1[-0.075] = +0.345;
584  f1[-0.070] = +0.342;
585  f1[-0.065] = +0.339;
586  f1[-0.060] = +0.336;
587  f1[-0.055] = +0.333;
588  f1[-0.050] = +0.330;
589  f1[-0.045] = +0.327;
590  f1[-0.040] = +0.325;
591  f1[-0.035] = +0.322;
592  f1[-0.030] = +0.319;
593  f1[-0.025] = +0.316;
594  f1[-0.020] = +0.313;
595  f1[-0.015] = +0.310;
596  f1[-0.010] = +0.308;
597  f1[-0.005] = +0.305;
598  f1[+0.000] = +0.302;
599  f1[+0.005] = +0.299;
600  f1[+0.010] = +0.296;
601  f1[+0.015] = +0.293;
602  f1[+0.020] = +0.289;
603  f1[+0.025] = +0.286;
604  f1[+0.030] = +0.283;
605  f1[+0.035] = +0.280;
606  f1[+0.040] = +0.277;
607  f1[+0.045] = +0.275;
608  f1[+0.050] = +0.272;
609  f1[+0.055] = +0.270;
610  f1[+0.060] = +0.268;
611  f1[+0.065] = +0.265;
612  f1[+0.070] = +0.262;
613  f1[+0.075] = +0.259;
614  f1[+0.080] = +0.256;
615  f1[+0.085] = +0.252;
616  f1[+0.090] = +0.249;
617  f1[+0.095] = +0.246;
618  f1[+0.100] = +0.244;
619  f1[+0.105] = +0.241;
620  f1[+0.110] = +0.239;
621  f1[+0.115] = +0.237;
622  f1[+0.120] = +0.234;
623  f1[+0.125] = +0.231;
624  f1[+0.130] = +0.229;
625  f1[+0.135] = +0.226;
626  f1[+0.140] = +0.224;
627  f1[+0.145] = +0.221;
628  f1[+0.150] = +0.218;
629  f1[+0.155] = +0.215;
630  f1[+0.160] = +0.212;
631  f1[+0.165] = +0.209;
632  f1[+0.170] = +0.206;
633  f1[+0.175] = +0.204;
634  f1[+0.180] = +0.202;
635  f1[+0.185] = +0.200;
636  f1[+0.190] = +0.199;
637  f1[+0.195] = +0.196;
638  f1[+0.200] = +0.194;
639  f1[+0.205] = +0.190;
640  f1[+0.210] = +0.187;
641  f1[+0.215] = +0.185;
642  f1[+0.220] = +0.182;
643  f1[+0.225] = +0.180;
644  f1[+0.230] = +0.177;
645  f1[+0.235] = +0.175;
646  f1[+0.240] = +0.173;
647  f1[+0.245] = +0.171;
648  f1[+0.250] = +0.169;
649  f1[+0.255] = +0.166;
650  f1[+0.260] = +0.164;
651  f1[+0.265] = +0.161;
652  f1[+0.270] = +0.159;
653  f1[+0.275] = +0.156;
654  f1[+0.280] = +0.154;
655  f1[+0.285] = +0.152;
656  f1[+0.290] = +0.150;
657  f1[+0.295] = +0.147;
658  f1[+0.300] = +0.145;
659  f1[+0.305] = +0.143;
660  f1[+0.310] = +0.141;
661  f1[+0.315] = +0.138;
662  f1[+0.320] = +0.136;
663  f1[+0.325] = +0.133;
664  f1[+0.330] = +0.131;
665  f1[+0.335] = +0.128;
666  f1[+0.340] = +0.126;
667  f1[+0.345] = +0.123;
668  f1[+0.350] = +0.121;
669  f1[+0.355] = +0.118;
670  f1[+0.360] = +0.115;
671  f1[+0.365] = +0.113;
672  f1[+0.370] = +0.110;
673  f1[+0.375] = +0.107;
674  f1[+0.380] = +0.103;
675  f1[+0.385] = +0.100;
676  f1[+0.390] = +0.097;
677  f1[+0.395] = +0.095;
678  f1[+0.400] = +0.092;
679  f1[+0.405] = +0.089;
680  f1[+0.410] = +0.086;
681  f1[+0.415] = +0.082;
682  f1[+0.420] = +0.079;
683  f1[+0.425] = +0.076;
684  f1[+0.430] = +0.074;
685  f1[+0.435] = +0.072;
686  f1[+0.440] = +0.069;
687  f1[+0.445] = +0.066;
688  f1[+0.450] = +0.063;
689  f1[+0.455] = +0.060;
690  f1[+0.460] = +0.057;
691  f1[+0.465] = +0.054;
692  f1[+0.470] = +0.051;
693  f1[+0.475] = +0.048;
694  f1[+0.480] = +0.046;
695  f1[+0.485] = +0.043;
696  f1[+0.490] = +0.040;
697  f1[+0.495] = +0.038;
698  f1[+0.500] = +0.035;
699  f1[+0.505] = +0.032;
700  f1[+0.510] = +0.029;
701  f1[+0.515] = +0.027;
702  f1[+0.520] = +0.024;
703  f1[+0.525] = +0.022;
704  f1[+0.530] = +0.019;
705  f1[+0.535] = +0.017;
706  f1[+0.540] = +0.015;
707  f1[+0.545] = +0.013;
708  f1[+0.550] = +0.011;
709  f1[+0.555] = +0.009;
710  f1[+0.560] = +0.008;
711  f1[+0.565] = +0.006;
712  f1[+0.570] = +0.005;
713  f1[+0.575] = +0.003;
714  f1[+0.580] = +0.002;
715  f1[+0.585] = +0.002;
716  f1[+0.590] = +0.001;
717  f1[+0.595] = +0.001;
718  f1[+0.600] = +0.001;
719  f1[+0.605] = +0.001;
720  f1[+0.610] = +0.001;
721  f1[+0.615] = +0.001;
722  f1[+0.620] = +0.001;
723  f1[+0.625] = +0.001;
724  f1[+0.630] = +0.001;
725  f1[+0.635] = +0.000;
726  f1[+0.640] = +0.000;
727  f1[+0.645] = +0.000;
728  f1[+0.650] = +0.000;
729  f1[+0.655] = +0.000;
730  f1[+0.660] = +0.000;
731  f1[+0.665] = +0.000;
732  f1[+0.670] = +0.000;
733  f1[+0.675] = +0.000;
734  f1[+0.680] = +0.000;
735  f1[+0.685] = +0.000;
736  f1[+0.690] = +0.000;
737  f1[+0.695] = +0.000;
738  f1[+0.700] = +0.000;
739  f1[+0.705] = +0.000;
740  f1[+0.710] = +0.000;
741  f1[+0.715] = +0.000;
742  f1[+0.720] = +0.000;
743  f1[+0.725] = +0.000;
744  f1[+0.730] = +0.000;
745  f1[+0.735] = +0.000;
746  f1[+0.740] = +0.000;
747  f1[+0.745] = +0.000;
748  f1[+0.750] = +0.000;
749  f1[+0.755] = +0.000;
750  f1[+0.760] = +0.000;
751  f1[+0.765] = +0.000;
752  f1[+0.770] = +0.000;
753  f1[+0.775] = +0.000;
754  f1[+0.780] = +0.000;
755  f1[+0.785] = +0.000;
756  f1[+0.790] = +0.000;
757  f1[+0.795] = +0.000;
758  f1[+0.800] = +0.000;
759  f1[+0.805] = +0.000;
760  f1[+0.810] = +0.000;
761  f1[+0.815] = +0.000;
762  f1[+0.820] = +0.000;
763  f1[+0.825] = +0.000;
764  f1[+0.830] = +0.000;
765  f1[+0.835] = +0.000;
766  f1[+0.840] = +0.000;
767  f1[+0.845] = +0.000;
768  f1[+0.850] = +0.000;
769  f1[+0.855] = +0.000;
770  f1[+0.860] = +0.000;
771  f1[+0.865] = +0.000;
772  f1[+0.870] = +0.000;
773  f1[+0.875] = +0.000;
774  f1[+0.880] = +0.000;
775  f1[+0.885] = +0.000;
776  f1[+0.890] = +0.000;
777  f1[+0.895] = +0.000;
778  f1[+0.900] = +0.000;
779  f1[+0.905] = +0.000;
780  f1[+0.910] = +0.000;
781  f1[+0.915] = +0.000;
782  f1[+0.920] = +0.000;
783  f1[+0.925] = +0.000;
784  f1[+0.930] = +0.000;
785  f1[+0.935] = +0.000;
786  f1[+0.940] = +0.000;
787  f1[+0.945] = +0.000;
788  f1[+0.950] = +0.000;
789  f1[+0.955] = +0.000;
790  f1[+0.960] = +0.000;
791  f1[+0.965] = +0.000;
792  f1[+0.970] = +0.000;
793  f1[+0.975] = +0.000;
794  f1[+0.980] = +0.000;
795  f1[+0.985] = +0.000;
796  f1[+0.990] = +0.000;
797  f1[+0.995] = +0.000;
798  f1[+1.000] = -0.000;
799 
800  f1.compile();
801 
802  f1.setExceptionHandler(new JGridPolint1Function1D_t::JDefaultResult(0.0));
803  }
804 
805 
806  return f1(x);
807  }
808 
809 
810  /**
811  * Quantum efficiency of 10-inch Hamamatsu PMT
812  *
813  * \param lambda wavelength of photon [nm]
814  * \param option include absorption in glass and gel (if true, otherwise not)
815  * \return quantum efficiency
816  */
817  inline double getQE(const double lambda, const bool option)
818  {
819  class tuple {
820  public:
821  tuple(const double __QE,
822  const double __l_gel,
823  const double __l_glass) :
824  QE (__QE),
825  l_gel (__l_gel),
826  l_glass(__l_glass)
827  {}
828 
829  double QE; // Quantum efficiiency
830  double l_gel; // gel absorption length [cm]
831  double l_glass; // glass absorption length [cm]
832  };
833 
834  static const tuple ntuple[] = {
835  tuple(0.000e-2, 0.00, 0.00),
836  tuple(1.988e-2, 100.81, 148.37),
837  tuple(2.714e-2, 99.94, 142.87),
838  tuple(3.496e-2, 99.89, 135.64),
839  tuple(4.347e-2, 96.90, 134.58),
840  tuple(5.166e-2, 96.42, 138.27),
841  tuple(6.004e-2, 94.36, 142.40),
842  tuple(6.885e-2, 89.09, 147.16),
843  tuple(8.105e-2, 90.10, 151.80),
844  tuple(10.13e-2, 86.95, 150.88),
845  tuple(13.03e-2, 85.88, 145.68),
846  tuple(15.29e-2, 84.49, 139.70),
847  tuple(16.37e-2, 81.08, 126.55),
848  tuple(17.11e-2, 78.18, 118.86),
849  tuple(17.86e-2, 76.48, 113.90),
850  tuple(18.95e-2, 74.55, 116.08),
851  tuple(20.22e-2, 72.31, 109.23),
852  tuple(21.26e-2, 68.05, 81.63),
853  tuple(22.10e-2, 66.91, 65.66),
854  tuple(22.65e-2, 64.48, 77.30),
855  tuple(23.07e-2, 62.53, 73.02),
856  tuple(23.14e-2, 59.38, 81.25),
857  tuple(23.34e-2, 56.64, 128.04),
858  tuple(22.95e-2, 53.29, 61.84),
859  tuple(22.95e-2, 48.96, 19.23),
860  tuple(22.74e-2, 45.71, 27.21),
861  tuple(23.48e-2, 41.88, 18.09),
862  tuple(22.59e-2, 37.14, 8.41),
863  tuple(20.61e-2, 30.49, 3.92),
864  tuple(17.68e-2, 23.08, 1.82),
865  tuple(13.18e-2, 15.60, 0.84),
866  tuple(7.443e-2, 8.00, 0.39),
867  tuple(2.526e-2, 0.00, 0.17),
868  tuple(0.000e-2, 0.00, 0.00)
869  };
870 
871  static const double cola = 0.9; // collection efficiency
872  static const double x_glass = 1.5; // glass thickness [cm]
873  static const double x_gel = 1.0; // gel thickness [cm]
874 
875  // ntuple
876  static const int N = sizeof(ntuple) / sizeof(ntuple[0]) - 1;
877 
878  static const double xmax = 620.0; // maximal wavelength [nm] (tuple[ 0 ])
879  static const double xmin = 290.0; // minimal wavelength [nm] (tuple[N-1])
880 
881  const double x = lambda;
882 
883  double y = 0.0;
884 
885  if (x > xmin && x < xmax) {
886 
887  const int i = (int) (N * (x - xmax) / (xmin - xmax));
888  const int j = (i == N ? i - 1 : i + 1);
889 
890  const double x1 = xmax + i * (xmin - xmax) / N;
891  const double x2 = xmax + j * (xmin - xmax) / N;
892 
893  const double dx = (x - x1) / (x1 - x2);
894 
895  const double QE = ntuple[i].QE + (ntuple[i].QE - ntuple[j].QE ) * dx;
896  const double l_gel = ntuple[i].l_gel + (ntuple[i].l_gel - ntuple[j].l_gel ) * dx;
897  const double l_glass = ntuple[i].l_glass + (ntuple[i].l_glass - ntuple[j].l_glass) * dx;
898 
899  y = cola * QE;
900 
901  if (option) {
902 
903  if (l_glass > 0.0 && l_gel > 0.0)
904  y *= exp(-x_glass/l_glass) * exp(-x_gel/l_gel);
905  else
906  y = 0.0;
907  }
908  }
909 
910  return y;
911  }
912 
913 
914  /**
915  * Quantum efficiency of 10-inch Hamamatsu PMT
916  *
917  * \param lambda wavelength of photon [nm]
918  * \return quantum efficiency
919  */
920  inline double getQE(const double lambda)
921  {
922  return getQE(lambda, true);
923  }
924 }
925 
926 
927 #endif
JTOOLS::PI
static const double PI
Constants.
Definition: JConstants.hh:20
JTOOLS::JPolint2Function1D_t
Type definition of a 2nd degree polynomial interpolation with result type double. ...
Definition: JFunction1D_t.hh:193
a
fi JEventTimesliceWriter a
Definition: JEventProcessor.sh:71
JTOOLS::JPolint3Function1D_t
Type definition of a 3rd degree polynomial interpolation with result type double. ...
Definition: JFunction1D_t.hh:201
JTOOLS::JGridPolint1Function1D_t
Type definition of a 1st degree polynomial interpolation based on a JGridCollection with result type ...
Definition: JFunction1D_t.hh:292
JFunction1D_t.hh
JConstants.hh
Constants.
TOULON::getQE
double getQE(const double lambda, const bool option)
Quantum efficiency of 10-inch Hamamatsu PMT.
Definition: Toulon.hh:817
TOULON::getPhotocathodeArea
const double getPhotocathodeArea()
Photo-cathode area 10 inch PMT.
Definition: Toulon.hh:32
TOULON::getAngularAcceptance
double getAngularAcceptance(const double x)
Angular acceptence of Antares PMT.
Definition: Toulon.hh:392
JTOOLS::JPolint1Function1D_t
Type definition of a 1st degree polynomial interpolation with result type double. ...
Definition: JFunction1D_t.hh:185
TOULON::getRayleighScatteringLength
double getRayleighScatteringLength(const double lambda)
Rayleigh scattering length.
Definition: Toulon.hh:199
TOULON::rayleigh
double rayleigh(const double a, const double x)
Auxiliary method to describe light scattering in water (Rayleigh)
Definition: Toulon.hh:277
JTOOLS::j
int j
Definition: JPolint.hh:634
TOULON::getMieScatteringLength
double getMieScatteringLength(const double lambda)
Mie scattering length.
Definition: Toulon.hh:122
N
then usage $script[input file[working directory[option]]] nWhere option can be N
Definition: JMuonPostfit.sh:37
TOULON::getAbsorptionLength
double getAbsorptionLength(const double lambda)
Absoption length.
Definition: Toulon.hh:44
TOULON::getPetzholdScatteringProbability
double getPetzholdScatteringProbability(const double ct)
Function to describe light scattering in water.
Definition: Toulon.hh:305
exp
then set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable OUTPUT_FILE histogram.root JHistogram1D -o $WORKDIR/$OUTPUT_FILE -F "$FORMULA" -
Generated by   doxygen 1.8.5