/*---------------------------------------------------------------------------*\

    geophase: A C++ template library for geometrical operations on non-convex
    polyhedra with non-planar faces
     
    Copyright (C) 2019 Tomislav Maric, TU Darmstadt 

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.

Description

    Perturb the geometry by a given tolerance using uniform random distribution
    within (-tolerance, tolerance).

Author
    Tomislav Maric maric@mma.tu-darmstadt.de

\*---------------------------------------------------------------------------*/

#ifndef Perturb_H 
#define Perturb_H 

#include <numeric>
#include <random>
#include "tag.hpp"
#include <iostream>
#include "tolerance.hpp"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace geophase { 

    static std::random_device randDev;  
    static std::mt19937 randGen(randDev()); 
    static const std::array<int, 2> signs = {-1, 1};
    static std::uniform_int_distribution<> randSign(0, 1);

    template<typename Vector> 
    enable_if_type_tagged_t<Vector, vector_tag>
    ulp_perturb(Vector& vec, std::size_t N)
    {
        std::uniform_int_distribution<> randInterv(1, N);
        for (auto& cmpt : vec)
            cmpt += signs[randSign(randGen)] * randInterv(randGen) * ulp(cmpt); 
    }

    template<typename Vector> 
    enable_if_type_tagged_t<Vector, vector_tag, Vector>
    ulp_perturbed(const Vector& vec, std::size_t N)
    {
        Vector copy(vec); 

        ulp_perturb(copy, N); 

        return copy; 
    }

    template<typename Vector> 
    enable_if_type_tagged_t<Vector, vector_tag>
    tolerance_perturb(Vector& vec, typename Vector::value_type epsilon) 
    {
        for (auto& cmpt : vec)
            cmpt = cmpt + signs[randSign(randGen)] * epsilon; 
    }

    template<typename Vector> 
    enable_if_type_tagged_t<Vector, vector_tag, Vector>
    tolerance_perturbed(const Vector& vec, typename Vector::value_type epsilon) 
    {
        Vector copy(vec); 

        tolerance_perturb(copy, epsilon); 

        return copy; 
    }

    //template<typename Vector> 
    //enable_if_type_tagged_t<Vector, vector_tag, Vector>
    //perturbation_vector(double tau)
    //{
        //static std::random_device randDevice;
        //static std::default_random_engine randEngine(randDevice());
        //static std::uniform_real_distribution<double> randDist(-tau, tau); 
        //return Vector(distribution(tau), distribution(tau), distribution(tau)); 
    //}

    //auto empty_ref_lambda = [](auto& x) {};

    //template
    //<
        //typename Point, 
        //typename UnaryFunction = decltype(empty_ref_lambda)
    //>
    //typename std::enable_if
    //<
        //Traits::tag_enabled<Point, point_tag>::value, 
        //void
    //>::type
    //perturb(
        //Point& g, 
        //double epsilon, 
        //UnaryFunction f = empty_ref_lambda
    //)
    //{
        //auto perturbationVector = perturbation_vector<Point>(epsilon);  

        //f(perturbationVector); 

        //g += perturbation_vector<Point>(epsilon); 
    //}

    //template
    //<
        //typename PointCollection, 
        //typename UnaryFunction = decltype(empty_ref_lambda)
    //>
    //typename std::enable_if
    //<
        //Traits::tag_enabled<PointCollection, point_collection_tag>::value ||
        //Traits::tag_enabled<PointCollection, triangle_tag>::value, 
        //void 
    //>::type
    //perturb(
        //PointCollection& polygon, 
        //double epsilon, 
        //UnaryFunction f = empty_ref_lambda
    //) 
    //{
        //for (auto& point : polygon)
        //{
            //perturb(point, epsilon, f); 
        //}
    //} 

    //template
    //<
        //typename PolygonCollection,
        //typename UnaryFunction = decltype(empty_ref_lambda)
    //>
    //typename std::enable_if
    //<
        //Traits::tag_enabled<PolygonCollection, polygon_collection_tag>::value ||
        //Traits::tag_enabled<PolygonCollection, tetrahedron_tag>::value,
        //void 
    //>::type
    //perturb(PolygonCollection& polyhedron, double epsilon, UnaryFunction f = empty_ref_lambda) 
    //{
        //for (auto& polygon : polyhedron)
        //{
            //perturb(polygon, epsilon, f);
        //}
    //} 

    //template
    //<
        //typename Triangulation,
        //typename UnaryFunction = decltype(empty_ref_lambda)
    //>
    //typename std::enable_if
    //<
        //Traits::tag_enabled<Triangulation, tetrahedron_collection_tag>::value,
        //void 
    //>::type
    //perturb(Triangulation& triangulation, double epsilon, UnaryFunction f = empty_ref_lambda) 
    //{
        //for (auto& tetrahedron : triangulation)
        //{
            //perturb(tetrahedron, epsilon, f);
        //}
    //} 

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace geophase 

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //
