Jpp  18.6.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JPDFTable.hh
Go to the documentation of this file.
1 #ifndef __JPHYSICS__JPDFTABLE__
2 #define __JPHYSICS__JPDFTABLE__
3 
4 #include <cmath>
5 
6 #include "JIO/JObjectBinaryIO.hh"
8 #include "JTools/JQuantiles.hh"
9 #include "JTools/JSet.hh"
10 #include "JTools/JRange.hh"
11 #include "JMath/JMathSupportkit.hh"
12 #include "JPhysics/JConstants.hh"
14 
15 
16 /**
17  * \author mdejong
18  */
19 
20 namespace JPHYSICS {}
21 namespace JPP { using namespace JPHYSICS; }
22 
23 namespace JPHYSICS {
24 
25  using JIO::JSerialisable;
26  using JIO::JReader;
27  using JIO::JWriter;
31  using JTOOLS::JRange;
32 
33 
34  /**
35  * Multi-dimensional PDF table for arrival time of Cherenkov light.
36  */
37  template<class JFunction1D_t,
38  class JMaplist_t,
40  class JPDFTable :
41  public JTransformableMultiFunction<JFunction1D_t, JMaplist_t, JDistance_t>,
42  public JSerialisable,
43  public JObjectBinaryIO< JPDFTable<JFunction1D_t, JMaplist_t, JDistance_t> >
44  {
45  public:
46 
48 
52 
55 
57 
61 
63 
64 
65  /**
66  * Default constructor.
67  */
70  {}
71 
72 
73  /**
74  * Constructor.
75  *
76  * \param input multi-dimensional function
77  */
78  template<class __JFunction_t, class __JMaplist_t, class __JDistance_t>
81  {}
82 
83 
84  /**
85  * Constructor.
86  *
87  * \param input multi-dimensional histogram
88  */
89  template<class JHistogram_t, class __JMaplist_t, class __JDistance_t>
92  {}
93 
94 
95  /**
96  * Blur PDF.
97  *
98  * The arrival times of Cherenkov light are smeared according to a Gaussian distribution
99  * with the specified width (i.e.\ TTS) using Gauss-Hermite integration.
100  * An exception is made when the time range according the specified quantile is
101  * smaller than the specified width (TTS) of the Gaussian distribution.
102  * In that case, the resulting PDF is a Gaussian distribution with the specified width (TTS)
103  * and normalisation according to the integral value of the input PDF.
104  * A smooth transition is imposed between the normal regime and this exeption.
105  *
106  * \param TTS TTS [ns]
107  * \param numberOfPoints number of points for Gauss-Hermite integration
108  * \param epsilon precision
109  * \param quantile quantile
110  */
111  void blur(const double TTS,
112  const int numberOfPoints = 25,
113  const double epsilon = 1.0e-10,
114  const double quantile = 0.99)
115  {
116  using namespace std;
117  using namespace JPP;
118 
119  typedef typename transformer_type::array_type array_type;
120 
121  const JGaussHermite engine(numberOfPoints, epsilon);
122 
123  for (super_iterator i = this->super_begin(); i != this->super_end(); ++i) {
124 
125  const array_type array = (*i).getKey();
126  function_type& f1 = (*i).getValue();
127 
128  if (!f1.empty()) {
129 
130  const typename function_type::supervisor_type& supervisor = f1.getSupervisor();
131 
132  const JMultiMapGetTransformer<NUMBER_OF_DIMENSIONS - 1, value_type> get(*(this->transformer), array);
133  const JMultiMapPutTransformer<NUMBER_OF_DIMENSIONS - 1, value_type> put(*(this->transformer), array);
134 
135  f1.transform(get);
136  f1.compile();
137 
138  const JQuantiles Q(f1, quantile);
139 
140  // abscissa
141 
142  JSet<double> X;
143 
144  for (JGaussHermite::const_iterator j = engine.begin(); j != engine.end(); ++j) {
145  X.insert(Q.getX() + TTS*sqrt(2.0)*j->getX());
146  }
147 
148  for (typename function_type::const_iterator j = f1.begin(); j != f1.end(); ++j) {
149 
150  if (j->getX() - TTS < X.getXmin()) {
151  X.insert(j->getX() - TTS);
152  }
153 
154  if (j->getX() + TTS > X.getXmax()) {
155  X.insert(j->getX() + TTS);
156  }
157  }
158 
159 
160  const double W = gauss(Q.getUpperLimit() - Q.getLowerLimit(), TTS);
161 
163 
164  for (JSet<double>::const_iterator x = X.begin(); x != X.end(); ++x) {
165 
166  double y = 0.0;
167 
168  for (JGaussHermite::const_iterator j = engine.begin(); j != engine.end(); ++j) {
169 
170  const double u = j->getX();
171  const double v = j->getY() / sqrt(PI);
172  const double w = get_value(f1(*x + u*TTS*sqrt(2.0)));
173 
174  y += v * w;
175  }
176 
177  buffer[*x] = W * Q.getIntegral() * Gauss(*x - Q.getX(), TTS) + (1.0 - W) * y;
178  }
179 
180  buffer.transform(put);
181  buffer.compile();
182 
183  f1 = buffer;
184 
185  f1.setExceptionHandler(supervisor);
186  }
187  }
188  }
189 
190 
191  /**
192  * Compresses PDF to given abscissa range.
193  *
194  * \param range abscissa range
195  */
197  {
198  for (super_iterator i = this->super_begin(); i != this->super_end(); ++i) {
199 
200  function_type& f1 = i.getValue();
201  function_type f2;
202 
203  typename function_type::iterator p = f1.lower_bound(range.getLowerLimit());
204  typename function_type::iterator q = f1.lower_bound(range.getUpperLimit());
205 
206  for (++q; p != q; ++p) {
207  f2.insert(*p);
208  }
209 
210  f2.swap(f1);
211  }
212 
213  this->compile();
214  }
215 
216 
217  /**
218  * Read from input.
219  *
220  * \param in reader
221  * \return reader
222  */
223  virtual JReader& read(JReader& in) override
224  {
225  if (in >> static_cast<transformablemultifunction_type&>(*this)) {
226 
227  // read optional transformer
228 
230 
231  if (buffer.read(in)) {
232 
233  this->transformer.reset(buffer.clone());
234 
235  } else {
236 
237  in.clear();
238 
240  }
241  }
242 
243  this->compile();
244 
245  return in;
246  }
247 
248 
249  /**
250  * Write from input.
251  *
252  * \param out writer
253  * \return writer
254  */
255  virtual JWriter& write(JWriter& out) const override
256  {
257  out << static_cast<const transformablemultifunction_type&>(*this);
258 
259  this->transformer->write(out);
260 
261  return out;
262  }
263  };
264 }
265 
266 #endif
function_type::argument_type argument_type
void insert(const JTransformableMultiFunction< __JFunction_t, __JMaplist_t, __JDistance_t > &input)
Insert multidimensional input.
void compile()
Compilation.
virtual JWriter & write(JWriter &out) const override
Write from input.
Definition: JPDFTable.hh:255
data_type w[N+1][M+1]
Definition: JPolint.hh:867
Interface for binary output.
Q(UTCMax_s-UTCMin_s)-livetime_s
multifunction_type::result_type result_type
T getLowerLimit() const
Get lower limit.
Definition: JRange.hh:202
multimap_type::super_iterator super_iterator
Auxiliary base class for storing and loading a single object to and from a binary file...
Auxiliary methods for mathematics.
transformablemultifunction_type::multimap_type multimap_type
Definition: JPDFTable.hh:53
Transformable multidimensional function.
JMultiMapTransformer< JMapLength< JMaplist_t >::value, argument_type > transformer_type
Template class for distance evaluation.
Definition: JDistance.hh:24
multifunction_type::argument_type argument_type
transformablemultifunction_type::function_type function_type
Definition: JPDFTable.hh:60
Template definition of transformer of the probability density function (PDF) of the time response of ...
void blur(const double TTS, const int numberOfPoints=25, const double epsilon=1.0e-10, const double quantile=0.99)
Blur PDF.
Definition: JPDFTable.hh:111
transformablemultifunction_type::value_type value_type
Definition: JPDFTable.hh:51
JTransformableMultiFunction< JFunction1D_t, JMaplist_t, JDistance_t > transformablemultifunction_type
Definition: JPDFTable.hh:47
JMultiMap< typename JFunction_t::argument_type, JFunction_t, JMaplist_t, JDistance_t > multimap_type
void compress(const JRange< typename function_type::abscissa_type > &range)
Compresses PDF to given abscissa range.
Definition: JPDFTable.hh:196
Forward declaration of binary output.
const JPolynome f1(1.0, 2.0, 3.0)
Function.
function_type::value_type value_type
JArray< NUMBER_OF_DIMENSIONS, argument_type > buffer
Transformable multidimensional histogram.
Multi-dimensional PDF table for arrival time of Cherenkov light.
Definition: JPDFTable.hh:40
JPDFTable()
Default constructor.
Definition: JPDFTable.hh:68
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
Physics constants.
static const double PI
Mathematical constants.
static JMultiMapTransformer * getClone()
Get clone of default transformer.
Interface for binary input.
One dimensional array of template objects with fixed length.
Definition: JArray.hh:40
Range of values.
Definition: JRange.hh:38
JPDFTable(const JTransformableMultiFunction< __JFunction_t, __JMaplist_t, __JDistance_t > &input)
Constructor.
Definition: JPDFTable.hh:79
transformablemultifunction_type::transformer_type transformer_type
Definition: JPDFTable.hh:54
transformablemultifunction_type::argument_type argument_type
Definition: JPDFTable.hh:49
z range($ZMAX-$ZMIN)< $MINIMAL_DZ." fi fi typeset -Z 4 STRING typeset -Z 2 FLOOR JPlot1D -f $
virtual JReader & read(JReader &in) override
Read from input.
Definition: JPDFTable.hh:223
Auxiliary class to define a range between two values.
then fatal The output file must have the wildcard in the e g root fi eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:48
then set_variable NUMBER_OF_TESTS else set_variable NUMBER_OF_TESTS fi function gauss()
int numberOfPoints
Definition: JResultPDF.cc:22
double Gauss(const double x, const double sigma)
Normalised Gauss function.
transformablemultifunction_type::super_const_iterator super_const_iterator
Definition: JPDFTable.hh:58
JPDFTable(const JTransformableMultiHistogram< JHistogram_t, __JMaplist_t, __JDistance_t > &input)
Constructor.
Definition: JPDFTable.hh:90
double getValue(const double x) const
Function value.
Definition: JMathlib.hh:1421
no fit printf nominal n $STRING awk v X
container_type::const_iterator const_iterator
Definition: JCollection.hh:91
transformablemultifunction_type::super_iterator super_iterator
Definition: JPDFTable.hh:59
int j
Definition: JPolint.hh:792
transformablemultifunction_type::result_type result_type
Definition: JPDFTable.hh:50
multifunction_type::super_iterator super_iterator
data_type v[N+1][M+1]
Definition: JPolint.hh:866
double u[N+1]
Definition: JPolint.hh:865
JLANG::JSharedPointer< transformer_type > transformer
multifunction_type::super_const_iterator super_const_iterator
virtual void clear()
Clear status of reader.
JResultEvaluator< JResult_t >::result_type get_value(const JResult_t &value)
Helper method to recursively evaluate a to function value.
Definition: JResult.hh:998
const double epsilon
Definition: JQuadrature.cc:21