Plotting FFT for audio WAV file?

Question

Avinash Kandalam 2020-11-4

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/636155-plotting-fft-for-audio-wav-file

评论： Walter Roberson 2024-8-27

采纳的回答： Mathieu NOE

myWAVaudiofile.zip

在 MATLAB Online 中打开

Dear all,

I tried to explain as clear as possible. I want to plot "Raw FFT" file for a "WAV" file. This WAV (audio) file is acquired from a microphone for a period of 1 minute. The goal is to plot frequency distribution (0 Hz - 20 kHz).

I want to acquire raw FFT (to see if there are any signal peaks at particular frequency) throughout 1 minute. The WAV (audio) file (only 1) is atttached to this question.
Please help me with the code and the output graph.

I tried to execute the following code (from previous answers here) and I think it is not the right way. I think what the code shows is basically amplitude vs frequency; but not a typical FFT spectrum.

Million Thanks,

Avinash.

CODE: I tried and most likely wrong. I think as said, it is just amp vs freq, which does not give me clear picture of frequencies which lies in different ranges.

        [y1,fs]=audioread('myWAVaudiofile.wav');
        t=linspace(0,length(y1)/fs,length(y1));
        Nfft=16777216; %power of 2 and I put a huge number so there are many data points
        f=linspace(0,fs,Nfft);
        X1=abs(fft(y1,Nfft));
        plot(f(1:Nfft/2),X1(1:Nfft/2))
        xlabel('Frequency'); 
        ylabel ('Power???');
        title ('FFT Spectrum');

OUTPUT: I only zoomed into 0-30 Hz using above code and WAV file attached (ofcourse the wole spectrum is until 20000 Hz)

6 个评论
显示 4更早的评论隐藏 4更早的评论

Avinash Kandalam 2020-11-4

编辑：Avinash Kandalam 2020-11-4

在 MATLAB Online 中打开

@Walter Roberson: Thanks for response. 1. The code which I used was freq vs amp? 2. Spectogram, meaning freq vs dBV? I'm confused as there are many possibilites how to represent a "WAV" file in FFT.. Is there a general/basic way which people represent WAV file?

I would like to know the difference-seems fundamental understnading issue between types of FFT's with respect to MATLAB coding, could you enlighten little more in detail about what is basic code/results people do under such cicumstances. Thanks.

Another code which I found is: I used the same WAV file as attached to the original question

 [y,Fs] = audioread('myWAVaudiofile.wav');            % Sampling frequency                    
T = 1/Fs;             % Sampling period       
L = 1500;             % Length of signal
t = (0:L-1)*T;        % Time vector
S = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
X = S + 2*randn(size(t));
plot(1000*t(1:50),X(1:50))
title('Signal Corrupted with Zero-Mean Random Noise')
xlabel('t (milliseconds)')
ylabel('X(t)')
Y = fft(S);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
plot(fs,P1) 
title('FFT spectrum')
xlabel('f (Hz)')
ylabel('|P1(f)|')

OUTPUT: Obviously, different output

Walter Roberson 2020-11-4

fft() is inherently a function for processing periodic signals. The assumption is that the signal continues on forever. Remember that you are decomposing into the sum of phased sine or cosine signals, and those have no endpoint.

fft() by itself is therefore not useful for processing speech, as speech consists mostly of changing information.

However, fft() can still be useful -- provided that you apply it to a window of data that is wide enough to be able to examine frequencies of interest, but which is short enough to not "drag down" the changing nature of the signal.

A typical way to handle this is to use Short-Time FFT (SFFT), or Spectrograms. Spectrograms use SFFT and some graphical representation of power, to give an idea of how power at different frequencies changes over time.

Avinash Kandalam 2020-11-4

@Walter: Thank you for the reply. The content is basically over the signal recorded, there are different sounds happening at different frequencies.

For example, at 30 Hz and 40 Hz, if I see a strong peak (in FFT) then I believe the frequencies produced are mostly at 30 & 50 Hz. Something like this, I want to see my WAV signal to be distinguished. Thanks

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Mathieu NOE 2020-11-4

3
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/636155-plotting-fft-for-audio-wav-file#answer_534105

在 MATLAB Online 中打开

dear friends, here my little contribution to wav file spectral analysis...

enjoy !

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FFT parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
NFFT = 8192;    % 
NOVERLAP = round(0.75*NFFT);
w = hanning(NFFT);
% spectrogram dB scale
spectrogram_dB_scale = 100;  % dB range scale (means , the lowest displayed level is XX dB below the max level)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% options 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% if you are dealing with acoustics, you may wish to have A weighted
% spectrums 
% option_w = 0 : linear spectrum (no weighting dB (L) )
% option_w = 1 : A weighted spectrum (dB (A) )
option_w = 0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% load signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% [signal,Fs]=wavread('myWAVaudiofile.wav');  %(older matlab)
% or 
[data,Fs]=audioread('myWAVaudiofile.wav'); %(newer matlab)
channel = 1;
signal = data(:,channel);
samples = length(signal);
dt = 1/Fs;
t = (0:dt:(samples-1)*dt);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 1 : averaged FFT spectrum
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[sensor_spectrum, freq] = pwelch(signal,w,NOVERLAP,NFFT,Fs);
% convert to dB scale (ref = 1)
sensor_spectrum_dB = 20*log10(sensor_spectrum);
% apply A weigthing if needed
if option_w == 1
    pondA_dB = pondA_function(freq);
    sensor_spectrum_dB = sensor_spectrum_dB+pondA_dB;
    my_ylabel = ('Amplitude (dB (A))');
else
    my_ylabel = ('Amplitude (dB (L))');
end
figure(1),semilogx(freq,sensor_spectrum_dB);grid
title(['Averaged FFT Spectrum  / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(freq(2)-freq(1)) ' Hz ']);
xlabel('Frequency (Hz)');ylabel(my_ylabel);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 2 : time / frequency analysis : spectrogram demo
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[sg,fsg,tsg] = specgram(signal,NFFT,Fs,w,NOVERLAP);  
% FFT normalisation and conversion amplitude from linear to dB (peak)
sg_dBpeak = 20*log10(abs(sg))+20*log10(2/length(fsg));     % NB : X=fft(x.*hanning(N))*4/N; % hanning only
% apply A weigthing if needed
if option_w == 1
    pondA_dB = pondA_function(fsg);
    sg_dBpeak = sg_dBpeak+(pondA_dB*ones(1,size(sg_dBpeak,2)));
    my_title = ('Spectrogram (dB (A))');
else
    my_title = ('Spectrogram (dB (L))');
end
% saturation of the dB range : 
% saturation_dB = 60;  % dB range scale (means , the lowest displayed level is XX dB below the max level)
min_disp_dB = round(max(max(sg_dBpeak))) - spectrogram_dB_scale;
sg_dBpeak(sg_dBpeak<min_disp_dB) = min_disp_dB;
% plots spectrogram
figure(2);
imagesc(tsg,fsg,sg_dBpeak);colormap('jet');
axis('xy');colorbar('vert');grid
title([my_title ' / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(fsg(2)-fsg(1)) ' Hz ']);
xlabel('Time (s)');ylabel('Frequency (Hz)');
function pondA_dB = pondA_function(f)
% dB (A) weighting curve
n = ((12200^2*f.^4)./((f.^2+20.6^2).*(f.^2+12200^2).*sqrt(f.^2+107.7^2).*sqrt(f.^2+737.9^2)));
r = ((12200^2*1000.^4)./((1000.^2+20.6^2).*(1000.^2+12200^2).*sqrt(1000.^2+107.7^2).*sqrt(1000.^2+737.9^2))) * ones(size(f));
pondA = n./r;
pondA_dB = 20*log10(pondA(:));
end

18 个评论
显示 16更早的评论隐藏 16更早的评论

George 2024-8-27

编辑：Walter Roberson 2024-8-27

在 MATLAB Online 中打开

The coded example from above doesnt work in 2024. Inline functions are not supported.

Moving forward, Anonymouse functions are supported.

Ive re coded the example, and provided automatic scaling for the FFT Y axis, and the ability to define the lower limits of the Y range.

In this case Ive allowed for 5dB of headroom on the FFT max and a 40 dB range of signal for the lower limit. Feel free to change these.

=======================================

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FFT parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%    FFT  spectral resolution
NFFT =  2048 ;  %  8192;   
% This changes temporal resolution and is used in Spectrogram 
NOVERLAP = round(0.50*NFFT);     
w = hanning(NFFT);
% spectrogram dB scale
Spectrogram_dB_scale = 100;  % eg  100 dB range scale (means , the lowest displayed level is XX dB below the max level)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% options 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% you may wish to have A weighted spectrum
% option_w = 0 : linear spectrum (no weighting dB (L) )
% option_w = 1 : A weighted spectrum (dB (A) )
option_w = 0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% load signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%(older matlab)
% [signal,Fs]=wavread('myWAVaudiofile.wav'); 
%(newer matlab)
[data,Fs]=audioread('myWAVaudiofile.wav'); 
channel = 1;
signal = data(:,channel);
samples = length(signal);
dt = 1/Fs;
t = (0:dt:(samples-1)*dt);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 1 : averaged FFT spectrum
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[sensor_spectrum, freq] = pwelch(signal,w,NOVERLAP,NFFT,Fs);
% convert to dB scale (ref = 1)
sensor_spectrum_dB = 20*log10(sensor_spectrum);
% apply A weighting if needed
if option_w == 1
    pondA_dB = pondA_function(freq);
    sensor_spectrum_dB = sensor_spectrum_dB+pondA_dB;
    my_ylabel = ('Amplitude (dB (A))');
else
    my_ylabel = ('Amplitude (dB (L))');
end
figure(1)
semilogx(freq, sensor_spectrum_dB)
grid on
% Find the peak value of the spectrum
peak_value = max(sensor_spectrum_dB);
% Set the Y-axis limits,  +5 sets  5db of headroom on peak 
% 40 sets Y range to 40dB down on peak signal 
ylim([peak_value  - 40, peak_value + 5])
title(['Averaged FFT Spectrum  with  Fs = ' num2str(Fs) ' Hz with delta f  = ' num2str(freq(2)-freq(1)) ' Hz '])
xlabel('Frequency (Hz)')
ylabel(my_ylabel)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 2 : time / frequency analysis : spectrogram 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[sg,fsg,tsg] = specgram(signal,NFFT,Fs,w,NOVERLAP);  
% FFT normalisation and conversion amplitude from linear to dB (peak)
sg_dBpeak = 20*log10(abs(sg))+20*log10(2/length(fsg));     % NB : X=fft(x.*hanning(N))*4/N; % hanning only
% apply A weigthing if needed
if option_w == 1
    pondA_dB = pondA_function(fsg);
    sg_dBpeak = sg_dBpeak+(pondA_dB*ones(1,size(sg_dBpeak,2)));
    my_title = ('Spectrogram (dB (A))');
else
    my_title = ('Spectrogram (dB (L))');
end
% saturation of the dB range : 
% saturation_dB = 60;  % dB range scale (means , the lowest displayed level is XX dB below the max level)
min_disp_dB = round(max(max(sg_dBpeak))) - spectrogram_dB_scale;
sg_dBpeak(sg_dBpeak<min_disp_dB) = min_disp_dB;
% plots spectrogram
figure(2);
imagesc(tsg,fsg,sg_dBpeak);colormap('jet');
axis('xy');colorbar('vert');grid
title([my_title ' / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(fsg(2)-fsg(1)) ' Hz ']);
xlabel('Time (s)');ylabel('Frequency (Hz)');
% dB (A) weighting curve
% Define the function for n
n = @(f) ((12200^2*f.^4)./((f.^2+20.6^2).*(f.^2+12200^2).*sqrt(f.^2+107.7^2).*sqrt(f.^2+737.9^2)));
% Define r (constant value)
r = ((12200^2*1000^4)./((1000^2+20.6^2).*(1000^2+12200^2)*sqrt(1000^2+107.7^2)*sqrt(1000^2+737.9^2)));
% Define pondA function
pondA = @(f) n(f) ./ r;
% Define PondA_dB function
PondA_dB = @(f) 20*log10(pondA(f));

Walter Roberson 2024-8-27

The coded example from above doesnt work in 2024. Inline functions are not supported.

I do not see anywhere that inline functions were used?

请先登录，再进行评论。

Plotting FFT for audio WAV file?

6 个评论
显示 4更早的评论隐藏 4更早的评论

采纳的回答

18 个评论
显示 16更早的评论隐藏 16更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

Plotting FFT for audio WAV file?

6 个评论 显示 4更早的评论隐藏 4更早的评论

采纳的回答

18 个评论 显示 16更早的评论隐藏 16更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

6 个评论
显示 4更早的评论隐藏 4更早的评论

18 个评论
显示 16更早的评论隐藏 16更早的评论