[Numpy-discussion] Unexpected reorganization of internal data
Matthew Brett
matthew.brett@gmail....
Tue Jan 31 08:40:59 CST 2012
Hi,
On Tue, Jan 31, 2012 at 8:29 AM, Mads Ipsen <madsipsen@gmail.com> wrote:
> Hi,
>
> I am confused. Here's the reason:
>
> The following structure is a representation of N points in 3D space:
>
> U = numpy.array([[x1,y1,z1], [x1,y1,z1],...,[xn,yn,zn]])
>
> So the array U has shape (N,3). This order makes sense to me since U[i] will
> give you the i'th point in the set. Now, I want to pass this array to a C++
> function that does some stuff with the points. Here's how I do that
>
> void Foo::doStuff(int n, PyObject * numpy_data)
> {
> // Get pointer to data
> double * const positions = (double *) PyArray_DATA(numpy_data);
>
> // Print positions
> for (int i=0; i<n; ++i)
> {
> float x = static_cast<float>(positions[3*i+0])
> float y = static_cast<float>(positions[3*i+1])
> float z = static_cast<float>(positions[3*i+2])
>
> printf("Pos[%d] = %f %f %f\n", x, y, z);
> }
> }
>
> When I call this routine, using a swig wrapped Python interface to the C++
> class, everything prints out nice.
>
> Now, I want to apply a rotation to all the positions. So I set up some
> rotation matrix R like this:
>
> R = numpy.array([[r11,r12,r13],
> [r21,r22,r23],
> [r31,r32,r33]])
>
> To apply the matrix to the data in one crunch, I do
>
> V = numpy.dot(R, U.transpose()).transpose()
>
> Now when I call my C++ function from the Python side, all the data in V is
> printed, but it has been transposed. So apparently the internal data
> structure handled by numpy has been reorganized, even though I called
> transpose() twice, which I would expect to cancel out each other.
>
> However, if I do:
>
> V = numpy.array(U.transpose()).transpose()
>
> and call the C++ routine, everything is perfectly fine, ie. the data
> structure is as expected.
>
> What went wrong?
The numpy array reserves the right to organize its data internally.
For example, a numpy array can be in Fortran order in memory, or C
order in memory, and many more complicated schemes. You might want to
have a look at:
http://docs.scipy.org/doc/numpy/reference/arrays.ndarray.html#internal-memory-layout-of-an-ndarray
If you depend on a particular order for your array memory, you might
want to look at:
http://docs.scipy.org/doc/numpy/reference/generated/numpy.ascontiguousarray.html
Best,
Matthew
More information about the NumPy-Discussion
mailing list