dvcm-2d-curvature.cpp

#include <iostream>
#include "DGtal/base/Common.h"
#include "DGtal/helpers/StdDefs.h"
#include "DGtal/topology/SetOfSurfels.h"
#include "DGtal/topology/DigitalSurface.h"
#include "DGtal/geometry/curves/GridCurve.h"
#include "DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.h"
#include "DGtal/geometry/volumes/estimation/VoronoiCovarianceMeasure.h"
#include "DGtal/geometry/surfaces/estimation/VoronoiCovarianceMeasureOnDigitalSurface.h"
#include "DGtal/geometry/surfaces/estimation/VCMDigitalSurfaceLocalEstimator.h"
#include "DGtal/io/colormaps/GradientColorMap.h"
#include "DGtal/io/boards/Board2D.h"
#include "ConfigExamples.h"



using namespace std;
using namespace DGtal;

int main( int /* argc */, char** /* argv */ )
{
  typedef Z2i::Space Space;
  typedef Z2i::KSpace KSpace;
  typedef Z2i::Curve Curve;
  typedef Z2i::Point Point;
  typedef ExactPredicateLpSeparableMetric<Space, 2> Metric; // L2-metric

  typedef Curve::SCellsRange LinelRange;
  typedef SetOfSurfels<KSpace> DigitalSurfaceContainer;
  typedef DigitalSurface<DigitalSurfaceContainer> Surface;
  typedef functors::BallConstantPointFunction<Point,double> KernelFunction;
  typedef VoronoiCovarianceMeasureOnDigitalSurface<DigitalSurfaceContainer,Metric,KernelFunction> VCMOnSurface;
  typedef functors::VCMAbsoluteCurvatureFunctor<VCMOnSurface> CurvatureFunctor;
  typedef VCMDigitalSurfaceLocalEstimator<DigitalSurfaceContainer,Metric,
                                          KernelFunction, CurvatureFunctor> VCMCurvatureEstimator;

  // Transform a 4-connected sequence of discrete points into a digital surface.
  //string inputSDP = examplesPath + "samples/flower-30-8-3.sdp";
  //string inputSDP = examplesPath + "samples/ellipse-20-7-0.4.sdp";
  //string inputSDP = examplesPath + "samples/accflower-20-5-5-0.1.sdp";
  string inputSDP = examplesPath + "samples/circle-43.sdp";
  trace.info() << "Reading input 2d discrete points file: " << inputSDP; 
  KSpace ks; ks.init( Point( -1000, -1000 ), Point( 1000, 1000 ), true );
  Curve curve( ks ); 
  fstream inputStream( inputSDP.c_str(), ios::in);
  curve.initFromVectorStream(inputStream);
  inputStream.close();
  DigitalSurfaceContainer* container 
    = new DigitalSurfaceContainer( ks,  SurfelAdjacency<KSpace::dimension>( true ) );
  LinelRange range = curve.getSCellsRange();
  for ( LinelRange::ConstIterator it = range.begin(), itE = range.end(); it != itE; ++it )
    container->surfelSet().insert( *it );
  CountedConstPtrOrConstPtr<Surface> ptrSurface( new Surface( container ) ); // acquired
  trace.info() << " [done] " << std::endl ; 
  const double R = 40;
  trace.info() << "Big radius   R = " << R << std::endl;
  const double r = 20;
  trace.info() << "Small radius r = " << r << std::endl;
  const double T = 0.2;
  trace.info() << "Curvature thres. T = " << T << std::endl; // threshold for displaying features as red.

  Metric l2;
  KernelFunction chi( 1.0, r );
  VCMCurvatureEstimator estimator( ptrSurface );
  estimator.setParams( Pointels, R, r, chi, 3.0, l2, true );
  estimator.init( 1.0, ptrSurface->begin(), ptrSurface->end() );

  // Flat zones are blue, more curved zones are yellow till red.
  GradientColorMap<double> colormap( 0, T );
  colormap.addColor( Color::Blue );
  colormap.addColor( Color::Yellow );
  colormap.addColor( Color::Red );
  Board2D board;
  for ( Surface::ConstIterator it = ptrSurface->begin(), itE = ptrSurface->end();
        it != itE; ++it )
    {
      double curv = estimator.eval( it );
      curv = std::min( T, curv );
      board << CustomStyle( it->className(), 
                            new CustomColors( colormap( curv ),  colormap( curv ) ) )
            << *it;
      std::cout << curv << std::endl;
    }      
  board.saveSVG("dvcm-curvature.svg");

  return 0;
}
//                                                                           //