图像特征提取方法-ZerNike程序

function [A_nm,zmlist,cidx,V_nm] = P_zernike(img,n,m) if nargin>0 if nargin==1 n = 0; end d = size(img); % dimension img = double(img); xstep = 2/(d(1)-1); ystep = 2/(d(2)-1); [x,y] = meshgrid(-1:xstep:1,-1:ystep:1); % 归一化区域 circle1 = x.^2 + y.^2; inside = find(circle1<=1); mask = zeros(d); mask(inside) = ones(size(inside)); [cimg,cidx] = clipimg(img,mask) z = clipimg(x+i*y,mask); % i 为虚数 p = 0.9*abs(z); %make sure that the p less than 0.9 according to the theory of p<=0.9; theta = angle(z); c = 1; for order=1:length(n) n1 = n(order) if nargin<3 m = 0:n1; end for r=1:length(m) V_nmt = zpoly(n1,m(r),p,theta); % pzpoly -> zpoly zprod = cimg.*conj(V_nmt); % 共轭 A_nm(c) = (n1+1)*sum(sum(zprod))/pi; zmlist(c,1:2) = [n1 m(r)]; if nargout==4 V_nm(:,c) = V_nmt; end c = c+1; end end else end function [cimg,cindex,dim] = clipimg(img,mask) dim = size(img); cindex = find(mask~=0); cimg = img(cindex); return; function [V_nm,mag,phase] = zpoly(n,m,p,theta) R_nm = zeros(size(p)); a = (n+abs(m))/2; b = (n-abs(m))/2; total = b; for s=0:total num = ((-1)^s)*fac(n-s)*(p.^(n-2*s)); den = fac(s)*fac(a-s)*fac(b-s); R_nm = R_nm + num/den; end mag = R_nm; phase = m*theta; V_nm = mag.*exp(i*phase); return; function [V_nm,mag,phase] = pzpoly(n,m,p,theta) R_nm = zeros(size(p)); a = n+abs(m)+1; b = n-abs(m); total = b; for s=0:total num = ((-1)^s)*fac(2*n+1-s)*(p.^(n-s)); den = fac(s)*fac(a-s)*fac(b-s); R_nm = R_nm + num/den; end mag = R_nm; phase = m*theta; V_nm = mag.*exp(i*phase); return; function [factorial] = fac(n) maxno = max(max(n)); zerosi = find(n<=0); n(zerosi) = ones(size(zerosi)); factorial = n; findex = n; for i=maxno:-1:2 cand = find(findex>2); candidates = findex(cand); findex(cand) = candidates-1; factorial(cand) = factorial(cand).*findex(cand); end return;