// Gnurbs - A curve and surface library
// Copyright (C) 2008-2026 Eric Bechet
//
// See the LICENSE file for contributions and license information.
// Please report all bugs and problems to <bechet@cadxfem.org>.
//

#include <iostream>
#include <vector>

#include "nvtkdisplay.h"
#include "nsurface.h"


class mysurface: public nparametricsurface
{
public:
  mysurface(): minu_(-1.0),maxu_(1.0),minv_(-1.0),maxv_(1.0),perm(false),perturb(0.) { transf[0]=npoint(1,0,0,0);transf[1]=npoint(0,1,0,0);transf[2]=npoint(0,0,1,0); transf[3]=npoint(0,0,0,1);}
  mysurface(npoint tr_[4]): minu_(-1.0),maxu_(1.0),minv_(-1.0),maxv_(1.0),perm(false),perturb(0.) {for (int i=0;i<4;++i) for (int j=0;j<4;++j) transf[i][j]=tr_[i][j];}
  mysurface(npoint tr1,npoint tr2,npoint tr3,npoint tr4): minu_(-1.0),maxu_(1.0),minv_(-1.0),maxv_(1.0),perm(false),perturb(0.) {for (int i=0;i<4;++i) {transf[0][i]=tr1[i];transf[1][i]=tr2[i];transf[2][i]=tr3[i];transf[3][i]=tr4[i];}}
  virtual int degree(void) const {return -1;}    // returns the degree in u and v or -1 if non polynomial
  virtual int degree_u(void) const {return -1;}    // returns the degree in u
  virtual int degree_v(void) const {return -1;}    // returnd the degree in v
  virtual int nb_CP_u(int iv=0) const {return 0;}    // number of control points in u
  virtual int nb_CP_v(int iu=0) const {return 0;}    // number of control points in v
  virtual double min_u() const {return minu_;} // minimal allowed value for u
  virtual double max_u() const {return maxu_;} // maximal allowed value for u
  virtual double min_v() const {return minv_;} // minimal allowed value for v
  virtual double max_v() const {return maxv_;} // maximal allowed value for v
  virtual void P(double u,double v, npoint& ret) const 
  { 
    double u_=2*(u-minu_)/(maxu_-minu_)-1.0;
    double v_=2*(v-minv_)/(maxv_-minv_)-1.0;
    if (perm) { double p=u_;u_=v_;v_=p; }
    double x=u_;
    double y=v_;
    double z=cos(sqrt(u_*u_+v_*v_)*2)*0.5+sin(u_*9)*0.05+sin(v_*10)*0.05;
//    double z=cos(2*u*u+2*v*v);
    double r=sqrt(x*x+y*y+z*z);
    if (r<0.5)
    {
      x=x+x*0.1*(r-0.5)*perturb;
      y=y+y*0.1*(r-0.5)*perturb;
      z=z+z*0.1*(r-0.5)*perturb;
    }
    npoint pt=npoint(x,y,z,1.0);

    for (int i=0;i<4;++i)
    {
      ret[i]=0.;
      for (int j=0;j<4;++j)
        ret[i]+=pt[j]*transf[j][i];
    }
  }    // point on the surface S(u,v)
  virtual npoint CP(int whichu,int whichv) const {return transf[3];}    // returns a control point (const)
  virtual int nb_CP(void) const {return 0;}    // returns the total number of control points
  virtual npoint CP(int which) const { return transf[3];}    // returns a control point
  virtual nsurface* clone() const {return new mysurface(*this);} // clones the surface (allocated on the heap)
private:
  npoint transf[4]; // transformation matrix.
public:
  double minu_;
  double maxu_;
  double minv_;
  double maxv_;
  bool perm;
  double perturb;
};


void find_reference_points(const nparametricsurface &surface, std::vector<npoint> &vecpts)
{
  // you have to implement this function
  double x=2.0,y=2.0,z=2.0;
  vecpts.push_back(npoint(x,y,z,1.0)); // example : adds a point at the end of the list.
}

int main(void)
{
  data_container data;
  nvtkdisplay display;
  mysurface surf1;
  mysurface surf2(npoint(-1,0,0,0),npoint(0,1,0,0),npoint( 0,0,-1,0),npoint(0,0,2,1));
  surf2.perturb=20.0;

  std::vector<npoint> pts1;
  std::vector<npoint> pts2;

  find_reference_points(surf1,pts1);
  find_reference_points(surf2,pts2);

  // displays the surface and the list of points
  surf1.Display(data);
  data.setcolorquads(color(255,0,0));
  surf2.Display(data);
  properties prop=data.getproppoints();
  prop.pointsize=3.0; // size of the points
  prop.c=color(0,255,0); // set color to green
  data.setproppoints(prop);
  for(int i=0;i<pts1.size();++i) 
    data.add_point(pts1[i]);
  for(int i=0;i<pts2.size();++i) 
    data.add_point(pts2[i]);

  display.init_data(data);
  display.display();
  return 0;
}