空间域中的图像卷积

Question

我正在尝试使用线性空间域卷积复制此链接的输出结果。

首先将图像转换为2D double数组，然后进行卷积。 图像和卷积核的大小相同。 在卷积之前对图像进行填充，然后在卷积后相应地裁剪。

与基于FFT的卷积相比，输出结果很奇怪且不正确。

我该如何解决这个问题？

请注意，我从Matlab获得了以下与C# FFT输出匹配的图像输出：

.

^{更新1：根据@Ben Voigt的评论，我将Rescale()函数替换为1而不是255.0，因此输出效果有了很大的改善。但是，输出结果仍与FFT输出（正确的结果）不匹配。}

.

^{更新2：根据@Cris Luengo的评论，我通过拼接填充了图像，然后执行了空间卷积。 结果如下：}

^{因此，输出结果比之前更糟糕。 但是，这与链接答案的第二个输出类似，这意味着循环卷积不是解决方案。}

.

^{更新3：我使用@Cris Luengo答案中提出的Sum()函数。 结果是**更新1**的更改版本：}

^{但是，它仍然与FFT版本不完全相同。}

.

^{更新4：根据@Cris Luengo的评论，我对两个输出结果进行了减法以查看差异：}
,

<sup>1. 空间域减频率域
2. 频率域减空间域</sup>

<

    public static double[,] LinearConvolutionSpatial(double[,] image, double[,] mask)
    {
        int maskWidth = mask.GetLength(0);
        int maskHeight = mask.GetLength(1);

        double[,] paddedImage = ImagePadder.Pad(image, maskWidth);

        double[,] conv = Convolution.ConvolutionSpatial(paddedImage, mask);

        int cropSize = (maskWidth/2);

        double[,] cropped = ImageCropper.Crop(conv, cropSize);

        return conv;
    } 
    static double[,] ConvolutionSpatial(double[,] paddedImage1, double[,] mask1)
    {
        int imageWidth = paddedImage1.GetLength(0);
        int imageHeight = paddedImage1.GetLength(1);

        int maskWidth = mask1.GetLength(0);
        int maskHeight = mask1.GetLength(1);

        int convWidth = imageWidth - ((maskWidth / 2) * 2);
        int convHeight = imageHeight - ((maskHeight / 2) * 2);

        double[,] convolve = new double[convWidth, convHeight];

        for (int y = 0; y < convHeight; y++)
        {
            for (int x = 0; x < convWidth; x++)
            {
                int startX = x;
                int startY = y;

                convolve[x, y] = Sum(paddedImage1, mask1, startX, startY);
            }
        }

        Rescale(convolve);

        return convolve;
    } 

    static double Sum(double[,] paddedImage1, double[,] mask1, int startX, int startY)
    {
        double sum = 0;

        int maskWidth = mask1.GetLength(0);
        int maskHeight = mask1.GetLength(1);

        for (int y = startY; y < (startY + maskHeight); y++)
        {
            for (int x = startX; x < (startX + maskWidth); x++)
            {
                double img = paddedImage1[x, y];
                double msk = mask1[x - startX, y - startY];
                sum = sum + (img * msk);
            }
        }

        return sum;
    }

    static void Rescale(double[,] convolve)
    {
        int imageWidth = convolve.GetLength(0);
        int imageHeight = convolve.GetLength(1);

        double maxAmp = 0.0;

        for (int j = 0; j < imageHeight; j++)
        {
            for (int i = 0; i < imageWidth; i++)
            {
                maxAmp = Math.Max(maxAmp, convolve[i, j]);
            }
        }

        double scale = 1.0 / maxAmp;

        for (int j = 0; j < imageHeight; j++)
        {
            for (int i = 0; i < imageWidth; i++)
            {
                double d = convolve[i, j] * scale;
                convolve[i, j] = d;
            }
        }
    } 

    public static Bitmap ConvolveInFrequencyDomain(Bitmap image1, Bitmap kernel1)
    {
        Bitmap outcome = null;

        Bitmap image = (Bitmap)image1.Clone();
        Bitmap kernel = (Bitmap)kernel1.Clone();

        //linear convolution: sum. 
        //circular convolution: max
        uint paddedWidth = Tools.ToNextPow2((uint)(image.Width + kernel.Width));
        uint paddedHeight = Tools.ToNextPow2((uint)(image.Height + kernel.Height));

        Bitmap paddedImage = ImagePadder.Pad(image, (int)paddedWidth, (int)paddedHeight);
        Bitmap paddedKernel = ImagePadder.Pad(kernel, (int)paddedWidth, (int)paddedHeight);

        Complex[,] cpxImage = ImageDataConverter.ToComplex(paddedImage);
        Complex[,] cpxKernel = ImageDataConverter.ToComplex(paddedKernel);

        // call the complex function
        Complex[,] convolve = Convolve(cpxImage, cpxKernel);

        outcome = ImageDataConverter.ToBitmap(convolve);

        outcome = ImageCropper.Crop(outcome, (kernel.Width/2)+1);

        return outcome;
    } 

- user366312

---

对于图像处理，“时域”意味着电影。如果您的意思是“静止图像上不是频率的领域”，那么应该使用“空间域”。 - Ben Voigt

我对C#一无所知，所以也许这是一个愚蠢的评论：双重循环中的string str = string.Empty;是做什么用的？你没有使用str，这可能会拖慢速度？-- 无论如何，在空间域中使用大内核进行卷积本质上是非常昂贵的。这就是为什么FFT方法如此重要的原因。 - Cris Luengo

@CrisLuengo，string.Empty 仅用于偶尔的调试目的。已删除。 - user366312

@CrisLuengo，我知道空间域卷积很耗费资源。这只是为了发展我的技能。我会移除问题中的性能部分。 - user366312

@CrisLuengo，convWidth是卷积/输出图像的大小。这里没有显示填充操作。填充在ConvolveTd()函数之前和之外完成。 - user366312

Answer 1

您当前的输出更像是自相关函数，而不是Lena与自身的卷积。我认为问题可能出在您的Sum函数上。

如果您查看卷积和的定义，您会发现核（或图像，无所谓）被镜像：

sum_m( f[n-m] g[m] )

对于一个函数，m 带有加号，而对于另一个函数，则带有减号。

您需要修改您的 Sum 函数以正确读取 mask1 图像的顺序：

static double Sum(double[,] paddedImage1, double[,] mask1, int startX, int startY)
{
    double sum = 0;

    int maskWidth = mask1.GetLength(0);
    int maskHeight = mask1.GetLength(1);

    for (int y = startY; y < (startY + maskHeight); y++)
    {
        for (int x = startX; x < (startX + maskWidth); x++)
        {
            double img = paddedImage1[x, y];
            double msk = mask1[maskWidth - x + startX - 1, maskHeight - y + startY - 1];
            sum = sum + (img * msk);
        }
    }

    return sum;
}

另一种选择是将 mask1 的镜像版本传递给此函数。

- Cris Luengo

---

谢谢您的回答。我认为修改Sum()不是一个好主意。您能告诉我这里的“镜像”是什么意思吗？是指水平翻转吗？ - user366312

我猜索引越界是由于一个偏移量错误导致的。这就是我发表了我甚至无法编译的代码所得到的结果。 :) - Cris Luengo

@anonymous：不，它是垂直和水平翻转。您需要翻转每个图像维度。 - Cris Luengo

好的，我尝试了你的Sum()函数。输出结果更加优化了。你还有什么窍门可以使它100%兼容吗？或者，我需要对FFT输出做些什么吗？ - user366312

@anonymous：周末我在SO上没有太多时间。我看到了你的更新，我没有发现任何问题。我现在唯一能想到的是让两者的填充相等——或者至少确保最后的裁剪产生相同大小的结果，这样输出就可以更好地进行比较。减法的结果很奇怪，但我认为这是因为图像大小不同造成的。 - Cris Luengo

Answer 2

我从这个链接找到了解决方案。主要的线索是引入了一个offset和一个factor。

• factor是内核中所有值的总和。
• offset是用于进一步修正输出的任意值。

@Cris Luengo的回答也提出了一个有效的观点。

下面的源代码在给定的链接中提供：

    private void SafeImageConvolution(Bitmap image, ConvMatrix fmat) 
    { 
        //Avoid division by 0 
        if (fmat.Factor == 0) 
            return; 

        Bitmap srcImage = (Bitmap)image.Clone(); 

        int x, y, filterx, filtery; 
        int s = fmat.Size / 2; 
        int r, g, b; 
        Color tempPix; 

        for (y = s; y < srcImage.Height - s; y++) 
        { 
            for (x = s; x < srcImage.Width - s; x++) 
            { 
                r = g = b = 0; 

                // Convolution 
                for (filtery = 0; filtery < fmat.Size; filtery++) 
                { 
                    for (filterx = 0; filterx < fmat.Size; filterx++) 
                    { 
                        tempPix = srcImage.GetPixel(x + filterx - s, y + filtery - s); 

                        r += fmat.Matrix[filtery, filterx] * tempPix.R; 
                        g += fmat.Matrix[filtery, filterx] * tempPix.G; 
                        b += fmat.Matrix[filtery, filterx] * tempPix.B; 
                    } 
                } 

                r = Math.Min(Math.Max((r / fmat.Factor) + fmat.Offset, 0), 255); 
                g = Math.Min(Math.Max((g / fmat.Factor) + fmat.Offset, 0), 255); 
                b = Math.Min(Math.Max((b / fmat.Factor) + fmat.Offset, 0), 255); 

                image.SetPixel(x, y, Color.FromArgb(r, g, b)); 
            } 
        } 
    } 

- user366312