interval.h

#pragma once
 
#include <cgv/type/traits/min.h>
#include <cgv/type/traits/max.h>
#include <iostream>
#include <algorithm>
 
namespace cgv {
        namespace math {
 
template <typename T>
class interval
{
protected:
        T lb, ub;
public:
        interval() : lb(1), ub(0) {}
        interval(bool) : lb(type::traits::min_fct<T>::get_value()), ub(type::traits::max_fct<T>::get_value()) {}
        interval(const interval<T>& I) : lb(I.get_lower_bound()), ub(I.get_upper_bound()) {}
        interval(const T& _lb, const T& _ub) : lb(std::min(_lb,_ub)), ub(std::max(_lb,_ub)) {}
        bool is_empty() const { return ub < lb; }
        void clear() { ub = 0; lb = 1; }
        bool contains(const T& v) const { return lb <= v && v <= ub; }
        void set_lower_bound(const T& _lb) { lb = _lb; }
        void set_upper_bound(const T& _ub) { ub = _ub; }
        const T& get_lower_bound() const { return lb; }
        const T& get_upper_bound() const { return ub; }
        T get_center() const { return (lb+ub)/2; }
        T get_size() const { return is_empty() ? 0 : ub-lb; }
        interval<T>& intersect(const interval<T>& I) {
                lb = std::max(lb, I.get_lower_bound());
                ub = std::min(ub, I.get_upper_bound());
                return *this;
        }
        interval<T> intersection(const interval<T>& I) const { interval J(*this); return J.intersect(I); } 
        interval<T>& extend(const T& v) { 
                if (is_empty()) 
                        lb=ub=v; 
                else { 
                        lb = std::min(lb,v);
                        ub = std::max(ub,v);
                }
                return *this;
        }
        interval<T>& extension(const T& v) { interval<T> I(*this); return I.extend(v); }
        interval<T>& extend(const interval<T>& I) { 
                if (!I.is_empty()) {
                        extend(I.get_lower_bound()); 
                        extend(I.get_upper_bound()); 
                }
                return *this;
        }
        interval<T>& extend(const T& v0, const T& v1, const T& v2, const T& v3) {
                extend(v0); extend(v1); extend(v2); extend(v3);
                return *this;
        }
        interval<T>& extension(const interval<T>& J) { interval<T> I(*this); return I.extend(J); }
 
        interval<T>& operator *= (const T& s) { lb *= s; ub *= s; return *this; } 
        interval<T> operator * (const T& s) const { interval I(*this); return I *= s; } 
        interval<T>& operator /= (const T& s) { return *this *= 1/s; } 
        interval<T> operator / (const T& s) const { return *this * (1/s); } 
        interval<T>& operator += (const T& s) { lb += s; ub += s; return *this; } 
        interval<T> operator + (const T& s) const { interval I(*this); return I += s; } 
        interval<T>& operator -= (const T& s) { lb -= s; ub -= s; return *this; } 
        interval<T> operator - (const T& s) const { interval I(*this); return I -= s; } 
 
        interval<T>& operator *= (const interval<T>& I) { 
                if (is_empty() || I.is_empty())
                        *this = interval<T>();
                else 
                        extend(lb*I.get_lower_bound(), lb*I.get_upper_bound(),
                                   ub*I.get_lower_bound(), ub*I.get_upper_bound());
                return *this; 
        } 
        interval<T> operator * (const interval<T>& J) const { interval I(*this); return I *= J; } 
        interval<T>& operator /= (const interval<T>& I) { 
                if (is_empty() || I.is_empty())
                        *this = interval<T>();
                else {
                        if (I.contains(0))
                                *this = interval<T>(true);
                        else
                                extend(lb/I.get_lower_bound(), lb/I.get_upper_bound(),
                                           ub/I.get_lower_bound(), ub/I.get_upper_bound());
                }
                return *this;
        } 
        interval<T> operator / (const interval<T>& J) const { interval I(*this); return I /= J; } 
        interval<T> operator - () const { if (is_empty()) return *this; else return interval<T>(-ub,-lb); }
        interval<T>& operator += (const interval<T>& I) { 
                if (is_empty() || I.is_empty())
                        *this = interval<T>();
                else {
                        extend(lb+I.get_lower_bound(), lb+I.get_upper_bound(),
                                   ub+I.get_lower_bound(), ub+I.get_upper_bound());
                }
                return *this; 
        } 
        interval<T> operator + (const interval<T>& J) const { interval I(*this); return I += J; } 
        interval<T>& operator -= (const interval<T>& I) { return *this += -I; }
        interval<T> operator - (const interval<T>& J) const { interval I(*this); return I -= J; } 
 
        bool operator == (const interval<T>& I) const {
                return is_empty() && I.is_empty() ||
                        lb == I.get_lower_bound() && ub == I.get_upper_bound();
        }
        bool operator != (const interval<T>& I) const {
                return !(*this == I);
        }
        bool operator < (const interval<T>& I) const {
                return !is_empty() && !I.is_empty() && ub < I.get_lower_bound();
        }
        bool operator > (const interval<T>& I) const {
                return !is_empty() && !I.is_empty() && lb > I.get_upper_bound();
        }
 
};
 
template <typename T> inline interval<T> operator + (const T& s, const interval<T>& I) { return  I+s; }
template <typename T> inline interval<T> operator - (const T& s, const interval<T>& I) { return -I+s; }
template <typename T> inline interval<T> operator * (const T& s, const interval<T>& I) { return  I*s; }
template <typename T> inline interval<T> operator / (const T& s, const interval<T>& I) { return interval<T>(s,s)/I; }
 
template <typename T>
std::ostream& operator << (std::ostream& os, const interval<T>& I) {
        return os << '[' << I.get_lower_bound() << ',' << I.get_upper_bound() << ']';
}
 
        }
}