Nurbs2D.cpp

/*
 * Nurbs2D.cpp
 *
 *  Created on: Jul 13, 2010
 *      Author: jordt
 */

#include "Nurbs2D.hh"
#include <Base/Math/Vector3.hh>
#include <Base/Math/Vector4.hh>
#include <Base/Geometry/HomgPoint2D.hh>
#include <Base/Geometry/HomgPoint3D.hh>
#include <Base/Math/Operators.hh>

#include "NurbsTemplateInit.hh"

using namespace std;
using namespace BIAS;



// This function allocates and initializes the control point- and weight-arrays
template <class ControlPointsT>
Nurbs2D<ControlPointsT>::Nurbs2D(int knotsPerDim, int degree ) : NurbsBase(knotsPerDim,degree,2){

  // Allocate and initialize
  controlPoints_ = new ControlPointsT*[numCtrlPoints_];
  weights_ = new double*[numCtrlPoints_];
  for (int x = 0; x < numCtrlPoints_; x++){
    controlPoints_[x] = new ControlPointsT[numCtrlPoints_];
    weights_[x] = new double[numCtrlPoints_];
    for (int y = 0; y < numCtrlPoints_; y++){
        TemplateInit::Init(controlPoints_[x][y]);
        weights_[x][y] = 1.0;
    }
  }

}

// The serialized not-inline implementation of the nurbs evaluation function
template <class ControlPointsT>
void Nurbs2D<ControlPointsT>::Eval(const Vector2<double> &nurbsParam, ControlPointsT &result) const{

  BIASASSERT(nurbsParam[0] >= 0.0 && nurbsParam[0] <= 1.0 &&
             nurbsParam[1] >= 0.0 && nurbsParam[1] <= 1.0)

  // Initialize the return value
  TemplateInit::Init(result);

  // Get the relevant control points
  double weightSum = 0.0;
  Vector2<int> relevantBegin;
  Vector2<int> relevantEnd;
  GetRelevantBaseCtrlPoints(nurbsParam, relevantBegin, relevantEnd);
  double* baseFuncCache = new double[degree_+1];

  // Evaluate the relevant base functions and use them to wheight the
  // control points
  for (int x = relevantBegin[0]; x <= relevantEnd[0]; x++)
    for (int y = relevantBegin[1]; y <= relevantEnd[1]; y++){
      double weight = BaseFunc(0,x,degree_,nurbsParam[0],baseFuncCache)*
                   BaseFunc(1,y,degree_,nurbsParam[1],baseFuncCache)*
                   weights_[x][y];
        weightSum += weight;

        result += (weight * controlPoints_[x][y]);
      }
  if (weightSum != 0.0)
    result = result * weightSum;

  delete[] baseFuncCache;
}



// The serialized not-inline implementation of the nurbs evaluation function
template <class ControlPointsT>
void Nurbs2D<ControlPointsT>::EvalCustom(const Vector2<double> &nurbsParam, ControlPointsT &result,
                                         const ControlPointsT** controlPoints) const{

  // Initialize the return value
  TemplateInit::Init(result);

  // Get the relevant control points
  double weightSum = 0.0;
  Vector2<int> relevantBegin;
  Vector2<int> relevantEnd;
  GetRelevantBaseCtrlPoints(nurbsParam, relevantBegin, relevantEnd);
  double* baseFuncCache = new double[degree_+1];

  // Evaluate the relevant base functions and use them to wheight the
  // control points
  for (int x = relevantBegin[0]; x <= relevantEnd[0]; x++)
    for (int y = relevantBegin[1]; y <= relevantEnd[1]; y++){
      double weight = BaseFunc(0,x,degree_,nurbsParam[0],baseFuncCache)*
                   BaseFunc(1,y,degree_,nurbsParam[1],baseFuncCache)*
                   weights_[x][y];
        weightSum += weight;

        result += (weight * controlPoints[x][y]);
      }

  // Warning: This has changed!
  if (weightSum != 0.0)
    result = result * (1.0 / weightSum);

  delete[] baseFuncCache;
}

template <class ControlPointsT>
void Nurbs2D<ControlPointsT>::Init(const ControlPointsT minX,const ControlPointsT maxX,
                                   const ControlPointsT minY,const ControlPointsT maxY){
  for (int x = 0; x < numCtrlPoints_; x++)
    for (int y = 0; y < numCtrlPoints_; y++){
      double ratioX = (double)x/(double)(numCtrlPoints_-1);
      double ratioY = (double)y/(double)(numCtrlPoints_-1);
      controlPoints_[x][y] = ratioX * maxX + (1.0 - ratioX) * minX
                           + ratioY * maxY + (1.0 - ratioY) * minY;
      }

}


namespace BIAS {

  // Implement NURBS2D for all common types of vectors:
  template class Nurbs2D<double >;
  template class Nurbs2D<Vector<double> >;
  template class Nurbs2D<Vector2<double> >;
  template class Nurbs2D<Vector3<double> >;
  template class Nurbs2D<Vector4<double> >;
  template class Nurbs2D<HomgPoint2D >;
  template class Nurbs2D<HomgPoint3D >;

}