频域滤波

1、基于傅里叶变换的同态滤波（高斯高通滤波器）

'''
基于傅里叶变换的同态滤波
高斯高通滤波器
参考：https://blog.csdn.net/m0_37661841/article/details/103716165

'''
import os
import cv2
from PIL import Image
import numpy as np


def homomorphic_filter(src, d0=10, r1=0.5, rh=2, c=0.5, h=2.0, l=0.5):#原程序c=4
    gray = src.copy()
    if len(src.shape) > 2:
        gray = cv2.cvtColor(src, cv2.COLOR_BGR2GRAY)
    gray = np.float64(gray) 
    rows, cols = gray.shape
    gray_fft = np.fft.fft2(gray)# 傅里叶变换
    gray_fftshift = np.fft.fftshift(gray_fft)# 中心化
    #print(len(gray_fftshift[0]))
    dst_fftshift = np.zeros_like(gray_fftshift)# 创建掩膜模板
    M, N = np.meshgrid(np.arange(-cols // 2, cols // 2), np.arange(-rows//2, rows//2))
    D = np.sqrt(M ** 2 + N ** 2)#到中心点距离
    Z = (rh - r1) * (1 - np.exp(-c * (D ** 2 / d0 ** 2))) + r1#定义高斯滤波器
    dst_fftshift = Z * gray_fftshift#滤波
    dst_fftshift = (h - l) * dst_fftshift + l#对滤波后的图线性变换
    dst_ifftshift = np.fft.ifftshift(dst_fftshift)#去中心化
    dst_ifft = np.fft.ifft2(dst_ifftshift)#傅里叶反变换
    dst = np.real(dst_ifft)#获取实部
    dst = np.uint8(np.clip(dst, 0, 255))
    return dst



if __name__ == '__main__':
    img = cv2.imread('./DIP3E_Original_Images_CH04/Fig0462(a)(PET_image).tif',0)
    img_new = homomorphic_filter(img)
    print("new img shape is {}".format(img_new.shape))
    cv2.imwrite('original.jpg',img)
    cv2.imwrite("3.png", img_new)

2、基于小波变换的同态滤波

'''
基于小波变换的同态滤波
读取本地图片：https://stackoom.com/question/3qAMe/%E5%A6%82%E4%BD%95%E4%BD%BF%E7%94%A8PyWavelets%E8%AF%BB%E5%8F%96%E5%9B%BE%E5%83%8F
源码：https://pywavelets.readthedocs.io/en/latest/
小波变换函数选择：https://blog.csdn.net/jbb0523/article/details/42586749
小波分析库中的一些基本使用方法:https://blog.csdn.net/qq_40587575/article/details/83154042?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522160594795319195264716818%2522%252C%2522scm%2522%253A%252220140713.130102334..%2522%257D&request_id=160594795319195264716818&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2~all~baidu_landing_v2~default-1-83154042.pc_search_result_no_baidu_js&utm_term=pywt.families&spm=1018.2118.3001.4449
小波变换库常用：https://blog.csdn.net/nanbei2463776506/article/details/64124841?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522160594711119725222412888%2522%252C%2522scm%2522%253A%252220140713.130102334..%2522%257D&request_id=160594711119725222412888&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2~all~sobaiduend~default-1-64124841.pc_search_result_no_baidu_js&utm_term=pywt.idwt2&spm=1018.2118.3001.4449
对分量的操作：https://www.cnblogs.com/xixixing/p/5834511.html
'''
import numpy as np
import matplotlib.pyplot as plt
import cv2
import pywt
import pywt.data

#print (pywt.families())
bior = pywt.Wavelet('bior1.3')   #创建一个小波对象
print(bior)
# Load image
original = pywt.data.camera()# 加载小波变换包中的示例图像
original = cv2.imread('E:\Test_code\light\original.jpg',0)
cv2.imshow('a0',original)
original = np.float32(original)
original /= 255
# Wavelet transform of image, and plot approximation and details
titles = ['Approximation', ' Horizontal detail',
          'Vertical detail', 'Diagonal detail']
coeffs2 = pywt.dwt2(original, 'bior1.3')#离散小波变换,多尺度变换pywt.wavedec2(data, wavelet, mode=‘symmetric’, level=多少级, axes=(-2, -1))
LL, (LH, HL, HH) = coeffs2
coeffs2 = 1*LL,(2*LH,2*HL,2*HH)
fig = plt.figure(figsize=(12, 3))
for i, a in enumerate([LL, LH, HL, HH]):
    ax = fig.add_subplot(1, 4, i + 1)
    ax.imshow(a, interpolation="nearest", cmap=plt.cm.gray)
    ax.set_title(titles[i], fontsize=10)
    ax.set_xticks([])
    ax.set_yticks([])
result = pywt.idwt2(coeffs2,'bior1.3')
result *= 255
result = np.uint8(result)
a = np.concatenate([result,original,],axis=1)
cv2.imshow('a',result)
cv2.imwrite('2.png',result)
cv2.waitKey(0)

fig.tight_layout()
plt.show()

以上两种方法均没有体现同态的思想，知识变换到频域和时频域进行高通、高/低同滤波。因为同态滤波是基于照明-反射模型的，明显的特征是取对数，但是上述代码并未体现。