matlab - ???索引超出矩阵尺寸PSD问题

Question

大家好，我试图找到我记录的.wav信号的功率谱密度，它实际上是一个浸入噪声中的正弦信号。我编写的函数应该记录长度为1024点的所有记录，并使用它来找到信号的Gxx，方法是为每条记录查找Gxx，然后将它们相加并除以记录数，这在下面的算法中有更好的解释:

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

b。在此记录上执行普通的PSD功能。

C。存储此 vector 。

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

我不知道如何解决它，我尝试了许多不同的方法，但仍然收到此错误。我对Matlab不太熟悉，但是对于19行，我真正想要做的就是:

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

有任何想法吗？？？谢谢

Answer 1

显然，这是家庭作业，所以我不会为您完成工作。但是您的代码有很多问题。首先修复所有这些问题:

使用更合适的函数名称，homework123并不是描述函数功能的好名字。

使用更合适的变量名。在这种情况下，更标准的是nfft而不是N和n_average而不是NumRec。我不在乎您使用的是什么，但是它应该准确描述变量的作用。

您的错误消息清楚地表明您正在尝试以某种非法方式为x编制索引。首先创建一个仅打印正确索引的循环(1..1024、1025..2048，...)，并确保它遵循您的指令E。只有在按预期方式工作时，才添加其余代码。

您使用三重嵌套的for循环。您只需要一个for循环或while循环即可解决此问题。