# Resample function is not working properly and damaging the signal

62 views (last 30 days)

Show older comments

Hi,

I have a question regarding the resample function (It is the documentation under "Resample Multichannel Signal"). In this function, I realized that the function damages the signal on both ends of the given signal. Here I am attaching the code to generate 100-sample sinusoidal signal.

p = 3;

q = 2;

tx = 0:p:300-p;

x = cos(2*pi*tx./(1:5)'/100);

plot(tx,x,'.:')

title('Original')

ylim([-1.5 1.5])

However, once you resample with the following give code, the signals are simply damaged at both ends. You can see the effects on both time points close to 0 and close to 300.

ty = 0:q:300-q;

y = resample(x,p,q,'Dimension',2);

plot(ty,y,'.:')

title('Upsampled')

Therefore, could you please help me to solve this issue? Because, I am currently trying to resample my multichannel motion capture signals, which are sampled at 60 Hz or 200 Hz, in "100 Hz". However, this function will damage my signal if I use it, and this will definitely affect the outcome of my research.

I will be waiting for your help and appreciate your assistance.

Thanks!

##### 6 Comments

Paul
on 17 Mar 2023

Bob S
on 22 Sep 2023

Edited: Bob S
on 22 Sep 2023

Hi,

I do have an additional comment regarding the post by @goksu avdan. There is an issue with using resample() for processing blocks of data. If I'm using resample to interpolate/decimate, I will get erroneous results at the edges. In the past, to get around this, I have used a "brute-force" approach with filter() as in the following:

while dataAvailable

% Say x is dimensioned 1 by Nsamples and I am upsampling by a factor

% Nup and downsampling by Ndown

% Zero fill between samples

xi = [x;zeros(Nup-1,Nsamples)];

xi = xi(:).';

% xi should now be dimensioned 1 by (Nup*Nsamples) with zeros inbetween the original samples

% Now low-pass filter and decimate

[xd, zFilt] = filter(bFilt, 1, xi, zFilt);

xd = xd(1:Ndown:end);

end

Assuming that I use the same bFilt that would come out of resamp(), the above result will match resample() if I ignore the first floor(length(bFilt)/2/Ndown) output samples. By using the above approach, I have to deal with a transient response for the first block, after that I have continuous correct results. In order to achieve that with resample, one would have to devise some kind of overlap-save approach.

The problem with the "brute force" approach is that it is inefficient and does not take advantage of all of the zeros that a proper polyphase filtering technique, as done in resample(), would handle.

What would be ideal would be something like:

[y, zFilt] = resmaple(x, p, q, zFilt);

I do realize that depending on what the ratio of p/q looks like, it might be near impossible to get an integer number of output points without "playing tricks", but I'm ignoring that right now.

Any thoughts on this? Am I missing something? Apologies for spelling errors and any code miscues.

Thanks in advance!

### Accepted Answer

Star Strider
on 16 Mar 2023

I had not noticed the ‘end effect’ problem before (although I use resample relatively frequently, and that is clearly visible in the signals in the resample documentation). Looking at the original resampled waveforms, it occurred to me that the waveforms that were zero at both ends didn’t have that probllem. I devised a single-line linear regression for another problem recently, and decided to use it here to first ‘detrend’ the signals so that they began and ended at zero, resample them, and then ‘retrend’ them to their original orientations.

This may not be perfect, since it uses resample on the detrended signals, however it does elliminate the ‘end effect’ problem.

I am not certain what other processing you’re doing, so experiment with it to see how it works with your signals in that context —

p = 3;

q = 2;

tx = 0:p:300-p;

x = cos(2*pi*tx./(1:5)'/100);

figure

plot(tx,x,'.:')

title('Original')

ylim([-1.5 1.5])

legend(compose('%d',1:5), 'Location','best')

ty = 0:q:300-q;

% y = resample(x,p,q,'Dimension',2);

Fs = 1/(ty(2)-ty(1));

LR = @(X,Y,x,y) [X(:) ones(size(X(:)))] * ([x(:) ones(size(x(:)))] \ y(:)); % Single-Line Linear Regression

for k = 1:size(x,1)

DeTr = LR(tx,x(k,:),tx([1 end]),x(k,[1 end])).';

[y(k,:),ty] = resample(x(k,:)-DeTr,tx,Fs,'Dimension',2);

ReTr = LR(ty,y(k,:),tx([end 1]),x(k,[end 1])).';

y(k,:) = y(k,:)+ReTr;

end

figure

plot(ty,y,'.:')

title('Upsampled')

legend(compose('%d',1:5), 'Location','best')

.

##### 14 Comments

Stephen23
on 18 Mar 2023

Edited: Stephen23
on 18 Mar 2023

"Using the inbuilt INTERP1 is even simpler than downloading a third-party function."

I did not want to imply that INTERP1 replaces RESAMPLE, only that INTERP1 is a simpler replacement for Star Strider's detrend and retrend function, exactly as I showed here:

I realize that my comment was ambiguous on this.

### More Answers (2)

Paul
on 17 Mar 2023

Edited: Paul
on 17 Mar 2023

Hi goksu,

The observed behavior is documented at resample specifcally here: "... resample assumes that the input sequence, x, is zero before and after the samples it is given. Large deviations from zero at the endpoints of x can result in unexpected values for y." See examples at that doc page.

If you're looking for theoretically pleasing approach, resample provides means to influence transients as shown here, for example. Don't know if that approach applies to the ends of the signal or just transients internal to the signal as shown in that example.

If you're looking for a heuristic approach, one option might be to extrapolate the signal off the ends, apply resample, then chop of the extensions. I used 20 sample extensions and linear extrapolation on the example data. The actual number and method will likely depend on the actual data and the resample filter parameters. I'm only thinking about extrapolations that minimize transients at the transitions to/from the endpoints of the original signal before resampling. spline or pchip might work better. Anyway, here's what I tried:

Original signals

p = 3;

q = 2;

tx = 0:p:300-p;

x = cos(2*pi*tx./(1:5)'/100);

plot(tx,x,'.:')

title('Original')

ylim([-1.5 1.5])

Extend the signals by 20 samples via linear extrapolation.

txext = [(-20:-1)*p tx (1:20)*p+tx(end)];

xext = (interp1(tx,x.',txext,'linear','extrap')).';

Resample

y = resample(xext,p,q,'Dimension',2);

fsx = 1/p;

fsy = fsx*p/q;

ty = txext(1) + (0:(size(y,2)-1))/fsy;

Plot the extended output, note the transients have been pushed away from the central region we care about 0 < t < 300.

figure;

plot(ty,y),grid

Lop off the extensions

y(:,ty<0) = [];

ty(ty<0) = [];

y(:,ty>tx(end)) = [];

ty(ty>tx(end)) = [];

figure

plot(ty,y),grid

##### 4 Comments

Paul
on 17 Mar 2023

Stephen23
on 18 Mar 2023

Edited: Stephen23
on 18 Mar 2023

It is easy to avoid the loop too, it is just a simple linear interpolation:

p = 3;

q = 2;

tx = (0:p:300-p).';

x = cos(2*pi*tx./(1:5)/100); % note orientation!

plot(tx,x,'.:')

title('Original')

ylim([-1.5,1.5])

ty = 0:q:300-q;

ax = interp1(tx([1,end]),x([1,end],:),tx);

ay = interp1(ty([1,end]),x([1,end],:),ty);

y = resample(x-ax,p,q)+ay;

plot(ty,y,'.:')

title('Upsampled')

##### 2 Comments

Stephen23
on 18 Mar 2023

Edited: Stephen23
on 18 Mar 2023

Note that you can easily avoid all of those loops and all of the superfluous data processsing.

Avoid fighting MATLAB with lots of loops and IFs and the like. Think in terms of arrays and indices.

S = load('Subject3.mat');

S = S.s;

% Select valid dataset:

X = strcmp({S.Data.Task},'Walking');

Y = cellfun(@(m)any(isnan(m(:))),{S.Data.Marker});

Z = find(X&~Y,1,'last');

old = S.Data(Z).Marker;

% Sample rate:

p = interp1([60,200],[5,1],S.KinFreq);

q = interp1([60,200],[3,2],S.KinFreq);

tx = p*(0:size(old,2)-1);

Fs = 1/q;

% Resample:

DeTr = interp1(tx([1,end]),old(:,[1,end]).',tx).';

[M,ty] = resample(old-DeTr,tx,Fs,'Dimension',2);

ReTr = interp1(ty([1,end]),old(:,[1,end]).',ty).';

new = M+ReTr;

% Plot original:

figure

plot(tx, old)

grid

xlabel('Time')

ylabel('Signal')

title('Original Vectors')

% Plot resampled:

figure

plot(ty, new)

grid

xlabel('Time')

ylabel('Signal')

title('Resampled Vectors')

### See Also

### Categories

### Community Treasure Hunt

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

Start Hunting!