1. 加性高斯白噪声
1.1 特点
加性高斯白噪声属于白噪声的一种,有如下两个特点:
- 功率谱密度为常数,即白
- 高斯则为其幅度的概率密度函数服从高斯分布
1.2 高斯分布的一维概率密度函数
- 为标准差
- 为均值
1.3 Python内置库说明
random.gauss(mu, sigma)其值即服从高斯分布,若想要是实现加性高斯白噪声,循环作加即可
1.4 Code
import numpy as np
import random
def gaussian_white_noise(intput_signal, mu, sigma):
'''
加性高斯白噪声(适用于灰度图)
:param intput_signal: 输入图像
:param mu: 均值
:param sigma: 标准差
:return:
'''
intput_signal_cp = np.copy(intput_signal) # 输入图像的副本
m, n = intput_signal_cp.shape # 输入图像尺寸(行、列)
# 添加高斯白噪声
for i in range(m):
for j in range(n):
intput_signal_cp[i, j] = intput_signal_cp[i, j] + random.gauss(mu, sigma)
return intput_signal_cp
1.5 测试对比
2. 逆滤波
实际上逆滤波是维纳滤波的一种理想情况,当不存在加性噪声时,维纳滤波与逆滤波等同。
2.1 基本原理
在时域内有
- 为输入信号
- 为系统函数
- 为输出信号
根据时域卷积定理,我们知道时域卷积等于频域乘积
- 分别为输出信号,输入信号,系统函数的傅里叶变换
则有
这意味着,当我们已知系统函数时,我们可以很简单的完成滤波。
2.2 Code
import numpy as np
def inverse_filtering(input_singal, h):
'''
逆滤波
:param input_singal: 输入信号
:param h: 退化函数(时域)
:return: 滤波后的信号(幅值)
'''
output_signal = [] #输出信号
output_signal_fft = [] #输出信号的傅里叶变换
h_fft = np.fft.fft2(h) # 退化函数的傅里叶变换
input_singal_fft = np.fft.fft2(input_singal) # 输入信号的傅里叶变换
output_signal_fft = input_singal_fft / h_fft
output_signal = np.abs(np.fft.ifft2(output_signal_fft))
return output_signal
3 维纳滤波
3.1 直观认识
理解了逆滤波的基本过程之后,实际上维纳滤波就不是太大问题了。实际上,逆滤波对于绝大多数情况滤波效果都不好,因为逆滤波是通过傅里叶变换将信号由时域转换到频域,再根据时域卷积定理,在频域作除法。对于乘性干扰这当然是没问题的,甚至是完美的。而如果存在加性噪声,例如:加性高斯白噪声。逆滤波效果就不好了,某些情况下几乎无法完成滤波情况。
3.2 数学解释
输入信号经过系统函数后
时域上
频域上
若存在加性噪声则为
时域上
- 其中
频域上
- 其中
则
于是,从上面对输入信号的估计表达式可以看出,多出了一项加性噪声的傅里叶变换与系统函数的比值。尤其当相对于很小时,滤波后的信号差距十分严重。
3.3 维纳滤波表达式
而我们又知道:白噪声的白为噪声的功率谱为常数,即为常数,于是,从直观上看,当相对于较大时,则较小,上式第一项则较小,而第二项较大从而保持相对平稳。
3.4 Code
import numpy as np
def wiener_filtering(input_signal, h, K):
'''
维纳滤波
:param input_signal: 输入信号
:param h: 退化函数(时域)
:param K: 参数K
:return: 维纳滤波后的信号(幅值)
'''
output_signal = [] # 输出信号
output_signal_fft = [] # 输出信号的傅里叶变换
input_signal_cp = np.copy(input_signal) # 输入信号的副本
input_signal_cp_fft = np.fft.fft2(input_signal_cp) # 输入信号的傅里叶变换
h_fft = np.fft.fft2(h) # 退化函数的傅里叶变换
h_abs_square = np.abs(h_fft)**2 # 退化函数模值的平方
# 维纳滤波
output_signal_fft = np.conj(h_fft) / (h_abs_square + K)
output_signal = np.abs(np.fft.ifft2(output_signal_fft * input_signal_cp_fft)) # 输出信号傅里叶反变换
return output_signal