Downsampling

Downsampling czyli zmniejszenie częstotliwości próbkowania sygnału o zadany całkowity współczynnik przy zachowaniu jego widma (poniżej częstotliwości do której zmniejszamy). Ponieważ tutaj rekonstrukcja odbywa się z częstotliwością mniejszą o całkowitą wartość w przeciwieństwie do rekonstrukcji wystarczy nam tylko jeden impuls sinc. Cały kod jest znacznie szybszy.

Funkcja:

function [ fout, tout ] = downsampling( fin, fsin, fsout )

scale = fsin / fsout;

if (fix(scale) ~= scale)
    error('scale should be integer');
end

if (scale <= 0)
    error('scale should be greater than zero');
end

% works only with odd window.
k = sinc_gen(fsin, fsout/2);
if (mod(size(k, 2), 2) == 0)
    k = sinc_gen(fsin, fsout/2, size(k, 2) + 1);
end

tout = linspace(0, (1/fsin) * (size(fin, 2)-1), ... 
                ceil(size(fin, 2)/scale));
fout = zeros(1, size(tout, 2));

for m=1:1:size(fout, 2)
    fout(m) = conv_point(fin, k, (m - 1) * scale + 1);
end

end

Przykład:

%% input data

% resampling
fsin = 30000;
fsout = 6000;
scale = fsin / fsout;

func = @(ti) sin(2*pi*650*ti-pi/4) + sin(2*pi*1500*ti+pi/2);
tend = 1/500*2.4; % czym dluzszy tym lepsze widmo

% build input func.
Tin = 1/fsin;
tin = (0:Tin:tend);
fin = func(tin);

% downsampling
[fout, tout] = downsampling(fin, fsin, fsout);

fout1 = decimate(fin, scale);

% calculate ideal reconstruction
forg = func(tout);

%% compare them

fig = figure;
hold on;
plot(tout, forg, 'r');
plot(tin, fin, 'g');
plot(tout, fout, 'b');
plot(tout, fout1, 'm');
legend('ideal','in', 'out', 'decimate');

%% compare fq

FIN = fftshift(abs(fft(fin)));
FOUT = fftshift(abs(fft(fout)));
FOUT1 = fftshift(abs(fft(fout1)));
s = size(FIN, 2) / size(FOUT, 2);
FORG = s*fftshift(abs(fft(forg)));
FOUT = s*FOUT;
FOUT1 = s*FOUT1;

fqout = fsout/2*linspace(-1, 1, size(FORG, 2));
fqin = (1/Tin)/2*linspace(-1, 1, size(FIN, 2));

hold on;
plot(fqout, FORG, 'r');
plot(fqin, FIN, 'g');
plot(fqout, FOUT, 'b');
plot(fqout, FOUT1, 'm');
legend('ideal','in', 'out', 'decimate');

%% error;

error_calc(fout, forg);
error_calc(fout1, forg);

Błąd mojej wersji jest porównywalny do tej wbudowanej w Matlab. Ilość próbek wejściowych jest równa wielokrotnością skalowania. Jak widać obie wersje rozmijają się, Matlab na początku, ja na końcu. Dla innych liczb próbek obie funkcję potrafią zwrócić bardzo porównywalne wyniki.