DensityScanND.hh

Go to the documentation of this file.

#ifndef NPSTAT_DENSITYSCANND_HH_
#define NPSTAT_DENSITYSCANND_HH_
 
/*!
// \file DensityScanND.hh
//
// \brief Fill multivariate arrays with multivariate density values
//
// Author: I. Volobouev
//
// December 2011
*/
 
#include <vector>
#include <cassert>
#include <stdexcept>
#include <memory>
#include <climits>
 
#include "npstat/stat/AbsDistributionND.hh"
 
namespace npstat {
    /**
    // The IndexMapper template parameter is a functor which takes
    // a double argument (representing unsigned int) and returns
    // a double (or something which can be automatically converted
    // to double). IndexMapper can, for example, be set to LinearMapper1d
    // or CircularMapper1d.
    */
    template<class IndexMapper>
    class DensityScanND
    {
    public:
        /**
        // A mapper for each coordinate in the "maps" constructor argument
        // will convert the array index into a proper argument for the scanned
        // density.
        //
        // This functor will NOT make a copy of the "fcn" parameter. This
        // parameter will be used by reference only (aliased). It is up to
        // the user of this class to ensure proper lifetime of this object.
        */
        inline DensityScanND(const AbsDistributionND& fcn,
                             const std::vector<IndexMapper>& maps,
                             double normfactor=1.0)
            : fcn_(fcn), mapping_(new std::vector<IndexMapper>(maps)),
              normfactor_(normfactor), dim_(fcn.dim())
        {
            if (!(dim_ && dim_ == mapping_->size())) throw std::invalid_argument(
                "In npstat::DensityScanND constructor: incompatible arguments");
            if (dim_ > CHAR_BIT*sizeof(unsigned long)) throw std::invalid_argument(
                "In npstat::DensityScanND constructor: dimensionality is too high");
        }
 
        inline void getCoords(const unsigned* index, const unsigned indexLen,
                              double* coords) const
        {
            if (dim_ != indexLen) throw std::invalid_argument(
                "In npstat::DensityScanND::getCoords: "
                "incompatible index dimensionality");
            assert(index);
            assert(coords);
            const IndexMapper* mappers = &(*mapping_)[0];
            for (unsigned i=0; i<indexLen; ++i)
                coords[i] = mappers[i](static_cast<double>(index[i]));
        }
 
        inline double operator()(const unsigned* index, unsigned indexLen) const
        {
            double x[CHAR_BIT*sizeof(unsigned long)];
            getCoords(index, indexLen, x);
            return normfactor_*fcn_.density(x, indexLen);
        }
 
    private:
        DensityScanND();
 
        // This functor will be recursively copied by the
        // scanning code, so it must be lightweight. This
        // is what std::shared_ptr is here for.
        const AbsDistributionND& fcn_;
        std::shared_ptr<std::vector<IndexMapper> > mapping_;
        double normfactor_;
        unsigned dim_;
 
#ifdef SWIG
    public:
        inline void mapIndex(const unsigned* index, unsigned indexLen,
                             double* coords, unsigned lenCoords)
        {
            if (dim_ != lenCoords) throw std::invalid_argument(
                "In npstat::DensityScanND::mapIndex: "
                "incompatible coordinate dimensionality");
            this->getCoords(index, indexLen, coords);
        }
#endif // SWIG
    };
 
    // Helper function for creating DensityScanND objects
    template<class IndexMapper>
    inline DensityScanND<IndexMapper> makeDensityScanND(
        const AbsDistributionND& fcn, const std::vector<IndexMapper>& axes,
        double normfactor=1.0)
    {
        return DensityScanND<IndexMapper>(fcn, axes, normfactor);
    }
}
 
#endif // NPSTAT_DENSITYSCANND_HH_