#pragma once

#ifndef T_GEOMETRY_INCLUDED

#define T_GEOMETRY_INCLUDED

#include "tutil.h"

#include <cmath></cmath>

#undef DVAPI

#undef DVVAR

#ifdef TGEOMETRY_EXPORTS

#define DVAPI DV_EXPORT_API

#define DVVAR DV_EXPORT_VAR

#else

#define DVAPI DV_IMPORT_API

#define DVVAR DV_IMPORT_VAR

#endif

//=============================================================================

inline double logNormalDistribuitionUnscaled(double x, double x0, double w)

  { return exp(-0.5*pow(log(x/x0)/w, 2.0))/x; }

inline double logNormalDistribuition(double x, double x0, double w)

  { return logNormalDistribuitionUnscaled(x, x0, w)/(w*sqrt(2.0*M_PI)); }

//=============================================================================

template <class t=""> class T3DPointT;</class>

template <class t=""> class T4DPointT;</class>

/*

* This is an example of how to use the TPointT, the TRectT and the TAffine

* classes.

*/

/*!  The template class TPointT defines the x- and y-coordinates of a point.

*/

template <class t=""></class>

class TPointT {

public:

  T x, y;

  inline TPointT() : x(0), y(0){};

  inline TPointT(T x, T y) : x(x), y(y){};

  inline explicit TPointT(const T3DPointT<t> &p);</t>

  inline explicit TPointT(const T4DPointT<t> &p);</t>

  inline TPointT& operator+=(const TPointT &a)

    { return x += a.x, y += a.y, *this; };

  inline TPointT& operator-=(const TPointT &a)

    { return x -= a.x, y -= a.y, *this; };

  inline TPointT operator+(const TPointT &a) const

    { return TPointT(x + a.x, y + a.y); };

  inline TPointT operator-(const TPointT &a) const

    { return TPointT(x - a.x, y - a.y); };

  inline TPointT operator-() const

    { return TPointT(-x, -y); };

  //! Scalar(dot) Product

  inline T operator*(const TPointT &a) const

    { return x*a.x + y*a.y; }

  inline TPointT operator*=(T a)

    { return x *= a, y *= a, *this; }

  inline TPointT operator*(T a) const

    { return TPointT(x*a, y*a); }

  friend inline TPointT operator*(T a, const TPointT &b)

    { return TPointT(a*b.x, a*b.y); }

  inline TPointT operator/(T a) const

    { return TPointT(x/a, y/a); }

  inline TPointT operator/=(T a)

    { return x /= a, y /= a, *this; }

  inline bool operator==(const TPointT &a) const

    { return x == a.x && y == a.y; }

  inline bool operator!=(const TPointT &a) const

    { return !(*this == a); }

  friend inline std::ostream &operator<<(std::ostream &out, const TPointT &p)

    { return out << "(" << p.x << ", " << p.y << ")"; }

  /*! Rotate a point 90 degrees (counterclockwise).

  \param p a point.

  \return the rotated point

  \sa rotate270 */

  friend inline TPointT rotate90(const TPointT &p) // 90 counterclockwise

    { return TPointT(-p.y, p.x); }

  /*! Rotate a point 270 degrees (clockwise).

  \param p a point.

  \return the rotated point

  \sa rotate90 */

  friend inline TPointT rotate270(const TPointT &p)  // 90 clockwise

    { return TPointT(p.y, -p.x); }

  //! This helper function returns the square of the absolute value of the point

  friend inline T norm2(const TPointT &a)

    { return a*a; }

  //! This helper function returns the square of the distance between two points

  friend inline T tdistance2(const TPointT &a, const TPointT &b)

    { return norm2(a - b); }

  //! the cross product

  friend inline T cross(const TPointT &a, const TPointT &b)

    { return a.x*b.y - a.y*b.x; }

};

template <>

inline bool TPointT<double>::operator==(const TPointT<double> &a) const</double></double>

  { return tdistance2(*this, a) <= TConsts::epsilon * TConsts::epsilon; }

//-----------------------------------------------------------------------------

typedef TPointT<int> TPoint, TPointI;</int>

typedef TPointT<double> TPointD;</double>

#ifdef _WIN32

template class DVAPI TPointT<int>;</int>

template class DVAPI TPointT<double>;</double>

#endif

//-----------------------------------------------------------------------------

/*!

\relates TPointT

This helper function converts a TPoint (TPointT<int>) into a TPointD</int>

*/

inline TPointD convert(const TPoint &p)

  { return TPointD(p.x, p.y); }

/*!

\relates TPointT

This helper function converts a TPointD (TPointT<double>) into a TPoint</double>

*/

inline TPoint convert(const TPointD &p)

  { return TPoint(tround(p.x), tround(p.y)); }

/*!

\relates TPointT

This helper function returns the absolute value of the specified point

*/

inline double norm(const TPointD &p) { return std::sqrt(norm2(p)); }

/*!

\relates TPointT

This helper function returns the normalized version of the specified point

*/

inline TPointD normalize(const TPointD &p) {

  double n = norm(p);

  assert(n);

  return p*(1/n);

}

/*!

\relates TPointT

This helper function returns the normalized version of the specified point

or zero if it is not possible

*/

inline TPointD normalizeOrZero(const TPointD &p) {

  double n = norm2(p);

  return fabs(n) > TConsts::epsilon*TConsts::epsilon ? p*(1/sqrt(n)) : TPointD();

}

/*!

\relates TPointT

This helper function returns the distance between two points

*/

inline double tdistance(const TPointD &p1, const TPointD &p2)

  { return norm(p2 - p1); }

/*!

returns the angle of the point p in polar coordinates

n.b atan(-y) = -pi/2, atan(x) = 0, atan(y) = pi/2, atan(-x) = pi

*/

inline double atan(const TPointD &p) { return atan2(p.y, p.x); }

//=============================================================================

template <class t=""></class>

class DVAPI T3DPointT {

public:

  T x, y, z;

  inline T3DPointT() : x(), y(), z() {}

  inline T3DPointT(T x, T y, T z) : x(x), y(y), z(z) {}

  inline T3DPointT(const TPointT<t> &p, T z) : x(p.x), y(p.y), z(z) {}</t>

  inline explicit T3DPointT(const T4DPointT<t> &p);</t>

  inline TPointT<t>& xy() { return *(TPointT<t>*)this; }</t></t>

  inline TPointT<t>& xy() const { return *(const TPointT<t>*)this; }</t></t>

  inline T3DPointT &operator+=(const T3DPointT &a)

    { return x += a.x, y += a.y, z += a.z, *this; }

  inline T3DPointT &operator-=(const T3DPointT &a)

    { return x -= a.x, y -= a.y, z -= a.z, *this; }

  inline T3DPointT operator+(const T3DPointT &a) const

    { return T3DPointT(x + a.x, y + a.y, z + a.z); }

  inline T3DPointT operator-(const T3DPointT &a) const

    { return T3DPointT(x - a.x, y - a.y, z - a.z); }

  inline T3DPointT operator-() const

    { return T3DPointT(-x, -y, -z); }

  //! Scalar(dot) Product

  inline T operator*(const T3DPointT &a) const

    { return x*a.x + y*a.y + z*a.z; }

  inline T3DPointT operator*=(T a)

    { return x *= a, y *= a, z *= a, *this; }

  inline T3DPointT operator*(T a) const

    { return T3DPointT(x*a, y*a, z*a); }

  friend inline T3DPointT operator*(T a, const T3DPointT &b)

    { return T3DPointT(a*b.x, a*b.y, a*b.z); }

  inline T3DPointT operator/(T a) const

    { return T3DPointT(x/a, y/a, z/a); }

  inline T3DPointT operator/=(T a)

    { return x /= a, y /= a, z /= a, *this; }

  inline bool operator==(const T3DPointT &a) const

    { return x == a.x && y == a.y && z == a.z; }

  inline bool operator!=(const T3DPointT &a) const

    { return !(*this == a); }

  friend inline std::ostream &operator<<(std::ostream &out, const T3DPointT &p)

    { return out << "(" << p.x << ", " << p.y << ", " << p.z << ")"; }

  //! This helper function returns the square of the absolute value of the point

  friend inline T norm2(const T3DPointT &a)

    { return a*a; }

  //! This helper function returns the square of the distance between two points

  friend inline T tdistance2(const T3DPointT &a, const T3DPointT &b)

    { return norm2(a - b); }

  //! the cross product

  friend inline T3DPointT cross(const T3DPointT &a, const T3DPointT &b) {

    return T3DPointT( a.y*b.z - b.y*a.z,

                      a.z*b.x - b.z*a.x,

                      a.x*b.y - b.x*a.y);

  }

};

template <>

inline bool T3DPointT<double>::operator==(const T3DPointT<double> &a) const</double></double>

  { return tdistance2(*this, a) <= TConsts::epsilon * TConsts::epsilon; }

template <class t=""></class>

inline TPointT<t>::TPointT(const T3DPointT<t> &p) : x(p.x), y(p.y) {};</t></t>

//=============================================================================

typedef T3DPointT<int> T3DPoint, T3DPointI;</int>

typedef T3DPointT<double> T3DPointD;</double>

#ifdef _WIN32

template class DVAPI T3DPointT<int>;</int>

template class DVAPI T3DPointT<double>;</double>

#endif

//-----------------------------------------------------------------------------

inline T3DPointD convert(const T3DPoint &p)

  { return T3DPointD(p.x, p.y, p.z); }

inline T3DPoint convert(const T3DPointD &p)

  { return T3DPoint(tround(p.x), tround(p.y), tround(p.z)); }

inline double norm(const T3DPointD &p)

  { return std::sqrt(norm2(p)); }

inline T3DPointD normalize(const T3DPointD &p) {

  double n = norm(p);

  assert(n);

  return p*(1/n);

}

inline T3DPointD normalizeOrZero(const T3DPointD &p) {

  double n = norm2(p);

  return fabs(n) > TConsts::epsilon*TConsts::epsilon ? p*(1/sqrt(n)) : T3DPointD();

}

inline double tdistance(const T3DPointD &p1, const T3DPointD &p2)

  { return norm(p2 - p1); }

//=============================================================================

template <class t=""></class>

class DVAPI T4DPointT {

public:

  T x, y, z, w;

  inline T4DPointT() : x(), y(), z(), w() {}

  inline T4DPointT(T x, T y, T z, T w) : x(x), y(y), z(z), w(w) {}

  inline T4DPointT(const TPointT<t> &p, T z, T w) : x(p.x), y(p.y), z(z), w(w) {}</t>

  inline T4DPointT(const T3DPointT<t> &p, T w) : x(p.x), y(p.y), z(p.z), w(w) {}</t>

  inline TPointT<t>& xy() { return *(TPointT<t>*)this; }</t></t>

  inline T3DPointT<t>& xyz() { return *(T3DPointT<t>*)this; }</t></t>

  inline TPointT<t>& xy() const { return *(const TPointT<t>*)this; }</t></t>

  inline T3DPointT<t>& xyz() const { return *(const T3DPointT<t>*)this; }</t></t>

  inline bool operator==(const T4DPointT &p) const

    { return x == p.x && y == p.y && z == p.z && w == p.w; }

  inline bool operator!=(const T4DPointT &p) const

    { return !(*this == p); }

  friend inline std::ostream &operator<<(std::ostream &out, const T4DPointT &p)

    { return out << "(" << p.x << ", " << p.y << ", " << p.z << ", " << p.w << ")"; }

};

template <>

inline bool T4DPointT<double>::operator==(const T4DPointT<double> &a) const {</double></double>

  T4DPointT<double> d(x - a.x, y - a.y, z - a.z, w - a.w);</double>

  return d.x*d.x + d.y*d.y + d.z*d.z + d.w*d.w

      <= TConsts::epsilon * TConsts::epsilon;

}

template <class t=""></class>

inline TPointT<t>::TPointT(const T4DPointT<t> &p) : x(p.x), y(p.y) {};</t></t>

template <class t=""></class>

inline T3DPointT<t>::T3DPointT(const T4DPointT<t> &p) : x(p.x), y(p.y), z(p.z) {};</t></t>

//=============================================================================

typedef T4DPointT<int> T4DPoint, T4DPointI;</int>

typedef T4DPointT<double> T4DPointD;</double>

#ifdef _WIN32

template class DVAPI T4DPointT<int>;</int>

template class DVAPI T4DPointT<double>;</double>

#endif

//-----------------------------------------------------------------------------

//!\relates T4DPointT

inline T4DPointD convert(const T4DPoint &p)

  { return T4DPointD(p.x, p.y, p.z, p.w); }

inline T4DPoint convert(const T4DPointD &p)

  { return T4DPoint(tround(p.x), tround(p.y), tround(p.z), tround(p.w)); }

//=============================================================================

/*!

TThickPoint describe a thick point.

\relates TThickQuadratic, TThickCubic

*/

class DVAPI TThickPoint final : public TPointD {

public:

  double thick;

  TThickPoint() : TPointD(), thick(0) {}

  TThickPoint(double _x, double _y, double _thick = 0)

      : TPointD(_x, _y), thick(_thick) {}

  TThickPoint(const TPointD &_p, double _thick = 0)

      : TPointD(_p.x, _p.y), thick(_thick) {}

  TThickPoint(const T3DPointD &_p) : TPointD(_p.x, _p.y), thick(_p.z) {}

  TThickPoint(const TThickPoint &_p) : TPointD(_p.x, _p.y), thick(_p.thick) {}

  inline TThickPoint &operator=(const TThickPoint &a) {

    x     = a.x;

    y     = a.y;

    thick = a.thick;

    return *this;

  }

  inline TThickPoint &operator+=(const TThickPoint &a) {

    x += a.x;

    y += a.y;

    thick += a.thick;

    return *this;

  }

  inline TThickPoint &operator-=(const TThickPoint &a) {

    x -= a.x;

    y -= a.y;

    thick -= a.thick;

    return *this;

  }

  inline TThickPoint operator+(const TThickPoint &a) const {

    return TThickPoint(x + a.x, y + a.y, thick + a.thick);

  }

  inline TThickPoint operator-(const TThickPoint &a) const {

    return TThickPoint(x - a.x, y - a.y, thick - a.thick);

  }

  inline TThickPoint operator-() const { return TThickPoint(-x, -y, -thick); }

  bool operator==(const TThickPoint &p) const {

    return x == p.x && y == p.y && thick == p.thick;

  }

  bool operator!=(const TThickPoint &p) const {

    return x != p.x || y != p.y || thick != p.thick;

  }

};

inline double operator*(const TThickPoint &a, const TThickPoint &b) {

  return a.x * b.x + a.y * b.y + a.thick * b.thick;

}

inline TThickPoint operator*(double a, const TThickPoint &p) {

  return TThickPoint(a * p.x, a * p.y, a * p.thick);

}

inline TThickPoint operator*(const TThickPoint &p, double a) {

  return TThickPoint(a * p.x, a * p.y, a * p.thick);

}

/*!

\relates TPointD

This helper function converts a TThickPoint into a TPointD

*/

inline TPointD convert(const TThickPoint &p) { return TPointD(p.x, p.y); }

/*!

\relates TThickPoint

This helper function returns the square of the distance between two thick points

(only x and y are used)

*/

inline double tdistance2(const TThickPoint &p1, const TThickPoint &p2) {

  return norm2(convert(p2 - p1));

}

/*!

\relates TThickPoint

This helper function returns the distance between two thick  points

(only x and y are used)

*/

inline double tdistance(const TThickPoint &p1, const TThickPoint &p2) {

  return norm(convert(p2 - p1));

}

inline std::ostream &operator<<(std::ostream &out, const TThickPoint &p) {

  return out << "(" << p.x << ", " << p.y << ", " << p.thick << ")";

}

//=============================================================================

//!	This is a template class representing a generic vector in a plane, i.e.

//! a point.

/*!

                It is a data structure with two objects in it representing

   coordinate of the point and

                the basic operations on it.

        */

template <class t=""></class>

class DVAPI TDimensionT {

public:

  T lx, ly;

  /*!

          Constructs a vector of two elements, i.e. a point in a plane.

  */

  TDimensionT() : lx(), ly() {}

  TDimensionT(T _lx, T _ly) : lx(_lx), ly(_ly) {}

  /*!

          Copy constructor.

  */

  TDimensionT(const TDimensionT &d) : lx(d.lx), ly(d.ly) {}

  /*!

          Vector addition.

  */

  TDimensionT &operator+=(TDimensionT a) {

    lx += a.lx;

    ly += a.ly;

    return *this;

  }

  /*!

          Difference of two vectors.

  */

  TDimensionT &operator-=(TDimensionT a) {

    lx -= a.lx;

    ly -= a.ly;

    return *this;

  }

  /*!

          Addition of two vectors.

  */

  TDimensionT operator+(TDimensionT a) const {

    TDimensionT ris(*this);

    return ris += a;

  }

  /*!

          Vector difference.

  */

  TDimensionT operator-(TDimensionT a) const {

    TDimensionT ris(*this);

    return ris -= a;

  }

  /*!

          Compare vectors and returns \e true if are equals element by element.

  */

  bool operator==(const TDimensionT &d) const {

    return lx == d.lx && ly == d.ly;

  }

  /*!

          Compare vectors and returns \e true if are not equals element by

     element.

  */

  bool operator!=(const TDimensionT &d) const { return !operator==(d); }

};

//=============================================================================

typedef TDimensionT<int> TDimension, TDimensionI;</int>

typedef TDimensionT<double> TDimensionD;</double>

//=============================================================================

template <class t=""></class>

inline std::ostream &operator<<(std::ostream &out, const TDimensionT<t> &p) {</t>

  return out << "(" << p.lx << ", " << p.ly << ")";

}

#ifdef _WIN32

template class DVAPI TDimensionT<int>;</int>

template class DVAPI TDimensionT<double>;</double>

#endif

//=============================================================================

//! Specifies the corners of a rectangle.

/*!

x0 specifies the x-coordinate of the bottom-left corner of a rectangle.

y0 specifies the y-coordinate of the bottom-left corner of a rectangle.

x1 specifies the x-coordinate of the upper-right corner of a rectangle.

y1 specifies the y-coordinate of the upper-right corner of a rectangle.

*/

template <class t=""></class>

class DVAPI TRectT {

public:

  /*! if x0>x1 || y0>y1 then rect is empty

if x0==y1 && y0==y1 and rect is a  TRectD then rect is empty */

  T x0, y0;

  T x1, y1;

  /*! makes an empty rect */

  TRectT();

  TRectT(T _x0, T _y0, T _x1, T _y1) : x0(_x0), y0(_y0), x1(_x1), y1(_y1){};

  TRectT(const TRectT &rect)

      : x0(rect.x0), y0(rect.y0), x1(rect.x1), y1(rect.y1){};

  TRectT(const TPointT<t> &p0, const TPointT<t> &p1)  // non importa l'ordine</t></t>

      : x0(std::min((T)p0.x, (T)p1.x)),

        y0(std::min((T)p0.y, (T)p1.y)),

        x1(std::max((T)p0.x, (T)p1.x)),

        y1(std::max((T)p0.y, (T)p1.y)){};

  TRectT(const TPointT<t> &bottomLeft, const TDimensionT<t> &d);</t></t>

  TRectT(const TDimensionT<t> &d);</t>

  void empty();

  /*! TRectD is empty if and only if (x0>x1 || y0>y1) || (x0==y1 && y0==y1);

TRectI  is empty if x0>x1 || y0>y1 */

  bool isEmpty() const;

  T getLx() const;

  T getLy() const;

  TDimensionT<t> getSize() const { return TDimensionT<t>(getLx(), getLy()); };</t></t>

  TPointT<t> getP00() const { return TPointT<t>(x0, y0); };</t></t>

  TPointT<t> getP10() const { return TPointT<t>(x1, y0); };</t></t>

  TPointT<t> getP01() const { return TPointT<t>(x0, y1); };</t></t>

  TPointT<t> getP11() const { return TPointT<t>(x1, y1); };</t></t>

  //! Returns the union of two source rectangles.

  //!The union is the smallest rectangle that contains both source rectangles.

  TRectT<t> operator+(const TRectT<t> &rect) const {  // unione</t></t>

    if (isEmpty())

      return rect;

    else if (rect.isEmpty())

      return *this;

    else

      return TRectT<t>(std::min((T)x0, (T)rect.x0), </t>

                       std::min((T)y0, (T)rect.y0),

                       std::max((T)x1, (T)rect.x1),

                       std::max((T)y1, (T)rect.y1));

  };

  TRectT<t> &operator+=(const TRectT<t> &rect) {  // unione</t></t>

    return *this = *this + rect;

  };

  TRectT<t> &operator*=(const TRectT<t> &rect) {  // intersezione</t></t>

    return *this = *this * rect;

  };

  //!Returns the intersection of two existing rectangles.

  //The intersection is the largest rectangle contained in both existing rectangles.

  TRectT<t> operator*(const TRectT<t> &rect) const {  // intersezione</t></t>

    if (isEmpty() || rect.isEmpty())

      return TRectT<t>();</t>

    else if (rect.x1 < x0 || x1 < rect.x0 || rect.y1 < y0 || y1 < rect.y0)

      return TRectT<t>();</t>

    else

      return TRectT<t>(std::max((T)x0, (T)rect.x0), std::max((T)y0, (T)rect.y0),</t>

                       std::min((T)x1, (T)rect.x1),

                       std::min((T)y1, (T)rect.y1));

  };

  TRectT<t> &operator+=(const TPointT<t> &p) {  // shift</t></t>

    x0 += p.x;

    y0 += p.y;

    x1 += p.x;

    y1 += p.y;

    return *this;

  };

  TRectT<t> &operator-=(const TPointT<t> &p) {</t></t>

    x0 -= p.x;

    y0 -= p.y;

    x1 -= p.x;

    y1 -= p.y;

    return *this;

  };

  TRectT<t> operator+(const TPointT<t> &p) const {</t></t>

    TRectT<t> ris(*this);</t>

    return ris += p;

  };

  TRectT<t> operator-(const TPointT<t> &p) const {</t></t>

    TRectT<t> ris(*this);</t>

    return ris -= p;

  };

  bool operator==(const TRectT<t> &r) const {</t>

    return x0 == r.x0 && y0 == r.y0 && x1 == r.x1 && y1 == r.y1;

  };

  bool operator!=(const TRectT<t> &r) const {</t>

    return x0 != r.x0 || y0 != r.y0 || x1 != r.x1 || y1 != r.y1;

  };

  bool contains(const TPointT<t> &p) const {</t>

    return x0 <= p.x && p.x <= x1 && y0 <= p.y && p.y <= y1;

  };

  bool contains(const TRectT<t> &b) const {</t>

    return x0 <= b.x0 && x1 >= b.x1 && y0 <= b.y0 && y1 >= b.y1;

  };

  bool overlaps(const TRectT<t> &b) const {</t>

    return x0 <= b.x1 && x1 >= b.x0 && y0 <= b.y1 && y1 >= b.y0;

  };

  TRectT<t> enlarge(T dx, T dy) const {</t>

    if (isEmpty()) return *this;

    return TRectT<t>(x0 - dx, y0 - dy, x1 + dx, y1 + dy);</t>

  };

  TRectT<t> enlarge(T d) const { return enlarge(d, d); };</t>

  TRectT<t> enlarge(TDimensionT<t> d) const { return enlarge(d.lx, d.ly); };</t></t>

};

//-----------------------------------------------------------------------------

typedef TRectT<int> TRect, TRectI;</int>

typedef TRectT<double> TRectD;</double>

#ifdef _WIN32

template class DVAPI TRectT<int>;</int>

template class DVAPI TRectT<double>;</double>

#endif

//=============================================================================

// check this, not final version

/*!

\relates TRectT

Convert a TRectD into a TRect

*/

inline TRect convert(const TRectD &r) {

  return TRect((int)(r.x0 + 0.5), (int)(r.y0 + 0.5), (int)(r.x1 + 0.5),

               (int)(r.y1 + 0.5));

}

/*!

\relates TRectT

Convert a TRect into a TRectD

*/

inline TRectD convert(const TRect &r) { return TRectD(r.x0, r.y0, r.x1, r.y1); }

// template?

/*!

\relates TRectT

\relates TPointT

*/

inline TRectD boundingBox(const TPointD &p0, const TPointD &p1) {

  return TRectD(std::min(p0.x, p1.x), std::min(p0.y, p1.y),  // bottom left

                std::max(p0.x, p1.x), std::max(p0.y, p1.y));  // top right

}

/*!

\relates TRectT

\relates TPointT

*/

inline TRectD boundingBox(const TPointD &p0, const TPointD &p1,

                          const TPointD &p2) {

  return TRectD(std::min({p0.x, p1.x, p2.x}), std::min({p0.y, p1.y, p2.y}),

                std::max({p0.x, p1.x, p2.x}), std::max({p0.y, p1.y, p2.y}));

}

/*!

\relates TRectT

\relates TPointT

*/

inline TRectD boundingBox(const TPointD &p0, const TPointD &p1,

                          const TPointD &p2, const TPointD &p3) {

  return TRectD(

      std::min({p0.x, p1.x, p2.x, p3.x}), std::min({p0.y, p1.y, p2.y, p3.y}),

      std::max({p0.x, p1.x, p2.x, p3.x}), std::max({p0.y, p1.y, p2.y, p3.y}));

}

//-----------------------------------------------------------------------------

// TRectT is a rectangle that uses thick points

template <>

inline TRectT<int>::TRectT() : x0(0), y0(0), x1(-1), y1(-1) {}</int>

template <>

inline TRectT<int>::TRectT(const TPointT<int> &bottomLeft,</int></int>

                           const TDimensionT<int> &d)</int>

    : x0(bottomLeft.x)

    , y0(bottomLeft.y)

    , x1(bottomLeft.x + d.lx - 1)

    , y1(bottomLeft.y + d.ly - 1){};

template <>

inline TRectT<int>::TRectT(const TDimensionT<int> &d)</int></int>

    : x0(0), y0(0), x1(d.lx - 1), y1(d.ly - 1){};

template <>

inline bool TRectT<int>::isEmpty() const {</int>

  return x0 > x1 || y0 > y1;

}

template <>

inline void TRectT<int>::empty() {</int>

  x0 = y0 = 0;

  x1 = y1 = -1;

}

// Is the adding of one here to account for the thickness?

template <>

inline int TRectT<int>::getLx() const {</int>

  return x1 >= x0 ? x1 - x0 + 1 : 0;

}

template <>

inline int TRectT<int>::getLy() const {</int>

  return y1 >= y0 ? y1 - y0 + 1 : 0;

}

template <>

inline TRectT<double>::TRectT() : x0(0), y0(0), x1(0), y1(0) {}</double>

template <>

inline TRectT<double>::TRectT(const TPointT<double> &bottomLeft,</double></double>

                              const TDimensionT<double> &d)</double>

    : x0(bottomLeft.x)

    , y0(bottomLeft.y)

    , x1(bottomLeft.x + d.lx)

    , y1(bottomLeft.y + d.ly){};

template <>

inline TRectT<double>::TRectT(const TDimensionT<double> &d)</double></double>

    : x0(0.0), y0(0.0), x1(d.lx), y1(d.ly){};

template <>

inline bool TRectT<double>::isEmpty() const {</double>

  return (x0 == x1 && y0 == y1) || x0 > x1 || y0 > y1;

}

template <>

inline void TRectT<double>::empty() {</double>

  x0 = y0 = x1 = y1 = 0;

}

template <>

inline double TRectT<double>::getLx() const {</double>

  return x1 >= x0 ? x1 - x0 : 0;

}

template <>

inline double TRectT<double>::getLy() const {</double>

  return y1 >= y0 ? y1 - y0 : 0;

}

//-----------------------------------------------------------------------------

inline TRectD &operator*=(TRectD &rect, double factor) {

  rect.x0 *= factor;

  rect.y0 *= factor;

  rect.x1 *= factor;

  rect.y1 *= factor;

  return rect;

}

//-----------------------------------------------------------------------------

inline TRectD operator*(const TRectD &rect, double factor) {

  TRectD result(rect);

  return result *= factor;

}

//-----------------------------------------------------------------------------

inline TRectD &operator/=(TRectD &rect, double factor) {

  assert(factor != 0.0);

  return rect *= (1.0 / factor);

}

//-----------------------------------------------------------------------------

inline TRectD operator/(const TRectD &rect, double factor) {

  assert(factor != 0.0);

  TRectD result(rect);

  return result *= 1.0 / factor;

}

//-----------------------------------------------------------------------------

template <class t=""></class>

inline std::ostream &operator<<(std::ostream &out, const TRectT<t> &r) {</t>

  return out << "(" << r.x0 << "," << r.y0 << ";" << r.x1 << "," << r.y1 << ")";

}

//=============================================================================