Fonte:
SIMG-782 Introduction to Digital Image Processing - Dr. Harvey Rhody
76-3144 Carlson Center for Imaging Science - 475-6215
rhody@cis.rit.edu
http://www.cis.rit.edu/class/simg782.old/talkMoments/momentEquations.html
Algoritmos:
PRO MESHGRID,N,M,X,Y ;+ ; PRO MESHGRID,N,M,X,Y calculates two N Column x M Row arrays, X ; and Y, that can be used to calculate and plot a 2D function. ; ; USAGE: ; MESHGRID,31,31,X,Y ; X=(X-15)/2. ; Y=(Y-11)/3. ; Z=EXP(-(X^2+Y^2)) ; SURFACE,Z ; ; HISTORY: ; Originated by Harvey Rhody September 17, 1997 ;- X=INDGEN(N)#(1+INTARR(M)) Y=INDGEN(M)##(1+INTARR(N)) END
FUNCTION Moments2D,A,p,q ;+ ;m=Moments2D(A,p,q) computes the Mpq moment of A. p and q may be ;vectors with a list of moment exponents. In that case, the result ;is a vector of values for the items in the p and q lists. The meshgrid ;procedure is used to construct X and Y arrays. The pq moment is ;computed by summing (X^p)*(Y^q)*A with the TOTAL function. ; ;USES meshgrid to compute the X and Y fields. ; ;HISTORY ;H. Rhody, Sept. 1999 ;- ;Construct the X and Y grids of the same ;size as the image.
aSiz=SIZE(A) meshGrid,aSiz[1],aSiz[2],X,Y X=FLOAT(X) & Y=FLOAT(Y) n=N_ELEMENTS(p)
;If p and q have one element then one calculation is needed.
IF n EQ 1 THEN RETURN,TOTAL(X^p*Y^q*A)
;If p and q are vectors, then compute a result for
;each value pair in the list.
IF n GT 1 THEN BEGIN m=FLTARR(n) FOR k=0,n-1 DO m[k]=TOTAL(X^p[k]*Y^q[k]*A) ENDIF RETURN,m END
FUNCTION CentralMoments2D,A,p,q ;+ ;m=CentralMoments2D(A,p,q) computes the Mpq central moment of A. ;p and q may be vectors with a list of moment exponents. In ;that case, the result is a vector of values for the items in ;the p and q lists. ; ;USES meshgrid to compute the X and Y fields. ; ; ;HISTORY ;H. Rhody, Sept. 1999 ;- ;Construct the X and Y arrays the same size as A. aSiz=SIZE(A) meshGrid,aSiz[1],aSiz[2],X,Y X=FLOAT(X) & Y=FLOAT(Y) ;Calculate the normalization factor and the centroid values. m00=TOTAL(A) xbar=TOTAL(X*A)/m00 ybar=TOTAL(Y*A)/m00 n=N_ELEMENTS(p) ;If p and q are scalar values, then return the related result IF n EQ 1 THEN RETURN,TOTAL((X-xbar)^p*(Y-ybar)^q*A) ;If p and q are vectors, then return a vector with values for ;each (p,q) in the list. IF n GT 1 THEN BEGIN m=FLTARR(n) FOR k=0,n-1 DO m[k]=TOTAL((X-xbar)^p[k]*(Y-ybar)^q[k]*A) ENDIF RETURN,m END
FUNCTION NormalMoments2D,A,p,q ;+ ;m=NormalMoments2D(A,p,q) computes the Mpq normalized central moment of A. ;p and q may be vectors with a list of moment exponents. In ;that case, the result is a vector of values for the items in ;the p and q lists. ; ;USES meshgrid to compute the X and Y fields. ; ; ;HISTORY ;H. Rhody, Sept. 1999 ;- ;Construct the X and Y arrays the same size as A.
aSiz=SIZE(A) meshGrid,aSiz[1],aSiz[2],X,Y X=FLOAT(X) & Y=FLOAT(Y)
;Calculate the normalization factor and the centroid values. m00=TOTAL(A) xbar=TOTAL(X*A)/m00 ybar=TOTAL(Y*A)/m00
;Calculate the exponent for the normalization term. n=N_ELEMENTS(p) gamma=FLOAT(p+q)/2+1
;If p and q are scalar values, then return the related result IF n EQ 1 THEN RETURN,TOTAL((X-xbar)^p*(Y-ybar)^q*A)/m00^gamma
;If p and q are vectors, then return a vector with values for ;each (p,q) in the list. IF n GT 1 THEN BEGIN m=FLTARR(n) FOR k=0,n-1 DO m[k]=TOTAL((X-xbar)^p[k]*(Y-ybar)^q[k]*A) ENDIF RETURN,m/m00^gamma END
FUNCTION InvariantMoments2D,A ;+ ;v=InvariantMoments2D(A) computes the seven invariant moment ;features described in Gonzalez and Woods, page 516. ; ;v[k]=phi[k+1], k=0,1,...,6 where phi[n] is element n in ;Gonzalez and Woods notation. ; ;USES NormalMoments2D. ; ;HISTORY ;H. Rhody, Sept. 1999 ;- ;Set up the p and q vectors to create the vector of normal moment ;results to support the calculation of the invariants. p=[0,0,1,1,2,2,3] q=[2,3,1,2,0,1,0] m=NormalMoments2D(A,p,q) ;Compute the invariants. Note that m[0]=m02, m[1]=m03, etc. from the ;order of coefficients in the p and q vectors. The invariants phi_1, ;phi_2,...,phi_7 of Gonzalez and Woods are the elements v[],v[1], etc. ;of the v vector. v=FLTARR(7) v[0]=m[0]+m[4] v[1]=(m[4]-m[0])^2+4*m[2]^2 v[2]=(m[6]-3*m[3])^2+(3*m[5]-m[1])^2 v[3]=(m[6]+m[3])^2+(m[5]+m[1])^2 v[4]=(m[6]-3*m[3])*(m[6]+m[3])*((m[6]+m[3])^2- $ 3*(m[5]+m[1])^2) +(3*m[5]-m[1])*(m[5]+m[1])* $ (3*(m[6]+m[3])^2-(m[5]+m[1])^2) v[5]=(m[4]-m[0])*((m[6]+m[3])^2-(m[5]+m[1])^2)+ $ 4*m[2]*(m[6]+m[3])*(m[5]+m[1]) v[6]=(3*m[5]-m[6])*(m[6]+m[3])*((m[6]+m[3])^2 $ -3*(m[5]+m[1])^2)+(3*m[3]-m[6])*(m[5]+m[1])* $ (3*(m[6]+m[3])^2-(m[3]+m[1])^2) RETURN,v END