basis_sample.cc

Go to the documentation of this file.

/************************************************************************/
/*                                                                      */
/*    vspline - a set of generic tools for creation and evaluation      */
/*              of uniform b-splines                                    */
/*                                                                      */
/*            Copyright 2017 - 2023 by Kay F. Jahnke                     */
/*                                                                      */
/*    Permission is hereby granted, free of charge, to any person       */
/*    obtaining a copy of this software and associated documentation    */
/*    files (the "Software"), to deal in the Software without           */
/*    restriction, including without limitation the rights to use,      */
/*    copy, modify, merge, publish, distribute, sublicense, and/or      */
/*    sell copies of the Software, and to permit persons to whom the    */
/*    Software is furnished to do so, subject to the following          */
/*    conditions:                                                       */
/*                                                                      */
/*    The above copyright notice and this permission notice shall be    */
/*    included in all copies or substantial portions of the             */
/*    Software.                                                         */
/*                                                                      */
/*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
/*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
/*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
/*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
/*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
/*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
/*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
/*    OTHER DEALINGS IN THE SOFTWARE.                                   */
/*                                                                      */
/************************************************************************/
 
/// basis_sample.cc - access to the b-spline basis functions
 
/// Calculate a set of equally-spaced samples of a b-spline basis function
/// This program takes four parameters on the command line:
/// - the spline's degree
/// - x0, the location of the initial sample
/// - step, the distance from one sample to the next
/// - normalize - if not zero, result will be normalized
/// These values are fed to 'sample_basis' (see below), which starts by
/// calculating the initial sample and then goes 'both ways' to add more
/// samples. The (optionally normalized) samples are echoed to std::cout.
/// example: basis_sample 3 .05 .25 1
 
/// vspline offers access to the b-spline basis functions via a helper
/// object 'basis_functor'. A good part of the calculation of the basis
/// function's value can be done beforehand and only depends on the
/// degree of the basis function. The basis_functor object performs
/// these calculations once when it's created, subsequent evaluations
/// are relatively fast and rely on the initially generated state.
/// For this example, we use long double values - this example is
/// not about speed, but rather about providing precise values.
/// The basis_functor object can also evaluate vectorized data,
/// but in this example we'll only use it's scalar form.
 
 
#include <deque>
#include <vector>
#include <iostream>
#include <iomanip>
#include <vspline/vspline.h>
 
/// 'sample_basis' returns a vector of samples of the basis function
/// if we have the basis function b(x), it returns these values:
/// bf ( x - x0 + step * k ) for all k E N,
/// wherever they are greater than zero. So for X0 = 0 and step = 1,
/// we get the ordinary unit-spaced sampling of the basis function,
/// with x0 != 0 and step == 1 a set of weights as it might be
/// applied to a set of coefficients to evaluate a b-spline, and
/// with step != 1, we get a set of values which follow the shape
/// of the basis function with more or less sample points than the
/// 'normal' unit-spaced sampling. What's surprising is that these
/// non-unit-step samplings (when multiplied with the step width)
/// also are (often very) close to a partition of unity.
 
typedef vspline::basis_functor < long double > bf_type ;
 
std::vector<long double> sample_basis ( bf_type bfunc ,
                                        long double x0 ,
                                        long double step ,
                                        bool normalize = false )
{
  auto y = bfunc ( x0 ) ;
  auto sum = y ;
 
  assert ( y > 0 ) ;
 
  std::deque<long double> accu ;
  accu.push_back ( y ) ;
  
  for ( int k = 1 ; ; k++ )
  {
    y = bfunc ( x0 - k * step ) ;
    if ( y <= 0 )
      break ;
    accu.push_front ( y ) ;
    sum += y ;
  }
  for ( int k = 1 ; ; k++ )
  {
    y = bfunc ( x0 + k * step ) ;
    if ( y <= 0 )
      break ;
    accu.push_back ( y ) ;
    sum += y ;
  }
 
  size_t sz = accu.size() ;
  std::vector < long double > result ( sz ) ;
  
  if ( normalize )
  {
    for ( size_t i = 0 ; i < sz ; i++ )
      result[i] = accu[i] / sum ;
  }
  else
  {
    for ( size_t i = 0 ; i < sz ; i++ )
      result[i] = accu[i] ;
  }
 
  return result ;
}
 
int main ( int argc, char * argv[] )
{
  if ( argc != 5 )
  {
    std::cerr
    << "use: basis_sample degree(int) x0(float) step(float) normalize(0/1)"
    << std::endl ;
    exit ( -1 ) ;
  }
 
  int degree = std::atoi ( argv[1] ) ;
  auto bfunc = bf_type ( degree ) ;
  auto x0 = std::atof ( argv[2] ) ;
  auto step = std::atof ( argv[3] ) ;
  bool normalize = std::atoi ( argv[4] ) ;
 
  auto sample = sample_basis ( bfunc , x0 , step , normalize ) ;
 
  std::cout << std::fixed << std::showpoint
            << std::setprecision
                (std::numeric_limits<long double>::max_digits10) ;
 
  for ( auto & e : sample )
    std::cout << e << std::endl ;
}