adafruit / Adafruit_ZeroPDM

Arduino Library for using ATSAMD processors with PDM microphones (w/I2S)
MIT License
15 stars 5 forks source link

Decimation Filter #4

Open mrmarkuese opened 4 years ago

mrmarkuese commented 4 years ago

We discovered that the sinc-filter used for the decimation is too short. It should have a length of at least 128 samples if using downsampling by factor 64. Using a windowed 64 tap sinc-filter doubles the cut-off frequency of the filter and causes heavy aliasing artifacts. Rather heavy aliasing will also occurr for 128 tap filters, and we would recommend sampling the sinc function approximately up to the 6th zero crossing. We've attached Matlab code for visualization. The sinc-filter is determined by the variables si.zeros and si.scale (sorry, no time for markdown and m-file attachement was not supported):

%% parameters ------------------------------------------------------------- close all; clear; clc

fs.PCM = 16e3; % sample rate of the final signal decimation = 128; fs.PDM = fs.PCM * decimation; % PDM sample rate

%% get PCM signal (reference) --------------------------------------------- f = [2e3 4e3 6e3 7e3 8.5e3]; % tested frequencies T = 16; t.PDM = (0:1/fs.PDM:T/max(f))'; t.PCM = (0:1/fs.PCM:T/max(f))'; x.analog = sin(2pif.*t.PDM) / 2; % test signals

%% get PDM signal (sigma delta modulation) --------------------------------

x.PDM = zeros(size(x.analog)); prev = zeros(size(f)); sigma = zeros(size(f)); thresh = 1;

for nn = 1:size(x.analog, 1) delta = x.analog(nn,:) + .5 - prev; sigma = sigma + delta;

x.PDM(nn,:) = sigma >= thresh;

prev         = x.PDM(nn,:);

end

clear prev sigma thresh delta nn

%% get the decimation filter ----------------------------------------------

% reasonable properties si.zeros = 6; % number of zero crossing up to wich the sinc function is evaluated si.scale = .85; % time scaling to decrease the low-pass cut-off frequency (by widening the sinc function) % ada fruit properties % si.zeros = .5; % si.scale = 1;

% generate the sinc function si.x = sinc(linspace(-si.zerossi.scale, si.zerossi.scale, round(2si.zerosdecimation))') ; si.x = si.x . hann(round(2si.zeros*decimation));

% dirty leveling (there is a theretical way, which I was to lazy to look up) si.level = abs(fft(si.x)); si.x = si.x / si.level(1,:);

% delay of the filter in samples @ fs.PDM, and fs.PCM si.delay = round(numel(si.x)/2/decimation) + 1;

%% convert PDM to PCM -----------------------------------------------------

% lowpass x.PCM = fftfilt(si.x, x.PDM);

% decimate x.PCM = x.PCM(1:decimation:end,:);

% remove DC x.PCM = x.PCM - 0.5;

%% plot -------------------------------------------------------------------

% plot time signals figure for ff = 1:numel(f) subplot(numel(f),1,ff) hold on plot(t.PDM/1e3, x.PDM(:,ff)-.5, 'color', [.6 .6 .6]) plot(t.PDM/1e3, x.analog(:,ff), 'k') plot(t.PCM(1:end-si.delay+1)/1e3, x.PCM(si.delay:end, ff), '.r') box on axis([0 t.PCM(end-si.delay+1)/1e3 -.6 .6]) legend('PDM', 'analogue input', 'PCM', 'Location', 'EastOutside') title(['Test signals (input: ' num2str(f(ff)), ' Hz, sampling rate: ' num2str(fs.PCM/1e3) ' kHz)']) xlabel 'Time in ms' ylabel 'Amplitude' end

% plot decimation filter figure subplot(2,1,1) plot(si.x, 'k') title({'Decimation filter' 'time signal'}) xlabel 'Samples @ fs_{PDM}' ylabel 'Amplitude' xlim([1 numel(si.x)]) box on; grid on subplot(2,1,2) df = 100; si_plot = [zeros(round(fs.PDM/df2),1); si.x; zeros(round(fs.PDM/df2),1)]; f_plot = (0:numel(si_plot)-1) fs.PDM/numel(si_plot); hold on plot(f_plot, 20log10(abs(fft(si_plot))), 'k') axis([100 40e3 -60 5]) plot([fs.PCM/2 fs.PCM/2], get(gca, 'Ylim'), '--r') set(gca, 'XScale', 'log') title({'Decimation filter' 'magnitude'}) xlabel 'Frequency in Hz' box on; grid on; ylabel 'Magnitude in dB'

ladyada commented 4 years ago

thanks, please submit a new filter as a PR and we'll try it out!

f-brinkmann commented 4 years ago

We attached filters for various sample rates that have an aliasing suppression of 40 dB at the cost of more filter coefficients. The filters (float numbers!) are stored inside the mat-files together with the PDM and PCM sampling rates.

We assume that a more efficient implementation could be realized using IIR filters at the cost of potentially inaudible phase distortions close to half the sampling frequency (e.g. Chebychev type 2). However, we did not try this.

2020 Adafruit_PDM_conversion.zip

ladyada commented 4 years ago

kk have you tried porting the filters to code? if we can't process a sample within a sample period, it doesnt really matter how good it is :)

f-brinkmann commented 4 years ago

We did not port them to code and they well exceed a sample period. However, I would look at it the other way around: if it's not good it doesn't matter if you can process the sample within a sample period :)

I'm not sure if there is a memory free method for doing this within sample time. From what I know, this can not be achieved with Sinc-filters. IIR filters would be shorter, but also require at least some memory.

ladyada commented 4 years ago

OK the current filter works pretty well - its not great but good enough for many purposes including voice recognition!

f-brinkmann commented 4 years ago

I don't doubt that there are cases in which the filter is sufficient. But I'm sure there are cases were it is not. I might be wrong, but I think the filter has an aliasing suppression of about 1 dB at half the sampling frequency, where it should have 40 dB and more. If I'm not mistaken, I would appreciate a comment to the code to let users now. This way everyone can evaluate if the filter is sufficient for their application.

nicolas-f commented 4 years ago

Hello, I'm working on it.

I measured a total harmonic distortion of 85%. It is too high for my use case (transfer data using audio) thd.zip