Python 2 D離散フーリエ変換

Python 2 D離散フーリエ変換
文書ディレクトリ

Python 2 D離散フーリエ変換

必要なライブラリ

は2枚のピクチャのPSNR

を計算する

二次元離散フーリエ変換

二次元離散フーリエ逆変換

周波数領域平行移動

周波数領域画像

を描画する
必要なライブラリ

import numpy as np
import cv2
import matplotlib.pyplot as plt

2枚の画像のPSNRを計算します

def PSNR(A, B):
    MSE = np.sum((A - B) ** 2) / A.shape[0] / A.shape[1]
    return 10 * np.log10((255 ** 2) / MSE)

にじげんりさんフーリエへんかん

def DFT_2(f, N, M):
    '''
    :param f: 2_dim marix
    :param N: even; denote number of rows
    :param M: even; denote number of cols
    :return: DFT results without shifting(complex number)
    '''
    if N % 2 != 0 or M % 2 != 0:
        print("param is wrong")
        return
    f = f.astype('float64')
    rows = f.shape[0]
    cols = f.shape[1]
    n = np.arange(0, N, 1).reshape((N, 1))
    row = np.arange(0, rows, 1).reshape((1, rows))
    left = np.exp(-1j*2*np.pi/N*(n @ row))
    col = np.arange(0, cols, 1).reshape((cols, 1))
    m = np.arange(0, M, 1).reshape((1, M))
    right = np.exp(-1j*2*np.pi/M*(col @ m))
    F = left @ f @ right
    return F

2 D離散フーリエ逆変換

def IDFT_2(F, rows, cols, option=True):
    '''
    :param F: 2-dim freqency spectrum without shifting
    :param rows: rows of original picture
    :param cols: cols of original picture
    :param option: whether make some process for results
    :return: picture in space domain
    '''
    if cols % 2 != 0 or rows % 2 != 0:
        print("IDFT params are wrong")
        return
    N = F.shape[0]
    M = F.shape[1]
    row = np.arange(0, rows, 1).reshape((rows, 1))
    n = np.arange(0, N, 1).reshape((1, N))
    left = np.exp(1j*2*np.pi/N*(row @ n))
    m = np.arange(0, M, 1).reshape((M, 1))
    col = np.arange(0, cols, 1).reshape((1, cols))
    right = np.exp(1j*2*np.pi/M*(m @ col))
    f = (left @ F @ right) / rows / cols
    if option:
        f = np.abs(f)
        f = np.around((f - np.min(f)) / (np.max(f) - np.min(f)) * 255)
    return f

ここで、optionがTrueの場合、戻り値が複素のモードであり、そうでない場合、戻り値が複素である.
しゅうはすうりょういきへいこう
DFTの結果は共役対称性を持つので,結果のこの対称性を見出すには結果を平行移動した後でなければならない.

def DFT_shifting(F):
    F_shift = 1j + np.zeros(F.shape)
    N = F.shape[0]
    M = F.shape[1]
    F_shift[N>>1:N, M>>1:M] = F[0:N>>1, 0:M>>1]
    F_shift[N>>1:N, 0:M>>1] = F[0:N>>1, M>>1:M]
    F_shift[0:N>>1, M>>1:M] = F[N>>1:N, 0:M>>1]
    F_shift[0:N>>1, 0:M>>1] = F[N>>1:N, M>>1:M]
    return F_shift

ここで、Fは上記DFT_を示す2()求めた周波数領域の結果.
周波数領域画像の描画
DFTで得られた結果は複素数であるため,そのモデリング値を求めるとともにlogをとる必要がある.

def DFT_2_pic(f, N, M):
    '''
    :param f: 2_dim marix
    :param N: even; denote number of rows
    :param M: even; denote number of cols
    :return: DFT results with shifting(complex number) and make log
    '''
    if N % 2 != 0 or M % 2 != 0:
        print("param is wrong")
        return
    F = DFT_2(f, N, M)
    shift = DFT_shifting(F)
    F_pic = np.log2(np.abs(shift))
    F_pic = np.round((F_pic-np.min(F_pic)) / (np.max(F_pic) - np.min(F_pic)) * 255)   # make pixel value in 0 - 255
    F_pic = F_pic.astype('uint8')
    return shift, F_pic

ここでF_picは対数を取った後のスペクトルを返す.shiftは,平行移動後の複素スペクトルを表す.