Distance Transformation Labels/Features

Distance transformation functions on most of the libraries (OpenCV, VXL etc.) only output the distance transformation. However, must of the time we also need the label/position of the shorted distance point along with the shortest distance for a given point, which is often referred as distance transformation feature in computer vision and machine learning. Following is the simplest implementation of distance transformation feature in C++  using OpenCV library.

   
 //============================================================================  
 // Name    : OpenCVTest.cpp  
 // Author   : Rudra Poudel  
 // Description : Example of distance transformation features. To find the shorted distance label we travel shorted distance path i.e. kind of gradient descent.  
 //============================================================================  
   
 #include "opencv2/core/core.hpp"  
 #include "opencv2/imgproc/imgproc.hpp"  
 #include "opencv2/calib3d/calib3d.hpp"  
 #include "opencv2/highgui/highgui.hpp"  
 #include <stdio.h>  
 #include <iostream>  
   
 #define VIDEO_FRAME_HEIGHT 480  
 #define VIDEO_FRAME_WIDTH 640  

 void VerifyDistTransformLabels(cv::Mat& depth_map, cv::Mat& labels, int row_start, int row_end, int col_start, int col_end) {  
  for( int y = row_start; y < row_end; y++ ) {  
   //src_data_row = dist_trans.ptr<float>(y);  
   for( int x = col_start; x < col_end; x++ ) {  
   
    cv::Vec2i alabel = labels.at<cv::Vec2i>(y,x);  
    assert(alabel[0] != -1);
    assert(alabel[1] != -1);
   
   }  
  }  
 }  

 cv::Vec2i GetLowestDistNeighbour(cv::Mat& dist_trans, int x, int y){  
  int x_start = x>0 ? x-1 : 0;  
  int y_start = y>0 ? y-1 : 0;  
   
  int x_end = x<(dist_trans.cols -1) ? x+1 : x;  
  int y_end = y<(dist_trans.rows -1) ? y+1 : y;  
   
  cv::Vec2i p; //p[0] = x; p[1] = y; assignment not necessary  
  float min_dist = dist_trans.at<float>(y, x);  
  float dist;  
  for( int row = y_start; row <= y_end; row++ ) {  
   for( int col = x_start; col <= x_end; col++ ) {  
    dist = dist_trans.at<float>(row, col);  
    if(dist<min_dist) {  
     min_dist = dist;  
     p[0] = col;  
     p[1] = row;  
    }  
   }  
  }  
   
  return p;  
 }  

 cv::Vec2i GetDistTransformLabel(cv::Mat& dist_trans, cv::Mat& labels, int x, int y) {  
   
  if(dist_trans.at<float>(y,x) == 0) {  
   cv::Vec2i alabel(x,y);  
   return alabel;  
  } else {  
   cv::Vec2i p = GetLowestDistNeighbour(dist_trans, x, y);  
   cv::Vec2i& alabel = labels.at<cv::Vec2i>(p[1], p[0]); // (y, x)  
   if(alabel[0] == -1 )  
    alabel = GetDistTransformLabel(dist_trans, labels, p[0], p[1]);  
   
   return alabel;  
  }  
 }  

 void GetDistTransformLabels(cv::Mat& dist_trans, cv::Mat& labels, int row_start, int row_end, int col_start, int col_end) {  
   
  assert(dist_trans.type() == CV_32FC1);  
  assert(labels.type() == CV_32SC2);  
   
  labels.setTo(-1); //float* src_data_row;  
  for( int y = row_start; y < row_end; y++ ) {  
   //src_data_row = dist_trans.ptr<float>(y);  
   for( int x = col_start; x < col_end; x++ ) {  
   
    cv::Vec2i& alabel = labels.at<cv::Vec2i>(y,x);  
    if(alabel[0] == -1 ) {  
     alabel = GetDistTransformLabel(dist_trans, labels, x, y);  
    }  
   
   }  
  }  
 }  

 int main( int argc, char** argv )  
 {  
  cv::Mat depth_map(VIDEO_FRAME_HEIGHT, VIDEO_FRAME_WIDTH, CV_8UC1, 1);  
  cv::Mat dist_trans(VIDEO_FRAME_HEIGHT, VIDEO_FRAME_WIDTH, CV_32FC1);  
  cv::Mat labels(VIDEO_FRAME_HEIGHT, VIDEO_FRAME_WIDTH, CV_32SC2);  
   
  cv::circle(depth_map, cv::Point2i(320, 240), 40, cvScalar(255,255,255), 10);  
  cv::MatExpr mask_fg = depth_map > 1;  
  cv::MatExpr mask_bg = depth_map == 1;  
  depth_map.setTo(0, mask_fg);  
  depth_map.setTo(255, mask_bg);  
   
  cv::distanceTransform(depth_map, dist_trans, CV_DIST_L2, CV_DIST_MASK_PRECISE);  
  GetDistTransformLabels(dist_trans, labels, 0, dist_trans.rows, 0, dist_trans.cols);  
  VerifyDistTransformLabels(depth_map, labels, 0, labels.rows, 0, labels.cols);  
   
  double min_value, max_value;  
  cv::minMaxLoc(dist_trans, &min_value, &max_value);  
  dist_trans *=(1/max_value);  
  std::cout<<"\nmin_value: "<<min_value<<"\tmax_value: "<<max_value;  
   
  std::cout<<"\n";  
  //cv::imshow("Dist transform for", depth_map);  
  cv::imshow("Distance Transformation", dist_trans);  
  cv::waitKey(40000);  
   
  return 0;  
 }