Fourier series to Fourier transform Masahiro Yamamoto September 9, 2016 OB (r j)j r (r i)i Figure 1: normal coordinate, projection, inner product 3 r

Size: px

Start display at page:

Download "Fourier series to Fourier transform Masahiro Yamamoto September 9, 2016 OB (r j)j r (r i)i Figure 1: normal coordinate, projection, inner product 3 r"

はなさどひら
5 years ago
Views:

1 Fourier series to Fourier transform Masahiro Yamamoto September 9, 2016 OB (r j)j r (r i)i Figure 1: normal coordinate, projection, inner product 3 r i.j, k x = r i, y = r j, z = r k r = xi + yj + zk N Fig.1 i, j, k i j = j k = k i = 0 i i = j j = k k = 1 1 f(x) g n (x) f(x) g n (x) (f, g n ) = b a dxf(x)g n(x) (1) f(x) [ l, l] 2l f(x) = n (f, g n )g n = 1 n= 2l n= c ne inπx/l (2) g n = einπx/l 2l (3) c n = (f, g n ) = 1 2l l l dxf(x)e inπx/l (4) 1

2 g n (g n, g m ) = 1 2l l l dxe inπx/l e imπx/l = { 1, (n = m) sin(n m)π 2(n m)π = 0, (n m) c n = c n/ 2l f(x) = (g n, g m ) = δ nm c n = 1 2l n= n= l l c n e inπx/l (5) dxf(x)e inπx/l (6) f(x + 2l) = n= n= c n e inπx/l e i2nπ = f(x) Fig.2 n 1.0 n= n=10 n=100 n= n=1 n=10 n=100 n= Figure 2: Examples of Fourier series Eq.(6) f(x) Eq.(5) c n = 1 2l = l l m= m= dx Eq.(5) c n Eq.(6) m= m= c m δ mn = c n c m e imπx/l e inπx/l f(x) = = n= n= l 1 l dx f(x )e inπx /l e inπx/l 2l l l dx f(x ) 1 2l n= e inπ(x x)/l n= 2

3 = = = l dx f(x ) 1 1 l 2l k l l dx f(x ) 1 2l l l l π n= e ik n(x x) k, n= dke ik(x x) dx f(x ) 1 2π 2πδ(x x) = f(x) ( k n = nπ l, k = π ) l Eq.(5)Eq.(6) L 2l 2 f(x) = c n = 1 L f(x + L) = + n= L/2 + L/2 n= c n e i2πnx/l (7) dxf(x)e i2πnx/l (8) c n e i2πnx/l e i2nπ = f(x) 2l k n nπ/l k, k = π/l 1 f(x) = lim l k f(k) c n / k f(x) =,f(x) f(k) f(x) = = n= n= c n e inπx/l k = lim n= l n= c n k e ik nx k (9) dkf(k)e ikx (10) f(k) = c n 1 l l = lim dxf(x)e iknx = 1 dxf(x)e ikx (11) k l 2l π l 2π = 1 2π dkf(k)e ikx dk 1 2π dx f(x ) 2.1 Gauss f(x) = A 2 πα f(k) = 1 2π = A π 2π dx f(x )e ikx e ikx dke ik(x x) = 1 dx f(x )2πδ(x x) = f(x) 2π f(x) = A exp( αx 2 ), α > 0 (12) α e k 2 dke k2 4α e ikx = kx α(x+i dxe 2α) 2 e k2 4α dxae αx2 e ikx = A 2π A 4α = 2 πα e k 2 4α (13) A 2 dke 1 4α (k 2iαx)2 e αx2 = πα A 2 πα 2 απe αx2 = Ae αx2 3

4 2.2 Gaussian wavepacket f(k ) = 1 2π = A [ 2π exp f(x) = = A 2π exp f(x) = Ae ikx e x2 /2 2 (14) dxae ikx e x2 /2 2 e ik x = A /2)(k k ) 2 ] 2π e[ ( (k k ) 2 = A [ exp 2 2π 2 k2 ] 2π 2 = A [ exp 2π 2 2 (k k ) 2 dk f(k )e ik x = A [ ] + dk exp 2 2π 2 (k k ) 2 [ ] [ k2 exp 2 (k + i ] [ + 2 x)2 dk exp ] [ 2 exp 2 (k + i ] 2π 2 x)2 2 = Ae ikx e x2 /(2 2 ) Ref(x) = A [ ] + dk exp 2 2π 2 (k k ) 2 cos(k x) = A [ ] (δk ) exp 2 2π 2 (k k ) 2 cos(k x) k = k k Ref(x) = A [ ] (δk ) exp 2 2π 2 k 2 cos[(k k )x] k e ik x 2 2 dxe [x i 2 (k k )] 2 ] { k (k + i } ] 2 x)2 k = k k = ndk, (n = N, N + 1,..., 2, 1, 0, 1, 2,..., N 1, N) (15) (16) 2.3 delta function dke ik(x x) L = lim dke ik(x x) 1 = lim x) L L L i(x x) [e il(x e il(x x) ] = lim L 2 sin[l(x x)] (x x) = Cδ(x x) C, X = x x 2 sin(lx)/x 2 sin(lx) 2 sin(lx) lim = 2L, = 0 at X = ± π X 0 X X L L (Fig.3) C Fig.4 C = 2π C dxδ(x) = C = = 2Im dx 2 sin(lx) X dz eilz z = 2Imiπ = 2π = 2Im dke ik(x x) = 2πδ(x x) dz eiz Z, (Z = zl) 4

5 !" )./01$/2/-%-&-&34561*/3 ))7)') """"" ))7)') """" ))7)') """ "!# "!" $"!# &'#()*' ()+*,' -&- """" $%" $ " " " %" x!",- #### $%&'(!# #!",- ##### #!# )* ) # * &'+ Figure 3: FT of Gaussian wave-packet 7 6 2*sin(4*x)/x 2*sin(2*x)/x 2*sin(x)/x Figure 4: functional form of 2 sin(lx)/x 5

6 Imz Rez Figure 5: Contour to calculate e iz /z 2.4 step function, step function ( step function θ(x x 0 ) = { 1, x > x0 0, x < x 0 (17) θ(k) θ(k) = 1 2π dxθ(x x 0 )e ikx = 1 dxe ikx (18) 2π x 0 x = + e ɛx ɛ [ ] 1 + θ(k) = lim dxe i(k iɛ)x 1 e i(k iɛ)x + 1 e i(k iɛ)x 0 = lim = lim (19) ɛ 0 2π x 0 ɛ 0 2π i(k iɛ) ɛ 0 2π i(k iɛ) = i 2π lim e i(k iɛ)x 0 ɛ 0 k iɛ lim ɛ 0 1 x ± iɛ x 0 (20) = P 1 iπδ(x) (21) x P 0-0+ step function θ(k) = θ(k; x 0 ) = i [P 1 ] 2π k + iπδ(k) e ikx 0 (22) θ(k) x 0 θ(k; x 0 ) 1 θ(x) 0 x 0 x Figure 6: step function θ(x x 0 ) θ(k; x 0 ) θ(x x 0 ) θ(x x 0 ) = dkθ(k; x 0 )e ikx = i 2π = i 2π P dk eik(x x 0) + i k = 1 2 i 2π P dk eik(x x 0) k 2π iπ dk [P 1 ] k + iπδ(k) e ik(x x 0) dkδ(k)e ik(x x 0) (23) (24) (25) 6

7 k k k = k +ik x x 0 > 0 e ik(x x 0) = e ik (x x 0 ) e k (x x 0 ) Fig.7 C x x 0 < 0 C x x 0 > 0 Fig.7 k = 0 0 = = P πi = P 0 dk eik(x x 0) k dk eik(x x 0) k dk eik(x x 0) k + dk C1 eik(x x 0) + dk eik(x x 0) C + dk eik(x x 0) k 0+ k k Fig.7 k = 0 (26) πie i0(x x 0) + 0 (27) (28) 2πie i0(x x 0) = P dk eik(x x 0) + πie i0(x x0) + 0 (29) k θ(x x 0 ) = 1 2 i 2π P dk eik(x x 0) k = 1 2 i 2π πi = 1, x x 0 > 0 (30) x x 0 < 0 Fig.7 k = 0 2πie i0(x x 0) πi = P = P dk eik(x x 0) k dk eik(x x 0) k πie i0(x x 0) + 0 (31) (32) k = 0 0 = P dk eik(x x 0) + πie i0(x x 0) k (33) θ(x x 0 ) = 1 2 i 2π P dk eik(x x 0) k = 1 2 i 2π ( πi) = 0, x x 0 < 0 (34) θ(k; x 0 ) θ(x x 0 ) = { 1, x > x0 0, x < x 0 (35) 2.5 Convolution: f(t) r(t) t r(t) = χ f(t) + φ(t t )f(t )dt (36) χ φ(t t ), t t t > t φ(t t ) = 0, t t < 0 r(t) = χ f(t) + φ(t t )f(t )dt (37) 7

8 C1 C1 C' C'' C' C2 C'' C2 Figure 7: contour to calculate θ(x x 0 ), x > x 0 (left) and θ(x x 0 ), x < x 0 (right) from θ(k), ω r(t)e iωt dt = χ f(t)e iωt dt + r(ω) = χ f(ω) + φ(t t )f(t )dt e iωt dt (38) φ(t t )f(t )e iω(t t ) e iωt dt dt (39) +, T = t t t χ(ω) = χ + φ(ω) r(ω) = χ f(ω) + φ(t )e iωt dt f(t )e iωt dt (40) = χ f(ω) + φ(ω)f(ω) (41) (convolution) φ(t) f(t) = r(ω) = χ(ω)f(ω) (42) φ(t t )f(t )dt (43) φ(t) f(t) φ(ω)f(ω) 2.6 (space-time) k Fourier e ikx (1D), e ikxx e ikyy e ikzz = e ik r (3D) Fourier e ikx (1D), e ikxx e ikyy e ikzz = e ik r (3D) e iωt, e iωt ω 1 F (k, ω) = (2π) 4 dr dtf(r, t) exp[ i(k r ωt)] (44) f(r, t) = dk dωf (k, ω) exp[i(k r ωt)] (45) 3 Fourier ρ(r) Poisson φ(r) [ɛ(r) φ(r)] = ρ(r) ɛ 0 (46) 8

9 if ɛ(r) is position-indepedent 2 φ(r) = ρ(r) ɛ 0 ɛ Fourier Fourier e ikx d n [f(x)/dx n ]dx = (ik) n e ikx f(x)dx k 2 φ(k) ρ(k) ɛ 0 ɛ (47) = 0 (48) φ(k) = ρ(k) ɛ 0 ɛk 2 (49) φ(r) Fourier Fourier 9

10 4 Fourier transform of the periodic system The structure of all crystal can be described in terms of lattice with a group of atoms attached to every lattice point. The group of atoms is called the basis. crystal structure = lattice{r} + basis{τ} (50) A lattice translation vector can be described by the primitive cell vectors R = n 1 a 1 + n 2 a 2 + n 3 a 3 (51) Any function f invariant under a lattice translation R may be expanded in a Fourier series of the form f(r) = G f(g) exp(ig r) (52) f(r + R) = f(r) (53) exp(ig R) = 1 (54) 1 f(g) = drf(r)e ig r = 1 drf(r)e ig r (55) Ω cell cell Ω f(g) = f(g 1 ) dre i(g G) r = f(g) (56) Ω G cell cell }{{} δ G,G Ω cell f(r) = 1 dr f(r )e ig r e ig r (57) Ω G cell cell 1 = dr f(r 1 ) e ig (r r ) (N cell Ω cell = Ω) (58) N cell Ω cell G }{{} Ωδ r,r e ig (r r ) G = f(r) (59) 1 = e ig (r r ) G x G y G z = Ω G x G y G z 8π 3 dge ig (r r ) = 8π3 Ω 8π 3 δ(r r ) If we define the reciprocal lattice( The primitive translation vectors of the reciprocal lattice are G a i b j = 2πδ ij (60) b 1 = 2π a 2 a 3 a 1 a 2 a 3 b 2 = 2π a 3 a 1 a 1 a 2 a 3 b 3 = 2π a 1 a 2 a 1 a 2 a 3 (61) Here a 1, a 2, a 3 are the primitive translation vectors of the crystal lattice. A reciprocal lattice vector has the form G = hb 1 + kb 2 + lb 3, (62) where h, k, l are integers or zero. This equation satisfies exp(ig R) = 1. If we consider a 1 /h, a 2 /k, a 3 /l plane (hkl Miller index: ), the reciprocal lattice vector G hkl = hb 1 + kb 2 + lb 3 is perpendicular to this plane. The distance between the two adjacent parallel planes of the lattice is d(hkl) = 2π/ G hkl. 1 1 proof: An arbitary vector t in the (hkl) plane can be described The inner product t = α(a 1 /h a 2 /k) + β(a 2 /k a 3 /l) (63) G hkl t = 2πα(h/h k/k) + 2πβ(k/k l/l) 10

11 5 Electrostatic Potential The Poisson equation is given in the SI unit if ɛ(r) is position-indepedent [ɛ(r) φ(r)] = ρ(r) ɛ 0 (66) 2 φ(r) = ρ(r) ɛ 0 ɛ The Fourier transform of the above equation becomes φ and ρ is [ G 2 G G 2 φ(g) ρ(g) ɛ 0 ɛ (67) φ(g) exp(ig r) = 1 ρ(g ) exp(ig r) (68) ɛ 0 ɛ G ] exp(ig r) = 0 (69) φ(g) = ρ(g) ɛ 0 ɛg 2 (70) Then we will back-fourier transform of the rhs of the above equation, we can get φ(r) In the 3-d FFT program f(g) = 1 Ω cell cell drf(r)e ig r is forward transformation and f(r) = f(g) exp(ig r) is backward transformation. G 6 Electrostatic Free Energy E elst = 1 drρ(r)φ(r) + drρ(r)v ext (r) (71) 2 = Ω ρ (G)φ(G) + Ω ρ (G)V ext (G) (72) 2 G G divd = ρ, D = ɛ 0 ɛ(r)e(r), E = gradφ(r) (73) ɛ(r) φ(r) = ρ(r) (74) ɛ 0 ρ(g ) exp(ig r) = ɛ 0 ɛ(g) exp(ig r) φ(g ) exp(ig r) (75) G G G = 0 (64) Thereby G hkl is perpendicular to the (hkl) plane. If we define the vector d is the one from the origin to the intersection of (hkl) plane and the vector G d = G hkl G hkl d(hkl) G hkl G hkl d = d(hkl) d = a1 h + tintersect G hkl d = G hkl a 1 + G hkl t intersect h = 2πh + 0 = 2π then h 2π d(hkl) = G hkl (65) 11

12 2D-case :MPA on Au(111) f(r) = f(r, z) (76) f(r + R, z) = f(r, z) (77) f(r, z) = f(g, z)e ig r G (78) f(g, z) = 1 Ω = 1 Ω dr f(r, z)e ig r (79) dr f(g, z)ei(g G ) r (80) φ(r, z) = = G G f(g, z) 1 Ω dr e i(g G ) r = f(g, z) (81) } {{ } Ω δ G,G ( i x + j y + k ) φ(g z, z)e i(g xx+g yy) G = (ig )φ(r, z) + k φ(g, z) e i(g xx+g y y) z G (82) (83) Then the PB equation becomes ɛ(z) φ(r, z) = ρ(r, z) ɛ 0 ρ(g, z) e ig r = (i ɛ G 0 x + j y + k [ z ) ɛ(z) (i(ig x ) + j(ig y ))φ(g, z)e i(gxx+gyy) G +k φ(g ], z) e i(g xx+g y y) z = G [ ɛ(z)g 2 φ(g, z) dɛ(z) dz φ(g, z) z (84) (85) ] ɛ(z) 2 φ(g, z) z 2 e ig r (86) ρ(g, z) = ɛ(z)g 2 ɛ φ(g, z) dɛ(z) φ(g, z) ɛ(z) 2 φ(g, z) 0 dz z z 2 (87) ρ(r) = i z i en i (r) + R,I ez I δ(r R r I ) + σ M δ(z) (88) ɛ 0 ɛ 2 φ(z + ) z D + z D z = σ 12 (89) ɛ 0 ɛ 2 E + z ɛ 0 ɛ 1 E z = σ 12 (90) + ɛ 0 ɛ 1 φ(z ) z = σ 12 (91) 12

13 7 FFT ( ) : ( modified with Numerical Redcipe by MY FFT FFT FFT 7.1 Fourier FFT Fourier Fourier Fourier FFT Fourier DFT Fourier DFT DFT Fourier DFT Fourier DFT m N DFT N m 1 N e CN N Fourier f max N 2f max ( ) N Fourier DFT DFT Fourier Fourier DFT Fourier DFT ( ) DFT Fourier DFT N DFT N N DFT DFT Fourier 2 H(ω) = = H(f) = h(t) = 1 2π = 2 dth(t)e iωt (92) dth(t)e i2πft, ω = 2πf (93) dωh(ω)e iωt (94) dfh(f)e i2πft (95) f(k, ω) = f(r, t) = 1 drdtf(r, t)e i(k r ωt) (2π) 4 dkdωf(k, ω)e i(k r ωt) (1/2π) e i.. e i.. 13

14 H(f) T/2 T/2 = t j=0 = t j=0 dth(t)e i2πft = t h(t)e i2πft (96) h( T/2 + j t )e i2πf( T/2+j t) (97) h(j t )e i2πfj t [t t + T/2] (98) t = T N, t(j) = j t, (j = 0, 1, 2, 3,..., N 2, N 1) (99) f min = 1 T = f, (f = 0, f, 2 f, 3 f,..., f max ) (100) f max = f max /f min = N 2T T = N 2 1 = N, Nyquist critical frequency. See Fig.2 below. (101) 2 t 2T (102) f(k) = k f = k 1 T, (k = N 2, N 2 + 1,...N 2 2, N 1) (103) 2 (N = even is assumed.) (104) H(f(k)) = t H k j=0 h(j t )e i2πjk t f = t j=0 h(j t )e ( i2π/n)jk (105) = t H k (106) j=0 h j e ( i2π/n)jk = j=0 h j (W N ) jk, W N e i2π/n (107) H k = H k+n, (W N ) k+n = (W N ) k (W N ) N = (W N ) k }{{} (108) =1 H k+n = H k, (W N ) k+n = (W N ) k (W N ) N = (W N ) k }{{} (109) =1 h(t(j)) = f H(f(k))e 2πif(k)t(j) (110) k=0 = f = 1 N k=0 k=0 t H k e 2πikj f t (111) H k e (2πi/N)kj = 1 N k=0 H k (W N ) kj (112) h j+n = h j, h j = h j+n (113) DFT DFT Fourier - Discrete Fourier Transform Fourier Fourier Fourier (a 0, a 1,..., a N 1 ) (A 0, A 1, A 2,..., A N 1 ) N DFT A k = j=0 a j (W N ) jk, W N = e 2πi/N (114) 14

15 cos(x) Figure 8: (left) Critical sampling of a cosine wave is two sample points per cycle. (peak to peak). Nyquist critical wavelength is 2 t. (right) DFT mesh a k = 1 N j=0 A j (W N ) jk (115) DFT (1/N, (1/N)) FFT DFT A k N W N N = 1 A k A N k A N+k A k N DFT 0 N 1 N/2 N/2 1 ( ) DFT A k = j=1 a j W (j+δ 1)(k+δ 2 ) N, W N = e 2πi/N (116) δ 1, δ 2 1/2 Odd DFT DFT DFT FFT DFT Fourier Fourier / a j A k A k = = j=1 j=1 a j(e 2πi/N ) jk (117) a j (e 2πi/N ) j( k) (118) = A k = A N k (119) 15

16 A 0, A N/2 a j A k a j (a j = a N j ) A k a j (a j = a N j ) A k 1 DFT Fourier DFT (Digital Convolution) a j h j N y k = j=0 a j h k j (120) a j h j y j a j h j NN Fourier DFT y j DFT Y k a j,h j DFT A k, H k Y k = A k H k (121) Y k DFT FFT Prime Factor FFT 8 Cooley-Tukey FFT FFT 1965 J.W.Cooley J.W.Tukey 3 FFT FFT N Fourier N 2 FFT N log(n) 10 9 N = 2 30 Fourier 30 3 FFT Cooley Tukey FFT 8.1 N DFT A k = j=0 a j (W N ) jk, W N = e 2πi/N (122) A 0 A N 1 N N 2 N 2 k N DFT N/2 DFT A 2k = A 2k+1 = N/2 1 j=0 N/2 1 j=0 (a j + a N/2+j )(W N/2 ) jk (123) [(a j a N/2+j )(W N ) j ](W N/2 ) jk (124) (125) N/2 DFT N 2 /4 2, 3,... 16

17 1/4, 1/8,... Cooley-Tukey FFT ( 2, Cooley-Tukey FFT) N = 8 (N 1 ) W 8 = e iπ/4 = U A 0 A 1 A 2 A 3 A 4 A 5 A 6 A 7 A k = = = j=0 (W N ) kj a j, W N = e 2πi/N W 0,0 8 W 0,1 8 W 0,2 8 W 0,3 8 W 0,4 8 W 0,5 8 W 0,6 8 W 0,7 8 W 1,0 8 W 1,1 8 W 1,2 8 W 1,3 8 W 1,4 8 W 1,5 8 W 1,6 8 W 1,7 8 W 2,0 8 W 2,1 8 W 2,2 8 W 2,3 8 W 2,4 8 W 2,5 8 W 2,6 8 W 2,7 8 W 3,0 8 W 3,1 8 W 3,2 8 W 3,3 8 W 3,4 8 W 3,5 8 W 3,6 8 W 3,7 8 W 4,0 8 W 4,1 8 W 4,2 8 W 4,3 8 W 4,4 8 W 4,5 8 W 4,6 8 W 4,7 8 W 5,0 8 W 5,1 8 W 5,2 8 W 5,3 8 W 5,4 8 W 5,5 8 W 5,6 8 W 5,7 8 W 6,0 8 W 6,1 8 W 6,2 8 W 6,3 8 W 6,4 8 W 6,5 8 W 6,6 8 W 6,7 8 W 7,0 8 W 7,1 8 W 7,2 8 W 7,3 8 W 7,4 8 W 7,5 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 1 U 2 U 3 U 4 U 5 U 6 U 7 U 0 U 2 U 4 U 6 U 8 U 10 U 12 U 14 U 0 U 3 U 6 U 9 U 12 U 15 U 18 U 21 U 0 U 4 U 8 U 12 U 16 U 20 U 24 U 28 U 0 U 5 U 10 U 15 U 20 U 25 U 30 U 35 U 0 U 6 U 12 U 18 U 24 U 30 U 36 U 42 U 0 U 7 U 14 U 21 U 28 U 35 U 42 U 49 (W N ) j = (W N ) j N = (W N ) j 2N... A 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 A 1 U 0 U 1 U 2 U 3 U 4 U 5 U 6 U 7 A 2 U 0 U 2 U 4 U 6 U 0 U 2 U 4 U 6 A 3 A = U 0 U 3 U 6 U 1 U 4 U 7 U 2 U 5 4 U 0 U 4 U 0 U 4 U 0 U 4 U 0 U 4 A 5 U 0 U 5 U 2 U 7 U 4 U 1 U 6 U 3 A 6 U 0 U 6 U 4 U 2 U 0 U 6 U 4 U 2 A 7 U 0 U 7 U 6 U 5 U 4 U 3 U 2 U 1 8 W 7,6 8 W 7, (= N N) a 0 A 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 a 1 a A 2 A 4 = U 0 U 2 U 4 U 6 U 0 U 2 U 4 U 6 2 U 0 U 4 U 0 U 4 U 0 U 4 U 0 U 4 a 3 a 4 A 6 U 0 U 6 U 4 U 2 U 0 U 6 U 4 U 2 a 5 a 6 a 7 U 1 A 0 U 0 U 0 U 0 U 0 a 0 + a 4 A 2 A 4 = U 0 U 2 U 4 U 6 a 1 + a 5 U 0 U 4 U 0 U 4 a 2 + a 6 A 6 U 0 U 6 U 4 U 2 a 3 + a 7 a 0 a 1 a 2 a 3 a 4 a 5 a 6 a 7 a 0 a 1 a 2 a 3 a 4 a 5 a 6 a 7 a 0 a 1 a 2 a 3 a 4 a 5 a 6 a 7 U 2 = V (= W 4 ) A 0 A 2 A 4 = A 6 V 0 V 0 V 0 V 0 V 0 V 1 V 2 V 3 V 0 V 2 V 0 V 2 V 0 V 3 V 2 V 1 a 0 + a 4 a 1 + a 5 a 2 + a 6 a 3 + a 7 17

18 A 0 A 2 A 4 A 6 = A 2k = N/2 1 j=0 W 0,0 4 W 0,1 4 W 0,2 4 W 0,3 4 W 1,0 4 W 1,1 4 W 1,2 4 W 1,3 4 W 2,0 4 W 2,1 4 W 2,2 4 W 2,3 4 W 3,0 4 W 3,1 4 W 3,2 4 W 3,3 4 a 0 + a 4 a 1 + a 5 a 2 + a 6 a 3 + a 7 (W N/2 ) kj (a j + a N/2+j ) (126), A 1 A 3 A 5 A 7 = = = U 0 U 1 U 2 U 3 U 4 U 5 U 6 U 7 U 0 U 3 U 6 U 1 U 4 U 7 U 2 U 5 U 0 U 5 U 2 U 7 U 4 U 1 U 6 U 3 U 0 U 7 U 6 U 5 U 4 U 3 U 2 U 1 U 0 U 0 U 0 U 0 U 4 U 4 U 4 U 4 U 0 U 2 U 4 U 6 U 4 U 6 U 0 U 2 U 0 U 4 U 0 U 4 U 4 U 0 U 4 U 0 U 0 U 6 U 4 U 2 U 4 U 2 U 0 U 6 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 0 U 2 U 4 U 6 U 0 U 2 U 4 U 6 U 0 U 4 U 0 U 4 U 0 U 4 U 0 U 4 U 0 U 6 U 4 U 2 U 0 U 6 U 4 U 2 a 0 a 1 a 2 a 3 a 4 a 5 a 6 a 7 U 0 a 0 U 1 a 1 U 2 a 2 U 3 a 3 U 0 a 4 U 1 a 5 U 2 a 6 U 3 a 7 U 0 a 0 U 1 a 1 U 2 a 2 U 3 a 3 U 0 a 4 U 1 a 5 U 2 a 6 U 3 a 7 U 0 = U 4, U 6 = U 2 1 A 1 U 0 U 0 U 0 U 0 U 0 (a 0 a 4 ) A 3 A 5 = U 0 U 2 U 4 U 6 U 1 (a 1 a 5 ) U 0 U 4 U 0 U 4 U 2 (a 2 a 6 ) A 7 U 0 U 6 U 4 U 2 U 3 (a 3 a 7 ) U A 1 A 3 A 5 A 7 = A 2k+1 = N/2 1 j=0 W 0,0 4 W 0,1 4 W 0,2 4 W 0,3 4 W 1,0 4 W 1,1 4 W 1,2 4 W 1,3 4 W 2,0 4 W 2,1 4 W 2,2 4 W 2,3 4 W 3,0 4 W 3,1 4 W 3,2 4 W 3,3 4 (W N/2 ) kj [(W N ) j (a j a N/2+j )] W 0 8 (a 0 a 4 ) W 1 8 (a 1 a 5 ) W 2 8 (a 2 a 6 ) W 3 8 (a 3 a 7 ) 18

19 8.2 2D FFT: example FFT f(x, y) = I A I exp[ {(x x I ) 2 + (y y I ) 2 }/σ] (127) A I = 255, σ = 5, L x = L y = 512 (128) (Eq.134-Eq.135) σ 1 f(k) = (2π) 2 dre ik r f(r) (129) f(r) = [ A I exp (r R I) 2 ] (130) σ I 1 [ f(k) = (2π) 2 A I dre ik r exp (r R I) 2 ] = 1 [ σ (2π) 2 A I e ik R I dre ik (r RI) exp (r R I) 2 ] σ I I πσ /4 = (2π) 2 e σk2 A I e ik R I I = σ 4π e σk2 /4 S(k) (131) S(k) = I A I e ik R I (structure factor) (132) f(k) = σ 4π e σk2 /4 S(k) (133) Gauss f(k) = 1 2π = A π 2π f(x) = A exp( αx 2 ), α > 0 (134) α e k 2 kx α(x+i dxe 2α) 2 e k2 4α dxae αx2 e ikx = A 2π A 4α = 2 πα e k 2 4α (135) f(x) = A + 2 dke k2 4α e ikx = πα A 2 dke 1 4α (k 2iαx)2 e αx2 = πα A 2 πα 2 απe αx2 = Ae αx2 19

20 For FFT cal. 2 n x 2 m size bmp file is required. Original B min-11.bmp (541x621 image) For 2D FFT calculation B min-11_512x512.bmp (512x512 image)

21 F(k) FFT image max is normalized by 1/500 of the peak value

22 Inverse-FFT image F(k) f ( r ) Inverse FT from the FT image we can get the almost same image as original one. So our FFT calculation may be right.

23 Real space f (r )" Gaussian peak with λ x,y = L" L F(k) k x,y = 2π / λ x,y = 2π / L pixcel =2π / L L L

24 Real space f (r )" Gaussian peak with λ x,y = L / 2" F(k) k x,y = 2π / λ x,y = 4π / L pixcel =2π / L L/2 L

25 Real space f (r )" Gaussian peak with λ x,y = L / 3" F(k) k x,y = 2π / λ x,y = 6π / L pixcel =2π / L

26 Real space f (r )" Gaussian peak with λ x,y = L / 4" F(k) k x,y = 2π / λ x,y = 8π / L pixcel =2π / L

27 Real space f (r )" Gaussian peak with λ x,y = L / 5" F(k) k x,y = 2π / λ x,y = 10π / L pixcel =2π / L

28 Real space f (r )" Gaussian peak with λ x,y = L / 8" F(k) k x,y = 2π / λ x,y = 16π / L pixcel =2π / L

29 Real space f (r )" Gaussian peak with λ x,y = L / 10" F(k) k x,y = 2π / λ x,y = 20π / L pixcel =2π / L

30 Real space f (r )" Gaussian peak with λ x,y = L / 20" F(k) k x,y = 2π / λ x,y = 40π / L pixcel =2π / L

31 Real space f (r )" Gaussian peak with λ x,y = L / 30" F(k) k x,y = 2π / λ x,y = 60π / L pixcel =2π / L

32 The other example original FFT FFT by our code IFFT in our code

Euler Appendix cos, sin 2π t = 0 kx = 0, 2π x = 0 (wavelength)λ kλ = 2π, k = 2π/λ k (wavenumber) x = 0 ωt = 0, 2π t = 0 (period)t T = 2π/ω ω = 2πν (fr

Euler Appendix cos, sin 2π t = 0 kx = 0, 2π x = 0 (wavelength)λ kλ = 2π, k = 2π/λ k (wavenumber) x = 0 ωt = 0, 2π t = 0 (period)t T = 2π/ω ω = 2πν (fr This manuscript is modified on March 26, 2012 3 : 53 pm [1] 1 ( ) Figure 1: (longitudinal wave) (transverse wave). P 7km S 4km P S P S x t x u(x, t) t = t 0 = 0 f(x) f(x) = u(x, 0) v +x (Fig.2) ( ) δt