Hi,
The problem occurs when there are any repeated values in X, GSL will
treat it as non-monotonically increasing data and cause GSL_ERROR
calling abort().
In GSL-1.12/interpolation/interp.c, line 79:
for (i = 1; i < size; i++)
{
if (!(x_array[i-1] < x_array[i]))
{
GSL_ERROR ("x values must be monotonically increasing",
GSL_EINVAL);
}
}
When there are any repeated values in X, x_array[i-1] will *not* be
smaller than x_array[i].
Regards,
Jason
On 10/27/2010 07:40 AM, Andres Jordan wrote:
> Hi Doug,
>
> That should not happen, the init method in GSL::INTERP is supposed to
> sort X by default (take a look at gsl_interp.pd). The error you
> mention should show up only if you pass a non-monotonically increasing
> set of X values *and* specify {Sort => 0} as an option. So I'm
> confused by the reported behaviour...
>
> -Andres
>
> ps. the relevant bit of code in gsl_interp.pd is
>
> if($$opt{Sort} != 0){
> my $idx = PDL::Ufunc::qsorti($x);
> $x = $x->index($idx);
> $y = $y->index($idx);
> }
>
>
>
> On Tue, Oct 26, 2010 at 1:46 PM, Doug Hunt<[email protected]> wrote:
>> Hi all: In the default interpolation interface to GSL, if one uses a set of
>> X values that are not monotonically increasing, GSL fails hard:
>>
>> gsl: interp.c:83: ERROR: x values must be monotonically increasing
>> Default GSL error handler invoked.
>>
>> this error cannot be caught by an exception handler, so is troublesome in
>> production code.
>>
>> My collegue Jason Lin has found a simple patch which fixes this problem
>> (attached). This patch just turns off the default GSL error handler which
>> just calls 'stop'.
>>
>> I've downloaded the latest GIT sources and applied it. It seems to work.
>> Any objections if I apply this patch? Here is a test case for it:
>>
>> --------------------------------------------------------------------------
>> my $nx = ($x)*($x<=3) + ($x-1)*($x>3); # x value not monotonically
>> increasing
>> my $i; eval { $i = PDL::GSL::INTERP->init('cspline',$nx, $y) };
>> like($@,qr/invalid argument supplied by user/,"module exception handling");
>> --------------------------------------------------------------------------
>>
>> I'm reluctant to add this to gsl_interp.t because it will cause the test
>> suite to fail hard if the patch has not been applied.
>>
>> If there are no objections, I'll apply this patch in a couple of days.
>>
>> Thanks,
>>
>> Doug
>>
>>
>> [email protected]
>> Software Engineer
>> UCAR - COSMIC, Tel. (303) 497-2611
>>
>>
>> _______________________________________________
>> Perldl mailing list
>> [email protected]
>> http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
>>
>>
--
Jason Lin
Software Engineer
COSMIC, UCAR
Tel. 303-497-2608
_______________________________________________
Perldl mailing list
[email protected]
http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
