Merge pull request #14 from jhlegarreta/FixExtractContourWithSnakesExampleStyle

STYLE: Conform to ITK Style Guide.

Author: thewtex

ImageSegmentation/ExtractContourWithSnakes.cxx

@@ -3,15 +3,15 @@
 #include "itkGradientRecursiveGaussianImageFilter.h"
 #include "itkGradientMagnitudeImageFilter.h"
 #include "itkImageFileReader.h"
- 
+
 #include "vnl/vnl_math.h"
 #include <vnl/vnl_matrix.h>
 #include "vnl/algo/vnl_determinant.h"
 #include "vnl/algo/vnl_matrix_inverse.h"
 #include <vnl/vnl_vector.h>
- 
+
 #include <iostream>
- 
+
 namespace
 {
 using ImageType = itk::Image< unsigned char, 2 >;
@@ -25,7 +25,7 @@ using GradMagfilterType = itk::GradientMagnitudeImageFilter<
  ImageType, FloatImageType >;
 using ReaderType = itk::ImageFileReader< ImageType >;
 }
- 
+
 vnl_vector<double> generateCircle(double cx, double cy,
  double rx, double ry, int n);
 void createImage(ImageType::Pointer image,
@@ -35,27 +35,27 @@ vnl_matrix<double> computeP(double alpha, double beta, double gamma,
 vnl_vector<double> sampleImage(vnl_vector<double> x, vnl_vector<double> y,
  OutputImageType::Pointer gradient,
  int position);
- 
+
 int main( int argc, char* argv[] )
 {
- //Image dimensions
+ // Image dimensions
  int w = 300;
  int h = 300;
  ImageType::Pointer image;
- if (argc < 7)
+ if( argc < 7 )
  {
  std::cout << "Usage " << argv[0]
  << " points alpha beta gamma sigma iterations [image]"
  << std::endl;
- return EXIT_SUCCESS;;
+ return EXIT_SUCCESS;
  }
- else if (argc < 8)
+ else if( argc < 8 )
  {
  //Synthesize the image
  image = ImageType::New();
  createImage(image, w, h, 150, 150, 50, 50);
  }
- else if (argc == 8)
+ else if( argc == 8 )
  {
  //Open the image
  ReaderType::Pointer reader = ReaderType::New();
@@ -73,200 +73,200 @@ int main( int argc, char* argv[] )
  return EXIT_FAILURE;
  }
  }
- 
- //Snake parameters
+
+ // Snake parameters
  double alpha = 0.001;
  double beta = 0.4;
  double gamma = 100;
  double iterations = 1;
  int nPoints = 20;
  double sigma;
- 
+
  nPoints = atoi(argv[1]);
  alpha = atof(argv[2]);
  beta = atof(argv[3]);
  gamma = atof(argv[4]);
  sigma = atof(argv[5]);
  iterations = atoi(argv[6]);
- 
- //Temporal variables
+
+ // Temporal variables
  vnl_matrix<double> P;
  vnl_vector<double> v;
  double N;
- 
- //Generate initial snake circle
+
+ // Generate initial snake circle
  v = generateCircle(130, 130, 50, 50, nPoints);
 
- //Computes P matrix.
+ // Compute P matrix.
  N = v.size()/2;
  try
  {
  P = computeP(alpha, beta, gamma, N);
  }
- catch (int n)
+ catch(...)
  {
- return EXIT_FAILURE;;
+ return EXIT_FAILURE;
  }
- 
- //Computes the magnitude gradient
+
+ // Compute the magnitude gradient
  GradMagfilterType::Pointer gradientMagnitudeFilter =
  GradMagfilterType::New();
  gradientMagnitudeFilter->SetInput( image );
  gradientMagnitudeFilter->Update();
- 
- //Computes the gradient of the gradient magnitude
+
+ // Compute the gradient of the gradient magnitude
  FilterType::Pointer gradientFilter = FilterType::New();
  gradientFilter->SetInput( gradientMagnitudeFilter->GetOutput() );
  gradientFilter->SetSigma( sigma );
  gradientFilter->Update();
- 
- //Loop
+
+ // Loop
  vnl_vector<double> x(N);
  vnl_vector<double> y(N);
- 
+
  std::cout << "Initial snake" << std::endl;
- for (int i = 0; i < N; i++)
+ for (int i = 0; i < N; ++i )
  {
  x[i] = v[2*i];
- y[i] = v[2*i+1];
+ y[i] = v[2*i + 1];
  std::cout << "(" << x[i] << ", " << y[i] << ")" << std::endl;
  }
- 
- for (int i = 0; i < iterations; i++)
+
+ for( int i = 0; i < iterations; ++i )
  {
  vnl_vector<double> fex;
  vnl_vector<double> fey;
  fex = sampleImage(x, y, gradientFilter->GetOutput(), 0);
  fey = sampleImage(x, y, gradientFilter->GetOutput(), 1);
- 
+
  x = (x+gamma*fex).post_multiply(P);
  y = (y+gamma*fey).post_multiply(P);
  }
- 
- //Display the answer
+
+ // Display the answer
  std::cout << "Final snake after " << iterations << " iterations" << std::endl;
  vnl_vector<double> v2(2*N);
- for (int i=0; i<N; i++)
+ for( int i = 0; i < N; ++i )
  {
  v2[2*i] = x[i];
- v2[2*i+1] = y[i];
+ v2[2*i + 1] = y[i];
  std::cout << "(" << x[i] << ", " << y[i] << ")" << std::endl;
  }
- 
- return EXIT_SUCCESS;;
+
+ return EXIT_SUCCESS;
 }
- 
+
 vnl_vector<double> generateCircle(
  double cx, double cy, double rx, double ry, int n)
 {
  vnl_vector<double> v(2*(n+1));
- 
- for (int i=0; i<n; i++)
+
+ for( int i = 0; i < n; ++i )
  {
  v[2*i] = cx + rx*cos(2*vnl_math::pi*i/n);
- v[2*i+1] = cy + ry*sin(2*vnl_math::pi*i/n);
+ v[2*i + 1] = cy + ry*sin(2*vnl_math::pi*i/n);
  }
- v[2*n]=v[0];
- v[2*n+1]=v[1];
+ v[2*n] = v[0];
+ v[2*n + 1] = v[1];
  return v;
 }
 
 void createImage(ImageType::Pointer image,
  int w, int h, double cx, double cy, double rx, double ry)
 {
- 
  itk::Size<2> size;
  size[0] = w;
  size[1] = h;
- 
+
  itk::RandomImageSource<ImageType>::Pointer randomImageSource = itk::RandomImageSource<ImageType>::New();
  randomImageSource->SetNumberOfThreads(1); // to produce non-random results
  randomImageSource->SetSize(size);
  randomImageSource->SetMin(200);
  randomImageSource->SetMax(255);
  randomImageSource->Update();
- 
+
  image->SetRegions(randomImageSource->GetOutput()->GetLargestPossibleRegion());
  image->Allocate();
- 
+
  IndexType index;
- 
- //Draw oval.
- for (int i=0; i<w; i++)
+
+ // Draw oval
+ for( int i = 0; i < w; ++i )
  {
- for (int j=0; j<h; j++)
+ for(int j = 0; j < h; ++j )
  {
- index[0] = i; index[1] = j;
- if ( ((i-cx)*(i-cx)/(rx*rx) + (j-cy)*(j-cy)/(ry*ry) ) < 1)
+ index[0] = i;
+ index[1] = j;
+ if( ((i - cx)*(i - cx)/(rx*rx) + (j - cy)*(j - cy)/(ry*ry) ) < 1 )
  {
  image->SetPixel(index, randomImageSource->GetOutput()->GetPixel(index)-100);
  }
  else
  {
  image->SetPixel(index, randomImageSource->GetOutput()->GetPixel(index));
  }
- 
+
  }
  }
 }
- 
+
 vnl_matrix<double> computeP(double alpha, double beta, double gamma, double N) throw (int)
 {
- 
- double a = gamma*(2*alpha+6*beta)+1;
- double b = gamma*(-alpha-4*beta);
+
+ double a = gamma*(2*alpha + 6*beta)+1;
+ double b = gamma*(-alpha - 4*beta);
  double c = gamma*beta;
- 
- vnl_matrix<double> P(N,N);
- 
+
+ vnl_matrix<double> P(N, N);
+
  P.fill(0);
- 
- //fill diagonal
+
+ // Fill diagonal
  P.fill_diagonal(a);
- 
- //fill next two diagonals
- for (int i=0; i<(N-1); i++)
+
+ // Fill next two diagonals
+ for( int i = 0; i < N - 1; ++i )
  {
- P(i+1,i) = b;
- P(i,i+1) = b;
+ P(i + 1, i) = b;
+ P(i, i + 1) = b;
  }
- //Moreover
- P(0, N-1)=b;
- P(N-1, 0)=b;
- 
- //fill next two diagonals
- for (int i=0; i<(N-2); i++)
+ // Moreover
+ P(0, N-1) = b;
+ P(N-1, 0) = b;
+
+ // Fill next two diagonals
+ for( int i = 0; i < N- 2; ++i )
  {
- P(i+2,i) = c;
- P(i,i+2) = c;
+ P(i + 2, i) = c;
+ P(i, i + 2) = c;
  }
- //Moreover
- P(0, N-2)=c;
- P(1, N-1)=c;
- P(N-2, 0)=c;
- P(N-1, 1)=c;
- 
- if ( vnl_determinant(P) == 0.0 )
+ // Moreover
+ P(0, N - 2) = c;
+ P(1, N - 1) = c;
+ P(N - 2, 0) = c;
+ P(N - 1, 1) = c;
+
+ if( vnl_determinant(P) == 0.0 )
  {
  std::cerr << "Singular matrix. Determinant is 0." << std::endl;
  throw 2;
  }
- 
- //Compute the inverse of the matrix P
+
+ // Compute the inverse of the matrix P
  vnl_matrix< double > Pinv;
  Pinv = vnl_matrix_inverse< double >(P);
- 
+
  return Pinv.transpose();
 }
- 
+
 vnl_vector<double> sampleImage(vnl_vector<double> x, vnl_vector<double> y, OutputImageType::Pointer gradient, int position)
 {
  int size;
  size = x.size();
  vnl_vector<double> ans(size);
- 
+
  IndexType index;
- for (int i=0; i<size; i++)
+ for( int i = 0; i < size; ++i )
  {
  index[0] = x[i];
  index[1] = y[i];