基于循环谱的水声通信信号特征分析和提取
- 归纳循环谱的几种方法，总结OFDM和BPSK是如何通过循环谱提取特征参数
- 实践仿真利用循环谱提取BPSK信号的特征参数
- 下周任务：利用循环谱提取OFDM信号的特征参数
目录
1.研究背景和发展现状
2.原理过程
补充上周所用的原理，这里总结了时域和频域两种角度的方法：
其中，时域平滑循环周期可以转换为FAM算法（快速FFT变换累积算法），FAM算法有两个不同点，一是FAM算法是在时域平滑循环周期法的基础上进行令两次FFT累计计算的，第二点石在进行离散傅里叶变换的时候，考虑令加权窗对循环谱密度估计的改进作用，从而使估计效果更准确。
频域平滑循环周期可以转黄为DFSM算法（直接频率平滑算法），DFSM算法是在频域平滑循环周期法的基础上直接计算谱分量之间的谱相关。算法框图与上图一致。
3.谱相关分析的特点
总结了几篇论文和书籍，谱相关分析方法主要有以下特点：
1.很多声呐信号都是循环平稳信号，而噪声却是平稳信号。当给回波信号做循环相关时，噪声大都集中在alpha=0上，而alpha=0相当于传统功率谱，如果采用传统功率谱方法，无法针对噪声进行抑制，而循环谱却可以在低信噪比的情况下很好的抑制平稳噪声。
2.由于循环谱包含幅度和相位，谱相关函数能为信号分析提供更多的信息（信号载频，带宽，符号速率，相位和时间）。所以采用循环谱方法，可以通过二维空间来提取更多的特征参数。
3.由于大多人造信号都是循环平稳信号，和传统的功率谱分析方法，谱相关分析更准确地揭示令循环平稳信号的本质。
拓展方向：解决频偏和循环谱泄露的问题，减轻计算量，还可以结合小波变换、高阶累积量等方法进行多种通信信号的调制识别，进行多种通信信号的特征提取
4.举例推导
4.1 BPSK
这周主要仿真了BPSK的循环谱，并利用循环谱进行特征分析。提取载频和符号速率。
| 内容 | 符号表示 | 数值 |
|---|---|---|
| 符号数 | N | 100 |
| 信号点数 | length(x) | 10000 |
| 循环谱检测点数 | K | 5000 |
| 平滑 | T | 25 |
| 载波频率 | fc | 125Hz |
| 符号速率 | Fb | 25个/s |
| 信噪比 | SNR | 10dBD |
循环谱：
仿真代码：
N=100; %%每行符号数/符号数
L=32; %%L行并行运输
M=2; %%多进制
SNR=10; %%单位dBD
fc=125; %%采样定理：小于25
fs=2500; %%采样频率
x=PSK(M,fc,SNR,N,fs); %%信号产生
K=N*50; %%循环谱检测采样长度，必须小于等于信号序列长度
T=25; %%平滑点数，检测采样长度K=T*频率采样个数
[alpha,f,S]=cyclic_spectrum(x, K, fs, T); %%做循环谱
S=data_sym(alpha,f,S); %%f轴的负半轴对称
%%% 循环谱作图 %%%
figure;
mesh(f, alpha, abs(S));
title('循环谱');xlabel('f');ylabel('a');
normal_data_alpha=S(floor(length(S(:,1))/2)+1,:); %alpha=0截面
normal_data_f=S(:,floor(length(S(1,:))/2)+1); %f=0截面
figure;
subplot(211);
plot(f,abs(normal_data_alpha));title('谱相关函数a=0的二维切面');
xlabel('频率 [Hz]');ylabel('S0(f)');grid on;
subplot(212);
plot(alpha,abs(normal_data_f));title('谱相关函数f=0的二维切面');
xlabel('循环频率 [Hz]');ylabel('S0(f)');grid on;
%% 载频估计
% normal_data_f =interp(normal_data_f,10);
% alphas = interp(alpha,10);
mid=fix(length(normal_data_f)/2);
for i=1:mid
if(abs(normal_data_f(i)) ==max(abs(normal_data_f(1:mid))))
f1=abs(alpha(i));
end
end
for ii=mid:length(normal_data_f)
if(abs(normal_data_f(ii)) ==max(abs(normal_data_f(mid:end))))
f2=abs(alpha(ii));
end
end
ff=(f1+f2)/4; %%估计的载波频率
data3=S(:,fix(ff)); %%找到f=fc截面
normal_data3=data3/max(data3);
figure;
plot(alpha,abs(normal_data3));title('谱相关函数f=fc的二维切面');
xlabel('频率 [Hz]');ylabel('S0(f)');grid on;
其中PSK产生函数：
function f=PSK(m,fc,SNR,N,fs)
%%产生所要求的PSK信号
% N=100; %%符号个数
% fc=2; %%载波频率
% m=2; %%m进制
% SNR=1; %%单位dBD
% fs=2500; %%采样频率
fb=100;
fc=fc*fb/fs; %%带入设定采样速率，计算实际载频
T=1; %%表达式中的参数,符号速率（一组载波运载fc个符号）
tb=1/fb;
t=0:tb:T-tb;
c=sqrt(2/T)*exp(j*2*pi*fc*t); %%载波信号生成两个正交载波
c1=sqrt(2/T)*cos(2*pi*fc*t);
c2=-sqrt(2/T)*sin(2*pi*fc*t);
snr=10.^(SNR/10); %%转换成单位线性值
msg_PSK=randi(m,1,N);
% msgmod_PSK=dpskmod(msg_PSK-1,m).'; %%dpsk调制
msgmod_PSK=pskmod(msg_PSK-1,m).'; %%psk调制
x=msgmod_PSK*c;
x1=reshape(x.',1,length(msgmod_PSK)*length(c));
s_pow=norm(x1).^2/(N); %%这里是符号的平均功率，而不是每一个采样点
sigma=sqrt(s_pow/(2*snr)); %%根据线性值的信噪比和信号功率求出噪声功率
f=x1+sigma*randn(1,length(x1)); %%加性噪声后的信号
循环谱函数（这里采样的是频域平滑周期图法）：
function [alpha,f,S] = cyclic_spectrum(x, K, fs, T)
%*******************************************************************
%%频域平滑周期图法求循环谱
% x ： 信号
% K ： 循环谱检测采样长度，必须小于等于信号序列长度
% fs ： 采样频率, 检测带宽为-fs/2至fs/2
% T ： 平滑点数, 时间分辨率*频率分辨率=k
%alpha f: alpha轴和f轴
win = 'hamming'; % 平滑窗类型
d_alpha = fs/(K); %% 1/时间分辨率=循环频率分辨率
alpha = -fs:d_alpha:fs; %% 循环频率, 分辨率=1/时间分辨率
a_len = length(alpha); %% 循环频率取样个数
f_len = floor(K/T-1)+1; %% 最大平滑窗个数, 即频率采样个数
f = -(fs/2-d_alpha*floor(T/2)) + d_alpha*T*(0:1:f_len-1);%% 频率采样点位置：-fs/2-fs/2
S = zeros(a_len, f_len); %%初始相关功率谱（负半轴）
i = 1; %%循环
%%% 信号fft变换 %%%
X = fftshift(fft(x(1:K)));
X = X';
%%% 遍历循环频率取值 %%%
for alfa = alpha
interval_f_N = round(abs(alfa)/d_alpha); % 循环频率所对应的频谱序列序号
f_N = floor((K-interval_f_N-T)/T)+1; % 平滑窗的个数
%%% 频域序列平滑模板 %%%
t = 1:T*f_N;
t = reshape(t, T, f_N);
%%% 计算X1,X2 %%%
X1 = X(t);
X2 = X(t+interval_f_N);
%%% 计算谱相关 %%%
St0 = conj(X1).*X2;
St = mean(St0, 1); % 平滑平均
S(i, floor((f_len-f_N)/2)+(1:f_N)) = St/K; % 将结果平移至序列中央以便作图
i = i+1;
end
补f负半轴函数：
function [SNEW]=data_sym(alpha,f,S)
N=length(alpha);
M=length(f);
mid=floor(M/2);
S_FU=zeros(N,M);
for i=1:N
S_FU(i,:)=fliplr(S(i,:));
end
SNEW=[S_FU(:,1:mid) S(:,mid+1:end)];
SNEW(:,mid)=2*S(:,mid);
4.2 OFDM
5.利用循环谱进行特征分析
5.1 BPSK
第一步：找到循环谱的f=0截面（alpha=0截面表示传统的功率谱，也可以找到载频，但是当信噪比过大可能就看不出来）
文献中给的方法是根据f=0截面（图2），峰值位置除以2就是载频fc，理论值是250Hz
这里x的位置有时候会是250HZ，有时候在250徘徊。目前猜测是平滑值采取不当和没有加窗导致的，后面会在改进中再看会不会更准确。通过alpha=0截面（图1），即功率谱密度也能判断出fc
第二步：
找到f=fc截面，我后面写的函数可以找到fc位置，只不过前面估计载频在理论值附件徘徊，所以不准确，有时候程序直接用fc代替fix(ff)。或者打开循环谱找到f=fc的次峰：
可以得到符号速率Fb=25个/s，这里的频率分辨率采取的不是很好，应该还是加窗平滑的问题
1041
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
