//#####################################################################
// Copyright 2002-2009, Robert Bridson, Ronald Fedkiw, Eran Guendelman, Geoffrey Irving, Sergey Koltakov, Nipun Kwatra, Neil Molino, Igor Neverov, Duc Nguyen, Avi Robinson-Mosher,
//     Craig Schroeder, Andrew Selle, Tamar Shinar, Jonathan Su, Jerry Talton, Joseph Teran.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
// Class RANGE
//#####################################################################
#ifndef __RANGE__
#define __RANGE__

#include <Core/Math_Tools/FIXED_NUMBER.h>
#include <Core/Math_Tools/minmax.h>
#include <Core/Matrices/SYMMETRIC_MATRIX.h>
#include <Core/Vectors/SCALAR_POLICY.h>
#include <Core/Vectors/VECTOR.h>
#include <cassert>
#include <cfloat>
#include <iostream>
#include <limits>
#include <string>
namespace PhysBAM{

template<class TV> class RANGE;
template<class TV> struct IS_SCALAR_BLOCK<RANGE<TV> >:public IS_SCALAR_BLOCK<TV>{};
template<class TV,class RW> struct IS_BINARY_IO_SAFE<RANGE<TV>,RW> {static const bool value=false;}; // required since memory format differs from disk format

template<class TV>
class RANGE
{
    typedef typename TV::SCALAR T;
public:
    template<class T2> struct REBIND{typedef RANGE<T2> TYPE;};
    typedef T SCALAR;
    typedef TV VECTOR_T;
    enum WORKAROUND {d=TV::m};

    union
    {
        struct
        {
            TV min_corner,max_corner;
        };
        TV corner[2];
    };

    // Default to an empty box.
    RANGE()
        :min_corner(TV::Constant_Vector(std::numeric_limits<T>::max())),max_corner(-TV::Constant_Vector(std::numeric_limits<T>::max()))
    {}

    RANGE(const TV& minimum_corner,const TV& maximum_corner)
        :min_corner(minimum_corner),max_corner(maximum_corner)
    {}

    template<class T2> explicit RANGE(const RANGE<T2>& box)
        :min_corner(TV(box.min_corner)),max_corner(TV(box.max_corner))
    {}

    explicit RANGE(const TV& point)
        :min_corner(point),max_corner(point)
    {}

    RANGE(const RANGE& r)
        :min_corner(r.min_corner),max_corner(r.max_corner)
    {}

    RANGE(const RANGE<TV>& range,const FRAME<VECTOR<T,1> >& frame) // allow 1d ranges to be used as oriented boxes
        :min_corner(frame*range.min_corner),max_corner(frame*range.max_corner)
    {}

    ~RANGE()
    {
        min_corner.~TV();
        max_corner.~TV();
    }

    RANGE& operator=(const RANGE& r)
    {
        min_corner=r.min_corner;
        max_corner=r.max_corner;
        return *this;
    }

    static RANGE<TV> Unit_Box()
    {return RANGE<TV>(TV(),TV::All_Ones_Vector());}

    static RANGE<TV> Centered_Box()
    {return RANGE<TV>(-TV::All_Ones_Vector(),TV::All_Ones_Vector());}

    static RANGE<TV> Empty_Box()
    {return RANGE<TV>();}

    static RANGE<TV> Full_Box()
    {return RANGE<TV>(TV::Constant_Vector(-std::numeric_limits<T>::max()),TV::Constant_Vector(std::numeric_limits<T>::max()));}

    bool Empty() const
    {return !min_corner.All_Less_Equal(max_corner);}

    bool Empty_Half_Open() const
    {return !min_corner.All_Less(max_corner);}

    RANGE<TV> To_Half_Opened() const
    {
        static_assert((is_same<typename TV::ELEMENT,int>::value));
        return RANGE<TV>(min_corner,max_corner+1);
    }

    RANGE<TV> To_Closed() const
    {
        static_assert((is_same<typename TV::ELEMENT,int>::value));
        return RANGE<TV>(min_corner,max_corner-1);
    }

    bool operator==(const RANGE<TV>& r) const
    {return min_corner==r.min_corner && max_corner==r.max_corner;}

    bool operator!=(const RANGE<TV>& r) const
    {return !(*this==r);}

    RANGE<TV> operator-() const
    {return RANGE<TV>(-max_corner,-min_corner);}

    RANGE<TV>& operator+=(const RANGE<TV>& r)
    {min_corner+=r.min_corner;max_corner+=r.max_corner;return *this;}

    RANGE<TV>& operator-=(const RANGE<TV>& r)
    {min_corner-=r.max_corner;max_corner-=r.min_corner;return *this;}

    RANGE<TV> operator+(const RANGE<TV>& r) const
    {return RANGE<TV>(min_corner+r.min_corner,max_corner+r.max_corner);}

    RANGE<TV> operator-(const RANGE<TV>& r) const
    {return RANGE<TV>(min_corner-r.max_corner,max_corner-r.min_corner);}

    RANGE<TV>& operator+=(const TV& v)
    {min_corner+=v;max_corner+=v;return *this;}

    RANGE<TV>& operator-=(const TV& v)
    {min_corner-=v;max_corner-=v;return *this;}

    RANGE<TV> operator*(const T a) const
    {return a>=0?RANGE<TV>(min_corner*a,max_corner*a):RANGE<TV>(max_corner*a,min_corner*a);}

    RANGE<TV> operator*(const TV& a) const
    {return RANGE<TV>(min_corner*a,max_corner*a);}

    RANGE<TV>& operator*=(const T a)
    {return *this=*this*a;}

    RANGE<TV>& operator*=(const TV& a)
    {return *this=*this*a;}

    RANGE<TV> operator/(const T a) const
    {assert(a!=0);return RANGE<TV>(min_corner/a,max_corner/a);}

    RANGE<TV>& operator/=(const T a)
    {return *this=*this/a;}

    TV Edge_Lengths() const
    {return max_corner-min_corner;}

    TV Center() const
    {return (min_corner+max_corner)/(T)2;} // Be nice to T=int

    TV Minimum_Corner() const
    {return min_corner;}

    TV Maximum_Corner() const
    {return max_corner;}

    void Corners(VECTOR<TV,1<<TV::m>& corners) const
    {TV s[2]={min_corner,max_corner};for(int i=0;i<(1<<TV::m);i++) for(int j=0;j<TV::m;j++) corners(i)(j)=s[(i>>j)&1](j);}

    T Size() const // assumes nonnegative entries
    {return Edge_Lengths().Product();}

    T Robust_Size() const
    {return Empty()?(T)0:Size();}

    void Reset_Bounds(const TV& point)
    {min_corner=max_corner=point;}

    void Enlarge_To_Include_Point(const TV& point)
    {min_corner=min_corner.Componentwise_Min(point);max_corner=max_corner.Componentwise_Max(point);}

    void Enlarge_Nonempty_Box_To_Include_Point(const TV& point)
    {assert(!Empty());for(int i=0;i<d;i++) if(point(i)<min_corner(i)) min_corner(i)=point(i);else if(point(i)>max_corner(i)) max_corner(i)=point(i);}

    void Enlarge_Nonempty_Box_To_Include_Points(const TV& p1,const TV& p2)
    {Enlarge_Nonempty_Box_To_Include_Point(p1);Enlarge_Nonempty_Box_To_Include_Point(p2);}

    void Enlarge_Nonempty_Box_To_Include_Points(const TV& p1,const TV& p2,const TV& p3)
    {Enlarge_Nonempty_Box_To_Include_Point(p1);Enlarge_Nonempty_Box_To_Include_Point(p2);Enlarge_Nonempty_Box_To_Include_Point(p3);}

    template<class T_ARRAY>
    void Enlarge_Nonempty_Box_To_Include_Points(const T_ARRAY& points)
    {static_assert((is_same<typename T_ARRAY::ELEMENT,TV>::value));
    for(int i=0;i<points.Size();i++) Enlarge_Nonempty_Box_To_Include_Point(points(i));}

    template<class T_ARRAY>
    void Enlarge_To_Include_Points(const T_ARRAY& points)
    {static_assert((is_same<typename T_ARRAY::ELEMENT,TV>::value));
    if(!points.Size()) return;
    Enlarge_To_Include_Point(points(0));
    for(int i=1;i<points.Size();i++) Enlarge_Nonempty_Box_To_Include_Point(points(i));}

    void Enlarge_To_Include_Box(const RANGE<TV>& box)
    {min_corner=min_corner.Componentwise_Min(box.min_corner);max_corner=max_corner.Componentwise_Max(box.max_corner);}

    void Change_Size(const T delta)
    {min_corner-=delta;max_corner+=delta;}

    void Change_Size(const TV& delta)
    {min_corner-=delta;max_corner+=delta;}

    void Shift(const TV& delta)
    {min_corner+=delta;max_corner+=delta;}

    void Shift(const T delta)
    {min_corner+=delta;max_corner+=delta;}

    RANGE<TV> Thickened(const T thickness_over_two) const
    {return RANGE<TV>(min_corner-thickness_over_two,max_corner+thickness_over_two);}

    RANGE<TV> Thickened(const TV& thickness_over_two) const
    {return RANGE<TV>(min_corner-thickness_over_two,max_corner+thickness_over_two);}

    static RANGE<TV> Combine(const RANGE<TV>& box1,const RANGE<TV>& box2)
    {return RANGE<TV>(box1.min_corner.Componentwise_Min(box2.min_corner),box1.max_corner.Componentwise_Max(box2.max_corner));}

    static RANGE<TV> Intersect(const RANGE<TV>& box1,const RANGE<TV>& box2) // assumes nonnegative entries
    {return RANGE<TV>(box1.min_corner.Componentwise_Max(box2.min_corner),box1.max_corner.Componentwise_Min(box2.max_corner));}

    RANGE<TV> Unite(const RANGE<TV>& box) const
    {return Combine(*this,box);}

    RANGE<TV> Intersect(const RANGE<TV>& box) const // assumes nonnegative entries
    {return Intersect(*this,box);}

    void Scale_About_Center(const T factor)
    {Scale_About_Point(Center(),factor);}

    void Scale_About_Center(const TV& factor)
    {Scale_About_Point(Center(),factor);}

    void Scale_About_Point(const TV& pt, const T factor)
    {TV length_over_two=factor*(T).5*(max_corner-min_corner);min_corner=pt-length_over_two;max_corner=pt+length_over_two;}

    void Scale_About_Point(const TV& pt, const TV& factor)
    {TV length_over_two=factor*(T).5*(max_corner-min_corner);min_corner=pt-length_over_two;max_corner=pt+length_over_two;}

    RANGE<TV> operator+(const TV& shift) const
    {return RANGE<TV>(min_corner+shift,max_corner+shift);}

    RANGE<TV> operator-(const TV& shift) const
    {return RANGE<TV>(min_corner-shift,max_corner-shift);}

    RANGE<TV> Translated(const TV& shift) const
    {return *this+shift;}

    bool Lazy_Inside(const TV& location) const
    {return location.All_Greater_Equal(min_corner) && location.All_Less_Equal(max_corner);}

    bool Lazy_Inside_Half_Open(const TV& location) const
    {return location.All_Greater_Equal(min_corner) && location.All_Less(max_corner);}

    bool Lazy_Inside(const RANGE<TV>& other) const
    {return other.min_corner.All_Greater_Equal(min_corner) && other.max_corner.All_Less_Equal(max_corner);}

    bool Inside(const TV& location,const T thickness_over_two) const
    {return Thickened(-thickness_over_two).Lazy_Inside(location);}

    bool Inside(const TV& location,const FIXED_NUMBER<T,0> thickness_over_two) const
    {return Lazy_Inside(location);}

    bool Lazy_Outside(const TV& location) const
    {return !Lazy_Inside(location);}

    bool Lazy_Outside_Half_Open(const TV& location) const
    {return !Lazy_Inside_Half_Open(location);}

    bool Outside(const TV& location,const T thickness_over_two) const
    {return Thickened(thickness_over_two).Lazy_Outside(location);}

    bool Outside(const TV& location,const FIXED_NUMBER<T,0> thickness_over_two) const
    {return Lazy_Outside(location);}

    bool Boundary(const TV& location,const T thickness_over_two) const
    {bool strict_inside=location.All_Greater(min_corner+thickness_over_two) && location.All_Less(max_corner-thickness_over_two);
    return !strict_inside && !Outside(location,thickness_over_two);}

    TV Clamp(const TV& location) const
    {return clamp(location,min_corner,max_corner);}

    TV Clamp_Half_Open(const TV& location) const
    {return clamp(location,min_corner,max_corner-1);}

    void Enlarge_By_Sign(const TV& v)
    {for(int i=0;i<d;i++) if(v(i)>0) max_corner(i)+=v(i);else min_corner(i)+=v(i);}

    TV Point_From_Normalized_Coordinates(const TV& weights) const
    {return min_corner+weights*(max_corner-min_corner);}

    bool Contains(const RANGE<TV>& box) const
    {return min_corner.All_Less_Equal(box.min_corner) && max_corner.All_Greater_Equal(box.max_corner);}

    bool Lazy_Intersection(const RANGE<TV>& box) const
    {return min_corner.All_Less_Equal(box.max_corner) && max_corner.All_Greater_Equal(box.min_corner);}

    bool Lazy_Intersection_Half_Open(const RANGE<TV>& box) const
    {return min_corner.All_Less(box.max_corner) && max_corner.All_Greater(box.min_corner);}

    bool Intersection(const RANGE<TV>& box,const T thickness_over_two) const
    {return Thickened(thickness_over_two).Lazy_Intersection(box);}

    bool Intersection(const RANGE<TV>& box,const FIXED_NUMBER<T,0> thickness_over_two) const
    {return Lazy_Intersection(box);}

    bool Intersection(const RANGE<TV>& box) const
    {return Lazy_Intersection(box);}

    T Intersection_Area(const RANGE<TV>& box) const
    {return Intersect(*this,box).Robust_Size();}

    void Project_Points_Onto_Line(const TV& direction,T& line_min,T& line_max) const
    {line_min=line_max=direction.Dot(min_corner);TV e=direction*(max_corner-min_corner);
    for(int i=0;i<d;i++) if(e(i)>0) line_max+=e(i);else line_min+=e(i);}

    static RANGE<TV> Bounding_Box(const TV& p1,const TV& p2)
    {RANGE<TV> box(p1);box.Enlarge_Nonempty_Box_To_Include_Point(p2);return box;}

    static RANGE<TV> Bounding_Box(const TV& p1,const TV& p2,const TV& p3)
    {RANGE<TV> box(p1);box.Enlarge_Nonempty_Box_To_Include_Points(p2,p3);return box;}

    static RANGE<TV> Bounding_Box(const TV& p1,const TV& p2,const TV& p3,const TV& p4)
    {RANGE<TV> box(p1);box.Enlarge_Nonempty_Box_To_Include_Points(p2,p3,p4);return box;}

    template<class T_ARRAY>
    static RANGE<TV> Bounding_Box(const T_ARRAY& points)
    {static_assert((is_same<typename T_ARRAY::ELEMENT,TV>::value));
    if(!points.Size()) return Empty_Box();
    RANGE<TV> box(points(0));
    for(int i=1;i<points.Size();i++) box.Enlarge_Nonempty_Box_To_Include_Point(points(i));
    return box;}

    RANGE<VECTOR<T,d-1> > Remove_Dimension(int dimension) const
    {return RANGE<VECTOR<T,d-1> >(min_corner.Remove_Index(dimension),max_corner.Remove_Index(dimension));}

    RANGE<VECTOR<T,1> > Dimension_Range(int dimension) const
    {return RANGE<VECTOR<T,1> >(VECTOR<T,1>(min_corner(dimension)),VECTOR<T,1>(max_corner(dimension)));}

    const RANGE<TV>& Bounding_Box() const
    {return *this;}

    RANGE<TV> Axis_Aligned_Bounding_Box() const
    {return *this;}

//#####################################################################
    TV Normal(const int aggregate) const;
    TV Normal(const TV& X) const;
    TV Surface(const TV& X) const;
    T Signed_Distance(const TV& X) const;
    SYMMETRIC_MATRIX<T,TV::m> Hessian(const TV& X) const;
    static std::string Name();
    template<class RW> void Read(UNTYPED_ISTREAM<RW>& input);
    template<class RW> void Write(UNTYPED_OSTREAM<RW>& output) const;
    void Fill_Uniform(RANDOM_NUMBERS<T>& rand,const RANGE<TV>& box);
//#####################################################################
};
template<class TV>
inline RANGE<TV> operator+(const TV& a,const RANGE<TV>& b)
{return RANGE<TV>(a+b.min_corner,a+b.max_corner);}

template<class TV>
inline RANGE<TV> operator-(const TV& a,const RANGE<TV>& b)
{return RANGE<TV>(a-b.max_corner,a-b.min_corner);}

template<class TV> inline RANGE<TV> operator*(const typename TV::SCALAR a,const RANGE<TV>& box)
{return box*a;}

template<class TV>
inline std::ostream& operator<<(std::ostream& output,const RANGE<TV>& box)
{output<<"("<<box.min_corner<<" "<<box.max_corner<<")";return output;}

template<class TV>
inline std::istream& operator>>(std::istream& input,RANGE<TV>& box)
{Ignore(input,'(');input>>box.min_corner>>box.max_corner;Ignore(input,')');return input;}
template<class TV> struct HASH_REDUCE<RANGE<TV> >
{static int H(const RANGE<TV>& key){return int_hash(HASH_REDUCE<TV>::H(key.min_corner),HASH_REDUCE<TV>::H(key.max_corner));}};
template<class T,class TV> RANGE<TV>
Scaled(const RANGE<TV>& b,const TV& center,T scale)
{
    return (b-center)*scale+center;
}
}
#include <Core/Math_Tools/RANGE_ITERATOR.h>
#endif
