//
// This file is part of MorphoGraphX - http://www.MorphoGraphX.org
// Copyright (C) 2012-2015 Richard S. Smith and collaborators.
//
// If you use MorphoGraphX in your work, please cite:
//   http://dx.doi.org/10.7554/eLife.05864
//
// MorphoGraphX is free software, and is licensed under under the terms of the 
// GNU General (GPL) Public License version 2.0, http://www.gnu.org/licenses.
//
#ifndef BOUNDINGBOX_HPP
#define BOUNDINGBOX_HPP

#include <Config.hpp>

#include <Vector.hpp>

#include <iostream>
#include <limits>

namespace mgx 
{
  template <size_t N, typename T> class BoundingBox {
  public:
    typedef Vector<N, T> Point;
    typedef std::numeric_limits<T> limits;
  
    CU_HOST_DEVICE
    BoundingBox() { reset(); }
  
    CU_HOST_DEVICE
    BoundingBox(const Point& vmin, const Point& vmax)
    {
      pts[0] = vmin;
      pts[1] = vmax;
    }
  
    CU_HOST_DEVICE
    BoundingBox(const BoundingBox& copy)
    {
      pts[0] = copy.pts[0];
      pts[1] = copy.pts[1];
    }
  
    CU_HOST_DEVICE
    BoundingBox(const Point& p)
    {
      pts[0] = pts[1] = p;
      if(limits::is_integer)
        pts[1] += T(1);       // What does this do?
    }
  
    CU_HOST_DEVICE
    void reset()
    {
      T v = limits::max();
      pts[0] = Point(v, v, v);
      if(limits::is_signed)
        pts[1] = -pts[0];
      else
        pts[1] = Point(0, 0, 0);
    }
  
    CU_HOST_DEVICE
    bool empty() const
    {
      for(size_t i = 0; i < N; ++i)
        if(pts[0][i] >= pts[1][i])
          return true;
      return false;
    }
  
    CU_HOST_DEVICE
    operator bool() const { return not empty(); }
  
    CU_HOST_DEVICE
    bool operator==(const BoundingBox& other) const 
   {
      return pts[0] == other.pts[0] and pts[1] == other.pts[1];
    }
  
    CU_HOST_DEVICE
    bool operator!=(const BoundingBox& other) const 
   {
      return pts[0] != other.pts[0] or pts[1] != other.pts[1];
    }
  
    CU_HOST_DEVICE
    Point size() const {
      return pts[1] - pts[0];
    }
  
    /**
     * Bounding box intersection
     */
    CU_HOST_DEVICE
    BoundingBox& operator&=(const BoundingBox& other)
    {
      for(size_t i = 0; i < N; ++i) {
        if(pts[0][i] < other.pts[0][i])
          pts[0][i] = other.pts[0][i];
        if(pts[1][i] > other.pts[1][i])
          pts[1][i] = other.pts[1][i];
      }
      return *this;
    }
  
    CU_HOST_DEVICE
    BoundingBox operator&(const BoundingBox& other) const
    {
      BoundingBox copy(*this);
      return copy &= other;
    }
  
    CU_HOST_DEVICE
    BoundingBox& operator*=(const BoundingBox& other) 
   {
      return (*this) &= (other);
    }
  
    CU_HOST_DEVICE
    BoundingBox operator*(const BoundingBox& other) const
    {
      BoundingBox copy(*this);
      return copy &= other;
    }
  
    /**
     * Bounding box union
     */
    CU_HOST_DEVICE
    BoundingBox& operator|=(const BoundingBox& other)
    {
      for(size_t i = 0; i < N; ++i) {
        if(pts[0][i] > other.pts[0][i])
          pts[0][i] = other.pts[0][i];
        if(pts[1][i] < other.pts[1][i])
          pts[1][i] = other.pts[1][i];
      }
      return *this;
    }
  
    CU_HOST_DEVICE
    BoundingBox operator|(const BoundingBox& other) const
    {
      BoundingBox copy(*this);
      return copy |= other;
    }
  
    CU_HOST_DEVICE
    BoundingBox& operator+=(const BoundingBox& other) 
   {
      return (*this) |= (other);
    }
  
    CU_HOST_DEVICE
    BoundingBox operator+(const BoundingBox& other) const
    {
      BoundingBox copy(*this);
      return copy |= other;
    }
  
    /**
     * Adding a point
     */
    CU_HOST_DEVICE
    BoundingBox& operator|=(const Point& p)
    {
      for(size_t i = 0; i < N; ++i) {
        if(pts[0][i] > p[i])
          pts[0][i] = p[i];
        if(limits::is_integer) {
          if(pts[1][i] < p[i] + 1)
            pts[1][i] = p[i] + 1;
        } else {
          if(pts[1][i] < p[i])
            pts[1][i] = p[i];
        }
      }
      return *this;
    }
  
    CU_HOST_DEVICE
    BoundingBox operator|(const Point& p) const
    {
      BoundingBox copy(*this);
      return copy |= p;
    }
  
    CU_HOST_DEVICE
    friend BoundingBox operator&(const Point& p, const BoundingBox& b)
    {
      BoundingBox copy(b);
      return copy |= p;
    }
  
    CU_HOST_DEVICE
    const Point& operator[](int i) const
    {
      if(i == 0)
        return pts[0];
      return pts[1];
    }
  
    CU_HOST_DEVICE
    Point& operator[](int i)
    {
      if(i == 0)
        return pts[0];
      return pts[1];
    }
  
    /**
     * Check if a point is in the BoundingBox
     */
    CU_HOST_DEVICE
    bool contains(const Point& p) const
    {
      return (p.x() >= pts[0].x() and p.y() >= pts[0].y() and p.z() >= pts[0].z() 
             and p.x() < pts[1].x() and p.y() < pts[1].y() and p.z() < pts[1].z());
    }
  
    CU_HOST_DEVICE
    Point* data() {
      return &pts[0];
    }
  
    CU_HOST_DEVICE
    Point pmin() {
      return pts[0];
    }
  
    CU_HOST_DEVICE
    Point pmin() const {
      return pts[0];
    }
  
    CU_HOST_DEVICE
    Point pmax() {
      return pts[1];
    }
  
    CU_HOST_DEVICE
    Point pmax() const {
      return pts[1];
    }
  
    Point pts[2];
  
    CU_HOST_DEVICE
    friend std::ostream& operator<<(std::ostream& s, const BoundingBox& bbox)
    {
      s << bbox.pts[0] << " " << bbox.pts[1];
      return s;
    }
  
    CU_HOST_DEVICE
    friend std::istream& operator>>(std::istream& s, BoundingBox& bbox)
    {
      s >> bbox.pts[0] >> bbox.pts[1];
      return s;
    }
  
  #ifndef COMPILE_CUDA
    friend QTextStream& operator<<(QTextStream& s, const BoundingBox& bbox)
    {
      s << bbox.pts[0] << " " << bbox.pts[1];
      return s;
    }
  
    friend QTextStream& operator>>(QTextStream& s, BoundingBox& bbox)
    {
      s >> bbox.pts[0] >> bbox.pts[1];
      return s;
    }
  #endif
  };
}
#endif // BOUNDINGBOX_HPP
