??? Индекс превышает размер матрицы PSD Proplem

Question

Эй, ребята, я пытаюсь найти спектральную плотность мощности сигнала .wav, который я записал, который по существу является синусом, погруженным в шум. Предполагается, что функция, которую я написал, записывает все 1024 точки по длине и использует ее для поиска Gxx сигнала путем поиска Gxx на запись, а затем добавления их и деления их на количество записей, лучше объясненных в алгоритме ниже :

a. Пройдите через wav-файл и извлеките первую длину записи (например, от 1 до 1024 точек). (Обратите внимание, что длина записи - это ваш новый «N», поэтому изменение частоты изменяется в соответствии с этим, а не общая длина wav-файла).

b. Выполните обычную функцию PSD на этой записи.

c. Сохраните этот вектор.

d. Извлеките следующие 1024 точки в wav-файле (например, 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 ow, чтобы исправить это, и я пробовал много разных методов, но все же я получаю эту ошибку. Я не очень хорошо знаком с 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-loop или while-loop.