7

我正在使用 CUDA 进行图像处理,但我对像素处理有疑问。

m x m应用卷积滤波器时,通常对图像的边界像素进行什么处理?

3 x 3卷积核中,忽略1图像的像素边界更容易处理,尤其是在使用共享内存改进代码时。实际上,在这种情况下,不需要检查给定像素是否具有所有可用的邻域(即坐标处的像素(0, 0)没有离开、左上、上邻居)。但是,删除1原始图像的像素边界可能会产生部分结果。

与此相反,我想处理图像中的所有像素,在使用共享内存改进时也是如此,例如,加载16 x 16像素,但计算内部14 x 14. 同样在这种情况下,忽略边界像素会生成更清晰的代码。

在这种情况下通常会做什么?

有没有人通常使用我的方法忽略边界像素?

当然,我知道答案取决于问题的类型,即按像素添加两个图像没有这个问题。

提前致谢。

4

3 回答 3

12

处理边框效果的常用方法是根据您的过滤器大小用额外的行和列填充原始图像。填充值的一些常见选择是:

  • 一个常数(例如零)
  • 根据需要多次复制第一行和最后一行/列
  • 在边界处反射图像(例如 column[-1] = column[1],column[-2] = column[2])
  • 包装图像值(例如 column[-1] = column[width-1], column[-2] = column[width-2])
于 2011-04-19T11:45:19.013 回答
5

tl; dr:这取决于您要解决的问题——没有适用于所有问题的解决方案。实际上,从数学上讲,我怀疑可能根本没有“解决方案”,因为我认为这是您被迫处理的不适定问题。

(为我鲁莽滥用数学提前道歉)

为了演示,让我们考虑一个假设所有像素分量和内核值都为正的情况。为了了解这些答案中的一些如何导致我们误入歧途,让我们进一步考虑一个简单的平均(“盒子”)过滤器。如果我们将图像边界外的值设置为零,那么这将明显降低边界 ceil(n/2)(曼哈顿距离)内每个像素的平均值。因此,您将在过滤后的图像上获得“暗”边框(假设单个强度分量或 RGB 颜色空间 - 您的结果将因颜色空间而异!)。请注意,如果我们将边界外的值设置为任意常数,则可以得出类似的论点——平均值将趋于该常数。如果您的典型图像的边缘无论如何都趋向于 0,则常数为零可能是合适的。像高斯这样的内核,但是问题不会那么明显,因为内核值会随着与中心的距离而迅速减小。

现在假设我们选择重复边缘值而不是使用常数。这与在图像周围制作边框并复制行、列或角足够多次以确保过滤器保持在新图像“内部”相同。您也可以将其视为对样本坐标进行钳位/饱和。这与我们的简单盒式过滤器存在问题,因为它过分强调了边缘像素的值。一组边缘像素会出现不止一次,但它们都获得相同的权重w=(1/(n*n))。假设我们对值为 K 的边缘像素进行 3 次采样。这意味着它对平均值的贡献是:

K*w + K*w + K*w  = K*3*w

如此有效,以至于一个像素的平均权重更高。请注意,由于这是一个平均过滤器,因此权重在内核上是一个常数。然而,这个论点也适用于权重因位置而异的内核(再次:想想高斯内核......)。

假设我们包装或反映采样坐标,以便我们仍然使用图像边界内的值。与使用常量相比,这具有一些有价值的优势,但也不一定是“正确的”。例如,您在上边框的对象与底部的对象相似的地方拍摄了多少张照片?除非你拍的是镜面光滑的湖泊,否则我怀疑这是真的。如果您正在拍摄岩石的照片以用作游戏中的纹理,则包裹或反射可能是合适的。我敢肯定,这里有很多重要的观点,关于包装和反射如何可能减少使用傅立叶变换产生的任何伪影。然而,这又回到了同样的想法:

那么,如果您要过滤蓝天下鲜红色岩石的照片,您能做些什么呢?显然,您不想在蓝天中添加橙色的薄雾,在红色的岩石上添加蓝色的绒毛。反射样本坐标是有效的,因为我们希望在反射坐标处找到与那些像素相似的颜色......除非,只是为了论证,我们想象过滤器内核太大以至于反射坐标会延伸到地平线之外。

让我们回到盒子过滤器的例子。使用此过滤器的另一种方法是停止考虑使用静态内核,并回想一下该内核的用途。平均/盒式滤波器旨在将像素分量相加,然后除以相加的像素数。这个想法是消除噪音。如果我们愿意在抑制边界附近的噪声方面降低效率,我们可以简单地将更少的像素相加并除以相应的更小的数字。这可以扩展到具有类似我将称之为“标准化”术语的过滤器 - 与过滤器的面积或体积相关的术语。对于“面积”术语,您计算边界内的内核权重的数量并忽略那些不在边界内的权重。然后将此计数用作“面积” (这可能涉及额外的乘法)。对于体积(再次:假设权重为正!)只需对内核权重求和。这个想法对于微分滤波器可能很糟糕,因为与噪声像素竞争的像素较少,而且差分对噪声非常敏感。此外,一些过滤器是通过数值优化和/或经验数据而非从头算/分析方法得出的,因此可能缺乏显而易见的“归一化”因子。

于 2014-10-07T21:24:25.577 回答
1

你的问题有点宽泛,我相信它混合了两个问题:

  1. 处理边界条件
  2. 处理晕区

遇到第一个问题(边界条件),例如,在计算图像和3 x 3内核之间的卷积时。当卷积窗口遇到边界时,就会出现将图像扩展到其边界之外的问题。

遇到第二个问题(晕区),例如,当16 x 16在共享内存中加载块时,必须处理内部块14 x 14以计算二阶导数。

对于第二个问题,我认为一个有用的问题是:Analyzing memory access coalescing of my CUDA kernel

关于将信号扩展到其边界之外,由于提供的不同寻址模式,纹理存储器在这种情况下提供了一个有用的工具,请参阅CUDA 纹理的不同寻址模式

下面,我将提供一个示例,说明如何使用纹理内存在周期性边界条件下实现中值滤波器。

#include <stdio.h>

#include "TimingGPU.cuh"
#include "Utilities.cuh"

texture<float, 1, cudaReadModeElementType> signal_texture;

#define BLOCKSIZE 32

/*************************************************/
/* KERNEL FUNCTION FOR MEDIAN FILTER CALCULATION */
/*************************************************/
__global__ void median_filter_periodic_boundary(float * __restrict__ d_vec, const unsigned int N){

    unsigned int tid = threadIdx.x + blockIdx.x * blockDim.x;

    if (tid < N) {

        float signal_center = tex1D(signal_texture, tid - 0);
        float signal_before = tex1D(signal_texture, tid - 1);
        float signal_after  = tex1D(signal_texture, tid + 1);

        printf("%i %f %f %f\n", tid, signal_before, signal_center, signal_after);

        d_vec[tid] = (signal_center + signal_before + signal_after) / 3.f;

    }
}


/********/
/* MAIN */
/********/
int main() {

    const int N = 10;

    // --- Input host array declaration and initialization
    float *h_arr = (float *)malloc(N * sizeof(float));
    for (int i = 0; i < N; i++) h_arr[i] = (float)i;

    // --- Output host and device array vectors
    float *h_vec = (float *)malloc(N * sizeof(float));
    float *d_vec;   gpuErrchk(cudaMalloc(&d_vec, N * sizeof(float)));

    // --- CUDA array declaration and texture memory binding; CUDA array initialization
    cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<float>();
    //Alternatively
    //cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc(32, 0, 0, 0, cudaChannelFormatKindFloat);

    cudaArray *d_arr;   gpuErrchk(cudaMallocArray(&d_arr, &channelDesc, N, 1));
    gpuErrchk(cudaMemcpyToArray(d_arr, 0, 0, h_arr, N * sizeof(float), cudaMemcpyHostToDevice));

    cudaBindTextureToArray(signal_texture, d_arr); 
    signal_texture.normalized = false; 
    signal_texture.addressMode[0] = cudaAddressModeWrap;

    // --- Kernel execution
    median_filter_periodic_boundary<<<iDivUp(N, BLOCKSIZE), BLOCKSIZE>>>(d_vec, N);
    gpuErrchk(cudaPeekAtLastError());
    gpuErrchk(cudaDeviceSynchronize());

    gpuErrchk(cudaMemcpy(h_vec, d_vec, N * sizeof(float), cudaMemcpyDeviceToHost));

    for (int i=0; i<N; i++) printf("h_vec[%i] = %f\n", i, h_vec[i]);

    printf("Test finished\n");

    return 0;
}
于 2015-09-02T08:46:01.067 回答