jailton_ex8

Autor: jailton
Data: 10/05/2009
 1 # 1. Estudar as funções da toolbox ia636: iacos, iadft, iadftview
 2 # 2. Estudar a demonstração iadftdecompose e procurar fazer a demonstração
 3 #    inversa, isto é, começar a reconstruir parcialmente a função original,
 4 #    adicionando-se as cossenóides, de menor frequência.
 5 # 3. Repetir parte do exercício que o Prof. Clésio demonstrou em classe usando as
 6 #    funções iadft
 7 #-------------------------------------------------------
 8 # f = iacos( s, t, theta, phi  )
 9 # s = matriz (imagem) de entrada [linha, coluna]
10 # t = período (número de pixels)
11 # theta = direção espacial da onda em radianos. theta = 0 , onda na direção horizontal.
12 # phi = fase
13 f1 = iacos([128,128],16, 0 , 0)
14 f2 = iacos([128,128],20, numpy.pi/2 , 0)
15 iashow(ianormalize(f1, [0,255]), title = " f1=iacos([128,128],16, 0 , 0)")
16 iashow(ianormalize(f2, [0,255]), title = " f2=iacos([128,128],20, numpy.pi/2 , 0)")
17 F = iadft(f1)
18 Fv = iadftview(F)
19 iashow(Fv, title = "iadft-Transformada de Fourier de f1")
20 F = iadft(f2)
21 Fv = iadftview(F)
22 iashow(Fv,title = "iadft-Transformada de Fourier de f2")
23 F = iadft(f1+f2)
24 Fv = iadftview(F)
25 iashow(Fv,title = "iadft-Transformada de Fourier de f1+f2")
26 def jcos(M,N,fase):
27     x, y = iameshgrid(range(M),range(N))
28     fs = cos(0.8*(x + y) + fase)
29     return fs
30 a=jcos(128,128,0)
31 iashow(ianormalize(a, [0,255]), title = "f(x)=cos(0.8(x+y))")
32 Fcos = iadft(a)
33 Fv = iadftview(Fcos)
34 iashow(Fv, title = "iadft-Transformada de Fourier de f(x)=cos(0.8*(x+y))")
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)

f1=iacos([128,128],16, 0 , 0)

f2=iacos([128,128],20, numpy.pi/2 , 0)

iadft-Transformada de Fourier de f1

iadft-Transformada de Fourier de f2

iadft-Transformada de Fourier de f1+f2

f(x)=cos(0.8(x+y))

iadft-Transformada de Fourier de f(x)=cos(0.8*(x+y))

 1 def fxyd(M,N):   # Discretiza f(x,y), no caso - cos(0.8*(x + y))
 2     x, y = iameshgrid(range(N),range(M))
 3     fxyd = cos(0.8*(x + y))
 4     return fxyd
 5 #---------------------------------------
 6 def jDFT2(fxyd):  # Calcula a transformada[F(u,v)] de Fourier
 7     M,N=fxyd.shape
 8     F=numpy.zeros((M,N),complex)
 9     for u in range(M):
10         for v in range(N):
11             F[u,v]=complex(0,0)
12             for x in range(M):
13                 for y in range(N):
14                     Fx=1.0/M*fxyd[x,y]*exp(-1j*2*pi*x*u/M) # Transformada-decomposta em duas partes
15                     Fy=1.0/N*exp(-1j*2*pi*y*v/N)
16                     F[u,v]=F[u,v]+Fx*Fy
17     return F
18 #---------------------------------------
19 def magnitudeF2(Fuv):   # Calcula a Magnitude da transformada|F(u,v)| de Fourier-2D
20     F=Fuv
21     M,N=Fuv.shape
22     Mag=numpy.zeros((M,N))
23     for u in range(M):
24         for v in range(N):
25             Mag[u,v]=(abs(Fuv[u,v]))
26     return Mag
27 #---------------------------------------
28 def jIDFT2(Fuv,fu1,fu2,fv1,fv2):   # Calcula a transformada inversa[f(x,y)] de Fourier-2D
29     F=Fuv
30     M,N=Fuv.shape
31     f=numpy.zeros((M,N))
32     for x in range(M):
33         for y in range(N):
34             f[x,y]=complex(0,0).real
35             for u in range(fu1,fu2):
36                 for v in range(fv1,fv2):
37                     Fu=Fuv[u,v]*exp(1j*2*pi*x*u/M)  # Transformada-decomposta em duas partes
38                     Fv=exp(1j*2*pi*y*v/N)
39                     f[x,y]=f[x,y]+(Fu*Fv).real
40     return f
41 #---------------------------------------
42 M=28
43 N=23
44 fxy1=fxyd(M,N)
45 Fuv=jDFT2(fxy1)
46 iashow(ianormalize(fxy1, [0,255]), title = "f(x,y)=cos(0.8*(x+y))-Entrada")
47 #---------------------------------------
48 fu1=0;fu2=M        # Inclui todas frequências
49 fv1=0;fv2=N
50 fxy=jIDFT2(Fuv,fu1,fu2,fv1,fv2)
51 iashow(ianormalize(fxy, [0,255]), title = "jIDFT2-f(x,y)-Não exclui frequência")
52 magnitude=magnitudeF2(Fuv)
53 iashow(ianormalize(magnitude, [0,255]), title = "jIDFT2-Magnitude de F(u,v)-Não exclui frequência")
54 Fv = iadft(fxy1)
55 Fv = iadftview(Fuv)
56 iashow(Fv, title = "Usando iadft-f(x)=cos(0.8(x+y))")
57 #---------------------------------------
58 fu1=0;fu2=0.15*M   # Exclui as frequências altas
59 fv1=0;fv2=0.15*N
60 fxy=jIDFT2(Fuv,fu1,fu2,fv1,fv2)
61 iashow(ianormalize(fxy, [0,255]),title = "jIDFT2-f(x,y)-Exclui freq. altas")
62 Fuv=jDFT2(fxy)
63 magnitude=magnitudeF2(Fuv)
64 iashow(ianormalize(magnitude, [0,255]),title = "jIDFT2-Magnitude de F(u,v)-Exclui freq. altas")
65 #---------------------------------------
66 fu1=0.9*M;fu2=M    # Exclui as frequências baixas
67 fv1=0.9*N;fv2=N
68 fxy=jIDFT2(Fuv,fu1,fu2,fv1,fv2)
69 iashow(ianormalize(fxy, [0,255]),title = "jIDFT2-f(x,y)-Exclui freq. baixas")
70 Fuv=jDFT2(fxy)
71 magnitude=magnitudeF2(Fuv)
72 iashow(ianormalize(magnitude, [0,255]),title = "jIDFT2-Magnitude de F(u,v)-Exclui freq. baixas")
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)
Warning: downcasting image from double to uint16 (may lose precision)

f(x,y)=cos(0.8*(x+y))-Entrada

jIDFT2-f(x,y)-Não exclui frequência

jIDFT2-Magnitude de F(u,v)-Não exclui frequência

Usando iadft-f(x)=cos(0.8(x+y))

jIDFT2-f(x,y)-Exclui freq. altas

jIDFT2-Magnitude de F(u,v)-Exclui freq. altas

jIDFT2-f(x,y)-Exclui freq. baixas

jIDFT2-Magnitude de F(u,v)-Exclui freq. baixas