Jpp  master_rocky-43-ge265d140c
the software that should make you happy
JDiffPD0.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <fstream>
4 #include <iomanip>
5 #include <cmath>
6 #include <algorithm>
7 
10 #include "JPhysics/JPDFTable.hh"
11 
12 #include "Jeep/JParser.hh"
13 #include "Jeep/JMessage.hh"
14 
15 namespace {
16  /**
17  * Compare two function values.
18  * If average function value less than one, absolute comparison; else relative comparison.
19  *
20  * \param u first value
21  * \param v second value
22  * \param precision precision
23  * \return true if comparable; else false
24  */
25  inline bool compare(const double u, const double v, const double precision)
26  {
27  using namespace std;
28 
29  return fabs(u - v) <= precision * max(1.0, 0.5*(u+v));
30  }
31 }
32 
33 /**
34  * \file
35  *
36  * Auxiliary program to compare PDF tables of the arrival time of the Cherenkov light from a bright point.
37  * \author mdejong
38  */
39 int main(int argc, char **argv)
40 {
41  using namespace std;
42  using namespace JPP;
43 
44  string inputFileA;
45  string inputFileB;
46  double precision;
47  int debug;
48 
49  try {
50 
51  JParser<> zap("Auxiliary program to compare PDF tables of the arrival time of the Cherenkov light from a bright point.");
52 
53  zap['a'] = make_field(inputFileA);
54  zap['b'] = make_field(inputFileB);
55  zap['p'] = make_field(precision) = 1.0e-6;
56  zap['d'] = make_field(debug) = 0;
57 
58  zap(argc, argv);
59  }
60  catch(const exception &error) {
61  FATAL(error.what() << endl);
62  }
63 
64  typedef JSplineFunction1D_t JFunction1D_t;
66  JPolint1FunctionalGridMap>::maplist JMapList_t;
68 
69  typedef JFunction1D_t::argument_type argument_type;
70  typedef JArray<2, argument_type> JArray_t;
71 
72  JPDF_t pdfA;
73  JPDF_t pdfB;
74 
75  const JFunction1D_t::JSupervisor supervisor(new JFunction1D_t::JDefaultResult(0.0));
76 
77  struct pair_type {
78  const string file_name;
79  JPDF_t& pdf;
80  } buffer[] = {
81  { inputFileA, pdfA },
82  { inputFileB, pdfB }
83  };
84 
85  for (const pair_type& i : buffer) {
86 
87  try {
88 
89  NOTICE("loading input from file " << i.file_name << "... " << flush);
90 
91  i.pdf.load(i.file_name.c_str());
92  i.pdf.compile();
93  i.pdf.setExceptionHandler(supervisor);
94 
95  NOTICE("OK" << endl);
96  }
97  catch(const JException& error) {
98  FATAL(error.what() << endl);
99  }
100  }
101 
102  bool ok = true;
103 
104  JPDF_t::super_const_iterator i = pdfA.super_begin();
105  JPDF_t::super_const_iterator j = pdfB.super_begin();
106 
107  for ( ; i != pdfA.super_end() &&
108  j != pdfB.super_end(); ++i, ++j) {
109 
110  const double Wa = pdfA.transformer->getWeight(JArray_t((*i).getKey()));
111  const double Wb = pdfB.transformer->getWeight(JArray_t((*j).getKey()));
112 
113  if (fabs(i->first - j->first) > precision ||
114  fabs(i->second->first - j->second->first) > precision ||
115  fabs(Wa - Wb) > precision) {
116 
117  DEBUG("a[] "
118  << i->first << ' '
119  << i->second->first << ' '
120  << Wa << endl);
121 
122  DEBUG("b[] "
123  << j->first << ' '
124  << j->second->first << ' '
125  << Wb << endl);
126 
127  ok = false;
128  }
129 
130  const JFunction1D_t& fa = (*i).getValue();
131  const JFunction1D_t& fb = (*j).getValue();
132 
133  JFunction1D_t::const_iterator p = fa.begin();
134  JFunction1D_t::const_iterator q = fb.begin();
135 
136  while (p != fa.end() &&
137  q != fb.end()) {
138 
139  if (fabs(p->getX() - q->getX()) > precision) {
140 
141  DEBUG("a.x " << p->getX() << endl);
142  DEBUG("b.x " << q->getX() << endl);
143 
144  ok = false;
145  }
146 
147  if (!compare(p->getY(), q->getY(), precision)) {
148 
149  DEBUG("a.y " << p->getX() << ' ' << p->getY() << endl);
150  DEBUG("b.y " << q->getX() << ' ' << q->getY() << endl);
151 
152  ok = false;
153  }
154 
155  if (p->getX() < q->getX())
156  ++p;
157  else if (q->getX() < p->getX())
158  ++q;
159  else {
160  ++p;
161  ++q;
162  }
163  }
164 
165  for ( ; p != fa.end(); ++p) {
166 
167  DEBUG("a.x " << p->getX() << endl);
168 
169  ok = false;
170  }
171 
172  for ( ; q != fb.end(); ++q) {
173 
174  DEBUG("b.x " << q->getX() << endl);
175 
176  ok = false;
177  }
178  }
179 
180  for ( ; i != pdfA.super_end(); ++i) {
181 
182  DEBUG("a[] "
183  << i->first << ' '
184  << i->second->first << endl);
185 
186  ok = false;
187  }
188 
189  for ( ; j != pdfB.super_end(); ++j) {
190 
191  DEBUG("b[] "
192  << j->first << ' '
193  << j->second->first << endl);
194 
195  ok = false;
196  }
197 
198  if (!ok) {
199  ERROR(inputFileA << " and " << inputFileB << " differ." << endl);
200  }
201 
202  return (ok ? 0 : 1);
203 }
int main(int argc, char **argv)
Definition: JDiffPD0.cc:39
Various implementations of functional maps.
General purpose messaging.
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
#define ERROR(A)
Definition: JMessage.hh:66
#define NOTICE(A)
Definition: JMessage.hh:64
#define FATAL(A)
Definition: JMessage.hh:67
int debug
debug level
Definition: JSirene.cc:69
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:2142
General exception.
Definition: JException.hh:24
virtual const char * what() const override
Get error message.
Definition: JException.hh:64
Utility class to parse command line options.
Definition: JParser.hh:1698
Multi-dimensional PDF table for arrival time of Cherenkov light.
Definition: JPDFTable.hh:44
One dimensional array of template objects with fixed length.
Definition: JArray.hh:43
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
double u[N+1]
Definition: JPolint.hh:865
data_type v[N+1][M+1]
Definition: JPolint.hh:866
int j
Definition: JPolint.hh:792
Definition: JSTDTypes.hh:14
Data structure for a pair of indices.
Auxiliary class for recursive map list generation.
Definition: JMapList.hh:109
Type definition of a 1st degree polynomial interpolation based on a JGridMap implementation.
Type definition of a 1st degree polynomial interpolation based on a JMap implementation.
Type definition of a spline interpolation method based on a JCollection with double result type.