special kind of interpolation

1 view (last 30 days)
Ingrid
Ingrid on 27 May 2015
Edited: Ingrid on 27 May 2015
I have some data shown in the example figure below for which I want to perform an interpolation to calculate the average value. However, I need a linear interpolation between the datapoints, except for points that are zero, there are values should be put at zero between the two surrounding points. How can I achieve this?

Accepted Answer

Thorsten
Thorsten on 27 May 2015
Edited: Thorsten on 27 May 2015
This "interpolation" is basically a plot with values added where y is zero:
x = x(:)'; y = y(:)'; % force row vectors
ind = find(diff(y == 0) == -1); % find the index where y changes from zero to non-zero
% add points before and after ind with zero y value
for i = 1:numel(ind)
indi = ind(i);
x = [x(1:indi-1) x(indi-1) x(indi) x(indi+1) x(indi+1:end)];
y = [y(1:indi-1) 0 y(indi) 0 y(indi+1:end)];
end
plot(x, y) % plot the curve
  1 Comment
Ingrid
Ingrid on 27 May 2015
this looks pretty good, but it is not the plot that I am interested in, I want to calculate the mean by using the interpolated values (I just use a x-vector with spacing 1 to do the interpolation). I thought it was not possible to use the interp1 function with two x-values equal to each other but with different y-values but I guess it is possible as my tests gives values as expected
I also have to this more than a million times so speed is an important factor so therefore I was looking for a vectorized form but I guess this comes pretty close so thanks

Sign in to comment.

More Answers (1)

Walter Roberson
Walter Roberson on 27 May 2015
Search the data for 0's. Suppose you find it at index J. Then introduce a new timepoint just after t(J-1), at time t(J-1)*(1+eps), with value 0, and introduce another just before t(J+1), at time t(J+1)*(1-eps), again with value 0. Now that I think of it, you might as well change the time t(J) to be t(J-1)*(1+eps) to move the 0 to right beside the previous point, and then you would only need to insert one new point.
Doing the stitching together in vectorized form could be a bit tricky. Mumble...
tidx = sort([(1:length(t)).', find(datapoints(:)==0)]);
newt = t(tidx);
newdata = datapoints(tidx);
new_at = find(diff(tidx)==0);
newt(new_at) = newt(new_at-1)*(1+eps);
newt(new_at+1) = newt(new_at+2)*(1-eps);
newdata(new_at) = 0;
newdata(new_at+1) = 0;
There, that should be pretty close.
  1 Comment
Ingrid
Ingrid on 27 May 2015
Edited: Ingrid on 27 May 2015
indeed I was having problems with vectorizing my code as with a for-loop it would be straightforward so thank you for your contribution. However, Thorsten's code was approximately two times faster which in my case where I have to repeat this a million times, is the crucial factor, I have accepted his answer
Elapsed time is 0.000208 seconds. -> walter
Elapsed time is 0.000109 seconds. -> thorsten
oh, and your code works and gives the same final answer, AFTER correction of the typo (just in case someone else want to try it):
tidx = sort([(1:length(t)).'; find(datapoints(:)==0)]);

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!