// 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>.
//

#ifndef __NBSPLINE_H
#define __NBSPLINE_H


#include "ncurve.h"
#include <vector>
#include <list>

// Forward declaration.
class Reparametrization;

/// A B-spline curve.
class nbspline : public nparametriccurve
{
protected:
  std::vector<npoint> val; //!< Control points.
  int nCP; //!< Number of control points.
  std::vector<double> knots; //!< Nodal sequence.
  int nknots; //!< Number of knots.
  int deg; //!< Curve degree.
public:

  /// \brief Default constructor.
  nbspline() : nCP(1),nknots(2),deg(0),knots(2),val(1)
  {
    knots[0]=0;
    knots[1]=1;
    val[0]=npoint();
  }

  /// \brief Construct a B-spline.
  /// \param[in] nCP_ number of control points.
  /// \param[in] degree_ curve degree.
  /// \warning The control points and the nodal sequence are unitialized.
  nbspline(int nCP_,int degree_) : val(nCP_),nCP(nCP_),deg(degree_),nknots(nCP_+degree_+1),knots(nCP_+degree_+1) {}

  /// \brief Construct a B-spline.
  /// \param[in] nCP_ number of control points.
  /// \param[in] degree_ curve degree.
  /// \param[in] knots_ knots vector.
  /// \param[in] pts_ control points vector.
  nbspline(int nCP_,int degree_, double knots_[], npoint pts_[]) : val(nCP_),nCP(nCP_),deg(degree_),nknots(nCP_+degree_+1),knots(nCP_+degree_+1) 
  {
    for (int i=0;i<nknots;++i) knots[i]=knots_[i];
    for (int i=0;i<nCP;++i) val[i]=pts_[i];
  }
  
  /// \brief Reset the B-spline to a new B-spline.
  /// \param[in] nCP_ new number of control points.
  /// \param[in] degree_ new curve degree.
  /// \warning If nCP_ > old number of control points, the new control points are not initialized.
  /// If nCP_ < old number of control points, the last old control points are lost. 
  /// The first control points are keept.
  /// Same for the nodal sequence.
  void reset(int nCP_,int degree_)
  {
    val.resize(nCP_);
    nCP=nCP_,deg=degree_,nknots= (nCP_+degree_+1);
    knots.resize(nCP_+degree_+1);
  }

  /// \brief Returns the curve's degree.
  /// \return Degree.
  virtual int degree(void) const
  {
    return deg ;
  }

  /// \brief Returns the number of control points.
  /// \return Number of control points.
  virtual int nb_CP(void) const
  {
    return nCP ;
  }

  /// \brief Returns the number of knots.
  /// \return Number of knots.
  virtual int nb_knots(void) const
  {
    return nknots ;
  }

  /// \brief Returns the minimum knot value.
  /// \return Minimum knot value.
  virtual double min_u(void) const
  {
    return knots[deg];
  }

  /// \brief Returns the maximum knot value.
  /// \return maximum knot value.
  virtual double max_u(void) const
  {
    return knots[nknots-deg-1];
  }

  /// \brief Evaluates basis functions.
  /// \param[in] which wanted basis function number.
  /// \param[in] u evaluation parameter.
  /// \return Value of the basis function number "which" at u.
  double Basis(int which,double u) const;

  /// \brief Position of a point on the curve.
  /// \param[in] u_ parameter.
  /// \param[out] ret point on the curve corresponding to u_.
  virtual void P(double u_,npoint& ret) const ;
  virtual void dP(double u_,npoint& p, npoint & dpdu) const ;
  virtual void d2P(double u_,npoint& p, npoint & dpdu, npoint & d2pdu2) const ;
  virtual void dPdu(double u_,npoint& ret) const ;    // derivative of a point on the curve
  virtual void d2Pdu2(double u_,npoint& ret) const ;    // second derivative of a point on the curve

  /// \brief Returns a control point (constant version).
  /// \param[in$ which wanted control point number.
  /// \return Control point.
  virtual npoint CP(int which) const
  {
    return val[which];
  }

  /// \brief Returns a control point (mutable version).
  /// \param[in$ which wanted control point number.
  /// \return Control point.
  virtual npoint& CP(int which)
  {
    return val[which];
  }

  /// \brief Set a control point.
  /// \param[in] which control point number to set.
  /// \param[in] pt new value of the control point.
  virtual void set_CP(int which,const npoint& pt)
  {
    val[which]=pt;
  } ;

  /// \brief Returns a knot (constant version).
  /// \param[in] which wanted knot number.
  /// \return Knot value.
  virtual double u(int which) const
  {
    return knots[which];
  }

  /// \brief Returns a knot (mutable version).
  /// \param[in] which wanted knot number.
  /// \return Knot value.
  virtual double& u(int which)
  {
    return knots[which];
  }

  /// \brief Insert a knot u r times in the nodal sequence.
  /// \param[in] u knot value.
  /// \param[in] r number of insertions.
  void insertknot(double u,int r);

  /// \brief Delete a knot. If the knot is repeated, all the repetitions are deleted.
  /// param[in] r knot number to delete.
  /// tolerance < 0 means always delete.
  /// \return Number of deletions.
  int deleteknot(int r,int s=1,double TOL=-1);

  /// \brief Removes unnecessary knots.
  int clean(double TOL=1e-9);

  /// \brief Extent B-spline to u.
  /// \param[in] u parameter.
  void extend(double u);

  /// \brief Saturate nodal sequence.
  void saturate();

  /// \brief Elevate degree by 1. The curve is unchanged.
  void degree_elevation();

  virtual nbspline* clone() const
  {
    return new nbspline(*this);
  }

  /// \brief Returns the knots number i in the nodal sequence such that
  /// \f$ u \in [U_i, U_{i+1}] \f$.
  /// \param[in] u parameter.
  /// \return Value of i.
  int findspan(double u) const ;
  void basisfuns(int i, double u, std::vector<double> &N) const ;
  void gradbasisfuns(int i, double u,int n,std::vector<std::vector<double> > &N) const ;
  int basisfuns(double u,std::vector<double> &N) const ;
  int gradbasisfuns(double u,int nu,std::vector<double> &N) const;
  double Nip(int p,int m,std::vector<double> &knots, int i, double u);
  /// \brief generates an (interpolating) uniform knot vector.
  void initknots(double start,double end);

private:
  void insertknot_(double u,int k,int s,int r);
  void curvepoint(double u,npoint &C) const ;
  void dercurvepoint(double u,int nu,npoint &C) const ;
  void dercurvepoint(double u,int nu,std::vector<npoint> &dp) const;
  void findspanmult(double u,int &k,int &s);
  void recursivedivide(nbspline &bs,int deep) const;
  friend class nbsplinesurface;
  friend class Reparametrization;
};

#endif // __NBSPLINE_H