PCL 点云库学习四（双边滤波器）

双边滤波器

论文：[Bilateral Filtering for Gray and Color Images]
临近采样点加权平均，同时剔除差异大的点，数据更加平滑去除噪点

1. 简介

• 双边滤波是一种简单的、非迭代的边缘保持平滑方法。
• 可以使用灰度值或者几何距离，采用值范围内的临近值
• 可以使用颜色的3个波段（RGB）或者CIE-LAB空间下使用
• 范围滤波：使用不同权重的的衰减来平均图像值，权重依赖于图像颜色或者光强度
• 阈值滤波+范围滤波

2. 思想

• 低通滤波器如下
  $c(\xi ,x):$x与临近点$\xi$之间的几何距离(贴近度)
  $$h(x)=k_d^{-1}(x){\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }f(\xi )c(\xi ,x)d\xi$$
  $k_d:$归一化
  $$k_d(x)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }c(\xi ,x)d\xi$$
  如果滤波器不随位移变化，那么向量差$\xi -x$和$k_d$是常量
• 范围滤波器
  $s(f(\xi ),f(x)):$临近点中心x与临近点$\xi$之间几何相似度
  $$h(x)=h_r^{-1}(x){\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }f(\xi )s(f(\xi ),f(x))d\xi$$
  $k_r(x):$归一化
  $$k_r(x)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }s(f(\xi ),f(x))d\xi$$
• 组合滤波
  $$h(x)=k^{(-1)(x)}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }f(\xi )c(\xi ,x)s(f(\xi ),f(x))d\xi$$
  归一化
  $$k(x)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }c(\xi ,x)s(f(\xi ),f(x))d\xi$$
• 它将x处的像素值替换为相似和邻近像素值的平均值
• 在平滑区域，小邻域中的像素值彼此相似，归一化相似函数k，$k^{(-1)(x)}$接近于1
• 双边滤波器基本上起到标准域滤波器的作用，并且平均掉由噪声引起的像素值之间的小的、弱相关的差异
• 作用：见上图a，假设对于明暗区域存在明显的边界，如图a所示，b为$k(x)$计算到的权重,最后得到c。

3. 典型双边滤波器 高斯滤波器

• $d(\xi ,x)=d(\xi -x)=||\xi -x||$
  $$c(\xi ,x)=exp(-\frac{1}{2}(\frac{d(\xi ,x)}{\sigma d}{)}^2)$$
• $\delta (\varphi ,f)=\delta (\varphi -f)=||\varphi -f||$:两个强度值之间的距离
  $$s(\xi ,x)=exp(-\frac{1}{2}(\frac{\delta (f(\xi ),f(x))}{{\sigma }_r}{)}^2)$$
• ${\sigma }_d:$该值越大，表示采用越多远方的信息，图像越模糊，如果图像放大或者缩小该值需要调整才能获得相应效果
• ${\sigma }_r:$表示实现实现组合所需的点的数量，距离小于该值的临近点将会采用，大于在计算时不采用，放大缩小时需要调整
• 该方法为无偏估计

4. pcl中使用(需要强度值即pcl::PointXYZI)

#include <pcl/kdtree/flann.h>
#include <pcl/io/pcd_io.h>
#include <pcl/filters/bilateral.h>

#include <pcl/kdtree/kdtree.h>
#include <pcl/filters/fast_bilateral.h>
#include <pcl/search/flann_search.h>
#include <pcl/point_cloud.h>
#include <pcl/point_types.h>
#include <pcl/filters/bilateral.h>  
using namespace std;
using namespace pcl;
int main(int argc, const char** argv) {
    pcl::PointCloud<pcl::PointXYZI>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZI>);
    pcl::io::loadPCDFile("/home/n1/notes/pcl/bilateral_filter/test2.pcd",*cloud);
    pcl::PointCloud<pcl::PointXYZI>::Ptr outcloud(new pcl::PointCloud<pcl::PointXYZI>);
    pcl::search::KdTree<pcl::PointXYZI>::Ptr tree1(new  pcl::search::KdTree<pcl::PointXYZI>);
    tree1->setInputCloud(cloud);
    pcl::BilateralFilter<pcl::PointXYZI> nf;
    //双边滤波器
    nf.setInputCloud(cloud);  //设置输入点云
    nf.setSearchMethod(tree1); //设置搜索方法
    nf.setHalfSize(1);      //高斯滤波器的一半窗口值
    nf.setStdDev(0.01);           //标准差
    nf.filter(*outcloud);

    pcl::io::savePCDFile("/home/n1/notes/pcl/bilateral_filter/output.pcd",*outcloud);
    return 0;
}

快速双边滤波器

• 提出了一种新的双边滤波器的信号处理分析方法,对空间关系和强度进行下采样进行加速
• 与典型双边滤波器相比，运行相同时间该方法更加准确
• 在信号处理框架解释了双边滤波器：证明在高维空间下，双边滤波器可以表示为简单的非线性卷积
• 证明在高维空间下，下采样对于卷积结果影响不大：这种近似技术可以在控制误差的同时使速度提高几个数量级。
• 使用FSNR（signal-to-noise tatio） 方法快速递进求解

1. 使用
   注：需要pcl数据为结构点云，点云类似于图像（点云宽高与位置相关），即point_cloud_ptr->height>1.
   不支持PointXYZI点云格式

#include <pcl/kdtree/flann.h>
#include <pcl/io/pcd_io.h>
#include <pcl/filters/bilateral.h>

#include <pcl/kdtree/kdtree.h>
#include <pcl/filters/fast_bilateral.h>
#include <pcl/filters/fast_bilateral_omp.h> 

#include <pcl/search/flann_search.h>
#include <pcl/point_cloud.h>
#include <pcl/point_types.h>
#include <pcl/filters/bilateral.h>  
using namespace std;
using namespace pcl;

void create_ellipse_pointcloud( pcl::PointCloud<pcl::PointXYZRGB>::Ptr& point_cloud_ptr ){
     // 以椭圆为边线沿z轴拉伸获取其点云，并赋予红绿蓝渐变色。
  uint8_t r(255), g(15), b(15);
  for (float z(-1.0); z <= 1.0; z += 0.05)
  {
    for (float angle(0.0); angle <= 360.0; angle += 5.0)
    {
      pcl::PointXYZ basic_point;
      basic_point.x = 0.5 * cosf (float(angle/180*M_PI));
      basic_point.y = sinf (float(angle/180*M_PI));
      basic_point.z = z;

      pcl::PointXYZRGB point;
      point.x = basic_point.x;
      point.y = basic_point.y;
      point.z = basic_point.z;
      uint32_t rgb = (static_cast<uint32_t>(r) << 16 |
              static_cast<uint32_t>(g) << 8 | static_cast<uint32_t>(b));
      point.rgb = *reinterpret_cast<float*>(&rgb);
      point_cloud_ptr->points.push_back (point);
    }
    if (z < 0.0)
    {
      r -= 12;
      g += 12;
    }
    else
    {
      g -= 12;
      b += 12;
    }
  }
  point_cloud_ptr->width = (int) point_cloud_ptr->points.size ()/20;
  point_cloud_ptr->height = 20;
}
int main(int argc, const char** argv) {
    pcl::PointCloud<pcl::PointXYZRGB>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZRGB>);
    create_ellipse_pointcloud(cloud);
    pcl::PointCloud<pcl::PointXYZRGB>::Ptr outcloud(new pcl::PointCloud<pcl::PointXYZRGB>);
    pcl::FastBilateralFilter<pcl::PointXYZRGB> fbf;
    //双边滤波器
    fbf.setInputCloud(cloud);  //设置输入点云
    //fbf.setSearchMethod(tree1);
    fbf.setSigmaS(1);      //高斯滤波器的一半窗口值
    fbf.setSigmaR(0.01);           //标准差
    fbf.applyFilter(*outcloud);

    pcl::io::savePCDFile("/home/n1/notes/pcl/fast_bilateral/output1.pcd",*outcloud);
    return 0;
}