单频信号的相位谱计算与误差修正-附Matlab代码一、问题描述二、原因分析三、解决方法四、Matlab程序获取与验证

71 阅读 0 评论 47 点赞

我是靠谱客的博主不安大炮，最近开发中收集的这篇文章主要介绍单频信号的相位谱计算与误差修正-附Matlab代码一、问题描述二、原因分析三、解决方法四、Matlab程序获取与验证，觉得挺不错的，现在分享给大家，希望可以做个参考。

概述

一、问题描述

我们在实际处理时经常遇到只有一个正弦信号的情况，其频率为 ${{f}_{0}}$ ，在谱分析以后，除了在频率为 ${{f}_{0}}$ 处有相位数值外，其他频率处都有相位数值，分析其他频谱出现相位值的原因。

例如，假设信号采样率为1000Hz，有一个余弦信号，其频率为f0=50Hz，幅值都为1，初始相角为pi/3=1.0472，信号长度为1000，该信号如下图所示：

对应Matlab代码如下：

fs=1000;                       % 采样频率
f0=50;                         % 信号频率
A=1;                           % 信号幅值
theta0=pi/3;                   % 信号初始相角
N=1000;                        % 信号长度
t=(0:N-1)/fs;                  % 设置时间序列
x=A*cos(2*pi*f0*t+theta0);     % 设置信号
figure,plot(x,'r')
xlabel('样点'); ylabel('幅值'); title('原始信号')

在这里插入图片描述

图1 原始信号

此时，对上述信号进行傅里叶变换，并获取其频谱图，如下图所示。从幅值谱图上可知，明显可看到在50Hz处有一个峰值，其他频率幅值都为0。

对应的Matlab代码如下所示：

X=fft(x);                      % FFT
n2=1:N/2+1;                    % 设置索引号序列
freq=(n2-1)*fs/N;              % 设置频率刻度
figure,plot(freq,abs(X(n2))*2/N,'r');title('幅值谱图-正频率');xlabel('频率/Hz'); ylabel('幅值')

在这里插入图片描述

图2 频谱图

绘制该信号的相位谱，如下图所示，可知该相位谱非常乱，除了50Hz频率外的其他频率都有一个初始相位。

THETA=angle(X(n2));            % 计算初始相角
Am=abs(X(n2));                 % 计算幅值
ph0=THETA(51);                 % 计算信号的初始相角
figure,plot(freq,THETA,'r')
xlabel('频率/Hz'); ylabel('初始角/弧度');title('相位谱图-正频率部分')

在这里插入图片描述

图3 相位谱图

二、原因分析

首先计算50Hz和其他频谱线的实部和虚部，以及各自的初始相位值。如下图所示：

在这里插入图片描述
可以发现，在50Hz频率出，计算的初相为1.0472与设置的初相角一致，然而，在频率为非50Hz以外，如1Hz等，它们的实部和虚部的值都是 ${{10}^{-14}}sim {{10}^{-13}}$ 量级，这是由计算误差生的，使用反三角函数后就能得到相应的初始相位数值，分布在-pi和pi之间，造成了相位值的混乱。而其他频率分量也差不多是这样的量级。

三、解决方法

要消除这种相角初始值的混乱，可在程序中设置一个阈值Th。在本例中有用的频率分量只有50Hz一条谱线，阈值比较容易选择；若频谱中有用的频率分量有多条谱线，则可以寻找有用的频率分量中的最小幅值Amin，而要把阈值Th设置成Th<Amin。

Matlab如下：

Th=0.1;                        % 设置阈值
thetadex=find(Am<Th);          % 寻找小于阈值的那线谱线的索引
THETA1=THETA;                  % 初始化THETA1
THETA1(thetadex)=0;            % 对于小于阈值的那线谱线初始相位都为0
figure,plot(freq,THETA1,'r')
xlabel('频率/Hz'); ylabel('初始角/弧度');title('相位谱图-正频率部分')