积分图像的应用（一）：局部标准差 分类： 图像处理 Matlab 2015-06-06 13:31 137人阅读 评论(0) 收藏

局部标准差在图像处理邻域具有广泛的应用，但是直接计算非常耗时，本文利用积分图像对局部标准差的计算进行加速。

局部标准差：

标准差定义如下（采用统计学中的定义，分母为

）：



其中

。

为了计算图像的局部标准差，首先设定局部区域的大小为

，则局部区域的像素点个数

。

对标准差的公式进行化简：



因

，故：



可以看出，局部标准差计算中需要对图像

与

进行局部求和操作，即

和

。

我们可以先分别计算出图像

和

的积分图像

、

，这样就能在常量时间计算出上述的局部和。

时间复杂度：

图像中共

个像素点，每个局部区域共有

个像素点，直接计算局部标准差的时间复杂度

，利用积分图像计算局部标准差的时间复杂度为

。（计算积分图像的时间为

、计算一个像素点的局部标准差需要常量时间

，共

个像素点，故最后的计算复杂度仍为

）

Matlab程序：

Localstd.m:直接计算图像的局部标准差

function Std=Localstd(Img,d)
[m,n]=size(Img);
Var=zeros(m,n);
Img = padarray(Img,[d d],'symmetric');%边界填充
N=(2*d+1)^2;
for i=1:m
    for j=1:n
        i1=i+d;
        j1=j+d;
        B=Img(i1-d:i1+d,j1-d:j1+d);%取局部区域
        Var(i,j)=sum(sum((B-mean(B(:))).^2))/(N-1);%标准差公式
    end
end
Std=sqrt(Var);

fastLocalstd.m:利用积分图像进行快速计算

function Std=fastLocalstd(Img,d)
[m,n]=size(Img);
Img = padarray(Img,[d+1 d+1],'symmetric');%边界填充
Img2=Img.^2;
Int=integral(Img);%计算Img的积分图像
Int2=integral(Img2);%计算Img^2的积分图像
Var=zeros(m,n);
N=(2*d+1)^2;%局部像素点个数
for i=1:m
for j=1:n
i1=i+d+1;
j1=j+d+1;
%利用积分图像求局部和
sumi2=Int2(i1+d,j1+d)+Int2(i1-d-1,j1-d-1)-Int2(i1+d,j1-d-1)-Int2(i1-d-1,j1+d);
sumi=Int(i1+d,j1+d)+Int(i1-d-1,j1-d-1)-Int(i1+d,j1-d-1)-Int(i1-d-1,j1+d);
Var(i,j)=(sumi2-sumi^2/N)/(N-1);
end
end
Std=sqrt(Var);

function I=Integral(Img)
Img=double(Img);
[m,n]=size(Img);
I=zeros(m,n);
for i=1:m
for j=1:n
if i==1 && j==1             %积分图像左上角
I(i,j)=Img(i,j);
elseif i==1 && j~=1         %积分图像第一行
I(i,j)=I(i,j-1)+Img(i,j);
elseif i~=1 && j==1         %积分图像第一列
I(i,j)=I(i-1,j)+Img(i,j);
else                        %积分图像其它像素
I(i,j)=Img(i,j)+I(i-1,j)+I(i,j-1)-I(i-1,j-1);
end
end
end

main.m:

clear all;
close all;
clc;
Img=double(imread('lena.tif'));
d=1;
tic
J1=Localstd(Img,d);
toc
tic
J2=fastLocalstd(Img,d);
toc
figure;
imshow([Img/max(Img(:)),J1/max(J1(:)),J2/max(J2(:))]);

lena.tif大小为256*256，Localstd运行时间为 0.984451s，而fastLocalstd运行时为0.025528s。

当然，对于局部标准差，Matlab本身提供了函数stdfilt，下面是它的核心代码：

function Std=Stdfilt(I,d)
I=double(I);
h=ones(2*d+1);
n=(2*d+1)^2;
conv1 = imfilter(I.^2,h,'symmetric') / (n-1);
conv2 = imfilter(I,h,'symmetric').^2 / (n*(n-1));
Std = sqrt(conv1-conv2);

可以看出，它按照化简后的公式直接对图像

和

进行卷积操作，因为imfilter函数是使用c实现的，且内部应该进行了优化，故速度较快。对上面的lena.tif运行时间为：0.064522s。