Skip to content

Commit 2f387ff

Browse files
authored
FIX - fixed bug related to amplitude scaling using the drop in upfirdn function (#2087)
1 parent 240fdf3 commit 2f387ff

File tree

3 files changed

+82
-3
lines changed

3 files changed

+82
-3
lines changed

external/signal/resample.m

Lines changed: 7 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -129,8 +129,13 @@
129129
hpad = [hpad; zeros(nz_post,1)];
130130

131131
% filtering
132-
xfilt = upfirdn(x,hpad,p,q);
133-
y = xfilt(offset+1:offset+Ly,:);
132+
if p==1 && q==1
133+
% no resampling needed.
134+
y = x;
135+
else
136+
xfilt = upfirdn(x,hpad,p,q);
137+
y = xfilt(offset+1:offset+Ly,:);
138+
end
134139

135140
if isrowvector
136141
y=y.';

external/signal/upfirdn.m

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -39,5 +39,8 @@
3939
end
4040

4141
yout = upsample(xin,p);
42-
yout = convn(yout, h); % original was filter(h, 1, yout);
42+
yout = convn(yout, h).*p; % original was filter(h, 1, yout);
43+
% the scaling with p is needed as per github issue 2085, causing the output
44+
% to be scaled by the value of p, with this change, the compat/matlab
45+
% versions will give an output that is about equal (scaled with about 0.9993)
4346
yout = downsample(yout,q);

test/test_issue2085.m

Lines changed: 71 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,71 @@
1+
function test_issue2085
2+
3+
% MEM 1gb
4+
% WALLTIME 00:10:00
5+
% DEPENDENCY ft_resampledata ft_specest_irasa
6+
7+
% test function to evaluate the reported scale difference when resampling
8+
% with the drop in compat replacement of the resample function (and lower
9+
% level dependencies).
10+
11+
global ft_default
12+
[ftver, ftpath] = ft_version;
13+
14+
15+
% simulate data
16+
t = (1:1000)/1000;
17+
fn = cumsum(randn(1,length(t)));
18+
in_dat = fn + cos(2*pi*10*t) + cos(2*pi*60*t);
19+
20+
in_dat = in_dat - mean(in_dat);
21+
22+
% resample, using a variety of p/q ratios
23+
p = [1:50 1:50]; q = [10*ones(1,50) 11*ones(1,50)];
24+
25+
restoredefaultpath
26+
addpath(ftpath);
27+
ft_default.toolbox.signal = 'compat';
28+
ft_defaults;
29+
fprintf(sprintf('resampling with: %s\n', which('resample')));
30+
for k = 1:numel(p)
31+
out_ft_rspu{k,1} = resample(in_dat, p(k), q(k));
32+
out_ft_rspd{k,1} = resample(in_dat, q(k), p(k));
33+
end
34+
35+
restoredefaultpath
36+
addpath(ftpath);
37+
ft_default.toolbox.signal = 'matlab';
38+
ft_defaults;
39+
fprintf(sprintf('resampling with: %s\n', which('resample')));
40+
for k = 1:numel(p)
41+
out_mat_rspu{k,1} = resample(in_dat, p(k), q(k));
42+
out_mat_rspd{k,1} = resample(in_dat, q(k), p(k));
43+
end
44+
45+
for k = 1:numel(p)
46+
bd(k,1) = out_ft_rspd{k}/out_mat_rspd{k};
47+
bu(k,1) = out_ft_rspu{k}/out_mat_rspu{k};
48+
end
49+
figure;plot(bu);hold on;plot(bd);
50+
51+
assert(all(bu>0.999) && all(bd>0.999));
52+
% there's a tiny scale difference between the legacy matlab and the compat
53+
% version, which could be titrated to be closer to 1, but that would not be
54+
% based on a mathematically informed choice, so I think we may accept this <0.1% amplitude diff.
55+
56+
57+
% figure(); hold on;
58+
% plot(in_dat, '-b');
59+
% plot(out_ft_rspu, '-r');
60+
% plot(out_mat_rspu, '-g');
61+
% legend({'original data', 'ft resampled', 'mat resampled'});
62+
% title('upsample')
63+
% hold off;
64+
%
65+
% figure(); hold on;
66+
% plot(in_dat, '-b');
67+
% plot(out_ft_rspd, '-r');
68+
% plot(out_mat_rspd, '-g');
69+
% legend({'original data', 'ft resampled', 'mat resampled'});
70+
% title('downsample');
71+
% hold off;

0 commit comments

Comments
 (0)