Mailman 3 November 2007 - NumPy-Discussion

PyArray_ContiguousFromObject & PyDECREF
by Michiel Jan Laurens de Hoon 16 May '11

16 May '11

Dear pythoneers, In Chapter 12 ("Writing a C extension") in the Numerical Python manual, in the first example PyArray_ContiguousFromObject is used together with PyDECREF, while in the last example PyArray_ContiguousFromObject is used but the arrays are not PyDECREF'ed before the function returns. From Chapter 13 ("C API Reference") it seems that a PyDECREF is needed in both cases. Or am I missing something here? Is PyDECREF not needed for some reason in the last example? Thanks in advance, Michiel de Hoon U Tokyo First example: ============== static PyObject * trace(PyObject *self, PyObject *args) { PyObject *input; PyArrayObject *array; double sum; int i, n; if (!PyArg_ParseTuple(args, "O", &input)) return NULL; array = (PyArrayObject *) PyArray_ContiguousFromObject(input, PyArray_DOUBLE, 2, 2); if (array == NULL) return NULL; n = array->dimensions[0]; if (n > array->dimensions[1]) n = array->dimensions[1]; sum = 0.; for (i = 0; i < n; i++) sum += *(double *)(array->data + i*array->strides[0] + i*array->strides[1]); Py_DECREF(array); return PyFloat_FromDouble(sum); } Last example: ============= static PyObject * matrix_vector(PyObject *self, PyObject *args) { PyObject *input1, *input2; PyArrayObject *matrix, *vector, *result; int dimensions[1]; double factor[1]; double real_zero[1] = {0.}; long int_one[1] = {1}; long dim0[1], dim1[1]; extern dgemv_(char *trans, long *m, long *n, double *alpha, double *a, long *lda, double *x, long *incx, double *beta, double *Y, long *incy); if (!PyArg_ParseTuple(args, "dOO", factor, &input1, &input2)) return NULL; matrix = (PyArrayObject *) PyArray_ContiguousFromObject(input1, PyArray_DOUBLE, 2, 2); if (matrix == NULL) return NULL; vector = (PyArrayObject *) PyArray_ContiguousFromObject(input2, PyArray_DOUBLE, 1, 1); if (vector == NULL) return NULL; if (matrix->dimensions[1] != vector->dimensions[0]) { PyErr_SetString(PyExc_ValueError, "array dimensions are not compatible"); return NULL; } dimensions[0] = matrix->dimensions[0]; result = (PyArrayObject *)PyArray_FromDims(1, dimensions, PyArray_DOUBLE); if (result == NULL) return NULL; dim0[0] = (long)matrix->dimensions[0]; dim1[0] = (long)matrix->dimensions[1]; dgemv_("T", dim1, dim0, factor, (double *)matrix->data, dim1, (double *)vector->data, int_one, real_zero, (double *)result->data, int_one); return PyArray_Return(result); } -- Michiel de Hoon, Assistant Professor University of Tokyo, Institute of Medical Science Human Genome Center 4-6-1 Shirokane-dai, Minato-ku Tokyo 108-8639 Japan http://bonsai.ims.u-tokyo.ac.jp/~mdehoon

3 2

Unpleasant behavior with poly1d and numpy scalar multiplication
by Fernando Perez 15 Feb '10

15 Feb '10

Hi all, consider this little script: from numpy import poly1d, float, float32 p=poly1d([1.,2.]) three=float(3) three32=float32(3) print 'three*p:',three*p print 'three32*p:',three32*p print 'p*three32:',p*three32 which produces when run: In [3]: run pol1d.py three*p: 3 x + 6 three32*p: [ 3. 6.] p*three32: 3 x + 6 The fact that multiplication between poly1d objects and numbers is: - non-commutative when the numbers are numpy scalars - different for the same number if it is a python float vs a numpy scalar is rather unpleasant, and I can see this causing hard to find bugs, depending on whether your code gets a parameter that came as a python float or a numpy one. This was found today by a colleague on numpy 1.0.4.dev3937. It feels like a bug to me, do others agree? Or is it consistent with a part of the zen of numpy I've missed thus far? Thanks, f

6 17

Compact way of performing array math with specified result type?
by Russell E. Owen 08 Feb '10

08 Feb '10

I often find myself doing simple math on sequences of numbers (which might or might not be numpy arrays) where I want the result (and thus the inputs) coerced to a particular data type. I'd like to be able to say: numpy.divide(seq1, seq2, dtype=float) but ufuncs don't allow on to specify a result type. So I do this instead: numpy.array(seq1, dtype=float) / numpy.array(seq2, dtype=float) Is there a more compact solution (without having to create the result array first and supply it as an argument)? -- Russell

4 3

Slicing/selection in multiple dimensions simultaneously
by Mike Ressler 17 Mar '09

17 Mar '09

The following seems to be a wart: is it expected? Set up a 10x10 array and some indexing arrays: a=arange(100) a.shape=(10,10) q=array([0,2,4,6,8]) r=array([0,5]) Suppose I want to extract only the "even" numbered rows from a - then print a[q,:] <works - output deleted> Every fifth column: print a[:,r] <works - output deleted> Only the even rows of every fifth column: print a[q,r] --------------------------------------------------------------------------- <type 'exceptions.ValueError'> Traceback (most recent call last) /.../.../.../<ipython console> in <module>() <type 'exceptions.ValueError'>: shape mismatch: objects cannot be broadcast to a single shape But, this works: print a[q,:][:,r] [[ 0 5] [20 25] [40 45] [60 65] [80 85]] So why does the a[q,r] form have problems? Thanks for your insights. Mike -- mike.ressler(a)alum.mit.edu

8 16

Request for advice: project to get NumPy working in IronPython
by Giles Thomas 26 Apr '08

26 Apr '08

Hi, At Resolver Systems, we have a product that is written in IronPython - the .NET Python implementation - and allows users to use that language to script a spreadsheet-like interface. Because they're using IronPython, they can access their existing .NET objects and libraries, which has worked out really well for us and for them. But there's an increasing number of users who would like access to CPython C extensions - in particular, NumPy. We *don't* want to just hack something together that would just allow users of our own product to access NumPy; what we'd really like to do is work on an OSS solution that will allow all IronPython users to access NumPy transparently, just as if they were using CPython. My question is, is this something that the readers of this list would like to get involved in, either in terms of contibuting, or even in terms on advice? We're under no illusions; this is a significant undertaking. And we're very keen on hearing from those who know NumPy better about how this should be structured, even if people don't have the time to write code for it. In particular: * How do people feel about a source-code-compatible solution, where (perhaps) we would maintain a project that basically provided an alternative makefile for NumPy - or, even better, we could work with the NumPy developers to contribute a .NET/IronPython package. * Would it be better to try for some kind of "binary compatibility", where we'd write some kind of "glue" that sat between the existing C extension .pyd files and the IronPython engine? Our gut feeling is that this would be much more work, but we might be missing something. * What should the work's relationship to the NumPy project be? Any thoughts, comments, or advice would be very much appreciated. Regards, Giles -- Giles Thomas MD & CTO, Resolver Systems Ltd. giles.thomas(a)resolversystems.com +44 (0) 20 7253 6372 We're hiring! http://www.resolversystems.com/jobs/ 17a Clerkenwell Road, London EC1M 5RD, UK VAT No.: GB 893 5643 79 Registered in England and Wales as company number 5467329. Registered address: 843 Finchley Road, London NW11 8NA, UK

12 16

Problem with numpy.linalg.eig?
by Peter Creasey 04 Apr '08

04 Apr '08

Hi all, The following code calling numpy v1.0.4 fails to terminate on my machine, which was not the case with v1.0.3.1 from numpy import arange, float64 from numpy.linalg import eig a = arange(13*13, dtype = float64) a.shape = (13,13) a = a%17 eig(a) Regards, Peter

6 19

Loading a > GB file into array
by Martin Spacek 21 Dec '07

21 Dec '07

I need to load a 1.3GB binary file entirely into a single numpy.uint8 array. I've been using numpy.fromfile(), but for files > 1.2GB on my win32 machine, I get a memory error. Actually, since I have several other python modules imported at the same time, including pygame, I get a "pygame parachute" and a segfault that dumps me out of python: data = numpy.fromfile(f, numpy.uint8) # where f is the open file 1382400000 items requested but only 0 read Fatal Python error: (pygame parachute) Segmentation Fault If I stick to just doing it at the interpreter with only numpy imported, I can open up files that are roughly 100MB bigger, but any more than that and I get a clean MemoryError. This machine has 2GB of RAM. I've tried setting the /3GB switch on winxp bootup, as well as all the registry suggestions at http://www.msfn.org/board/storage-process-command-t62001.html. No luck. I get the same error in (32bit) ubuntu for a sufficiently big file. I find that if I load the file in two pieces into two arrays, say 1GB and 0.3GB respectively, I can avoid the memory error. So it seems that it's not that windows can't allocate the memory, just that it can't allocate enough contiguous memory. I'm OK with this, but for indexing convenience, I'd like to be able to treat the two arrays as if they were one. Specifically, this file is movie data, and the array I'd like to get out of this is of shape (nframes, height, width). Right now I'm getting two arrays that are something like (0.8*nframes, height, width) and (0.2*nframes, height, width). Later in the code, I only need to index over the 0th dimension, i.e. the frame index. I'd like to access all the data using a single range of frame indices. Is there any way to combine these two arrays into what looks like a single array, without having to do any copying within memory? I've tried using numpy.concatenate(), but that gives me a MemoryError because, I presume, it's doing a copy. Would it be better to load the file one frame at a time, generating nframes arrays of shape (height, width), and sticking them consecutively in a python list? I'm using numpy 1.0.4 (compiled from source tarball with Intel's MKL library) on python 2.5.1 in winxp. Thanks for any advice, Martin

13 29

Re: [Numpy-discussion] display numpy array as image
by Joe Harrington 03 Dec '07

03 Dec '07

If you want to explore the array interactively, blink images, mess with colormaps using the mouse, rescale the image values, mark regions, add labels, look at dynamic plots of rows and columns, etc., get the ds9 image viewer and the xpa programs that come with it that allow it to communicate with other programs: ftp://sao-ftp.harvard.edu/pub/rd/ds9 http://hea-www.harvard.edu/RD/ds9/index.html Then get the Python numdisplay package, which uses xpa. You have to get numdisplay from inside the stsci_python package: http://www.stsci.edu/resources/software_hardware/pyraf/stsci_python/current… Just grab the numdisplay directory from within that. Older versions of numdisplay are standalone but don't work perfectly. Beware, there are outdated web sites about numdisplay on the stsci site. Don't google! Run ds9 before you load numdisplay. Then you can send your python arrays to a real interactive data viewer at will. There are even mechanisms to define physical coordinates mapped from the image coordinates. --jh--

4 4

swig numpy2carray converters
by Georg Holzmann 01 Dec '07

01 Dec '07

Hallo! Because I had some troubles in wrapping my C++ library in python/numpy, I did (another) numpy2carray.i file for SWIG. With that interface file it is possible to input/output arrays with or without copying data (I also included an example for an interface to fortran style arrays). I am still a numpy newbie, so this might not be perferct, but however, I finally found some time to put it online - maybe it is helpfull for other beginners: http://grh.mur.at/software/numpy2carray.html LG Georg

7 31

Re: [Numpy-discussion] display numpy array as image (Giorgio F. Gilestro)
by Joe Harrington 01 Dec '07

01 Dec '07

I was misinformed about the status of numdisplay's pages. The package is available as both part of stsci_python and independently, and its (up-to-date) home page is here: http://stsdas.stsci.edu/numdisplay/ Googling numdisplay finds that page. My apologies to those inconvenienced by my prior post. Note that there is practical documentation and examples of how to use it for astronomical image display in the Data Analysis Tutorial: http://www.scipy.org/wikis/topical_software/Tutorial --jh--

2 1