c++ - C 和 C++ 中的快速 7x7 2D 中值滤波器

Question

我正在尝试将以下代码从 MATLAB 转换为 C++:

function data = process(data)
    data = medfilt2(data, [7 7], 'symmetric');
    mask = fspecial('gaussian', [35 35], 12);
    data = imfilter(data, mask, 'replicate', 'same');
    maximum = max(data(:));
    data = 1 ./ (data/maximum);
    data(data > 10) = 16;
end

我的问题出在 medfilt2 上，它是一个二维中值滤波器。我需要它支持每像素 10 位和更 multimap 像。

1. 我调查了OpenCV它有一个支持 16 位的 5x5 中值滤波器，但 7x7 只支持字节。

   medianBlur

2. 我还研究了英特尔 IPP ，但我只能看到一维中值滤波器。 https://software.intel.com/en-us/node/502283

有二维滤镜的快速实现吗？

我正在寻找类似的东西:

1. Fast Median Search: An ANSI C Implementation 使用并行编程和向量化 ( AVX/SSE ) 操作...
2. 二维数字信号处理 II。转换和中值滤波器。 由 T.S.Huang 编辑。施普林格出版社。 1981 年。

中有更多代码示例Fast median filtering with implementations in C/C++/C#/VB.NET/Delphi .

我还找到了 Median Filtering in Constant Time .

Answer 1

由于 OpenCV 没有为大内核大小(大于 5)实现 16 位中值滤波器这一事实，我尝试了三种不同的策略。

都是基于Huang的[2]滑动窗口算法。也就是说，当窗口从左向右滑动时，通过删除和插入像素条目来更新直方图。这对于 8 位图像来说非常简单，并且已经在 OpenCV 中实现了。然而，一个大的 65536 bin 直方图使得计算有点困难。

...The algorithm still remains O(log r), but storage considerations render it impractical for 16-bit images and impossible for floating-point images. [3]

我在适用的地方使用了算法 C++ 标准库，并且没有实现 Weiss 的额外优化策略。

1) 一个朴素的排序实现。我认为这是任意像素类型(尤其是 float )的最佳起点。

// copy pixels in the sliding window to a temporary vec and
// compute the median value (size is always odd)
memcpy( &v[0], &window[0], window.size() * sizeof(_Type) );
std::vector< _Type >::iterator it = v.begin() + v.size()/2;
std::nth_element( v.begin(), it, v.end() );
return *it;

2) 稀疏直方图。我们不想跨越 65536 个 bin 来找到每个像素的中值，那么如何存储稀疏直方图呢？同样，这适用于所有像素类型，但如果窗口中的所有像素都不同(例如 float )，则没有意义。

typedef std::map< _Type, int > Map;
//...
// inside the sliding window, update the histogram as follows
for ( /* pixels to remove */ )
{
    // _Type px
    Map::iterator it = map.find( px );
    if ( it->second > 1 )
        it->second -= 1;
    else
        map.erase( it );
}
// ...
for ( /* pixels to add */ )
{
    // _Type px
    Map::iterator lower = map.lower_bound( px );
    if ( lower != map.end() && lower->first == px )
        lower->second += 1;
    else
        map.insert( lower, std::pair<_Type,int>( px, 1 ) );
}
//... and compute the median by integrating from the one end until
// until the appropriate sum is reached ..

3) 密集的直方图。所以这是密集的直方图，但不是简单的 65536 数组，我们通过将其分成子箱来使搜索更容易一些，例如:

[0...65535] <- px
[0...4095] <- px / 16
[0...255] <- px / 256
[0...15] <- px / 4096

这会使插入速度变慢(按常数时间)，但搜索速度会快很多。我发现 16 是个不错的数字。

图 我测试了方法 (1) 红色、(2) 蓝色和 (3) 黑色以及 8bpp OpenCV(绿色)。对于除 OpenCV 之外的所有图像，输入图像都是 16-bpp 灰度。虚线在动态范围 [0,255] 处被截断，平滑线在 [0, 8020] 处被截断(通过乘以 16 和平滑以增加像素值的更多方差)。

有趣的是稀疏直方图随着像素值方差的增加而发散。第 N 个元素始终是安全的赌注，OpenCV 是最快的(如果 8bpp 可以)并且密集直方图落后。

我使用 Windows 7、8 x 3.4 GHz 和 Visual Studio v.10。我的运行多线程，OpenCV 实现是单线程的。输入图像大小 2136x3201(http://i.imgur.com/gg9Z2aB.jpg，来自 Vogue)。

[2]: Huang, T: 《二维信号处理II:变换 和中值滤波器”，1981

[3]:Weiss，B:“快速中值和双边过滤”，2006 年