clean-up, code adaptation towards matlab-suitable data sizes, code commentations, stable version

segment-graph.h

@@ -35,67 +35,59 @@ bool operator<(const edge &a, const edge &b) {
   return a.w < b.w;
 }
 
-/*
- * Segment a graph
- *
- * Returns a disjoint-set forest representing the segmentation.
- *
- * num_vertices: number of vertices in graph.
- * num_edges: number of edges in graph
- * edges: array of edges.
- * c: constant for treshold function.
+ /**
+ * @brief Segment a graph, a disjoint-set forest representing the segmentation
+ * @author Pedro Felzenszwalb, Alexander Freytag
+ * @date 27-03-2014 ( dd-mm-yyyy, last updated)
+ * 
+ * @param[in] i_numVertices: number of vertices in graph
+ * @param[in] i_numEdges: number of edges in graph
+ * @param[in] edges: array of edges
+ * @param[in] c: constant for treshold function
+ * 
+ * @param[out] u universe (disjoint-set forest representing the segmentation)
  */
-universe *segment_graph(int num_vertices, int num_edges, edge *edges, 
-			float c) { 
+universe *segment_graph( const int i_numVertices, 
+                         const int i_numEdges, 
+                         edge *edges, 
+                         const float c
+                       )
+{ 
   // sort edges by weight
   // note: if weights occure more then once, this might lead to non-deterministic (non-reproducable) behaviour, since
   // "Equivalent elements are not guaranteed to keep their original relative order"
-//   std::sort(edges, edges + num_edges);
-  
-/*   for (int i = 0; i < 5; i++) {
-      int a = edges[i].a;
-      int b = edges[i].b;
-      mexPrintf("edges before %d: %d <-> %d with weight %f\n",i, a, b, edges[i].w);      
-    } */  
+  //  std::sort(edges, edges + i_numEdges);
   
   // adaptation: use stable_sort instead, which will keep the relative position of equal elements :)
-  std::stable_sort(edges, edges + num_edges);
+  std::stable_sort(edges, edges + i_numEdges);
 
   // make a disjoint-set forest
-  universe *u = new universe(num_vertices);
+  universe *u = new universe(i_numVertices);
 
   // init thresholds
-  float *threshold = new float[num_vertices];
-  for (int i = 0; i < num_vertices; i++)
+  float *threshold = new float[i_numVertices];
+  for (int i = 0; i < i_numVertices; i++)
+  {
     threshold[i] = THRESHOLD(1,c);
+  }
   
-    //that's the same as with find, and also non-reproducable
-//    for (int i = 0; i < 5; i++) {
-//       int a = edges[i].a;
-//       int b = edges[i].b;
-//       mexPrintf("edges after %d: %d <-> %d with weight %f\n",i, a, b, edges[i].w);    
-//     } 
-  
-//      for (int i = 0; i < 5; i++) {
-//       int a = u->find(edges[i].a);
-//       int b = u->find(edges[i].b);
-//       mexPrintf("components %d: %d <-> %d\n",i, a, b);      
-//     } 
-  // this is already non-reproducable
+
 
   // for each edge, in non-decreasing weight order...
-  for (int i = 0; i < num_edges; i++) {
+  for (int i = 0; i < i_numEdges; i++)
+  {
     edge *pedge = &edges[i];
     
     // components conected by this edge
     int a = u->find(pedge->a);
     int b = u->find(pedge->b);
-    if (a != b) {
-      if ((pedge->w <= threshold[a]) &&
-	  (pedge->w <= threshold[b])) {
-	u->join(a, b);
-	a = u->find(a);
-	threshold[a] = pedge->w + THRESHOLD(u->size(a), c);
+    if (a != b)
+    {
+      if ( (pedge->w <= threshold[a]) && (pedge->w <= threshold[b]) )
+      {
+        u->join(a, b);
+        a = u->find(a);
+        threshold[a] = pedge->w + THRESHOLD(u->size(a), c);
       }
     }
   }

segment-image-labelOutput.h

@@ -13,24 +13,24 @@
 #include <mex.h>
 
  /**
- * @brief Segment an image, returning an int-image giving the index of the corresponding segment for every pixel
+ * @brief Segment an image, returning an ushort-image giving the index of the corresponding segment for every pixel
  * @author Alexander Freytag
  * @date 27-03-2014 ( dd-mm-yyyy, last updated)
  * 
  * @param[in] im: image to segment
- * @param[in] sigma: to smooth the image
+ * @param[in] d_sigma: to smooth the image
  * @param[in] c: constant for treshold function
- * @param[in] min_size: minimum component size (enforced by post-processing stage)
+ * @param[in] i_minSize: minimum component size (enforced by post-processing stage)
  * @param[in] num_ccs: number of connected components in the segmentation
  * 
- * @param[out] output int-image giving the index of the corresponding segment for every pixel
+ * @param[out] output ushort-image giving the index of the corresponding segment for every pixel
  */
-image<int> *segment_image_labelOutput(   image<rgb> *im,
-                                         float sigma,
-                                         float c,
-                                         int min_size,
-                                         int *num_ccs
-                                     )
+image<unsigned short> *segment_image_labelOutput(   image<rgb> *im,
+                                                    const float d_sigma,
+                                                    const float c,
+                                                    const int i_minSize,
+                                                    int *num_ccs
+                                                )
 {
   int width = im->width();
   int height = im->height();
@@ -39,7 +39,7 @@ image<int> *segment_image_labelOutput(   image<rgb> *im,
   image<float> *g = new image<float>(width, height);
   image<float> *b = new image<float>(width, height);
 
-  // smooth each color channel  
+  // pre-processing: smooth each color channel  
   for (int y = 0; y < height; y++)
   {
     for (int x = 0; x < width; x++)
@@ -49,16 +49,17 @@ image<int> *segment_image_labelOutput(   image<rgb> *im,
       imRef(b, x, y) = imRef(im, x, y).b;
     }
   }
-  image<float> *smooth_r = smooth(r, sigma);
-  image<float> *smooth_g = smooth(g, sigma);
-  image<float> *smooth_b = smooth(b, sigma);
+  image<float> *smooth_r = smooth(r, d_sigma);
+  image<float> *smooth_g = smooth(g, d_sigma);
+  image<float> *smooth_b = smooth(b, d_sigma);
   delete r;
   delete g;
   delete b;
  
   // build graph
   edge *edges = new edge[width*height*4];
-  int num = 0;
+  int num ( 0 );
+  
   for (int y = 0; y < height; y++)
   {
     for (int x = 0; x < width; x++)
@@ -108,7 +109,7 @@ image<int> *segment_image_labelOutput(   image<rgb> *im,
   {
     int a = u->find(edges[i].a);
     int b = u->find(edges[i].b);
-    if ( (a != b) && ((u->size(a) < min_size) || (u->size(b) < min_size))  )
+    if ( (a != b) && ((u->size(a) < i_minSize) || (u->size(b) < i_minSize))  )
       u->join(a, b);
   }
   delete [] edges;
@@ -117,11 +118,11 @@ image<int> *segment_image_labelOutput(   image<rgb> *im,
   // max_int segments ( and therefore also not imgs with more pixels)
   *num_ccs = u->num_sets();
   
-  image<int> *output = new image<int>(width, height);
+  image<unsigned short> *output = new image<unsigned short>(width, height);
  
   //how many different regions do we finally have, and which of them are the corresponding indices?
-  std::map<int,int> regionLabels;
-  int idx ( 0 );
+  std::map<int,unsigned short> regionLabels;
+  unsigned short idx ( 0 );
   
   for (int y = 0; y < height; y++)
   {
@@ -130,7 +131,7 @@ image<int> *segment_image_labelOutput(   image<rgb> *im,
       int comp = u->find(y * width + x);
       if ( regionLabels.find( comp ) == regionLabels.end() )
       {
-        regionLabels.insert( std::pair<int,int>(comp,idx) );
+        regionLabels.insert( std::pair<int,unsigned short>(comp,idx) );
         idx++;
       }
     }

segment-image.h

@@ -38,26 +38,43 @@ rgb random_rgb(){
 }
 
 // dissimilarity measure between pixels
-static inline float diff(image<float> *r, image<float> *g, image<float> *b,
-			 int x1, int y1, int x2, int y2) {
+static inline float diff(   image<float> *r, 
+                            image<float> *g, 
+                            image<float> *b,
+                            int x1, 
+                            int y1, 
+                            int x2, 
+                            int y2
+                        )
+{
   return sqrt(square(imRef(r, x1, y1)-imRef(r, x2, y2)) +
-	      square(imRef(g, x1, y1)-imRef(g, x2, y2)) +
-	      square(imRef(b, x1, y1)-imRef(b, x2, y2)));
+              square(imRef(g, x1, y1)-imRef(g, x2, y2)) +
+              square(imRef(b, x1, y1)-imRef(b, x2, y2))
+             );
 }
 
-/*
- * Segment an image
- *
- * Returns a color image representing the segmentation.
- *
- * im: image to segment.
- * sigma: to smooth the image.
- * c: constant for treshold function.
- * min_size: minimum component size (enforced by post-processing stage).
- * num_ccs: number of connected components in the segmentation.
+ /**
+ * @brief Segment an image, returning a color image representing the segmentation
+ * @author Pedro Felzenszwalb, Alexander Freytag
+ * @date 27-03-2014 ( dd-mm-yyyy, last updated)
+ * 
+ * @param[in] im: image to segment
+ * @param[in] d_sigma: to smooth the image
+ * @param[in] c: constant for treshold function
+ * @param[in] i_minSize: minimum component size (enforced by post-processing stage)
+ * @param[in] num_ccs: number of connected components in the segmentation
+ * 
+ * @param[out] output rgb-image (ushort per channel) visualizing the resulting segments
  */
-image<rgb> *segment_image(image<rgb> *im, float sigma, float c, int min_size,
-			  int *num_ccs) {
+ 
+
+image<rgb> *segment_image(  image<rgb> *im, 
+                            const float d_sigma, 
+                            const float c, 
+                            const int i_minSize,
+                            int *num_ccs
+                         )
+{
   int width = im->width();
   int height = im->height();
 
@@ -66,16 +83,19 @@ image<rgb> *segment_image(image<rgb> *im, float sigma, float c, int min_size,
   image<float> *b = new image<float>(width, height);
 
   // smooth each color channel  
-  for (int y = 0; y < height; y++) {
-    for (int x = 0; x < width; x++) {
+  for (int y = 0; y < height; y++)
+  {
+    for (int x = 0; x < width; x++)
+    {
       imRef(r, x, y) = imRef(im, x, y).r;
       imRef(g, x, y) = imRef(im, x, y).g;
       imRef(b, x, y) = imRef(im, x, y).b;
     }
   }
-  image<float> *smooth_r = smooth(r, sigma);
-  image<float> *smooth_g = smooth(g, sigma);
-  image<float> *smooth_b = smooth(b, sigma);
+  
+  image<float> *smooth_r = smooth(r, d_sigma);
+  image<float> *smooth_g = smooth(g, d_sigma);
+  image<float> *smooth_b = smooth(b, d_sigma);
   delete r;
   delete g;
   delete b;
@@ -83,34 +103,40 @@ image<rgb> *segment_image(image<rgb> *im, float sigma, float c, int min_size,
   // build graph
   edge *edges = new edge[width*height*4];
   int num = 0;
-  for (int y = 0; y < height; y++) {
-    for (int x = 0; x < width; x++) {
-      if (x < width-1) {
-	edges[num].a = y * width + x;
-	edges[num].b = y * width + (x+1);
-	edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y);
-	num++;
+  for (int y = 0; y < height; y++)
+  {
+    for (int x = 0; x < width; x++)
+    {
+      if (x < width-1)
+      {
+        edges[num].a = y * width + x;
+        edges[num].b = y * width + (x+1);
+        edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y);
+        num++;
       }
 
-      if (y < height-1) {
-	edges[num].a = y * width + x;
-	edges[num].b = (y+1) * width + x;
-	edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x, y+1);
-	num++;
+      if (y < height-1)
+      {
+        edges[num].a = y * width + x;
+        edges[num].b = (y+1) * width + x;
+        edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x, y+1);
+        num++;
       }
 
-      if ((x < width-1) && (y < height-1)) {
-	edges[num].a = y * width + x;
-	edges[num].b = (y+1) * width + (x+1);
-	edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y+1);
-	num++;
+      if ((x < width-1) && (y < height-1)) 
+      {
+        edges[num].a = y * width + x;
+        edges[num].b = (y+1) * width + (x+1);
+        edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y+1);
+        num++;
       }
 
-      if ((x < width-1) && (y > 0)) {
-	edges[num].a = y * width + x;
-	edges[num].b = (y-1) * width + (x+1);
-	edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y-1);
-	num++;
+      if ((x < width-1) && (y > 0))
+      {
+        edges[num].a = y * width + x;
+        edges[num].b = (y-1) * width + (x+1);
+        edges[num].w = diff(smooth_r, smooth_g, smooth_b, x, y, x+1, y-1);
+        num++;
       }
     }
   }
@@ -122,10 +148,11 @@ image<rgb> *segment_image(image<rgb> *im, float sigma, float c, int min_size,
   universe *u = segment_graph(width*height, num, edges, c);
   
   // post process small components
-  for (int i = 0; i < num; i++) {
+  for (int i = 0; i < num; i++)
+  {
     int a = u->find(edges[i].a);
     int b = u->find(edges[i].b);
-    if ((a != b) && ((u->size(a) < min_size) || (u->size(b) < min_size)))
+    if ((a != b) && ((u->size(a) < i_minSize) || (u->size(b) < i_minSize)))
       u->join(a, b);
   }
   delete [] edges;
@@ -136,10 +163,14 @@ image<rgb> *segment_image(image<rgb> *im, float sigma, float c, int min_size,
   // pick random colors for each component
   rgb *colors = new rgb[width*height];
   for (int i = 0; i < width*height; i++)
+  {
     colors[i] = random_rgb();
+  }
   
-  for (int y = 0; y < height; y++) {
-    for (int x = 0; x < width; x++) {
+  for (int y = 0; y < height; y++)
+  {
+    for (int x = 0; x < width; x++)
+    {
       int comp = u->find(y * width + x);
       imRef(output, x, y) = colors[comp];
     }

segmentFelzenszwalb.cpp

@@ -3,6 +3,7 @@
 #include <stdlib.h>
 #include <mex.h>
 
+
 #include "./image.h"
 #include "./misc.h"
 #include "./pnmfile.h"
@@ -45,15 +46,15 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
   *                      Get the arguments
   * -------------------------------------------------------------- */  
   
-  bool verbose ( false );
+  bool b_verbose ( false );
   if ( nInput >= 7)
   {
     if ( mxIsLogicalScalar( pInput[6] ) && mxIsLogicalScalarTrue( pInput[6] ) )
     {
-      verbose = true;
+      b_verbose = true;
     }
     else
-      verbose = false; 
+      b_verbose = false; 
   }
    
   /* Get string of input image*/  
@@ -65,7 +66,7 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
     char *input;     
     input = (char *) mxCalloc( mxGetN(pInput[0] )+1, sizeof(char) );
     mxGetString(pInput[0], input, mxGetN( pInput[0] )+1);
-    if ( verbose ) 
+    if ( b_verbose ) 
       mexPrintf("The input string is:  %s\n", input);    
     
     imgInput = loadPPM( input );
@@ -102,57 +103,57 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
   char buf [1024] ;
   
   /* Get sigma */  
-  double sigma;
+  double d_sigma;
   if ( nInput < 2) 
   {
-    sigma = 0.5; 
+    d_sigma = 0.5; 
   }
   else
-    sigma = mxGetScalar( pInput[1] ); 
-  if ( verbose ) 
-    mexPrintf("The sigma value is:  %f\n", sigma);
+    d_sigma = mxGetScalar( pInput[1] ); 
+  if ( b_verbose ) 
+    mexPrintf("The sigma value is:  %f\n", d_sigma);
     
   /* Get k */  
-  int k;
+  int i_k;
   if ( nInput < 3) 
   {
-    k = 500;
+    i_k = 500;
   }
   else
-    k = mxGetScalar( pInput[2] );
+    i_k = mxGetScalar( pInput[2] );
   
-  if ( verbose ) 
-    mexPrintf("The k is:  %i\n", k);
+  if ( b_verbose ) 
+    mexPrintf("The k is:  %i\n", i_k);
    
   /* Get minSize*/  
-  int minSize;  
+  int i_minSize;  
   if ( nInput < 4)
   {
-    minSize = 50;
+    i_minSize = 50;
   }
   else
-    minSize = mxGetScalar(  pInput[3] );
+    i_minSize = mxGetScalar(  pInput[3] );
   
-  if ( verbose ) 
-    mexPrintf("The minSize is:  %i\n", minSize);
+  if ( b_verbose ) 
+    mexPrintf("The minSize is:  %i\n", i_minSize);
    
   /* Get bool whether to compute the label img (int) or colored img (rgb)*/  
-  bool computeColorOutput;
+  bool b_computeColorOutput;
   if (nInput >=  5)
   {
     if ( mxIsLogicalScalar( pInput[4] ) && mxIsLogicalScalarTrue( pInput[4] ) )
     {
-      computeColorOutput = true;
+      b_computeColorOutput = true;
     }
     else
-      computeColorOutput = false;   
+      b_computeColorOutput = false;   
   }
   else
   {
-    computeColorOutput = false;
+    b_computeColorOutput = false;
   }
-  if ( verbose ) 
-    mexPrintf("To we compute RGB colored segmentation result?  :  %i\n", computeColorOutput );
+  if ( b_verbose ) 
+    mexPrintf("To we compute RGB colored segmentation result?  :  %i\n", b_computeColorOutput );
    
   /* Get string of output image if given*/  
   char * output;
@@ -161,11 +162,11 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
     output = (char *) mxCalloc(mxGetN(pInput[5])+1, sizeof(char));
     mxGetString(pInput[5], output, mxGetN(pInput[5])+1);
     
-    if ( verbose ) 
+    if ( b_verbose ) 
       mexPrintf("The output string is:  %s\n", output);
   }
     
-  if ( verbose )
+  if ( b_verbose )
     mexPrintf("image loaded, now start segmentation\n");  
   
    
@@ -177,18 +178,18 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
    
   int num_ccs; 
   image<rgb> * imgResultRGB = NULL;
-  image<int> * imgResult  = NULL;
+  image<unsigned short> * imgResult  = NULL;
    
-  if ( computeColorOutput )
+  if ( b_computeColorOutput )
   {
-    imgResultRGB = segment_image(imgInput, sigma, k, minSize, &num_ccs); 
+    imgResultRGB = segment_image(imgInput, d_sigma, i_k, i_minSize, &num_ccs); 
   }
   else
   {
-    imgResult = segment_image_labelOutput(imgInput, sigma, k, minSize, &num_ccs); 
+    imgResult = segment_image_labelOutput(imgInput, d_sigma, i_k, i_minSize, &num_ccs); 
   }
        
-  if ( verbose )
+  if ( b_verbose )
     mexPrintf("segmentation done\n");   
    
    
@@ -201,7 +202,7 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
     //was a filename given to store the image in?
     if (nInput >=  6) 
     {
-      if ( computeColorOutput )
+      if ( b_computeColorOutput )
         savePPM( imgResultRGB, output );
       else
         save_image( imgResult, output );
@@ -214,24 +215,24 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
   }
   else
   {
-    if ( verbose )
+    if ( b_verbose )
       mexPrintf("convert to matlab structure and hand result back to main program\n"); 
      
     /* convert segmentation result to matlab matrix*/
         
-    if ( computeColorOutput )
+    if ( b_computeColorOutput )
     {
       int width ( imgResultRGB->width() );
       int height (imgResultRGB->height() );       
       /* keep in mind that matlab stores images  height × width × color, whereas C does it the other way round*/
       int dims[] = {height, width,3};
-      pOutput[0] = mxCreateNumericArray (3, dims, mxINT32_CLASS, mxREAL);
+      pOutput[0] = mxCreateNumericArray (3, dims, mxUINT8_CLASS, mxREAL);
       //unsigned char *out1; /* pointer to output 1 */
       //out1 = (unsigned char *)mxGetPr( pOutput[0] ); /* pointer to output 1 */
       //unsigned char *out1; /* pointer to output 1 */
       //out1 = (unsigned char *)mxGetPr( pOutput[0] ); /* pointer to output 1 */
-      int *out1; /* pointer to output 1 */
-      out1 = (int *)mxGetData( pOutput[0] ); /* pointer to output 1 *	
+      unsigned char *out1; /* pointer to output 1 */
+      out1 = (unsigned char *)mxGetData( pOutput[0] ); /* pointer to output 1 *	
 
               
       /* start with RED channel*/     
@@ -241,7 +242,7 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
         uint rowOffset ( x*height );
         for ( uint y = 0; y < height; y++)
         {
-          out1[rowOffset + y ] = (int) (imRef(imgResultRGB, x, y)).r;
+          out1[rowOffset + y ] = (unsigned char) (imRef(imgResultRGB, x, y)).r;
         }
       }
         
@@ -252,7 +253,7 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
         uint rowOffset ( x*height );
           for ( uint y = 0; y < height; y++)
           {
-            out1[channelOffsetG + rowOffset + y ] = (int) (imRef(imgResultRGB, x, y)).g;
+            out1[channelOffsetG + rowOffset + y ] = (unsigned char) (imRef(imgResultRGB, x, y)).g;
           }
       }
         
@@ -263,11 +264,11 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
         uint rowOffset ( x*height );
           for ( uint y = 0; y < height; y++)
           {
-            out1[channelOffsetB + rowOffset + y ] = (int) (imRef(imgResultRGB, x, y)).b;
+            out1[channelOffsetB + rowOffset + y ] = (unsigned char) (imRef(imgResultRGB, x, y)).b;
           }
       }
     }
-    else /* do not compute colored rgb segmentation image, but only an int img*/
+    else /* do not compute colored rgb segmentation image, but only an unsigned short-int img*/
     {
       int width ( imgResult->width() );
       int height (imgResult->height() ); 
@@ -285,15 +286,15 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
           uint rowOffset ( x*height );
           for ( uint y = 0; y < height; y++)
           {
-        out1[rowOffset + y ] = (unsigned char ) (imRef(imgResult, x, y));
+            out1[rowOffset + y ] = (unsigned char) (imRef(imgResult, x, y));
           }
         }
       }
       else
       {
-        pOutput[0] = mxCreateNumericArray (2, dims, mxINT32_CLASS, mxREAL);
-        int *out1; /* pointer to output 1 */
-        out1 = (int *)mxGetData( pOutput[0] ); /* pointer to output 1 *
+        pOutput[0] = mxCreateNumericArray (2, dims, mxUINT16_CLASS, mxREAL);
+        unsigned short *out1; /* pointer to output 1 */
+        out1 = (unsigned short *)mxGetData( pOutput[0] ); /* pointer to output 1 *
         
         /* keep in mind that matlab stores images col by width, whereas C does it the other way round*/
         for ( uint x = 0; x < width; x++)
@@ -301,7 +302,7 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
           uint rowOffset ( x*height );
           for ( uint y = 0; y < height; y++)
           {
-            out1[rowOffset + y ] = (int) (imRef(imgResult, x, y));
+            out1[rowOffset + y ] = (unsigned short) (imRef(imgResult, x, y));
           }
         }
       }
@@ -321,17 +322,17 @@ void mexFunction(int nOutput, mxArray *pOutput[], /* Output variables */
       }
       else
       {
-        pOutput[1] = mxCreateNumericArray(1,dims,mxINT32_CLASS,mxREAL);
+        pOutput[1] = mxCreateNumericArray(1,dims,mxINT16_CLASS,mxREAL);
         //unsigned char *out2; /* pointer to output 2 */
         //out2 = (unsigned char *)mxGetPr( pOutput[1] ); /* pointer to output 2 */
-        int *out2; /* pointer to output 2 */
-        out2 = (int *) mxGetData( pOutput[1] ); /* pointer to output 2 */
+        unsigned short *out2; /* pointer to output 2 */
+        out2 = (unsigned short *) mxGetData( pOutput[1] ); /* pointer to output 2 */
           
         out2[0] = num_ccs;
 
       }
        
-      if ( verbose ) 
+      if ( b_verbose ) 
         mexPrintf("number of components: %i", num_ccs);
     }