Jpp  18.3.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JPhysicsSupportkit.hh
Go to the documentation of this file.
1 #ifndef __JPHYSICS_JPHYSICSSUPPORTKIT__
2 #define __JPHYSICS_JPHYSICSSUPPORTKIT__
3 
4 #include <cmath>
5 
7 
8 #include "JPhysics/JConstants.hh"
9 
10 
11 /**
12  * \author mdejong
13  */
14 
15 namespace JPHYSICS {}
16 namespace JPP { using namespace JPHYSICS; }
17 
18 /**
19  * Auxiliary methods for light properties of deep-sea water.
20  */
21 namespace JPHYSICS {
22 
25 
26 
27  /**
28  * Get multiple Coulomb scattering angle.
29  *
30  * The formula is taken from reference:
31  * Particle Data Book, formula 27.14.
32  *
33  * \param E Energy [GeV]
34  * \param x distance [m]
35  * \param X0 radiation length [m]
36  * \param M mass [GeV]
37  * \param Q charge [unit]
38  * \return angle [rad]
39  */
40  inline double getThetaMCS(const double E,
41  const double x,
42  const double X0,
43  const double M,
44  const double Q)
45  {
46  if (E > M) {
47 
48  const double p = sqrt((E + M) * (E - M));
49  const double beta = p / E;
50 
51  return THETA_MCS * Q * sqrt(x/X0) * (1.0 + 0.038*log(x/X0)) / (beta*p);
52  }
53 
54  return 0.0;
55  }
56 
57 
58  /**
59  * Get multiple Coulomb scattering angle for muon.
60  *
61  * \param E Energy [GeV]
62  * \param x distance [m]
63  * \return angle [rad]
64  */
65  inline double getThetaMCS(const double E, const double x)
66  {
67  return getThetaMCS(E, x, X0_WATER_M, MASS_MUON, 1.0);
68  }
69 
70 
71  /**
72  * Auxiliary method to describe light scattering in water (Henyey-Greenstein).
73  *
74  * \param g angular dependence parameter
75  * \param x cosine scattering angle
76  * \return probability
77  */
78  inline double henyey_greenstein(const double g,
79  const double x)
80 
81  {
82  const double a0 = (1.0 - g*g) / (4*PI);
83  const double y = 1.0 + g*g - 2.0*g*x;
84 
85  return a0 / (y*sqrt(y));
86  }
87 
88 
89  /**
90  * Auxiliary method to describe light scattering in water (Heneyey-Greenstein).
91  *
92  * \param x cosine scattering angle
93  * \return probability
94  */
95  inline double henyey_greenstein(const double x)
96  {
97  static const double g = 0.924;
98 
99  return henyey_greenstein(g, x);
100  }
101 
102 
103  /**
104  * Auxiliary method to describe light scattering in water (Rayleigh).
105  *
106  * \param a angular dependence parameter
107  * \param x cosine scattering angle
108  * \return probability
109  */
110  inline double rayleigh(const double a,
111  const double x)
112 
113  {
114  const double a0 = 1.0 / (1.0 + a/3.0) / (4*PI);
115 
116  return a0 * (1.0 + a*x*x);
117  }
118 
119 
120  /**
121  * Auxiliary method to describe light scattering in water (Rayleigh).
122  *
123  * \param x cosine scattering angle
124  * \return probability
125  */
126  inline double rayleigh(const double x)
127  {
128  return rayleigh(0.835, x);
129  }
130 
131 
132  /**
133  * Model specific function to describe light scattering in water (f4).
134  *
135  * \param x cosine scattering angle
136  * \return probability
137  */
138  inline double f4(const double x)
139  {
140  static const double g1 = 0.77;
141  static const double g2 = 0.00;
142  static const double f = 1.00;
143 
144  return f * henyey_greenstein(g1,x) + (1.0 - f) * henyey_greenstein(g2,x);
145  }
146 
147 
148  /**
149  * Model specific function to describe light scattering in water (p00075).
150  *
151  * \param x cosine scattering angle
152  * \return probability
153  */
154  inline double p00075(const double x)
155  {
156  static const double g = 0.924;
157  static const double f = 0.17;
158 
159  return f * rayleigh(x) + (1.0 - f) * henyey_greenstein(g,x);
160  }
161 
162 
163  /**
164  * Measurement of light scattering in water.
165  *
166  * Values are taken from reference C.D. Mobley "Light and Water", ISBN 0-12-502750-8, pag. 111, Table 3.10 (right column).
167  */
168  class JPetzhold :
169  public JPolint1Function1D_t
170  {
171  public:
172  /**
173  * Default constructor.
174  */
176  {
177  // angle probability
178  // [deg]
179  (*this)[ 0.100] = 1.767e+03;
180  (*this)[ 0.126] = 1.296e+03;
181  (*this)[ 0.158] = 9.502e+02;
182  (*this)[ 0.200] = 6.991e+02;
183  (*this)[ 0.251] = 5.140e+02;
184  (*this)[ 0.316] = 3.764e+02;
185  (*this)[ 0.398] = 2.763e+02;
186  (*this)[ 0.501] = 2.188e+02;
187  (*this)[ 0.631] = 1.444e+02;
188  (*this)[ 0.794] = 1.022e+02;
189  (*this)[ 1.000] = 7.161e+01;
190  (*this)[ 1.259] = 4.958e+01;
191  (*this)[ 1.585] = 3.395e+01;
192  (*this)[ 1.995] = 2.281e+01;
193  (*this)[ 2.512] = 1.516e+01;
194  (*this)[ 3.162] = 1.002e+01;
195  (*this)[ 3.981] = 6.580e+00;
196  (*this)[ 5.012] = 4.295e+00;
197  (*this)[ 6.310] = 2.807e+00;
198  (*this)[ 7.943] = 1.819e+00;
199  (*this)[ 10.000] = 1.153e+00;
200  (*this)[ 15.000] = 4.893e-01;
201  (*this)[ 20.000] = 2.444e-01;
202  (*this)[ 25.000] = 1.472e-01;
203  (*this)[ 30.000] = 8.609e-02;
204  (*this)[ 35.000] = 5.931e-02;
205  (*this)[ 40.000] = 4.210e-02;
206  (*this)[ 45.000] = 3.067e-02;
207  (*this)[ 50.000] = 2.275e-02;
208  (*this)[ 55.000] = 1.699e-02;
209  (*this)[ 60.000] = 1.313e-02;
210  (*this)[ 65.000] = 1.046e-02;
211  (*this)[ 70.000] = 8.488e-03;
212  (*this)[ 75.000] = 6.976e-03;
213  (*this)[ 80.000] = 5.842e-03;
214  (*this)[ 85.000] = 4.953e-03;
215  (*this)[ 90.000] = 4.292e-03;
216  (*this)[ 95.000] = 3.782e-03;
217  (*this)[100.000] = 3.404e-03;
218  (*this)[105.000] = 3.116e-03;
219  (*this)[110.000] = 2.912e-03;
220  (*this)[115.000] = 2.797e-03;
221  (*this)[120.000] = 2.686e-03;
222  (*this)[125.000] = 2.571e-03;
223  (*this)[130.000] = 2.476e-03;
224  (*this)[135.000] = 2.377e-03;
225  (*this)[140.000] = 2.329e-03;
226  (*this)[145.000] = 2.313e-03;
227  (*this)[150.000] = 2.365e-03;
228  (*this)[155.000] = 2.506e-03;
229  (*this)[160.000] = 2.662e-03;
230  (*this)[165.000] = 2.835e-03;
231  (*this)[170.000] = 3.031e-03;
232  (*this)[175.000] = 3.092e-03;
233  (*this)[180.000] = 3.154e-03;
234 
235  this->mul(1.166); // normalisation
236  this->compile();
237  }
238 
239 
240  /**
241  * Measurement of light scattering in water.
242  *
243  * \param x cosine scattering angle
244  * \return probability
245  */
246  double operator()(const double x) const
247  {
248  double a;
249 
250  if (x >= +1.0)
251  a = 0.0;
252  else if (x <= -1.0)
253  a = 180.0;
254  else
255  a = acos(x) * 180.0 / PI;
256 
257  if (a >= this->getXmax())
258  return this->rbegin()->getY();
259  else if (a <= this->getXmin())
260  return this-> begin()->getY();
261  else
263  }
264  };
265 
266 
267  /**
268  * Function object for measurement of light scattering in water.
269  */
270  static const JPetzhold petzhold;
271 
272 
273  /**
274  * Auxiliary data structure for scattering lengths of deep-sea water.
275  *
276  * Use the Kopelevich model for spectral volume scattering functions.\n
277  * The model separates the contributions by:
278  * - pure: pure sea water;
279  * - small: 'small' particles (size < 1 um);
280  * - large: 'large' particles (size > 1 um).
281  *
282  * Values are taken from reference C.D. Mobley "Light and Water", ISBN 0-12-502750-8, pag. 119.
283  */
284  struct JMobley {
285  /**
286  * Constructor.
287  *
288  * \param lambda wavelength of light [nm]
289  */
290  JMobley(const double lambda)
291  {
292  static const double Vs = 0.0075;
293  static const double Vl = 0.0075;
294  static const double bw = 0.0017;
295  static const double bs = 1.340;
296  static const double bl = 0.312;
297 
298  const double x = 550.0/lambda;
299 
300  Lis_pure = bw * pow(x, 4.3);
301  Lis_small = bs * Vs * pow(x, 1.7);
302  Lis_large = bl * Vl * pow(x, 0.3);
303  }
304 
305  double Lis_pure; //!< inverse scattering length due to pure water [m^-1]
306  double Lis_small; //!< inverse scattering length due to small particles [m^-1]
307  double Lis_large; //!< inverse scattering length due to large particles [m^-1]
308  };
309 
310 
311  /**
312  * Absorption length of pure water.
313  *
314  * CITATION:
315  * Jonasz M. 2007. Absorption coefficient of water: Data sources (www.tpdsci.com/Tpc/AbsCfOfWaterDat.php).
316  * In: Top. Part. Disp. Sci. (www.tpdsci.com).
317  *
318  * DATA FROM:
319  * Wozniak B., Wozniak S. B., Tyszka K., Ostrowska M., Majchrowski R., Ficek D., Dera J. 2005.
320  * Modelling the light absorption properties of particulate matter forming organic particles suspended in seawater. Part 2.
321  * Modelling results. Oceanologia 47, 621-662.
322  * see also
323  * Wozniak B., Dera J. 2007.
324  * Light absorption in sea water. Springer, Berlin, 456 pp. (see p. 62)
325  *
326  * NOTES:
327  * As stated by the data authors, the data are based on measurement results obtained by various authors
328  * (interpolated by a linear approximation where applicable):
329  * Wavelength Reference
330  * - 200-335 nm Smith R.C., Baker K.S. 1981. Optical properties of the clearest natural waters (200-800 nm). Appl. Opt. 20, 177-184.
331  * - 340-370 nm Sogandares F.M., Fry, E.S. 1997. Absorption spectrum (340 -640 nm) of pure water. I. Photothermal measurements Appl. Opt. 36, 8699-8709.
332  * - 380-700 nm Pope R.M., Fry E.S. 1997. Absorption spectrum (380 -700 nm) of pure water. II. Integrating cavity measurements. Appl. Opt. 36, 8710-8723
333  */
336  {
337  public:
338  /**
339  * Default constructor.
340  */
342  {
343  // wave- absorption
344  // length coefficient
345  // [um] [1/m]
346  (*this)[0.200e3] = 3.07;
347  (*this)[0.205e3] = 2.53;
348  (*this)[0.210e3] = 1.99;
349  (*this)[0.215e3] = 1.65;
350  (*this)[0.220e3] = 1.31;
351  (*this)[0.225e3] = 1.1185;
352  (*this)[0.230e3] = 0.927;
353  (*this)[0.235e3] = 0.8235;
354  (*this)[0.240e3] = 0.72;
355  (*this)[0.245e3] = 0.6395;
356  (*this)[0.250e3] = 0.559;
357  (*this)[0.255e3] = 0.508;
358  (*this)[0.260e3] = 0.457;
359  (*this)[0.265e3] = 0.415;
360  (*this)[0.270e3] = 0.373;
361  (*this)[0.275e3] = 0.3305;
362  (*this)[0.280e3] = 0.288;
363  (*this)[0.285e3] = 0.2515;
364  (*this)[0.290e3] = 0.215;
365  (*this)[0.295e3] = 0.178;
366  (*this)[0.300e3] = 0.141;
367  (*this)[0.305e3] = 0.123;
368  (*this)[0.310e3] = 0.105;
369  (*this)[0.315e3] = 0.0947;
370  (*this)[0.320e3] = 0.0844;
371  (*this)[0.325e3] = 0.0761;
372  (*this)[0.330e3] = 0.0678;
373  (*this)[0.335e3] = 0.06195;
374  (*this)[0.340e3] = 0.0325;
375  (*this)[0.345e3] = 0.02645;
376  (*this)[0.350e3] = 0.0204;
377  (*this)[0.355e3] = 0.018;
378  (*this)[0.360e3] = 0.0156;
379  (*this)[0.365e3] = 0.0135;
380  (*this)[0.370e3] = 0.0114;
381  (*this)[0.375e3] = 0.011385;
382  (*this)[0.380e3] = 0.01137;
383  (*this)[0.385e3] = 0.00941;
384  (*this)[0.390e3] = 0.00851;
385  (*this)[0.395e3] = 0.00813;
386  (*this)[0.400e3] = 0.00663;
387  (*this)[0.405e3] = 0.0053;
388  (*this)[0.410e3] = 0.00473;
389  (*this)[0.415e3] = 0.00444;
390  (*this)[0.420e3] = 0.00454;
391  (*this)[0.425e3] = 0.00478;
392  (*this)[0.430e3] = 0.00495;
393  (*this)[0.435e3] = 0.0053;
394  (*this)[0.440e3] = 0.00635;
395  (*this)[0.445e3] = 0.00751;
396  (*this)[0.450e3] = 0.00922;
397  (*this)[0.455e3] = 0.00962;
398  (*this)[0.460e3] = 0.00979;
399  (*this)[0.465e3] = 0.01011;
400  (*this)[0.470e3] = 0.0106;
401  (*this)[0.475e3] = 0.0114;
402  (*this)[0.480e3] = 0.0127;
403  (*this)[0.485e3] = 0.0136;
404  (*this)[0.490e3] = 0.015;
405  (*this)[0.495e3] = 0.0173;
406  (*this)[0.500e3] = 0.0204;
407  (*this)[0.505e3] = 0.0256;
408  (*this)[0.510e3] = 0.0325;
409  (*this)[0.515e3] = 0.0396;
410  (*this)[0.520e3] = 0.0409;
411  (*this)[0.525e3] = 0.0417;
412  (*this)[0.530e3] = 0.0434;
413  (*this)[0.535e3] = 0.0452;
414  (*this)[0.540e3] = 0.0474;
415  (*this)[0.545e3] = 0.0511;
416  (*this)[0.550e3] = 0.0565;
417  (*this)[0.555e3] = 0.0596;
418  (*this)[0.560e3] = 0.0619;
419  (*this)[0.565e3] = 0.0642;
420  (*this)[0.570e3] = 0.0695;
421  (*this)[0.575e3] = 0.0772;
422  (*this)[0.580e3] = 0.0896;
423  (*this)[0.585e3] = 0.11;
424  (*this)[0.590e3] = 0.1351;
425  (*this)[0.595e3] = 0.1672;
426  (*this)[0.600e3] = 0.2224;
427  (*this)[0.605e3] = 0.2577;
428  (*this)[0.610e3] = 0.2644;
429  (*this)[0.615e3] = 0.2678;
430  (*this)[0.620e3] = 0.2755;
431  (*this)[0.625e3] = 0.2834;
432  (*this)[0.630e3] = 0.2916;
433  (*this)[0.635e3] = 0.3012;
434  (*this)[0.640e3] = 0.3108;
435  (*this)[0.645e3] = 0.325;
436  (*this)[0.650e3] = 0.34;
437  (*this)[0.655e3] = 0.371;
438  (*this)[0.660e3] = 0.41;
439  (*this)[0.665e3] = 0.429;
440  (*this)[0.670e3] = 0.439;
441  (*this)[0.675e3] = 0.448;
442  (*this)[0.680e3] = 0.465;
443  (*this)[0.685e3] = 0.486;
444  (*this)[0.690e3] = 0.516;
445  (*this)[0.695e3] = 0.559;
446  (*this)[0.700e3] = 0.624;
447 
448  compile();
449  }
450 
451 
452  /**
453  * Absorption length of pure water.
454  *
455  * \param lambda wavelength of light [nm]
456  * \return absorption length [m]
457  */
458  double operator()(const double lambda) const
459  {
460  const double y = JGridSplineFunction1D_t::operator()(lambda);
461 
462  return 1.0 / y;
463  }
464  };
465 
466 
467  /**
468  * Function object for absorption length of pure water.
469  */
471 }
472 
473 #endif
double f4(const double x)
Model specific function to describe light scattering in water (f4).
Q(UTCMax_s-UTCMin_s)-livetime_s
then usage $script< input file >[option[primary[working directory]]] nWhere option can be E
Definition: JMuonPostfit.sh:40
result_type operator()(const argument_type x) const
Function value evaluation.
Definition: JFunctional.hh:348
JMobley(const double lambda)
Constructor.
then set_variable DIR else fatal Wrong number of arguments fi for INPUT_FILE in ls rt $DIR stage * log
static const JAbsorptionLengthOfPureWater getAbsorptionLengthOfPureWater
Function object for absorption length of pure water.
double operator()(const double lambda) const
Absorption length of pure water.
double operator()(const double x) const
Measurement of light scattering in water.
double Lis_large
inverse scattering length due to large particles [m^-1]
static const double MASS_MUON
muon mass [GeV]
o $QUALITY_ROOT d $DEBUG!CHECK_EXIT_CODE JPlot1D f
Definition: JDataQuality.sh:76
double rayleigh(const double a, const double x)
Auxiliary method to describe light scattering in water (Rayleigh).
Type definition of a spline interpolation based on a JGridCollection with result type double...
double Lis_small
inverse scattering length due to small particles [m^-1]
double getThetaMCS(const double E, const double x, const double X0, const double M, const double Q)
Get multiple Coulomb scattering angle.
Absorption length of pure water.
static const JPetzhold petzhold
Function object for measurement of light scattering in water.
then JCalibrateToT a
Definition: JTuneHV.sh:113
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:97
Physics constants.
static const double PI
Mathematical constants.
static const double THETA_MCS
Multiple Coulomb scattering constant [GeV].
JAbsorptionLengthOfPureWater()
Default constructor.
Measurement of light scattering in water.
Type definition of a 1st degree polynomial interpolation with result type double. ...
Auxiliary data structure for scattering lengths of deep-sea water.
double Lis_pure
inverse scattering length due to pure water [m^-1]
double p00075(const double x)
Model specific function to describe light scattering in water (p00075).
double henyey_greenstein(const double g, const double x)
Auxiliary method to describe light scattering in water (Henyey-Greenstein).
static const double X0_WATER_M
Radiation length pure water [m].
JPetzhold()
Default constructor.
Double_t g1(const Double_t x)
Function.
Definition: JQuantiles.cc:25