#ifndef __JMULTIFUNCTION__
#define __JMULTIFUNCTION__

#include <memory>
#include <stdarg.h>
#include <functional>
#include <cmath>

#include "JTools/JCollection.hh"
#include "JTools/JMixedMultiMap.hh"
#include "JTools/JArray.hh"
#include "JTools/JFunctional.hh"
#include "JTools/JFunctionalMap.hh"
#include "JIO/JSerialisable.hh"


namespace JTOOLS {

  namespace {

    using JIO::JSerialisable;
    using JIO::JReader;
    using JIO::JWriter;


    /**
     * Template definition of auxiliary class to evaluate result type of multi-dimensional function.
     */
    template<class JFunction1D_t, class JMapList_t>
    class JResultType;
    
    
    /**
     * Template specialisation of auxiliary class to evaluate result type of multi-dimensional function.
     */
    template<class JFunction1D_t, template<class, class> class JMap_t, class JTail_t>
    class JResultType<JFunction1D_t, JMapList<JMap_t, JTail_t> >
    {
    public:
      
      typedef typename JFunction1D_t::argument_type                                key_type;
      typedef JMap_t<key_type, JMixedMultiMap<key_type, JFunction1D_t, JTail_t> >  JMultiMap_t;
      
      typedef typename JMultiMap_t::result_type                                    result_type;
    };
    
    
    /**
     * Terminator class to evaluate result type of multi-dimensional function.
     */
    template<class JFunction1D_t, template<class, class> class JMap_t>
    class JResultType<JFunction1D_t, JMapList<JMap_t> >
    {
    public:
      
      typedef typename JFunction1D_t::argument_type                                key_type;
      typedef JMap_t<key_type, JFunction1D_t>                                      JMultiMap_t;
      
      typedef typename JMultiMap_t::result_type                                    result_type;
    };
  }

 
  /**
   * Multi-dimensional interpolation method.
   *
   * The first template parameter refers to one-dimensional function object,
   * the second to the list of maps used and  
   * the third to the return type of the multi-dimensional functio.
   */
  template<class JFunction1D_t, class JMapList_t>
  class JMultiFunction :
    public JMixedMultiMap<typename JFunction1D_t::argument_type, JFunction1D_t, JMapList_t>
  {
  public:
    
    typedef JMixedMultiMap<typename JFunction1D_t::argument_type, JFunction1D_t, JMapList_t> JMultiMap_t;

    using JMultiMap_t::NUMBER_OF_DIMENSIONS;

    typedef typename JFunction1D_t::argument_type                         argument_type;
    typedef typename JResultType<JFunction1D_t, JMapList_t>::result_type  result_type;

    typedef typename JMultiMap_t::key_type                                key_type;
    typedef typename JMultiMap_t::mapped_type                             mapped_type;

    typedef typename JMultiMap_t::iterator                                iterator;
    typedef typename JMultiMap_t::reverse_iterator                        reverse_iterator;
    typedef typename JMultiMap_t::const_iterator                          const_iterator;
    typedef typename JMultiMap_t::const_reverse_iterator                  const_reverse_iterator;

    typedef typename JMultiMap_t::super_iterator                          super_iterator;
    typedef typename JMultiMap_t::super_const_iterator                    super_const_iterator;

    typedef JArray<NUMBER_OF_DIMENSIONS + 1, argument_type>               array_type;


    /**
     * Default constructor.
     */
    JMultiFunction() :
      JMultiMap_t()
    {}


    /**
     * Multi-dimensional interpolation method call.
     *
     * \param  x      comma seperated argument list
     * \return        function value
     */
    result_type operator()(const argument_type& x, ...) const
    {
      va_start(ap, x);

      typename array_type::iterator p = buffer.begin();

      *p = x;                            // first argument

      while (++p != buffer.end())
        *p = va_arg(ap, argument_type);  // next argument

      va_end(ap);

      return this->evaluate(buffer.data());
    }


  protected:
    mutable va_list    ap;
    mutable array_type buffer;
  };


  /**
   * Interface for weight application and coordinate transformation of function.
   */
  template<unsigned int N, class JArgument_t>
  class JFunctionTransformerInterface : 
    public JSerialisable 
  {
  public:

    typedef JArgument_t                     argument_type;
    typedef JArray<N, const argument_type>  array_type;

    /**
     * Default constructor.
     */
    JFunctionTransformerInterface() :
      JSerialisable()
    {}


    /**
     * Clone object.
     * 
     * \return            pointer to newly created JFunctionTransformerInterface
     */
    virtual JFunctionTransformerInterface* clone() const 
    {
      return new JFunctionTransformerInterface(*this);
    }


    /**
     * Evaluate xn value as a function of {x0, ..., xn-1}
     *
     * \param  buffer     x0 - xn-1 values
     * \param  xn         xn value
     * \return            new xn value
     */
    virtual argument_type putXn(const array_type& buffer, const argument_type& xn) const
    {
      return xn;
    }


    /**
     * Evaluate xn value as a function of {x0, ..., xn-1}
     *
     * \param  buffer     x0 - xn-1 values
     * \param  xn         xn value
     * \return            new xn value
     */
    virtual argument_type getXn(const array_type& buffer, const argument_type& xn) const
    {
      return xn;
    }


    /**
     * Weight function.
     *
     * \param  buffer     x0 - xn-1 values
     * \return            weight
     */
    virtual double getWeight(const array_type& buffer) const
    {
      return 1.0;
    }


    /**
     * Weight function.
     *
     * \param  x          comma seperated argument list
     * \return            weight
     */
    double getWeight(const argument_type& x, ...) const
    {
      va_list  ap;

      va_start(ap, x);

      JArray<N, argument_type> buffer;

      typename JArray<N, argument_type>::iterator p = buffer.begin();

      *p = x;                            // first argument

      while (++p != buffer.end())
        *p = va_arg(ap, argument_type);  // next argument

      va_end(ap);

      return getWeight(buffer);
    }


    /**
     * Read from input.
     *
     * \param  in     JReader
     * \return        JReader
     */
    virtual JReader& read(JReader& in) { return in; }


    /**
     * Write to output.
     *
     * \param  out    JWriter
     * \return        JWriter
     */
    virtual JWriter& write(JWriter& out) const { return out; }
  };


  /**
   * Auxiliary class to convert JFunctionTransformerInterface to JCollectionElementTransformer.
   */
  template<unsigned int N, class JElement_t>
  class JFunctionPutTransformer : 
    public JCollectionElementTransformer<JElement_t> 
  {
  public:

    typedef typename JElement_t::key_type                key_type;
    typedef typename JElement_t::mapped_type             mapped_type;
    typedef JFunctionTransformerInterface<N, key_type>   transformer_type;
    typedef typename transformer_type::array_type        array_type;


    /**
     * Constructor.
     *
     * \param  __transformer  JFunctionTransformerInterface<> object
     * \param  __buffer       x0 - xn-1 values
     */
    JFunctionPutTransformer(const transformer_type& __transformer,
			    const array_type&       __buffer) :
      transformer(__transformer),
      buffer     (__buffer),
      W(transformer.getWeight(buffer))
    {}


    /**
     * Transform element.
     *
     * \param  element  input  element
     * \return          output element
     */
    virtual JElement_t operator()(const JElement_t& element) const
    {
      const key_type&    xn = element.getX();
      const mapped_type& y  = element.getY();

      return JElement_t(transformer.putXn(buffer, xn), y / W);
    }
    
  private:    
    const transformer_type& transformer;
    const array_type&       buffer;
    const double W;
  };


  /**
   * Auxiliary class to convert JFunctionTransformerInterface to JCollectionElementTransformer.
   */
  template<unsigned int N, class JElement_t>
  class JFunctionGetTransformer : 
    public JCollectionElementTransformer<JElement_t> 
  {
  public:

    typedef typename JElement_t::key_type                key_type;
    typedef typename JElement_t::mapped_type             mapped_type;
    typedef JFunctionTransformerInterface<N, key_type>   transformer_type;
    typedef typename transformer_type::array_type        array_type;


    /**
     * Constructor.
     *
     * \param  __transformer  JFunctionTransformerInterface<> object
     * \param  __buffer       x0 - xn-1 values
     */
    JFunctionGetTransformer(const transformer_type& __transformer,
			    const array_type&       __buffer) :
      transformer(__transformer),
      buffer     (__buffer),
      W(transformer.getWeight(buffer))
    {}


    /**
     * Transform element.
     *
     * \param  element  input  element
     * \return          output element
     */
    virtual JElement_t operator()(const JElement_t& element) const
    {
      const key_type&    xn = element.getX();
      const mapped_type& y  = element.getY();

      return JElement_t(transformer.getXn(buffer, xn), y * W);
    }
    
  private:
    const transformer_type& transformer;
    const array_type&       buffer;
    const double W;
  };


  /**
   * Multi-dimensional interpolation method.
   *
   * A weight function is used to improve the accuracy of the linear interpolation.
   * The handling of the weight function is done in the method transform(..) and the function operator ()(..).
   */
  template<class JFunction1D_t, class JMapList_t>
  class JTransformableMultiFunction :
    public JMultiFunction<JFunction1D_t, JMapList_t>,
    public JSerialisable
  {
  public:

    typedef JMultiFunction<JFunction1D_t, JMapList_t>            JMultiMap_t;

    using JMultiMap_t::NUMBER_OF_DIMENSIONS;

    typedef typename JFunction1D_t::argument_type                argument_type;
    typedef typename JFunction1D_t::result_type                  result_type;
    typedef typename JFunction1D_t::value_type                   value_type;

    typedef typename JMultiMap_t::key_type                       key_type;
    typedef typename JMultiMap_t::mapped_type                    mapped_type;

    typedef typename JMultiMap_t::iterator                       iterator;
    typedef typename JMultiMap_t::reverse_iterator               reverse_iterator;
    typedef typename JMultiMap_t::const_iterator                 const_iterator;
    typedef typename JMultiMap_t::const_reverse_iterator         const_reverse_iterator;

    typedef typename JMultiMap_t::super_iterator                 super_iterator;
    typedef typename JMultiMap_t::super_const_iterator           super_const_iterator;

    typedef JFunctionTransformerInterface<NUMBER_OF_DIMENSIONS, argument_type>      transformer_type;
    typedef JFunction1D_t                                                           function_type;


    /**
     * Default constructor.
     */
    JTransformableMultiFunction() :
      JMultiMap_t(),
      transformer(new transformer_type())
    {}


    /**
     * Multi-dimensional interpolation method call.
     *
     * \param  x      comma seperated argument list
     * \return        function value
     */
    result_type operator()(const argument_type& x, ...) const
    {
      va_start(this->ap, x);

      typename JArray<NUMBER_OF_DIMENSIONS + 1, argument_type>::iterator p = this->buffer.begin();

      *p = x;                                  // first argument

      while (++p != this->buffer.end())
        *p = va_arg(this->ap, argument_type);  // next argument

      va_end(this->ap);

      this->buffer[NUMBER_OF_DIMENSIONS] = transformer->putXn(this->buffer, this->buffer[NUMBER_OF_DIMENSIONS]);

      const result_type y = this->evaluate(this->buffer.data());
      const double      W = transformer->getWeight(this->buffer);

      return y * W;
    }


    /**
     * Application of weight function and coordinate transformation.
     *
     * \param  __transformer      JFunctionTransformerInterface object
     */
    template<class JFunctionTransformer_t>
    void transform(const JFunctionTransformer_t& __transformer)
    {
      typedef JArray<NUMBER_OF_DIMENSIONS, argument_type>  array_type;

      for (super_iterator i = this->super_begin(); i != this->super_end(); ++i) {

	const array_type  array    = (*i).getKey();
	function_type&    function = (*i).getValue();

	const JFunctionGetTransformer<NUMBER_OF_DIMENSIONS, value_type> get( *transformer, array);
	const JFunctionPutTransformer<NUMBER_OF_DIMENSIONS, value_type> put(__transformer, array);
	
	function.transform(get);
	function.transform(put);
      }

      transformer.reset(__transformer.clone());
    }


    /**
     * Read from input.
     *
     * \param  in     JReader
     * \return        JReader
     */
    virtual JReader& read(JReader& in)
    {
      in >> (JMultiMap_t&) *this;

      return transformer->read(in);
    }


    /**
     * Write to output.
     *
     * \param  out     JWriter
     * \return         JWriter
     */
    virtual JWriter& write(JWriter& out) const
    {
      out << (JMultiMap_t&) *this;

      return transformer->write(out);
    }


    /**
     * Read JTransformableMultiFunction from input.
     *
     * \param  in       JReader
     * \param  buffer   JTransformableMultiFunction
     * \return          JReader
     */
    friend inline JReader& operator>>(JReader& in, JTransformableMultiFunction& buffer)
    {
      return buffer.read(in);
    }


    /**
     * Write JTransformableMultiFunction to output.
     *
     * \param  out      JWriter
     * \param  buffer   JTransformableMultiFunction
     * \return          JWriter
     */
    friend inline JWriter& operator<<(JWriter& out, const JTransformableMultiFunction& buffer)
    {
      return buffer.write(out);
    }


    std::auto_ptr<transformer_type> transformer;
  };
}

#endif