#include <cmath>
#include <fstream>
#include "Contour.hpp"

#include <QFile>
#include <QTextStream>
#include <stdio.h>
#include <QString>
#include <QStringList>

static QTextStream err(stderr);


#if defined(__APPLE__) || defined(linux)
#  include <sys/file.h>
#endif

#define REPORT_ERROR(filename, line_nb, error) \
  err << "Error in file " << filename << " on line " << line_nb << ":" << error << endl

namespace mgx
{
  typedef Vector<3,double> Point3d;

  /** @brief Default constructor. */
  Contour::Contour(): closed(true), regular(true) {}

  /** @brief Constructor from file.
    @bug   If the specified contains an incomplete or improper
    specification, the result of construction will be
    unknown.  Most probably, the object will unuseable;
    however, this is not currently reported.
    @param filename The contour file to load.
    */
  Contour::Contour(const QString &filename): fileName(filename), closed(true), regular(true)
  {
    reread();
  }

  Contour::Contour(const Contour& copy) : fileName(copy.fileName), pts(copy.pts), max(copy.max), 
                           min(copy.min), closed(copy.closed), regular(copy.regular) {}

  void Contour::reread() 
  {
    closed = false;
    QFile f(fileName);
    if(!f.open(QIODevice::ReadOnly))
    {
      err << "Error opening file " << fileName << ": " << f.errorString() << endl;
      return;
    }
#if defined(__APPLE__) || defined(linux)
    flock(f.handle(), LOCK_SH);
#endif
    QTextStream ts(&f);
    int ptid = 0, line_nb = 1;

    pts.clear();

    QString line = ts.readLine().trimmed();
    QStringList ver = line.split(" ");
    if(ver.size() != 3 or ver[0] != "cver")
    {
      err << "The file '" << fileName << "' doesn't contain a valid contour" << endl;
      return;
    }

    bool ok;
    int vmaj = ver[1].toInt(&ok);
    if(!ok)
    {
      REPORT_ERROR(fileName, line_nb, "The major version number is not a valid integer");
      return;
    }

    int vmin = ver[2].toInt(&ok);
    if(!ok)
    {
      REPORT_ERROR(fileName, line_nb, "The minor version number is not a valid integer");
      return;
    }

    int version = 100*vmaj+vmin;
    if(version == 101)
    {
      while(!ts.atEnd())
      {
        ++line_nb;
        QString line = ts.readLine().trimmed();
        if(line.startsWith("name:"))
        {
          // Do nothing of the name
        }
        else if(line.startsWith("points:"))
        {
          QStringList nbpts = line.mid(7).trimmed().split(" ");
          if(nbpts.size() != 2)
          {
            REPORT_ERROR(fileName, line_nb, "Error, there should be two values for the points number, not " << nbpts.size());
            return;
          }
          bool ok;
          int nb_pts = nbpts[0].toInt(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, number of control points is not a number");
            return;
          }
          pts.resize(nb_pts);
        }
        else if(line.startsWith("type:"))
        {
          QString type = line.mid(5).trimmed().toLower();
          if(type == "closed")
            closed = true;
          else if(type == "open")
            closed = false;
          else
          {
            REPORT_ERROR(fileName, line_nb, "Unknown contour type: " << type);
            return;
          }
        }
        else
        {
          if(ptid >= (int)pts.size())
          {
            REPORT_ERROR(fileName, line_nb, "There are more control points defined than declared");
            return;
          }
          QStringList pt = line.split(" ");
          if(pt.size() != 4)
          {
            REPORT_ERROR(fileName, line_nb, "A points should be described with four values");
            return;
          }
          bool ok;
          double x,y,z; //,m;
          x = pt[0].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, x coordinate is not a valid floating point value");
            return;
          }
          y = pt[1].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, y coordinate is not a valid floating point value");
            return;
          }
          z = pt[2].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, z coordinate is not a valid floating point value");
            return;
          }
          //m = pt[3].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, m coordinate is not a valid floating point value");
            return;
          }
          Point3d v(x,y,z);
          pts[ptid++] = v;
        }
      }
    }
    else if(version == 102)
    {
      while(!ts.atEnd())
      {
        ++line_nb;
        QString line = ts.readLine().trimmed();
        if(line.startsWith("name:"))
        {
          // Do nothing of the name
        }
        else if(line.startsWith("points:"))
        {
          QStringList nbpts = line.mid(7).trimmed().split(" ");
          if(nbpts.size() != 2)
          {
            REPORT_ERROR(fileName, line_nb, "Error, there should be two values for the points number, not " << nbpts.size());
            return;
          }
          bool ok;
          int nb_pts = nbpts[1].toInt(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, number of control points is not a number");
            return;
          }
          pts.resize(nb_pts);
        }
        else if(line.startsWith("type:"))
        {
          QString type = line.mid(5).trimmed().toLower();
          if(type.size() != 2)
          {
            REPORT_ERROR(fileName, line_nb, "Error reading contour type: '" << type << "'");
            return;
          }
          if(type[0] == 'c')
            closed = true;
          else if(type[0] == 'o')
            closed = false;
          else
          {
            REPORT_ERROR(fileName, line_nb, "Unknown contour type: " << type);
            return;
          }
          if(type[1] == 'e')
          {
            regular = false;
            err << "Warning, this contour reader doesn't handle contours that are not regular" << endl;
          }
          else if(type[1] == 'r')
          {
            regular = true;
          }
          else
          {
            REPORT_ERROR(fileName, line_nb, "Unknown contour type: " << type);
            return;
          }
        }
        else
        {
          if(ptid >= (int)pts.size())
          {
            REPORT_ERROR(fileName, line_nb, "There are more control points defined than declared");
            return;
          }
          QStringList pt = line.split(" ");
          if(pt.size() != 4)
          {
            REPORT_ERROR(fileName, line_nb, "A points should be described with four values");
            return;
          }
          bool ok;
          double x,y,z;//,m;
          x = pt[0].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, x coordinate is not a valid floating point value");
            return;
          }
          y = pt[1].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, y coordinate is not a valid floating point value");
            return;
          }
          z = pt[2].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, z coordinate is not a valid floating point value");
            return;
          }
          //m = pt[3].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, m coordinate is not a valid floating point value");
            return;
          }
          Point3d v(x,y,z);
          pts[ptid++] = v;
        }
      }
    }
    else if(version == 103)
    {
      while(!ts.atEnd())
      {
        ++line_nb;
        QString line = ts.readLine().trimmed();
        if(line.startsWith("name:"))
        {
          // Do nothing of the name
        }
        else if(line.startsWith("points:"))
        {
          QStringList nbpts = line.mid(7).trimmed().split(" ");
          if(nbpts.size() != 2)
          {
            REPORT_ERROR(fileName, line_nb, "Error, there should be two values for the points number, not " << nbpts.size());
            return;
          }
          bool ok;
          int nb_pts = nbpts[1].toInt(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, number of control points is not a number");
            return;
          }
          pts.resize(nb_pts);
        }
        else if(line.startsWith("type:"))
        {
          QString type = line.mid(5).trimmed().toLower();
          if(type.size() != 2)
          {
            REPORT_ERROR(fileName, line_nb, "Error reading contour type: '" << type << "'");
            return;
          }
          if(type[0] == 'c')
            closed = true;
          else if(type[0] == 'o')
            closed = false;
          else
          {
            REPORT_ERROR(fileName, line_nb, "Unknown contour type: " << type);
            return;
          }
          if(type[1] == 'e')
          {
            regular = false;
            err << "Warning, this contour reader doesn't handle contours that are not regular" << endl;
          }
          else if(type[1] == 'r')
          {
            regular = true;
          }
          else
          {
            REPORT_ERROR(fileName, line_nb, "Unknown contour type: " << type);
            return;
          }
        }
        else if(line.startsWith("samples:"))
        {
          // ignore that
        }
        else
        {
          if(ptid >= (int)pts.size())
          {
            REPORT_ERROR(fileName, line_nb, "There are more control points defined than declared");
            return;
          }
          QStringList pt = line.split(" ");
          if(pt.size() != 4)
          {
            REPORT_ERROR(fileName, line_nb, "A points should be described with four values");
            return;
          }
          bool ok;
          double x,y,z,m;
          x = pt[0].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, x coordinate is not a valid floating point value");
            return;
          }
          y = pt[1].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, y coordinate is not a valid floating point value");
            return;
          }
          z = pt[2].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, z coordinate is not a valid floating point value");
            return;
          }
          m = pt[3].toDouble(&ok);
          if(!ok)
          {
            REPORT_ERROR(fileName, line_nb, "Error, m coordinate is not a valid floating point value");
            return;
          }
          Point3d v(x,y,z);
          pts[ptid++] = v;
        }
      }
    }
    else
    {
      REPORT_ERROR(fileName, line_nb, "Cannot handle version " << vmaj << "  " << vmin);
      return;
    }

    if(ptid < (int)pts.size())
    {
      REPORT_ERROR(fileName, line_nb, "There are less control points defined than declared");
      return;
    }

    if (closed) {
      pts.push_back(pts[0]);
      pts.push_back(pts[1]);
      pts.push_back(pts[2]);
    }
#if defined(__APPLE__) || defined(linux)
    flock(f.handle(), LOCK_UN);
#endif
    f.close();
  }

  /** @brief Copy constructor. */
  Contour& Contour::operator=(const Contour& c) {
    pts = c.pts;
    max = c.max;
    min = c.min;
    closed = c.closed;

    return *this;
  }

  /** @brief Function operator.
    @param t Time along the curve, between zero and one.

    This calculates the 3-dimensional coordinates of the curve
    at time t, where t is between zero and one.  If t is less then
    zero, it is set to zero.  If it is larger than one, it is set to one.
    */
  Vector<3,double> Contour::operator()(double t) const {
    if (t < 0.0) t = 0.0;
    else if (t > 1.0) t = 1.0;

    double k = 1.0 + double(pts.size() - 3) * t;
    int i = int(k) + 2;
    double p = k - floor(k);

    return
      pts[i - 3] * Basis0(p)
      + pts[i - 2] * Basis1(p)
      + pts[i - 1] * Basis2(p)
      + pts[i]     * Basis3(p);
  }

  /** @brief Return the maximum x, y and z coordinates. */
  const Vector<3,double>& Contour::getMax() const {
    return max;
  }

  /** @brief Return the minimum x, y and z coordinates. */
  const Vector<3,double>& Contour::getMin() const {
    return min;
  }

  /** @brief Return the arc length for a parameter interval.
    @param a Parameter for the start of the interval.
    @param b Parameter for the end of the interval.
    @param dt The step used for forward differencing.

    It is expected that a and b are both in the interval of [0, 1] and
    that dt is less than b - a.
    */
  double Contour::length(double a, double b, double dt) {
    if (a < 0.0) a = 0.0;
    if (a > 1.0) a = 1.0;
    if (b > 1.0) b = 1.0;
    if (b < a + dt) return 0.0;

    double l = 0.0;
    Vector<3,double> p = (*this)(a), q;
    while (a <= b) {
      a += dt;
      q = (*this)(a);
      l += norm(p-q);
      p = q;
    }
    return l;
  }

  /** @brief Return the paramter after traveling a given arc length.
    @param t Parameter for the start of the interval.
    @param l The arc length to travel.
    @param dt The step used for forward differencing.

    It is expected that t is in the interval of [0, 1] and that l is
    positive.
    */
  double Contour::travel(double t, double l, double dt) {
    if (t < 0.0) t = 0.0;
    if (l < 0.0) return 0.0;
    if (l < dt) return t;

    //Vector<3,double> start = (*this)(t);
    double t_length = 0.0;
    double u = t + dt;
    while (t_length < l && t < 1.0) {
      t_length += this->length(t, u, dt);
      t = u;
      u += dt;
    }
    return t;
  }

  /** @brief Finds the tangent vector at a given parameter value.
    @param t Parameter value for where to take the tangent.
    @param dt The interval size used for the numerical evaluation.
    */
  Vector<3,double> Contour::tangent(double t, double dt) {
    Vector<3,double> r = (*this)(t + dt) - (*this)(t - dt);
    r.normalize();
    return r;
  }

  /** @brief Finds the normal vector at a given parameter value.
    @param t Parameter value for where to take the normal.
    @param dt The interval size used for the numerical evaluation.
    */
  Vector<3,double> Contour::normal(double t, double dt) {
    Vector<3,double> tvec = tangent(t, dt);
    return Vector<3,double>(tvec.y(), tvec.x(), 0);
  }

  /** @brief A basis function. */
  double Contour::Basis0(double t) const {
    return 1.0/6.0 * (-pow(t, 3) + 3 * pow(t, 2) - 3 * t + 1);
  }

  /** @brief A basis function. */
  double Contour::Basis1(double t) const {
    return 1.0/6.0 * (3 * pow(t, 3) - 6 * pow(t, 2) + 4);
  }

  /** @brief A basis function. */
  double Contour::Basis2(double t) const {
    return 1.0/6.0 * (-3 * pow(t, 3) + 3 * pow(t, 2) + 3 * t + 1);
  }

  /** @brief A basis function. */
  double Contour::Basis3(double t) const {
    return 1.0/6.0 * (pow(t, 3));
  }

}
