Vectors must be the same length

13 Ansichten (letzte 30 Tage)
Nathan Espindola
Nathan Espindola am 18 Nov. 2021
Kommentiert: Star Strider am 19 Nov. 2021
Hello, im having the error Vectors must be the same length for the code below:
%% An example of noise removal from an audio file
clearvars;
close all;
%% Read in the file
[f, Fs] = audioread('handel_J.wav');
T = 1 / Fs; % Sampling period
L = length(f); % Length of signal
t = (0:L-1) * T; % Time vector
%% Plot audio channels
N = size(f, 1);
figure; stem(t, f);
title('Original: Time-domain'); xlabel('time(seconds)');
%% Plot the spectrum
df = Fs / N;
w = (-(N/2):(N/2)-1)*df;
y = fft(f) / N; % For normalizing
y2 = fftshift(y);
figure; plot(w, abs(y2));
title('Original: Amplitude Spectrum'); xlabel('Frequency(Hz)');
%% plot the single-sided amplitude spectrum
figure;
plot(Fs*(0:(L/2))/L, abs(y2(N/2:end)))
title('Original: Single-Sided Amplitude Spectrum')
xlabel('f (Hz)')
%% find f1
% [fmax1, f1index] = max(abs(y2(N/2:end)));
% fmax1, f1index = f1index*df
%% filtering f1: 1043Hz
n = 2;
beginFreq = 1040 / (Fs/2);
endFreq = 1050 / (Fs/2);
[b,a] = butter(n, [beginFreq, endFreq], 'stop');
fOut = filter(b, a, f);
% extract f1
[b,a] = butter(n, [beginFreq, endFreq], 'bandpass');
f1 = filter(b, a, f);
figure; plot(w,abs(fftshift(fft(f1)/N))); title('f1: 1043Hz');
%% find f2
% [fmax2, f2index] = max(abs(fftshift(fft(fOut)/N)));
% fmax2, f2index = Fs/2 - f2index*df
%% filtering f2: 1956Hz
n = 2;
beginFreq = 1950 / (Fs/2);
endFreq = 1960 / (Fs/2);
[b,a] = butter(n, [beginFreq, endFreq], 'stop');
fOut = filter(b, a, fOut);
% extract f2
[b,a] = butter(n, [beginFreq, endFreq], 'bandpass');
f2 = filter(b, a, f);
figure; plot(w,abs(fftshift(fft(f2)/N))); title('f2: 1956Hz');
%% For normalizing
Ym = max(max(max(fOut)),max(abs(min(fOut))));
fOut = fOut ./ Ym;
%% After processing
figure; stem(t, fOut);
title('After processing: Time-domain'); xlabel('time(seconds)');
figure; plot(w,abs(fftshift(fft(fOut)/N)));
title('After processing: Amplitude Spectrum'); xlabel('Frequency(Hz)');
%% Create object for playing audio
pOrig = audioplayer(f,Fs); % original signal
% pOrig.play;
p = audioplayer(fOut, Fs); % filtered signal
% p.play;
pf1 = audioplayer(f1, Fs); % noise: f1: 1043Hz
% pf1.play;
pf2 = audioplayer(f2, Fs); % noise: f2: 1956Hz
% pf2.play;
%% Write out to audio file
audiowrite('handel_J_processed.wav', fOut, Fs);
audiowrite('f1.wav',f1,Fs);
audiowrite('f2.wav',f2,Fs);

Antworten (1)

Star Strider
Star Strider am 18 Nov. 2021
I can’t run the code, however I believe the problem is with the legth of ‘w’ (since the call that threw that error is not stated). The time vectors appear to be correct, so I doubt those plots are throuwing the error.
For two-sided fft plots using fftshift, I usually define the frequency vector as going from the negative Nyquist frerquency to the positive Nyquist frequency with the vector length equal to the length of the fft argument. That ’s just easier.
w = linspace(-Fs/2, Fs/2, L);
That should work, and should solve the unequal vector lengths problem.
.
  2 Kommentare
Nathan Espindola
Nathan Espindola am 18 Nov. 2021
The error is at line 27 and the audio file is attach.
I still coudnt find the problem.
Star Strider
Star Strider am 19 Nov. 2021
The frequency vector and amplitude length mismatch were the problem.
Try this —
uz = unzip('https://www.mathworks.com/matlabcentral/answers/uploaded_files/805604/handel_J.zip')
uz = 1×1 cell array
{'handel_J.wav'}
uzv = uz{1}
uzv = 'handel_J.wav'
[f,Fs] = audioread(uzv)
f = 354221×2
1.0e+00 * 0 0 0 -0.0000 0.0000 0.0001 -0.0000 -0.0000 0.0000 0.0000 0 0 0 0 0.0000 0 -0.0000 0.0000 0.0001 0
Fs = 44100
T = 1 / Fs; % Sampling period
L = length(f); % Length of signal
t = (0:L-1) * T; % Time vector
%% Plot audio channels
N = size(f, 1);
figure; stem(t, f);
title('Original: Time-domain'); xlabel('time(seconds)');
%% Plot the spectrum
df = Fs / N;
% w = (-(N/2):(N/2)-1)*df;
w = linspace(-Fs/2, Fs/2, L);
y = fft(f) / N; % For normalizing
y2 = fftshift(y);
figure; plot(w, abs(y2));
title('Original: Amplitude Spectrum'); xlabel('Frequency(Hz)');
%% plot the single-sided amplitude spectrum
figure;
% plot(Fs*(0:(L/2))/L, abs(y2(N/2:end)))
w1 = linspace(0, 1, fix(L/2)+1);
Iv = 1:numel(w1);
plot(w1, abs(y(Iv))*2)
title('Original: Single-Sided Amplitude Spectrum')
xlabel('f (Hz)')
%% find f1
% [fmax1, f1index] = max(abs(y2(N/2:end)));
% fmax1, f1index = f1index*df
%% filtering f1: 1043Hz
n = 2;
beginFreq = 1040 / (Fs/2);
endFreq = 1050 / (Fs/2);
[b,a] = butter(n, [beginFreq, endFreq], 'stop');
fOut = filter(b, a, f);
% extract f1
[b,a] = butter(n, [beginFreq, endFreq], 'bandpass');
f1 = filter(b, a, f);
figure; plot(w,abs(fftshift(fft(f1)/N))); title('f1: 1043Hz');
%% find f2
% [fmax2, f2index] = max(abs(fftshift(fft(fOut)/N)));
% fmax2, f2index = Fs/2 - f2index*df
%% filtering f2: 1956Hz
n = 2;
beginFreq = 1950 / (Fs/2);
endFreq = 1960 / (Fs/2);
[b,a] = butter(n, [beginFreq, endFreq], 'stop');
fOut = filter(b, a, fOut);
% extract f2
[b,a] = butter(n, [beginFreq, endFreq], 'bandpass');
f2 = filter(b, a, f);
figure; plot(w,abs(fftshift(fft(f2)/N))); title('f2: 1956Hz');
%% For normalizing
Ym = max(max(max(fOut)),max(abs(min(fOut))));
fOut = fOut ./ Ym;
%% After processing
figure; stem(t, fOut);
title('After processing: Time-domain'); xlabel('time(seconds)');
figure; plot(w,abs(fftshift(fft(fOut)/N)));
title('After processing: Amplitude Spectrum'); xlabel('Frequency(Hz)');
Success!
.

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Time-Frequency Analysis finden Sie in Help Center und File Exchange

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by