2013-10-03 173 views
-1

嘿家伙我试图找到我记录的.wav信号的功率谱密度,这实际上是一个沉浸在噪声中的正弦波。我写的函数应该记录全部1024个点的记录,并通过查找每个记录的Gxx,然后添加它们并将它们除以在下面的算法中更好地解释的记录数来找到信号的Gxx :???索引超过矩阵尺寸PSD Proplem

a。浏览wav文件并提取第一个记录长度(例如1到1024点)。 (请注意,记录长度是新的“N”,因此频率间隔会根据此而变化,而不是wav文件的总长度)。

b。在此记录上执行正常的PSD功能。

c。存储这个矢量。 d)。提取wav文件中的下一个1024点(例如1025:2048)并在该记录上执行PSD。

e。将其添加到以前存储的记录中,并继续执行步骤c到e,直到达到wav文件的结尾或所需的记录总数。 (请记住,总记录数*记录长度必须小于wavfile的总长度!)

F。用平均数(或记录数)除以PSD。 这是你的平均PSD

我创建的功能如下:

%Function to plot PSD 
function[f1, GxxAv] = HW3_A_Fn_811003472_RCT(x,fs,NumRec) 
Gxx = 0; 
GxxAv = 0; 

N = 1024; 
df = fs/N; 
f1 = 0:df:fs/2; 
dt = 1/fs; 
T = N*dt; 

q = 0; 
e = 1; 

for i = 1:NumRec; 
    for r = (1+q):(N*e); 

     L = x(1+q:N*e); 
     M = length(L); 
     Xm = fft(L).*dt; 
     aXm = abs(Xm); 

     Gxx(1)=(1/T).*(aXm(1).*aXm(1)); 
     for k = 2:(M/2); 
      Gxx(k) = (2/T) *(aXm(k).*(aXm(k))); 
      %Gxx = Gxx + Gxx1(k); 
     end 
     Gxx((M/2)+1)= (1/T)*(aXm((M/2)+1)).*(aXm((M/2)+1)); 
     q = q+1024; 
     e = e+1; 
     %Gxx = Gxx + Gxx1((M/2)+1); 
    end 
    GxxAv = GxxAv + Gxx; 
    %Gxx = Gxx + Gxx1; 
end 
GxxAv = GxxAv/NumRec; 

我用来调用这个函数的代码如下:

[x,fs] = wavread('F:\Final\sem1Y3\Acoustics\sinenoise5s.wav'); 

[f1,GxxAv] = HW3_A_Fn_811003472_RCT(x,fs,100); %where 100 is the number of records to generated 

plot(f1,GxxAv) 

xlabel ('Frequency/Hz', 'fontsize', 18) 
ylabel ('Amplitude Squared per Frequency/WU^2/Hz', 'fontsize', 18) 
title ('Plot of the single sided PSD, using Averaging', 'fontsize', 18) 
grid on 

当试图绘制该图表观察到以下错误:

??? Index exceeds matrix dimensions. 

Error in ==> HW3_A_Fn_811003472_RCT at 19 
     L = x(1+q:N*e); 

Error in ==> HW3_A_3_811003472_RCT at 3 
[f1,GxxAv] = HW3_A_Fn_811003472_RCT(x,fs,100); %where 100 is the number of records to generated 

我不确定h解决它,我已经尝试了许多不同的方法,但仍然得到这个错误。我对Matlab不太熟悉,但是我真正想为第19行做的所有事情就像:

x(1:1024),x(1025:2048),x(2049:3072),x(3072 :4096)...等等到100条记录

任何想法???由于

+0

L = x(i *(N-1):i * N); – durasm

+0

它只是意味着你的输入'x'没有那么多元素。检查函数第一行中的'size(x)'并将其发布到此处。 –

回答

1

这显然是家庭作业,所以我不会那样做你的工作你。但是你的代码有很多错误。通过固定所有这些首先开始:

  • 使用更合适的函数名,homework123不是一个好名字来描述函数的功能。

  • 使用更合适的变量名。在这方面更多的标准是nfft而不是Nn_average而不是NumRec。我不关心你使用的确切的东西,但它应该准确地描述变量的作用。

  • 你的错误信息明确暗示你正试图指数x在一些非法的方式。首先创建一个只打印正确索引(1..1024,1025..2048,...)的循环,并确保它遵循指令E.只有当它按预期工作时,才会添加其余的代码。

  • 您使用三重嵌套的循环。您只需要一个for-loop或while循环来解决这个问题。

+0

谢谢你的建设性批评:)我自己想出来了,不过谢谢你的提示 – user2643355