JToolsToolkit.hh

Go to the documentation of this file.

#ifndef __JTOOLS__JTOOLSTOOLKIT__
#define __JTOOLS__JTOOLSTOOLKIT__
 
#include <type_traits>
 
 
#include "JLang/JException.hh"
#include "JMath/JZero.hh"
#include "JTools/JAssembler.hh"
#include "JTools/JElement.hh"
#include "JTools/JCollection.hh"
#include "JTools/JGarbageCollection.hh"
#include "JTools/JRange.hh"
#include "JTools/JMultiFunction.hh"
#include "JTools/JMultiHistogram.hh"
#include "JTools/JFunctionalMap.hh"
 
 
/**
 * \file
 *
 * This include file contains various recursive methods to operate on multi-dimensional collections.
 * \author mdejong
 */
namespace JTOOLS {}
namespace JPP { using namespace JTOOLS; }
 
/**
 * Auxiliary classes and methods for multi-dimensional interpolations and histograms.
 */
namespace JTOOLS {
 
  using JLANG::JException;
  using JLANG::JEmptyCollection;
  using JLANG::JDivisionByZero;
 
  
  /**
   * Conversion of data points to cumulative probability distribition (CDF).
   *
   * Note that the data type of the input container should be preserved 
   * so that the corresponding method <tt>integrate()</tt> is used.
   *
   * \param  input             collection   
   * \param  output            mappable collection
   * \param  eps               minimal step size
   * \return                   integral value
   */
  template<class JContainer_t,
           class JKey_t,
           class JValue_t>
  inline typename JContainer_t::ordinate_type
  makeCDF(const JContainer_t&                        input,
          JMappableCollection<JKey_t, JValue_t>&     output, 
          const typename JContainer_t::ordinate_type eps = JMATH::zero)
  {
    typedef typename JContainer_t::ordinate_type       ordinate_type;
    typedef typename JContainer_t::abscissa_type       abscissa_type;
    typedef JElement2D<abscissa_type, ordinate_type>   element_type;
    typedef JCollection<element_type>                  buffer_type;
 
 
    if (input.getSize() > 1) {
      
      buffer_type buffer;
 
      const ordinate_type V = integrate(input, buffer);
      
      if (V == ordinate_type(JMATH::zero)) {
        THROW(JDivisionByZero, "Method makeCDF(): integral equals zero.");
      }
      
      output.clear();
      
      typename buffer_type::const_iterator i = buffer.begin();
            
      for ( ; i != buffer.end() && i->getY() <= 0.5 * eps * V; ++i) {}
      
      if (i != buffer.end()) {
        
        // force first point at (0, ymin)
      
        JKey_t   x = 0.0;
        JValue_t y = i->getX();
        
        output.put(x, y);
      
        JKey_t   xmax = x;
        JValue_t ymax = y;
      
        for (++i; i != buffer.end(); ++i) {
          
          x = i->getY() / V;
          y = i->getX();
          
          if (x > xmax) {
          
            ymax = y;
            
            if (x > xmax + eps) {
              
              output.put(x, y);
              
              xmax = x;
            }
          }
        }
        
        // force last point at (1, ymax)
        
        x = 1.0;
        y = ymax;
        
        output.put(x,y);
        
      } else {
        THROW(JDivisionByZero, "Method makeCDF(): no remaining data.");
      }
 
      return V;
 
    } else {
      THROW(JEmptyCollection, "Method makeCDF(): not sufficient input data.");
    }
  }
 
  
  /**
   * Get integral of input data points.
   *
   * \param  input           input data
   * \return                 integral value
   */
  template<class JContainer_t>
  inline typename JContainer_t::ordinate_type getIntegral(const JContainer_t& input)
  { 
    JGarbageCollection<typename JContainer_t::abscissa_type, typename JContainer_t::ordinate_type> garbage;
    
    return integrate(input, garbage);
  }
 
 
  /**
   * Get integral of input data points.
   *
   * \param  input           input data
   * \param  range           range
   * \return                 integral value
   */
  template<class JContainer_t>
  inline typename JContainer_t::ordinate_type getIntegral(const JContainer_t& input,
                                                          const JRange<typename JContainer_t::argument_type>& range)
  {
    JContainer_t output;
 
    integrate(input, output);
 
    if (!output.empty())
      return ((output.in_range(range.getUpperLimit()) ? get_value(output(range.getUpperLimit())) : output.rbegin()->getY()) -
              (output.in_range(range.getLowerLimit()) ? get_value(output(range.getLowerLimit())) : output. begin()->getY()));
    else
      return 0.0;
  }
 
 
  /**
   * Auxiliary method to get integral of input data points.
   *
   * \param  input           input data
   * \return                 integral value
   */
  template<class JFunction_t,
           template<class, class, class> class JMap_t,
           class JTail_t>
  inline typename JFunction_t::ordinate_type
  getIntegral(const JMultiMap<typename JFunction_t::abscissa_type, JFunction_t, JMapList<JMap_t, JTail_t>, typename JFunction_t::distance_type>& input)
  { 
    typedef typename JFunction_t::abscissa_type                                               abscissa_type;
    typedef typename JFunction_t::ordinate_type                                               ordinate_type;
    typedef typename JFunction_t::distance_type                                               distance_type;
    typedef typename JFunctionalMap<JMap_t>::template function_type<abscissa_type,
                                                                    ordinate_type,
                                                                    ordinate_type,
                                                                    distance_type>            buffer_type;
    typedef JMultiMap<abscissa_type, JFunction_t, JMapList<JMap_t, JTail_t>, distance_type>   multimap_type;
    
    static buffer_type buffer;
    
    buffer.configure(input);
    
    typename buffer_type::iterator out = buffer.begin();
    
    for (typename multimap_type::const_iterator in = input.begin(); in != input.end(); ++in, ++out) {
      out->getY() = getIntegral(in->getY());
    }
    
    buffer.compile();
    
    return getIntegral(buffer);
  }
  
 
  /**
   * Get integral of input data points.
   *
   * \param  input           input data
   * \return                 integral value
   */
  template<class JFunction_t,
           template<class, class, class> class JMap_t,
           class JTail_t>
  inline typename JFunction_t::ordinate_type
  getIntegral(const JMultiFunction<JFunction_t, JMapList<JMap_t, JTail_t>, typename JFunction_t::distance_type>& input)
  { 
    typedef typename JFunction_t::abscissa_type                                      abscissa_type;
    typedef typename JFunction_t::ordinate_type                                      ordinate_type;
    typedef typename JFunction_t::distance_type                                      distance_type;
    typedef typename JFunctionalMap<JMap_t>::template function_type<abscissa_type,
                                                                    ordinate_type,
                                                                    ordinate_type,
                                                                    distance_type>   buffer_type;
    typedef JMultiFunction<JFunction_t, JMapList<JMap_t, JTail_t>, distance_type>    multifunction_type;
    
    buffer_type buffer;
    
    for (typename multifunction_type::const_iterator i = input.begin(); i != input.end(); ++i) {      
      buffer.put(i->getX(), getIntegral(i->getY()));
    }
    
    buffer.compile();
    
    return getIntegral(buffer);
  }
 
 
  /**
   * Reset value.
   *
   * The value is set to (the equivalent of) zero, see method JMATH::getZero.
   *
   * \param  value           value
   */
  template<class T>
  inline void reset(T& value)
  {
    value = JMATH::getZero<T>();
  }
  
 
  /**
   * Recursive reset of collection.
   *
   * \param  collection      collection
   */
  template<class JElement_t, class JDistance_t>
  inline void reset(JCollection<JElement_t, JDistance_t>& collection)
  {
    typedef typename JCollection<JElement_t, JDistance_t>::iterator  iterator;
 
    for (iterator i = collection.begin(); i != collection.end(); ++i) {
      reset(i->getY());
    }
  }
 
  
  /**
   * Copy of input to output.
   *
   * The output value is set to the input value.
   *
   * \param  input           input
   * \param  output          output
   */
  template<class T>
  inline void copy(const T& input, T& output)
  {
    output = input;
  }
 
 
  /**
   * Recursive copy of input collection to output collection.
   *
   * \param  input           input
   * \param  output          output
   */
  template<class JElement_t, class JDistance_t, class JKey_t, class JValue_t>
  inline void copy(const JCollection<JElement_t, JDistance_t>& input, JMappableCollection<JKey_t, JValue_t>& output)
  {
    typedef typename JCollection<JElement_t, JDistance_t>::const_iterator  const_iterator;
 
    output.clear();
 
    for (const_iterator i = input.begin(); i != input.end(); ++i) {
      copy(i->getY(), output.get(i->getX()));
    }
  }
 
 
  /**
   * Configuration of value.
   *
   * This is a fall-back method; it does nothing.
   *
   * \param  value           value
   * \param  bounds          bounds
   * \param  option          false
   */
  template<class T, class JAbscissa_t>
  inline void configure(const T& value, const JAbstractCollection<JAbscissa_t>& bounds, std::false_type option)
  {}
 
 
  /**
   * Recursive configuration of collection.
   *
   * \param  collection      collection
   * \param  bounds          bounds
   * \param  option          true
   */
  template<class JElement_t, class JDistance_t>
  inline void configure(JCollection<JElement_t, JDistance_t>&                          collection,
                        const JAbstractCollection<typename JElement_t::abscissa_type>& bounds,
                        std::true_type                                                    option = std::true_type())
  {
    typedef typename JCollection<JElement_t, JDistance_t>::iterator       iterator;
    typedef typename JCollection<JElement_t, JDistance_t>::ordinate_type  ordinate_type;
 
    collection.configure(bounds);
 
    for (iterator i = collection.begin(); i != collection.end(); ++i) {
      configure(i->getY(), bounds, std::bool_constant<JAssembler<ordinate_type>::is_collection>());
    }
  }
 
 
  /**
   * Accumulation of value.
   *
   * This is a fall-back method; it does nothing.
   *
   * \param  value           value
   * \param  option          false
   */
  template<class T>
  inline void accumulate(T& value, std::false_type option)
  {}
 
 
  /**
   * Recursive accumulation of collection.
   *
   * \param  collection      collection
   * \param  option          true
   */
  template<class JElement_t, class JDistance_t>
  inline void accumulate(JCollection<JElement_t, JDistance_t>& collection, std::true_type option = std::true_type())
  {
    typedef typename JCollection<JElement_t, JDistance_t>::iterator       iterator;
    typedef typename JCollection<JElement_t, JDistance_t>::ordinate_type  ordinate_type;
 
    for (iterator i = collection.begin(); i != collection.end(); ++i) {
      accumulate(i->getY(), std::bool_constant<JAssembler<ordinate_type>::is_collection>());
    }
 
    if (collection.getSize() > 1) {
 
      for (iterator j = collection.begin(), i = j++; j != collection.end(); ++i, ++j) {
        j->getY() += i->getY();
      }
 
      reset(collection.begin()->getY());
    }
  }
}
 
#endif