Visualization of the NARF descriptor and descriptor distances

This tutorial is about the visualization of how the NARF descriptor is calculated and to test how the descriptor distances between certain points in a range image behave. Compared to the other tuturials, this one is not really about the code, but about trying the program and looking at the visualization. Of course, nothing keeps you from having a look at it anyway.

The code

First, create a file called, let’s say, narf_descriptor_visualization.cpp in your favorite editor, and place the following code inside it:

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
/* \author Bastian Steder */

#include <iostream>

#include <pcl/point_cloud.h>
#include <pcl/io/pcd_io.h>
#include <pcl/visualization/range_image_visualizer.h>
#include <pcl/range_image/range_image.h>
#include <pcl/features/narf.h>
#include <pcl/console/parse.h>
#include <pcl/common/file_io.h> // for getFilenameWithoutExtension

float angular_resolution = 0.5f;
int rotation_invariant = 0;
float support_size = 0.3f;
int descriptor_size = 36;
pcl::RangeImage::CoordinateFrame coordinate_frame = pcl::RangeImage::CAMERA_FRAME;
bool setUnseenToMaxRange = false;

typedef pcl::PointXYZ PointType;

void 
printUsage (const char* progName)
{
  std::cout << "\n\nUsage: "<<progName<<" [options] <scene.pcd>\n\n"
            << "Options:\n"
            << "-------------------------------------------\n"
            << "-r <float>   angular resolution in degrees (default "<<angular_resolution<<")\n"
            << "-s <float>   support size for the interest points (diameter of the used sphere - "
            <<                                                     "default "<<support_size<<")\n"
            << "-d <int>     descriptor size (default "<<descriptor_size<<")\n"
            << "-c <int>     coordinate frame of the input point cloud (default "<< (int)coordinate_frame<<")\n"
            << "-o <0/1>     switch rotational invariant version of the feature on/off"
            <<               " (default "<< (int)rotation_invariant<<")\n"
            << "-m           set unseen pixels to max range\n"
            << "-h           this help\n"
            << "\n\n";
}

int 
main (int argc, char** argv)
{
  // --------------------------------------
  // -----Parse Command Line Arguments-----
  // --------------------------------------
  if (pcl::console::find_argument (argc, argv, "-h") >= 0)
  {
    printUsage (argv[0]);
    return 0;
  }
  if (pcl::console::find_argument (argc, argv, "-m") >= 0)
  {
    setUnseenToMaxRange = true;
    std::cout << "Setting unseen values in range image to maximum range readings.\n";
  }
  if (pcl::console::parse (argc, argv, "-o", rotation_invariant) >= 0)
    std::cout << "Switching rotation invariant feature version "<< (rotation_invariant ? "on" : "off")<<".\n";
  int tmp_coordinate_frame;
  if (pcl::console::parse (argc, argv, "-c", tmp_coordinate_frame) >= 0)
  {
    coordinate_frame = pcl::RangeImage::CoordinateFrame (tmp_coordinate_frame);
    std::cout << "Using coordinate frame "<< (int)coordinate_frame<<".\n";
  }
  if (pcl::console::parse (argc, argv, "-s", support_size) >= 0)
    std::cout << "Setting support size to "<<support_size<<".\n";
  if (pcl::console::parse (argc, argv, "-d", descriptor_size) >= 0)
    std::cout << "Setting descriptor size to "<<descriptor_size<<".\n";
  if (pcl::console::parse (argc, argv, "-r", angular_resolution) >= 0)
    std::cout << "Setting angular resolution to "<<angular_resolution<<"deg.\n";
  angular_resolution = pcl::deg2rad (angular_resolution);
  

  // -----------------------
  // -----Read pcd file-----
  // -----------------------
  pcl::PointCloud<PointType>::Ptr point_cloud_ptr (new pcl::PointCloud<PointType>);
  pcl::PointCloud<PointType>& point_cloud = *point_cloud_ptr;
  pcl::PointCloud<pcl::PointWithViewpoint> far_ranges;
  Eigen::Affine3f scene_sensor_pose (Eigen::Affine3f::Identity ());
  std::vector<int> pcd_filename_indices = pcl::console::parse_file_extension_argument (argc, argv, "pcd");
  if (!pcd_filename_indices.empty ())
  {
    std::string filename = argv[pcd_filename_indices[0]];
    if (pcl::io::loadPCDFile (filename, point_cloud) == -1)
    {
      std::cout << "Was not able to open file \""<<filename<<"\".\n";
      printUsage (argv[0]);
      return 0;
    }
    scene_sensor_pose = Eigen::Affine3f (Eigen::Translation3f (point_cloud.sensor_origin_[0],
                                                               point_cloud.sensor_origin_[1],
                                                               point_cloud.sensor_origin_[2])) *
                        Eigen::Affine3f (point_cloud.sensor_orientation_);
    std::string far_ranges_filename = pcl::getFilenameWithoutExtension (filename)+"_far_ranges.pcd";
    if (pcl::io::loadPCDFile (far_ranges_filename.c_str (), far_ranges) == -1)
      std::cout << "Far ranges file \""<<far_ranges_filename<<"\" does not exists.\n";
  }
  else
  {
    std::cout << "\nNo *.pcd file for scene given.\n\n";
    printUsage (argv[0]);
    return 1;
  }
  
  // -----------------------------------------------
  // -----Create RangeImage from the PointCloud-----
  // -----------------------------------------------
  float noise_level = 0.0;
  float min_range = 0.0f;
  int border_size = 1;
  pcl::RangeImage::Ptr range_image_ptr (new pcl::RangeImage);
  pcl::RangeImage& range_image = *range_image_ptr;
  range_image.createFromPointCloud (point_cloud, angular_resolution, pcl::deg2rad (360.0f), pcl::deg2rad (180.0f),
                                   scene_sensor_pose, coordinate_frame, noise_level, min_range, border_size);
  range_image.integrateFarRanges (far_ranges);
  if (setUnseenToMaxRange)
    range_image.setUnseenToMaxRange ();
  
  // Extract NARF features:
  std::cout << "Now extracting NARFs in every image point.\n";
  std::vector<std::vector<pcl::Narf*> > narfs;
  narfs.resize (range_image.size ());
  int last_percentage=-1;
  for (unsigned int y=0; y<range_image.height; ++y)
  {
    for (unsigned int x=0; x<range_image.width; ++x)
    {
      const auto index = y*range_image.width+x;
      const auto percentage = ((100*index) / range_image.size ());
      if (percentage > last_percentage)
      {
        std::cout << percentage<<"% "<<std::flush;
        last_percentage = percentage;
      }
      pcl::Narf::extractFromRangeImageAndAddToList (range_image, x, y, descriptor_size,
                                                    support_size, rotation_invariant != 0, narfs[index]);
      //std::cout << "Extracted "<<narfs[index].size ()<<" features for pixel "<<x<<","<<y<<".\n";
    }
  }
  std::cout << "100%\n";
  std::cout << "Done.\n\n Now you can click on points in the image to visualize how the descriptor is "
       << "extracted and see the descriptor distances to every other point..\n";
  
  //---------------------
  // -----Show image-----
  // --------------------
  pcl::visualization::RangeImageVisualizer range_image_widget ("Scene range image"),
                                           surface_patch_widget("Descriptor's surface patch"),
                                           descriptor_widget("Descriptor"),
                                           descriptor_distances_widget("descriptor distances");
  range_image_widget.showRangeImage (range_image);
  //range_image_widget.visualize_selected_point = true;

  //--------------------
  // -----Main loop-----
  //--------------------
  while (true) 
  {
    range_image_widget.spinOnce ();  // process GUI events
    surface_patch_widget.spinOnce ();  // process GUI events
    descriptor_widget.spinOnce ();  // process GUI events
    pcl_sleep(0.01);
    
    //if (!range_image_widget.mouse_click_happened)
      continue;
    //range_image_widget.mouse_click_happened = false;
    //float clicked_pixel_x_f = range_image_widget.last_clicked_point_x,
          //clicked_pixel_y_f = range_image_widget.last_clicked_point_y;
    int clicked_pixel_x, clicked_pixel_y;
    //range_image.real2DToInt2D (clicked_pixel_x_f, clicked_pixel_y_f, clicked_pixel_x, clicked_pixel_y);
    if (!range_image.isValid (clicked_pixel_x, clicked_pixel_y))
      continue;
      //Vector3f clicked_3d_point;
      //range_image.getPoint (clicked_pixel_x, clicked_pixel_y, clicked_3d_point);
    
    //surface_patch_widget.show (false);
    //descriptor_widget.show (false);"
    
    int selected_index = clicked_pixel_y*range_image.width + clicked_pixel_x;
    pcl::Narf narf;
    if (!narf.extractFromRangeImage (range_image, clicked_pixel_x, clicked_pixel_y,
                                                                         descriptor_size, support_size))
    {
      continue;
    }
    
    int surface_patch_pixel_size = narf.getSurfacePatchPixelSize ();
    float surface_patch_world_size = narf.getSurfacePatchWorldSize ();
    surface_patch_widget.showFloatImage (narf.getSurfacePatch (), surface_patch_pixel_size, surface_patch_pixel_size,
                                         -0.5f*surface_patch_world_size, 0.5f*surface_patch_world_size, true);
    float surface_patch_rotation = narf.getSurfacePatchRotation ();
    float patch_middle = 0.5f* (float (surface_patch_pixel_size-1));
    float angle_step_size = pcl::deg2rad (360.0f)/narf.getDescriptorSize ();
    float cell_size = surface_patch_world_size/float (surface_patch_pixel_size),
          cell_factor = 1.0f/cell_size,
          max_dist = 0.5f*surface_patch_world_size,
          line_length = cell_factor* (max_dist-0.5f*cell_size);
    for (int descriptor_value_idx=0; descriptor_value_idx<narf.getDescriptorSize (); ++descriptor_value_idx)
    {
      float angle = descriptor_value_idx*angle_step_size + surface_patch_rotation;
      //surface_patch_widget.markLine (patch_middle, patch_middle, patch_middle+line_length*sinf (angle),
                                     //patch_middle+line_length*-std::cos (angle), pcl::visualization::Vector3ub (0,255,0));
    }
    std::vector<float> rotations, strengths;
    narf.getRotations (rotations, strengths);
    float radius = 0.5f*surface_patch_pixel_size;
    for (unsigned int i=0; i<rotations.size (); ++i)
    {
      //surface_patch_widget.markLine (radius-0.5, radius-0.5, radius-0.5f + 2.0f*radius*sinf (rotations[i]),
                                                //radius-0.5f - 2.0f*radius*std::cos (rotations[i]), pcl::visualization::Vector3ub (255,0,0));
    }
    
    descriptor_widget.showFloatImage (narf.getDescriptor (), narf.getDescriptorSize (), 1, -0.1f, 0.3f, true);

    //===================================
    //=====Compare with all features=====
    //===================================
    const std::vector<pcl::Narf*>& narfs_of_selected_point = narfs[selected_index];
    if (narfs_of_selected_point.empty ())
      continue;
    
    //descriptor_distances_widget.show (false);
    float* descriptor_distance_image = new float[range_image.size ()];
    for (unsigned int point_index=0; point_index<range_image.size (); ++point_index)
    {
      float& descriptor_distance = descriptor_distance_image[point_index];
      descriptor_distance = std::numeric_limits<float>::infinity ();
      std::vector<pcl::Narf*>& narfs_of_current_point = narfs[point_index];
      if (narfs_of_current_point.empty ())
        continue;
      for (unsigned int i=0; i<narfs_of_selected_point.size (); ++i)
      {
        for (unsigned int j=0; j<narfs_of_current_point.size (); ++j)
        {
          descriptor_distance = (std::min)(descriptor_distance,
                                           narfs_of_selected_point[i]->getDescriptorDistance (*narfs_of_current_point[j]));
        }
      }
    }
    descriptor_distances_widget.showFloatImage (descriptor_distance_image, range_image.width, range_image.height,
                                               -std::numeric_limits<float>::infinity (), std::numeric_limits<float>::infinity (), true);
    delete[] descriptor_distance_image;
  }
}

Compiling and running the program

Add the following lines to your CMakeLists.txt file:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
cmake_minimum_required(VERSION 2.6 FATAL_ERROR)

project(narf_descriptor_visualization)

find_package(PCL 1.3 REQUIRED)

include_directories(${PCL_INCLUDE_DIRS})
link_directories(${PCL_LIBRARY_DIRS})
add_definitions(${PCL_DEFINITIONS})

add_executable (narf_descriptor_visualization narf_descriptor_visualization.cpp)
target_link_libraries (narf_descriptor_visualization ${PCL_LIBRARIES})

You can now try it with a point cloud file from your hard drive:

$ ./narf_descriptor_visualization <point_cloud.pcd>

It will take a few second, during which you will see the status in the terminal. During this time, a NARF feature is extracted in every point of the range image created from the given point cloud. When it is done, a widget showing the range image pops up. Now click on a point in the range image. If it is a valid image point, three additional widgets will pop up. One visualizing the actual descriptor as a row of gray values, one showing a local range image patch of the area on which you clicked, overlaid with a star shaped pattern. Each beam corresponds to one of the cells in the descriptor. The one facing upwards to the first cell and then going clockwise. The basic intuition is, that the more the surface changes under the beam, the higher (brighter) the value of the corresponding descriptor cell. There is also one or more red beams, which mark the extracted dominant orientations of the image patch, which, together with the normal, is used to create a unique orientation for the feature coordinate frame. The last image visualizes the descriptor distances to every other point in the scene. The darker the value, the more similar the point is to the clicked image point.

The result should look similar to this:

_images/narf_descriptor_visualization.png

Also have a look at:

$ ./narf_descriptor_visualization -h

for a list of parameters.