A-weighting filter to lookup table

20 views (last 30 days)
Davide
Davide on 28 Aug 2013
Commented: Daniel Shub on 16 Jan 2014
hello,
i'm working at an embedded SPL meter, it has to filter out the audio using A-weighting.
I've build a script to test A-weighting filters, but here comes the problem.
My idea is to build a lookup table, that contains the weight value for a specific frequency, then use this table to weighting the signal using fft.
I'm investigating three filters, one build upong the fdesign function, one two using filter coefficients based functions.
I've two question and one problem:
1) is my method correct?
2) what filter table i've to use? if i test it seems that the two coef based filters differs from an offset.
download sample code here: https://www.dropbox.com/s/62ueiv3rmc8l3tc/test.zip
thank you!
  1 Comment
Daniel Shub
Daniel Shub on 8 Jan 2014
It is not clear what you are doing, what you are getting, and what you expect. Without this information, it is hard to answer your question.

Sign in to comment.

Answers (4)

Wayne King
Wayne King on 2 Jan 2014
How about using
fdesign.audioweighting
if you have the DSP System Toolbox? That will give you A-weighting.

Al Hasan
Al Hasan on 2 Jan 2014
Hey Davide,
Was your method correct ?
I have attached an image is that A weighted filtering
By any chance do you know what they are doing to the above image to get the below image in my attachment, I have data to make the above image dunno how to make the below one from the above
Thanks, Hasan

Davide
Davide on 8 Jan 2014
no-one reply me, so i'm not whure my method worked. one of my questions, is about the use of fdesign.audioweighting, i've found other implementation of the A-weighting, that differs from DSP matlab one. i'm asking in my question why it differs. i mean A-weighting should always the same.
  1 Comment
Wayne King
Wayne King on 8 Jan 2014
You should provide more detail about your above statement.

Sign in to comment.


Davide
Davide on 11 Jan 2014
i'm using this script
%clear matlab workspace
clear;
%build the filter using fdesign
Fsf = 22050; %frequency used for filter creation
HawfA = fdesign.audioweighting('WT,Class','A',1,Fsf);
Afilter = design(HawfA);
%plot filter 1
%get filter freq resp
[AH, AW] = freqz(Afilter);
%plot filter
figure(1)
semilogx(AW*1000,20*log10(abs(AH)));
xlabel('Frequency [Hz]'); ylabel('Gain [dB]');
title(['A-weighting, Fs = ',int2str(Fsf),' Hz']);
legend('Filter');
grid on
%gaintable size
gaintablesize = 32;
%build filter lookup table
Flookup = 22050; %frequency used for lookuptable creation
WindowSize = gaintablesize*2; %samples for lookuptable creation / window length
filter_gaintable = freqz(Afilter, WindowSize, 'whole'); %built for the sized window
%build test signal
Fs = 22050; %sampling frequency of signal
T = 1/Fs; %sample time
L = WindowSize*20; %length of signal (20 * window)
t = (0:L-1)*T; %time vector
y = 0.7*sin(2*pi*50*t) + sin(2*pi*2000*t) + 2*randn(size(t)); %sum of a 50Hz and a 2000 Hz sinusoid plus noise
signal = y';
%estimate signal level (within each windowed segment).
windowIntervals = [1:WindowSize:(length(signal)-WindowSize)];
windowTime = t(windowIntervals+round((WindowSize-1)/2));
%will contain the full signal
signalT = [];
filteredsignalT = [];
filteredsignalfreqrespT = [];
%will contain the dBA computations
dBA1 = zeros(length(windowIntervals),1);
dBA2 = zeros(length(windowIntervals),1);
dBA3 = zeros(length(windowIntervals),1);
for i = [1:length(windowIntervals)]
%take a part of the signal
signalP = signal(windowIntervals(i)-1+[1:WindowSize]);
%filter signal part by filter function
filteredsignalP = filter(Afilter, signalP);
%filter signal part by frequency response
NFFTP = length(signalP);
fftsignalP = fft(signalP);
fftfilteredsignalP = fftsignalP .* filter_gaintable;
filteredsignalfreqrespP = ifft(fftfilteredsignalP, NFFTP);
%built the full signal
signalT = cat(1, signalT, signalP);
filteredsignalT = cat(1, filteredsignalT, filteredsignalP);
filteredsignalfreqrespT = cat(1, filteredsignalfreqrespT, filteredsignalfreqrespP);
%calcuatate rms of the signal part
rmssignalP = sqrt(mean(signalP.^2));
rmsfilteredsignalP = sqrt(mean(filteredsignalP.^2));
rmsfilteredsignalfreqrespP = sqrt(mean(filteredsignalfreqrespP.^2));
% Estimate dBA value using Parseval's relation.
dBA1(i) = 20*log10(rmssignalP);
dBA2(i) = 20*log10(rmsfilteredsignalP);
dBA3(i) = 20*log10(rmsfilteredsignalfreqrespP);
end
% %plot full signal, and filtered signal
% figure(2)
% plot([1:length(signalT)], signalT, '-', [1:length(filteredsignalT)], filteredsignalT, '-', [1:length(filteredsignalfreqrespT)], filteredsignalfreqrespT);
% title(['compare signal (1 window), to filtered signal']);
% legend('signal', 'filtered (using filter function)', 'filtered (using response table)');
% grid on
%plot out results
figure(3)
plot([1:length(dBA1)],dBA1,'-',[1:length(dBA2)],dBA2,'-',[1:length(dBA3)],dBA3);
legend('dBA signal','dBA filtered (using filter function)', 'dBA filtered (using weight table)',3);
grid on
%save gaintable to xls
C(:,1) = filter_gaintable(1:length(filter_gaintable)/2); %retain frequencies below Nyquist rate
M = xlswrite('aweightingtable.xls', C);
  1 Comment
Daniel Shub
Daniel Shub on 16 Jan 2014
Okay, and what happens when you run the code? What are you expecting to happen? What would you like to happen? Ideally you should strip the code down to be as minimal as possible. Right now it looks like it should plot 3 figures. If those figures are not related to the problem, don't post that part of the code. You right to an xls file, presumably this isn't the problem. The simpiler the code, the more easily we can help you.

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!