39int main(
int argc,
char **argv)
58 JParser<> zap(
"Program to plot PDF of Cherenkov light from EM-shower using interpolation tables.");
63 zap[
'D'] =
make_field(dir,
"(theta, phi) of PMT [rad]");
66 zap[
'E'] =
make_field(E,
"Energy [GeV]") = 1.0;
72 catch(
const exception &error) {
73 FATAL(error.what() << endl);
77 const double theta = dir.first;
78 const double phi = dir.second;
98 NOTICE(
"loading input from file " << pdfFile <<
"... " << flush);
100 pdf.load(pdfFile.c_str());
102 pdf.setExceptionHandler(
new JFunction1D_t::JDefaultResult(
JMATH::zero));
108 NOTICE(
"loading input from file " << cdfFile <<
"... " << flush);
110 cdf.load(cdfFile.c_str());
115 FATAL(error.what() << endl);
118 if (!x.is_valid()) { x =
JHistogram_t(250, 0.0, 250.0); }
119 if (!y.is_valid()) { y =
JHistogram_t(200, -1.0, +1.0); }
121 TH2D h0(
"pdf", NULL, x.getNumberOfBins(), x.getLowerLimit(), x.getUpperLimit(), y.getNumberOfBins(), y.getLowerLimit(), y.getUpperLimit());
122 TH2D h1(
"npe", NULL, x.getNumberOfBins(), x.getLowerLimit(), x.getUpperLimit(), y.getNumberOfBins(), y.getLowerLimit(), y.getUpperLimit());
123 TH2D h2(
"cdf", NULL, x.getNumberOfBins(), x.getLowerLimit(), x.getUpperLimit(), y.getNumberOfBins(), y.getLowerLimit(), y.getUpperLimit());
125 for (
int i = 1; i <= h0.GetNbinsX(); ++i) {
126 for (
int j = 1; j <= h0.GetNbinsY(); ++j) {
128 const double D = h0.GetXaxis()->GetBinCenter(i);
129 const double cd = h0.GetYaxis()->GetBinCenter(j);
133 const double y0 = get_integral(pdf(D, cd, theta, phi, 1.0e3));
134 const double y1 = npe (D, cd, theta, phi);
135 const double y2 = cdf.getNPE (D, cd, theta, phi);
138 <<
FIXED(5,1) << D <<
' '
139 <<
FIXED(5,1) << cd <<
' '
144 h0.SetBinContent(i, j, E * y0);
145 h1.SetBinContent(i, j, E * y1);
146 h2.SetBinContent(i, j, E * y2);
148 catch(
const exception& error) {
158 out << h0 << h1 << h2;