# bug #41838: interp2 needs input validation for strictly monotonic property

 Submitted by: Muhali Submitted on: Tue 11 Mar 2014 08:43:51 PM UTC Category: Octave Function Severity: 3 - Normal Priority: 5 - Normal Item Group: Segfault, Bus Error, etc. Status: Fixed Assigned to: None Originator Name: Open/Closed: Closed Release: dev Operating System: GNU/Linux

You are not logged in

## Discussion

Wed 09 Apr 2014 03:03:49 AM UTC, comment #11:

I added a test for strictly monotonic x and y to the default branch in this changeset (http://hg.savannah.gnu.org/hgweb/octave/rev/4792a115c735).  Closing report.

Rik <rik5>
Mon 31 Mar 2014 03:05:37 PM UTC, comment #10:

no, the one that you get from ordering y and yi increasingly, which is actually only a different way of axis labeling, isn't it?

Muhali <muhali>
Mon 31 Mar 2014 02:49:50 PM UTC, comment #9:

Was the same result an error message about the data not being monotonic?

Rik <rik5>
Mon 31 Mar 2014 02:40:06 PM UTC, comment #8:

as ML wants to have at least 3 points I did

Z = randn(3,3);
interp2(1:3, 3:-1:1, Z, (1:0.5:3)', 3:-0.5:1, 'cubic')

[X, Y] = meshgrid ([1 2 3], [3 2 1]);
interp2(X, Y, Z, (1:0.5:3)', 3:-0.5:1, 'cubic')

and obtained in both cases the exact same result.

Muhali <muhali>
Mon 31 Mar 2014 01:56:51 PM UTC, comment #7:

Matlab uses 'cubic' for an interpolation method with C1 continuity rather than 'pchip'.

Could you try changing the method and repeating the examples below.

interp2(1:2, 2:-1:1, randn(2,2), (1:0.5:2)', 2:-0.5:1, 'cubic')

[X, Y] = meshgrid ([1 2], [2 1]);
interp2(X, Y, 2:-1:1, randn(2,2), (1:0.5:2)', 2:-0.5:1, 'cubic')

I think the solution is going to come back to doing better input validation as is done in interp1.  I overhauled interp1, interp2, interp3 on the development branch in this changeset (http://hg.savannah.gnu.org/hgweb/octave/rev/0ede4dbb37f1).  Any improvements should start after that changeset so that the patches apply cleanly.

Rik <rik5>
Mon 31 Mar 2014 06:48:06 AM UTC, comment #6:

There is no 2d pchip method in Matlab.

The error does not occur for interp1 because there are at least two checks/repairs for decreasing data:

## interp1 line 184
## check whether x is sorted; sort if not.
if (! issorted (x, "either"))
[x, p] = sort (x);
y = y(p,:);
endif

## pchip line 93
h = diff (x);
if (all (h < 0))
x = fliplr (x);
h = diff (x);
y = fliplr (y);
elseif (any (h <= 0))
error ("pchip: X must be strictly monotonic");
endif

Without that, the slatec call to _pchip_deriv_ returns the error seen in the segfault message.

Muhali <muhali>
Sat 29 Mar 2014 10:58:40 PM UTC, comment #5:

The test for strict monotonic is there, but written in fortran. If it is added to C++ code it could be removed in fortran. I don't have Matlab to test, and pchip is not a valid method according to online documentation. However, we had a discussion on IRC yesterday about this issue. Seems Matlab 2014a returns some results.

Stefan Mahr <dac922>
Sat 29 Mar 2014 07:14:44 PM UTC, comment #4:

I applied Stefan's patch to the stable branch here (http://hg.savannah.gnu.org/hgweb/octave/rev/e7369c43ebe0).  This at least stops the segmentation fault.

Back to the root problem, there should be input validation of X and Y.  According to the documentation, when X and Y are grid vectors (rather than matrices from meshgrid()), the must be strictly monotonic, but there is no test for that.  We could add one though.

What does Matlab do for these two test cases?

interp2(1:2, 2:-1:1, randn(2,2), (1:0.5:2)', 2:-0.5:1, 'pchip')

[X, Y] = meshgrid ([1 2], [2 1]);
interp2(X, Y, 2:-1:1, randn(2,2), (1:0.5:2)', 2:-0.5:1, 'pchip')

Rik <rik5>
Sun 16 Mar 2014 07:08:39 AM UTC, comment #3:

I found the error is due to this line in file interp2.m:
DY = _pchip_deriv_ (Y, Z, 1);
and in the input Y is decreasing.
So I think it's the problem with function _pchip_deriv_
But I don't know how to find the definition of this function.

Ma Jun <yefllower>
Wed 12 Mar 2014 01:54:27 AM UTC, comment #2:

Attached patch fixes this bug.

Stefan Mahr <dac922>
Tue 11 Mar 2014 09:54:31 PM UTC, comment #1:

Confirmed on Windows XP (3.8.1) and Linux (Mageia-2 32b; 4.1.0+ and 3.8.0).

Tried a backtrace from 3.8.0, but got, apart from Muhali's uppercase messages, only:

Program received signal SIGSEGV, Segmentation fault.
#0 0x8b93?? ()
Backtrace stopped: Not enough registers or memory available to unwind further

Philip Nienhuis <philipnienhuis>
Tue 11 Mar 2014 08:43:51 PM UTC, original submission:

interp2(1:2, 2:-1:1, randn(2,2), (1:0.5:2)', 2:-0.5:1, 'pchip')

causes a crash in my current octave-dev:

***MESSAGE FROM ROUTINE DPCHIM IN LIBRARY SLATEC.
***POTENTIALLY RECOVERABLE ERROR, PROG ABORTED, TRACEBACK REQUESTED

• X-ARRAY NOT STRICTLY INCREASING
• ERROR NUMBER = -3
•

***END OF MESSAGE

***JOB ABORT DUE TO UNRECOVERED ERROR.

1.          ERROR MESSAGE SUMMARY

LIBRARY    SUBROUTINE MESSAGE START             NERR     LEVEL     COUNT
SLATEC     DPCHIM     X-ARRAY NOT STRICTLY        -3         1         1

Segmentation fault (core dumped)

It happens each time one of the arguments is decreasing (consistent with the interp2 help which only requires them to be monotonic).

Muhali <muhali>

## Attached Files

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:

Comment:

Attached Files
file #30869:  interp2_pchip.diff added by dac922 (2KiB - text/x-patch)

## Dependencies

Depends on the following items: None found

Items that depend on this one: None found

Carbon-Copy List
• -email is unavailable- added by rik5 (Posted a comment)
• -email is unavailable- added by mtmiller (Updated the item)
• -email is unavailable- added by yefllower (Posted a comment)
• -email is unavailable- added by dac922 (Updated the item)
• -email is unavailable- added by muhali (Submitted the item)
•

Do you think this task is very important?
This task has 0 encouragements so far.

Only project members can vote.