Jpp  17.3.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JLED.hh
Go to the documentation of this file.
1 #ifndef __JPHYSICS__JLED__
2 #define __JPHYSICS__JLED__
3 
8 #include "JPhysics/JConstants.hh"
9 #include "JTools/JElement.hh"
10 #include "JTools/JQuadrature.hh"
11 
12 
13 /**
14  * \author mdejong
15  */
16 
17 namespace JPHYSICS {}
18 namespace JPP { using namespace JPHYSICS; }
19 
20 namespace JPHYSICS {
21 
22 
23  using JTOOLS::JQuadrature;
25  using JTOOLS::JCotangent;
26 
29 
30 
31  /**
32  * Interface for emission profile from LED.
33  */
34  class JAbstractLED {
35  public:
36 
37  /**
38  * Virtual destructor.
39  */
40  virtual ~JAbstractLED()
41  {}
42 
43 
44  /**
45  * Light yield from LED (number of p.e. per unit solid angle per unit time).
46  *
47  * \param ct zenith angle of emission
48  * \param phi azimuth angle of emission
49  * \param dt time of emission [ns]
50  * \return d^2P / dOmega dt [npe/ns/sr]
51  */
52  virtual double getLightFromLED(const double ct,
53  const double phi,
54  const double dt) const = 0;
55  };
56 
57 
58  /**
59  * Probability Density Functions of the time response of a PMT.
60  */
61  class JLED :
62  public virtual JDispersionInterface,
63  public JAbstractPMT,
64  public JAbstractLED,
65  public JAbstractMedium
66  {
67  public:
68  /**
69  * Constructor.
70  *
71  * \param lambda wavelength photon [nm]
72  * \param Tmin_ns minimal time of emmision [ns]
73  * \param Tmax_ns minimal time of emmision [ns]
74  * \param engine scattering angle integrator
75  * \param numberOfPoints number of points for integration
76  * \param epsilon precision of points for integration
77  */
78  JLED(const double lambda,
79  const double Tmin_ns,
80  const double Tmax_ns,
81  const JQuadrature& engine = JCotangent(20),
82  const int numberOfPoints = 20,
83  const double epsilon = 1e-12) :
84  wavelength(lambda),
85  tmin(Tmin_ns),
86  tmax(Tmax_ns),
88  beta_engine(engine)
89  {
90  const double dp = PI / numberOfPoints;
91 
92  for (double x = 0.5*dp; x < PI; x += dp) {
93  phi_engine.push_back(JElement3D_t(cos(x),sin(x), dp));
94  }
95  }
96 
97 
98  /**
99  * Probability density function for direct light from LED.
100  *
101  * \param D_m distance between LED and PMT [m]
102  * \param cd cosine angle LED orientation and LED - PMT position
103  * \param theta zenith angle orientation PMT
104  * \param phi azimuth angle orientation PMT
105  * \param t_ns time difference relative to direct light [ns]
106  * \return dP/dt [npe/ns]
107  */
108  double getDirectLightFromLED(const double D_m,
109  const double cd,
110  const double theta,
111  const double phi,
112  const double t_ns) const
113  {
114  const double sd = sqrt((1.0 + cd)*(1.0 - cd));
115  const double d = D_m; // photon path [m]
116 
117  const double A = getPhotocathodeArea();
118 
119  const double px = sin(theta)*cos(phi);
120  //const double py = sin(theta)*sin(phi);
121  const double pz = cos(theta);
122 
123  const double ct = sd*px + cd*pz; // cosine angle of incidence on PMT
124 
125  const double qe = getQE(wavelength);
126  const double l_abs = getAbsorptionLength(wavelength);
127  const double ls = getScatteringLength(wavelength);
128 
129  const double U = getAngularAcceptance(ct);
130  const double V = exp(-d/l_abs) * exp(-d/ls); // absorption & scattering
131  const double W = A/(d*d); // solid angle
132 
133  return qe * getLightFromLED(cd, 0.0, t_ns) * U * V * W;
134  }
135 
136 
137  /**
138  * Probability density function for scattered light from LED.
139  *
140  * \param D_m distance between LED and PMT [m]
141  * \param cd cosine angle LED orientation and LED - PMT position
142  * \param theta zenith angle orientation PMT
143  * \param phi azimuth angle orientation PMT
144  * \param t_ns time difference relative to direct light [ns]
145  * \return dP/dt [npe/ns]
146  */
147  double getScatteredLightFromLED(const double D_m,
148  const double cd,
149  const double theta,
150  const double phi,
151  const double t_ns) const
152  {
153  using namespace std;
154 
155  double value = 0;
156 
157  if (t_ns >= tmin) {
158 
159  const double ng = getIndexOfRefractionGroup(wavelength);
160 
161  const double sd = sqrt((1.0 + cd)*(1.0 - cd));
162  const double l = D_m; // distance [m]
163 
164  const double A = getPhotocathodeArea();
165  //const double sr = sqrt(A/PI) / D_m;
166  //const double cr = sqrt((1.0 + sr)*(1.0 - sr));
167 
168  const double px = sin(theta)*cos(phi);
169  const double py = sin(theta)*sin(phi);
170  const double pz = cos(theta);
171 
172  const double qx = cd*px + 0 - sd*pz;
173  const double qy = 1*py;
174  const double qz = sd*px + 0 + cd*pz;
175 
176  const double qe = getQE(wavelength);
177  const double l_abs = getAbsorptionLength(wavelength);
178  const double ls = getScatteringLength(wavelength);
179 
180  const double Jc = 1.0 / ls; // dN/dx
181 
182  const double xmin = tmin;
183  const double xmax = std::min(tmax, t_ns);
184 
185  JQuadrature buffer;
186 
187  if (xmax > xmin) {
188 
189  for (JQuadrature::const_iterator i = main_engine.begin(); i != main_engine.end(); ++i) {
190 
191  const double t = 0.5 * (xmax + xmin) + i->getX() * 0.5 * (xmax - xmin);
192  const double dt = i->getY() * 0.5 * (xmax - xmin);
193 
194  buffer[t] = dt;
195  }
196 
197  } else {
198 
199  buffer[tmin] = 1.0;
200  }
201 
202  for (JQuadrature::const_iterator i = buffer.begin(); i != buffer.end(); ++i) {
203 
204  const double t = i->getX();
205  const double dt = i->getY();
206 
207  const double d = l + C*(t_ns - t)/ng; // photon path
208  const double V = exp(-d/l_abs); // absorption
209 
210  for (JQuadrature::const_iterator j = beta_engine.begin(); j != beta_engine.end(); ++j) {
211 
212  const double cb = j->getX();
213  const double dc = j->getY();
214 
215  const double sb = sqrt((1.0 + cb)*(1.0 - cb));
216 
217  const double v = 0.5 * (d + l) * (d - l) / (d - l*cb);
218  const double u = d - v;
219 
220  if (u <= 0) continue;
221  if (v <= 0) continue;
222 
223  const double cts = (l*cb - v) / u; // cosine scattering angle
224 
225  const double W = min(A/(v*v), 2.0*PI); // solid angle
226  const double Ja = getScatteringProbability(cts); // P(cos(),phi)
227  const double Jd = ng * (1.0 - cts) / C; // dt/du
228 
229  const double ca = (l - v*cb) / u;
230  const double sa = v*sb / u;
231 
232  for (std::vector<JElement3D_t>::const_iterator k = phi_engine.begin(); k != phi_engine.end(); ++k) {
233 
234  const double cp = k->getX();
235  const double sp = k->getY();
236  const double dp = k->getZ();
237 
238  const double dom = dc * dp * v*v / (u*u);
239 
240  const double ct0 = cd*ca - sd*sa*cp;
241  const double ph0 = atan2(sa*sp, cd*sa*cp + sd*ca);
242 
243  const double vx = -sb*cp * qx;
244  const double vy = -sb*sp * qy;
245  const double vz = cb * qz;
246 
247  const double U =
248  getAngularAcceptance(vx + vy + vz) * getLightFromLED(ct0, +ph0, t) +
249  getAngularAcceptance(vx - vy + vz) * getLightFromLED(ct0, -ph0, t);
250 
251  value += qe * dt * dom * Ja * Jc * U * V * W / Jd;
252  }
253  }
254  }
255  }
256 
257  return value;
258  }
259 
260  protected:
261  double wavelength; // wavelength photon [nm]
262 
263  double tmin; // minimal time of emmision [ns]
264  double tmax; // maximal time of emmision [ns]
265 
266  JQuadrature main_engine; // numerical integration
267  JQuadrature beta_engine; // numerical integration of scattering angle
268  std::vector<JElement3D_t> phi_engine; // numerical integration of phi angle
269  };
270 
271 
272  /**
273  * Probability Density Functions of the time response of a PMT (C-like interface)
274  */
275  class JLED_C :
276  public JLED,
277  public JDispersion
278  {
279  public:
280  /**
281  * Constructor.
282  *
283  * \param Area photo-cathode area [m^2]
284  * \param LED pointer to interface for emission profile from LED
285  * \param QE pointer to function for quantum efficiency of PMT
286  * \param Pmt pointer to function for angular acceptance of PMT
287  * \param L_abs absorption length [m]
288  * \param L_s scattering length [m]
289  * \param Km3 pointer to model specific function to describe light scattering in water
290  * \param P_atm ambient pressure [atm]
291  * \param wavelength wavelength of light [nm]
292  * \param Tmin_ns minimal time of emmision [ns]
293  * \param Tmax_ns minimal time of emmision [ns]
294  * \param engine scattering angle integrator
295  * \param numberOfPoints number of points for integration
296  * \param epsilon precision of points for integration
297  */
298  JLED_C(const double Area,
299  //
300  // pointers to static functions
301  //
302  const JAbstractLED* LED,
303  double (*QE) (const double),
304  double (*Pmt) (const double),
305  const double L_abs,
306  const double L_s,
307  double (*Km3) (const double),
308  //
309  // parameters for base class
310  //
311  const double P_atm,
312  const double wavelength,
313  const double Tmin_ns,
314  const double Tmax_ns,
315  const JQuadrature& engine = JCotangent(20),
316  const int numberOfPoints = 20,
317  const double epsilon = 1e-12) :
318  JLED(wavelength, Tmin_ns, Tmax_ns, engine, numberOfPoints, epsilon),
319  JDispersion(P_atm),
320  A (Area ),
321  led (LED ),
322  qe (QE ),
323  l_abs(L_abs),
324  ls (L_s ),
325  pmt (Pmt ),
326  km3 (Km3 )
327  {}
328 
329 
330  /**
331  * Light yield from LED (number of p.e. per unit solid angle per unit time).
332  *
333  * \param ct zenith angle of emission
334  * \param phi azimuth angle of emission
335  * \param dt time of emission [ns]
336  * \return d^2P / dOmega dt [npe/ns/sr]
337  */
338  virtual double getLightFromLED(const double ct,
339  const double phi,
340  const double dt) const
341  {
342  return led->getLightFromLED(ct, phi, dt);
343  }
344 
345 
346  /**
347  * Photo-cathode area of PMT.
348  *
349  * \return photo-cathode area [m^2]
350  */
351  virtual double getPhotocathodeArea() const override
352  {
353  return A;
354  }
355 
356 
357  /**
358  * Quantum efficiency of PMT (incl. absorption in glass, gel, etc.).
359  *
360  * \param lambda wavelenth [nm]
361  * \return QE
362  */
363  virtual double getQE(const double lambda) const override
364  {
365  return qe(lambda);
366  }
367 
368 
369  /**
370  * Angular acceptance of PMT.
371  *
372  * \param ct cosine angle of incidence
373  * \return relative efficiency
374  */
375  virtual double getAngularAcceptance(const double ct) const override
376  {
377  return pmt(ct);
378  }
379 
380 
381  /**
382  * Absorption length.
383  *
384  * \param lambda wavelenth [nm]
385  * \return absorption length [m]
386  */
387  virtual double getAbsorptionLength(const double lambda) const override
388  {
389  return l_abs;
390  }
391 
392 
393  /**
394  * Scattering length.
395  *
396  * \param lambda wavelenth [nm]
397  * \return scattering length [m]
398  */
399  virtual double getScatteringLength(const double lambda) const override
400  {
401  return ls;
402  }
403 
404 
405  /**
406  * Model specific function to describe light scattering in water
407  * (integral over full solid angle normalised to one).
408  *
409  * \param ct cosine scattering angle
410  * \return probability
411  */
412  virtual double getScatteringProbability(const double ct) const override
413  {
414  return km3(ct);
415  }
416 
417  protected:
418  /**
419  * photo-cathode area [m2]
420  */
421  const double A;
422 
423  /**
424  * Pointer to interface for emission profile from LED
425  */
427 
428  /**
429  * Quantum efficiency of PMT (incl. absorption in glass, gel, etc.)
430  *
431  * \param lambda wavelenth [nm]
432  * \return QE
433  */
434  double (*qe)(const double lambda);
435 
436  /**
437  * Absorption length
438  */
439  double l_abs;
440 
441  /**
442  * Scattering length
443  */
444  double ls;
445 
446  /**
447  * Angular acceptance of PMT
448  *
449  * \param ct cosine angle of incidence
450  * \return relative efficiency
451  */
452  double (*pmt)(const double ct);
453 
454  /**
455  * Model specific function to describe light scattering in water
456  *
457  * \param ct cosine scattering angle
458  * \return probability
459  */
460  double (*km3)(const double ct);
461  };
462 }
463 
464 #endif
const double xmax
Definition: JQuadrature.cc:24
JTOOLS::JElement3D< double, double > JElement3D_t
Definition: JLED.hh:28
double tmin
Definition: JLED.hh:263
then fatal No hydrophone data file $HYDROPHONE_TXT fi sort gr k
virtual double getAbsorptionLength(const double lambda) const override
Absorption length.
Definition: JLED.hh:387
Numerical integrator for .
Definition: JQuadrature.hh:111
virtual ~JAbstractLED()
Virtual destructor.
Definition: JLED.hh:40
The elements in a collection are sorted according to their abscissa values and a given distance opera...
virtual double getQE(const double lambda) const override
Quantum efficiency of PMT (incl.
Definition: JLED.hh:363
const JAbstractLED * led
Pointer to interface for emission profile from LED.
Definition: JLED.hh:426
Implementation of dispersion for water in deep sea.
Definition: JDispersion.hh:26
JQuadrature main_engine
Definition: JLED.hh:266
Numerical integrator for .
Definition: JQuadrature.hh:482
double tmax
Definition: JLED.hh:264
3D Element.
Definition: JElement.hh:494
Interface for emission profile from LED.
Definition: JLED.hh:34
V(JDAQEvent-JTriggerReprocessor)*1.0/(JDAQEvent+1.0e-10)
JLED_C(const double Area, const JAbstractLED *LED, double(*QE)(const double), double(*Pmt)(const double), const double L_abs, const double L_s, double(*Km3)(const double), const double P_atm, const double wavelength, const double Tmin_ns, const double Tmax_ns, const JQuadrature &engine=JCotangent(20), const int numberOfPoints=20, const double epsilon=1e-12)
Constructor.
Definition: JLED.hh:298
Light yield from LED (number of p.e.
Definition: JDrawLED.cc:31
PMT interface.
Definition: JAbstractPMT.hh:17
then usage set_variable ACOUSTICS_WORKDIR $WORKDIR set_variable FORMULA sin([0]+2 *$PI *([1]+[2]*x)*x)" set_variable DY 1.0e-8 mkdir $WORKDIR for DETECTOR in $DETECTORS[*]
static const double C
Physics constants.
virtual double getLightFromLED(const double ct, const double phi, const double dt) const
Light yield from LED (number of p.e.
Definition: JLED.hh:338
virtual double getScatteringLength(const double lambda) const override
Scattering length.
Definition: JLED.hh:399
std::vector< JElement3D_t > phi_engine
Definition: JLED.hh:268
virtual double getScatteringLength(const double lambda) const =0
Scattering length.
double l_abs
Absorption length.
Definition: JLED.hh:439
Probability Density Functions of the time response of a PMT (C-like interface)
Definition: JLED.hh:275
JQuadrature beta_engine
Definition: JLED.hh:267
virtual double getAbsorptionLength(const double lambda) const =0
Absorption length.
Physics constants.
static const double PI
Mathematical constants.
virtual double getPhotocathodeArea() const =0
Photo-cathode area of PMT.
Type definition for numerical integration.
Definition: JQuadrature.hh:37
double ls
Scattering length.
Definition: JLED.hh:444
const double A
photo-cathode area [m2]
Definition: JLED.hh:421
const double xmin
Definition: JQuadrature.cc:23
virtual double getScatteringProbability(const double ct) const override
Model specific function to describe light scattering in water (integral over full solid angle normali...
Definition: JLED.hh:412
Auxiliary classes for numerical integration.
virtual double getIndexOfRefractionGroup(const double lambda) const
Index of refraction for group velocity.
then JMuonMCEvt f $INPUT_FILE o $INTERMEDIATE_FILE d
Definition: JMuonPath.sh:47
virtual abscissa_type getX(int index) const override
Get abscissa value.
Definition: JCollection.hh:209
double getScatteredLightFromLED(const double D_m, const double cd, const double theta, const double phi, const double t_ns) const
Probability density function for scattered light from LED.
Definition: JLED.hh:147
int numberOfPoints
Definition: JResultPDF.cc:22
Probability Density Functions of the time response of a PMT.
Definition: JLED.hh:61
virtual double getAngularAcceptance(const double ct) const =0
Angular acceptence of PMT.
virtual double getAngularAcceptance(const double ct) const override
Angular acceptance of PMT.
Definition: JLED.hh:375
virtual double getLightFromLED(const double ct, const double phi, const double dt) const =0
Light yield from LED (number of p.e.
2D Element.
Definition: JElement.hh:46
then cp
JTOOLS::JElement2D< double, double > JElement2D_t
Definition: JLED.hh:27
container_type::const_iterator const_iterator
Definition: JCollection.hh:91
virtual double getPhotocathodeArea() const override
Photo-cathode area of PMT.
Definition: JLED.hh:351
int j
Definition: JPolint.hh:703
virtual double getQE(const double lambda) const =0
Quantum efficiency of PMT (incl.
data_type v[N+1][M+1]
Definition: JPolint.hh:777
double u[N+1]
Definition: JPolint.hh:776
Light dispersion inteface.
double getDirectLightFromLED(const double D_m, const double cd, const double theta, const double phi, const double t_ns) const
Probability density function for direct light from LED.
Definition: JLED.hh:108
source $JPP_DIR setenv csh $JPP_DIR &dev null eval JShellParser o a A
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
const double epsilon
Definition: JQuadrature.cc:21
virtual double getScatteringProbability(const double ct) const =0
Model specific function to describe light scattering in water (integral over full solid angle normali...
JLED(const double lambda, const double Tmin_ns, const double Tmax_ns, const JQuadrature &engine=JCotangent(20), const int numberOfPoints=20, const double epsilon=1e-12)
Constructor.
Definition: JLED.hh:78
double wavelength
Definition: JLED.hh:261