声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 6760|回复: 15

[其他] 请帮我看看广义互相关时延估计程序

[复制链接]
发表于 2006-12-25 19:30 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
用多个正弦产生仿真信号,两个阵元接受到信号后,估计两信号的延迟

SN1=20; %信噪比
nn=4096; %信号长度
tt=1:nn;
fs=1000;
S=3*sin(2*pi*17*tt/fs)+2.5*sin(2*pi*35*tt/fs)+1.5*sin(2*pi*51*tt/fs)+sin(2*pi*102*tt/fs); %采样频率为1000

Ps=S*S'/nn;
ps=diag(Ps);
refp=2*10.^(SN1/10);
tmp=sqrt(refp./ps);
SS=diag(tmp)*S;
%SS=SS+randn(size(SS)); 。。。。。。。。。1

for j=1:nn
S2(j+549)=1/dist(2)*SS(j); %传送到阵2的时候延迟了549个点
end
%S2=S2+randn(size(S2));。。。。。。。。。。2

for j=1:nn
S3(j+580)=1/dist(3)*SS(j); %信号传送到阵3为延迟了580个点
end%
%S3=S3+randn(size(S3));。。。。。。。。。。。。3

n=200; %取信号的长度不超过200ms
m=580-549;
SS2=[S2(tao(2)+1:tao(2)+1+n),zeros(1,m)];
SS3=[zeros(1,m),S3(tao(3)+1:tao(3)+1+n)]; %做相关的信号

SS2=[SS2,zeros(1,length(SS2))];
SS3=[zeros(1,length(SS3)),SS3];

fx2=fft(SS2);fx3=fft(SS3);
Pxy=fx2.*conj(fx3)./(abs(fx2).*abs(fx3)+eps); 。。。。。。4

R=ifft(Pxy);

R=real(R);
c=max(abs(R));
a=find(abs(R)==c) %找出最大点
m=(a-length(R)/2)*1/fs %换算出信号的延迟时间
figure(3)
plot(R)

信号2与3的时延为0。0306秒,我现在如果把噪声加在 1 处的时候m可以得比较接近的值,但是如果把噪声加在2 和 3 处就不行,而且每次估计的 m值差别变化很大,但我认为在 2 3 处加噪声比较符合逻辑,因为每个阵的噪声应该是不相关的 ,请问是什么原因? 噪声的影响能有这样大吗? 另外我的采样频率是否有问题?

另外 4 处求的是两信号的互功率谱,请问互功率谱和功率谱密度的区别,如果去掉 4 处的分母,可以吗?

本帖被以下淘专辑推荐:

回复
分享到:

使用道具 举报

发表于 2007-1-18 17:35 | 显示全部楼层

我也在做这个东西,你的tao是多少?没写清楚啊?

我也在做这个东西,你的tao是多少?没写清楚啊?
咱们可以交流交流
发表于 2007-1-19 22:34 | 显示全部楼层
我的QQ是51492475,加我交流交流啊,呵呵
发表于 2007-4-1 15:08 | 显示全部楼层
我的QQ:49312483  加我一个吧
发表于 2007-4-10 14:25 | 显示全部楼层

有兴趣

我也在做这个,qq191174105希望能探讨
发表于 2008-5-11 19:42 | 显示全部楼层

毕设相关

我也在做这个,还不是很清楚,我的qq361141239
希望共同探讨!!
发表于 2008-5-14 09:44 | 显示全部楼层
本帖最后由 VibInfo 于 2016-11-9 15:07 编辑
原帖由 pjj 于 2006-12-25 19:30 发表
用多个正弦产生仿真信号,两个阵元接受到信号后,估计两信号的延迟
...
信号2与3的时延为0。0306秒,我现在如果把噪声加在 1 处的时候m可以得比较接近的值,但是如果把噪声加在2 和 3 处就不行,而且每次估计的 m值差别变化很大,但我认为在 2 3 处加噪声比较符合逻辑,因为每个阵的噪声应该是不相关的 ,请问是什么原因? 噪声的影响能有这样大吗? 另外我的采样频率是否有问题?

另外 4 处求的是两信号的互功率谱,请问互功率谱和功率谱密度的区别,如果去掉 4 处的分母,可以吗?

我对于dist函数不熟悉,用了它以后,在SS2和SS3中前500多个数都为0。我把这部分取消了,程序为
SN1=20; %信噪比
nn=4096; %信号长度
tt=1:nn;
fs=1000;
S=3*sin(2*pi*17*tt/fs)+2.5*sin(2*pi*35*tt/fs)+1.5*sin(2*pi*51*tt/fs)+sin(2*pi*102*tt/fs); %采样频率为1000
Ps=S*S'/nn;
ps=diag(Ps);
refp=2*10.^(SN1/10);
tmp=sqrt(refp./ps);
SS=diag(tmp)*S;
SS=SS+randn(size(SS)); %。。。。。。。。。1
S2=SS+randn(size(SS)); %。。。。。。。。。。2
S3=SS+randn(size(SS)); %。。。。。。。。。。。。3

n=200; %取信号的长度不超过200ms
m=580-549;
tao(2)=0; tao(3)=0;
SS2=[S2(tao(2)+1:tao(2)+1+n),zeros(1,m)];
SS3=[zeros(1,m),S3(tao(3)+1:tao(3)+1+n)]; %做相关的信号
SS2=[SS2,zeros(1,length(SS2)-1)];
SS3=[zeros(1,length(SS3)-1),SS3];
fx2=fft(SS2);fx3=fft(SS3);
Pxy=fx2.*conj(fx3)./(abs(fx2).*abs(fx3)+eps); %。。。。。。4
R=ifft(Pxy);
R=real(R);
c=max(abs(R));
a=find(abs(R)==c) %找出最大点
Rl2=fix(length(R)/2);
m=(a-(Rl2+2))/fs %换算出信号的延迟时间
figure(3);
time=(-Rl2:Rl2)/fs;
plot(time,R); grid;
对楼主谈到的,在1,2,3中都增加了噪声,一样能得到很好的峰值。见下图,延迟时间是-0.031
同时对4是这样的,若单独求互谱,应为
Pxy=fx2.*conj(fx3)
可楼主是用SCOT方法求广义相关,SCOT方法为
Rxy=Pxy./sqrt(Pxx.*Pyy)
楼主的关系4实际上是把上两式合并在一起。所以分母不能去掉,去掉了就不是求广义相关函数了。

[ 本帖最后由 songzy41 于 2008-5-14 10:21 编辑 ]
pj21a.jpg
发表于 2009-9-14 15:02 | 显示全部楼层
功率谱密度也简称为功率谱
回复 支持 1 反对 0

使用道具 举报

发表于 2011-3-26 23:59 | 显示全部楼层
大家有人用别的方法算时延的吗》?如高阶谱,小波
发表于 2011-9-23 20:04 | 显示全部楼层
我写了一个根据互功率谱来求时延的MATALB程序,如果我先加噪声,再延时得到两路信号,那得出的结果比较贴切;如果先延时,再加噪声,就得不到像样的结果。我感觉应该先延时再加噪声,比较符合实际,可是这又得不到理想的结果,很是费解啊!我看过你在这方面的帖子的评价,所以来请教!期待您的回复
clear;
fs=65536;
c=340;
S0=5*cos(2*pi*17*(0:1324)/fs)+cos(2*pi*117*(0:1324)/fs);
%S0=awgn(S0,10,'measured'); 先加噪声,后延时;
S1=S0(1:1024);
S2=S0(200:1223);
S1=awgn(S1,20,'measured');%先延时,在给1路信号加噪声
S2=awgn(S2,20,'measured');%先延时,在给2路信号加噪声
figure(1);
subplot(2,1,1);
plot(S1);
subplot(2,1,2);
plot(S2);
h=hamming(1024)';
S1=S1.*h;
S2=S2.*h;
figure(2);
subplot(2,1,1);
plot(S1);
subplot(2,1,2);
plot(S2);
SK1=fft(S1,1024);
SK2=fft(S2,1024);
SK1j=conj(SK1);
SK2j=conj(SK2);
fai1=1./(abs(SK2.*SK1j)+eps);
%fai2=1./(sqrt(SK2.*SK2j.*SK1.*SK1j)+eps);
R12=fai1.*SK2.*SK1j;
r12=real(ifft(R12,1024));
find(r12==max(r12))
figure(3);
plot(r12);
发表于 2013-1-12 21:33 | 显示全部楼层

也遇到同样问题 也解决不了! scot phat 基本互相关可以得到正确延迟,加权后都得不到正确延迟!
发表于 2017-5-6 10:47 | 显示全部楼层
xiaoaq 发表于 2011-3-26 23:59
大家有人用别的方法算时延的吗》?如高阶谱,小波

我用小波
发表于 2017-5-6 14:43 | 显示全部楼层
songzy41 发表于 2008-5-14 09:44
我对于dist函数不熟悉,用了它以后,在SS2和SS3中前500多个数都为0。我把这部分取消了,程序为
SN1=20;  ...

R=ifft(Pxy);经过逆傅里叶变换后, R已经是实数不需要R=real(R);
发表于 2017-5-6 19:33 | 显示全部楼层
由于ifft计算中有误差存在,使ifft有时会有虚部,只是很小。所以在这种情况下一般 取实数。
发表于 2017-5-20 13:40 | 显示全部楼层

我也用小波, 能加一下QQ吗?我312411529
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-4-19 15:34 , Processed in 0.082836 second(s), 24 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表