Numpy和Tensorflow傅立叶变换示例:幅值和相位

2019-04-14 20:31发布

一、numpy傅立叶变换和反变换

import cv2 import numpy as np import matplotlib.pyplot as plt img = cv2.imread('Set5/001_HR.png') (r, g, b)=cv2.split(img) img=cv2.merge([b,g,r]) #cv2.imshow("test",img) #cv2.waitKey(0) #plt.imshow(img) #img = cv2.imread('flower.jpg',0) #直接读为灰度图像 f = np.fft.fft2(img) fshift = np.fft.fftshift(f) #取绝对值:将复数变化成实数 #取对数的目的为了将数据变化到0-255 s1 = np.log(np.abs(fshift)) plt.subplot(221),plt.imshow(img),plt.title('original') plt.xticks([]),plt.yticks([]) #--------------------------------------------- # 逆变换--取绝对值就是振幅 f1shift = np.fft.ifftshift(np.abs(fshift)) img_back = np.fft.ifft2(f1shift) #出来的是复数,无法显示 img_back = np.abs(img_back) #调整大小范围便于显示 img_back = (img_back-np.amin(img_back))/(np.amax(img_back)-np.amin(img_back)) plt.subplot(222),plt.imshow(img_back,'gray'),plt.title('only Amplitude') plt.xticks([]),plt.yticks([]) #--------------------------------------------- # 逆变换--取相位 f2shift = np.fft.ifftshift(np.angle(fshift)) img_back = np.fft.ifft2(f2shift) #出来的是复数,无法显示 img_back = np.abs(img_back) #调整大小范围便于显示 img_back = (img_back-np.amin(img_back))/(np.amax(img_back)-np.amin(img_back)) plt.subplot(223),plt.imshow(img_back,'gray'),plt.title('only phase') plt.xticks([]),plt.yticks([]) #--------------------------------------------- # 逆变换--将两者合成看看 s1 = np.abs(fshift) #取振幅 s1_angle = np.angle(fshift) #取相位 s1_real = s1*np.cos(s1_angle) #取实部 s1_imag = s1*np.sin(s1_angle) #取虚部 s2 = np.zeros(img.shape,dtype=complex) s2.real = np.array(s1_real) #重新赋值给s2 s2.imag = np.array(s1_imag) f2shift = np.fft.ifftshift(s2) #对新的进行逆变换 img_back = np.fft.ifft2(f2shift) #出来的是复数,无法显示 img_back = np.abs(img_back) #调整大小范围便于显示 img_back = (img_back-np.amin(img_back))/(np.amax(img_back)-np.amin(img_back)) plt.subplot(224),plt.imshow(img_back,'gray'),plt.title('another way') plt.xticks([]),plt.yticks([]) plt.show()

二、tensorflow傅立叶变换示例:幅值和相位

#coding=utf-8 import cv2 import numpy as np import matplotlib.pyplot as plt import tensorflow as tf from skimage import measure from skimage import img_as_ubyte # plt显示fft的幅值图 """ img = cv2.imread('Set5/001_HR.png',0) f = np.fft.fft2(img) fshift = np.fft.fftshift(f) s1 = np.log(np.abs(f)) s2 = np.log(np.abs(fshift)) plt.subplot(121),plt.imshow(s1,'gray'),plt.title('original') #plt.subplot(122),plt.imshow(s2,'gray'),plt.title('center') sess = tf.InteractiveSession() xs = tf.placeholder(dtype=tf.complex64,shape=(img.shape[0],img.shape[1])) f1=tf.fft2d(xs) f1_dis=f1.eval(feed_dict={xs:img}) f1shift = np.fft.fftshift(f1_dis) ss1 = np.log(np.abs(f1_dis)) ss1 = np.float64(ss1) plt.subplot(122),plt.imshow(ss1,'gray'),plt.title('fft_original') plt.show() """ # plt显示fft的相位图 """ Ihr = cv2.imread('Set5/003_HR.png',1) Ilr = cv2.imread('Set5/003_EnhanceNet-PAT.png',1) ph_out=[] for i,img in enumerate([Ihr,Ilr]): f = np.fft.fft2(img) fshift = np.fft.fftshift(f) #取绝对值:将复数变化成实数 #取对数的目的为了将数据变化到较小的范围(比如0-255) ph_f = np.angle(f) ph_out.append(ph_f) ph_fshift = np.angle(fshift) index='22'+ str(i+1) plt.subplot(int(index)),plt.imshow(ph_f,'gray'),plt.title('original') #index = ['22' + str(i * 2 + 1), '22' + str(i * 2 + 2)] #plt.subplot(int(index[0])),plt.imshow(ph_f,'gray'),plt.title('original') #plt.subplot(int(index[1])),plt.imshow(ph_fshift,'gray'),plt.title('center') shape = Ihr.shape sess = tf.InteractiveSession() xs = tf.placeholder(dtype=tf.complex64, shape=([2, shape[1], shape[0], shape[2]])) fhr = tf.fft2d(xs[0]) flr = tf.fft2d(xs[1]) #fshift_hr, fshift_lr = np.fft.fftshift([fhr_dis,flr_dis]) ph_fhr = tf.cast(tf.angle(fhr),dtype=tf.float64) ph_flr = tf.cast(tf.angle(flr),dtype=tf.float64) #ph_fhr = tf.atan2(tf.imag(fhr),tf.real(fhr)) #ph_flr = tf.atan2(tf.imag(flr),tf.real(flr)) ph_fhr, ph_flr = sess.run([ph_fhr,ph_flr],feed_dict={xs:[Ihr, Ilr]}) print ("MSE of ph_fhr and ph_out[0]:",measure.compare_mse(ph_fhr,ph_out[0])) print ("MSE of ph_fhr and ph_flr:",measure.compare_mse(ph_fhr,ph_flr)) plt.show() """ # plt显示fft的赋值大小 """ Ihr = cv2.imread('Set5/003_HR.png',1) Ilr = cv2.imread('Set5/003_EnhanceNet-PAT.png',1) nhr = np.fft.fft2(Ihr) nlr = np.fft.fft2(Ilr) amp_nhr = np.log(np.abs(nhr)) amp_nlr = np.log(np.abs(nlr)) shape = Ihr.shape sess = tf.InteractiveSession() xs = tf.placeholder(dtype=tf.complex64, shape=([2, shape[1], shape[0], shape[2]])) fhr = tf.fft2d(xs[0]) flr = tf.fft2d(xs[1]) amp_fhr = tf.cast(tf.log(tf.abs(fhr)),dtype=tf.float64) amp_flr = tf.cast(tf.log(tf.abs(flr)),dtype=tf.float64) amp_fhr, amp_flr = sess.run([amp_fhr,amp_flr],feed_dict={xs:[Ihr, Ilr]}) print ("MSE of ph_fhr and ph_out[0]:",measure.compare_mse(amp_fhr,amp_flr)) """ sess = tf.InteractiveSession() theta0 = tf.atan2(1.0, 1.0) # pi/4 x0 = tf.constant(1.0) y = tf.Variable(0.4) c = tf.complex(x0, y) c = tf.expand_dims(c,0) c = tf.fft(c) theta = tf.angle(c) err = (theta - theta0)**2 optimize = tf.train.GradientDescentOptimizer(1e-1).minimize(err) sess.run(tf.global_variables_initializer()) for i in range(10000): sess.run(optimize) if i % 10 == 0: print(sess.run(c)) # should go to 1 + 1j a=tf.constant(1e-8,tf.complex64,tf.shape())

三、numpy相位交换,幅值不变

import cv2 import numpy as np import matplotlib.pyplot as plt from skimage import measure from skimage import img_as_ubyte #img_flower = cv2.imread('Set5/001_HR.png',1) #直接读为灰度图像 #img_man = cv2.imread('Set5/001_EnhanceNet-E.png',1) #直接读为灰度图像 #img_man = cv2.imread('Set5/004_EnhanceNet-PAT.png',1) #直接读为灰度图像 img_flower = cv2.imread('../../datasets/SuperResolution/EnhanceNet_Urban100/Urban100/039_HR.png',1) #直接读为灰度图像 img_man = cv2.imread('../../datasets/SuperResolution/EnhanceNet_Urban100/Urban100/039_EnhanceNet-E.png',1) #直接读为灰度图像 #img_man = cv2.imread('../../datasets/SuperResolution/EnhanceNet_BSD100BSD100/043_EnhanceNet-PAT.png',1) #cv2.imshow("test",img_flower) #cv2.waitKey(0) print ("PSNR of 0 and 1:",measure.compare_psnr(img_flower,img_man)) (r, g, b)=cv2.split(img_flower) img_flower=cv2.merge([b,g,r]) (r, g, b)=cv2.split(img_man) img_man=cv2.merge([b,g,r]) print ("PSNR of 0 and 1:",measure.compare_mse(img_flower,img_man)) plt.subplot(221),plt.imshow(img_flower,'gray'),plt.title('001_HR') plt.xticks([]),plt.yticks([]) plt.subplot(222),plt.imshow(img_man,'gray'),plt.title('001_EnhanceNet-PAT') plt.xticks([]),plt.yticks([]) #-------------------------------- f1 = np.fft.fft2(img_flower) f1shift = np.fft.fftshift(f1) f1_A = np.abs(f1shift) #取振幅 f1_P = np.angle(f1shift) #取相位 #-------------------------------- f2 = np.fft.fft2(img_man) f2shift = np.fft.fftshift(f2) f2_A = np.abs(f2shift) #取振幅 f2_P = np.angle(f2shift) #取相位 #---图1的振幅--图2的相位-------------------- img_new1_f = np.zeros(img_flower.shape,dtype=complex) img1_real = f1_A*np.cos(f2_P) #取实部 img1_imag = f1_A*np.sin(f2_P) #取虚部 img_new1_f.real = np.array(img1_real) img_new1_f.imag = np.array(img1_imag) f3shift = np.fft.ifftshift(img_new1_f) #对新的进行逆变换 img_new1 = np.fft.ifft2(f3shift) #出来的是复数,无法显示 img_new1 = np.abs(img_new1) #调整大小范围便于显示 img_new1 = (img_new1-np.amin(img_new1))/(np.amax(img_new1)-np.amin(img_new1)) plt.subplot(223),plt.imshow(img_new1,'gray'),plt.title('another way') plt.xticks([]),plt.yticks([]) #---图2的振幅--图1的相位-------------------- img_new2_f = np.zeros(img_flower.shape,dtype=complex) img2_real = f2_A*np.cos(f1_P) #取实部 img2_imag = f2_A*np.sin(f1_P) #取虚部 img_new2_f.real = np.array(img2_real) img_new2_f.imag = np.array(img2_imag) f4shift = np.fft.ifftshift(img_new2_f) #对新的进行逆变换 img_new2 = np.fft.ifft2(f4shift) #出来的是复数,无法显示 img_new2 = np.abs(img_new2) #调整大小范围便于显示 img_new2 = (img_new2-np.amin(img_new2))/(np.amax(img_new2)-np.amin(img_new2)) plt.subplot(224),plt.imshow(img_as_ubyte(img_new2),'gray'),plt.title('another way') plt.xticks([]),plt.yticks([]) plt._imsave("test.jpg",img_new2) print ("PSNR of 0 and 4:",measure.compare_mse(img_flower,img_as_ubyte(img_new2))) plt.show()

四、几个简单的画图函数:

#!/usr/bin/env python # -*- coding: utf-8 -*- import numpy as np from PIL import Image import matplotlib import matplotlib.path as path import matplotlib.pyplot as plt import matplotlib.patches as patches from matplotlib.ticker import MultipleLocator fig = plt.figure(figsize=(8,6), dpi=72,facecolor="white") axes = plt.subplot(211, aspect=1, axisbelow=True) lena = np.array(Image.open("output/001_HR.png")) im = plt.imshow(lena, origin = 'upper', interpolation='nearest') axes.set_xticks([]) axes.set_yticks([]) axes = plt.axes([0.05,0.00,0.40,0.50], frameon=False) patch1 = patches.Circle((.5,.5), radius=.45, transform=axes.transAxes, linewidth=2, facecolor='none', edgecolor='black',zorder=2) patch2 = patches.Circle((.5,.5), radius=.46, transform=axes.transAxes, linewidth=2, facecolor='white', edgecolor='white',zorder=1) im = plt.imshow(lena[215:315,285:385], origin = 'upper', zorder=3, interpolation='nearest', clip_path=patch1) axes.set_xticks([]) axes.set_yticks([]) axes.add_patch(patch1) axes.add_patch(patch2) plt.show() #coding=utf-8 import numpy as np import matplotlib.pyplot as plt img001_HR = plt.imread('output/001_HR.png',0) img001_En = plt.imread('output/001_EnhanceNet-PAT.png',0) img002_HR = plt.imread('output/002_HR.png',0) img002_En = plt.imread('output/002_EnhanceNet-PAT.png',0) img003_HR = plt.imread('output/003_HR.png',0) img003_En = plt.imread('output/003_EnhanceNet-PAT.png',0) img004_HR = plt.imread('output/004_HR.png',0) img004_En = plt.imread('output/004_EnhanceNet-PAT.png',0) img005_HR = plt.imread('output/001_HR.png',0) img005_En = plt.imread('output/001_EnhanceNet-PAT.png',0) plt.subplot(521),plt.imshow(img001_HR,'gray'),plt.title('img001_HR') plt.subplot(522),plt.imshow(img001_En,'gray'),plt.title('img001_En') plt.subplot(523),plt.imshow(img002_HR,'gray'),plt.title('img002_HR') plt.subplot(524),plt.imshow(img002_En,'gray'),plt.title('img002_En') plt.subplot(525),plt.imshow(img003_HR,'gray'),plt.title('img003_HR') plt.subplot(526),plt.imshow(img003_En,'gray'),plt.title('img003_En') plt.subplot(527),plt.imshow(img004_HR,'gray'),plt.title('img004_HR') plt.subplot(528),plt.imshow(img004_En,'gray'),plt.title('img004_En') """ plt.subplot(529),plt.imshow(img005_HR,'gray'),plt.title('img005_HR') #plt.subplot(5210),plt.imshow(img005_En,'gray'),plt.title('img005_En') """ plt.show()