(摘自冈萨雷斯的《数字图像处理(使用Matlab)》)

：

l function [P,npix]=histroi(f,c,r)

% HISTROI Computes the histogram of an ROI

in an image.

% [P,npix]=histroi(f,c,r) computes the

histogram,P,of a

% polygonal region of interest (ROI) in

image F.The polygonal

% region is defined by the colunm and row

coordinates of its

% vertices, which are by the column and

row coordinates of its

% respectively. All pixels of F must be

>=0. Parameter NPIX is

% number of pixels in the polygonal

region.

% ****from ggbondg****

% Generate the binary mask

image.

B=roipoly(f,c,r);

% Compute the histogram of the pixels in

the ROI.

P=imhist(f(B));

%Obtain the number of pixels in ROI if

requested in the output.

if nargout>1

npix=sum(B(:));

end

%****from ggbondg****

l function [v,unv]=statmoments(p,n)

% STATMOMENTS Computes statistical central

moments of image histogram.

% [v,unv]=statmoments(p,n) computes up to

the Nth statistical

% central moment of a histogram whose

components are in vector

% P. The length of must equal 256 or

65536.

% ****from ggbondg****

% The program outputs a V with V(1)=mean,

V(2) = variance.

% V(3) = 3rd moment,...V(N)=Nth central

moment. The random

% variable values are normalized to the

range [0,1], so all

% moments also are in this

range.

% ****from ggbondg****

% The program also outputs a Vector UNV

containing the same moments

% as V,but using un-normalized random

variable values (e.g., 0 to

% 255 if length(P)=2^8). For example, if

length(P)=256 and V(1)

% = 0.5, then UNV(1) would have the value

UNV(1)=127.5 (half of

% the [0 255] range).

%****from ggbondg****

Lp=length(p);

if (Lp~=256)&(Lp~=65536)

error('P

must be a 256- or 65536- elements vector.');

end

G=Lp-1;

% Make sure the histogram has unit area,

and convert it to a

% column vector.

p=p/sum(p);

p=p(:);

% Form a vector of all the possible values

of the

% random variable.

z=0:G;

%****from ggbondg****

% Now normalize the z's to the range

[0,1].

z=z./G;

%****from ggbondg****

% The mean.

m=z*p;

%****from ggbondg****

% Cencter random variable about the

mean.

z=z-m;

%****from ggbondg****

% Compute the central moments.

v=zeros(1,n);

v(1)=m;

for j=2:n

v(j)=(z.^j)*p;

end

%****from ggbondg****

if nargout>1

% Compute the uncentralized

moments.

unv=zeros(1,n);

unv(1)=m.*G;

for

j=2:n

unv(j)=((z*G).^j)*p;

end

end

%****from ggbondg****

l function [H]=lpfilter(type,M,N,D0,n)

% LPFILTER Computes frequency domain

lowpass filters.

% H=lpfilter(type,M,N,D0,n) creates the

transfer function of

% a lowpass filter ,H,of the specified

TYPE and size (M-by-N).To

% view the filter as an image or mesh

plot,it should be centered

% using H=fftshift(H).

% %****from ggbondg****

% Valid values for TYPE,D0,and n

are:

% %****from ggbondg****

%

'ideal' Ideal lowpass filter with cutoff frequency D0. n need

% not be supplied. D0 must be positive.

% %****from ggbondg****

%

'btw' Butterworth lowpass filter of order n, and cutoff

% D0.The default value for n is 1.0. D0 must be positive.

%

% 'gasussian' Gaussian

lowpass filter with cutoff (standard

% deviation) D0. n need not be supplied. D0 must be

% positive.

% %****from ggbondg****

% use function dftuv to set up the

meshgrid arrays needed for

% computing the required

distances.

[U,V]=dftuv(M,N);

%****from ggbondg****

% Compute the distances D(U,V).

D=sqrt(U.^2+V.^2);

%****from ggbondg****

% Begin filter computations.

switch type

case

'ideal'

H=double(D<=D0);

case

'btw'

if nargin==4

n=1;

end

H=1./(1+(D./D0).^(2*n));

case

'gaussian'

H=exp(-(D.^2)./(2*(D0^2)));

otherwise

error('Unknown filter type.')

end

%****from ggbondg****

l function [U,V]=dftuv(M,N)

% DFTUV Computes meshgrid frequency

matrices.

% [U,V]=dftuv(M,N) computes meshgrid

frequency matrices U and V.

% U and V are useful for computing

frequency-domain filter

% functions.U and V are both

M-by-N.

% %****from ggbondg****

% Set up range of variables.

u=0:M-1;

v=0:N-1;

%****from ggbondg****

% Compute the incides for use in

meshgrid.

idx=find(u>M/2);

u(idx)=u(idx)-M;

idy=find(v>N/2);

v(idy)=v(idy)-N;

%****from ggbondg****

% Compute the meshgrid arrays.

[V,U]=meshgrid(v,u);

%****from ggbondg****