Re: Should interp1.m function allow jumps in X-values?

[Ben Abbott]

On Aug 21, 2012, at 10:53 AM, Ben Abbott wrote:

> On Aug 21, 2012, at 10:35 AM, Juan Pablo Carbajal wrote:
> 
>> On Tue, Aug 21, 2012 at 3:26 PM, vivek dogra <address@hidden> wrote:
>> 
>>> On Tue, Aug 21, 2012 at 6:30 PM, Michael D Godfrey
>>> <address@hidden> wrote:
>>>> 
>>>> On 08/21/2012 08:57 AM, vivek dogra wrote:
>>>> 
>>>> Although if it has identical values then it is always favoring increasing
>>>> order
>>>> 
>>>> I think you need to read the code again taking into account:
>>>> 
>>>> octave:1> X = [ 1 2 2 3];
>>>> octave:2> issorted(X)
>>>> ans =  1
>>>> 
>>> 
>>> Sorry i am not able to make myself clear. My understanding of the code is as
>>> follows:
>>> 
>>> ## check whether x is sorted; sort if not.
>>> if (! issorted (x, "either"))
>>> ## Check if X is sorted either in ascending or descending order.
>>>   [x, p] = sort (x);
>>> ##  If it is not sorted, sort it  ("X is sorted in ascending order by
>>> default")
>>>   y = y(p,:);
>>> endif
>>> 
>>> My point:
>>> 1st case: X having distinct values, everything is alright!
>>> 
>>> 2nd case: X is increasing or decreasing (sorted) and X has at most one
>>> identical pair. interp1 gives answer as mentioned in the documentation.
>>> 
>>> 3rd case: X is neither increasing nor decreasing (unsorted) and have
>>> identical values ==> Octave will always sort it in ascending order.
>>> 
>>> In the following example which is 3rd case (X is unsorted and and have
>>> identical values):
>>> 
>>> 
>>> X = [ 2 1 3 2];
>>> Y = [ 2 1 3 10];                               # i get 10 because of some
>>> error in an experiment
>>> 
>>> Now because Octave will always sort it in ascending order, i.e. After
>>> sorting:
>>> 
>>> X = [ 1 2 2 3];
>>> Y = [ 1 2 10 3];
>>> 
>>> Now as X is increasing now it is right continuous therefore to interpolate
>>> Octave will use value of 2 as 10, which is incorrect. (Octave has decided
>>> for the user)
>>> 
>>> Please reply to this situation. I would like to have a discussion on what
>>> will happen in some practical situation like this?
>> 
>> I repeat.
>> If X there are repeated values of X with different values of Y. You
>> can't do interpolation. It is not defined, no continuous function can
>> go through your values. You may want to fix you data or do a
>> regression.
>> 
>> interp1 calculates an interpolant, so it should either ignore the
>> repeated values or give an error. If it does some other calculation,
>> and therefore not returning an interpolant, it should give a warning.
> 
> 
> I agree that linear piece-wise continuous interpolation is not possible for 
> repeated values of X.  However, discontinuous interpolation is not only 
> possible but has is already implemented by interp1.m, and is documented in 
> its doc-string.
> 
> If I understand Vivek's concern, the issue is how to tell interp1() that the 
> interpolant is to produce a left-continuous or right-continuous output.  
> Perhaps, we should support an option? What about "-right" and "-left", with 
> "-right" being the default when X is unsorted?
> 
>       X = [ 2 1 3 2];
>       Y = [ 2 1 3 10]; 
>       x = 1:0.5:3
>       y_right = interp1 (X, Y, x, "-right")
>       y_left = interp1 (X, Y, x, "-left")
> 
> Ben

I've attached a changeset which reverts the recent change, adds the "-left" and 
"-right" options, and modifies the documentation in the hope of clarifying the 
implementation.

I've also added two tests, and a demo, to clarify how the discontinuous 
interpolation was implemented.

This changeset has also been attached to the bug report below.

        https://savannah.gnu.org/bugs/?36732

Vivek/others, does this address your concerns?

Ben

changeset.patch
Description: Binary data