#include "low_discrepancy.hpp"
#include <vector>
#include <gsl/gsl_cdf.h>
#include "var_alea.hpp"


double frac_part(double x);

//fonctions de compositions de monte_carlo.hpp

template <class _Result>
struct generator 
{
	typedef _Result result_type;
};

template <typename Fct, typename VA>
struct compose_t : public generator< typename Fct::result_type >
{
	compose_t(Fct f, VA X) : f(f), X(X) {};
	typename Fct::result_type operator()() { 
		return f(X()); 
	};
	private:
		Fct f; VA X;
};

template <typename Fct, typename VA>
inline compose_t<Fct, VA>
compose(Fct f, VA X) {
    return compose_t<Fct, VA>(f, X);
};

//

template <typename L>
std::vector<double> monte_carlo(int n, L X)
{
	std::vector<double> result(3,0);
	double x;
	for (int j = 0; j < n; j++) {
		x = X();
		result[0] += x;
		result[1] += x*x;
	}
	result[0] /= (double) n;
	result[1] = (result[1] - n*result[0]*result[0])/(double)(n-1);
	result[2] = 1.96*sqrt(result[1]/(double) n);
	return result;
};

template <typename L, typename Fct>
std::vector<double> monte_carlo(int n, Fct f, L X)
{
	return monte_carlo(n, compose(f, X)); 
};

struct quasi_gaussian : public var_alea<std::vector<double> >
{
    quasi_gaussian (int d, double mean = 0, double std =1)
        : d(d), mean(mean), std(std), sob(d) {};
    std::vector<double> operator ()(){
        std::vector<double> result = sob();
        for(int i=0;i<d; i++){
            result[i] = mean + gsl_cdf_gaussian_Pinv(result[i], std);
        }
        return value = result;
    }
    private: 
    int d;
    double mean, std;
    sobol sob;
};

struct gaussian_d : public var_alea<std::vector<double> >
{
    gaussian_d(int d, double mean = 0, double std =1)
        : d(d), mean(mean), std(std), G(mean, std) {};
    std::vector<double> operator ()(){
        std::vector<double> result(d);
        for(int i=0;i<d; i++){
            result[i] = G();
        }
        return value = result;
    }
    private: 
    int d;
    double mean, std;
    gaussian G;
};