31 class JEMShowerCorrection
41 const double P_Atm = 240.0;
47 const double xmin = 1.0 / wmax;
48 const double xmax = 1.0 / wmin;
53 for (
double x = xmin, dx = (xmax - xmin) / (
double) nx;
x <=
xmax;
x += dx) {
55 const double w = 1.0 /
x;
56 const double dw = dx * w*w;
58 const double n = dispersion.getIndexOfRefractionPhase(w);
66 f1[MASS_ELECTRON] = 0.00;
67 f1[ 1.0e-3] = 90.96 / 1.0e-3;
68 f1[ 2.0e-3] = 277.36 / 2.0e-3;
69 f1[ 3.0e-3] = 485.82 / 3.0e-3;
70 f1[ 4.0e-3] = 692.83 / 4.0e-3;
71 f1[ 5.0e-3] = 890.01 / 5.0e-3;
72 f1[ 6.0e-3] = 1098.53 / 6.0e-3;
73 f1[ 7.0e-3] = 1285.47 / 7.0e-3;
74 f1[ 8.0e-3] = 1502.86 / 8.0e-3;
75 f1[ 9.0e-3] = 1687.15 / 9.0e-3;
76 f1[10.0e-3] = 1891.00 / 10.0e-3;
81 f2[-2.0] = log(1.891e3 / 1.0e-2);
82 f2[-1.0] = log(1.905e4 / 1.0e-1);
83 f2[ 0.0] = log(1.889e5 / 1.0e+0);
84 f2[+1.0] = log(1.875e6 / 1.0e+1);
85 f2[+2.0] = log(1.881e7 / 1.0e+2);
98 double operator()(
const double E)
const
100 if (E <=
f1.getXmin()) {
104 }
else if (E <=
f1.getXmax()) {
110 const double x = log10(E);
112 if (x <= f2.getXmin()) {
114 return exp(f2.begin()->getY());
116 }
else if (x <= f2.getXmax()) {
122 return exp(f2.rbegin()->getY());
137 static const JEMShowerCorrection getEMShowerCorrection;
158 JParser<> zap(
"Auxiliary program to draw npe as a function of EM-energy.");
167 catch(
const exception &error) {
168 FATAL(error.what() << endl);
175 pdf(NAMESPACE::getPhotocathodeArea(),
177 NAMESPACE::getAngularAcceptance,
178 NAMESPACE::getAbsorptionLength,
179 NAMESPACE::getScatteringLength,
180 NAMESPACE::getScatteringProbability,
181 NAMESPACE::getAmbientPressure(),
182 getMinimalWavelength(),
183 getMaximalWavelength(),
191 const double precision = 1.0e-10;
193 while (fabs(xmax - xmin) > precision) {
195 const double x = 0.5 * (xmin + xmax);
196 const double E = pow(10.0, x);
198 const double y = getEMShowerCorrection(E);
206 const double EMIN = pow(10.0, 0.5*(xmax + xmin));
208 NOTICE(
"Threshold kinetic energy [GeV] " << sqrt((EMIN + MASS_ELECTRON) * (EMIN - MASS_ELECTRON)) << endl);
213 TH1D h0(
"h0", NULL, 10000, -4.0, +4.0);
214 TH1D h1(
"h1", NULL, 10000, -4.0, +4.0);
216 for (
int i = 1; i <= h0.GetNbinsX(); ++i) {
218 const Double_t x = h0.GetBinCenter(i);
219 const Double_t E = pow(10.0, x);
222 h1.SetBinContent(i, getEMShowerCorrection(E));
int main(int argc, char **argv)
General purpose messaging.
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
I/O formatting auxiliaries.
Properties of KM3NeT PMT and deep-sea water.
Utility class to parse command line options.
Implementation of dispersion for water in deep sea.
Probability Density Functions of the time response of a PMT with an implementation of the JAbstractPM...
double getNumberOfPhotons() const
Number of Cherenkov photons per unit track length.
const JPolynome f1(1.0, 2.0, 3.0)
Function.
double cherenkov(const double lambda, const double n)
Number of Cherenkov photons per unit track length and per unit wavelength.
double getMinimalWavelength()
Get minimal wavelength for PDF evaluations.
double geanc()
Equivalent muon track length per unit shower energy.
double getMaximalWavelength()
Get maximal wavelength for PDF evaluations.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
JF1_t & div(const double factor)
Scale function.