clc

clear all

NFFT = 1024;

NOVERLAP = round(0.75*NFFT);

w = hanning(NFFT);

spectrogram_dB_scale = 80;

option_w = 0;

[data,Fs]=audioread('myWAVaudiofile.wav');

channel = 1;

signal = data(:,channel);

samples = length(signal);

dt = 1/Fs;

[sensor_spectrum, freq] = pwelch(signal,w,NOVERLAP,NFFT,Fs);

sensor_spectrum_dB = 20*log10(sensor_spectrum);

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

[sg,fsg,tsg] = specgram(signal,NFFT,Fs,w,NOVERLAP);

sg_dBpeak = 20*log10(abs(sg))+20*log10(2/length(fsg));

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

min_disp_dB = round(max(max(sg_dBpeak))) - spectrogram_dB_scale;

sg_dBpeak(sg_dBpeak<min_disp_dB) = min_disp_dB;

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)

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

## Best Answer