[vspline] 40/72: added two examples, gsm.cc and gsm2.cc, which calculate the gradient squared magnitude of an image, the first implementation using dirct evaluation of the gradients in a loop as a quick shot, the second implementation creating a dedicated gsm evaluator and producing the result by using vspline::index_remap

Kay F. Jahnke kfj-guest at moszumanska.debian.org
Sun Jul 2 09:02:41 UTC 2017


This is an automated email from the git hooks/post-receive script.

kfj-guest pushed a commit to branch master
in repository vspline.

commit f36607f6266ebbfd6cb6cc3bdcc1255b50b87916
Author: Kay F. Jahnke <kfjahnke at gmail.com>
Date:   Sat Mar 25 10:18:41 2017 +0100

    added two examples, gsm.cc and gsm2.cc, which calculate the gradient squared magnitude
    of an image, the first implementation using dirct evaluation of the gradients in a loop
    as a quick shot, the second implementation creating a dedicated gsm evaluator and
    producing the result by using vspline::index_remap
---
 bspline.h                         |  71 ++++++++++++-
 eval.h                            |   4 +-
 example/{slice.cc => channels.cc} |  65 ++++++++----
 example/gsm.cc                    | 136 ++++++++++++++++++++++++
 example/gsm2.cc                   | 216 ++++++++++++++++++++++++++++++++++++++
 example/slice.cc                  |   4 +
 example/splinus.cc                |   2 +
 unary_functor.h                   |  36 +++++++
 8 files changed, 510 insertions(+), 24 deletions(-)

diff --git a/bspline.h b/bspline.h
index d0c1f91..4c3fc72 100644
--- a/bspline.h
+++ b/bspline.h
@@ -173,8 +173,11 @@ struct bspline
   /// nD type for one boundary condition per axis
   typedef typename vigra::TinyVector < bc_code , dimension > bcv_type ;
 
-  /// elementary type of value_type, float or double
+  /// elementary type of value_type, like float or double
   typedef typename ExpandElementResult < value_type >::type real_type ;
+  enum { channels = ExpandElementResult < value_type >::size } ;
+
+  typedef bspline < real_type , dimension > channel_view_type ;
   
   // for internal calculations in the filter, we use the elementary type of value_type.
   // Note how in class bspline we make very specific choices about the
@@ -399,6 +402,72 @@ public:
     core = coeffs.subarray ( left_brace , left_brace + core_shape ) ;
   } ;
 
+  /// copy constructor. This will result in a view to the same data and is therefore
+  /// lightweight. But the data viewed by this spline must remain accessible.
+
+  bspline ( const bspline& other )
+  : strategy ( other.strategy ) ,
+    boundary_conditions ( other.boundary_conditions ) ,
+    container ( other.container ) ,
+    coeffs ( other.coeffs ) ,
+    core ( other.core ) ,
+    spline_degree ( other.spline_degree ) ,
+    bcv ( other.bcv ) ,
+    tolerance ( other.tolerance ) ,
+    smoothing ( other.smoothing ) ,
+    braced ( other.braced ) ,
+    horizon ( other.horizon ) ,
+    left_brace ( other.left_brace ) ,
+    right_brace ( other.right_brace ) ,
+    left_frame ( other.left_frame ) ,
+    right_frame ( other.right_frame ) ,
+    container_shape ( other.container_shape ) ,
+    core_shape ( other.core_shape ) ,
+    braced_shape ( other.braced_shape )
+  { } ;
+  
+  bspline operator= ( const bspline& other )
+  {
+    return bspline ( *this ) ;
+  }
+  
+  /// get a bspline object for a single channel of the data. This is also lightweight
+  /// and requires the viewed data to remain present as long as the channel view is used.
+  /// This is useful because, when working with multichannel data, often they will come
+  /// interleaved in memory, which assures good locality of access and therefore efficient
+  /// prefiltering and evaluation. If the channels require separate treatment, the
+  /// channel-specific operations are best performed at the same location with the
+  /// channel-specific operations in sequence, since this accesses memory in the same
+  /// vicinity - rather than processing each channel by itself, which traverses the
+  /// memory as many times as there are channels. Of course, if only one specific channel
+  /// needs to be operated upon, getting a view to the channel first and constructing
+  /// a b-spline over it is more efficient.
+  /// the channel view inherits all metrics from it's parent, only the MultiArrayViews
+  /// to the data are different.
+  
+  channel_view_type get_channel_view ( const int & channel )
+  {
+    assert ( channel < channels ) ;
+    
+    real_type * base = (real_type*) ( container.data() ) ;
+    base += channel ;
+    auto stride = container.stride() ;
+    stride *= channels ;
+    
+    MultiArrayView < dimension , real_type >
+      channel_container ( container.shape() , stride , base ) ;
+
+    return channel_view_type ( core_shape , 
+                               spline_degree ,
+                               bcv ,
+                               strategy ,
+                               horizon ,
+                               channel_container , // coefficient storage to 'adopt'
+                               tolerance ,
+                               smoothing ,
+                               0 ) ;
+  } ;
+
   /// prefilter converts the knot point data in the 'core' area into b-spline
   /// coefficients. Depending on the strategy chosen in the b-spline object's
   /// constructor, bracing/framing may be applied. Even if the degree of the
diff --git a/eval.h b/eval.h
index 54d4e60..c7be940 100644
--- a/eval.h
+++ b/eval.h
@@ -1350,8 +1350,8 @@ public:
   /// splines, albeit les efficiently.
   
   void eval ( const ic_v& select ,       // offsets to coefficient windows
-                const in_v& tune ,      // fractional parts of the coordinates
-                out_v & result ) const  // target
+              const in_v& tune ,      // fractional parts of the coordinates
+              out_v & result ) const  // target
   {
     if ( specialize == 0 )
     {
diff --git a/example/slice.cc b/example/channels.cc
similarity index 73%
copy from example/slice.cc
copy to example/channels.cc
index 6047db0..e07c4af 100644
--- a/example/slice.cc
+++ b/example/channels.cc
@@ -29,13 +29,17 @@
 /*                                                                      */
 /************************************************************************/
 
-/// slice.cc
+/// channels.cc
 ///
-/// build a 3D volume from samples of the RGB colour space
-/// build a spline over it and extract a 2D slice
+/// demonstrates the use of 'channel views'
+/// This example is derived from 'slice.cc', we use the same volume
+/// as source data. But instead of producing an image output, we create
+/// three separate colour channels of the bspline object and assert that
+/// the evaluation of the channel views is identical with the evaluation
+/// of the 'mother' spline.
 ///
 /// compile with:
-/// clang++ -std=c++11 -march=native -o slice -O3 -pthread -DUSE_VC=1 slice.cc -lvigraimpex -lVc
+/// clang++ -std=c++11 -march=native -o channels -O3 -pthread -DUSE_VC=1 channels.cc -lvigraimpex -lVc
 /// g++ also works.
 
 #include <vspline/vspline.h>
@@ -67,8 +71,17 @@ int main ( int argc , char * argv[] )
   // create quintic 3D b-spline object containing voxels
   vspline::bspline < voxel_type , 3 >
     space ( vigra::Shape3 ( 10 , 10 , 10 ) , 5 , bcv ) ;
-  
+
+  // here we create the channel view. Since these are merely views
+  // to the same data, no data will be copied, and it doesn't matter
+  // whether we create these views before or after prefiltering.
+
+  auto red_channel = space.get_channel_view ( 0 ) ;
+  auto green_channel = space.get_channel_view ( 1 ) ;
+  auto blue_channel = space.get_channel_view ( 2 ) ;
+
   // fill the b-spline's core with a three-way gradient
+
   for ( int z = 0 ; z < 10 ; z++ )
   {
     for ( int y = 0 ; y < 10 ; y++ )
@@ -86,10 +99,6 @@ int main ( int argc , char * argv[] )
   // prefilter the b-spline
   space.prefilter ( true ) ;
   
-  // get an evaluator for the b-spline
-  typedef vspline::evaluator < coordinate_type , voxel_type > ev_type ;
-  ev_type ev ( space ) ;
-  
   // now make a warp array with 1920X1080 3D coordinates
   warp_type warp ( vigra::Shape2 ( 1920 , 1080 ) ) ;
   
@@ -108,20 +117,34 @@ int main ( int argc , char * argv[] )
     }
   }
   
-  // this is where the result should go:
-  target_type target ( vigra::Shape2 ( 1920 , 1080 ) ) ;
-
-  // now we perform the remap, yielding the result
-  vspline::remap < ev_type , 2 > ( ev , warp , target ) ;
+  // get an evaluator for the b-spline
 
-  // store the result with vigra impex
-  vigra::ImageExportInfo imageInfo ( "slice.tif" );
+  typedef vspline::evaluator < coordinate_type , voxel_type > ev_type ;
+  ev_type ev ( space ) ;
   
-  vigra::exportImage ( target ,
-                      imageInfo
-                      .setPixelType("UINT8")
-                      .setCompression("100")
-                      .setForcedRangeMapping ( 0 , 255 , 0 , 255 ) ) ;
+  // the evaluators of the channel views have their own type:
   
+  typedef vspline::evaluator < coordinate_type , float > ch_ev_type ;
+  
+  // we create the three evaluators for the three channel views
+
+  ch_ev_type red_ev ( red_channel ) ;
+  ch_ev_type green_ev ( green_channel ) ;
+  ch_ev_type blue_ev ( blue_channel ) ;
+
+  // and make sure the evaluation results match
+
+  for ( int y = 0 ; y < 1080 ; y++ )
+  {
+    for ( int x = 0 ; x < 1920 ; x++ )
+    {
+      coordinate_type & c ( warp [ vigra::Shape2 ( x , y ) ] ) ;
+      assert ( ev ( c ) [ 0 ] == red_ev ( c ) ) ;
+      assert ( ev ( c ) [ 1 ] == green_ev ( c ) ) ;
+      assert ( ev ( c ) [ 2 ] == blue_ev ( c ) ) ;
+    }
+  }
+
+  std::cout << "success" << std::endl ;
   exit ( 0 ) ;
 }
diff --git a/example/gsm.cc b/example/gsm.cc
new file mode 100644
index 0000000..8a6718f
--- /dev/null
+++ b/example/gsm.cc
@@ -0,0 +1,136 @@
+/************************************************************************/
+/*                                                                      */
+/*    vspline - a set of generic tools for creation and evaluation      */
+/*              of uniform b-splines                                    */
+/*                                                                      */
+/*            Copyright 2015, 2016 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.                                   */
+/*                                                                      */
+/************************************************************************/
+
+/// gsm.cc
+///
+/// implementation of gsm.cc, performing the calculation of the
+/// gradient squared magnitude in a loop using two evaluators for
+/// the two derivatives
+///
+/// compile with:
+/// clang++ -std=c++11 -march=native -o gsm -O3 -pthread -DUSE_VC=1 gsm.cc -lvigraimpex -lVc
+/// g++ also works.
+
+#include <vspline/vspline.h>
+
+#include <vigra/stdimage.hxx>
+#include <vigra/imageinfo.hxx>
+#include <vigra/impex.hxx>
+
+#ifdef USE_VC
+const int VSIZE = vspline::vector_traits < float > :: size ;
+#else
+const int VSIZE = 1 ;
+#endif
+
+// we silently assume we have a colour image
+typedef vigra::RGBValue<float,0,1,2> pixel_type; 
+
+// coordinate_type has a 2D coordinate
+typedef vigra::TinyVector < float , 2 > coordinate_type ;
+
+// target_type is a 2D array of pixels  
+typedef vigra::MultiArray < 2 , pixel_type > target_type ;
+
+// b-spline evaluator producing float pixels
+typedef vspline::evaluator < coordinate_type , // incoming coordinate's type
+                             pixel_type ,      // singular result data type
+                             -1 , -1 ,         // no specializations
+                             VSIZE             // vector size
+                           > ev_type ;
+
+int main ( int argc , char * argv[] )
+{
+  vigra::ImageImportInfo imageInfo ( argv[1] ) ;
+
+  // we want a b-spline with natural boundary conditions
+  vigra::TinyVector < vspline::bc_code , 2 > bcv ( vspline::NATURAL ) ;
+  
+  // create cubic 2D b-spline object containing the image data
+  vspline::bspline < pixel_type , 2 > bspl ( imageInfo.shape() , 3 , bcv ) ;
+  
+  // load the image data into the b-spline's core
+  vigra::importImage ( imageInfo , bspl.core ) ;
+  
+  // prefilter the b-spline
+  bspl.prefilter() ;
+  
+  // we create two evaluators for the b-spline, one for the horizontal and
+  // one for the vertical gradient. The derivatives for a b-spline are requested
+  // by passing a TinyVector with as many elements as the spline's dimension
+  // with the desired derivative degree for each dimension. Here we want the
+  // first derivative in x and y direction:
+  
+  vigra::TinyVector < float , 2 > dx1_spec { 1 , 0 } ;
+  vigra::TinyVector < float , 2 > dy1_spec { 0 , 1 } ;
+  
+  // now we can construct the evaluators
+  
+  ev_type xev ( bspl , dx1_spec ) ;
+  ev_type yev ( bspl , dy1_spec ) ;
+  
+  // this is where the result should go:
+  target_type target ( imageInfo.shape() ) ;
+
+  // quick-shot solution, iterating in a loop, not vectorized
+  
+  auto start = vigra::createCoupledIterator ( target ) ;
+  auto end = start.getEndIterator() ;
+  
+  for ( auto it = start ; it < end ; ++it )
+  {
+    // we fetch the discrete coordinate from the coupled iterator
+    // and instantiate a coordinate_type from it. Note that we can't pass
+    // the discrete coordinate directly to the evaluator's operator()
+    // because this fails to be disambiguated.
+    
+    coordinate_type crd ( it.get<0>() ) ;
+    
+    // now we get the two gradients by evaluating the gradient evaluators
+    // at the given coordinate
+    
+    pixel_type dx = xev ( crd ) ;
+    pixel_type dy = yev ( crd ) ;
+    
+    // and conclude by writing the sum of the squared gradients to target
+
+    it.get<1>() = dx * dx + dy * dy ;
+  }
+  
+  // store the result with vigra impex
+  vigra::ImageExportInfo eximageInfo ( "gsm.tif" );
+  
+  vigra::exportImage ( target ,
+                       eximageInfo
+                       .setPixelType("UINT8") ) ;
+  
+  exit ( 0 ) ;
+}
diff --git a/example/gsm2.cc b/example/gsm2.cc
new file mode 100644
index 0000000..6612d80
--- /dev/null
+++ b/example/gsm2.cc
@@ -0,0 +1,216 @@
+/************************************************************************/
+/*                                                                      */
+/*    vspline - a set of generic tools for creation and evaluation      */
+/*              of uniform b-splines                                    */
+/*                                                                      */
+/*            Copyright 2015, 2016 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.                                   */
+/*                                                                      */
+/************************************************************************/
+
+/// gsm2.cc
+///
+/// alternative implementation of gsm.cc, performing the calculation of the
+/// gradient squared magnitude with a functor and index_remap, which is faster since
+/// the whole operation is multithreaded and potentially vectorized.
+///
+/// compile with:
+/// clang++ -std=c++11 -march=native -o gsm -O3 -pthread -DUSE_VC=1 gsm.cc -lvigraimpex -lVc
+/// g++ also works. And, as ever, you can compile without USE_VC, if you haven't got Vc
+/// or you don't want to use it.
+
+#include <vspline/vspline.h>
+
+#include <vigra/stdimage.hxx>
+#include <vigra/imageinfo.hxx>
+#include <vigra/impex.hxx>
+
+#ifdef USE_VC
+const int VSIZE = vspline::vector_traits < float > :: size ;
+#else
+const int VSIZE = 1 ;
+#endif
+
+// we silently assume we have a colour image
+typedef vigra::RGBValue<float,0,1,2> pixel_type; 
+
+// coordinate_type has a 2D coordinate
+typedef vigra::TinyVector < float , 2 > coordinate_type ;
+
+// type of b-spline object
+typedef vspline::bspline < pixel_type , 2 > spline_type ;
+
+// target_type is a 2D array of pixels  
+typedef vigra::MultiArray < 2 , pixel_type > target_type ;
+
+// b-spline evaluator producing float pixels
+typedef vspline::evaluator < coordinate_type , // incoming coordinate's type
+                             pixel_type ,      // singular result data type
+                             -1 , -1 ,         // no specializations
+                             VSIZE             // vector size
+                           > ev_type ;
+
+/// we build a vspline::unary_functor which calculates the gradient squared
+/// magnitude. The code here isn't generic, for a general-purpose gsm evaluator
+/// one would need a bit more of a coding effort, but we want to demonstrate
+/// the principle here.
+/// Note how the 'compound evaluator' we construct follows a pattern of
+/// - derive from vspline::unary_functor
+/// - keep const references to 'inner' types
+/// - pass these in the constructor, yielding a 'pure' functor
+/// - if the vector code is identical to the unvectorized code, implement
+///   eval() with a template                      
+
+struct ev_gsm
+: public vspline::unary_functor < coordinate_type , // incoming coordinate
+                                  pixel_type ,      // desired result type
+                                  VSIZE >
+{
+  // pull in vspline::unary_functor's type system
+  typedef vspline::unary_functor < coordinate_type ,
+                                   pixel_type ,
+                                   ev_type::vsize > base_type ;
+                                   
+  using_unary_functor_types(base_type) ;
+  
+  // we create two evaluators for the b-spline, one for the horizontal and
+  // one for the vertical gradient. The derivatives for a b-spline are requested
+  // by passing a TinyVector with as many elements as the spline's dimension
+  // with the desired derivative degree for each dimension. Here we want the
+  // first derivative in x and y direction:
+
+  const vigra::TinyVector < float , 2 > dx1_spec { 1 , 0 } ;
+  const vigra::TinyVector < float , 2 > dy1_spec { 0 , 1 } ;
+
+  // we keep two 'inner' evaluators, one for each direction
+  
+  const ev_type xev , yev ;
+  
+  // which are initialized in the constructor, using the bspline and the
+  // derivative specifiers
+  
+  ev_gsm ( const spline_type & bspl )
+  : xev ( bspl , dx1_spec ) ,
+    yev ( bspl , dy1_spec )
+  { } ;
+  
+  // since the code is the same for vectorized and unvectorized
+  // operation, we can write a template:
+  
+  template < class IN , class OUT = base_type::out_type_of<IN> >
+  void _eval ( const IN & c ,
+                     OUT & result ) const
+  {
+    auto dx = xev ( c ) ; // get the gradient in x direction
+    auto dy = yev ( c ) ; // get the gradient in y direction
+    
+    // TODO: really, we'd like to write:
+    // result = dx * dx + dy * dy ;
+    // but fail due to problems with type inference, so we need to be
+    // a bit more explicit:
+    
+    dx *= dx ;            // square the gradients
+    dy *= dy ;
+    dx += dy ;
+    
+    result = dx ;         // assign to result
+  } 
+  
+  // and implement the pure virtual method(s) like this:
+  void eval ( const in_type & c ,
+                    out_type & result ) const
+  {
+    _eval ( c , result ) ;
+  } ;
+
+#ifdef USE_VC
+  
+  void eval ( const in_v & c ,
+                    out_v & result ) const
+  {
+    _eval ( c , result ) ;
+  } ;
+
+#endif
+} ;
+
+int main ( int argc , char * argv[] )
+{
+  // get the image file name
+  
+  vigra::ImageImportInfo imageInfo ( argv[1] ) ;
+
+  // we want a b-spline with natural boundary conditions
+  
+  vigra::TinyVector < vspline::bc_code , 2 > bcv ( vspline::NATURAL ) ;
+  
+  // create cubic 2D b-spline object containing the image data
+  
+  spline_type bspl ( imageInfo.shape() , // the shape of the data for the spline
+                     3 ,                 // degree 3 == cubic spline
+                     bcv                 // specifies natural BCs along both axes
+                   ) ;
+  
+  // load the image data into the b-spline's core. This is a common idiom:
+  // the spline's 'core' is a MultiArrayView to that part of the spline's
+  // data container which is meant to hold the input data. This saves loading
+  // the image to some memory first and then transferring the data into
+  // the spline. Since the core is a vigra::MultiarrayView, we can pass it
+  // to importImage as the desired target for loading the image from disk.
+                   
+  vigra::importImage ( imageInfo , bspl.core ) ;
+  
+  // prefilter the b-spline
+
+  bspl.prefilter() ;
+  
+  // now we can construct the gsm evaluator
+  
+  ev_gsm ev ( bspl ) ;
+  
+  // this is where the result should go:
+  
+  target_type target ( imageInfo.shape() ) ;
+
+  // now we obtain the result by performing an index_remap. index_remap
+  // successively passes discrete indices into the target to the evaluator
+  // it's invoked with, storing the result of the evaluator's evaluation
+  // at the self-same coordinates. This is done multithreaded and vectorized
+  // automatically, so it's very convenient, if the evaluator is at hand.
+  // So here we have invested moderately more coding effort in the evaluator
+  // and are rewarded with being able to use the evaluator with vspline's
+  // high-level code for a very fast implementation of our gsm problem.
+  
+  vspline::index_remap < ev_gsm , 2 > ( ev , target ) ;
+
+  // store the result with vigra impex
+
+  vigra::ImageExportInfo eximageInfo ( "gsm2.tif" );
+  
+  vigra::exportImage ( target ,
+                       eximageInfo
+                       .setPixelType("UINT8") ) ;
+  
+  exit ( 0 ) ;
+}
diff --git a/example/slice.cc b/example/slice.cc
index 6047db0..d3e19de 100644
--- a/example/slice.cc
+++ b/example/slice.cc
@@ -68,6 +68,10 @@ int main ( int argc , char * argv[] )
   vspline::bspline < voxel_type , 3 >
     space ( vigra::Shape3 ( 10 , 10 , 10 ) , 5 , bcv ) ;
   
+  auto red_channel = space.get_channel_view ( 0 ) ;
+  auto green_channel = space.get_channel_view ( 1 ) ;
+  auto blue_channel = space.get_channel_view ( 2 ) ;
+
   // fill the b-spline's core with a three-way gradient
   for ( int z = 0 ; z < 10 ; z++ )
   {
diff --git a/example/splinus.cc b/example/splinus.cc
index 360dd6e..e6b6a78 100644
--- a/example/splinus.cc
+++ b/example/splinus.cc
@@ -50,6 +50,8 @@
 
 int main ( int argc , char * argv[] )
 {
+  assert ( argc > 1 ) ;
+  
   int degree = std::atoi ( argv[1] ) ;
   
   assert ( degree >= 0 && degree < 25 ) ;
diff --git a/unary_functor.h b/unary_functor.h
index 5ec29ef..b2f4a89 100644
--- a/unary_functor.h
+++ b/unary_functor.h
@@ -139,6 +139,42 @@ struct unary_functor
    virtual void eval ( const in_v & c ,             // simdized nD input
                        out_v & result ) const = 0 ; // simdized mD output
 
+  /// for single-channel data, we want a variant evaluating to an out_ele_v
+  /// instead of an out_v, because in this case out_v is a TinyVector
+  /// with one element, which is awkward and unexpected.
+
+  void eval ( const in_v & input ,        // number of dimensions * coordinate vectors
+              out_ele_v & result )  const // single value vector
+  {
+    static_assert ( dim_out == 1 , "this evaluation routine is for single-channel data only" ) ;
+    out_v helper ;
+    eval ( input , helper ) ;
+    result = helper[0] ;
+  }
+
+  /// ditto, for 1D coordinates
+
+  void eval ( const in_ele_v & input , // single coordinate vector
+              out_v & result ) const   // simdized mD output
+  {
+    static_assert ( dim_in == 1 , "this evaluation routine is for 1D data only" ) ;
+    in_v helper ;
+    helper[0] = input ;
+    eval ( helper , result ) ;
+  }
+
+  void eval ( const in_ele_v & input ,    // single coordinate vector
+              out_ele_v & result )  const // single value vector
+  {
+    static_assert ( dim_in == 1 && dim_out == 1  ,
+                    "this evaluation routine is for single-channel, 1D data only") ;
+    in_v in_helper ;
+    in_helper[0] = input ;
+    out_v out_helper ;
+    eval ( in_helper , out_helper ) ;
+    result = out_helper[0] ;
+  }
+
   /// here I provide an (ineffecient) method to evaluate vectorized data by
   /// delegating to the unvectorized evaluation routine. This will work, but
   /// it's slow compared to proper vector code.

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/vspline.git



More information about the debian-science-commits mailing list