voxel_grid.h

Go to the documentation of this file.

/*
 * Software License Agreement (BSD License)
 *
 *  Point Cloud Library (PCL) - www.pointclouds.org
 *  Copyright (c) 2010-2011, Willow Garage, Inc.
 *
 *  All rights reserved.
 *
 *  Redistribution and use in source and binary forms, with or without
 *  modification, are permitted provided that the following conditions
 *  are met:
 *
 *   * Redistributions of source code must retain the above copyright
 *     notice, this list of conditions and the following disclaimer.
 *   * Redistributions in binary form must reproduce the above
 *     copyright notice, this list of conditions and the following
 *     disclaimer in the documentation and/or other materials provided
 *     with the distribution.
 *   * Neither the name of Willow Garage, Inc. nor the names of its
 *     contributors may be used to endorse or promote products derived
 *     from this software without specific prior written permission.
 *
 *  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 *  "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 *  LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
 *  FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
 *  COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
 *  INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
 *  BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
 *  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
 *  CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
 *  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
 *  ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 *  POSSIBILITY OF SUCH DAMAGE.
 *
 * $Id: voxel_grid.h 5654 2012-05-01 05:32:03Z rusu $
 *
 */

#ifndef PCL_FILTERS_VOXEL_GRID_MAP_H_
#define PCL_FILTERS_VOXEL_GRID_MAP_H_

#include <pcl/filters/filter.h>
#include <map>
#include <boost/unordered_map.hpp>
#include <boost/fusion/sequence/intrinsic/at_key.hpp>

namespace pcl
{
  PCL_EXPORTS void 
  getMinMax3D (const sensor_msgs::PointCloud2ConstPtr &cloud, int x_idx, int y_idx, int z_idx, 
               Eigen::Vector4f &min_pt, Eigen::Vector4f &max_pt);

  PCL_EXPORTS void 
  getMinMax3D (const sensor_msgs::PointCloud2ConstPtr &cloud, int x_idx, int y_idx, int z_idx, 
               const std::string &distance_field_name, float min_distance, float max_distance, 
               Eigen::Vector4f &min_pt, Eigen::Vector4f &max_pt, bool limit_negative = false);

  inline Eigen::MatrixXi
  getHalfNeighborCellIndices ()
  {
    Eigen::MatrixXi relative_coordinates (3, 13);
    int idx = 0;

    // 0 - 8
    for (int i = -1; i < 2; i++) 
    {
      for (int j = -1; j < 2; j++) 
      {
        relative_coordinates (0, idx) = i;
        relative_coordinates (1, idx) = j;
        relative_coordinates (2, idx) = -1;
        idx++;
      }
    }
    // 9 - 11
    for (int i = -1; i < 2; i++) 
    {
      relative_coordinates (0, idx) = i;
      relative_coordinates (1, idx) = -1;
      relative_coordinates (2, idx) = 0;
      idx++;
    }
    // 12
    relative_coordinates (0, idx) = -1;
    relative_coordinates (1, idx) = 0;
    relative_coordinates (2, idx) = 0;

    return (relative_coordinates);
  }

  inline Eigen::MatrixXi
  getAllNeighborCellIndices ()
  {
    Eigen::MatrixXi relative_coordinates = getHalfNeighborCellIndices ();
    Eigen::MatrixXi relative_coordinates_all( 3, 26);
    relative_coordinates_all.block<3, 13> (0, 0) = relative_coordinates;
    relative_coordinates_all.block<3, 13> (0, 13) = -relative_coordinates;
    return (relative_coordinates_all);
  }

  template <typename PointT> void 
  getMinMax3D (const typename pcl::PointCloud<PointT>::ConstPtr &cloud, 
               const std::string &distance_field_name, float min_distance, float max_distance,
               Eigen::Vector4f &min_pt, Eigen::Vector4f &max_pt, bool limit_negative = false);

  template <typename PointT>
  class VoxelGrid: public Filter<PointT>
  {
    protected:
      using Filter<PointT>::filter_name_;
      using Filter<PointT>::getClassName;
      using Filter<PointT>::input_;
      using Filter<PointT>::indices_;

      typedef typename Filter<PointT>::PointCloud PointCloud;
      typedef typename PointCloud::Ptr PointCloudPtr;
      typedef typename PointCloud::ConstPtr PointCloudConstPtr;

    public:
      VoxelGrid () : 
        leaf_size_ (Eigen::Vector4f::Zero ()),
        inverse_leaf_size_ (Eigen::Array4f::Zero ()),
        downsample_all_data_ (true), 
        save_leaf_layout_ (false),
        leaf_layout_ (),
        min_b_ (Eigen::Vector4i::Zero ()),
        max_b_ (Eigen::Vector4i::Zero ()),
        div_b_ (Eigen::Vector4i::Zero ()),
        divb_mul_ (Eigen::Vector4i::Zero ()),
        filter_field_name_ (""), 
        filter_limit_min_ (-FLT_MAX), 
        filter_limit_max_ (FLT_MAX),
        filter_limit_negative_ (false)
      {
        filter_name_ = "VoxelGrid";
      }

      virtual ~VoxelGrid ()
      {
      }

      inline void 
      setLeafSize (const Eigen::Vector4f &leaf_size) 
      { 
        leaf_size_ = leaf_size; 
        // Avoid division errors
        if (leaf_size_[3] == 0)
          leaf_size_[3] = 1;
        // Use multiplications instead of divisions
        inverse_leaf_size_ = Eigen::Array4f::Ones () / leaf_size_.array ();
      }

      inline void
      setLeafSize (float lx, float ly, float lz)
      {
        leaf_size_[0] = lx; leaf_size_[1] = ly; leaf_size_[2] = lz;
        // Avoid division errors
        if (leaf_size_[3] == 0)
          leaf_size_[3] = 1;
        // Use multiplications instead of divisions
        inverse_leaf_size_ = Eigen::Array4f::Ones () / leaf_size_.array ();
      }

      inline Eigen::Vector3f 
      getLeafSize () { return (leaf_size_.head<3> ()); }

      inline void 
      setDownsampleAllData (bool downsample) { downsample_all_data_ = downsample; }

      inline bool 
      getDownsampleAllData () { return (downsample_all_data_); }

      inline void 
      setSaveLeafLayout (bool save_leaf_layout) { save_leaf_layout_ = save_leaf_layout; }

      inline bool 
      getSaveLeafLayout () { return (save_leaf_layout_); }

      inline Eigen::Vector3i 
      getMinBoxCoordinates () { return (min_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getMaxBoxCoordinates () { return (max_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getNrDivisions () { return (div_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getDivisionMultiplier () { return (divb_mul_.head<3> ()); }

      inline int 
      getCentroidIndex (const PointT &p)
      {
        return (leaf_layout_.at ((Eigen::Vector4i (static_cast<int> (floor (p.x * inverse_leaf_size_[0])), 
                                                   static_cast<int> (floor (p.y * inverse_leaf_size_[1])), 
                                                   static_cast<int> (floor (p.z * inverse_leaf_size_[2])), 0) - min_b_).dot (divb_mul_)));
      }

      inline std::vector<int> 
      getNeighborCentroidIndices (const PointT &reference_point, const Eigen::MatrixXi &relative_coordinates)
      {
        Eigen::Vector4i ijk (static_cast<int> (floor (reference_point.x * inverse_leaf_size_[0])), 
                             static_cast<int> (floor (reference_point.y * inverse_leaf_size_[1])), 
                             static_cast<int> (floor (reference_point.z * inverse_leaf_size_[2])), 0);
        Eigen::Array4i diff2min = min_b_ - ijk;
        Eigen::Array4i diff2max = max_b_ - ijk;
        std::vector<int> neighbors (relative_coordinates.cols());
        for (int ni = 0; ni < relative_coordinates.cols (); ni++)
        {
          Eigen::Vector4i displacement = (Eigen::Vector4i() << relative_coordinates.col(ni), 0).finished();
          // checking if the specified cell is in the grid
          if ((diff2min <= displacement.array()).all() && (diff2max >= displacement.array()).all())
            neighbors[ni] = leaf_layout_[((ijk + displacement - min_b_).dot (divb_mul_))]; // .at() can be omitted
          else
            neighbors[ni] = -1; // cell is out of bounds, consider it empty
        }
        return (neighbors);
      }

      inline std::vector<int> 
      getLeafLayout () { return (leaf_layout_); }

      inline Eigen::Vector3i 
      getGridCoordinates (float x, float y, float z) 
      { 
        return (Eigen::Vector3i (static_cast<int> (floor (x * inverse_leaf_size_[0])), 
                                 static_cast<int> (floor (y * inverse_leaf_size_[1])), 
                                 static_cast<int> (floor (z * inverse_leaf_size_[2])))); 
      }

      inline int 
      getCentroidIndexAt (const Eigen::Vector3i &ijk)
      {
        int idx = ((Eigen::Vector4i() << ijk, 0).finished() - min_b_).dot (divb_mul_);
        if (idx < 0 || idx >= static_cast<int> (leaf_layout_.size ())) // this checks also if leaf_layout_.size () == 0 i.e. everything was computed as needed
        {
          //if (verbose)
          //  PCL_ERROR ("[pcl::%s::getCentroidIndexAt] Specified coordinate is outside grid bounds, or leaf layout is not saved, make sure to call setSaveLeafLayout(true) and filter(output) first!\n", getClassName ().c_str ());
          return (-1);
        }
        return (leaf_layout_[idx]);
      }

      inline void
      setFilterFieldName (const std::string &field_name)
      {
        filter_field_name_ = field_name;
      }

      inline std::string const
      getFilterFieldName ()
      {
        return (filter_field_name_);
      }

      inline void
      setFilterLimits (const double &limit_min, const double &limit_max)
      {
        filter_limit_min_ = limit_min;
        filter_limit_max_ = limit_max;
      }

      inline void
      getFilterLimits (double &limit_min, double &limit_max)
      {
        limit_min = filter_limit_min_;
        limit_max = filter_limit_max_;
      }

      inline void
      setFilterLimitsNegative (const bool limit_negative)
      {
        filter_limit_negative_ = limit_negative;
      }

      inline void
      getFilterLimitsNegative (bool &limit_negative)
      {
        limit_negative = filter_limit_negative_;
      }

      inline bool
      getFilterLimitsNegative ()
      {
        return (filter_limit_negative_);
      }

    protected:
      Eigen::Vector4f leaf_size_;

      Eigen::Array4f inverse_leaf_size_;

      bool downsample_all_data_;

      bool save_leaf_layout_;

      std::vector<int> leaf_layout_;

      Eigen::Vector4i min_b_, max_b_, div_b_, divb_mul_;

      std::string filter_field_name_;

      double filter_limit_min_;

      double filter_limit_max_;

      bool filter_limit_negative_;

      typedef typename pcl::traits::fieldList<PointT>::type FieldList;

      void 
      applyFilter (PointCloud &output);
  };

  template <>
  class PCL_EXPORTS VoxelGrid<sensor_msgs::PointCloud2> : public Filter<sensor_msgs::PointCloud2>
  {
    using Filter<sensor_msgs::PointCloud2>::filter_name_;
    using Filter<sensor_msgs::PointCloud2>::getClassName;

    typedef sensor_msgs::PointCloud2 PointCloud2;
    typedef PointCloud2::Ptr PointCloud2Ptr;
    typedef PointCloud2::ConstPtr PointCloud2ConstPtr;

    public:
      VoxelGrid () : 
        leaf_size_ (Eigen::Vector4f::Zero ()),
        inverse_leaf_size_ (Eigen::Array4f::Zero ()),
        downsample_all_data_ (true), 
        save_leaf_layout_ (false),
        leaf_layout_ (),
        min_b_ (Eigen::Vector4i::Zero ()),
        max_b_ (Eigen::Vector4i::Zero ()),
        div_b_ (Eigen::Vector4i::Zero ()),
        divb_mul_ (Eigen::Vector4i::Zero ()),
        filter_field_name_ (""), 
        filter_limit_min_ (-FLT_MAX), 
        filter_limit_max_ (FLT_MAX),
        filter_limit_negative_ (false)
      {
        filter_name_ = "VoxelGrid";
      }

      virtual ~VoxelGrid ()
      {
      }

      inline void 
      setLeafSize (const Eigen::Vector4f &leaf_size) 
      { 
        leaf_size_ = leaf_size; 
        // Avoid division errors
        if (leaf_size_[3] == 0)
          leaf_size_[3] = 1;
        // Use multiplications instead of divisions
        inverse_leaf_size_ = Eigen::Array4f::Ones () / leaf_size_.array ();
      }

      inline void
      setLeafSize (float lx, float ly, float lz)
      {
        leaf_size_[0] = lx; leaf_size_[1] = ly; leaf_size_[2] = lz;
        // Avoid division errors
        if (leaf_size_[3] == 0)
          leaf_size_[3] = 1;
        // Use multiplications instead of divisions
        inverse_leaf_size_ = Eigen::Array4f::Ones () / leaf_size_.array ();
      }

      inline Eigen::Vector3f 
      getLeafSize () { return (leaf_size_.head<3> ()); }

      inline void 
      setDownsampleAllData (bool downsample) { downsample_all_data_ = downsample; }

      inline bool 
      getDownsampleAllData () { return (downsample_all_data_); }

      inline void 
      setSaveLeafLayout (bool save_leaf_layout) { save_leaf_layout_ = save_leaf_layout; }

      inline bool 
      getSaveLeafLayout () { return (save_leaf_layout_); }

      inline Eigen::Vector3i 
      getMinBoxCoordinates () { return (min_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getMaxBoxCoordinates () { return (max_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getNrDivisions () { return (div_b_.head<3> ()); }

      inline Eigen::Vector3i 
      getDivisionMultiplier () { return (divb_mul_.head<3> ()); }

      inline int 
      getCentroidIndex (float x, float y, float z)
      {
        return (leaf_layout_.at ((Eigen::Vector4i (static_cast<int> (floor (x * inverse_leaf_size_[0])), 
                                                   static_cast<int> (floor (y * inverse_leaf_size_[1])), 
                                                   static_cast<int> (floor (z * inverse_leaf_size_[2])), 
                                                   0) 
                - min_b_).dot (divb_mul_)));
      }

      inline std::vector<int> 
      getNeighborCentroidIndices (float x, float y, float z, const Eigen::MatrixXi &relative_coordinates)
      {
        Eigen::Vector4i ijk (static_cast<int> (floor (x * inverse_leaf_size_[0])), 
                             static_cast<int> (floor (y * inverse_leaf_size_[1])), 
                             static_cast<int> (floor (z * inverse_leaf_size_[2])), 0);
        Eigen::Array4i diff2min = min_b_ - ijk;
        Eigen::Array4i diff2max = max_b_ - ijk;
        std::vector<int> neighbors (relative_coordinates.cols());
        for (int ni = 0; ni < relative_coordinates.cols (); ni++)
        {
          Eigen::Vector4i displacement = (Eigen::Vector4i() << relative_coordinates.col(ni), 0).finished();
          // checking if the specified cell is in the grid
          if ((diff2min <= displacement.array()).all() && (diff2max >= displacement.array()).all())
            neighbors[ni] = leaf_layout_[((ijk + displacement - min_b_).dot (divb_mul_))]; // .at() can be omitted
          else
            neighbors[ni] = -1; // cell is out of bounds, consider it empty
        }
        return (neighbors);
      }
      
      inline std::vector<int> 
      getNeighborCentroidIndices (float x, float y, float z, const std::vector<Eigen::Vector3i> &relative_coordinates)
      {
        Eigen::Vector4i ijk (static_cast<int> (floorf (x * inverse_leaf_size_[0])), static_cast<int> (floorf (y * inverse_leaf_size_[1])), static_cast<int> (floorf (z * inverse_leaf_size_[2])), 0);
        std::vector<int> neighbors;
        neighbors.reserve (relative_coordinates.size ());
        for (std::vector<Eigen::Vector3i>::const_iterator it = relative_coordinates.begin (); it != relative_coordinates.end (); it++)
          neighbors.push_back (leaf_layout_[(ijk + (Eigen::Vector4i() << *it, 0).finished() - min_b_).dot (divb_mul_)]);
        return (neighbors);
      }

      inline std::vector<int> 
      getLeafLayout () { return (leaf_layout_); }

      inline Eigen::Vector3i 
      getGridCoordinates (float x, float y, float z) 
      { 
        return (Eigen::Vector3i (static_cast<int> (floor (x * inverse_leaf_size_[0])), 
                                 static_cast<int> (floor (y * inverse_leaf_size_[1])), 
                                 static_cast<int> (floor (z * inverse_leaf_size_[2])))); 
      }

      inline int 
      getCentroidIndexAt (const Eigen::Vector3i &ijk)
      {
        int idx = ((Eigen::Vector4i() << ijk, 0).finished() - min_b_).dot (divb_mul_);
        if (idx < 0 || idx >= static_cast<int> (leaf_layout_.size ())) // this checks also if leaf_layout_.size () == 0 i.e. everything was computed as needed
        {
          //if (verbose)
          //  PCL_ERROR ("[pcl::%s::getCentroidIndexAt] Specified coordinate is outside grid bounds, or leaf layout is not saved, make sure to call setSaveLeafLayout(true) and filter(output) first!\n", getClassName ().c_str ());
          return (-1);
        }
        return (leaf_layout_[idx]);
      }

      inline void
      setFilterFieldName (const std::string &field_name)
      {
        filter_field_name_ = field_name;
      }

      inline std::string const
      getFilterFieldName ()
      {
        return (filter_field_name_);
      }

      inline void
      setFilterLimits (const double &limit_min, const double &limit_max)
      {
        filter_limit_min_ = limit_min;
        filter_limit_max_ = limit_max;
      }

      inline void
      getFilterLimits (double &limit_min, double &limit_max)
      {
        limit_min = filter_limit_min_;
        limit_max = filter_limit_max_;
      }

      inline void
      setFilterLimitsNegative (const bool limit_negative)
      {
        filter_limit_negative_ = limit_negative;
      }

      inline void
      getFilterLimitsNegative (bool &limit_negative)
      {
        limit_negative = filter_limit_negative_;
      }

      inline bool
      getFilterLimitsNegative ()
      {
        return (filter_limit_negative_);
      }

    protected:
      Eigen::Vector4f leaf_size_;

      Eigen::Array4f inverse_leaf_size_;

      bool downsample_all_data_;

      bool save_leaf_layout_;

      std::vector<int> leaf_layout_;

      Eigen::Vector4i min_b_, max_b_, div_b_, divb_mul_;

      std::string filter_field_name_;

      double filter_limit_min_;

      double filter_limit_max_;

      bool filter_limit_negative_;

      void 
      applyFilter (PointCloud2 &output);
  };
}

#endif  //#ifndef PCL_FILTERS_VOXEL_GRID_MAP_H_