Filling numpy array with C routine

classic Classic list List threaded Threaded
1 message Options
Reply | Threaded
Open this post in threaded view
|

Filling numpy array with C routine

jeinarsson
This post has NOT been accepted by the mailing list yet.
Dear list,

I need an opinion on good practice.

I'd like to write simple programs that
1) (In Python) allocates numpy array,
2) (In C/C++) fills said numpy array with data.

To this end I use Boost.Python to compile an extension module. I use the (possibly obsolete?) boost/python/numeric.hpp to allow passing an ndarray to my C-functions. Then I use the numpy C API directly to extract a pointer to the underlying data.

This seemingly works very well, and I can check for correct dimensions and data types, etcetera.

As documentation is scarce, I ask you if this is an acceptable procedure? Any pitfalls nearby?

Sample code: C++
void fill_array(numeric::array& y)
{
	const int ndims = 2;

	// Get pointer to np array
	PyArrayObject* a = (PyArrayObject*)PyArray_FROM_O(y.ptr());
	if (a == NULL) {
                throw std::exception("Could not get NP array.");
        }
	if (a->descr->elsize != sizeof(double)) 
	{
		throw std::exception("Must be double ndarray");
	}
	if (a->nd != ndims) 
	{
		throw std::exception("Wrong dimension on array.");
	}
	int rows = *(a->dimensions);
	int cols = *(a->dimensions+1);
	double* data = (double*)a->data;

	for (int i = 0; i < rows; i++)
	{
		for (int j = 0; j < cols; j++) 
		{
			*(data + i*cols + j) = really_cool_function(i,j);
		}
	}
}

BOOST_PYTHON_MODULE(Practical01)
{
	import_array();
	boost::python::numeric::array::set_module_and_type("numpy", "ndarray");

	def("fill_array",&fill_array);
}


And in python this could be used such as:
import Practical01
import numpy
import matplotlib.pyplot as plt
import matplotlib.cm as colormaps
import time


w=500
h=500
large_array = numpy.ones( (h,w) );

t1 = time.time()
Practical01.fill_array(large_array)
t2 = time.time()
print 'Horrible calculation took %0.3f ms' % ((t2-t1)*1000.0)

plt.imshow(large_array,cmap=colormaps.gray)
plt.show()


Simplicity is a major factor for me. I don't want a complete wrapper for ndarrays, I just want to compute and shuffle data to Python for further processing. Letting Python handle allocation and garbage collection also seems like a good idea.

Sincerely,
Jonas Einarsson