From bug-request at octave dot org Thu Jan  5 22:47:08 2006
Subject: Re: interp1.m (octave-forge) *linear extrapolation
From: Paul Kienzle <pkienzle at users dot sourceforge dot net>
To: Olaf Till <olaf dot till at uni-jena dot de>
Cc: bug at octave dot org
Date: Thu, 5 Jan 2006 23:45:07 -0500

Please provide an example where it fails.

Thanks,

	- Paul
On Jan 5, 2006, at 4:57 AM, Olaf Till wrote:

> Sorry, but there is a further issue:
>
> in the same lines of code cited below, line 134
>
> idx = floor(t);
>
> seems not to allow extrapolation (I get 'invalid row index 0' from
> octave when it evaluates line 139.
>
> Adding (in line 135)
>
> idx = max (idx, 1); idx = min (idx, nx);
>
> seems to solve this problem for me.
>
> But a remark: for my data, then, '*linear' seems to last almost twice
> as long as 'linear'.
>
> Olaf
>
> On Thu, Jan 05, 2006 at 09:02:41AM +0000, olaf wrote:
>> This is a bug in a component of octave-forge, please forward if this
>> is not the correct mailing list. If someone wants to respond --- I am
>> not subscribed to bug at octave dot org dot 
>>
>> Function: interp1.m
>> Author: Paul Kienzle
>> Version: Date of last change documented in the file: 2002-01-23
>>
>> Description:
>> In method '*linear' the results are different from method 'linear'
>> (despite of even spacing of x) when the difference between neighboured
>> x's is not equal to 1.
>>
>> Suggested fix:
>>
>> in lines 131 to 139:
>>
>>   elseif strcmp(method, '*linear')
>>     ## find the interval containing the test point
>>     t = (xi - minx)/dx + 1;
>>     idx = floor(t);
>>
>>     ## use the endpoints of the interval to define a line
>>     dy = [y(2:ny,:) - y(1:ny-1,:); zeros(1,nc)];
>>     s = (t - idx)./dx;
>>     yi(range,:) = s(:,ones(1,nc)).*dy(idx,:) + y(idx,:);
>>
>> there is a doubled division by `dx' (in line 133 and in line 138).
>> Changing line 138
>>
>>     s = (t - idx)./dx;
>>
>> to
>>
>>     s = t - idx;
>>
>> makes method '*linear' return the same results as 'linear' even with
>> non integer x's for me.
>>
>> Regards, Olaf
>
>
>
> -------------------------------------------------------------
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web:  http://www.octave.org
> How to fund new projects:  http://www.octave.org/funding.html
> Subscription information:  http://www.octave.org/archive.html
> -------------------------------------------------------------
>



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------