[Matrix-SIG] patch - generalize arrayfns.interp to mulitple dimensions

[Joe Van Andel]

This is a multipart MIME message.

--==_Exmh_19846673670
Content-Type: text/plain; charset=us-ascii

I just stumbled across the arrayfns.interp function:
   /* interp (y, x, z) treats (x, y) as a piecewise linear function
    * whose value is y [0] for x < x [0] and y [len (y) -1] for x >
    * x [len (y) -1]. An array of floats the same length as z is
    * returned, whose values are ordinates for the corresponding z
    * abscissae interpolated into the piecewise linear function.         */

It lives in

LLNLDistribution11/Graphics/Arrayfcns/Src/arrayfnsmodule.c

(Because I found this function, I didn't need to write my
function in 'C', and could write a fuzzy logic recognizer in about 40 lines of 
Python.)

As distributed, the input and output arrays can only have a single 
dimension.  I needed a more general version, and rather than reshaping 
my input and output arrays, I modified this routine to accept arrays 
with up to 'MAX_INTERP_DIMS (currently 6) dimensions.

With the (limited) testing I've done, the new version works just fine.

I've attached a context diff, and hope this change is accepted for the 
next release from LLNL.



Joe VanAndel  		Internet: vanandel@ucar.edu
National Center for	http://www.atd.ucar.edu/~vanandel/home.html
Atmospheric Research


--==_Exmh_19846673670
Content-Type: text/plain ; name="arrayfns.diff"; charset=us-ascii
Content-Description: arrayfns.diff

*** 1.1	1999/04/14 22:58:32
--- arrayfnsmodule.c	1999/04/15 14:14:16
***************
*** 6,11 ****
--- 6,13 ----
  #include <stdio.h>
  #include <stdlib.h>
  
+ #define MAX_INTERP_DIMS 6
+ 
  static PyObject *ErrorObject;
  
  /* Define 2 macros for error handling:
***************
*** 34,39 ****
--- 36,43 ----
  #define A_DIM(a,i) (((PyArrayObject *)a)->dimensions[i])
  #define GET_ARR(ap,op,type,dim) \
    Py_Try(ap=(PyArrayObject *)PyArray_ContiguousFromObject(op,type,dim,dim))
+ #define GET_ARR2(ap,op,type,min,max) \
+   Py_Try(ap=(PyArrayObject *)PyArray_ContiguousFromObject(op,type,min,max))
  #define ERRSS(s) ((PyObject *)(PyErr_SetString(ErrorObject,s),0))
  #define SETERR(s) if(!PyErr_Occurred()) ERRSS(errstr ? errstr : s)
  #define DECREF_AND_ZERO(p) do{Py_XDECREF(p);p=0;}while(0)
***************
*** 597,609 ****
         Py_DECREF(ay);
         Py_DECREF(ax);
         return NULL ;}
!    GET_ARR(az,oz,PyArray_DOUBLE,1);
     lenz = A_SIZE (az);
     dy = (double *) A_DATA (ay);
     dx = (double *) A_DATA (ax);
     dz = (double *) A_DATA (az);
!    Py_Try (_interp = (PyArrayObject *) PyArray_FromDims (1, &lenz, 
!       PyArray_DOUBLE));
     dres = (double *) A_DATA (_interp) ;
     slopes = (double *) malloc ( (leny - 1) * sizeof (double)) ;
     for (i = 0 ; i < leny - 1; i++) {
--- 601,614 ----
         Py_DECREF(ay);
         Py_DECREF(ax);
         return NULL ;}
!    GET_ARR2(az,oz,PyArray_DOUBLE,1,MAX_INTERP_DIMS);
     lenz = A_SIZE (az);
     dy = (double *) A_DATA (ay);
     dx = (double *) A_DATA (ax);
     dz = (double *) A_DATA (az);
!    /* create output array with same size as 'Z' input array */
!    Py_Try (_interp = (PyArrayObject *) PyArray_FromDims
!            (A_NDIM(az), az->dimensions, PyArray_DOUBLE));
     dres = (double *) A_DATA (_interp) ;
     slopes = (double *) malloc ( (leny - 1) * sizeof (double)) ;
     for (i = 0 ; i < leny - 1; i++) {

--==_Exmh_19846673670--