Delphi图像处理 -- 图像高斯模糊处理（改进版）

我在文章《Delphi图像处理 -- 图像卷积及高斯模糊》中，介绍了利用通用的图像卷积过程对图像进行高斯模糊处理，其处理效果还不错，处理小型图像时感觉也还行，但是处理较大图像时的速度还是嫌慢，在我的P4 2.8G、1G内存的机器上对千万像素图像进行Q=3，R=5的高斯模糊处理，不包括图像装载和前期数据转换，耗时达8600ms以上，虽经几次修改，其处理速度始终得不到明显提高，主要原因还是采用通用卷积过程处理的问题：用R=5得到的卷积模板为11*11像素，一个像素有4个分量（32位ARGB），对每个象素必须作11*11*4=484个乘法、484个加法及4个除法，最后还得作4个分量是否超界的判断处理，想不慢都难啦！如果不是采用BASM定点数处理代码，其处理速度更是难以想象。

我在网上多次查找图像高斯模糊的优化算法，不少算法和处理方式，包括代码优化还不如我的那个高斯模糊处理过程，使我很失望。前天查找其它资料时，在国外某个网站上发现介绍图像高斯模糊处理方法时似乎与常规的算法有所不同，但都没有详细的资料（因为不懂外语，很少上国外网站，但看些公式、伪代码还是行的），经过反复琢磨，可以将其处理流程归纳如下：

1、用给定的确定Q和长度size，计算一个size+1长的高斯分布权数数据weights：

radius = size / 2 // 计算初始数据 for (i = -radius; i <= radius; i ++) { x = i / Q; weights[i+radius] = exp(-x * x / 2) } // 求和 sum = 0 for (i = -radius; i <= radius; i ++) { sum += weights[i+radius] } // 数据归一，即归一后的数据之和等于1 for (i = -radius; i <= radius; i ++) { weights[i+radius] /= sum }

2、使用weights对原图像作垂直的模糊运算，即以像素(x, y)为中心，对(x, y - radius)和(x, y + radius)的像素点与weights对应的值相乘后求和得到新的像素，并写入到一个临时的图像上相应的点上（因为数据进行了归一处理，乘积和不必再作除法运算）；

3、使用weights对临时图像作水平的模糊运算，即以像素(x, y)为中心，对(x- radius, y)和(x + radius, y)的像素点与weights对应的相乘后求和得到新的像素，并写入到目标图像上相应的点上。

处理过程结束。

由于上面的处理流程只是对图像每个象素作了一个“十”字型的运算，使得对每个象素点的运算大大减少，模糊长度越大，减少的越多。如前面所说的Q=3、R=5的模糊运算只需要11*2*4=88个乘法、88个加法即可。

我还是采用BASM按上面的流程作定点数运算，改进后的高斯模糊过程代码如下：

过程定义： // 高斯模糊图像。参数： // Dest输出图，Source原图，Data自身操作图像 // Q为标准差，Radius模糊半径，如果为0，则自动按Q计算 // Callback回调函数，返回True终止操作，CallbackData回调函数参数地址 procedure ImageGaussiabBlur(var Dest: TImageData; const Source: TImageData; Q: double; Radius: Integer = 0); overload; {$IF RTLVersion >= 17.00}inline;{$IFEND} procedure ImageGaussiabBlur(var Data: TImageData; Q: double; Radius: Integer = 0); overload; {$IF RTLVersion >= 17.00}inline;{$IFEND} // 无效参数或者被回调函数终止操作返回False。 function ImageGaussiabBlur(var Dest: TImageData; const Source: TImageData; Q: double; Radius: Integer; Callback: TImageAbort; CallbackData: Pointer): Boolean; overload; 代码实现： procedure CrossBlur(Weights: Pointer; Size, dstWidth, srcStride: Integer; var x); pascal; var height: Integer; dstOffset, srcOffset: Integer; reds, greens, blues: Integer; asm mov height, edx mov dstOffset, ebx push esi push edi push edx push ecx push eax // blur col // x是在分块模糊时，判断是否为列的首块。 // 如果是，按列模糊时，宽度加size-1，否则扫描线右移size-1 mov edi, x mov edx, [edi] add edx, ecx // edx = *x + width cmp edx, dstWidth jb @@1 xor edx, edx // if (edx >= dst.width) edx = 0 @@1: mov ebx, size dec ebx cmp dword ptr [edi], 0 jne @@2 // if (*px == 0) add ecx, ebx // { shl ebx, 2 // width += size - 1 sub eax, ebx // srcOffset -= ((size - 1) * 4) jmp @@3 // } @@2: shl ebx, 2 // else esi += (size - 1) * 4 add esi, ebx @@3: mov [edi], edx // *x = edx mov srcOffset, eax mov edi, weights // edi = weights @@cyLoop: push ecx @@cxLoop: push ecx push esi push edi xor ebx, ebx mov reds, ebx mov greens, ebx mov blues, ebx mov ecx, Size @@cblurLoop: movzx eax, [esi].TARGBQuad.Blue movzx edx, [esi].TARGBQuad.Green imul eax, [edi] imul edx, [edi] add blues, eax add greens, edx movzx eax, [esi].TARGBQuad.Red movzx edx, [esi].TARGBQuad.Alpha imul eax, [edi] imul edx, [edi] add reds, eax add ebx, edx add edi, 4 add esi, srcStride loop @@cblurLoop pop edi pop esi mov eax, blues mov edx, greens mov ecx, reds shr eax, 16 shr edx, 16 shr ecx, 16 shr ebx, 16 mov [esi].TARGBQuad.Blue, al mov [esi].TARGBQuad.Green, dl mov [esi].TARGBQuad.Red, cl mov [esi].TARGBQuad.Alpha, bl add esi, 4 pop ecx loop @@cxLoop add esi, srcOffset pop ecx dec height jnz @@cyLoop pop srcOffset pop ecx pop height pop edi pop esi // blur row @@ryLoop: push ecx @@rxLoop: push ecx push esi push edi xor ebx, ebx mov reds, ebx mov greens, ebx mov blues, ebx mov ecx, Size mov edi, Weights @@rblurLoop: movzx eax, [esi].TARGBQuad.Blue movzx edx, [esi].TARGBQuad.Green imul eax, [edi] imul edx, [edi] add blues, eax add greens, edx movzx eax, [esi].TARGBQuad.Red movzx edx, [esi].TARGBQuad.Alpha imul eax, [edi] imul edx, [edi] add reds, eax add ebx, edx add edi, 4 add esi, 4 loop @@rblurLoop pop edi pop esi mov eax, blues mov edx, greens mov ecx, reds shr eax, 16 shr edx, 16 shr ecx, 16 shr ebx, 16 mov [edi].TARGBQuad.Blue, al mov [edi].TARGBQuad.Green, dl mov [edi].TARGBQuad.Red, cl mov [edi].TARGBQuad.Alpha, bl add esi, 4 add edi, 4 pop ecx loop @@rxLoop add esi, srcOffset add edi, dstOffset pop ecx dec height jnz @@ryLoop end; function ImageGaussiabBlur(var Dest: TImageData; const Source: TImageData; Q: double; Radius: Integer; Callback: TImageAbort; CallbackData: Pointer): Boolean; var dst, src: TImageData; fweights: array of Single; weights: array of Integer; i, x, size: Integer; fx: Double; begin Result := False; if ImageEmpty(Dest) or ImageEmpty(Source) or (Q = 0.0) then Exit; if Radius <= 0 then begin if Abs(Q) < 1.0 then Radius := 1 else Radius := Round(Abs(Q)) + 2; end; size := Radius shl 1 + 1; SetLength(fweights, size); for i := 1 to Radius do begin fx := i / Q; fweights[Radius + i] := exp(-fx * fx / 2); fweights[Radius - i] := fweights[Radius + i]; end; fweights[Radius] := 1.0; fx := 0.0; for i := 0 to size - 1 do fx := fx + fweights[i]; SetLength(weights, size); for i := 0 to size - 1 do weights[i] := Round(fweights[i] / fx * 65536.0); SetLength(fweights, 0); GetExpandData(Dest, Source, Radius, dst, src); try x := 0; if Assigned(Callback) then Result := ExecuteAbort(dst, src, @CrossBlur, [weights, size, dst.Width, src.Stride, @x], Callback, CallbackData) else Result := ExecuteProc(dst, src, @CrossBlur, [weights, size, dst.Width, src.Stride, @x]); finally FreeImageData(src); end; end; procedure ImageGaussiabBlur(var Dest: TImageData; const Source: TImageData; Q: double; Radius: Integer); begin ImageGaussiabBlur(Dest, Source, Q, Radius, nil, nil); end; procedure ImageGaussiabBlur(var Data: TImageData; Q: double; Radius: Integer); begin ImageGaussiabBlur(Data, Data, Q, Radius, nil, nil); end;

用改进后的高斯模糊处理过程在我的机器上对千万像素图像进行Q=3，R=5的高斯模糊处理，不包括图像装载和前期数据转换，耗时为1390ms，处理速度确实得到了大幅度的提高。我是按32位ARGB颜色处理图像像素的，如果改为24位RGB颜色处理图像像素，耗时还可以减少，不过，RGB颜色没法处理PNG等32位像素格式的图像。

不用模板卷积方式，而采用“十”字运算进行高斯模糊处理，效果如何呢？请看下面的简单例子代码及处理效果图：

例子代码：

var data: TImageData; begin data := GetImageData(TGpBitmap.Create('../../Media/56-3.jpg', False), True); ImageGaussiabBlur(data, 3, 6); DrawImage(Canvas, 0, 0, data, 1); FreeImageData(data); end;

处理原图：

处理效果与Photoshop高斯模糊处理对比图：

左上是Photoshop半径3.0高斯模糊效果图，右上是本文过程Q=3.0,R=6高斯模糊效果图。

左下是Photoshop半径5.0高斯模糊效果图，右下是本文过程Q=5.0,R=9高斯模糊效果图。

怎么样，效果还不错吧！

遗憾的是我没能找到按照Q自动计算模糊半径的方法，所以处理过程给出了2个参数Q和Radius。

文章中所用数据类型及一些过程见《Delphi图像处理 -- 数据类型及内部过程》和《Delphi图像处理 -- 图像像素结构与图像数据转换》。

例子使用GDI+版本下载地址和说明见《GDI+ for VCL基础 -- GDI+ 与 VCL》。

尽管我十分努力，但水平有限，错误在所难免，欢迎指正和指导。邮箱地址：

maozefa@hotmail.com

本文代码于2010.5.20重新修订过。增加了拷贝形式的调整过程和响应回调函数的调整过程。代码中的ExecuteAbort过程和ExecuteProc过程见《Delphi图像处理 -- 图像像素结构与图像数据转换》。