参考:https://www.cnblogs.com/oloroso/p/10758029.html
Streaming Maximum-Minimum Filter Using No More than Three Comparisons per Element
A fast algorithm for local minimum and maximum filters on rectangular and octagonal kernels
首先,对于一个多维的数据,都可以逐个维度进行处理。比如说一个图片,也就是二维数据,可以先对每一行进行处理,然后再对每一列进行处理,这样得到的结果与行列同时处理是一样的。
这个算法的过程大概是这样的:
1、首先遍历一行数据中最左边的r*2+1
个数据,获取最小值和最小值的位置。然后对左边边界部分的处理,直接赋最小值。
2、从r+1
位置开始向后遍历,一直到右边界部分。
3、遍历的时候,判断上一次获取的最小值索引minIndex
,是否在当前位置的领域r
以内。
如果不在,则遍历当前位置的领域r
范围,找出最小值的位置。也可以先与当前位置领域r
内最右边的比较,如果最右边的小于minIndex
位置的值,则minIndex
就是这最右边的这个,否则就需要遍历当前位置领域r
范围内。
如果在,则说明当前位置领域r
内,除了最右边的元素,肯定都小于minIndex
处的值。因为minIndex
是当前位置上一个的领域r
内的最小值,而上一个位置的领域r
范围与当前位置的领域r
范围只偏移了一个位置。
MinFilterOneRow 单行滤波代码
/**********************************************************************//**
* @brief 对一行数据进行滤波,每个值用邻域 r 内的最小值替换.
* @author solym@sohu.com/ymwh@foxmail.com
* @date 2019/4/23
* @param srcData 待滤波数据地址.
* @param srcChanelCount 待滤波数据每个像素的通道数.
* @param srcChanelIndex 待滤波数据要进行滤波的通道[0,srcChanelCount).
* @param dstData 滤波后输出数据地址.
* @param dstChanelCount 滤波后输出数据每个像素的通道数.
* @param dstChanelIndex 滤波后数据要输出的通道索引[0,srcChanelCount).
* @param colnumCount 该行数据要滤波的像素数.
* @param radius 滤波的半径大小.
*************************************************************************/
template<typename PixelDataType>
void MinFilterOneRow(
PixelDataType* srcData, const size_t srcChanelCount, const size_t srcChanelIndex,
PixelDataType* dstData, const size_t dstChanelCount, const size_t dstChanelIndex,
const size_t colnumCount, const size_t radius)
{
PixelDataType* pSrc = srcData;
PixelDataType* pDst = dstData;
size_t minIndex = 0; // 记录最小值的下标
size_t blockSize = radius * 2 + 1; // 块大小,以当前点为中心,左右各radius的宽度
PixelDataType minValue = pSrc[srcChanelIndex]; // 比较中获取最小值进行记录
// 对第一个块进行处理(i从1开始,比较(i-1,i)位置像素值)
// 找出最小值(第一个块内的最小值,就是r位置(块中心)处的输出值)
for(size_t iPixel=1; iPixel < blockSize; ++iPixel){
PixelDataType value = pSrc[iPixel*srcChanelCount + srcChanelIndex];
// 使用 >= 比 > 更快推进minIndex向前走
if(minValue >= value){
minValue = value;
minIndex = iPixel;
}
}
// 输出到第一个块中心(r)位置处的值。
// 它已经是第一个块内的最小值,也就是该块左边都只能是这个值
for(size_t i=0;i<=radius;++i){
pDst[i*dstChanelCount + dstChanelIndex] = minValue;
}
// 开始处理r+1位置之后的值
for (size_t iPixel = radius + 1; iPixel < colnumCount - radius; ++iPixel) {
/* i-r i i+r
* |____________|___________|_
* └min
* 当前最小的索引在当前位置为中心的块的内(一定位于当前块内或前一个)
* iPixel是当前块的中心,下面说的当前位置都指iPixel
*/
if(minIndex >= (iPixel - radius)) {
// 当前最小索引位置值与当前位置为中心的块的最后一个值比较
// 根据下面的代码可知,如果mIndex在块的内部,它所在位置的值一定是最小的
// 进入本次循环时,minIndex是上次比较的值,而上一个块与当前块等长,位置差一位
// 所以可以直接和当前块最后一个像素值进行比较了,当前块也就完全比较完了
size_t nextBlockFirstIndex = iPixel + radius;
if(pSrc[minIndex*srcChanelCount + srcChanelIndex] >
pSrc[nextBlockFirstIndex*srcChanelCount + srcChanelIndex]){
// 赋值当前最小值索引和值
minIndex = nextBlockFirstIndex;
minValue = pSrc[nextBlockFirstIndex*srcChanelCount + srcChanelIndex];
}
}else{
// 如果不在当前位置为中心的块内,则对当前块进行查找最小值
// 则将minIndex设置该块的最左边位置
minIndex = iPixel - radius;
// 获取当前位置为中心的块的最小值和索引
minValue = pSrc[minIndex*srcChanelCount + srcChanelIndex];
size_t blockEnd = minIndex + blockSize;
for (size_t iBPixel = minIndex; iBPixel < blockEnd; ++iBPixel) {
PixelDataType value = pSrc[iBPixel*srcChanelCount + srcChanelIndex];
if(minValue >= value){
minIndex = iBPixel;
minValue = value;
}
}
} // end if minIndex > ...
pDst[iPixel*dstChanelCount + dstChanelIndex] = minValue;
} // end for iPixel
// 最后一个块中心位置的右边,一定都是和它中心位置的值是一样的
for (size_t i = colnumCount-radius; i < colnumCount; ++i) {
pDst[i*dstChanelCount + dstChanelIndex] = minValue;
}
}
template<typename PixelDataType>
void MinFilterOneMatrix(
PixelDataType* srcData, const size_t srcBytePerRow,
const size_t srcChanelCount, const size_t srcChanelIndex,
PixelDataType* dstData, const size_t dstBytePerRow,
const size_t dstChanelCount, const size_t dstChanelIndex,
const size_t rowCount, const size_t colCount,
const size_t radius)
{
unsigned char* pSrc = reinterpret_cast<unsigned char*>(srcData);
unsigned char* pDst = reinterpret_cast<unsigned char*>(dstData);
// 保存中间结果
std::vector<PixelDataType> tmpData(rowCount * colCount);
// 逐行进行滤波
for (size_t row = 0; row < rowCount; ++row) {
// 获取输入和输出每行的行首位置
PixelDataType* pSrcRowFirst = (PixelDataType*)(pSrc + row * srcBytePerRow);
PixelDataType* pDstRowFirst = tmpData.data() + row * colCount;
// 对当前行进行滤波
MinFilterOneRow<PixelDataType>(pSrcRowFirst, srcChanelCount, srcChanelIndex,
pDstRowFirst, 1, 0,
colCount, radius);
}
// 将行滤波后的结果进行 行列转置(进行列滤波)
std::vector<PixelDataType> tmpDataT(rowCount * colCount);
for (size_t row = 0; row < rowCount; ++row) {
for (size_t col = 0; col < colCount; ++col) {
tmpDataT[col * rowCount + row] = tmpData[row * colCount + col];
}
}
// 对转置后的矩阵进行 逐行滤波(就是原行滤波后结果进行列滤波)
for (size_t col = 0; col < colCount; ++col) {
PixelDataType* pSrcColFirst = tmpDataT.data() + col * rowCount;
PixelDataType* pDstColFirst = tmpData.data() + col * rowCount;
// 对当前行进行滤波
MinFilterOneRow<PixelDataType>(pSrcColFirst, 1, 0,
pDstColFirst, 1, 0,
rowCount, radius);
}
// 将行列滤波后的结果输出
for (size_t row = 0; row < rowCount; ++row) {
PixelDataType* pDstRowFirst = (PixelDataType*)(pDst + row * dstBytePerRow);
for (size_t col = 0; col < colCount; ++col) {
pDstRowFirst[col * dstChanelCount + dstChanelIndex] = tmpData[col * rowCount + row];
}
}
}