??? Index exceeds matrix dimensions PSD Proplem

Question

Hey guys I am trying to find the Power Spectral Density of a .wav signal i recorded which is essentially a sine immersed in noise. The function that i have written is supposed to take records all of 1024 points in length and use it to find the Gxx of the signal by finding Gxx per record and then adding them and dividing them by the number of records better explained in the algorithm below:

a. Step through the wav file and extract the first record length (e.g. 1 to 1024 points). (Note that the record length is your new “N”, hence the frequency spacing changes in accordance with this, NOT the total length of the wav file).

b. Perform the normal PSD function on this record.

c. Store this vector.

d. Extract the next 1024 points in the wav file (e.g. 1025:2048) and perform PSD on this record.

e. Add this to the previously stored record and continue through steps c to e until you reach the end of your wav file or the total number of records you want. (Remember that total records*record length must be less than the total length of the wavfile!)

f. Divide the PSD by the number of averages (or number of records). This is your averaged PSD

The function I created is as follows:

%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;

And the code I used to call this function is as follows:

[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

When Trying to plot this graph the following error was observed:

??? 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

I am not sure how to fix it and i have tried many different methods but still i get this error. I am not too familiar with Matlab but all I really want to do for line 19 is to go like:

x(1:1024), x(1025:2048), x(2049:3072), x(3072:4096)...etc to 100 records

Any ideas??? Thanks

Answer 1

This is obviously homework, so I am not gonna do your work for you. But there quite some things wrong with your code. Start by fixing all of those first:

• Use more appropriate function names, homework123 is not a good name to describe what the function does.

• Use more appropriate variable names. More standard in this context would be nfft instead of N and n_average instead of NumRec. I don't care about the exact thing you use, but it should describe exactly what the variable does.

• Your error message clearly hints that you are trying to index x in some illegal way. Start with making a loop that just prints the right indices (1..1024, 1025..2048, ...) and make sure it follows your instruction E. Only when this works as expected add the rest of the code.

• you use a triple-nested for-loop. You only need a single for-loop or while-loop to solve this problem.