傅里叶全息变换

 1 # -*- coding:utf-8 -*-
 2 import numpy as np
 3 import matplotlib.pyplot as plt
 4 import cv2
 5 import random,cmath
 6 
 7 '''设置反馈系数、迭代次数'''
 8 alpha=0.1
 9 iteration=int(input('请输入迭代次数'))
10 
11 NRMSD=np.zeros([1,iteration],float)
12 
13 '''读取图像数组，并归一化'''
14 imgGray=cv2.imread('C:\\Users\hanli\Desktop\study\Einstein.tif',0)
15 [r,c]=imgGray.shape
16 Iorg=imgGray/imgGray.max()
17 
18 Iorg=np.fft.ifftshift(Iorg)
19 Iang=2*cmath.pi*random.random()
20 I=Iorg*cmath.exp(1j*Iang)
21 for i in range(iteration):
22     h=np.fft.ifft2(I)
23     hphi=np.angle(h)
24     for m in range(r):
25         for n in range(c):
26             h[m,n]=cmath.exp(1j*hphi[m,n])
27     I=np.fft.fft2(h)
28     Iabs=abs(I)
29     Iabs=(Iabs-Iabs.min())/(Iabs.max()-Iabs.min())
30     Iang=np.angle(I)
31     resid=Iorg-Iabs
32     Inew=Iorg+resid*alpha
33     Iexp=np.ones((r,c),dtype=complex)
34     for m in range(r):
35         for n in range(c):
36             Iexp[m,n]=cmath.exp(1j*Iang[m,n])
37     I=Inew*Iexp
38 
39 I=np.fft.fft2(h)
40 plt.subplot(121)
41 plt.imshow(hphi,cmap='gray')
42 plt.subplot(122)
43 plt.imshow(abs(np.fft.fftshift(I)),cmap='gray')
44 plt.show()

上图是迭代60次的结果

转载于:https://www.cnblogs.com/haliQ/p/6637092.html