Python/C API, Numeric Python, Type Conversion

[ronysk]

Hi,

I am posting here to seek for help on type conversion between Python
(Numeric Python) and C.

Attachment A is a math function written in C, which is called by a Python
program. I had studied SWIG and Python/C API a bit. I was able to pass
numeric array (Numeric Python) into the C function and access/change these
arrays.

Problem I am having is that I couldn't quite get the array converted back to
some PyObject or something that Python Interpreter can understand. So,
alternatively I passed an empty numeric array as well as other arguments,
and had return value(result) stored in that chunk of place. It sort of
worked. However, if such function was called few times, it would result in
segmentation fault.

This is a lousy way to get around the return type matter, and probably is
causing seg fault. Besides a correct type conversion, a proper SWIG
interface probably could've saved all the fuss, maybe? ( attachment B is the
SWIG interface I am using for this function)

I am trying to study more on C API section of NumPy document, and to gain
more knowledge off SWIG(which I only understand very basic way of it). Any
feedback will be great. Thanks for your time, friends. Truly appreicate
your precious time.

Ronnie


P.S. Those PyObject pointers are created with Numeric.Array( data,
Numeric,Float64).


Attachment A

/*
*******************************************************************************
*
*  Function Name : 
*     cressman()
*  
*  Usage :
*     PyObject* cressman
*     (int nx,int ny,PyObject* result,PyObject* us,PyObject* rlats,
PyObject* 
*      rlons, int nst,float disinf,float slon,float elon,float slat,float
elat)
*
*  Description :
*     Objective analysis by Cressman method
*
*  Aumument :
*     nx(I)    : x-dimension of grid data 
*     ny(I)    : y-dimension of grid data
*     result(I): an empty result set
*     us(I)    : input station data                
*     rlons(I) : x-coordinate of station
*     rlats(I) : y-coordinate of station
*     nst(I)   : number of station data
*     disinf(I): influence radius distance(km)
*     slon(I)  : x-coordinate of start grid longitude
*     elon(I)  : x-coordinate of end grid longitude
*     slat(I)  : y-coordinate of start grid latitude
*     elat(I)  : y-coordinate of end grid latitude
*     
*   Data Files : None
*  
*   Called Function :
*     ardist()    
*
*   Return Value : 
*      data    : grid data value    
*
*   Linking Lib : None
*
***************************************************************************                         
*/
PyObject* cressman
(int nx,int ny,PyObject* result,PyObject* us,PyObject* rlats, PyObject*
rlons, int nst,float disinf,float slon,float elon,float slat,float elat)
{
   if( nst < 1)
       return;

   // initialization    
   float xdelta,ydelta,disinf2,wsum,usum,xx,yy;
   float **data;
   float data_ave, sum;
   int i,j,k,count;
    
   // PyObject(python-C API) related 
   PyArrayObject *p_us, *p_rlats, *p_rlons, *p_result;
   double *array_us, *array_rlats, *array_rlons, *array_result;

   p_us= (PyArrayObject *) us;
   p_rlats= (PyArrayObject *) rlats;
   p_rlons= (PyArrayObject *) rlons;
   p_result= (PyArrayObject *) result;
   
   if ( nst != p_us->dimensions[0])
       return;

   array_us = (double*) p_us->data;
   array_rlats = (double*) p_rlats->data;
   array_rlons = (double*) p_rlons->data;
   array_result = (double*) p_result->data;
   
   // find out data avg
   sum=0;
   count=0;
   for (i=0;i<nst;i++)
   { 
       if (array_us[i]>spv)
       {
           sum+=array_us[i];
           count++;
       }
   }
   if (count!=0)
       data_ave= sum / count;

   xdelta = (elon - slon)/(nx -1);  
   ydelta = (elat- slat)/(ny -1);
   disinf2 = pow(disinf,2);

   for( i = 0; i < ny; i++)
   {
        for( j = 0; j < nx; j++)
        {
             wsum = 0.;
             usum = 0.;
             xx = slon + j*xdelta;
             yy = slat + i*ydelta;
             
             for( k = 0; k < nst; k++)
             {
                  if( array_us[k] <= spv )
                      continue;
                  
                  //ardist just a function that calculates distance between
2 points
                  float adist = ardist(array_rlons[k],array_rlats[k],xx,yy);
                  if( adist > disinf )
                      continue;
                  float adist2 =  pow(adist,2);
                  float w = (disinf2 - adist2)/(disinf2+adist2);
                 
                  wsum = wsum + w;
                  usum = usum + (w*array_us[k]);       
             }
             if( wsum != 0 )
                 //data[i][j] = usum/wsum;
                 array_result[i*nx+j] = usum/wsum;
             else
                 //data[i][j] = data_ave;
                 array_result[i*nx+j] = data_ave;
             
        }
   }
   return result;       
}

Attachment B - SWIG interface

 /* objmod.i */
 %module objmod
 %{
 /* Put header files here or function declarations like below */
 PyObject* cressman(int nx,int ny,PyObject* result,PyObject* us,
                     PyObject* rlats,PyObject* rlons,int nst,float disinf,
                     float slon,float elon,float slat,float elat);
 float ardist(float rlon,float rlat,float xx,float yy);
 %}
 
 extern PyObject* cressman(int nx,int ny,PyObject* result,PyObject* us,
                            PyObject* rlats,PyObject* rlons,int nst,float
disinf,
                            float slon,float elon,float slat,float elat);
 extern float ardist(float rlon,float rlat,float xx,float yy);