一种基于交叉验证的傅里叶变换去噪方法与流程

本发明公开了一种基于交叉验证的傅里叶变换去噪方法，用于滤除光谱中的高频噪声，属于谱处理技术领域。

背景技术：

x射线荧光分析在冶金、考古、探矿等领域有着广泛应用，具有快速、准确、不破坏样本的优点。根据荧光光谱中x射线的荧光强度可以对待测样品进行定性和定量分析，但由于统计涨落、电子学噪声等的影响，荧光光谱含有高频噪声，影响光谱的基线校正、峰的识别等。因此，光谱的去噪技术在光谱分析中具有重要作用。

以往采用的去噪方法多为savitzky-golay滤波器法，对多个连续数据点进行多项式拟合，用拟合后的数据作为去噪后的数据。该方法具有简单、快速的优点，但去噪效果一般。近年来，基于小波变换的去噪方法被广泛应用于光谱去噪领域，由于小波变换具有时频域局部化特性和多分辨率分析特性，其去噪效果较好。但小波变换去噪方法需要确定的参数过多，比如小波函数的选取、分解层数的确定、阈值的计算等，对使用者的经验要求高，不利于光谱数据的自动化处理。傅里叶变换去噪方法是基于噪声频率和有用信号频率的分布差异的去噪方法，具有去噪效果好、使用简单的优点，但去噪阈值的选取往往依靠使用者的经验。

技术实现要素：

本发明针对傅里叶变换去噪方法中阈值不能自动确定问题，提出了一种交叉验证和傅里叶变换相结合确定最优阈值，再利用该最优阈值基于傅里叶变换去噪。该方法能比较准确地确定最优阈值，使傅里叶变换去噪能自动进行，且能有效滤除噪声。

为解决上述技术问题，本发明提供了一种基于交叉验证的傅里叶变换去噪方法，包括如下步骤：

s1、将采样频率为fs、长度为n的含噪信号x(n)按照奇偶序号分为两组，原偶数序号序列记为xe(n)，原奇数序号序列记为xo(n)，每组长度均为n/2，n为偶数；

s2、对奇数序号序列xo(n)进行插值，得到偶数序号序列的估计序列插值公式为：

其中，规定

s3、对估计序列做快速傅里叶变换fft，得到频域信号

s4、用阈值t对频域信号进行阈值处理，滤除频率大于阈值t的高频噪声，得到去噪后的频域信号

式中，表示频域信号第m个点对应的频率；

s5、对频域信号做快速傅里叶逆变换ifft，得到去噪后的时域信号

s6、将xe(n)作为参考值，计算xe(n)和的误差平方和me(t)：

s7、对偶数序号序列xe(n)应用步骤s2)-步骤s6)，得到估计的奇数序号序列的误差平方和mo(t)，则总的误差平方和为

m(t)＝me(t)+mo(t)(6)

s8、计算m(t)取最小值时的阈值tm；

s9、对原始信号x(n)进行快速傅里叶变换，再滤除频域信号中频率大于tm的高频噪声，对去噪后的频域信号进行快速傅里叶逆变换，得到去噪后的光滑信号。

进一步的，偶数序号序列xe(n)和奇数序号序列xo(n)的采样频率均为fs/2，所述阈值t的取值范围为[0,fs/4]。

有益效果：本发明与现有技术相比，本发明具有以下优点：

1)本发明针对有用信号频率和噪声频率不重合的光谱，根据噪声频率的分布特点，提出一种交叉验证和傅里叶变换相结合确定最优阈值，再利用该最优阈值基于傅里叶变换去噪的方法。

2)本发明的核心去噪方法为傅里叶变换去噪法，当高频的噪声频率和低频的有用信号频率分离时，可以有效滤出高频噪声，具有快速、准确、简单的特点。

3)本发明使用交叉验证和傅里叶变换相结合确定最优阈值，确定的最优阈值和人工确定的最优阈值相接近，表明该阈值去噪效果较佳，且使得傅里叶变换去噪可以自动化进行，提高光谱处理的效率和便捷性。

附图说明

图1为交叉验证和和傅里叶变换相结合确定最优阈值的流程图；

图2为信噪比为20db的仿真光谱及其频谱图；

图3为信噪比为20db的仿真光谱的奇数序号序列和偶数序号序列的时域图及其频谱图，其中：(a)奇数序号序列的时域图；(b)奇数序号序列的频谱图；(c)偶数序号序列的时域图；(d)偶数序号序列的频谱图；

图4为不同阈值下估计的奇数序号序列和偶数序号序列的总误差平方和m(t)曲线；

图5为不同阈值下原含噪信号的傅里叶变换去噪后光滑信号的误差平方和曲线；

图6为基于交叉验证的傅里叶变换去噪法和savitzky-golay滤波器法对不同信噪比的仿真光谱去噪后的信噪比曲线；

图7为基于交叉验证的傅里叶变换去噪法和savitzky-golay滤波器法对信噪比为20db的仿真光谱的去噪效果对比图；

图8为对实测能量色散x射线荧光光谱去噪效果图。

具体实施方式

下面以能量色散x射线荧光光谱为例，详细说明本发明的实施例，所述实施例的示例在附图中示出。

图1为本发明具体实施时的流程图，图2(a)为信噪比为20db的仿真光谱。具体步骤如下：

s1、将信噪比为20db、长度为1024的仿真能量色散x射线荧光光谱x(n)按照奇偶序号分为两组，原偶数序号序列记为xe(n)，原奇数序号序列记为xo(n)，每组长度均为512。含噪信号x(n)如图2(a)所示，奇数序号序列xo(n)如图3(a)所示，偶数序号序列xe(n)如图3(c)所示。假设含噪信号x(n)的采样频率fs＝1，则奇数序号序列xo(n)和偶数序号序列xe(n)的采样频率为fs/2＝0.5，对x(n)、xo(n)和xe(n)做快速傅里叶变换fft，得到它们相应的频谱图。含噪信号x(n)的频谱图如图2(b)所示，奇数序号序列xo(n)的频谱图如图3(b)所示，偶数序号序列xe(n)的频谱图如图3(d)所示，频谱图的横坐标为频率，纵坐标为频率的幅度。

s2、对奇数序号序列xo(n)进行插值，得到偶数序号序列的估计序列插值公式为

其中，规定xo(513)＝xo(1)。

s3、估计序列做快速傅里叶变换fft，得到频域信号

s4、用阈值t对频域信号进行阈值处理，滤除频率大于t的高频噪声，得到去噪后的频域信号

式中，表示频域信号第m个点对应的频率。根据奈奎斯特采样定理，信号的最大频率应小于等于采样频率的一半，所以阈值t的取值范围为[0,0.25]。以0.001为间隔获得阈值序列，用不同的阈值处理频域信号，得到不同的去噪后的频域信号。对于不同应用，阈值序列的间隔可以增加或缩小，如果对去噪效果要求高，可以适当减小间隔，如果对算法的执行时间要求高，可以适当增加间隔。

s5、对频域信号做快速傅里叶逆变换ifft，得到去噪后的时域信号

s6、将xe(n)作为参考值，计算xe(n)和的误差平方和me(t)，

s7、对偶数序号序列xe(n)应用步骤s2)-步骤s6)，得到估计的奇数序号序列的误差平方和mo(t)，则总的误差平方和为

m(t)＝me(t)+mo(t)(6)

图4为不同阈值下估计的奇数序号序列和偶数序号序列的总误差平方和m(t)曲线；

s8、计算m(t)取最小值时的阈值tm，对于本示例的仿真信号，tm＝0.063。为了判断交叉验证和傅里叶变换相结合确定的最优阈值是否合适，直接用不同阈值对含噪信号x(n)进行傅里叶变换去噪，以无噪信号为参考值计算光滑信号的误差平方和，误差平方和曲线如图5所示。由图5得，人工确定的最优阈值为0.060，两者相差不大，所以交叉验证和傅里叶变换相结合确定的最优阈值是合适的。

s9、对原始信号x(n)进行快速傅里叶变换，再滤除频域信号中频率大于tm的高频噪声，对去噪后的频域信号进行快速傅里叶逆变换，得到去噪后的光滑信号，其信噪比为28.78。而基于人工确定的最优阈值的傅里叶变换去噪，得到去噪后的光滑信号的信噪比为28.93。

为了证明本发明提出的方法的优越性和适应性，用本方法和savitzky-golay滤波器法对不同信噪比的仿真信号进行去噪处理，计算光滑信号的信噪比。图6为基于交叉验证的傅里叶变换去噪法和savitzky-golay滤波器法对不同信噪比的仿真光谱去噪后的信噪比曲线，横坐标为原始含噪信号的信噪比，纵坐标为去噪后光滑信号的信噪比。由图6得，不同信噪比下，基于交叉验证的傅里叶变换去噪法比savitzky-golay滤波器法对仿真光谱的去噪效果更好。为了更直观的对比两种方法的去噪效果，在图7中显示两种方法对信噪比为20db的仿真光谱的去噪效果，且显示时相对原始含噪信号均向下偏移。为了证明本方法的实用性，对一个实测的能量色散x射线荧光光谱进行去噪处理，结果如图8所示。由图8得，本方法去噪后的信号更加光滑，且无明显失真，表明本方法可以用于实际光谱的去噪处理。