——————————=head1 NAME
PDL::Interpolate - provide a consistent interface to the interpolation routines available in PDL
=head1 SYNOPSIS
use PDL::Interpolate;
my $i = new PDL::Interpolate( x => $x, y = $y );
my $y = $i->interpolate( $xi );
=head1 DESCRIPTION
This module aims to provide a relatively-uniform interface
to the various interpolation methods available to PDL.
The idea is that a different interpolation scheme
can be used just by changing the C<new> call.
At present, PDL::Interpolate itself just provides
a somewhat-convoluted interface to the C<interpolate>
function of L<PDL::Primitive|PDL::Primitive/interpolate>.
However, it is expected that derived classes,
such as
L<PDL::Interpolate::Slatec|PDL::Interpolate::Slatec>,
will actually be used in real-world situations.
To use, create a PDL::Interpolate (or a derived class)
object, supplying it with its required attributes.
=head1 LIBRARIES
Currently, the avaliable classes are
=over 4
=item PDL::Interpolate
Provides an interface to the interpolation routines of PDL.
At present this is the linear interpolation routine
L<PDL::Primitive::interpol|PDL::Primitive/interpol>.
=item PDL::Interpolate::Slatec
The SLATEC library contains several approaches to interpolation:
piecewise cubic Hermite functions and B-splines.
At present, only the former method is available.
=back
It should be relatively easy to provide an interface to other
interpolation routines, such as those provided by the
Gnu Scientific Library (GSL).
=head1 ATTRIBUTES
The attributes (or options) of an object are as follows;
derived classes may modify this list.
Attribute Flag Description
x sgr x positions of data
y sgr function values at x positions
bc g boundary conditions
err g error flag
type g type of interpolation
A flag of C<s> means that a user can set this attribute
with the L<new|/new> or L<set|/set> methods,
a flag of C<g> means that the user can obtain the
value of this attribute using L<get|/get>,
and a flag of C<r> means that the attribute is required
when an object is created (see the L<new|/new> method).
Attribute Default value
bc "none"
type "linear"
If a routine is sent an attribute it does not understand, then
it ignores that attribute, except for L<get|/get>, which
returns C<undef> for that value.
=head1 METHODS
The default methods are described below. However, defined classes
may extend them as they see fit, and add new methods.
Throughout this documentation, C<$x> and C<$y> refer to the function
to be interpolated whilst C<$xi> and C<$yi> are the interpolated values.
=head1 THREADING
The class will thread properly if the routines it calls do so.
See the SYNOPSIS section of L<PDL::Interpolate::Slatec>
(if available) for an example.
=cut
package
PDL::Interpolate;
use
Carp;
use
strict;
####################################################################
## Public routines:
=head2 new
=for usage
$obj = new PDL::Interpolate( x => $x, y => $y );
=for ref
Create a PDL::Interpolate object.
The required L<attributes|/attributes> are
C<x> and C<y>.
At present the only available interpolation method
is C<"linear"> - which just uses
L<PDL::Primitive::interpolate|PDL::Primitive::interpolate> - and
there are no options for boundary conditions, which is why
the C<type> and C<bc> attributes can not be changed.
=cut
# meaning of types:
# required - required, if this attr is changed, we need to re-initialise
# settable - can be changed with a new() or set() command
# gettable - can be read with a get() command
#
sub
new {
my
$this
=
shift
;
my
$class
=
ref
(
$this
) ||
$this
;
# class structure
my
$self
= {
attributes
=> {},
values
=> {},
types
=> {
required
=> 0,
settable
=> 0,
gettable
=> 0 },
flags
=> {
library
=>
"PDL"
,
status
=> 1,
routine
=>
"none"
,
changed
=> 1 },
};
# make $self into an object
bless
$self
,
$class
;
# set up default attributes
#
$self
->_add_attr(
x
=> {
required
=> 1,
settable
=> 1,
gettable
=> 1 },
y
=> {
required
=> 1,
settable
=> 1,
gettable
=> 1 },
bc
=> {
gettable
=> 1 },
err
=> {
gettable
=> 1 },
type
=> {
gettable
=> 1 },
);
$self
->_set_value(
bc
=>
"none"
,
type
=>
"linear"
,
);
# set variables
# - expect sub-classes to call this new with no variables, so $#_ == -1
$self
->set(
@_
)
if
(
@_
);
# return the object
return
$self
;
}
# sub: new()
#####################################################################
# in _add_attr(), _change_attr() and _add_attr_type()
# we set flags->changed to 1 when something changes. It's
# a bit over the top to do this, as these should only be called when
# creating the object, when the changed flag should be set to 1 anyway
# add attributes to the object and sets value to undef
#
# supply a hash array, keys == variable name,
# values are a hash array with keys matching
# $self->{values}, which also gives the default value
# for the type
#
# this can only be used to create an attribute -
# see _change_attr() to change an already exsiting attribute.
#
# the fields are set to the default values, then filled in with the supplied values
# any value that is unknown is ignored
#
sub
_add_attr {
my
$self
=
shift
;
my
%attrs
= (
@_
);
foreach
my
$attr
(
keys
%attrs
) {
croak
"ERROR: adding an attribute ($attr) which is already known.\n"
if
defined
$self
->{attributes}->{
$attr
};
# set default values
foreach
my
$type
(
keys
%{
$self
->{types}} ) {
$self
->{attributes}->{
$attr
}->{
$type
} =
$self
->{types}->{
$type
};
}
# change the values to those supplied
foreach
my
$type
(
keys
%{
$attrs
{
$attr
}} ) {
$self
->{attributes}->{
$attr
}->{
$type
} =
$attrs
{
$attr
}->{
$type
}
if
exists
$self
->{types}->{
$type
};
}
# set value to undef
$self
->{
values
}->{
$attr
} =
undef
;
}
$self
->{flags}->{changed} = 1;
}
# sub: _add_attr()
# changes attributes of the object
#
# the given attributes MUST already exist
#
sub
_change_attr {
my
$self
=
shift
;
my
%attrs
= (
@_
);
foreach
my
$attr
(
keys
%attrs
) {
croak
"ERROR: changing an attribute ($attr) which is not known.\n"
unless
defined
$self
->{attributes}->{
$attr
};
# change the values to those supplied
foreach
my
$type
(
keys
%{
$attrs
{
$attr
}} ) {
if
(
exists
$self
->{types}->{
$type
} ) {
$self
->{attributes}->{
$attr
}->{
$type
} =
$attrs
{
$attr
}->{
$type
};
$self
->{flags}->{changed} = 1;
}
}
}
}
# sub: _change_attr()
# adds the given types to the allowed list, and
# updates all attributes to contain the default value
#
# Useful for sub-classes which add new types
#
sub
_add_attr_type {
my
$self
=
shift
;
my
%types
= (
@_
);
foreach
my
$type
(
keys
%types
) {
croak
"ERROR: adding type ($type) that is already known.\n"
if
exists
$self
->{types}->{
$type
};
$self
->{types}->{
$type
} =
$types
{
$type
};
# loop through each attribute, adding this type
foreach
my
$attr
(
keys
%{
$self
->{attributes}} ) {
$self
->{attributes}->{
$attr
}->{
$type
} =
$types
{
$type
};
}
$self
->{flags}->{changed} = 1;
}
}
# sub: _add_attr_type()
####################################################################
# if an attribute has changed, check all required attributes
# still exist and re-initialise the object (for PDL::Interpolate
# this is a nop)
#
sub
_check_attr {
my
$self
=
shift
;
return
unless
$self
->{flags}->{changed};
my
@emsg
;
foreach
my
$name
(
keys
%{
$self
->{attributes} } ) {
if
(
$self
->{attributes}->{
$name
}->{required} ) {
push
@emsg
,
$name
unless
defined
(
$self
->{
values
}->{
$name
});
}
}
croak
"ERROR - the following attributes must be supplied:\n [ @emsg ]\n"
unless
$#emsg
== -1;
$self
->{flags}->{routine} =
"none"
;
$self
->{flags}->{status} = 1;
$self
->_initialise;
$self
->{flags}->{new} = 0;
}
# sub: check_attr()
####################################################################
#
# method to be over-ridden by sub-classes
# PDL::Interpolate needs no initialisation
#
sub
_initialise {}
####################################################################
# a version of set that ignores the settable flag
# - for use by the class, not by the public
#
# it still ignores unknown attributes
#
sub
_set_value {
my
$self
=
shift
;
my
%attrs
= (
@_
);
foreach
my
$attr
(
keys
%attrs
) {
if
(
exists
(
$self
->{
values
}->{
$attr
}) ) {
$self
->{
values
}->{
$attr
} =
$attrs
{
$attr
};
$self
->{flags}->{changed} = 1;
}
}
}
# sub: _set_value()
# a version of get that ignores the gettable flag
# - for use by the class, not by the public
#
# an unknown attribute returns an undef
#
sub
_get_value {
my
$self
=
shift
;
my
@ret
;
foreach
my
$name
(
@_
) {
if
(
exists
$self
->{
values
}->{
$name
} ) {
push
@ret
,
$self
->{
values
}->{
$name
};
}
else
{
push
@ret
,
undef
;
}
}
return
wantarray
?
@ret
:
$ret
[0];
}
# sub: _get_value()
####################################################################
=head2 set
=for usage
my $nset = $obj->set( x = $newx, $y => $newy );
=for ref
Set attributes for a PDL::Interpolate object.
The return value gives the number of the supplied attributes
which were actually set.
=cut
sub
set {
my
$self
=
shift
;
my
%vals
= (
@_
);
my
$ctr
= 0;
foreach
my
$name
(
keys
%vals
) {
if
(
exists
$self
->{attributes}->{
$name
}->{settable} ) {
$self
->{
values
}->{
$name
} =
$vals
{
$name
};
$ctr
++;
}
}
$self
->{flags}->{changed} = 1
if
$ctr
;
return
$ctr
;
}
# sub: set()
####################################################################
=head2 get
=for usage
my $x = $obj->get( x );
my ( $x, $y ) = $obj->get( qw( x y ) );
=for ref
Get attributes from a PDL::Interpolate object.
Given a list of attribute names, return a list of
their values; in scalar mode return a scalar value.
If the supplied list contains an unknown attribute,
C<get> returns a value of C<undef> for that
attribute.
=cut
sub
get {
my
$self
=
shift
;
my
@ret
;
foreach
my
$name
(
@_
) {
if
(
exists
$self
->{attributes}->{
$name
}->{gettable} ) {
push
@ret
,
$self
->{
values
}->{
$name
};
}
else
{
push
@ret
,
undef
;
}
}
return
wantarray
?
@ret
:
$ret
[0];
}
# sub: get()
####################################################################
=head2 interpolate
=for usage
my $yi = $obj->interpolate( $xi );
=for ref
Returns the interpolated function at a given set of points.
A status value of -1, as returned by the C<status> method,
means that some of the C<$xi> points lay outside the
range of the data. The values for these points
were calculated using linear extrapolation.
=cut
sub
interpolate {
my
$self
=
shift
;
my
$xi
=
shift
;
croak
'Usage: $obj->interpolate( $xi )'
.
"\n"
unless
defined
$xi
;
# check everything is fine
$self
->_check_attr();
# get values in one go
my
(
$x
,
$y
) =
$self
->_get_value(
qw( x y )
);
my
(
$yi
,
$err
) = PDL::Primitive::interpolate(
$xi
,
$x
,
$y
);
if
( any
$err
) {
$self
->{flags}->{status} = -1;
}
else
{
$self
->{flags}->{status} = 1;
}
$self
->_set_value(
err
=>
$err
);
$self
->{flags}->{routine} =
"interpolate"
;
return
$yi
;
}
####################################################################
#
# access to flags - have individual methods for these
=head2 status
=for usage
my $status = $obj->status;
=for ref
Returns the status of a PDL::Interpolate object
Returns B<1> if everything is okay, B<0> if
there has been a serious error since the last time
C<status> was called, and B<-1> if there
was a problem which was not serious.
In the latter case, C<$obj-E<gt>get("err")> may
provide more information, depending on the
particular class.
=cut
sub
status {
my
$self
=
shift
;
return
$self
->{flags}->{status}; }
=head2 library
=for usage
my $name = $obj->library;
=for ref
Returns the name of the library used by a PDL::Interpolate object
For PDL::Interpolate, the library name is C<"PDL">.
=cut
sub
library {
my
$self
=
shift
;
return
$self
->{flags}->{library}; }
=head2 routine
=for usage
my $name = $obj->routine;
=for ref
Returns the name of the last routine called by a PDL::Interpolate object.
For PDL::Interpolate, the only routine used is C<"interpolate">.
This will be more useful when calling derived classes,
in particular when trying to decode the values stored in the
C<err> attribute.
=cut
sub
routine {
my
$self
=
shift
;
return
$self
->{flags}->{routine}; }
=head2 attributes
=for usage
$obj->attributes;
PDL::Interpolate::attributes;
=for ref
Print out the flags for the attributes of an object.
Useful in case the documentation is just too opaque!
=for example
PDL::Interpolate->attributes;
Flags Attribute
SGR x
SGR y
G err
G type
G bc
=cut
# note, can be called with the class, rather than just
# an object
#
# to allow this, I've used a horrible hack - we actually
# create an object and then print out the attributes from that
# Ugh!
#
sub
attributes {
my
$self
=
shift
;
# ugh
$self
=
$self
->new
unless
ref
(
$self
);
"Flags Attribute\n"
;
while
(
my
(
$attr
,
$hashref
) =
each
%{
$self
->{attributes}} ) {
my
$flag
=
""
;
$flag
.=
"S"
if
$hashref
->{settable};
$flag
.=
"G"
if
$hashref
->{gettable};
$flag
.=
"R"
if
$hashref
->{required};
printf
" %-3s %s\n"
,
$flag
,
$attr
;
}
return
;
}
####################################################################
=head1 AUTHOR
Copyright (C) 2000 Doug Burke (burke@ifa.hawaii.edu).
All rights reserved. There is no warranty.
You are allowed to redistribute this software / documentation as
described in the file COPYING in the PDL distribution.
=head1 SEE ALSO
L<PDL>, perltoot(1).
=cut
####################################################################
# End with a true
1;