Jpp  master_rocky
the software that should make you happy
JSpline.hh
Go to the documentation of this file.
1 #ifndef __JTOOLS__JSPLINE__
2 #define __JTOOLS__JSPLINE__
3 
4 #include <utility>
5 #include <sstream>
6 
7 #include "JMath/JZero.hh"
8 #include "JLang/JException.hh"
9 #include "JLang/JClass.hh"
11 #include "JTools/JFunctional.hh"
12 #include "JTools/JDistance.hh"
13 #include "JTools/JResult.hh"
15 #include "JTools/JMapCollection.hh"
16 
17 
18 /**
19  * \author mdejong
20  */
21 
22 namespace JTOOLS {}
23 namespace JPP { using namespace JTOOLS; }
24 
25 namespace JTOOLS {
26 
27  using JLANG::JNoValue;
31 
32 
33  /**
34  * Auxiliary class to define first derivates of the spline function at the two extrema.
35  */
36  template<class JOrdinate_t>
37  class JSplineBounds {
38  public:
39 
40  typedef JOrdinate_t ordinate_type;
42 
43 
44  /**
45  * Default constructor.
46  */
48  fp_at_xmin(false, ordinate_type()),
49  fp_at_xmax(false, ordinate_type())
50  {}
51 
52 
53  /**
54  * Constructor.
55  *
56  * \param fpAtXmin 1st derivative at minimal abscissa value
57  * \param fpAtXmax 1st derivative at maximal abscissa value
58  */
60  argument_type fpAtXmax) :
61  fp_at_xmin(true, fpAtXmin),
62  fp_at_xmax(true, fpAtXmax)
63  {}
64 
65 
66  /**
67  * Set first derivative of function at minimal abscissa value.
68  *
69  * \param fp 1st derivative
70  */
72  {
73  fp_at_xmin.first = true;
74  fp_at_xmin.second = fp;
75  }
76 
77 
78  /**
79  * Set first derivative of function at maximal abscissa value.
80  *
81  * \param fp 1st derivative
82  */
84  {
85  fp_at_xmax.first = true;
86  fp_at_xmax.second = fp;
87  }
88 
89 
90  /**
91  * Has first derivative of function at minimal abscissa value.
92  *
93  * \return true if 1st derivative is set; else false
94  */
95  const bool& hasFirstDerivativeAtXmin() const
96  {
97  return fp_at_xmin.first;
98  }
99 
100 
101  /**
102  * Has first derivative of function at maximal abscissa value.
103  *
104  * \return true if 1st derivative is set; else false
105  */
106  const bool& hasFirstDerivativeAtXmax() const
107  {
108  return fp_at_xmax.first;
109  }
110 
111 
112  /**
113  * Get first derivative of function at minimal abscissa value.
114  *
115  * \return 1st derivative
116  */
118  {
119  if (fp_at_xmin.first)
120  return fp_at_xmin.second;
121  else
122  throw JNoValue("JSplineBounds: missing 1st derivative.");
123  }
124 
125 
126  /**
127  * Get first derivative of function at maximal abscissa value.
128  *
129  * \return 1st derivative
130  */
132  {
133  if (fp_at_xmax.first)
134  return fp_at_xmax.second;
135  else
136  throw JNoValue("JSplineBounds: missing 1st derivative.");
137  }
138 
139  protected:
142  };
143 
144 
145  /**
146  * Helper method for JSplineBounds.
147  *
148  * \param fpAtXmin 1st derivative at minimal abscissa value
149  * \param fpAtXmax 1st derivative at maximal abscissa value
150  * \return spline bounds
151  */
152  template<class JOrdinate_t>
153  inline JSplineBounds<JOrdinate_t> make_spline_bounds(const JOrdinate_t fpAtXmin,
154  const JOrdinate_t fpAtXmax)
155  {
156  return JSplineBounds<JOrdinate_t>(fpAtXmin, fpAtXmax);
157  }
158 
159 
160  /**
161  * Template base class for spline interpolations.
162  *
163  * This class partially implements the JFunctional interface.
164  *
165  * Note that the data structure of the elements in the collection should have the additional methods:
166  * <pre>
167  * ordinate_type getU() const;
168  * void setU(ordinate_type u);
169  * </pre>
170  * to get and set the second derivatives, respectively.
171  *
172  * Spline interpolation code is taken from reference:
173  * Numerical Recipes in C++, W.H. Press, S.A. Teukolsky, W.T. Vetterling and B.P. Flannery,
174  * Cambridge University Press.
175  */
176  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
178  public JCollection_t<JElement_t, JDistance_t>,
179  public virtual JFunctional<>
180  {
181  public:
182 
183  typedef JCollection_t<JElement_t, JDistance_t> collection_type;
184 
185  typedef typename collection_type::abscissa_type abscissa_type;
186  typedef typename collection_type::ordinate_type ordinate_type;
187  typedef typename collection_type::value_type value_type;
188 
189  typedef typename collection_type::const_iterator const_iterator;
190  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
191  typedef typename collection_type::iterator iterator;
192  typedef typename collection_type::reverse_iterator reverse_iterator;
193 
195 
196 
197  /**
198  * Determination of second derivatives with specified bounds.
199  *
200  * \param bounds 1st derivatives at two extrema.
201  */
203  {
204  const int numberOfElements = this->size();
205 
206  using namespace std;
207 
208  if (numberOfElements > 2) {
209 
210  std::vector<double> buffer(numberOfElements);
211 
212  if (bounds.hasFirstDerivativeAtXmin()) {
213 
214  iterator j = this->begin();
215  iterator i = j++;
216 
217  const double dx = this->getDistance(i->getX(), j->getX());
218  const ordinate_type dy = (j->getY() - i->getY());
219 
220  buffer[0] = -0.5;
221 
222  i->setU((3.0/dx) * (dy/dx - bounds.getFirstDerivativeAtXmin()));
223 
224  } else {
225 
226  buffer[0] = 0.0;
227 
228  this->begin()->setU(JMATH::zero);
229  }
230 
231  int index = 1;
232 
233  for (iterator k = this->begin(), i = k++, j = k++; k != this->end(); ++i, ++j, ++k, ++index) {
234 
235  const abscissa_type x1 = i->getX();
236  const abscissa_type x2 = j->getX();
237  const abscissa_type x3 = k->getX();
238 
239  const ordinate_type& y1 = i->getY();
240  const ordinate_type& y2 = j->getY();
241  const ordinate_type& y3 = k->getY();
242 
243  const double sig = this->getDistance(x1, x2) / this->getDistance(x1, x3);
244  const double h = sig * buffer[index-1] + 2.0;
245 
246  buffer[index] = (sig - 1.0) / h;
247 
248  j->setU((y3 - y2) / this->getDistance(x2, x3) -
249  (y2 - y1) / this->getDistance(x1, x2));
250 
251  j->setU((6.0 * j->getU() / this->getDistance(x1, x3) - sig * i->getU()) / h);
252  }
253 
254 
255  if (bounds.hasFirstDerivativeAtXmax()) {
256 
257  reverse_iterator j = this->rbegin();
258  reverse_iterator i = j++;
259 
260  index = numberOfElements - 2;
261 
262  const double dx = this->getDistance(i->getX(), j->getX());
263  const ordinate_type dy = (j->getY() - i->getY());
264 
265  i->setU((3.0/dx) * (bounds.getFirstDerivativeAtXmax() - dy/dx));
266 
267  i->setU((i->getU() - 0.5*j->getU()) / (0.5*buffer[index] + 1.0));
268 
269  } else {
270 
271  this->rbegin()->setU(JMATH::zero);
272  }
273 
274 
275  reverse_iterator j = this->rbegin();
276  reverse_iterator i = j++;
277 
278  index = numberOfElements - 2;
279 
280  for ( ; j != this->rend(); ++i, ++j, --index) {
281  j->setU(j->getU() + i->getU() * buffer[index]);
282  }
283 
284  } else {
285 
286  for (iterator i = this->begin(); i != this->end(); ++i) {
287  i->setU(JMATH::zero);
288  }
289  }
290  }
291 
292 
293  protected:
294  /**
295  * Default constructor.
296  */
298  {}
299 
300 
301  /**
302  * Determination of second derivatives with no bounds.
303  */
304  virtual void do_compile() override
305  {
307  }
308  };
309 
310 
311  /**
312  * Template definition for functional collection with spline interpolation.
313  */
314  template<class JElement_t,
315  template<class, class> class JCollection_t,
316  class JResult_t,
317  class JDistance_t>
319 
320 
321  /**
322  * Template specialisation for functional collection with spline interpolation.
323  */
324  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
325  class JSplineFunction<JElement_t,
326  JCollection_t,
327  typename JResultType<typename JElement_t::ordinate_type>::result_type,
328  JDistance_t> :
329  public JSplineCollection<JElement_t, JCollection_t, JDistance_t>,
330  public virtual JFunction<typename JElement_t::abscissa_type,
331  typename JResultType<typename JElement_t::ordinate_type>::result_type>
332  {
333  public:
334 
336 
337  typedef typename collection_type::abscissa_type abscissa_type;
338  typedef typename collection_type::ordinate_type ordinate_type;
339  typedef typename collection_type::value_type value_type;
340  typedef typename collection_type::distance_type distance_type;
341 
342  typedef typename collection_type::const_iterator const_iterator;
343  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
344  typedef typename collection_type::iterator iterator;
345  typedef typename collection_type::reverse_iterator reverse_iterator;
346 
349 
353 
354 
355  /**
356  * Default constructor.
357  */
359  {}
360 
361 
362  /**
363  * Recursive interpolation method implementation.
364  *
365  * \param pX pointer to abscissa values
366  * \return function value
367  */
368  virtual result_type evaluate(const argument_type* pX) const override
369  {
370  const argument_type x = *pX;
371 
372  if (this->size() > 1u) {
373 
374  const_iterator p = this->lower_bound(x);
375 
376  if ((p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) ||
377  (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision)) {
378 
379  try {
380  return this->getExceptionHandler().action();
381  }
382  catch (const JException& error) {
383 
384  std::ostringstream os;
385 
386  os << __FILE__ << ':' << __LINE__ << " abscissa out of range "
387  << STREAM("?") << x << " <> "
388  << STREAM("?") << this->begin() ->getX() << ' '
389  << STREAM("?") << this->rbegin()->getX();
390 
391  throw JValueOutOfRange(os.str());
392  }
393  }
394 
395  const_iterator q = p--;
396 
397  const double dx = this->getDistance(p->getX(), q->getX());
398  const double a = this->getDistance(x, q->getX()) / dx;
399  const double b = 1.0 - a;
400 
401  return (a * p->getY() + b * q->getY()
402  - a*b * ((a + 1.0)*p->getU() + (b + 1.0)*q->getU()) * dx*dx/6);
403 
404  } else if (this->size() == 1u && this->getDistance(x, this->begin()->getX()) <= distance_type::precision) {
405 
406  return this->begin()->getY();
407 
408  } else {
409 
410  try {
411  return this->getExceptionHandler().action();
412  }
413  catch (const JException& error) {
414 
415  std::ostringstream os;
416 
417  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
418 
419  throw JFunctionalException(os.str());
420  }
421  }
422  }
423  };
424 
425 
426  /**
427  * Template specialisation for spline interpolation method with returning JResultDerivative data structure.
428  */
429  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
430  class JSplineFunction<JElement_t,
431  JCollection_t,
432  JResultDerivative<typename JResultType<typename JElement_t::ordinate_type>::result_type>,
433  JDistance_t> :
434  public JSplineCollection<JElement_t, JCollection_t, JDistance_t>,
435  public virtual JFunction<typename JElement_t::abscissa_type,
436  JResultDerivative<typename JResultType<typename JElement_t::ordinate_type>::result_type> >
437  {
438  public:
439 
441 
442  typedef typename collection_type::abscissa_type abscissa_type;
443  typedef typename collection_type::ordinate_type ordinate_type;
444  typedef typename collection_type::value_type value_type;
445  typedef typename collection_type::distance_type distance_type;
446 
447  typedef typename collection_type::const_iterator const_iterator;
448  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
449  typedef typename collection_type::iterator iterator;
450  typedef typename collection_type::reverse_iterator reverse_iterator;
451 
454 
458 
460 
461 
462  /**
463  * Default constructor.
464  */
466  {}
467 
468 
469  /**
470  * Recursive interpolation method implementation.
471  *
472  * \param pX pointer to abscissa values
473  * \return function value
474  */
475  virtual result_type evaluate(const argument_type* pX) const override
476  {
477  const argument_type x = *pX;
478 
479  if (this->size() <= 1u) {
480 
481  try {
482  return this->getExceptionHandler().action();
483  }
484  catch (const JException& error) {
485 
486  std::ostringstream os;
487 
488  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
489 
490  throw JFunctionalException(os.str());
491  }
492  }
493 
494  const_iterator p = this->lower_bound(x);
495 
496 
497  if ((p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) ||
498  (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision)) {
499 
500  try {
501  return this->getExceptionHandler().action();
502  }
503  catch (const JException& error) {
504 
505  std::ostringstream os;
506 
507  os << __FILE__ << ':' << __LINE__ << " abscissa out of range "
508  << STREAM("?") << x << " <> "
509  << STREAM("?") << this->begin() ->getX() << ' '
510  << STREAM("?") << this->rbegin()->getX();
511 
512  throw JValueOutOfRange(os.str());
513  }
514  }
515 
516  const_iterator q = p--;
517 
518  const double dx = this->getDistance(p->getX(), q->getX());
519  const double a = this->getDistance(x, q->getX()) / dx;
520  const double b = 1.0 - a;
521 
522  result.f = a * p->getY() + b * q->getY()
523  - a*b * ((a + 1.0)*p->getU() + (b + 1.0)*q->getU()) * dx*dx/6;
524 
525  result.fp = (q->getY() - p->getY() + (p->getU()*(1.0 - 3*a*a) -
526  q->getU()*(1.0 - 3*b*b)) * dx*dx/6) / dx;
527 
528  return result;
529  }
530 
531 
532  private:
534  };
535 
536 
537  /**
538  * Template specialisation for spline interpolation method with returning JResultPDF data structure.
539  *
540  * Note that the data structure of the elements in the collection should have the additional methods:
541  * <pre>
542  * ordinate_type getIntegral() const;
543  * void setIntegral(ordinate_type v);
544  * </pre>
545  * to get and set the integral values, respectively.
546  */
547  template<class JElement_t, template<class, class> class JCollection_t, class JDistance_t>
548  class JSplineFunction<JElement_t,
549  JCollection_t,
550  JResultPDF<typename JResultType<typename JElement_t::ordinate_type>::result_type>,
551  JDistance_t> :
552  public JSplineCollection<JElement_t, JCollection_t, JDistance_t>,
553  public virtual JFunction<typename JElement_t::abscissa_type,
554  JResultPDF<typename JResultType<typename JElement_t::ordinate_type>::result_type> >
555  {
556  public:
557 
559 
560  typedef typename collection_type::abscissa_type abscissa_type;
561  typedef typename collection_type::ordinate_type ordinate_type;
562  typedef typename collection_type::value_type value_type;
563  typedef typename collection_type::distance_type distance_type;
564 
565  typedef typename collection_type::const_iterator const_iterator;
566  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
567  typedef typename collection_type::iterator iterator;
568  typedef typename collection_type::reverse_iterator reverse_iterator;
569 
572 
576 
578 
579 
580  /**
581  * Default constructor.
582  */
584  {}
585 
586 
587  /**
588  * Determination of second derivatives with specified bounds.
589  *
590  * \param bounds 1st derivatives at two extrema.
591  */
593  {
594  if (this->size() >= 2u) {
595 
596  collection_type::compile(bounds);
597 
598  this->begin()->setIntegral(JMATH::zero);
599 
600  for (iterator j = this->begin(), i = j++; j != this->end(); ++i, ++j) {
601 
602  const double dx = this->getDistance(i->getX(), j->getX());
603  const ordinate_type y = i->getY() + j->getY();
604  const ordinate_type z = i->getU() + j->getU();
605 
606  const ordinate_type v = dx * 0.50 * y;
607  const ordinate_type w = dx * 0.25 * z*dx*dx/6;
608 
609  j->setIntegral(i->getIntegral() + v - w);
610  }
611  }
612  }
613 
614 
615  /**
616  * Recursive interpolation method implementation.
617  *
618  * \param pX pointer to abscissa values
619  * \return function value
620  */
621  virtual result_type evaluate(const argument_type* pX) const override
622  {
623  const argument_type x = *pX;
624 
625  if (this->size() <= 1u) {
626 
627  try {
628  return this->getExceptionHandler().action();
629  }
630  catch (const JException& error) {
631 
632  std::ostringstream os;
633 
634  os << __FILE__ << ':' << __LINE__ << " not enough data " << STREAM("?") << x;
635 
636  throw JFunctionalException(os.str());
637  }
638  }
639 
640  const_iterator p = this->lower_bound(x);
641 
642  if (p == this->begin() && this->getDistance(x, (p++)->getX()) > distance_type::precision) {
643 
644  try {
645 
646  result = this->getExceptionHandler().action();
647 
648  // overwrite integral values
649 
650  result.v = 0;
651  result.V = this->rbegin()->getIntegral();
652 
653  } catch(const JException& error) {
654 
655  std::ostringstream os;
656 
657  os << __FILE__ << ':' << __LINE__ << " abscissa out of range " << STREAM("?") << x << " < " << STREAM("?") << this->begin() ->getX();
658 
659  throw JValueOutOfRange(os.str());
660  }
661 
662  return result;
663 
664  } else if (p == this->end() && this->getDistance((--p)->getX(), x) > distance_type::precision) {
665 
666  try {
667 
668  result = this->getExceptionHandler().action();
669 
670  // overwrite integral values
671 
672  result.v = this->rbegin()->getIntegral();
673  result.V = this->rbegin()->getIntegral();
674 
675  } catch(const JException& error) {
676 
677  std::ostringstream os;
678 
679  os << __FILE__ << ':' << __LINE__ << " abscissa out of range " << STREAM("?") << x << " > " << STREAM("?") << this->rbegin()->getX();
680 
681  throw JValueOutOfRange(os.str());
682  }
683 
684  return result;
685  }
686 
687  const_iterator q = p--;
688 
689  const double dx = this->getDistance(p->getX(), q->getX());
690  const double a = this->getDistance(x, q->getX()) / dx;
691  const double b = 1.0 - a;
692 
693  result.f = a * p->getY() + b * q->getY()
694  - a*b * ((a + 1.0)*p->getU() + (b + 1.0)*q->getU()) * dx*dx/6;
695 
696  result.fp = (q->getY() - p->getY() + (p->getU()*(1.0 - 3*a*a) -
697  q->getU()*(1.0 - 3*b*b)) * dx*dx/6) / dx;
698 
699  result.v = p->getIntegral()
700  + 0.5*dx * (p->getY() - 0.5*p->getU()*dx*dx/6)
701  - 0.5*dx * ((a*a*p->getY() - b*b*q->getY()) +
702  (p->getU() * a*a*(0.5*a*a - 1.0) -
703  q->getU() * b*b*(0.5*b*b - 1.0)) * dx*dx/6);
704 
705  result.V = this->rbegin()->getIntegral();
706 
707  return result;
708  }
709 
710 
711  protected:
712  /**
713  * Determination of second derivatives with no bounds.
714  */
715  virtual void do_compile() override
716  {
718  }
719 
720 
721  private:
723  };
724 
725 
726  /**
727  * Template class for spline interpolation in 1D
728  *
729  * This class implements the JFunction1D interface.
730  */
731  template<class JElement_t,
732  template<class, class> class JCollection_t,
733  class JResult_t = typename JElement_t::ordinate_type,
736  public JSplineFunction<JElement_t, JCollection_t, JResult_t, JDistance_t>,
737  public JFunction1D<typename JElement_t::abscissa_type, JResult_t>
738  {
739  public:
740 
742 
746  typedef typename collection_type::distance_type distance_type;
747 
752 
754 
758 
759 
760  /**
761  * Default contructor.
762  */
764  {}
765  };
766 
767 
768  /**
769  * \cond NEVER
770  * Forward declarations.
771  * \endcond
772  */
773  template<class JAbscissa_t, class JOrdinate_t>
774  struct JSplineElement2D;
775 
776  template<template<class, class, class> class JMap_t>
777  struct JMapCollection;
778 
779 
780  /**
781  * Functional map with spline interpolation.
782  */
783  template<class JKey_t,
784  class JValue_t,
785  template<class, class, class> class JMap_t,
786  class JResult_t,
787  class JDistance_t = JDistance<JKey_t> >
788  class JSplineMap :
789  public JMap_t<JKey_t, JValue_t, JDistance_t>,
790  public JFunction<JKey_t, JResult_t>
791  {
792  public:
793 
794  typedef JMap_t<JKey_t, JValue_t, JDistance_t> collection_type;
796 
797  typedef typename collection_type::abscissa_type abscissa_type;
798  typedef typename collection_type::ordinate_type ordinate_type;
799  typedef typename collection_type::value_type value_type;
800  typedef typename collection_type::distance_type distance_type;
801 
802  typedef typename collection_type::const_iterator const_iterator;
803  typedef typename collection_type::const_reverse_iterator const_reverse_iterator;
804  typedef typename collection_type::iterator iterator;
805  typedef typename collection_type::reverse_iterator reverse_iterator;
806 
810 
815 
816 
817  /**
818  * Default constructor.
819  */
821  {}
822 
823 
824  /**
825  * Recursive interpolation method implementation.
826  *
827  * \param pX pointer to abscissa values
828  * \return function value
829  */
830  virtual result_type evaluate(const argument_type* pX) const override
831  {
832  const argument_type x = *pX;
833 
834  ++pX; // next argument value
835 
836  const_iterator p = this->begin();
837 
838  for (typename JSplineFunction1D_t::iterator q = buffer.begin(); q != buffer.end(); ++q, ++p) {
839  q->getY() = JFunction<argument_type, data_type>::getValue(p->getY(), pX);
840  }
841 
842  buffer.compile();
843 
844  return buffer(x);
845  }
846 
847 
848  private:
849  /**
850  * Function compilation.
851  */
852  virtual void do_compile() override
853  {
854  buffer.clear();
855 
856  for (iterator i = this->begin(); i != this->end(); ++i) {
857  buffer.put(i->getX(), data_type());
858  }
859  }
860 
861 
863  };
864 
865 
866  /**
867  * Conversion of data points to integral values.
868  *
869  * The integration includes the use of 2nd derivatives of the data points of the input spline interpolating function.
870  *
871  * \param input collection
872  * \param output mappable collection
873  * \return integral
874  */
875  template<class JElement_t,
876  template<class, class> class JCollection_t,
877  class JResult_t,
878  class JDistance_t>
879  inline typename JElement_t::ordinate_type
881  typename JMappable<JElement_t>::map_type& output)
882  {
883  typedef typename JElement_t::ordinate_type ordinate_type;
885 
887 
888  if (input.getSize() > 1) {
889 
890  output.put(input.begin()->getX(), V);
891 
892  for (const_iterator j = input.begin(), i = j++; j != input.end(); ++i, ++j) {
893 
894  const double dx = input.getDistance(i->getX(), j->getX());
895  const ordinate_type y = i->getY() + j->getY();
896  const ordinate_type z = i->getU() + j->getU();
897  const ordinate_type v = dx * 0.50 * y;
898  const ordinate_type w = dx * 0.25 * z*dx*dx/6;
899 
900  V += v - w;
901 
902  output.put(j->getX(), V);
903  }
904  }
905 
906  return V;
907  }
908 
909 
910  /**
911  * Conversion of data points to integral values.
912  *
913  * The integration directly uses the integral values of the input spline interpolating function.
914  *
915  * \param input collection
916  * \param output mappable collection
917  * \return integral
918  */
919  template<class JElement_t,
920  template<class, class> class JCollection_t,
921  class JDistance_t>
922  inline typename JElement_t::ordinate_type
923  integrate(const JSplineFunction1D<JElement_t, JCollection_t, JResultPDF<typename JElement_t::ordinate_type>, JDistance_t>& input,
924  typename JMappable<JElement_t>::map_type& output)
925  {
926  typedef typename JElement_t::ordinate_type ordinate_type;
929 
930  if (input.getSize() > 1) {
931 
932  for (const_iterator i = input.begin(); i != input.end(); ++i) {
933  output.put(i->getX(), i->getIntegral());
934  }
935 
936  return input.rbegin()->getIntegral();
937  }
938 
939  return JMATH::zero;
940  }
941 }
942 
943 #endif
Exceptions.
This include file containes various data structures that can be used as specific return types for the...
Definition of zero value for any class.
Exception for division by zero.
Definition: JException.hh:288
General exception.
Definition: JException.hh:24
Exception for a functional operation.
Definition: JException.hh:144
Exception for missing value.
Definition: JException.hh:216
Exception for accessing a value in a collection that is outside of its range.
Definition: JException.hh:180
Template definition of function object interface.
Definition: JFunctional.hh:72
JArgument_t argument_type
Definition: JFunctional.hh:84
const JExceptionHandler & getExceptionHandler() const
Get exception handler.
Definition: JFunctional.hh:292
JResult_t result_type
Definition: JFunctional.hh:87
static result_type getValue(const JFunctional &function, const argument_type *pX)
Recursive function value evaluation.
Definition: JFunctional.hh:107
Auxiliary class to define first derivates of the spline function at the two extrema.
Definition: JSpline.hh:37
JLANG::JClass< ordinate_type >::argument_type argument_type
Definition: JSpline.hh:41
std::pair< bool, ordinate_type > fp_at_xmin
Definition: JSpline.hh:140
JSplineBounds()
Default constructor.
Definition: JSpline.hh:47
const bool & hasFirstDerivativeAtXmin() const
Has first derivative of function at minimal abscissa value.
Definition: JSpline.hh:95
void setFirstDerivativeAtXmax(argument_type fp)
Set first derivative of function at maximal abscissa value.
Definition: JSpline.hh:83
JOrdinate_t ordinate_type
Definition: JSpline.hh:40
void setFirstDerivativeAtXmin(argument_type fp)
Set first derivative of function at minimal abscissa value.
Definition: JSpline.hh:71
std::pair< bool, ordinate_type > fp_at_xmax
Definition: JSpline.hh:141
ordinate_type getFirstDerivativeAtXmin() const
Get first derivative of function at minimal abscissa value.
Definition: JSpline.hh:117
ordinate_type getFirstDerivativeAtXmax() const
Get first derivative of function at maximal abscissa value.
Definition: JSpline.hh:131
const bool & hasFirstDerivativeAtXmax() const
Has first derivative of function at maximal abscissa value.
Definition: JSpline.hh:106
JSplineBounds(argument_type fpAtXmin, argument_type fpAtXmax)
Constructor.
Definition: JSpline.hh:59
Template base class for spline interpolations.
Definition: JSpline.hh:180
JSplineCollection()
Default constructor.
Definition: JSpline.hh:297
collection_type::const_iterator const_iterator
Definition: JSpline.hh:189
collection_type::abscissa_type abscissa_type
Definition: JSpline.hh:185
void compile(const JSplineBounds< ordinate_type > &bounds)
Determination of second derivatives with specified bounds.
Definition: JSpline.hh:202
JCollection_t< JElement_t, JDistance_t > collection_type
Definition: JSpline.hh:183
collection_type::const_reverse_iterator const_reverse_iterator
Definition: JSpline.hh:190
collection_type::value_type value_type
Definition: JSpline.hh:187
virtual void do_compile() override
Determination of second derivatives with no bounds.
Definition: JSpline.hh:304
collection_type::reverse_iterator reverse_iterator
Definition: JSpline.hh:192
collection_type::ordinate_type ordinate_type
Definition: JSpline.hh:186
collection_type::iterator iterator
Definition: JSpline.hh:191
Template class for spline interpolation in 1D.
Definition: JSpline.hh:738
collection_type::const_iterator const_iterator
Definition: JSpline.hh:748
collection_type::reverse_iterator reverse_iterator
Definition: JSpline.hh:751
collection_type::ordinate_type ordinate_type
Definition: JSpline.hh:744
JFunction1D< abscissa_type, JResult_t > function_type
Definition: JSpline.hh:753
function_type::argument_type argument_type
Definition: JSpline.hh:755
function_type::JExceptionHandler exceptionhandler_type
Definition: JSpline.hh:757
function_type::result_type result_type
Definition: JSpline.hh:756
collection_type::iterator iterator
Definition: JSpline.hh:750
collection_type::distance_type distance_type
Definition: JSpline.hh:746
collection_type::value_type value_type
Definition: JSpline.hh:745
collection_type::const_reverse_iterator const_reverse_iterator
Definition: JSpline.hh:749
collection_type::abscissa_type abscissa_type
Definition: JSpline.hh:743
JSplineFunction1D()
Default contructor.
Definition: JSpline.hh:763
JSplineCollection< JElement_t, JCollection_t, JDistance_t > collection_type
Definition: JSpline.hh:741
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
Definition: JSpline.hh:475
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
Definition: JSpline.hh:621
void compile(const JSplineBounds< ordinate_type > &bounds)
Determination of second derivatives with specified bounds.
Definition: JSpline.hh:592
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
Definition: JSpline.hh:368
Template definition for functional collection with spline interpolation.
Definition: JSpline.hh:318
Functional map with spline interpolation.
Definition: JSpline.hh:791
function_type::result_type result_type
Definition: JSpline.hh:808
JSplineFunction1D< JSplineElement2D< argument_type, data_type >, JMapCollection< JMap_t >::template collection_type, result_type > JSplineFunction1D_t
Definition: JSpline.hh:814
collection_type::value_type value_type
Definition: JSpline.hh:799
collection_type::const_reverse_iterator const_reverse_iterator
Definition: JSpline.hh:803
JSplineMap()
Default constructor.
Definition: JSpline.hh:820
function_type::JExceptionHandler exceptionhandler_type
Definition: JSpline.hh:809
JMap_t< JKey_t, JValue_t, JDistance_t > collection_type
Definition: JSpline.hh:794
function_type::argument_type argument_type
Definition: JSpline.hh:807
JFunction< JKey_t, JResult_t > function_type
Definition: JSpline.hh:795
JResultType< ordinate_type >::result_type data_type
Definition: JSpline.hh:811
collection_type::const_iterator const_iterator
Definition: JSpline.hh:802
collection_type::abscissa_type abscissa_type
Definition: JSpline.hh:797
collection_type::distance_type distance_type
Definition: JSpline.hh:800
collection_type::reverse_iterator reverse_iterator
Definition: JSpline.hh:805
JSplineFunction1D_t buffer
Definition: JSpline.hh:862
collection_type::iterator iterator
Definition: JSpline.hh:804
virtual void do_compile() override
Function compilation.
Definition: JSpline.hh:852
virtual result_type evaluate(const argument_type *pX) const override
Recursive interpolation method implementation.
Definition: JSpline.hh:830
collection_type::ordinate_type ordinate_type
Definition: JSpline.hh:798
const double a
Definition: JQuadrature.cc:42
static const JZero zero
Function object to assign zero value.
Definition: JZero.hh:105
double getDistance(const JFirst_t &first, const JSecond_t &second)
Get distance between objects.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
std::vector< event_type > data_type
Definition: JPerth.cc:82
Auxiliary classes and methods for multi-dimensional interpolations and histograms.
data_type w[N+1][M+1]
Definition: JPolint.hh:867
data_type v[N+1][M+1]
Definition: JPolint.hh:866
JElement_t::ordinate_type integrate(const JSplineFunction1D< JElement_t, JCollection_t, JResultPDF< typename JElement_t::ordinate_type >, JDistance_t > &input, typename JMappable< JElement_t >::map_type &output)
Conversion of data points to integral values.
Definition: JSpline.hh:923
int j
Definition: JPolint.hh:792
JSplineBounds< JOrdinate_t > make_spline_bounds(const JOrdinate_t fpAtXmin, const JOrdinate_t fpAtXmax)
Helper method for JSplineBounds.
Definition: JSpline.hh:153
Definition: JSTDTypes.hh:14
JArgument< T >::argument_type argument_type
Definition: JClass.hh:82
Template definition of function object interface in one dimension.
Definition: JFunctional.hh:334
functional_type::result_type result_type
Definition: JFunctional.hh:339
functional_type::argument_type argument_type
Definition: JFunctional.hh:338
Template definition of function object interface in multidimensions.
Definition: JFunctional.hh:320
functional_type::argument_type argument_type
Definition: JFunctional.hh:322
functional_type::result_type result_type
Definition: JFunctional.hh:323
Exception handler for functional object.
Definition: JFunctional.hh:133
Template class to define the corresponding JCollection for a given template JMap.
Template interface definition for associative collection of elements.
void put(typename JClass< key_type > ::argument_type key, typename JClass< mapped_type >::argument_type value)
Put pair-wise element (key,value) into collection.
Data structure for result including value and first derivative of function.
Definition: JResult.hh:45
Data structure for result including value, first derivative and integrals of function.
Definition: JResult.hh:339
Auxiliary class to evaluate result type.
Definition: JFunctional.hh:392
2D Element for spline interpolations.
Definition: JElement.hh:155
Auxiliary data structure for handling std::ostream.