---

gunu_bezier_derivatives.cpp

Go to the documentation of this file.

/* LIBGUL - Geometry Utility Library
 * Copyright (C) 1998-1999 Norbert Irmer
 *
 * This library is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Library General Public
 * License as published by the Free Software Foundation; either
 * version 2 of the License, or (at your option) any later version.
 *
 * This library is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Library General Public License for more details.
 *
 * You should have received a copy of the GNU Library General Public
 * License along with this library; if not, write to the
 * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
 * Boston, MA 02111-1307, USA.
 */
 
#include "stdafx.h"

#include <iostream>

#include "gul_types.h"
#include "gul_vector.h"
#include "guar_bincoeff.h"
#include "gunu_basics.h"
#include "gunu_derivatives.h"


namespace gunu {

using gul::rtr;
using gul::point;
using gul::point1;
using gul::point2;
using gul::hpoint;
using gul::hpoint1;
using gul::hpoint2;
using gul::set;
using gul::euclid;
using gul::ortho;
using gul::Max;

/************************************************************************
  The following functions are just simplifications for Bezier functions.
  The only difference to the corresponding Nurbs functions is, that
  the knot span searching is missing, since Bezier knot vectors have
  only one knot span
*************************************************************************/

/*-------------------------------------------------------------------------
  Calculates the first 'd' derivatives of a curve at a point 'u'. The results
  are homogeneous points (like the control points), returned in 'CK[]'
  (see "The NURBS Book")
--------------------------------------------------------------------------- */  
template< class T, class EP >
void BezierBSPCurveDerivatives( 
           const T u, const int p, const Ptr< T >& U,
           const Ptr< EP >& Pw, const int d, Ptr< EP >& CK )
{
  int du,k,j,span;
  EP v1;
  Ptr< Ptr< T > > ders;
    
  du = Min( d, p );

  ders.reserve_place( reserve_stack(Ptr< T >,du+1), du+1 );
  for( j = 0; j < p+1; j++ )
    ders[j].reserve_place( reserve_stack(T,p+1), p+1 );

  for( k = p+1; k <= d; k++ )
    set( CK[k], (T)0.0 );   /* Ableitungen > Grad = 0 */

  span = p;
  CalcDersBasisFuns( u, span, p, du, U, ders );
  
  for( k = 0; k <= du; k++ )
  {
    set( CK[k], (T)0.0 );
    for( j = 0; j <= p; j++ )
    {
      v1 = ders[k][j] * Pw[span-p+j];
      CK[k] = CK[k] + v1;
    }
  }
}  
// template instantiation
template void BezierBSPCurveDerivatives( 
     const float u, const int p, const Ptr< float >& U,
     const Ptr< point<float> >& Pw, const int d, Ptr< point<float> >& CK );
template
void BezierBSPCurveDerivatives( 
      const double u, const int p, const Ptr< double >& U,
      const Ptr< point<double> >& Pw, const int d, Ptr< point<double> >& CK );
template
void BezierBSPCurveDerivatives( 
       const float u, const int p, const Ptr< float >& U,
       const Ptr< point2<float> >& Pw, const int d, Ptr< point2<float> >& CK );
template
void BezierBSPCurveDerivatives( 
      const double u, const int p, const Ptr< double >& U,
      const Ptr< point2<double> >& Pw, const int d, Ptr< point2<double> >& CK );

/*---------------------------------------------------------------------------
  Calculates mixed partial derivatives of a NURBS surface. Let 'k' and 'l' be
  the number of derivatives in 'u' and 'v' direction, then all mixed ders with
  'k+l <= d' are calculated. They are returned in 'SKL[k][l]'.  
  (see "The NURBS Book")
------------------------------------------------------------------------- */  
template< class T, class EP >
void BezierBSPSurfaceDerivatives(
                const T u, const T v,
                const int pu, const Ptr< T >& U,
                const int pv, const Ptr< T >& V,
                const Ptr< Ptr< EP > >& Pw,
                const int d, Ptr< Ptr< EP > >& SKL )
{                
  int du,dv,k,l,uspan,vspan,s,r,dd;
  EP v1;
  Ptr< EP > temp;
  Ptr< Ptr< T > > Nu,Nv;

  du = Min( d, pu );
  dv = Min( d, pv );

  Nu.reserve_place( reserve_stack(Ptr< T >,du+1), du+1 );
  for( k = 0; k < pu+1; k++ )
    Nu[k].reserve_place( reserve_stack(T,pu+1), pu+1 );

  Nv.reserve_place( reserve_stack(Ptr< T >,dv+1), dv+1 );
  for( k = 0; k < pv+1; k++ )
    Nv[k].reserve_place( reserve_stack(T,pv+1), pv+1 );

  temp.reserve_place( reserve_stack(EP,pv+1), pv+1 );

  for( k = pu+1; k <= d; k++ )
    for( l = 0; l <= d - k; l++ )
      set( SKL[k][l], (T)0.0 );
      
  for( l = pv+1; l <= d; l++ )
    for( k = 0; k <= d - l; k++ )
      set( SKL[k][l], (T)0.0 );
      
  uspan = pu;
  CalcDersBasisFuns( u, uspan, pu, du, U, Nu );
 
  vspan = pv;
  CalcDersBasisFuns( v, vspan, pv, dv, V, Nv );
  
  for( k = 0; k <= du; k++ )
  {
    for( s = 0; s <= pv; s++ )
    {
      set( temp[s], (T)0.0 );
      for( r = 0; r <= pu; r++ )
      {
        v1 = Nu[k][r] * Pw[vspan-pv+s][uspan-pu+r];
        temp[s] = temp[s] + v1;
      }  
    }
    dd = Min( d - k, dv );
    for( l = 0; l <= dd; l++ )
    {
      set( SKL[k][l], (T)0.0 );
      for( s = 0; s <= pv; s++ )
      {
        v1 = Nv[l][s] * temp[s];
        SKL[k][l] = SKL[k][l] + v1;
      }
    }
  }
}  
// template instantiation
template void BezierBSPSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr< point<float> > >& Pw,
                const int d, Ptr< Ptr< point<float> > >& SKL );           
template void BezierBSPSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr< point<double> > >& Pw,
                const int d, Ptr< Ptr< point<double> > >& SKL );

template void BezierBSPSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr< hpoint<float> > >& Pw,
                const int d, Ptr< Ptr< hpoint<float> > >& SKL );           
template void BezierBSPSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr< hpoint<double> > >& Pw,
                const int d, Ptr< Ptr< hpoint<double> > >& SKL );

/*-------------------------------------------------------------------------
  Calculates the first 'd' partial derivatives of the projection of a NURBS
  curve into euclidian space (3-dimensions), so the calculated ders are points
  in 3-dimensional space   
  (see "The NURBS Book")
------------------------------------------------------------------------- */  
template< class T, class HP, class EP >
void BezierCurveDerivatives(
                    const T u, const int p,
                    const Ptr< T >& U, const Ptr< HP >& Pw, const int d,
                    Ptr< EP >& CK )             
{
  Ptr< HP > CKh;
  EP v,v1;
  T w0;
  int k,i;

  CKh.reserve_place( reserve_stack(HP,d+1), d+1 );
     
  BezierBSPCurveDerivatives<T,HP>( u, p, U, Pw, d, CKh );

  w0 = CKh[0].w;
  
  for( k = 0; k <= d; k++ )
  {
    v = ortho( CKh[k] );
    
    for( i = 1; i <= k; i++ )
    {
      v1 = (rtr<T>::BinCoeff(k,i) * CKh[i].w) * CK[k-i];
      v = v - v1;
    }
    CK[k] = ((T)1.0 / w0) * v;       
  }
}
// template instantiation
template void BezierCurveDerivatives(
          const float u, const int p,
          const Ptr< float >& U, const Ptr< hpoint<float> >& Pw, const int d,
          Ptr< point<float> >& CK );
template void BezierCurveDerivatives(
          const double u, const int p,
          const Ptr< double >& U, const Ptr< hpoint<double> >& Pw, const int d,
          Ptr< point<double> >& CK );

template void BezierCurveDerivatives(
         const float u, const int p,
         const Ptr< float >& U, const Ptr< hpoint2<float> >& Pw, const int d,
         Ptr< point2<float> >& CK );
template void BezierCurveDerivatives(
         const double u, const int p,
         const Ptr< double >& U, const Ptr< hpoint2<double> >& Pw, const int d,
         Ptr< point2<double> >& CK );

/*-------------------------------------------------------------------------
  Calculates the first 'd' mixed partial derivatives of the projection of a
  NURBS surface into euclidian space (3-dimensions), so the calculated ders
  are points in 3-dimensional space. They are returned in 'SKL[k][l]'.  
  (see "The NURBS Book")
------------------------------------------------------------------------- */
template< class T, class HP, class EP >          
void BezierSurfaceDerivatives(
                const T u, const T v,
                const int pu, const Ptr< T >& U,
                const int pv, const Ptr< T >& V,
                const Ptr< Ptr < HP > >& Pw,
                const int d, Ptr< Ptr < EP > >& SKL )
{
  int i,j,k,l;
  EP v1,v2,vh;
  T w00;
  Ptr< Ptr < HP > > SKLh;
  
  SKLh.reserve_place( reserve_stack(Ptr < HP >,d+1), d+1 );
  for( i = 0; i < d+1; i++ )
    SKLh[i].reserve_place( reserve_stack(HP,d+1), d+1 );
    
  BezierBSPSurfaceDerivatives<T,HP>( u, v, pu, U, pv, V, Pw, d, SKLh );

  w00 = SKLh[0][0].weight();
  
  for( k = 0; k <= d; k++ )
  {
    for( l = 0; l <= d-k; l++ )
    {
      v1 = ortho( SKLh[k][l] );

      for( j = 1; j <= l; j++ )
      {
        vh = (rtr<T>::BinCoeff(l,j) * SKLh[0][j].weight()) * SKL[k][l-j]; 
        v1 = v1 - vh;
      }  
      for( i = 1; i <= k; i++ )
      {
        vh = (rtr<T>::BinCoeff(k,i) * SKLh[i][0].weight()) * SKL[k-i][l];
        v1 = v1 - vh;

        set( v2, (T)0.0 );
        for( j = 1; j <= l; j++ )
        {
          vh = (rtr<T>::BinCoeff(l,j) * SKLh[i][j].weight()) * SKL[k-i][l-j];
          v2 = v2 + vh;
        }
        vh =  rtr<T>::BinCoeff(k,i) * v2;
        v1 = v1 - vh;
      }
      SKL[k][l] = ((T)1.0 / w00) * v1;
    }    
  }
}

// template instantiation
template void BezierSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr < hpoint<float> > >& Pw,
                const int d, Ptr< Ptr < point<float> > >& SKL );
template void BezierSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr < hpoint<double> > >& Pw,
                const int d, Ptr< Ptr < point<double> > >& SKL );

template void BezierSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr < hpoint1<float> > >& Pw,
                const int d, Ptr< Ptr < point1<float> > >& SKL );
template void BezierSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr < hpoint1<double> > >& Pw,
                const int d, Ptr< Ptr < point1<double> > >& SKL );

#ifdef _MSC_VER
template void BezierSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr < point<float> > >& Pw,
                const int d, Ptr< Ptr < point<float> > >& SKL );
template void BezierSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr < point<double> > >& Pw,
                const int d, Ptr< Ptr < point<double> > >& SKL );

template void BezierSurfaceDerivatives(
                const float u, const float v,
                const int pu, const Ptr< float >& U,
                const int pv, const Ptr< float >& V,
                const Ptr< Ptr < point1<float> > >& Pw,
                const int d, Ptr< Ptr < point1<float> > >& SKL );
template void BezierSurfaceDerivatives(
                const double u, const double v,
                const int pu, const Ptr< double >& U,
                const int pv, const Ptr< double >& V,
                const Ptr< Ptr < point1<double> > >& Pw,
                const int d, Ptr< Ptr < point1<double> > >& SKL );
#endif

}

---