Mailman 3 August 2008 - NumPy-Discussion

mixed mode arithmetic
by Neal Becker July 11, 2023

July 11, 2023

I've been browsing the numpy source. I'm wondering about mixed-mode arithmetic on arrays. I believe the way numpy handles this is that it never does mixed arithmetic, but instead converts arrays to a common type. Arguably, that might be efficient for a mix of say, double and float. Maybe not. But for a mix of complex and a scalar type (say, CDouble * Double), it's clearly suboptimal in efficiency. So, do I understand this correctly? If so, is that something we should improve?

4 6

PyArray_ContiguousFromObject & PyDECREF
by Michiel Jan Laurens de Hoon May 16, 2011

May 16, 2011

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 Feb. 15, 2010

Feb. 15, 2010

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 Feb. 8, 2010

Feb. 8, 2010

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

Is anyone knowledgeable about dll deployment on windows ?
by David Cournapeau Dec. 1, 2009

Dec. 1, 2009

Hi, My head hurt trying to understand dll management with windows. I wanted to find a sane way to use dll for numscons, but I can't see how to do this, so I was wondering if anyone on this ML had any deep knowledge on how to install dll, and reuse them with python extensions ? The problem is the following: - you have some dll (say MKL) in a path (say C:\program files\MKL). - how do you tell another dll (more exactly a python extension, a .pyd) to look for mkl dll into the MKL path ? It seems that there is no rpath-like capability on windows, and the only way to load dll is to put them into some windows directory, or to put the dll path into PATH. Is it really the only way ? cheers, David

8 25

Slicing/selection in multiple dimensions simultaneously
by Mike Ressler March 17, 2009

March 17, 2009

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

access ndarray in C++
by Thomas Hrabe Dec. 10, 2008

Dec. 10, 2008

Hi all! I am currently developing a python module in C/C++ which is supposed to access nd arrays as defined by the following command in python a = numpy.array([1,1,1]) I want to access the array the following way and use the nd array data for further processing in C. mymod.doSthg(a) The example code on http://numpy.sourceforge.net/numdoc/HTML/numdoc.htm#pgfId-36721 (!PyArg_ParseTuple(args, "O!",&PyArray_Type, &array)) does not work for nd arrays. I always get TypeError: argument 1 must be array, not numpy.ndarray I assume the error is the constant provided as the third parameter, saying that the imput is of PyArray_Type and no nd array. So here are my questions: 1. is there any good tutorial / example code for acessing nd arrays in C? 2. what is the difference between both (arrays and nd arrays? - I am new to python and heard there were different approaches to arrays and that nd arrays work better for multi dimensional applications. Is this right?) 3. which one of both will be used in the future? Thank you in advance for your help, Thomas

15 31

error importing a f2py compiled module.
by Fabrice Silva Nov. 30, 2008

Nov. 30, 2008

Dear all I've tried to run f2py on a fortran file which used to be usable from python some months ago. Following command lines are applied with success (no errors raised) : f2py -m modulename -h tmpo.pyf --overwrite-signature tmpo.f f2py -m modulename -c --f90exec=/usr/bin/f95 tmpo.f The output of these commands is available here: http://paste.debian.net/7307 When importing in Python with "import modulename", I have an ImportError: Traceback (most recent call last): File "Solveur.py", line 44, in <module> import modulename as Modele ImportError: modulename.so: failed to map segment from shared object: Operation not permitted How can that be fixed ? Any suggestion ? Thanks -- Fabrice Silva LMA UPR CNRS 7051 - équipe S2M

3 6

Bilateral filter
by Nadav Horesh Nov. 25, 2008

Nov. 25, 2008

Attached here my cython implementation of the bilateral filter, which is my first cython program. I would ask for the following: 1. Is there any way to speed up the code just by "cosmetic" modifications? 2. In particular I use the unportable gcc function __builtin_abs: Is there any way to access this in a portable way? 3. I would like to donate this code to scipy or any other suitable project. What can/should I do to realise it? Remarks: The code contains 3 end-user routines: 1. A pure python implementation: Easy to read and modify --- it can be cut out into a python source code. 2. A straight forward cython implementation: About 4 times as fast as the python implementation. 3. An optimised cython implementation earning another factor of 2 in speed (depends on the parameters used). I see this code as a "research grade" that would evolve in the near future, as I am working on a related project, and hopefully following your comments. Nadav

6 7

Masked arrays and pickle/unpickle
by Anthony Floyd Nov. 7, 2008

Nov. 7, 2008

We have an application that has previously used masked arrays from Numpy 1.0.3. Part of saving files from that application involved pickling data types that contained these masked arrays. In the latest round of library updates, we've decided to move to the most recent version of matplotlib, which requires Numpy 1.1. Unfortunately, when we try to unpickle the data saved with Numpy 1.0.3 in the new code using Numpy 1.1.0, it chokes because it can't import numpy.core.ma for the masked arrays. A check of Numpy 1.1.0 shows that this is now numpy.ma.core. Does anyone have any advice on how we can unpickle the old data files and update the references to the new classes? Thanks, Anthony. -- Anthony Floyd, PhD Convergent Manufacturing Technologies Inc. 6190 Agronomy Rd, Suite 403 Vancouver BC V6T 1Z3 CANADA Email: Anthony.Floyd(a)convergent.ca | Tel: 604-822-9682 x102 WWW: http://www.convergent.ca | Fax: 604-822-9659 CMT is hiring: See http://www.convergent.ca for details

5 10

patch for new mgrid / ogrid functionality
by David M. Kaplan Sept. 13, 2008

Sept. 13, 2008

Hi, A while back, I sent some changes to index_tricks.py that would allow mgrid and ogrid to mesh things other than slices. For example: >>> mgrid[['a','b'],[float,int],:3] [array([[['a', 'a', 'a'], ['a', 'a', 'a']], [['b', 'b', 'b'], ['b', 'b', 'b']]], dtype='|S1'), array([[[<type 'float'>, <type 'float'>, <type 'float'>], [<type 'int'>, <type 'int'>, <type 'int'>]], [[<type 'float'>, <type 'float'>, <type 'float'>], [<type 'int'>, <type 'int'>, <type 'int'>]]], dtype=object), array([[[0, 1, 2], [0, 1, 2]], [[0, 1, 2], [0, 1, 2]]])] At the time, there wasn't much follow-up, but I am hoping that there is still interest in this functionality, as I have gone ahead and finished the patch including documentation changes and updates to test_index_tricks.py. Attached is a patch set to the latest subversion of the numpy trunk. I don't think I am allowed to commit the changes myself - correct me if I am wrong. This functionality seems like a nice addition to me as it allows one to mesh things that are not uniformly spaced and potentially not even numbers. The changes don't affect functionality that existed previously except for one minor change - instead of returning a numpy array of arrays, mgrid/ogrid now return a list of arrays. However, this is unlikely to be a problem as the majority of users generally unpack the results of mgrid/ogrid so that each matrix can be used individually. Comments welcome. Cheers, David -- ********************************** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement Centre de Recherche Halieutique Mediterraneenne et Tropicale av. Jean Monnet B.P. 171 34203 Sete cedex France Phone: +33 (0)4 99 57 32 27 Fax: +33 (0)4 99 57 32 95 http://www.ur097.ird.fr/team/dkaplan/index.html **********************************

3 4

Updated Numpy reference guide
by Pauli Virtanen Sept. 4, 2008

Sept. 4, 2008

Hi all, I finished the first iteration of incorporating material from Travis Oliphant's "Guide to Numpy" to the Sphinxy reference guide we were constructing in the Doc marathon. Result is here: (the PDF is a bit ugly, though, some content is almost randomly scattered there) http://www.iki.fi/pav/tmp/numpy-refguide/index.xhtml http://www.iki.fi/pav/tmp/numpy-refguide/NumPy.pdf Source is here: (Stéfan, if it looks ok to you, could you pull and check if it builds for you when you have time?) https://code.launchpad.net/~pauli-virtanen/scipy/numpy-refguide What I did with the "Guide to Numpy" material was: - Collapsed each of the reference Chapters 3, 6, 8, 9 (ndarrays, scalars, dtypes, ufuncs) with the more introductory material in Chapter 2. - As this was supposed to be a reference guide, I tried to compress the text from Chapter 2 as much as possible, by sticking to definitions and dropping some more tutorial-oriented parts. This may have reduced readability at some points... - I added some small bits or rewrote parts in the above sections in places where I thought it would improve the result. - I did not include material that I thought was better to be put into appropriate docstrings in Numpy. What to do with class docstrings and obscure __xxx__ attributes was not so clear a decision, so what I did for these varies. - The sections about Ufuncs and array indexing are taken almost verbatim from the "Guide to Numpy". The ndarray, scalar and dtype sections somewhat follow the structure of the Guide, but the text is more heavily edited from the original. Some things to do: - Descriptions about constructing items with __new__ methods should probably still be clarified; I just replaced references to __new__ with references to the corresponding classes. - What to do with the material from numpy.doc.* should be decided, as the text there doesn't look like it should go into a reference manual. Some questions: - Is this good enough to go into Numpy SVN at some point? Or should we redo it and base the work closer to the original "Guide to Numpy"? - Does it build for you? (I'd recommend using the development 0.5 version of Sphinx, so that you get the nifty Inter-Sphinx links to the Python documentation.) We are unfortunately beating the Sphinx with a big stick to make it place the documentation of each function or class into a separate file, and to convert the Numpy docstring format to something the Sphinx can fathom. There's also some magic in place to make toctrees:: of function listings more pleasant to the eye. Any comments of what should be improved are welcome. (Even better: clone the bzr branch, make the changes yourself, and put the result somewhere available! E.g. as a bzr bundle or a branch on the launchpad.) -- Pauli Virtanen

2 1

at my wits end over an error message...
by Alan Jackson Sept. 3, 2008

Sept. 3, 2008

I been beating myself up over this bit of code for far too long now - I know I must be missing something really simple, but what is it? TYPICALLY_UINT_COLUMNS = ['Track', 'Bin', 'code', 'horizon'] ...... dtypes = [ ] for i in range(0, len(self.var_list)) : if TYPICALLY_UINT_COLUMNS.count(self.var_list[i]) > 0: dtypes.append((self.var_list[i], '<i8')) else : dtypes.append((self.var_list[i], '<f8')) print "dtypes = ", dtypes print outdata[0] outdata = np.array(outdata, dtype=dtypes) dtypes = [('Track', '<i8'), ('Bin', '<i8'), ('Depth', '<f8'), ('extrema', '<f8'), ('code', '<i8'), ('horizon', '<i8')] (345, 615, 4826.0, 470.588012695, 47, 636) Exception occurred in traits notification handler. Please check the log file for details. Exception occurred in traits notification handler for object: <__main__.MaplEdit object at 0x8bb54dc>, trait: Run, old value: <undefined>, new value: 0 Traceback (most recent call last): File "/usr/local/lib/python2.5/site-packages/enthought.traits-2.0.4-py2.5-linux-i686.egg/enthought/traits/trait_notifiers.py", line 325, in call_1 self.handler( object ) File "TrimMapl_1.py", line 98, in _Run_fired outdata = np.array(outdata, dtype=dtypes) TypeError: expected a readable buffer object Thanks, Alan -- ----------------------------------------------------------------------- | Alan K. Jackson | To see a World in a Grain of Sand | | alan(a)ajackson.org | And a Heaven in a Wild Flower, | | www.ajackson.org | Hold Infinity in the palm of your hand | | Houston, Texas | And Eternity in an hour. - Blake | -----------------------------------------------------------------------

3 3

floating point char - bug?
by Christian K. Aug. 31, 2008

Aug. 31, 2008

Hi, I just came across somethin I never noticed before. I cannot say whether this is due to an update of numpy but it is possible - I am running 1.1.1 on __german__ windows. Here is the observation: a = N.linspace(0,1,5) a array([ 0. , 0.25, 0.5 , 0.75, 1. ]) a.astype(float) array([ 0. , 0.25, 0.5 , 0.75, 1. ]) a[0].astype(float) 0.0 a[1].astype(float) 0,25 As you see in the last line, suddenly numpy picks up the german locale setting and converts the floating point into a comma. It does not affect the '0.0' in the line above and I believe to have seen some other numbers ending with '.0' where the point has not been replaced. I guess this is a bug. In fact I do not like the idea that repr() of a numpy float honours the locale settings. Reagrds, Christian

2 2

sort documentation
by Alan G Isaac Aug. 31, 2008

Aug. 31, 2008

I find this confusing: numpy.sort(a, axis=-1, kind='quicksort', order=None) Return copy of 'a' sorted along the given axis. Perform an inplace sort along the given axis using the algorithm specified by the kind keyword. I suppose the last bit is supposed to refer to the ``sort`` method rather than the function, but I do not see any signal that this is the case. Cheers, Alan Isaac

5 7

Advice on converting iterator into array efficiently
by Alan Jackson Aug. 31, 2008

Aug. 31, 2008

Looking for advice on a good way to handle this problem. I'm dealing with large tables (Gigabyte large). I would like to efficiently subset values from one column based on the values in another column, and get arrays out of the operation. For example, say I have 2 columns, "energy" and "collection". Collection is basically an index that flags values that go together, so all the energy values with a collection value of 18 belong together. I'd like to be able to set up an iterator on collection that would hand me an array of energy on each iteration : if table is all my data, then something like for c in table['collection'] : e = c['energy'] ... do array operations on e I've been playing with pytables, and they help, but I can't quite seem to get there. I can get an iterator for energy within a collection, but I can't figure out an efficient way to get an array out. What I have so far is for h in np.unique(table.col('collection')) : rows = table.where('collection == c') for row in rows : print c,' : ', row['energy'] but I really want to convert rows['energy'] to an array. I've thought about building a nasty set of pointers and whatnot - I did it once in perl - but I'm hoping to avoid that. -- ----------------------------------------------------------------------- | Alan K. Jackson | To see a World in a Grain of Sand | | alan(a)ajackson.org | And a Heaven in a Wild Flower, | | www.ajackson.org | Hold Infinity in the palm of your hand | | Houston, Texas | And Eternity in an hour. - Blake | -----------------------------------------------------------------------

3 5

huge array calculation speed
by Dan Lussier Aug. 30, 2008

Aug. 30, 2008

Hello, I am relatively new to numpy and am having trouble with the speed of a specific array based calculation that I'm trying to do. What I'm trying to do is to calculate the total total potential energy and coordination number of each atom within a relatively large simulation. Each atom is at a position (x,y,z) given by a row in a large array (approximately 1e6 by 3) and presently I have no information about its nearest neighbours so each its position must be checked against all others before cutting the list down prior to calculating the energy. My current numpy code is below and in it I have tried to push as much of the work for this computation into compiled C (ideally) of numpy. However it is still taking a very long time at approx 1 second per row. At this speed even some simple speed ups like weave won't really help. Are there any other numpy ways that it would be possible to calculate this, or should I start looking at going to C/C++? I am also left wondering if there is something wrong with my installation of numpy (i.e. not making use of lapack/ATLAS). Other than that if it is relevant - I am running 32 bit x86 Centos 5.1 linux on a dual Xeon 3.0 GHz Dell tower with 4 GB of memory. Any suggestions of things to try would be great. Dan ======================================================================== ======================= def epairANDcoord(xyz,cutoff=4.0): box_length = 160.0 cut2 = cutoff**2 ljTotal = numpy.zeros(xyz.shape[0]) coord = numpy.zeros(xyz.shape[0]) i = 0 for r0 in xyz: r2 = r0-xyz # application of minimum image convention so each atom interacts with nearest periodic image w/i box_length/2 r2 -= box_length*numpy.rint(r2/box_length) r2 = numpy.power(r2,2).sum(axis=1) r2 = numpy.extract(r2<cut2, r2) coord[i] = r2.shape[0]-1 r2 = numpy.extract(r2 != 0, r2) # avoid problems of inf when considering own atom. r2 = numpy.divide(1.0,r2) ljTotal[i] = numpy.multiply(2,(numpy.power(r2,6)- numpy.power(r2,3))).sum(axis=0) # ljTotal[i] = ljPerAtom.sum(axis=0) i += 1 del r2 # if i>25: break # used to cut down the length of calculation as needed return ljTotal,coord

11 13

How do I do this?
by Christopher Barker Aug. 30, 2008

Aug. 30, 2008

HI all, I need to do something I thought would be simple -- set all the values of an array to some minimum. so I did this: >>> min_value = 2 >>> a = np.array((1, 2, 3, 4, 5,)) >>> np.maximum(a, min_value) array([2, 2, 3, 4, 5]) all was well... then I realized that a could have negative numbers in in, and I really wanted the absolute value to be greater than than minimum: >>> a = np.array((1, 2, 3, 4, -5,)) >>> np.maximum(a, min_value) array([2, 2, 3, 4, 2]) oops! so I added a sign() and abs(): >>> np.sign(a) * np.maximum(np.abs(a), min_value) array([ 2, 2, 3, 4, -5]) all was well. However it turns out a could contain a zero: >>> a = np.array((0, 1, 2, 3, 4, -5,)) >>> np.sign(a) * np.maximum(np.abs(a), min_value) array([ 0, 2, 2, 3, 4, -5]) Darn! I want that zero to become a 2, but sign(0) = 0, so that doesn't work. How can I do this without another line of code special casing the 0, which isn't that big I deal, but it seems kind of ugly... >>> a[a==0] = min_value >>> np.sign(a) * np.maximum(np.abs(a), min_value) array([ 2, 2, 2, 3, 4, -5]) thanks, -Chris -- Christopher Barker, Ph.D. Oceanographer Emergency Response Division NOAA/NOS/OR&R (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception Chris.Barker(a)noaa.gov

4 16

default index inconsistency?
by Alan G Isaac Aug. 30, 2008

Aug. 30, 2008

I think the question of the timing of removing the remaining default index inconsistencies in NumPy was raised by not answered? Anyway, what is the answer? I just ran into this again, this time with ``argsort``, which still has ``axis=-1`` as the default. But maybe it is intentional to keep this? E.g., a.sort(axis=None) is refused. I would have expected it to return the same as np.sort(a.flat), possibly (?) reshaped. Alan Isaac

1 0

np.test() Bus error on OSX. Reopened Ticket 816
by Christopher Burns Aug. 29, 2008

Aug. 29, 2008

I reopened Ticket 816: http://scipy.org/scipy/numpy/ticket/816 Running numpy.test() causes a bus error on OSX. -- Christopher Burns Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

2 1

isn't it a bug in array.fill()?
by dmitrey Aug. 29, 2008

Aug. 29, 2008

hi all, isn't it a bug (latest numpy from svn, as well as my older version) from numpy import array print array((1,2,3)).fill(10) None Regards, D.

5 10

Re: [Numpy-discussion] isn't it a bug in array.fill()
by dmitrey Aug. 29, 2008

Aug. 29, 2008

sorry, it isn't a bug, it's my fault, fill() returns None and do in-place modification. D.

1 0

Problem with numpy in pyinstaller
by Jeffrey Starn Aug. 29, 2008

Aug. 29, 2008

I have a program that imports numpy and scipy. When I try to use PyInstaller, I get the following error: File "<string>", line 21, in <module> NameError: name 'io' is not defined The code up to line 21 is (line numbers differ becasue of comments): from numpy import * from numpy import random from scipy import * from scipy.linalg import cholesky from scipy.io import write_array model=raw_input('MF2K or UCode? ') model=model.upper() fname=raw_input('Root name? ') NumRan=raw_input('Number of random realizations? ') NumRan=int(NumRan) outfile=fname+".ranvar" MvName=fname+"._mv" Mv=io.read_array(MvName,lines=(1,-1) I have seen related psots, but I really don't understand what the fix is. Thanks. --------------------------------------------------------------------------------------------------------- Jeff Starn, Ground-water Specialist U.S. Geological Survey 101 Pitkin Street East Hartford, CT 06108 (860) 291-6746 jjstarn(a)usgs.gov ---------------------------------------------------------------------------------------------------------

1 0

#882: PyArray_ConvertToCommon and 0 sized sequences
by David Cournapeau Aug. 28, 2008

Aug. 28, 2008

Hi, The error in #882 looks pretty straightforward: from numpy import numarray a = numarray.array(1) a.choose([]) core dumps because a.choose calls PyArray_ConvertCommonType with a 0-sized array, which is not handled correctly. Should the 0-sized case be handled in this function, or in PyArray_Choose (the caller of PyArray_ConvertCommonType) ? cheers, David

1 0

Where does _assert_func come from.
by Charles R Harris Aug. 28, 2008

Aug. 28, 2008

Where does _assert_func come from in this fragment from a test class? class _GenericTest(object): def _test_equal(self, a, b): self._assert_func(a, b) I hope it isn't monkey patched in. Chuck

2 5

numpy and py2exe
by Christopher Barker Aug. 28, 2008

Aug. 28, 2008

Hi all, An issue just came up with numpy and py2exe on the wxPython list. A solution has been found, but I thought I'd post here, so that it could be on the numpy developers' radar when considering the structure of numpy for the future: > I'm compile a wxPython script using GUI2Exe, and am getting the > following when running the produced exe: > > $ ./dist/svm.exe > Traceback (most recent call last): > File "svm.py", line 28, in <module> > import graphic_svm_data > File "graphic_svm_data.pyo", line 29, in <module> > File "floatcanvas\NavCanvas.pyo", line 7, in <module> > File "floatcanvas\FloatCanvas.pyo", line 7, in <module> > File "numpy\__init__.pyo", line 93, in <module> > File "numpy\add_newdocs.pyo", line 9, in <module> > File "numpy\lib\__init__.pyo", line 19, in <module> > File "numpy\lib\financial.pyo", line 78, in <module> > TypeError: unsupported operand type(s) for +=: 'NoneType' and 'str' > I tried to add numpy.lib.financial to the Excludes and Ignores list, but am still getting this error. > > Googling found that this problem is a known bug: http://article.gmane.org/gmane.comp.python.py2exe/2920 The solution was to: > 2) Use Optimize=1 or Optimize=0 in GUI2Exe and re-build your > executable. Uing Optimize=2 strips all docstrings in the executable, > which may cause troubles (as in your case); > [Ron Barak] YES!!! > With Optimize=1 my script finally compiles to exe. While this worked, I have a couple observations: 1) This broke on numpy.lib.financial. that had nothing to do with the issue, but it would be nice if people could use numpy with py2exe and the like without having to include all those unused packages. I think this ties in with the discussion about speeding up with importing process, and may be addressed in future versions if I understood the results of some of the discussion at SciPy right. 2) it would also be nice if one could use numpy with the docstrings removed, like in this case. I have no idea what the issues are with that. -Chris -- Christopher Barker, Ph.D. Oceanographer Emergency Response Division NOAA/NOS/OR&R (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception Chris.Barker(a)noaa.gov

2 3

ImportError: No module named my_module
by Nils Wagner Aug. 28, 2008

Aug. 28, 2008

====================================================================== ERROR: Failure: ImportError (No module named my_module) ---------------------------------------------------------------------- Traceback (most recent call last): File "/data/home/nwagner/local/lib/python2.5/site-packages/nose-0.10.1-py2.5.egg/nose/loader.py", line 364, in loadTestsFromName addr.filename, addr.module) File "/data/home/nwagner/local/lib/python2.5/site-packages/nose-0.10.1-py2.5.egg/nose/importer.py", line 39, in importFromPath return self.importFromDir(dir_path, fqname) File "/data/home/nwagner/local/lib/python2.5/site-packages/nose-0.10.1-py2.5.egg/nose/importer.py", line 84, in importFromDir mod = load_module(part_fqname, fh, filename, desc) File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/doc/__init__.py", line 10, in <module> __import__(__name__ + '.' + f) File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/doc/example.py", line 35, in <module> from my_module import my_func, other_func ImportError: No module named my_module

3 5

Re: [Numpy-discussion] rc1 update
by Travis E. Oliphant Aug. 28, 2008

Aug. 28, 2008

Andrew Straw wrote: > Joining the party late, here, so I don't know how much has already been > said... Particularly in favor of renaming NPY_VERSION. What's the benefit? > > I'm -1 on removing the name NPY_VERSION. It would cause unnecessary > breakage. > > I'm -0 on defining an additional NPY_ABI_VERSION to equal NPY_VERSION. > It just seems unnecessary to me. Won't a comment do, e.g. /* We wanted > to name this NPY_ABI_VERSION, but it was already NPY_VERSION since pre > 1.0, so we're sticking with it. (Until 1.2, it also served a dual role > as NPY_API_VERSION, when we added that, below.) */ ? > > I'm +0 on renaming the recently added NPY_FEATURE_VERSION to > NPY_API_VERSION. I don't care what it's called, really. > > Since A) this is a C developer thing only (aka people who should know > what they're doing), B) there is already IMO fairly good documentation > immediately above the NPY_VERSION and NPY_FEATURE_VERSION in the header > file, and C) very little impetus for anyone to actually directly use > NPY_VERSION anyway (there's a runtime ABI check during import_array() > which is the most important thing -- it checks at runtime if the compile > time ABI version of the caller and numpy itself are the same), I think > it's best to leave NPY_VERSION as NPY_VERSION. > > The only code that I know of that currently explicitly uses NPY_VERSION > is C preprocessor tests for which C API to use -- from now, new code > being written will check the new NPY_FEATURE_VERSION (or NPY_API_VERSION > or whatever it ends up being called) for API features. Anyone updating > old code will already have a bit of additional complexity to deal with, > and I want to minimize the amount of additional complexity. We don't > want a rat's nest of C preprossor #ifdefs as people track the numpy > version checking defines. Given the addition of NPY_FEATURE_VERSION, > code that's already been written with explicit NPY_VERSION checks for > API information is already going to become outdated, but as long as we > keep the NPY_VERSION name, at least it won't break or require another > level of complexity to maintain forward compatibility. > > +1 on what Andrew said. -Travis

4 3

Array subclassing - question for Travis
by Matthew Brett Aug. 28, 2008

Aug. 28, 2008

Hi Travis and team, I am just writing some docs for subclassing, and ran into some behavior I didn't understand: In [143]: class A(np.ndarray): pass In [144]: arr = np.arange(5) In [145]: obj = arr.copy().view(A) In [146]: type(obj) Out[146]: <class '__main__.A'> In [147]: obj.__array_priority__ # shouldn't this be 1.0 by default (from the numpy book)? Out[147]: 0.0 In [148]: type(np.multiply(arr, obj)) # this is what I expected Out[148]: <class '__main__.A'> In [149]: type(np.multiply.outer(arr, obj)) # this is not - I expected class A again Out[149]: <type 'numpy.ndarray'> In [151]: class A(np.ndarray): __array_priority__ = 20.0 # setting array priority does not affect behavior In [152]: obj = arr.copy().view(A) In [153]: type(np.multiply.outer(arr, obj)) Out[153]: <type 'numpy.ndarray'> Is this what you would expect? Thanks for clarification, and apologies if I misunderstood... Matthew

2 3

numpy trunk build broken on 32-bit ubuntu
by Fernando Perez Aug. 28, 2008

Aug. 28, 2008

Howdy, building numpy from trunk right now (r5708) on a 32-bit ubuntu box gives: compiling C sources C compiler: gcc -pthread -fno-strict-aliasing -DNDEBUG -g -fwrapv -O2 -Wall -Wstrict-prototypes -fPIC creating build/temp.linux-i686-2.5/build creating build/temp.linux-i686-2.5/build/src.linux-i686-2.5 creating build/temp.linux-i686-2.5/build/src.linux-i686-2.5/numpy creating build/temp.linux-i686-2.5/build/src.linux-i686-2.5/numpy/core creating build/temp.linux-i686-2.5/build/src.linux-i686-2.5/numpy/core/src compile options: '-Ibuild/src.linux-i686-2.5/numpy/core/src -Inumpy/core/include -Ibuild/src.linux-i686-2.5/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.5 -c' gcc: build/src.linux-i686-2.5/numpy/core/src/umathmodule.c numpy/core/src/umathmodule.c.src:332: error: static declaration of 'trunc' follows non-static declaration numpy/core/src/umathmodule.c.src: In function 'trunc': numpy/core/src/umathmodule.c.src:333: warning: unused variable 'r' numpy/core/src/umathmodule.c.src:333: warning: unused variable 'y' numpy/core/src/umathmodule.c.src:332: error: static declaration of 'trunc' follows non-static declaration numpy/core/src/umathmodule.c.src: In function 'trunc': numpy/core/src/umathmodule.c.src:333: warning: unused variable 'r' numpy/core/src/umathmodule.c.src:333: warning: unused variable 'y' error: Command "gcc -pthread -fno-strict-aliasing -DNDEBUG -g -fwrapv -O2 -Wall -Wstrict-prototypes -fPIC -Ibuild/src.linux-i686-2.5/numpy/core/src -Inumpy/core/include -Ibuild/src.linux-i686-2.5/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.5 -c build/src.linux-i686-2.5/numpy/core/src/umathmodule.c -o build/temp.linux-i686-2.5/build/src.linux-i686-2.5/numpy/core/src/umathmodule.o" failed with exit status 1 Notes (Matthew Brett found all of this): - on 64-bit fedora it builds fine - calling a similar test on 32-bit ubuntu with "gcc -std=c99" works. Perhaps the build is using something similar... - Not that I want to blame anyone, but I suspect this has something to do with David C's recent checkins: ------------------------------------------------------------------------ r5699 | cdavid | 2008-08-25 15:05:29 -0700 (Mon, 25 Aug 2008) | 2 lines Add a trunc function in umath module. ------------------------------------------------------------------------ r5698 | cdavid | 2008-08-25 15:05:21 -0700 (Mon, 25 Aug 2008) | 2 lines Add pure C trunc function implementations for platform which do not have it. as that is where the errror is coming from. Cheers, f

5 4

subsampling array without loops
by Catherine Moroney Aug. 28, 2008

Aug. 28, 2008

I'm looking for a way to acccomplish the following task without lots of loops involved, which are really slowing down my code. I have a 128x512 array which I want to break down into 2x2 squares. Then, for each 2x2 square I want to do some simple calculations such as finding the maximum value, which I then store in a 64x256 array. Here is the actual code involved. It's only slightly more complex than what I described above, since it also involves doing a bit of masking on the 2x2 sub-arrays. Any hints are much appreciated! An inordinate amount of time is being spent in this function and another one like it. Catherine def calc_sdcm_at_rlra(self,iblock): npixl = self.nlineht/self.nlinerl npixs = self.nsmpht/self.nsmprl sdcm_out = numpy.array([Constants.CLOUDMASK_NR] *self.nlinerl*self.nsmprl,'int8') \ .reshape(self.nlinerl,self.nsmprl) for ilinerl in range(0,self.nlinerl): for ismprl in range(0,self.nsmprl): height = self.data[iblock].height[ilinerl*2:ilinerl*2 +2, ismprl*2:ismprl*2+2] sdcm = self.data[iblock].sdcm[ilinerl*2:ilinerl*2 +2, ismprl*2:ismprl*2+2] source = self.data[iblock].heightsrc [ilinerl*2:ilinerl*2+2, ismprl*2:ismprl*2+2] mask1 = (source == Constants.HEIGHT_STEREO) mask2 = ( (source == Constants.HEIGHT_SURFACE) | \ (source == Constants.HEIGHT_DEFAULT) ) if (mask1.any()): loc = height[mask1].argmax() sdcm_out[ilinerl,ismprl] = sdcm[mask1].ravel()[loc] elif (mask2.any()): loc = height[mask2].argmax() sdcm_out[ilinerl,ismprl] = sdcm[mask2].ravel()[loc] return sdcm_out

4 3

Re: [Numpy-discussion] [SciPy08] Documentation BoF
by Neil Crighton Aug. 28, 2008

Aug. 28, 2008

> > - Should we have a separate User manual and a Reference manual, or only > a single manual? > Are there still plans to write a 10 page 'Getting started with NumPy' document? I think this would be very useful. Ideally a 'getting started' document, the docstrings, and a reference manual is all the documentation you'd need. I try to avoid reading long user manuals unless I'm forced to, but I don't know what other people do. Neil

2 1

normalizing a vector so it has magnitude 1
by SimonPalmer Aug. 27, 2008

Aug. 27, 2008

bit of a newb question, is there a method for normalising a 1D vector so it ends up with magnitude 1? I can do it manually but I was hoping there was a neat numpy - or scipy - trick. I've been web surfing but nothing really leaps out

2 2

About asarray (list)
by Claude Gouedard Aug. 27, 2008

Aug. 27, 2008

Ok , The same for asarray(1) .. The problem is that aa=asarray(1) is an numpy.array (right ? ) with a size 1 and a shape ( ) ! No surprising ? Claude

3 2

About asarray (list)
by Claude Gouedard Aug. 27, 2008

Aug. 27, 2008

Hi , I'm just surprised by the behaviour of numpy.asarray on lists. Can someone comment this : ===================== a=(1) aa=asarray(a) print aa.size , aa.shape >> 1 ( ) ===================== The shape doesnot reflect the actual size. If a=(1,2) there is no problem . ===================== aa=asarray((1,2)) print aa.size , aa.shape >> 2 (2, ) ===================== ... same for sequences: ===================== aa=asarray( [1] ) print aa.size , aa.shape >> 1 (1, ) ===================== Thanks Claude Gouedard

4 3

converting array of timedelta to array of integers.
by Ryan Neve Aug. 26, 2008

Aug. 26, 2008

Apologies in advance if this is an obvious answer, but I'm new to most of this. My overall goal is to produce a contour plot of some irregular time series data. I've imported the data from mySQL into three arrays x,y,and z where x is an array of datetime.timedelta objects. I need an array of the values of x converted into an array of integers representing the number of minutes in the timedelta. I've tried a few things, and spent hours searching for examples, but to no avail. Once I get the x array as integers, I can use griddata()

2 2

Compile error w/ rev5703
by Pierre GM Aug. 26, 2008

Aug. 26, 2008

Here's the interesting part: -------------------- creating build/temp.linux-x86_64-2.4/build/src.linux-x86_64-2.4/numpy/core/src compile options: '-Ibuild/src.linux-x86_64-2.4/numpy/core/src -Inumpy/core/include -Ibuild/src.linux-x86_64-2.4/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.4 -c' x86_64-pc-linux-gnu-gcc: build/src.linux-x86_64-2.4/numpy/core/src/umathmodule.c numpy/core/src/umathmodule.c.src:332: error: static declaration of âtruncâ follows non-static declaration ---------------------- Here's some generic info about the system. gcc (GCC) 4.1.2 (Gentoo 4.1.2 p1.0.1) Copyright (C) 2006 Free Software Foundation, Inc. C compiler: x86_64-pc-linux-gnu-gcc -pthread -fno-strict-aliasing -DNDEBUG -fPIC compile options: '-c' x86_64-pc-linux-gnu-gcc: _configtest.c x86_64-pc-linux-gnu-gcc -pthread _configtest.o -L/usr/lib64 -lf77blas -lcblas -latlas -o _configtest ATLAS version 3.8.2 built by root on Thu Jun 19 14:40:27 EDT 2008: UNAME : Linux back2tux 2.6.22-suspend2-r1 #8 Fri Feb 15 17:00:17 EST 2008 x86_64 Mobile AMD Athlon(tm) 64 Processor 3400+ AuthenticAMD GNU/Linux INSTFLG : -1 0 -a 1 ARCHDEFS : -DATL_OS_Linux -DATL_ARCH_HAMMER -DATL_CPUMHZ=2205 -DATL_SSE2 -DATL_SSE1 -DATL_3DNow -DATL_USE64BITS -DATL_GAS_x8664 F2CDEFS : -DAdd_ -DF77_INTEGER=int -DStringSunStyle CACHEEDGE: 262144 F77 : gfortran, version GNU Fortran 95 (GCC) 4.1.2 (Gentoo 4.1.2 p1.0.1) F77FLAGS : -O2 -m64 SMC : x86_64-pc-linux-gnu-gcc, version x86_64-pc-linux-gnu-gcc (GCC) 4.1.2 (Gentoo 4.1.2 p1.0.1) SMCFLAGS : -march=k8 -O2 -pipe -m64 SKC : x86_64-pc-linux-gnu-gcc, version x86_64-pc-linux-gnu-gcc (GCC) 4.1.2 (Gentoo 4.1.2 p1.0.1) SKCFLAGS : -march=k8 -O2 -pipe -m64

1 0

[SciPy08] Documentation BoF
by Stéfan van der Walt Aug. 26, 2008

Aug. 26, 2008

Hi all, This is my personal recollection of the documentation BoF. Feel free to comment or correct the text below. Regards Stéfan Summary of the Documentation Birds-of-a-Feather Session ======================================================= Topics proposed for discussion ------------------------------ - Path towards getting SciPy documented using the online documentation editor - Intersphinx: linking Sphinx-using projects - Culture / methodology: what do we expect from code accepted into NumPy/SciPy? - Tools for reading documentation (Help 2.0) - The complete SciPy book -- how to get there - NumPy Documentation Standard: is it useful beyond NumPy? Conclusions reached ------------------- Documentation format ```````````````````` - Question: could we standardise the descriptions of data-types between docs of different projects? They should overlap to some extent. - If we could agree on some common standard (or at least a common subset of a format) across projects, IPython could parse and display docstrings with special markup. - Matplotlib should be involved in this conversation. They have already written a vast amount of documentation, and their needs may be different to those of the SciPy projects. Formal requirements ``````````````````` For functions going into NumPy and SciPy, we have a certain minimum documentation requirements. However, we'd rather take a somewhat liberal, open approach, and accept code with somwhat inferior docstrings, working with the author to improve them. This way, we remain accessible as a community, while striving to improve the quality of the code base. Marketing SciPy ``````````````` The current entry point for SciPy on the web (http://www.scipy.org) may be confusing to beginners. The web page should be redesigned, in a similar fashion as the new Sage and code.enthought.com pages, linking to the messy wiki behind it in a consistent fashion. Joe Harrington may be able to hire someone to work on this important aspect of SciPy marketing. Cross-project documentation searching ````````````````````````````````````` Ideally, one should be able to perform a documentation search across several projects. In order to achieve this, we should negotiate a common convention for installing Sphinx-generated documentation into a standard location. A javascript-backed webpage can then be provided to do a search, using the json indices, or some other means.

6 7

Can't get to numpy.org or scipy.org today...
by SimonPalmer Aug. 26, 2008

Aug. 26, 2008

...are they down?

2 1

numpy.core.numerictypes questions
by William Reade Aug. 26, 2008

Aug. 26, 2008

Hi I hope this is the right place for this question... if not, I apologise, and ask only that you please direct me somewhere more appropriate. Line 532 of numerictypes.py reads "_unicodesize = array('u','U1').itemsize", but _unicodesize itself does not appear to be referenced anywhere else. My questions are: * Does this line have any use? * If not, would it be possible to remove it from a future release? I ask because I'm working on a project to get CPython extensions (specifically numpy) working with IronPython, and this line is blocking me. If I comment it out, I can actually import numpy and use some features, but I'm not sure how important unicode arrays are: if they're not widely used, I'd be keen to ignore them for a while and work on, say, floating-point arrays. Best William

3 4

Converting audio data to array for further processing
by Dan Colesworthy Aug. 26, 2008

Aug. 26, 2008

Hello all, I need an efficient way to convert 24 bit signed audio data to a numpy array for further processing. The data will be in a .wav file, and can be recovered via the python wave module. At that point it is a byte string - likely in little endian order. Somewhere in that data will be a 1 second length tone of prescribed but varying amplitude and frequency which must be measured for amplitude, frequency, duration, thd% etc. Also need to make the conversion quickly - over relatively large amounts of data. Something like 5 seconds sampled at 32000 samples/second. This same job will have to be repeated > 300 times as part of a test module. Thanks in advance for any and all ideas. Dan Colesworthy Otologics LLC

3 2

Buildbot failure
by Alan McIntyre Aug. 25, 2008

Aug. 25, 2008

In case it hasn't been noted yet, three of the buildbots (Windows_XP_x86_64_MSVC, Linux_SPARC_64_Debian, Linux_SPARC_64_Debian_gcc4) are failing to build. The Linux builds have the error: gcc: build/src.linux-sparc64-2.4/numpy/core/src/umathmodule.c numpy/core/src/umathmodule.c.src:332: error: static declaration of 'trunc' follows non-static declaration The MSVC build appears to fail with build\src.win32-2.5\numpy\core\include\numpy\__umath_generated.c(580) : error C2065: 'truncf' : undeclared identifier

2 1

Re: [Numpy-discussion] [SciPy08] Documentation BoF
by Sebastien Binet Aug. 25, 2008

Aug. 25, 2008

On Monday 25 August 2008 10:36:03 Travis E. Oliphant wrote: > Pauli Virtanen wrote: > > - I vaguely remember someone mentioned having done numpybook -> RST > > conversion. If so, is the result somewhere available? > > Was something done towards this in the Scipy'08 sprints? > > Yes, Gabriel Gellner made progress on this during the sprints, and I > asked him to post his scripts for doing it. I don't have his email > address currently, but perhaps he will respond. Note that a set of latex->reST converters seem to already exist there: http://svn.python.org/view/doctools/converter/ (haven't tried myself, so, YMMV) cheers, sebastien. -- ################################### # Dr. Sebastien Binet # Lawrence Berkeley National Lab. # 1 Cyclotron Road # Berkeley, CA 94720 ###################################

1 0

Re: [Numpy-discussion] subsampling array without loops
by Catherine Moroney Aug. 25, 2008

Aug. 25, 2008

This almost works. Is there a way to do some masking on tiles, for instance taking the maximum height of each 2x2 square that is an odd number? I've tried playing around with masking and where, but they don't return an array of the original size and shape of "tiles" below. Catherine > Perhaps you could take the relevant slices of the arrays with > appropriate striding: > > a = numpy.arange(128*512).reshape((128, 512)) > top_left = a[::2, ::2] > top_right = a[::2, 1::2] > bottom_left = a[1::2, ::2] > bottom_right = a[1::2, 1::2] > > tiles = numpy.array([top_left, top_right, bottom_left, bottom_right]) > maxes = tiles.max(axis=0) > > Similarly if you want overlapping tiles, you could leave out the > final :2 in the slice specifications above. > > Zach > > > >> I'm looking for a way to acccomplish the following task without lots >> of loops involved, which are really slowing down my code. >> >> I have a 128x512 array which I want to break down into 2x2 squares. >> Then, for each 2x2 square I want to do some simple calculations >> such as finding the maximum value, which I then store in a 64x256 >> array. >> >> Here is the actual code involved. It's only slightly more complex >> than what I described above, since it also involves doing a bit of >> masking on the 2x2 sub-arrays. >> >> Any hints are much appreciated! An inordinate amount of time is >> being spent in this function and another one like it. >> >> Catherine >> >> def calc_sdcm_at_rlra(self,iblock): >> >> npixl = self.nlineht/self.nlinerl >> npixs = self.nsmpht/self.nsmprl >> >> sdcm_out = numpy.array([Constants.CLOUDMASK_NR] \ >> *self.nlinerl*self.nsmprl,'int8') \ >> .reshape(self.nlinerl,self.nsmprl) >> >> for ilinerl in range(0,self.nlinerl): >> for ismprl in range(0,self.nsmprl): >> >> height = self.data[iblock].height >> [ilinerl*2:ilinerl*2+2, \ >> ismprl*2:ismprl*2+2] >> sdcm = self.data[iblock].sdcm[ilinerl*2:ilinerl*2 >> +2, \ >> ismprl*2:ismprl*2+2] >> source = self.data[iblock].heightsrc >> [ilinerl*2:ilinerl*2+2, \ >> >> ismprl*2:ismprl*2+2] >> >> mask1 = (source == Constants.HEIGHT_STEREO) >> mask2 = ( (source == Constants.HEIGHT_SURFACE) | \ >> (source == Constants.HEIGHT_DEFAULT) ) >> >> if (mask1.any()): >> loc = height[mask1].argmax() >> sdcm_out[ilinerl,ismprl] = sdcm[mask1].ravel() >> [loc] >> elif (mask2.any()): >> loc = height[mask2].argmax() >> sdcm_out[ilinerl,ismprl] = sdcm[mask2].ravel() >> [loc] >> >> return sdcm_out >

3 2

Inplace dot and blas access
by Dan Goodman Aug. 25, 2008

Aug. 25, 2008

Hi all, I'm sorry for what is probably a very simple question but despite looking I didn't manage to find anything. I have some code (solving linear differential equation) which is doing the following operation lots of times: S[:] = dot(A,S)+C here S is a (M,N) matrix, A is (M,M) and C is (M,1). M is usually smallish and N is largish. My question is, can this be turned into an inplace operation to save memory and time? As far as I can see, although some array functions allow you to specify the destination of the output, dot always generates a new array. So I would like to have something like this: dot(A,S, out=S) S+=C Ideally, I could do this in one line. I believe the BLAS DGEMM routine does what I want (if I converted C to a matrix rather than a vector). Is there any way to access this using Python and NumPy (e.g. something like numpy.blas.DGEMM) or would I need to use C code to do this? Many thanks in advance, Dan Goodman

1 0

doing zillions of 3x3 determinants fast
by Daniel Lenski Aug. 25, 2008

Aug. 25, 2008

Hi all, I need to take the determinants of a large number of 3x3 matrices, in order to determine for each of N points, in which of M tetrahedral cells they lie. I arrange the matrices in an ndarray of shape (N,M,5,3,3). As far as I can tell, Numpy doesn't have a function to do determinants over a specific set of axes... so I can't say: dets = np.linalg.det(a, axis1=-1, axis2=-1) Using an explicit Python loop is very slow... doing the determinants of 100,000 random 3x3 matrices in a loop and inserting the results into a preallocated results array takes about 5.1 seconds. So instead I coded up my own routine to do the determinants over the last 2 axes, based on the naive textbook formula for a 3x3 determinant: def det3(ar): a=ar[...,0,0]; b=ar[...,0,1]; c=ar[...,0,2] d=ar[...,1,0]; e=ar[...,1,1]; f=ar[...,1,2] g=ar[...,2,0]; h=ar[...,2,1]; i=ar[...,2,2] return (a*e*i+b*f*g+c*d*h)-(g*e*c+h*f*a+i*d*b) This is *much much* faster... it does the same 100,000 determinants in 0.07 seconds. And the results differ from linalg.det's results by less than 1 part in 10^15 on average (w/ a std dev of about 1 part in 10^13). Does anyone know of any well-written routines to do lots and lots of determinants of a fixed size? My routine has a couple problems I can see: * it's fast enough for 100,000 determinants, but it bogs due to all the temporary arrays when I try to do 1,000,000 determinants (=72 MB array) * I've heard the numerical stability of the naive determinant algorithms is fairly poor, so I'm reluctant to use this function on real data. Any advice will be appreciated! Dan

4 4

Report from SciPy
by Travis E. Oliphant Aug. 25, 2008

Aug. 25, 2008

Hi everybody, Robert K, Chuck H, Stefan VdW, Jarrod M, David C, and I had a nice discussion about the future directions of NumPy. We resolved some things and would like community feedback on them if there are opinions. * we will be moving to time-based releases (at least 2 times a year -- November / May) with major changes not accepted about 4 weeks before the release. * The releases will be numbered major.minor.bugfix * There will be no ABI changes in minor releases * There will be no API changes in bugfix releases * Any API changes in minor releases will be done in a backward compatible fashion (possibly with deprecations). * Thus 1.2 will not break ABI compatibility but will add new API features. * We will push the generalized ufuncs off to 1.3 * NumPy 2.0 will be a library and will not automagically import numpy.fft * We will suggest that other libraries use from numpy import fft instead of import numpy as np; np.fft * We will apply most of Andrew Dalke's speed up patches but will keep ctypeslib import * We will remove automatic import of numpy.doc from trunk * NumPy 1.2 will come out shortly after the conference (rc1 on Sunday). * SciPy 0.7b1 will come out after the sprints. If there is anything else I missed from those who were there, please let us all know. By the way, as promised, the NumPy book is now available for download and the source to the book is checked in to the numpy SVN tree: http://svn.scipy.org/svn/numpy/trunk/numpy/doc/numpybook Best regards, -Travis

10 10

Remove ParametricTestCase from numpy.testing
by Stéfan van der Walt Aug. 25, 2008

Aug. 25, 2008

Hey all, Now that we have switched to Nose as a test framework, we should look into deprecating ParametricTestCase. Alan, have you looked into this before? Cheers Stéfan

3 4

weights parameter of np.average() doesn't work (BUG?)
by Dan Lenski Aug. 25, 2008

Aug. 25, 2008

Hi all, Is there a good reason why the weights parameter of np.average() doesn't broadcast properly? This is with the Ubuntu Hardy x86_64 numpy package, version 1.0.4. In [293]: a=arange(100).reshape(10,10) # Things work fine when weights have the exact same shape as a In [297]: average(a, axis=1, weights=ones((10,10))) Out[297]: array([ 4.5, 14.5, 24.5, 34.5, 44.5, 54.5, 64.5, 74.5, 84.5, 94.5]) # Bizarre and incorrect result with length-10 weight array In [298]: average(a, axis=1, weights=ones(10)) Out[298]: array([[[[[[[[[ 0., 1., 2., 3., 4., 5., 6., 7., 8., 9.], [ 10., 11., 12., 13., 14., 15., 16., 17., 18., 19.], [ 20., 21., 22., 23., 24., 25., 26., 27., 28., 29.], [ 30., 31., 32., 33., 34., 35., 36., 37., 38., 39.], [ 40., 41., 42., 43., 44., 45., 46., 47., 48., 49.], [ 50., 51., 52., 53., 54., 55., 56., 57., 58., 59.], [ 60., 61., 62., 63., 64., 65., 66., 67., 68., 69.], [ 70., 71., 72., 73., 74., 75., 76., 77., 78., 79.], [ 80., 81., 82., 83., 84., 85., 86., 87., 88., 89.], [ 90., 91., 92., 93., 94., 95., 96., 97., 98., 99.] ]]]]]]]]) Doing the weighted-sum explicitly works fine for me: In [311]: sum(a*ones(10), axis=-1)/sum(ones(10)) Out[311]: array([ 4.5, 14.5, 24.5, 34.5, 44.5, 54.5, 64.5, 74.5, 84.5, 94.5]) This seems like a bug, especially since average.__doc__ states that: If weights are given, result is: sum(a * weights,axis) / sum(weights,axis), where the weights must have a's shape or be 1D with length the size of a in the given axis. Integer weights are converted to Float. Not specifying weights is equivalent to specifying weights that are all 1. Frankly, I don't even see why weights is constrained to be 1D or the same shape as a... why not anything that's broadcastable to the same shape as a? Dan

3 3

How to replace a loop for a ufunc (to test optimizations)
by Travis E. Oliphant Aug. 23, 2008

Aug. 23, 2008

Yesterday at the SciPy conference I suggested it would be a good idea to add a function to replace one of the inner loops in a ufunc so that an extension module could test an optimized inner loop more easily. It turns out that this is such a good idea that it's actually already in the code. The UFUNC_API call that does this is called: PyUFunc_ReplaceLoopBySignature You would use it something like this static void _new_add_loopfunc(char **args, intp *dimensions, intp *steps, void *func) { /* my new faster code */ } Then in the extension module (either in the init or in a method call) something like this: PyUFuncObject *ufunc; PyUFuncGenericFunction oldloopfunc, dummy; int signature[3] /* This should be the number of arguments to the ufunc you are replacing */ PyObject *module, *dict; module = PyImport_ImportModule("umath"); dict = PyModule_GetDict(module); ufunc = PyDict_GetItemString(dict, "add") /* get the ufunc where you want to replace one of the loops */ signature[0] = NPY_DOUBLE; signature[1] = NPY_DOUBLE; signature[2] = NPY_DOUBLE; if (PyUFunc_ReplaceLoopBySignature(ufunc, _new_add_loopfunc, signature, &old1dloop) < 0) { PyErr_SetString(PyExc_RuntimeError, "problem replacing loop"); /* error return i.e. return NULL or return or return -1*/ } Then to restore the original function: if (PyUFunc_ReplaceLoopBySignature(ufunc, old1dloop, signature, &dummy) < 0) { /* error */ } Remember to use import_umath() in the module you are using the UFUNC_API in. Now, lot's of folks should be able to test optimizations. -Travis

1 0

sprint: IRC chat?
by Damian Eads Aug. 23, 2008

Aug. 23, 2008

Friends, Are we meeting over IRC chat? I'd like to help with the sprint but remotely. I have to leave LA today, unfortunately. Thanks! Damian

2 1

problems generating negative binomials
by Chris Fonnesbeck Aug. 23, 2008

Aug. 23, 2008

I'm trying to generate negative binomial random numbers where the n parameter is non integer (<1 in fact). This should be possible, as n only needs to be real. However, I get the following message: ValueError: n <= 0 I assume this is because some rounding is going on behind the scenes. The reason I am doing this is because I'm trying to use the alternative parameterisation of the NB, with parameters k,m for which p=k/(k+m). I noticed that this parameterisation does not exist over in the scipy distributions either. Has anyone else run into this problem? The latter parameterisation is far more common in ecology, because it is an overdispersed Poisson. Thanks, cf

1 0

generating float sequences
by Chris Fonnesbeck Aug. 22, 2008

Aug. 22, 2008

Is there any prospect for easily getting ranges of floats in numpy, rather than just integers? In R, for example, you can specify a decimal value for the step size. It would be great to be able to do the following: arange(0, 100, 0.1) for a sequence from 0 to 100 in steps of 0.1. It seems rather unpythonic to have to do this: arange(0, 1000)/10.

3 4

Revisiting numpy/scipy on 64 bit OSX
by Michael Abshoff Aug. 22, 2008

Aug. 22, 2008

Hi, I have been playing with 64 bit numpy/scipy on OSX 10.5 Intel again. I thought everything worked after the last time we discussed this, but I noticed that I had Scipy import failures in for example optimize since the gfortran I used creates 32 bit code per default. So I build a gcc 4.2.4 on OSX, then build python in 64 bit mode (fixing configure.in since Python assumes flags only available with the Apple gcc on Darwin), added a wrapper around gfortran that injects a "-m64" compile flag. Then I checkout out and build numpy (r5671), scipy (r4661) and the current nose release and ran the tests. I had some more failures than I expected (details below). Any thoughts on the failures? I am about to build ATLAS and a couple other dependencies on that box. William Stein is at Scipy 08 and there seems to be at least some interest in 64 bit OSX binaries. If anyone wants one let me know and I can put one together once Scipy 0.7 and Numpy 1.1.2 are out. Cheers, Michael numpy failed tests: **************************************************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_polynomial.py", line 29, in test_polynomial Failed example: p / q Expected: (poly1d([ 0.33333333]), poly1d([ 1.33333333, 2.66666667])) Got: (poly1d([ 0.333]), poly1d([ 1.333, 2.667])) ********************************************************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_polynomial.py", line 51, in test_polynomial Failed example: p.integ() Expected: poly1d([ 0.33333333, 1. , 3. , 0. ]) Got: poly1d([ 0.333, 1. , 3. , 0. ]) **********************************************************************File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_polynomial.py", line 53, in test_polynomialFailed example: p.integ(1) Expected: poly1d([ 0.33333333, 1. , 3. , 0. ]) Got: poly1d([ 0.333, 1. , 3. , 0. ]) ********************************************************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_polynomial.py", line 55, in test_polynomial Failed example: p.integ(5) Expected: poly1d([ 0.00039683, 0.00277778, 0.025 , 0. , 0. , 0. , 0. , 0. ]) Got: poly1d([ 0. , 0.003, 0.025, 0. , 0. , 0. , 0. , 0. ]) ******************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 48, in test_ufunclike Failed example: U.log2(a) Expected: array([ 2.169925 , 1.20163386, 2.70043972]) Got: array([ 2.17 , 1.202, 2.7 ]) ********************************************************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 53, in test_ufunclike Failed example: U.log2(a, y) Expected: array([ 2.169925 , 1.20163386, 2.70043972]) Got: array([ 2.17 , 1.202, 2.7 ]) ********************************************************************** File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/lib/tests/test_ufunclike.py", line 55, in test_ufunclike Failed example: y Expected: array([ 2.169925 , 1.20163386, 2.70043972]) Got: array([ 2.17 , 1.202, 2.7 ]) ====================================================================== FAIL: test_umath.TestComplexFunctions.test_it ---------------------------------------------------------------------- Traceback (most recent call last): File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/nose/case.py", line 182, in runTest self.test(*self.arg) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/core/tests/test_umath.py", line 196, in test_it assert_almost_equal(fz.real, fr, err_msg='real part %s'%f) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/testing/utils.py", line 207, in assert_almost_equal assert round(abs(desired - actual),decimal) == 0, msg AssertionError: Items are not equal: real part <ufunc 'log10'> ACTUAL: -0.3010299956639812 DESIRED: 0.0 scipy: ====================================================================== ERROR: test_nonsymmetric_modes (test_arpack.TestEigenNonSymmetric) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/tests/test_arpack.py", line 204, in test_nonsymmetric_modes self.eval_evec(m,typ,k,which) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/tests/test_arpack.py", line 186, in eval_evec eval,evec=eigen(a,k,which=which,**kwds) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/arpack.py", line 220, in eigen raise RuntimeError("Error info=%d in arpack"%info) RuntimeError: Error info=-8 in arpack ====================================================================== ERROR: test_starting_vector (test_arpack.TestEigenNonSymmetric) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/tests/test_arpack.py", line 214, in test_starting_vector self.eval_evec(self.symmetric[0],typ,k,which='LM',v0=v0) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/tests/test_arpack.py", line 186, in eval_evec eval,evec=eigen(a,k,which=which,**kwds) File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/sparse/linalg/eigen/arpack/arpack.py", line 220, in eigen raise RuntimeError("Error info=%d in arpack"%info) RuntimeError: Error info=-8 in arpack ====================================================================== FAIL: test_implicit (test_odr.TestODR) ---------------------------------------------------------------------- Traceback (most recent call last): File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/scipy/odr/tests/test_odr.py", line 118, in test_implicit 2.0778813389755596e-03]]), File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/testing/utils.py", line 304, in assert_array_almost_equal header='Arrays are not almost equal') File "/Users/mabshoff/64bit_numpy_scipy/lib/python2.5/site-packages/numpy/testing/utils.py", line 289, in assert_array_compare assert cond, msg AssertionError: Arrays are not almost equal (mismatch 4.0%) x: array([[ 2.10892628e+00, -1.94376722e+00, 7.02635228e-02, -4.71752493e-02, 5.25155371e-02], [ -1.94376722e+00, 2.04814914e+00, -6.16004809e-02,... y: array([[ 2.10892746e+00, -1.94376864e+00, 7.02635509e-02, -4.71752674e-02, 5.25155759e-02], [ -1.94376864e+00, 2.04815092e+00, -6.16005159e-02,...

4 8

generalized eigenvector problem
by Manuel Metz Aug. 22, 2008

Aug. 22, 2008

Hi list, are there any plans to implement a routine to solve the "generalized eigenvector problem" as is done in matlab ? see http://www.mathworks.com/access/helpdesk/help/techdoc/ref/eig.html manuel

2 2

Problem with correlate?
by Hanno Klemm Aug. 22, 2008

Aug. 22, 2008

Hi All, after the discussion on numpy.correlate some time ago, regarding complex conjugation, etc. I today was pointed to yet another oddity, which I hope somebody could explain to me, as to why that's a feature, rather than a bug. I'm thoroughly confused by the following behaviour: In [29]: x = array([0.,0.,1, 0, 0]) In [30]: y = array([0, 0, 0, 1., 0, 0, 0]) In [31]: correlate(x,y) Out[31]: array([ 0., 1., 0.]) In [32]: correlate(x,y, mode='full') Out[32]: array([ 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0.]) In [33]: y = array([0,0,1.,0,0]) In [34]: correlate(x,y,mode='full') Out[34]: array([ 0., 0., 0., 0., 1., 0., 0., 0., 0.]) So far, everything is allright and as expected. Now: In [35]: y1 = array([1,0,0,0,0]) In [36]: correlate(x,y1,mode='full') Out[36]: array([ 0., 0., 0., 0., 0., 0., 1., 0., 0.]) In [37]: y2 = array([1,0,0,0,0,0,0]) In [38]: correlate(x,y2,mode='full') Out[38]: array([ 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0.]) Could somebody please help me understanding, why the one switches places? Furthermore, this behaviour does not seem to be very consistent, if we switch the x and y argument: In [39]: correlate(y2,x,mode='full') Out[39]: array([ 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0.]) In [40]: correlate(y1,x,mode='full') Out[40]: array([ 0., 0., 1., 0., 0., 0., 0., 0., 0.]) The answer remains the same. Any help would be appreciated, the numpy version is 1.0.4. Best regards, Hanno -- Hanno Klemm klemm(a)phys.ethz.ch

2 4

The "NumPy" Cython release + tutorial
by Dag Sverre Seljebotn Aug. 22, 2008

Aug. 22, 2008

Cython just had a release, and amongst the new features are efficient NumPy array indexing for integers, real floats and Python objects. You can get it at http://cython.org For those new to Cython, I've written a tutorial specifically targeted for NumPy users here: http://wiki.cython.org/tutorials/numpy ...documenting the process of getting from 1.8 seconds to 6 milliseconds for the same data by adding typing in a sample for-loop-intensive code snippet. The last section of that tutorial (Efficient indexing) should be useful for those who know Cython as well. -- Dag Sverre

4 3

first post, simple questions
by Prashant Saxena Aug. 21, 2008

Aug. 21, 2008

Hi, numpy rocks!!! import numpy linsp = numpy.linspace red = linsp(0, 255, 50) green = linsp(125, 150, 50) blue = linsp(175, 255, 50) array's elements are float. How do I convert them into integer? I need to build a new array from red, green, blue. like this: [[ red[0], green[0], blue[0] ], [ red[1], green[1], blue[1] ], [ red[2], green[3], blue[3] ], .... .... [ red[49], green[49], blue[49] ], ] how do I create such an array? Are there any issues regarding compilation of simple python scripts using numpy functions and cython? Thanks

5 4

numpy build fail under cygwin (Vista)
by Charles Doutriaux Aug. 20, 2008

Aug. 20, 2008

Hi both trunk and 1.1.1 fail to build under cygwin vista (latest version) I'm copy/pasting the end of the log C. copying numpy/doc/reference/__init__.py -> build/lib.cygwin-1.5.25-i686-2.5/numpy/doc/reference running build_ext customize UnixCCompiler customize UnixCCompiler using build_ext customize GnuFCompiler gnu: no Fortran 90 compiler found gnu: no Fortran 90 compiler found customize GnuFCompiler gnu: no Fortran 90 compiler found gnu: no Fortran 90 compiler found customize GnuFCompiler using build_ext building 'numpy.core.multiarray' extension compiling C sources C compiler: gcc -fno-strict-aliasing -DNDEBUG -g -O3 -Wall -Wstrict-prototypes creating build/temp.cygwin-1.5.25-i686-2.5 creating build/temp.cygwin-1.5.25-i686-2.5/numpy creating build/temp.cygwin-1.5.25-i686-2.5/numpy/core creating build/temp.cygwin-1.5.25-i686-2.5/numpy/core/src compile options: '-g -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/src -Inumpy/core/include -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.5 -I/usr/include/python2.5 -c' gcc: numpy/core/src/multiarraymodule.c gcc -shared -Wl,--enable-auto-image-base -g build/temp.cygwin-1.5.25-i686-2.5/numpy/core/src/multiarraymodule.o -L/usr/lib/python2.5/config -lpython2.5 -o build/lib.cygwin-1.5.25-i686-2.5/numpy/core/multiarray.dll building 'numpy.core.umath' extension compiling C sources C compiler: gcc -fno-strict-aliasing -DNDEBUG -g -O3 -Wall -Wstrict-prototypes creating build/temp.cygwin-1.5.25-i686-2.5/build creating build/temp.cygwin-1.5.25-i686-2.5/build/src.cygwin-1.5.25-i686-2.5 creating build/temp.cygwin-1.5.25-i686-2.5/build/src.cygwin-1.5.25-i686-2.5/numpy creating build/temp.cygwin-1.5.25-i686-2.5/build/src.cygwin-1.5.25-i686-2.5/numpy/core creating build/temp.cygwin-1.5.25-i686-2.5/build/src.cygwin-1.5.25-i686-2.5/numpy/core/src compile options: '-g -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/src -Inumpy/core/include -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.5 -I/usr/include/python2.5 -c' gcc: build/src.cygwin-1.5.25-i686-2.5/numpy/core/src/umathmodule.c /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s: Assembler messages: /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s:72843: Error: suffix or operands invalid for `fnstsw' /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s:73098: Error: suffix or operands invalid for `fnstsw' /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s: Assembler messages: /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s:72843: Error: suffix or operands invalid for `fnstsw' /cygdrive/c/Users/CHARLE~1/AppData/Local/Temp/ccib5Stl.s:73098: Error: suffix or operands invalid for `fnstsw' error: Command "gcc -fno-strict-aliasing -DNDEBUG -g -O3 -Wall -Wstrict-prototypes -g -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/src -Inumpy/core/include -Ibuild/src.cygwin-1.5.25-i686-2.5/numpy/core/include/numpy -Inumpy/core/src -Inumpy/core/include -I/usr/include/python2.5 -I/usr/include/python2.5 -c build/src.cygwin-1.5.25-i686-2.5/numpy/core/src/umathmodule.c -o build/temp.cygwin-1.5.25-i686-2.5/build/src.cygwin-1.5.25-i686-2.5/numpy/core/src/umathmodule.o" failed with exit status 1

2 4

Possible new multiplication operators for Python
by Fernando Perez Aug. 20, 2008

Aug. 20, 2008

Hi all, [ please keep all replies to this only on the numpy list. I'm cc'ing the scipy ones to make others aware of the topic, but do NOT reply on those lists so we can have an organized thread for future reference] In the Python-dev mailing lists, there were recently two threads regarding the possibility of adding to the language new multiplication operators (amongst others). This would allow one to define things like an element-wise and a matrix product for numpy arrays, for example: http://mail.python.org/pipermail/python-dev/2008-July/081508.html http://mail.python.org/pipermail/python-dev/2008-July/081551.html It turns out that there's an old pep on this issue: http://www.python.org/dev/peps/pep-0225/ which hasn't been ruled out, simply postponed. At this point it seems that there is room for some discussion, and obviously the input of the numpy/scipy crowd would be very welcome. I volunteered to host a BOF next week at scipy so we could collect feedback from those present, but it's important that those NOT present at the conference can equally voice their ideas/opinions. So I wanted to open this thread here to collect feedback. We'll then try to have the bof next week at the conference, and I'll summarize everything for python-dev. Obviously this doesn't mean that we'll get any changes in, but at least there's interest in discussing a topic that has been dear to everyone here. Cheers, f

12 33

global overloading of 1+1 -> MyClass(1, 1)
by Ondrej Certik Aug. 20, 2008

Aug. 20, 2008

Hi, with Andrew's permission, I am starting a new thread, where our discussion is ontopic. :) My original question was, that I would like to override 1+1 to return MyClass(1, 1) or something. Robert said it would break other libraries and Andrew said this: On Mon, Aug 18, 2008 at 9:23 PM, Andrew Dalke <dalke(a)dalkescientific.com> wrote: >> There are basically two possible options: >> >> 1) global overloading of + >> 2) allow to tell Python that "1" is not its int, but ours Integer. >> >> BTW, Sage just preparses Python for exactly this reason and >> substitutes all numbers like 1 with Integer(1). This really sucks >> imho. > > How would you like to do it? Any solution I can think of > would cause huge problems. For example, being able to change > int.__add__ would cause systemic problems throughout the > code base and would prevent, or very much hinder, migration > of Python code to C. Agree. > > Changes to support a different constructor for basic types, > on a per-module basis has its own problems. When are the > new types specified? In the module itself or by the code > which imports the module? Can the definitions be changed > on a scope-by-scope level? > > Now, I can imagine a Python with builtin multimethod > dispatch defined via static scoping, so that > > def __add__(left:int, right:int): > return __builtin__.__add__(__builtin__.__add__(left, right), 2) > > def __int__(s): > return Decimal(s) > > might work, but even here there's a problem because the "2" in > __add__ gets replaced by a Decimal, when I really want it to > be an integer. Not getting why is "2" replaced by Decimal if you don't want to? If you don't want it, you can just reimplement __int__(s). This would of course be only active in the module where it was defined. > > So I don't see how you can do it in the context of Python. > In the context of a very different programming language, > sure, it's possible. > > >> But that's what Sage is doing, they just preparse the input. But in >> SymPy we'd like to use it as a regular Python library in users >> programs, so imho preparsing, or custom compiling using your package >> is not really an option, is it? > > Use a different extension, like ".sympy" and an import hook > which does the .sympy -> .pyc conversion. > > But it's not one I recommend for real code, at least not > without a lot of personal experience to see if it's worthwhile > in the face of the downsides. Yes, this is a mess, this is just like preparsing. Well, not like -- this is preparsing. Ok, in the current state, you don't know either what's going to happen. If you write In [1]: x/2*3/4 you have no idea what the result is going to be, you need to analyze x.__div__() and start from there. But if you write In [2]: 1/2*3/4 currently you know it will be 0. But imho you could as well analyze the global __mul__ (or global __int__, depending on how this would be technically implemented) to see what's going to happen. I mean what is the difference between [1] and [2]? Ondrej

6 12

C-API change for 1.2
by Travis E. Oliphant Aug. 19, 2008

Aug. 19, 2008

Hi all, The 1.2 version of NumPy is going to be tagged. There is at least one change I'd like to add: The hasobject member of the PyArray_Descr structure should be renamed to "flags" and converted to a 32-bit integer. What does everybody think about this change? It should have minimal affect except to require a re-compile of extension modules using NumPy. The only people requiring code changes would be those making intimate use of the PyArray_Descr structure instead of using the macros. It's a simple change if there is no major opposition. -Travis

10 28

Possible new multiplication operators for Python
by Sebastien Loisel Aug. 19, 2008

Aug. 19, 2008

Hi, Just a quick note to let you guys know that I'm the one who most recently revived the discussion on python-dev. My needs are research and teaching in mathematics. So ideally, the programming notation should require only minimal explanations, and (equally importantly), it should be easy to set up your programming environment so that you can type your programs in. I like the PEP225 very much. I don't like @ operators because the @ is too large and I find it hard to read. ~ operators are fine (the dotted operators are fine too). I don't like unicode operators because they are hard to set up (I have no idea how), and also I'm wary of APL. Most of the unicode operators that were proposed on this list got rendered as indistinct blobs in my firefox. (Screenshot attached, I don't know if the mailing list will allow it -- if not, the picture is here: http://www.math.temple.edu/~loisel/unicode-operators.png) Sincerely, -- Sébastien Loisel

1 0

a good polygon class?
by mark Aug. 19, 2008

Aug. 19, 2008

Hello List - I am looking for a good polygon class. My main interest it to be able to figure out if a point is inside or outside the polygon, which can have any shape (well, as long as it is a polygon). Any suggestions? Thanks, Mark

4 3

indexing (compared to matlab)
by Brian Blais Aug. 19, 2008

Aug. 19, 2008

Hello, I wanted to do the following thing that I do in Matlab (on a bigger problem), setting the values of a part of a matrix with indexing: >> a=floor(rand(5,5)*10) % make an example matrix to work on a = 2 4 7 9 8 6 9 2 5 2 6 3 5 1 8 1 5 6 1 2 1 2 8 2 9 >> ind=[2,4] ind = 2 4 >> a(ind,ind)=a(ind,ind)+100 a = 2 4 7 9 8 6 109 2 105 2 6 3 5 1 8 1 105 6 101 2 1 2 8 2 9 =========================== In numpy, the same gives: In [11]:a=floor(random.rand(5,5)*10) In [14]:a Out[14]: array([[ 7., 7., 8., 1., 9.], [ 0., 4., 9., 0., 5.], [ 4., 3., 7., 8., 3.], [ 2., 0., 4., 2., 4.], [ 9., 5., 0., 9., 9.]]) In [15]:ind=[1,3] In [20]:a[ind,ind]+=100 In [21]:a Out[21]: array([[ 7., 7., 8., 1., 9.], [ 0., 104., 9., 0., 5.], [ 4., 3., 7., 8., 3.], [ 2., 0., 4., 102., 4.], [ 9., 5., 0., 9., 9.]]) which only replaces 2 values, not all the values in the row,col combinations of [1,1],[1,2],etc...[3,3] like matlab. Is there a preferred way to do this, which I think should be fairly common. If I know that the indices are regular (like a slice) is there a way to do this? thanks, Brian Blais -- Brian Blais bblais(a)bryant.edu http://web.bryant.edu/~bblais

4 7

NumPy 1.2.0b2 released
by Jarrod Millman Aug. 18, 2008

Aug. 18, 2008

Hey, NumPy 1.2.0b2 is now available. Please test this so that we can uncover any problems ASAP. SVN tag: http://svn.scipy.org/svn/numpy/tags/1.2.0b2 Mac binary: https://cirl.berkeley.edu/numpy/numpy-1.2.0b2-py2.5-macosx10.5.dmg Windows binary: http://www.enthought.com/~gvaroquaux/numpy-1.2.0b2-win32.zip Source tarball: https://cirl.berkeley.edu/numpy/numpy-1.2.0b2.tar.gz Thanks, -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

16 39

Re: [Numpy-discussion] Generalized ufuncs?
by Engel, Hans-Andreas Aug. 18, 2008

Aug. 18, 2008

"Charles R Harris" <charlesr.harris(a)gmail.com> wrote in message news:<e06186140808171929n705f347dj7b158162b580baef(a)mail.gmail.com>... On Sun, Aug 17, 2008 at 7:56 PM, Stéfan van der Walt <stefan(a)sun.ac.za>wrote: > 2008/8/17 Robert Kern <robert.kern(a)gmail.com>: > > I suggested that we move it to a branch for the time being so we can > > play with it and come up with examples of its use. > > That branch is here: > > http://stefan@svn.scipy.org/svn/numpy/branches/gen_ufuncs > > For an earlier thread about using vector valued ufuncs for sorts and such -- and negative reactions to the suggestion -- go here<http://thread.gmane.org/gmane.comp.python.numeric.general/20552/focus=20560>. One of the major objections was how to call such functions with the ufunc machinery and the needed machinery for type promotions, sub classes, and all that nonsense. Are these dealt with in the patch? The current numpy code for all that is a bit of a mess anyway, and it would be nice to figure out some unified interface to call through and clean up the current code in the process. In fact, I've been making some preliminary efforts in that direction by formatting the current code and working through it. Also, do we also want reduce methods and such? I suspect there is still a lot of work to do to get this whole thing up and running. Chuck ---------- The good news is that the patch just uses of the existing code to deal with all the tricky issues (this is why the patch is so short). By the way, sort could be implemented with the proposed specifications, its signature would be "(i)->(i)". I agree that it would be nice if that code could be made somewhat clearer; however, I think that this task is orthogonal to the generalized ufuncs patch, because there is no code overlap. The way the suggested implementation basically works is to remove the "core dimensions" from the input/output arrays, and then have the existing code handle all the intricacies over the "loop" dimensions. Reduce methods are currently not supported (an error is raised). Therefore, the current patch does not forestall anything and the desired functionality can be added whenever it is clear what would be best. I do not think that it would makes sense to specify/implement all possible extensions, optimizations, concrete ufuncs, morphing of existing numpy functions to ufuncs, etc. at once; presumably it is much better to start with a small but extremely flexible specification of generalized ufuncs first. Best, Hansres

3 2

Re: [Numpy-discussion] Generalized ufuncs?
by Engel, Hans-Andreas Aug. 18, 2008

Aug. 18, 2008

"Robert Kern" <robert.kern(a)gmail.com> wrote in message news:<3d375d730808172015t62739266hbee205b97bdc86c3(a)mail.gmail.com>... > On Sun, Aug 17, 2008 at 21:55, Anne Archibald <peridot.faceted(a)gmail.com> wrote: > > 2008/8/17 Robert Kern <robert.kern(a)gmail.com>: > >> > >> I suggested that we move it to a branch for the time being so we can > >> play with it and come up with examples of its use. If you have > >> examples that you have already written, I would love to see them. I, > >> for one, am amenable to seeing this in 1.2.0, but only if we push back > >> the release by at least a few weeks. > > > > This is not a worked example, but this is exactly what is needed to > > make possible the "arrays of matrices" functions that were discussed > > some time ago. For example, users frequently want to do something like > > multiply an array of matrices by an array of matrices; that is not > > currently feasible without a very large temporary array (or a loop). > > > > But I think you were looking for examples of code using the interface, > > to see whether it worked well. > > I'll take what I can get. In order of increasing utility: > > 1. Descriptions of use cases (as above). > 2. Mockups of the Python code demonstrating the use case (e.g. > nonfunctional Python code showing a potential generalized ufunc > with inputs and outputs). > 3. The C code implementing the actual functionality for the use case. > > -- > Robert Kern > > "I have come to believe that the whole world is an enigma, a harmless > enigma that is made terrible by our own mad attempt to interpret it as > though it had an underlying truth." > -- Umberto Eco Please find an example on "inner1d" in the following. 1. One use case for an inner function is described on http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions and http://thread.gmane.org/gmane.comp.python.numeric.general/17694. (Another one would be the "array of matrices" usage mentioned above; we have called this "dot2d" with signature "(m,n),(n,p)->(m,p)" on http://scipy.org/scipy/numpy/ticket/887: here the matrix multiplication would occur with respect to the last two dimensions.) 2. The mockup python code would be: >>> from numpy import * >>> N = 10 >>> a = random.randn(3, 5, N) >>> b = random.randn(5, N) >>> # standard inner function ... inner(a, b).shape (3, 5, 5) >>> # new ufunc "inner1d" with signature "(i),(i)->()", satisfying GeneralLoopingFunctions use case ... inner1d(a, b).shape (3, 5) 3. Concrete implementation of inner1d in C: /* * This implements the function out = inner1d(in1, in2) with * out[K] = sum_i { in1[K, i] * in2[K, i] } * and multi-index K, as described on * http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions * and on http://projects.scipy.org/scipy/numpy/ticket/887. */ static void DOUBLE_inner1d(char **args, intp *dimensions, intp *steps, void *func) { /* Standard ufunc loop length and strides. */ intp dn = dimensions[0]; intp s0 = steps[0]; intp s1 = steps[1]; intp s2 = steps[2]; intp n; /* Additional loop length and strides for dimension "i" in elementary function. */ intp di = dimensions[1]; intp i_s1 = steps[3]; intp i_s2 = steps[4]; intp i; /* Outer loop: equivalent to standard ufuncs */ for (n = 0; n < dn; n++, args[0] += s0, args[1] += s1, args[2] += s2) { char *ip1 = args[0], *ip2 = args[1], *op = args[2]; /* Implement elementary function: out = sum_i { in1[i] * in2[i] } */ npy_double sum = 0; for (i = 0; i < di; i++) { sum += (*(npy_double *)ip1) * (*(npy_double *)ip2); ip1 += i_s1; ip2 += i_s2; } *(npy_double *)op = sum; } } /* Actually create the ufunc object */ static PyUFuncGenericFunction inner1d_functions[] = { DOUBLE_inner1d }; static void *inner1d_data[] = { (void *)NULL }; static char inner1d_sigs[] = { PyArray_DOUBLE, PyArray_DOUBLE, PyArray_DOUBLE }; static void addInner1d(PyObject *dictionary) { PyObject *f = PyUFunc_FromFuncAndDataAndSignature(inner1d_functions, inner1d_data, inner1d_sigs, 1, 2, 1, PyUFunc_None, "inner1d", "inner on the last dimension and broadcast on the other dimensions", 0, "(i),(i)->()"); PyDict_SetItemString(dictionary, "inner1d", f); }

1 0

Please volunteer a Mac OS X buildbot slave
by Stéfan van der Walt Aug. 18, 2008

Aug. 18, 2008

Hi all, Due to hardware failure, we recently lost another of our build slaves. OSX is one of our primary deployment platforms, so I'd like to replace this one. If you have a mac available with good network connectivity, please consider volunteering it as a buildslave (setup should not take more than 15 minutes of your time). Thanks, Stéfan

3 4

Interrupting long running calculations ....
by Sebastian Haase Aug. 18, 2008

Aug. 18, 2008

Hi, Could someone remind me of the current state of numpy with regards to honoring KeyboardInterrupts !? I think KeyboardInterrupt is the exception that would have to be used to do this kind of thing - right !? E.g. by pressing Ctrl-C.... Another question, and not quite numpy specific, is how to generate this exception from another thread -- effectively "injecting" it into the calling stack !? I ask this, because I'm using wxPython, which makes things even more complicated .... :-( Thanks for any help, Sebastian Haase

2 1

Re: [Numpy-discussion] Generalized ufuncs?
by Engel, Hans-Andreas Aug. 18, 2008

Aug. 18, 2008

I am sorry that our submission http://projects.scipy.org/scipy/numpy/ticket/887 has created some annoyance; presumably we have taken the "Make contributions (e.g. code patches), (...) by submitting a 'ticket' on the Trac pages linked below" on http://scipy.org/Developer_Zone somewhat too literally. It was great to receive very positive responses to our patch and to receive a very timely review; this is encouraging for submitting code to the numpy repository. I would like to add that the proposed change is not that arbitrary; it is a well-established concept -- it is outlined on the GeneralLoopingFunctions wiki page, and it is a prominent concept of perl's PDL vector library. Of course, there is still room for arguing about details. The fact that no explicit "generalized" ufuncs are provided, should in my opinion not be an argument why not to include the change in the 1.2.0 release. Writing extension libraries that implement such generalized ufuncs, while being able to use the standard numpy distribution, would certainly be very valuable. Furthermore, the risk for including the proposed patch in 1.2.0 is very low: existing functionality is not touched. (Except the glitch we had by declaring variables in a gcc-way.) For standard ufuncs, it should be straightforward to see that the patch does not modify the behavior at all. I wish everyone a great SciPy'08 conference and very much hope that you can include the proposed functionality in the forthcoming NumPy release. Best regards, Hansres

5 7

[ANNOUNCE] EPD with Py2.5 v4.0.3001 Beta1 now available
by Dave Peterson Aug. 17, 2008

Aug. 17, 2008

Hello, Thanks to heroic efforts by Chris Galvan this weekend, and significant efforts by the team that finalized ETS 3.0.0 this week, we've been able to publish public beta releases of EPD with Py2.5 v4.0.30001 Beta1 for Windows and Mac OS X today. I've uploaded them to the downloads website and updated the EPD product pages to provide download links for the public. You can find the link to the betas here: http://www.enthought.com/products/epddownload.php Please give them a try and report any bugs to the EPD Trac site at https://svn.enthought.com/epd. In this release, EPD has been updated to include ETS 3.0.0, NumPy 1.1.1, IPython 0.9.beta, Matplotlib 0.98.1, Sphinx 0.4.2, pyhdf 0.8, VTK 5.0.4, wxPython 2.8.8.1, and many more updated projects. There are a few issues known at this time, but remember these are our first beta release of this version: * The included documentation hasn't been updated to the current versions of the third-party libraries. * Some of the product branding is not up-to-date with regard to the product name change to "EPD with Py2.5", nor with the version number of 4.0.30001 Beta 1. -- Dave

1 0

Re: [Numpy-discussion] Generalised ufuncs branch
by Engel, Hans-Andreas Aug. 17, 2008

Aug. 17, 2008

"Stéfan van der Walt" <stefan(a)sun.ac.za> wrote in message news:<9457e7c80808151645v457e155apa192bdbcf3f91d47(a)mail.gmail.com>... Hi all, I have moved the generalised ufuncs functionality off to a branch: http://svn.scipy.org/svn/numpy/branches/gen_ufuncs Please try it out and give us your feedback. We shall also pound on it at the sprint during SciPy'08, and thereafter decide how and when to best integrate it into NumPy. Thanks to Wenjie Fu and Hans-Andreas Engel for taking the time to think this issue through and to submit such a high-quality patch. Regards Stéfan ---------- Motivated by Stefan's and Chuck's off-line comments, I have attached an example implementation of a generalized ufunc: /* * This implements the function out = inner1d(in1, in2) with * out[K] = sum_i { in1[K, i] * in2[K, i] } * and multi-index K, as described on * http://scipy.org/scipy/numpy/wiki/GeneralLoopingFunctions * and on http://projects.scipy.org/scipy/numpy/ticket/887. */ static void DOUBLE_inner1d(char **args, intp *dimensions, intp *steps, void *func) { /* Standard ufunc loop length and strides. */ intp dn = dimensions[0]; intp s0 = steps[0]; intp s1 = steps[1]; intp s2 = steps[2]; intp n; /* Additional loop length and strides for dimension "i" in elementary function. */ intp di = dimensions[1]; intp i_s1 = steps[3]; intp i_s2 = steps[4]; intp i; /* Outer loop: equivalent to standard ufuncs */ for (n = 0; n < dn; n++, args[0] += s0, args[1] += s1, args[2] += s2) { char *ip1 = args[0], *ip2 = args[1], *op = args[2]; /* Implement elementary function: out = sum_i { in1[i] * in2[i] } */ npy_double sum = 0; for (i = 0; i < di; i++) { sum += (*(npy_double *)ip1) * (*(npy_double *)ip2); ip1 += i_s1; ip2 += i_s2; } *(npy_double *)op = sum; } } /* Actually create the ufunc object */ static PyUFuncGenericFunction inner1d_functions[] = { DOUBLE_inner1d }; static void *inner1d_data[] = { (void *)NULL }; static char inner1d_sigs[] = { PyArray_DOUBLE, PyArray_DOUBLE, PyArray_DOUBLE }; static void addInner1d(PyObject *dictionary) { PyObject *f = PyUFunc_FromFuncAndDataAndSignature(inner1d_functions, inner1d_data, inner1d_sigs, 1, 2, 1, PyUFunc_None, "inner1d", "inner on the last dimension and broadcast on the other dimensions", 0, "(i),(i)->()"); PyDict_SetItemString(dictionary, "inner1d", f); }

1 0

Re: [Numpy-discussion] Problem with the mailing list?
by Engel, Hans-Andreas Aug. 17, 2008

Aug. 17, 2008

"Charles R Harris" <charlesr.harris(a)gmail.com> wrote in message news:<e06186140808171654o468f2a8co94832ed6266dc6b4(a)mail.gmail.com>... > Hi All, > > I received an email from Hans Andreas -- the gen-ufuncs guy -- and he is > unable to post to the list even though subscribed. Anyone know what might be > the problem? Please cc Hans-Andreas.Engel(a)deshaw.com if you reply to this > post. > > Chuck > Thanks Chuck for following up on this -- it works now! (Initially I signed up a different address than the one I used for posting.)

1 0

Problem with the mailing list?
by Charles R Harris Aug. 17, 2008

Aug. 17, 2008

Hi All, I received an email from Hans Andreas -- the gen-ufuncs guy -- and he is unable to post to the list even though subscribed. Anyone know what might be the problem? Please cc Hans-Andreas.Engel(a)deshaw.com if you reply to this post. Chuck

1 0

long double woes on win32
by David Cournapeau Aug. 16, 2008

Aug. 16, 2008

Hi, The test failures on windows with 1.2b2 are due to buggy long double behavior of mingw. My understanding is that on windows, long double is 8 bytes (that's the sizeof value returned by VS 2003), but mingw says 12 bytes. One solution would be forcing numpy to configure itself to handle long double as double. But doing this breaks the configuration in many ways, and this would require relatively important changes which is quite fragile in numpy ATM (the whole math function things in umathmodule.c.src). I don't see this happening for 1.2. Another one is to use VS 2003 for the binaries: but then we are back to the problem of non buildable numpy with VS 2003. The last one is to ignore this for 1.2 (long double was already broken before in numpy win32 binaries). What do people think ? cheers, David

3 18

ANNOUNCE: ETS 3.0.0 released!
by Dave Peterson Aug. 16, 2008

Aug. 16, 2008

Hello, I'm pleased to announce that ETS 3.0.0 has just been tagged and released! Source distributions have been pushed to PyPi and over the next couple hours, Win32 and OSX binaries will also be uploaded to PyPi. This means you can install ETS, assuming you have the prereq software installed, via the simple command: easy_install ETS[nonets] Please see the Install page on our wiki for more detailed installation instructions: https://svn.enthought.com/enthought/wiki/Install Developers of ETS will find that the projects' trunks have already been bumped up to the next version numbers so a simple "ets up" (or svn up) should bring you up to date. Others may wish to grab a complete new checkout via a "ets co ETS". The release branches that had been created are now removed. The next release is currently expected to be ETS 3.0.1 -- Dave Enthought Tool Suite --------------------------- The Enthought Tool Suite (ETS) is a collection of components developed by Enthought and open source participants, which we use every day to construct custom scientific applications. It includes a wide variety of components, including: * an extensible application framework * application building blocks * 2-D and 3-D graphics libraries * scientific and math libraries * developer tools The cornerstone on which these tools rest is the Traits package, which provides explicit type declarations in Python; its features include initialization, validation, delegation, notification, and visualization of typed attributes. More information is available for all the packages within ETS from the Enthought Tool Suite development home page at http://code.enthought.com/projects/tool-suite.php. Testimonials ---------------- "I set out to rebuild an application in one week that had been developed over the last seven years (in C by generations of post-docs). Pyface and Traits were my cornerstones and I knew nothing about Pyface or Wx. It has been a hectic week. But here ... sits in front of me a nice application that does most of what it should. I think this has been a huge success. ... Thanks to the tools Enthought built, and thanks to the friendly support from people on the [enthought-dev] list, I have been able to build what I think is the best application so far. I have built similar applications (controlling cameras for imaging Bose-Einstein condensate) in C+MFC, Matlab, and C+labWindows, each time it has taken me at least four times longer to get to a result I regard as inferior. So I just wanted to say a big "thank you". Thank you to Enthought for providing this great software open-source. Thank you for everybody on the list for your replies." — Gaël Varoquaux, Laboratoire Charles Fabry, Institut d’Optique, Palaiseau, France "I'm currently writing a realtime data acquisition/display application … I'm using Enthought Tool Suite and Traits, and Chaco for display. IMHO, I think that in five years ETS/Traits will be the most comonly used framework for scientific applications." — Gary Pajer, Department of Chemistry, Biochemistry and Physics, Rider University, Lawrenceville NJ

1 0

C99 on windows
by David Cournapeau Aug. 16, 2008

Aug. 16, 2008

Hi, The current trunk has 14 failures on windows (with mingw). 12 of them are related to C99 (see ticket 869). Can the people involved in recent changes to complex functions take a look at it ? I think this is high priority for 1.2.0 thanks, David

4 5

[REVIEW] Update NumPy API format to support updates that don't break binary compatibility
by Stéfan van der Walt Aug. 16, 2008

Aug. 16, 2008

The current NumPy API number, stored as NPY_VERSION in the header files, needs to be incremented every time the NumPy C-API changes. The counter tells developers with exactly which revision of the API they are dealing. NumPy does some checking to make sure that it does not run against an old version of the API. Currently, we have no way of distinguishing between changes that break binary compatibility and those that don't. The proposed fix breaks the version number up into two counters -- one that gets increased when binary compatibility is broken, and another when the API is changed without breaking compatibility. Backward compatibility with packages such as Matplotlib is maintained by renaming NPY_VERSION to NPY_BINARY_VERSION. Please review the proposed change at http://codereview.appspot.com/2946 Regards Stéfan

2 2

Addition of a dict object to all NumPy objects
by Travis E. Oliphant Aug. 16, 2008

Aug. 16, 2008

Hello all, While we are on the subject of C-API changes, I've noticed that quite a few of the sub-classes of ndarray are constructed to basically add meta-information to the array. What if the base-class ndarray grew a dict object at it's end to hold meta information. Naturally, several questions arise: 1) How would you access the dictionary? (i.e. __dict__?) 2) Would attribute setting and getting retrieve from this dictionary (how are conflicts managed). * I think I would prefer a dict attribute on the numpy array that gets and sets into the dictionary. 3) Are the additional 4-8 bytes too expensive 4) Should there instead be a C-level array sub-class with the dict attached. 5) How should dict's be propagated through views? My preference is to not propagate them at all. Thoughts? -Travis

4 9

Different results from repeated calculation, part 2
by Keith Goodman Aug. 15, 2008

Aug. 15, 2008

I get slightly different results when I repeat a calculation. I've seen this problem before (it went away but has returned): http://projects.scipy.org/pipermail/numpy-discussion/2007-January/025724.ht… A unit test is attached. It contains three tests: In test1, I construct matrices x and y and then repeatedly calculate z = calc(x,y). The result z is the same every time. So this test passes. In test2, I construct matrices x and y each time before calculating z = calc(x,y). Sometimes z is slightly different. But the x's test to be equal and so do the y's. This test fails (on Debian Lenny, Core 2 Duo, with libatlas3gf-sse2 but not with libatlas3gf-sse). test3 is the same as test2 but I calculate z like this: z = calc(100*x,y) / (100 * 100). This test passes. I get: ====================================================================== FAIL: repeatability #2 ---------------------------------------------------------------------- Traceback (most recent call last): File "/home/[snip]/test/repeat_test.py", line 73, in test_repeat_2 self.assert_(result, msg) AssertionError: Max difference = 2.04946e-16 ---------------------------------------------------------------------- Should a unit test like this be added to numpy?

6 7

Generalised ufuncs branch
by Stéfan van der Walt Aug. 15, 2008

Aug. 15, 2008

Hi all, I have moved the generalised ufuncs functionality off to a branch: http://svn.scipy.org/svn/numpy/branches/gen_ufuncs Please try it out and give us your feedback. We shall also pound on it at the sprint during SciPy'08, and thereafter decide how and when to best integrate it into NumPy. Thanks to Wenjie Fu and Hans-Andreas Engel for taking the time to think this issue through and to submit such a high-quality patch. Regards Stéfan

1 0

Questions about some of the random functions
by Alan Jackson Aug. 15, 2008

Aug. 15, 2008

Just a question - I'm gradually working through the distributions for the documentation marathon and I realised that there is a whole nest of them named "standard-xxxx". For several (e.g., normal) they are just the regular distribution with all the parameters except size set to "standard" values. So my first question is - why? They seem very redundant. Second question - why is there is a "standard-t" for Student's T-test (or the distribution associated with it) but no corresponding "t" distribution. I probably have too much time to think about it this weekend. 8-) -- ----------------------------------------------------------------------- | Alan K. Jackson | To see a World in a Grain of Sand | | alan(a)ajackson.org | And a Heaven in a Wild Flower, | | www.ajackson.org | Hold Infinity in the palm of your hand | | Houston, Texas | And Eternity in an hour. - Blake | -----------------------------------------------------------------------

2 1

Re: [Numpy-discussion] min() of array containing NaN
by Joe Harrington Aug. 15, 2008

Aug. 15, 2008

> I'm doing nothing. Someone else must volunteer. Fair enough. Would the code be accepted if contributed? > There is a > reasonable design rule that if you have a boolean argument which you > expect to only be passed literal Trues and Falses, you should instead > just have two different functions. Robert, can you list some reasons to favor this design rule? Here are some reasons to favor richly functional routines: User's code is more readable because subtle differences affect args, not functions Easier learning for new users Much briefer and more readable docs Similar behavior across languages Smaller number of functions in the core package (a recent list topic) Many fewer routines to maintain, particularly if multiple switches exist Availability of the NaN functionality in a method of ndarray The last point is key. The NaN behavior is central to analyzing real data containing unavoidable bad values, which is the bread and butter of a substantial fraction of the user base. In the languages they're switching from, handling NaNs is just part of doing business, and is an option of every relevant routine; there's no need for redundant sets of routines. In contrast, numpy appears to consider data analysis to be secondary, somehow, to pure math, and takes the NaN functionality out of routines like min() and std(). This means it's not possible to use many ndarray methods. If we're ready to handle a NaN by returning it, why not enable the more useful behavior of ignoring it, at user discretion? --jh--

3 2

VS 2003 problems with cython-generated code
by David Cournapeau Aug. 15, 2008

Aug. 15, 2008

Hi, I noticed this morning that numpy 1.2 is not buildable with VS 2003 (the one you have to use for official python releases for at least python 2.5, and maybe 2.4). When we generate C code, both with internal code (numpy/core/code_generator) and with external tools (cython/pyrex for mtrand), the string literals generated for docstrings are too long for visual studio. We have to break them (e.g. "foo bar" becomes "foo""bar"), but doing so with cython-generated code is only doable by changing cython itself. So I did patch cython to break those, and regenerate the mtrand.c. This is done in the vs_longstring branch. Is it ok to put this for 1.2 ? Without it, I don't see a way to have numpy 1.2 buildable with VS. cheers, David P.S: I attached the necessary patches to cython bug tracker too, so that hopefully, the problem can be solved for a future version of cython.

2 2

Generalized ufuncs?
by Charles R Harris Aug. 15, 2008

Aug. 15, 2008

Stefan, I notice that you have merged some new ufunc infrastructure. I think these sort of things should be discussed and reviewed on the list before being committed. Could you explain what the purpose of these patches is? The commit messages are rather skimpy. Chuck

7 21

Re: [Numpy-discussion] min() of array containing NaN
by Joe Harrington Aug. 15, 2008

Aug. 15, 2008

> If you're willing to do arithmetic you might even be able to > pull it off, since NaNs tend to propagate: > if (new<min) min -= (min-new); > Whether the speed of this is worth its impenetrability I couldn't say. Code comments cure impenetrability, and have no cost in speed. One could write a paragraph explaining it (if it really needed that much). The comments could even reference the current discussion. --jh--

1 0

reading *big* inhomogenous text matrices *fast*?
by Dan Lenski Aug. 14, 2008

Aug. 14, 2008

Hi all, I'm using NumPy to read and process data from ASCII UCD files. This is a file format for describing unstructured finite-element meshes. Most of the file consists of rectangular, numerical text matrices, easily and efficiently read with loadtxt(). But there is one particularly nasty section that consists of matrices with variable numbers of columns, like this: # index property type nodes 1 1 tet 620 583 1578 1792 2 1 tet 656 551 553 566 3 1 tet 1565 766 1600 1646 4 1 tet 1545 631 1566 1665 5 1 hex 1531 1512 1559 1647 1648 1732 6 1 hex 777 1536 1556 1599 1601 1701 7 1 quad 296 1568 1535 1604 8 1 quad 54 711 285 666 As you might guess, the "type" label in the third column does indicate the number of following columns. Some of my files contain sections like this of *more than 1 million lines*, so I need to be able to read them fast. I have not yet come up with a good way to do this. What I do right now is I split them up into separate arrays based on the "type" label: lines = [f.next() for i in range(n)] lines = [l.split(None, 3) for l in lines] id, prop, types, nodes = apply(zip, lines) # THIS TAKES /FOREVER/ id = array(id, dtype=uint) prop = array(id, dtype=uint) types = array(types) cells = {} for t in N.unique(types): these = N.nonzero(types==t) # THIS NEXT LINE TAKES FOREVER these_nodes = array([nodes[ii].split() for ii in these], dtype=uint).T cells[t] = N.row_stack(( id[these], prop[these], these_nodes )) This is really pretty slow and sub-optimal. Has anyone developed a more efficient way to read arrays with variable numbers of columns??? Dan

7 12

min() of array containing NaN
by Thomas J. Duck Aug. 14, 2008

Aug. 14, 2008

Determining the minimum value of an array that contains NaN produces a surprising result: >>> x = numpy.array([0,1,2,numpy.nan,4,5,6]) >>> x.min() 4.0 I expected 0.0. Is this the intended behaviour or a bug? I am using numpy 1.1.1. Thanks, Tom -- Thomas J. Duck <tom.duck(a)dal.ca> Associate Professor, Department of Physics and Atmospheric Science, Dalhousie University, Halifax, Nova Scotia, Canada, B3H 3J5. Tel: (902)494-1456 | Fax: (902)494-5191 | Lab: (902)494-3813 Web: http://aolab.phys.dal.ca/

15 21

Mac OSX 4-way universal
by Kees, Christopher E Aug. 14, 2008

Aug. 14, 2008

Hi, I'm trying to install numpy 1.1.1 into a 4-way universal build of python 2.6. It builds fine, but I get an error in install saying I can't install when cross-compiling. Does anybody know how to get around this? Chris

7 9

unique1d returning indices
by Robert Cimrman Aug. 14, 2008

Aug. 14, 2008

Hi, due to popular demand, I have updated unique1d() to optionally return both kinds of indices: In [3]: b, i, j = nm.unique1d( a, return_index=True, return_inverse=True ) In [4]: a Out[4]: array([1, 1, 8, 3, 3, 5, 4]) In [6]: b Out[6]: array([1, 3, 4, 5, 8]) In [7]: a[i] Out[7]: array([1, 3, 4, 5, 8]) In [8]: b[j] Out[8]: array([1, 1, 8, 3, 3, 5, 4]) In [9]: a Out[9]: array([1, 1, 8, 3, 3, 5, 4]) Could someone needing this try the attached code? If it is ok I will update the tests, the docstring and commit it. Note also that the order of outputs has changed (previously unique1d() returned (i, b) for return_index=True). cheers, r.

2 7

non-linear array manipulation
by Gong, Shawn (Contractor) Aug. 13, 2008

Aug. 13, 2008

hi list, The following array manipulation takes long time because I can't find ways to do in row/column, and have to do cell by cell. Would you check to see if there is a nicer/faster way for this non-linear operation? for i in range(rows): for j in range(columns): a[i][j] = math.sqrt( max(0.0, a[i][j]*a[i][j] - b[j]*c[j]) ) thanks, Shawn

4 5

test errors with numpy-1.1.0
by Mathew Yeates Aug. 13, 2008

Aug. 13, 2008

On an AMD x86_64 with ATLAS installed I am getting errors like ValueError: On entry to DLASD0 parameter number 9 had an illegal value ValueError: On entry to ILAENV parameter number 2 had an illegal value Anybody seen this before? Mathew

2 1

Re: [Numpy-discussion] min() of array containing NaN
by Joe Harrington Aug. 13, 2008

Aug. 13, 2008

>On Tue, Aug 12, 2008 at 19:28, Charles R Harris ><charlesr.harris(a)gmail.com> wrote: >> >> >> On Tue, Aug 12, 2008 at 5:13 PM, Andrew Dalke <dalke(a)dalkescientific.com> >> wrote: >>> >>> On Aug 12, 2008, at 9:54 AM, Anne Archibald wrote: >>> > Er, is this actually a bug? I would instead consider the fact that >>> > np.min([]) raises an exception a bug of sorts - the identity of min is >>> > inf. >> >> <snip> >> >>> >>> Personally, I expect that if my array 'x' has a NaN then >>> min(x) must be a NaN. >> >> I suppose you could use >> >> min(a,b) = (abs(a - b) + a + b)/2 >> >> which would have that effect. >Or we could implement the inner loop of the minimum ufunc to return >NaN if there is a NaN. Currently it just compares the two values >(which causes the unpredictable results since having a NaN on either >side of the < is always False). I would be amenable to that provided >that the C isnan() call does not cause too much slowdown in the normal >case. While you're doing that, can you do it so that if keyword nan=False it returns NaN if NaNs exist, and if keyword nan=True it ignores NaNs? We can argue which should be the default (see my prior post). Both are compatible with the current undefined behavior. I assume that the fastest way to do it is two separate loops for the separate cases, but it might be fast enough straight (with a conditional in the inner loop), or with some other trick (macro magic, function pointer, whatever). Thanks, --jh--

3 2

C99 math funcs and windows
by David Cournapeau Aug. 13, 2008

Aug. 13, 2008

Hi, Many tests (14) are failing on windows because of the recent changes in some maths functions (C99 conformance). I believe the problem is caused by mingw compiler (which is stuck at 3.4.5), but I still have to check that. To solve this, I see three solutions: - disabling the C99 tests on windows. Easy - using a more recent version of mingw: mingw does not officially support gcc 4, but we can built it ourselves. Not too difficult (I have some makefiles to build a toolchain from scratch on linux/mac os x), but I don't like the idea of using our own compiler. - we don't care, tests fail on windows. Do some people have a strong opinion on that ? Disabling sounds like the way to go for 1.2.0, but at some point, we will have to deal with this C99 problem on windows if people add more C99 code in it. cheers, David

3 3

ravel() in ma/core.py
by Christoph T. Weidemann Aug. 13, 2008

Aug. 13, 2008

Hi! I'm working with a subclass of ndarray and ran into an issue that seems to be caused by a line in numpy/ma/core.py The offending line is no. 1837 in version 1.1.0 or 2053 in the latest SVN version (revision 5635): r = ndarray.ravel(self._data).view(type(self)) The problem is that my subclass of ndarray has its own implementation of the ravel method that takes care of some important things before calling the ndarray.ravel class. Calling ndarray.ravel directly on this data structure doesn't allow for the necessary house-keeping to take place and my data structure complains by throwing an exception.I believe this issue would be solved if the above line instead read: r = self.ravel() or alternatively r = np.ravel(self) or something along those lines -- can anybody here please take a look and see if these (or something else) would be a viable alternative that would cause the subclass's ravel function to be called instead of ndarray.ravel()? Thanks, CTW PS: I ran into this issue while trying to plot my data with the imshow function in matplotlib. If vmin and/or vmax are not specified, a masked array version of the input array is generated in matplotlib/colors.py in this line: val = ma.asarray(value).astype(np.float) Then ma.minimum() and ma.maximum() are called and eventually the line highlighted above gets executed. I don't know if/why this conversion to masked arrays is necessary in the matplotlib code, but it seems pretty clear that the source of the problem is the application of the ndarray.ravel() method rather than that for the subclass in the line highlighted above.

3 2

SSEPlus + Framewave
by Holger Rapp Aug. 13, 2008

Aug. 13, 2008

Hi everyone, I have a question concerning performance of numpy. I'm using it for heavy image processing tasks and often need more speed than the stock numpy delivers. Especially in numeric calculations with big arrays (in my current case shape = (8,4,1200,1600), 'float'). So I often rely on self written special modules which do use the IPP (Intel Performance Primitives) to speed up simple tasks like dividing a big array by a number. I realize that it was not really feasible to support a proprietary library like the IPP in a beautifully crafted Open Source Project, but quite recently, AMD came up with two very interesting projects (SSEPlus and Framewave, links provided below) which are more or less a direct response to intels IPP. And the best: They are OpenSource (under a Apache license, afaik). My question is now: Is it intended/Is there interest to get this performance gain into numpy? Are their any political restrictions (license/project identity)? Is there already work underway? I for one would consider helping in a effort like that, because it would probably safe me time in the long run. (Sidenote: I'm aware that this optimization would only help INTEL/AMD boxes, but hardware acceleration is so common these days that it is a shame NOT to use it in a numbercrunching library. Implementing such a library on one architecture might make it easier to implement something similar on others too with other libs. Maybe we see CUDA support in numpy in the future....) Here are the links: http://sseplus.sourceforge.net/ http://framewave.sourceforge.net/ Greetings, Holger

4 11

has f2py changed from 1.0.4 to 1.1.1 ????
by mark Aug. 13, 2008

Aug. 13, 2008

Hello all - I have been banging my head against the wall. f2py worked fine on windows under 1.0.4 (Python 2.5) After upgrading to numpy 1.1.1, f2py keeps complaining it cannot find msvccompiler, yet it is still in the same place (under distutils) as it was under 1.0.4. Did something change? I read something about some f2py branch that was ended in 1.1.1. Any suggestions on how to solve this? Thanks, Mark

3 4

Re: [Numpy-discussion] min() of array containing NaN
by Thomas J. Duck Aug. 12, 2008

Aug. 12, 2008

Christopher Barker wrote: > well, it's not a bug because the result if there is a NaN is > undefined. > However, it sure could trip people up. If you know there is likely > to be > a NaN in there, then you could use nanmin() or masked arrays. The > problem comes up when you have no idea there might be a NaN in > there, in > which case you get a bogus answer -- this is very bad. This is exactly what happened to me. I was getting crazy results when contour plotting with matplotlib, although the pcolor plots looked fine. In particular, the colorscale had incorrect limits. This led me to check the min() and max() values in my array, which were clearly wrong as illustrated by the pcolor plot. Further investigation revealed unexpected NaNs in my array. > Is there an error state that will trigger an error or warning in these > situations? Otherwise, I'd have to say that the default should be to > test for NaN's, and either raise an error or return NaN. If that > really > does slow things down too much, there could be a flag that lets you > turn > it off. It is quite often the case that NaNs are unexpected, so it would be helpful to raise an Exception. Thanks for all of the helpful discussion on this issue. -- Thomas J. Duck <tom.duck(a)dal.ca> Associate Professor, Department of Physics and Atmospheric Science, Dalhousie University, Halifax, Nova Scotia, Canada, B3H 3J5. Tel: (902)494-1456 | Fax: (902)494-5191 | Lab: (902)494-3813 Web: http://aolab.phys.dal.ca/

2 1

Re: [Numpy-discussion] min() of array containing NaN
by Joe Harrington Aug. 12, 2008

Aug. 12, 2008

> It really isn't very hard to replace > np.sum(A) > with > np.sum(A[~isnan(A)]) > if you want to ignore NaNs instead of propagating them. So I don't > feel a need for special code in sum() that treats NaN as 0. That's all well and good, until you want to set the axis= keyword. Then you're stuck with looping. As doing stats for each pixel column in a stack of astronomical images with bad pixels and cosmic-ray hits is one of the most common actions in astronomical data analysis, this is an issue for a significant number of current and future users. >>> a=np.arange(9, dtype=float) >>> a.shape=(3,3) >>> a[1,1]=np.nan >>> a array([[ 0. , 1. , 2. ], [ 3. , nan, 5. ], [ 6. , 7. , 8. ]]) >>> np.sum(a) nan >>> np.sum(a[~np.isnan(a)]) 32.0 Good, but... >>> np.sum(a[~np.isnan(a)], axis=1) Traceback (most recent call last): File "<stdin>", line 1, in <module> File "/usr/local/lib/python2.5/site-packages/numpy/core/fromnumeric.py", line 634, in sum return sum(axis, dtype, out) ValueError: axis(=1) out of bounds Uh-oh... >>> np.sum(a[~np.isnan(a)], axis=0) 32.0 Worse: wrong answer but not an exception, since >>> a[~np.isnan(a)] array([ 0., 1., 2., 3., 5., 6., 7., 8.]) has the undesired side effect of irreversibly flattening the array. --jh--

1 0

has f2py changed from 1.0.4 to 1.1.1 ????
by mark Aug. 12, 2008

Aug. 12, 2008

Hello all - I have been banging my head against the wall. f2py worked fine on windows under 1.0.4 (Python 2.5) After upgrading to numpy 1.1.1, f2py keeps complaining it cannot find msvccompiler, yet it is still in the same place (under distutils) as it was under 1.0.4. Did something change? I read something about some f2py branch that was ended in 1.1.1. Any suggestions on how to solve this? Thanks, Mark

1 0

masked_equal not commutative
by Charles Doutriaux Aug. 12, 2008

Aug. 12, 2008

As always as i clicked send i realized my error it is indeed not commutative and that makes sense but i'm not sure the case: numpy.ma.masked_equal(1,a) should have worked, since we don't really know how to do this comparison, the only thing that could make sense would be commutation but i thinks it's probably dangerous to assume what the user wants to do. C.

1 0

masked_equal not commutative?
by Charles Doutriaux Aug. 12, 2008

Aug. 12, 2008

Hi I'm using 1.1.1 and found that numpy.ma.masked_equal is not commutative! I would expect it to be in this case. Or raise an error for uncompatible shape in the first case, no ? >>> a = numpy.ma.arange(100) >>> a.shape=(10,10) >>> b=numpy.ma.masked_equal(1,a) >>> b Traceback (most recent call last): File "<stdin>", line 1, in <module> File "/lgm/cdat/latest/lib/python2.5/site-packages/numpy/ma/core.py", line 1691, in __repr__ 'data': str(self), File "/lgm/cdat/latest/lib/python2.5/site-packages/numpy/ma/core.py", line 1665, in __str__ res[m] = f IndexError: 0-d arrays can't be indexed. >>> b=numpy.ma.masked_equal(a,1) >>> b masked_array(data = [[0 -- 2 3 4 5 6 7 8 9] [10 11 12 13 14 15 16 17 18 19] [20 21 22 23 24 25 26 27 28 29] [30 31 32 33 34 35 36 37 38 39] [40 41 42 43 44 45 46 47 48 49] [50 51 52 53 54 55 56 57 58 59] [60 61 62 63 64 65 66 67 68 69] [70 71 72 73 74 75 76 77 78 79] [80 81 82 83 84 85 86 87 88 89] [90 91 92 93 94 95 96 97 98 99]], mask = [[False True False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False] [False False False False False False False False False False]], fill_value=999999)

1 0

Ticket #750 (Swig interface - Bug in FORTRAN-ordering support)
by Stéfan van der Walt Aug. 12, 2008

Aug. 12, 2008

Hi all, I've given ticket #750 a shot. Since I was not previously familiar with SWIG, I would appreciate it if those of you who know it could take a look at the proposed patch: http://codereview.appspot.com/2890 Thanks for your time, Stéfan

2 1

Re: [Numpy-discussion] min() of array containing NaN
by Joe Harrington Aug. 12, 2008

Aug. 12, 2008

Masked arrays are a bit clunky for something as simple and standard as NaN handling. They also have the inverse of the standard truth sense, at least as used in my field. 1 (or True) usually means the item is allowed, not denied, so that you can multiply the mask by the data to zero all bad values, add and subtract masks in sensible ways and get what's expected, etc. For example, in the "stacked, masked mean" image processing algorithm, you sum the data along an axis, sum the masks along that axis, and divide the results to get the mean image without bad pixels. This is much more accurate than taking a median, and admits to error analysis, which the median does not (easily). While the regular behavior is "just a ~ away", as Stefan pointed out to me once, that's not acceptable if the image cube is large and memory or speed are at issue, and it's also very prone to bugs if you're negating everything all the time. Further, with ma you have to convert to using an entirely different and redundant set of routines instead of having the very standard handling of NaNs found in our competitor programs, such as IDL. The issue of not having an in-place method in ma was also raised earlier. I'll add the difficulty of converting code if a standard thing like NaN handling has to be simulated in multiple calls. So, I endorse extending min() and all other statistical routines to handle NaNs, possibly with a switch to turn it on if a suitably fast algorithm cannot be found (which is competitor IDL's solution). Certainly without a switch the default behavior should be to return NaN, not to return some random value, if a NaN is present. Otherwise the user may never know a NaN is present, and therefore has to check every use for NaNs. That constand manual NaN checking is slower and more error-prone than any numerical speed advantage. So to sum, proposed for statistical routnes: if NaN is not present, return value if NaN is present, return NaN if NaN is present and nan=True, return value ignoring all NaNs OR: if NaN is not present, return value if NaN is present, return value ignoring all NaNs if NaN is present and nan=True, return NaN I'd prefer the latter. IDL does the former and it is a pain to do /nan all the time. However, the latter might trip up the unwary, whereas the former never does. This would apply at least to: min max sum prod mean median std and possibly many others. --jh--

3 2

ATTENTION: 1.2.0b1 tagged tomorrow
by Jarrod Millman Aug. 12, 2008

Aug. 12, 2008

Hello, The C-API change was committed in revision 5626: http://scipy.org/scipy/numpy/changeset/5626 Please note that this change requires recompiling extension modules. A draft of the release notes can be found here: http://scipy.org/scipy/numpy/milestone/1.2.0 Please take a look and let me know if I am missing anything. I will be tagging the 1.2.0b1 tomorrow. The beta release will be followed by a release candidate at the beginning of the SciPy conference with the final release planned for the end of the conference. Thanks, -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

1 1

are fortran extensions numpy version dependent?
by mark Aug. 11, 2008

Aug. 11, 2008

Hello list - Easy question: do fortran extension created with numpy 1.0.4 also work with numpy 1.1.1? And vice versa? Just wondering before I start updating my version (and maybe I am just too worried, but those fortran extesions are still magic to me), Mark

2 1

can I install superpack numpy 1.1.1 over an egg install of 1.0.4 egg?
by mark Aug. 11, 2008

Aug. 11, 2008

Hello - I downloaded the latest Enthought distribution, where everything is an egg. I needed to upgrade to the new numpy 1.1.1, but eggs don't seem available on sourceforge. Can I use the numpy 1.1.1 super pack installer (for windows) to install over the already available egg of 1.0.4? Any reason why there are no numpy eggs? I thought eggs are the future, Thanks, Mark

1 0

compress() signature
by Christoph T. Weidemann Aug. 10, 2008

Aug. 10, 2008

The signature for numpy.compress() seems somewhat inconsistent with other functions, because the first argument it expects is not the array but the condition. This tends to throw me off when I use this function, because the common order for other functions seems to be: 1) input array 2) function specific argument(s) if applicable 3) axis 4) out (if applicable) 5) optional function specific argument(s) if applicable I think it would be nice to change the signature of compress() to bring it in line with the other functions, but obviously such a change would break some code without added functionality. What do people here think?

2 1

ndarray subclassing bug?
by Christoph T. Weidemann Aug. 10, 2008

Aug. 10, 2008

I've come across a strange behavior for classes subclassed from ndarray. Here's a minimal example that illustrates the problem: import numpy as np class TestArray(np.ndarray): def __new__(cls, data, info=None, dtype=None, copy=False): subarr = np.array(data, dtype=dtype, copy=copy) subarr = subarr.view(cls) return subarr def sort(self,*args,**kwargs): print type(self) print type(self.base) Now consider this: In [1]: tst = TestArray(np.random.rand(2,3)) In [2]: tst.sort() <class '__main__.TestArray'> <type 'numpy.ndarray'> In [3]: np.sort(tst) <class '__main__.TestArray'> <type 'NoneType'> Out[3]: TestArray([[ 0.90489484, 0.950291 , 0.80753772], [ 0.49020689, 0.84582283, 0.61532922]]) Why whould tst.sort() show the correct base class and np.sort show NoneType as base class for tst? I'd appreciate any insights ...

2 6

New python module to simulate arbitrary fixed and infinite precision binary floating point
by Rob Clewley Aug. 10, 2008

Aug. 10, 2008

Dear Pythonistas, How many times have we seen posts recently along the lines of "why is it that 0.1 appears as 0.10000000000000001 in python?" that lead to posters being sent to the definition of the IEEE 754 standard and the decimal.py module? I am teaching an introductory numerical analysis class this fall, and I realized that the best way to teach this stuff is to be able to play with the representations directly, in particular to be able to see it in action on a simpler system than full 64-bit precision, especially when str(f) or repr(f) won't show *all* of the significant digits stored in a float. The decimal class deliberately avoids binary representation issues, and I can't find what I want online. Consequently, I have written a module to simulate the machine representation of binary floating point numbers and their arithmetic. Values can be of arbitrary fixed precision or infinite precision, along the same lines as python's in-built decimal class. The code is here: http://www2.gsu.edu/~matrhc/binary.html The design is loosely based on that decimal module, although it doesn't get in to threads, for instance. You can play with different IEEE 754 representations with different precisions and rounding modes, and compare with infinite precision Binary numbers. For instance, it is easy to learn about machine epsilon, representation/rounding error using a much simpler format such as a 4-bit exponent and 6-bit mantissa. Such a format is easily defined in the new module and can be manipulated easily: >>> context = define_context(4, 6, ROUND_DOWN) >>> zero = context(0) Binary("0", (4, 6, ROUND_DOWN)) >>> print zero # sign, characteristic, significand bits 0 0000 000000 >>> zero.next() Binary("0.001E-9", (4, 6, ROUND_DOWN)) >>> print zero.next() 0 0000 000001 >>> largest_denormalized = context('0 0000 111111') # direct spec of the sign, characteristic, and significand bits >>> largest_denormalized Binary("0.111111E-6", (4, 6, ROUND_DOWN)) >>> largest_denormalized.as_decimal() Decimal("0.015380859375") >>> n01 = context(0.1) # nearest representable is actually stored >>> print n01, " rounded to ", n01.as_decimal() 0 0011 100110 rounded to 0.099609375 >>> Binary('-10111.0000001').as_decimal() Decimal("-23.0078125") >>> Binary('-10111.0000001', context).as_decimal() # not enough precision in this context Decimal("-23.0000") >>> diff = abs(Binary('-10111.0000001') - Binary('-10111.0000001', context)) >>> diff Binary("0.1E-6", (4, 6, ROUND_DOWN)) The usual arithmetic operations are permitted on these objects, as well as representations of their values in decimal or binary form. Default contexts for half, single, double, and quadruple IEEE 754 precision floats are provided. Binary integer classes are also provided, and some other utility functions for converting between decimal and binary string representations. The module is compatible with the numpy float classes and requires numpy to be installed. The source code is released under the BSD license, but I am amenable to other licensing ideas if there is interest in adapting the code for some other purpose. Full details of the functionality and known issues are in the module's docstring, and many examples of usage are in the accompanying file binary_tests.py (which also acts to validate the common representations against the built-in floating point types). I look forward to hearing feedback, especially in case of bugs or suggestions for improvements. -Rob -- Robert H. Clewley, Ph. D. Assistant Professor Department of Mathematics and Statistics Georgia State University 720 COE, 30 Pryor St Atlanta, GA 30303, USA tel: 404-413-6420 fax: 404-413-6403 http://www2.gsu.edu/~matrhc http://brainsbehavior.gsu.edu/

1 0

numpy 1.1.rc2: win32 binaries
by David Cournapeau Aug. 10, 2008

Aug. 10, 2008

Hi, After some delay, here are the win32 binaries for numpy 1.1.1rc2: http://www.ar.media.kyoto-u.ac.jp/members/david/numpy-1.1.1.dev5559-win32-s… Notes on those binaries: - Based on Atlas 3.8.2 (the 1.1.0 was built against 3.8.0, which had a serious bug wrt dgemm, which is used for numpy.dot. It should solve #844 (problem with numpy.inner) - It is not stricly based on 1.1.rc2, but on 1.1.x trunk. The only difference, though, is a small fix in MANIFEST.in which was broken wrt sdist target and the version of course. If those work out, I will also prepare 2.4 binaries. I am sorry for the delay, cheers, David

7 16

packbits / unpackbits bugs or misfeatures?
by Zachary Pincus Aug. 9, 2008

Aug. 9, 2008

Hello all, numpy.unpackbits has a docstring that states that it returns a boolean array, but the function instead returns a uint8 array. Should I enter this in trac as a documentation bug or a functionality bug? Also, numpy.packbits will not accept a bool-typed array as input (only scalar-typed arrays) -- should it have this ability, for symmetry with unpackbits? (Assuming that unpackbits should indeed unpack into bool arrays...) Zach

2 2

Horizontal lines in diffraction image (NumPy FFT)
by Matthias Hillenbrand Aug. 9, 2008

Aug. 9, 2008

Hello, I am trying to calculate the propagation of a Gaussian beam with an angular spectrum propagation algorithm. My program does the following steps: 1. Generation and multiplication of the Gaussian beam, an aperture, and a lens -> u 2. FFT(u) 3. Multiplication of FFT(u) with the angular spectrum kernel H 4. IFFT(FU*H) Steps 2-4 are repeated 1000 times to show the propagation of the beam in z-direction Unfortunately the calculation results in a lot of horizontal lines that should not be there. The program produces reasonable results besides this problem. It is especially interesting that an equivalent calculation with Matlab or Octave has the same results without these artifacts. Both calculations take approximately 90 seconds. I am not sure whether I am doing something wrong or there is a difference in the FFT codes. I assume that the problem has something to do with the propagation kernel H but I cannot see a difference between both programs. Thank you very much for your help! Matthias --------------------------------------------------------------------------------------------------- Python code with comments (version without comments below) --------------------------------------------------------------------------------------------------- from pylab import * from numpy.lib import scimath #Complex roots import os ######################################### ########## Defintion of the input variables ######### ######################################### ##General variables lam = 532.e-9 #wavelength in m k = 2*pi/lam #wave number ##Computational array arr_size = 80.e-3 #size of the computation array in m N = 2**17 #points of the 1D array ##Aperture ap = 40.e-3 #aperture size of the 1D array in m ##Gaussian beam w0 =10.e-3 #waist radius of the gaussian beam in m ##Lenses (definition of the focal length) n = 1.56 #refractive index of the glass f = .1 #focal length in m Phi = 1/f #focal power of the lens r2 = 1.e18 #radius of the second lens surface r1 = 1 / ( Phi/(n-1) + 1/r2 ) #computation of the radius of the first lens surface ##z distances zmin = 0.99*f #initial z position zmax = 1.01*f #final z position zsteps = 1000 #number of computated positions ROI = 1000 #Region of interest in the diffraction image ############################################## ############### Program execution 1D ################ ############################################### x = linspace(-arr_size/2,arr_size/2,N) A = ones(N) A[where(ap/2<=abs(x))] = 0 #Definition of the aperture G = exp(- x**2/w0**2) #Generation of the Gaussian beam delta = -r1*(1 - scimath.sqrt(1 - x**2/r1**2)) + r2*(1 - scimath.sqrt(1 - x**2/r2**2)) m = (r1**2 <= x**2 ) delta[m] = 0 #correction in case of negative roots m = (r2**2 <= x**2 ) delta[m] = 0 #correction in case of negative roots lens = exp(1j * 2 * pi / lam * (n-1) * delta) #Calculation of the lens phase function u = A*G*lens #Complex amplitude before beam propagation ############################################ ########### Start angular spectrum method ########### deltaf = 1/arr_size #spacing in frequency domain fx = r_[-N/2*deltaf:(N/2)*deltaf:deltaf] #whole frequency domain FU = fftshift(fft(u)) #1D FFT U = zeros((zsteps,ROI)) #Generation of the image array z = linspace(zmin,zmax,zsteps) #Calculation of the different z positions for i in range(zsteps): H = exp(1j * 2 * pi / lam * z[i] * scimath.sqrt(1-(lam*fx)**2)) U[i] = (ifft(fftshift(FU*H)))[N/2 - ROI/2 : N/2 + ROI/2] if i%10 == 0: t = 'z position: %4i' % i print t ############ End angular spectrum method ############ ############################################# res = abs(U) imshow(res) show() ---------------------------------------------------------- Python code without comments ---------------------------------------------------------- from pylab import * from numpy.lib import scimath #Complex roots import os lam = 532.e-9 k = 2*pi/lam arr_size = 80.e-3 N = 2**17 ap = 40.e-3 w0 =10.e-3 n = 1.56 f = .1 Phi = 1/f r2 = 1.e18 r1 = 1 / ( Phi/(n-1) + 1/r2 ) zmin = 0.99*f zmax = 1.01*f zsteps = 1000 ROI = 1000 x = linspace(-arr_size/2,arr_size/2,N) A = ones(N) A[where(ap/2<=abs(x))] = 0 G = exp(- x**2/w0**2) delta = -r1*(1 - scimath.sqrt(1 - x**2/r1**2)) + r2*(1 - scimath.sqrt(1 - x**2/r2**2)) m = (r1**2 <= x**2 ) delta[m] = 0 m = (r2**2 <= x**2 ) delta[m] = 0 lens = exp(1j * 2 * pi / lam * (n-1) * delta) u = A*G*lens deltaf = 1/arr_size fx = r_[-N/2*deltaf:(N/2)*deltaf:deltaf] FU = fftshift(fft(u)) U = zeros((zsteps,ROI)) z = linspace(zmin,zmax,zsteps) for i in range(zsteps): H = exp(1j * 2 * pi / lam * z[i] * scimath.sqrt(1-(lam*fx)**2)) U[i] = (ifft(fftshift(FU*H)))[N/2 - ROI/2 : N/2 + ROI/2] if i%10 == 0: t = 'z position: %4i' % i print t res = abs(U) imshow(res) show() ---------------------------------------------------------- Matlab, Octave code ---------------------------------------------------------- lam = 532.e-9 k = 2*pi/lam arr_size = 80.e-3 N = 2^17 ap = 40.e-3 w0 =10.e-3 n = 1.56 f = .1 Phi = 1/f r2 = 1.e18 r1 = 1 / ( Phi/(n-1) + 1/r2 ) zmin = 0.99*f zmax = 1.01*f zsteps = 1000 ROI = 1000 x=linspace(-arr_size/2,arr_size/2,N); A = ones(1,N); A(find(ap/2<=abs(x)))=0; G = exp(- x.^2/w0^2); delta = -r1*(1 - sqrt(1 - x.^2/r1^2)) + r2*(1 - sqrt(1 - x.^2/r2^2)); delta(find(r1.^2 <= x.^2 ))=0; delta(find(r2.^2 <= x.^2 ))=0; lens = exp(1j * 2 * pi / lam * (n-1) * delta); u = A.*G.*lens; deltaf = 1/arr_size; fx = [-N/2*deltaf:deltaf:(N/2-1)*deltaf]; FU = fftshift(fft(u)); U = zeros(zsteps,ROI); z = linspace(zmin,zmax,zsteps); for i=1:zsteps H = exp(1j * 2 * pi / lam * z(i) * sqrt(1-lam.^2*fx.^2)); U(i,:) = ifft(fftshift(FU.*H))(N/2 - ROI/2 : N/2 + ROI/2-1); end imagesc(abs(U))

5 9

Appending data to a big ndarray
by oc-spam66 Aug. 9, 2008

Aug. 9, 2008

Hello, I would like to build a big ndarray by adding rows progressively. I considered the following functions : append, concatenate, vstack and the like. It appears to me that they all create a new array (which requires twice the memory). Is there a method for just adding a row to a ndarray without duplicating the data ? Regards, O.C. Créez votre adresse électronique prenom.nom(a)laposte.net 1 Go d'espace de stockage, anti-spam et anti-virus intégrés.

3 2

Re: [Numpy-discussion] numpy.fromstring() : error handling ?
by oc-spam66 Aug. 9, 2008

Aug. 9, 2008

> > Shouldn't it raise an exception ValueError ? (because "abcd" is not a float) > > I don't think so, but it shouldn't return a zero either. > > That call should mean: scan this whitespace separated string for as many > floating point numbers as it has. There are none, so it should return > and empty array: Is it really how it should be ? I thought this method should be restrictive in what it accepts (because the goal is to fill an array with precise dimensions). When floating point numbers are mixed with something else, it should be parsed with a regexp. Créez votre adresse électronique prenom.nom(a)laposte.net 1 Go d'espace de stockage, anti-spam et anti-virus intégrés.

3 2

np.bincount() won't accept list or array from np.random.beta and acts like len
by Michael Aug. 9, 2008

Aug. 9, 2008

Hi list, I have found np.bincount() not to behave as expected when provided with data from np.random.beta. It works fine with lists, but not the type returned by np.random.beta... there seems to be a block on casting too. When it does give output on np.array it just reports the number of items, like len() ?? x = np.array([1,1,1,1,2,2,4,4,5,6,6,6]) print np.bincount(x) prior = np.random.beta(2,7,1000) a = prior.tolist() print type(a) print np.bincount(a) b = [ round(i,2) for i in a ] print b; print type(b) print np.bincount(b) #c = np.array(a) #print np.bincount(c); print type(c) i know its easily solved but bincount is probably faster than making a pass over the probabilities and hashing them as strings etc have got a class that subclasses dict and accepts np.random.beta as 'seq': ... for prob in seq: self.bincount( prob ) def bincount( self, pr, increment=1 ): strPr = str(round(pr,self.dp)) self[strPr] = increment + self.get( strPr, 0 ) would prefer the one-liner: self.bincount = np.bincount(seq) any ideas? thanks. Michael Nandris (btw apologies for the cross posting)

2 1

Links for installing Numpy 1.1.1
by Matthieu Brucher Aug. 9, 2008

Aug. 9, 2008

Hi, There seems to be missing parameters when I try to download numpy 1.1.1 for Python 2.5 (David's superpack). And when I try to go to the numpy project on SF, the project is invalid as well. Is there something I'm missing ? Matthieu -- French PhD student Website : http://matthieu-brucher.developpez.com/ Blogs : http://matt.eifelle.com and http://blog.developpez.com/?blog=92 LinkedIn : http://www.linkedin.com/in/matthieubrucher

2 2

Re: [Numpy-discussion] numpy.fromstring() : error handling ?
by oc-spam66 Aug. 8, 2008

Aug. 8, 2008

Thank you for the answers, I am now disturbed by this result : > In [1]: import numpy > In [2]: numpy.fromstring("abcd", dtype = float, sep = ' ') > Out[2]: array([ 0.]) Shouldn't it raise an exception ValueError ? (because "abcd" is not a float) Regards, O.C. Créez votre adresse électronique prenom.nom(a)laposte.net 1 Go d'espace de stockage, anti-spam et anti-virus intégrés.

2 1

ATTENTION: 1.2.0b1 tagged tomorrow
by Jarrod Millman Aug. 8, 2008

Aug. 8, 2008

Just a reminder: I will be tagging 1.2.0b1 tomorrow evening. -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

1 0

numpy.fromstring() : error handling ?
by oc-spam66 Aug. 7, 2008

Aug. 7, 2008

Hello, the following code drives python into an endless loop : >>> import numpy >>> numpy.fromstring("abcd", dtype = float, sep = ' ') I think the function numpy.fromstring is lacking an adequate error handling for that case. Is it a bug ? Regards, -- O.C. Python 2.5.2 Debian Lenny Créez votre adresse électronique prenom.nom(a)laposte.net 1 Go d'espace de stockage, anti-spam et anti-virus intégrés.

3 2

Bug? Numpy crash with bad __array_interface__
by Zachary Pincus Aug. 7, 2008

Aug. 7, 2008

Hello all, As per my previous email, I encountered a strange sometimes-segfault when using 'numpy.array(thing)' to convert 'thing' (which provides an __array_interface__) to a numpy array. The offending __array_interface__ has a 'data' item that is a python string (not, as specified in the protocol, a pointer to a memory area) which is too small for the provided array size. Which, in some cases, causes a segfault. Is this a bug? That is, should numpy check the size of the 'data' (when 'data' is in fact a sized python object and not just a bare pointer) against the size of the array and provided typestr? Or is this too much of a corner case to deal with... Zach

1 0

Error using rec.fromrecords with specific dtype
by Sameer DCosta Aug. 7, 2008

Aug. 7, 2008

Hi, I am unable to specify a dtype when using numpy.rec.fromrecords. This problem only occurs when one of the dtypes is an object. (Please see below) Thanks in advance for your help. Sameer In [7]: import datetime, numpy as np In [8]: print np.__version__ 1.2.0.dev5253 In [9]: dtype=np.dtype([('spam', '|O4'),('ham', '<i4')] ) In [10]: data = [ (datetime.date(2001,1,1), 2), (datetime.date(2001,1,2), 4), (datetime.date(2001,1,3), 6) ] In [11]: n = np.rec.fromrecords(data, dtype=dtype) --------------------------------------------------------------------------- TypeError Traceback (most recent call last) /home/titan/sameer/<ipython console> /home/titan/sameer/local/mins/lib/python2.4/site-packages/numpy/core/records.pyc in fromrecords(recList, dtype, shape, formats, names, titles, aligned, byteorder) 469 res = retval.view(recarray) 470 --> 471 res.dtype = sb.dtype((record, res.dtype)) 472 return res 473 /home/titan/sameer/local/mins/lib/python2.4/site-packages/numpy/core/records.pyc in __setattr__(self, attr, val) 287 if attr not in fielddict: 288 exctype, value = sys.exc_info()[:2] --> 289 raise exctype, value 290 else: 291 fielddict = ndarray.__getattribute__(self,'dtype').fields or {} TypeError: Cannot change data-type for object array. In [12]: names = "spam,ham" In [13]: n = np.rec.fromrecords(data, names=names) In [14]: print n [(datetime.date(2001, 1, 1), 2) (datetime.date(2001, 1, 2), 4) (datetime.date(2001, 1, 3), 6)]

1 0

numpy / ipython heisenbug crasher help?
by Zachary Pincus Aug. 7, 2008

Aug. 7, 2008

Hello all, I have a bizarre bug that causes numpy to segfault, but: - only when run under ipython - only when numpy is imported *after* an other library (that does not import numpy) Here's what the code looks like. Crashes (only in ipython): import celltool.utility.pil_lite.Image as Image, numpy i = Image.open('untreated_2_logical.png') a = numpy.array(i) Does not crash: import numpy, celltool.utility.pil_lite.Image as Image i = Image.open('untreated_2_logical.png') a = numpy.array(i) The 'i' object exposes an __array_interface__, and is produced from my personal fork of the PIL. The crash is thusly reported: Exception Type: EXC_BAD_ACCESS (SIGSEGV) Exception Codes: KERN_INVALID_ADDRESS at 0x000000000251f000 0 libSystem.B.dylib 0xffff0a13 __memcpy + 627 1 multiarray.so 0x022f7af8 _array_copy_into + 2840 (arrayobject.c:1113) 2 multiarray.so 0x022f88af PyArray_FromArray + 495 (arrayobject.c:8334) 3 multiarray.so 0x022ee314 PyArray_FromAny + 308 (arrayobject.c: 8814) 4 multiarray.so 0x022fb01c PyArray_CheckFromAny + 92 (arrayobject.c:8995) 5 multiarray.so 0x022fb2e6 _array_fromobject + 310 (multiarraymodule.c:5787) [And then the python stack frames] Can anyone suggest what a good next step in reducing this patently insane bug down to a more minimal test case is? I had this bug with numpy 1.1 as well as the current SVN head. Zach

1 1

FFTW and config vs system_info
by jordan＠math.ucsb.edu Aug. 7, 2008

Aug. 7, 2008

Hi all, I'm trying to get scipy to recognize my FFTW install. I know there is another thread on this, but I didn't understand the discussion fully. Which file do I modify to let scipy know where my FFTW install is? As a possibly related note, running scipy.config() and scipy.distutil.system_info.py give conflicting answers as to whether I have BLAS and LAPACK etc. It'd be nice to know for sure what scipy is using. Any light on this matter would be greatly appreciated. Cheers, Jordan

2 2

[ANNOUNCE] EPD 2.5.2001 for OS X released!
by Dave Peterson Aug. 7, 2008

Aug. 7, 2008

I'm pleased to announce that Enthought has released the Enthought Python Distribution (EPD) 2.5.2001 for OS X! EPD is a Distribution of the Python Programming Language (currently version 2.5.2) that includes over 60 additional libraries, including ETS 2.7.1. Please visit the EPD website (http://www.enthought.com/epd) to get the OS X release, or just to find out more information about EPD 2.5.2001 So what’s the big deal? In addition to making everyones’ life easier with installation, EPD also represents a common suite of functionality deployed across platforms such as Windows XP, RedHat Linux, and now OS X 10.4 and above. The cross-platform promise of Python is better realized because it’s trivial for everyone to get substantially the same set of libraries installed on their system with a single-click install. What’s the catch? You knew it was coming, huh? If you’d like to use EPD in a Commercial or Governmental entity, we do ask you to pay for an annual subscription to download and update EPD. For academics and non-profit, private-sector organizations, EPD is and will remain free. Let’s be clear, though. EPD is the bundle of software. People pay for the subscription to download the bundle. The included libraries are, of course, freely available separately under the terms of license for each individual package (this should sound familiar). The terms for the bundle subscription are available at http://www.enthought.com/products/epdlicense.php. BTW, anyone can try it out for free for 30 days. If you have questions, check out the FAQ at http://www.enthought.com/products/epdfaq.php or drop us a line at epd-users(a)enthought.com. And just one more note: Enthought is deeply grateful to all those who have contributed to these libraries over the years. We’ve built a business around the things they allow us to do and we appreciate such a nice set of tools and the privilege of being part of the community that created them.

1 0

FAIL: Test corrcoef 1 1D variable w/missing values
by Nils Wagner Aug. 7, 2008

Aug. 7, 2008

Hi all, Can someone reproduce the failure with >>> numpy.__version__ '1.2.0.dev5618' ====================================================================== FAIL: Test corrcoef 1 1D variable w/missing values ---------------------------------------------------------------------- Traceback (most recent call last): File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/ma/tests/test_extras.py", line 463, in test_1d_w_missing corrcoef(x, x[::-1], rowvar=False)) File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/ma/testutils.py", line 120, in assert_equal return assert_array_equal(actual, desired, err_msg) File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/ma/testutils.py", line 186, in assert_array_equal header='Arrays are not equal') File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/ma/testutils.py", line 180, in assert_array_compare verbose=verbose, header=header) File "/data/home/nwagner/local/lib/python2.5/site-packages/numpy/testing/utils.py", line 289, in assert_array_compare assert cond, msg AssertionError: Arrays are not equal (mismatch 50.0%) x: array([[ 1. , -0.38605094], [-0.38605094, 1. ]]) y: array([[ 1. , -0.38605094], [-0.38605094, 1. ]]) Nils

2 1

confusion on importing numpy
by mark Aug. 6, 2008

Aug. 6, 2008

Hello list. I am confused about importing numpy. When I do from numpy import * and I try the min function, I get the default Python min function. On the other hand, when I do import numpy as np and use the np.min function, I get the numpy min function (which is obviously what I want). I know, the latter is preferred, but the former is so dang easy to type. I always thought both imports work the same, but I am obviously wrong. What's the difference? Thanks, Mark

6 7

FFT usage / consistency
by Felix Richter Aug. 6, 2008

Aug. 6, 2008

Hi all, I found myself busy today trying to understand what went wrong in my FFT code. I wrote a minimal example/testing code to check the FFT output against an analytic result and also tried to reverse the transformation to get the original function back. Most curiously, the results depend on whether I first do the fft and then the ifft or the other way round. For the test function, I use the Lorentz function 1/(x^2+1). The exact FT is exp(-|t|)*sqrt(pi/2), the IFT yields the same. 1) First FFT and then IFFT: The real part of FFT oscillates, the imaginary part is not zero, and the magnitudes do not match. All this should not be, but the IFFT reproduces the original function just fine. 2) First IFFT and then FFT: The IFFT is perfect, but the FFT does not reproduce the original function. Could someone please have a look and tell me what I am doing wrong? The code is attached, it also plots the results nicely. Maybe the (corrected) code could be used as an example in the documentation (correct use of (i)fftshift and (i)fftfreq is not trivial!) or as an additional test case. The existing test cases for numpy only seem to check that the fft function can be called, but they do not check consistency of results. I'm using NumPy versions 1.0.4 and 1.1.0 on Linux with fftpack_lite.so (even though fftw3 is installed and configured, but I'll probably ask for that later...) Thanks a lot, Felix

5 5

Re: [Numpy-discussion] Calculating roots with negative numbers
by Matthias Hillenbrand Aug. 6, 2008

Aug. 6, 2008

Hello, > When you convolve two signals, of lengths N and M, you need to pad the > FFTs to length (N+M-1) before multiplication. > You can take a look at my linear position-invariant filtering code at: > http://mentat.za.net/hg/filter I understand your comments that I need zero padding when doing FFT filtering. I will need some time to understand Stéfan's program as I am a NumPy/Python newbie and not used to object oriented programming. What I don't understand is that Octave/Matlab and Numpy generate different results with nearly the same code. Why are the results with Octave/Matlab OK although I have not applied additional zero padding. In my opinion I should get exactly the same results. Regards, Matthias --------------------------------------------------------------------------------------------------- Python code with comments (version without comments below) --------------------------------------------------------------------------------------------------- from pylab import * from numpy.lib import scimath #Complex roots import os ######################################### ########## Defintion of the input variables ######### ######################################### ##General variables lam = 532.e-9 #wavelength in m k = 2*pi/lam #wave number ##Computational array arr_size = 80.e-3 #size of the computation array in m N = 2**17 #points of the 1D array ##Aperture ap = 40.e-3 #aperture size of the 1D array in m ##Gaussian beam w0 =10.e-3 #waist radius of the gaussian beam in m ##Lenses (definition of the focal length) n = 1.56 #refractive index of the glass f = .1 #focal length in m Phi = 1/f #focal power of the lens r2 = 1.e18 #radius of the second lens surface r1 = 1 / ( Phi/(n-1) + 1/r2 ) #computation of the radius of the first lens surface ##z distances zmin = 0.99*f #initial z position zmax = 1.01*f #final z position zsteps = 1000 #number of computated positions ROI = 1000 #Region of interest in the diffraction image ############################################## ############### Program execution 1D ################ ############################################### x = linspace(-arr_size/2,arr_size/2,N) A = ones(N) A[where(ap/2<=abs(x))] = 0 #Definition of the aperture G = exp(- x**2/w0**2) #Generation of the Gaussian beam delta = -r1*(1 - scimath.sqrt(1 - x**2/r1**2)) + r2*(1 - scimath.sqrt(1 - x**2/r2**2)) m = (r1**2 <= x**2 ) delta[m] = 0 #correction in case of negative roots m = (r2**2 <= x**2 ) delta[m] = 0 #correction in case of negative roots lens = exp(1j * 2 * pi / lam * (n-1) * delta) #Calculation of the lens phase function u = A*G*lens #Complex amplitude before beam propagation ############################################ ########### Start angular spectrum method ########### deltaf = 1/arr_size #spacing in frequency domain fx = r_[-N/2*deltaf:(N/2)*deltaf:deltaf] #whole frequency domain FU = fftshift(fft(u)) #1D FFT U = zeros((zsteps,ROI)) #Generation of the image array z = linspace(zmin,zmax,zsteps) #Calculation of the different z positions for i in range(zsteps): H = exp(1j * 2 * pi / lam * z[i] * scimath.sqrt(1-(lam*fx)**2)) U[i] = (ifft(fftshift(FU*H)))[N/2 - ROI/2 : N/2 + ROI/2] if i%10 == 0: t = 'z position: %4i' % i print t ############ End angular spectrum method ############ ############################################# res = abs(U) imshow(res) show() ---------------------------------------------------------- Python code without comments ---------------------------------------------------------- from pylab import * from numpy.lib import scimath #Complex roots import os lam = 532.e-9 k = 2*pi/lam arr_size = 80.e-3 N = 2**17 ap = 40.e-3 w0 =10.e-3 n = 1.56 f = .1 Phi = 1/f r2 = 1.e18 r1 = 1 / ( Phi/(n-1) + 1/r2 ) zmin = 0.99*f zmax = 1.01*f zsteps = 1000 ROI = 1000 x = linspace(-arr_size/2,arr_size/2,N) A = ones(N) A[where(ap/2<=abs(x))] = 0 G = exp(- x**2/w0**2) delta = -r1*(1 - scimath.sqrt(1 - x**2/r1**2)) + r2*(1 - scimath.sqrt(1 - x**2/r2**2)) m = (r1**2 <= x**2 ) delta[m] = 0 m = (r2**2 <= x**2 ) delta[m] = 0 lens = exp(1j * 2 * pi / lam * (n-1) * delta) u = A*G*lens deltaf = 1/arr_size fx = r_[-N/2*deltaf:(N/2)*deltaf:deltaf] FU = fftshift(fft(u)) U = zeros((zsteps,ROI)) z = linspace(zmin,zmax,zsteps) for i in range(zsteps): H = exp(1j * 2 * pi / lam * z[i] * scimath.sqrt(1-(lam*fx)**2)) U[i] = (ifft(fftshift(FU*H)))[N/2 - ROI/2 : N/2 + ROI/2] if i%10 == 0: t = 'z position: %4i' % i print t res = abs(U) imshow(res) show() ---------------------------------------------------------- Matlab, Octave code ---------------------------------------------------------- lam = 532.e-9 k = 2*pi/lam arr_size = 80.e-3 N = 2^17 ap = 40.e-3 w0 =10.e-3 n = 1.56 f = .1 Phi = 1/f r2 = 1.e18 r1 = 1 / ( Phi/(n-1) + 1/r2 ) zmin = 0.99*f zmax = 1.01*f zsteps = 1000 ROI = 1000 x=linspace(-arr_size/2,arr_size/2,N); A = ones(1,N); A(find(ap/2<=abs(x)))=0; G = exp(- x.^2/w0^2); delta = -r1*(1 - sqrt(1 - x.^2/r1^2)) + r2*(1 - sqrt(1 - x.^2/r2^2)); delta(find(r1.^2 <= x.^2 ))=0; delta(find(r2.^2 <= x.^2 ))=0; lens = exp(1j * 2 * pi / lam * (n-1) * delta); u = A.*G.*lens; deltaf = 1/arr_size; fx = [-N/2*deltaf:deltaf:(N/2-1)*deltaf]; FU = fftshift(fft(u)); U = zeros(zsteps,ROI); z = linspace(zmin,zmax,zsteps); for i=1:zsteps H = exp(1j * 2 * pi / lam * z(i) * sqrt(1-lam.^2*fx.^2)); U(i,:) = ifft(fftshift(FU.*H))(N/2 - ROI/2 : N/2 + ROI/2-1); end imagesc(abs(U))

1 1

Functions that stores results without using the out keyword
by Bruce Southey Aug. 6, 2008

Aug. 6, 2008

Hi, The following functions support a optional second argument that stores the result: isfinite, isnan, isinf, isposinf and isneginf. There may be other functions with the same behavior. Is this usage equivalent to 'out' keyword that appears in other functions like mean? If so, should these functions have the 'out' keyword? Also, I am correct that the 'out' keyword must always be an array and not a scalar? Numpy does returns the error 'TypeError: output must be an array' if a scalar is used. I would like to be explicit in the documentation by stating that a scalar in any form can not be used. Thanks Bruce import numpy as np a = np.array([[1,2, np.nan],[3,4, np.inf]]) b=np.array([[0.0],[0.0]]) c = np.array([[2,2,2],[2,2,2]]) np.mean(a, axis=1, out=b) #returns: array([[ NaN], # [ Inf]]) np.isinf(a, c) #returns: array([[0, 0, 0], # [0, 0, 1]]) np.isnan(a, c) #returns: array([[0, 0, 1], # [0, 0, 0]]) d=a.mean() np.mean(a, out=d) #Traceback (most recent call last): # File "<stdin>", line 1, in <module> # File "/usr/lib64/python2.5/site-packages/numpy/core/fromnumeric.py", line 1806, in mean # return mean(axis, dtype, out) #TypeError: output must be an array

1 0

cubic spline function in numarray or numpy
by Gong, Shawn (Contractor) Aug. 6, 2008

Aug. 6, 2008

hi list, I am trying to find 1-D cubic spline function. But I am not able to call it. Error message can't find "ndimage" Would someone point me to 1-D cubic spline functions either in numarray or numpy? and their help pages? I use python 2.3.5 and numpy 1.0.1. Are they too old to use cubic spline? thanks, Shawn

6 6

Cython/NumPy syntax
by Dag Sverre Seljebotn Aug. 6, 2008

Aug. 6, 2008

I'd like input from anyone interested in the following syntax questions. The current experimental Cython syntax for efficient ndarray indexing looks like this: cdef numpy.ndarray[numpy.int64, 2] arr Issue 1: One can argue that "2" above looks like a length specifier to newcomers. Issue 2: This clashes with some forms of C array notation (you can do "sizeof(int[3])" or "extern ... def myfunc(int[])" in Cython, and though such syntax is redundant it will stay in order not to break code). Solving issue #1 helps with #2 as well because without a single integer literal in it it should look a lot less like an array declaration. So, some quick proposals: - Require a "D" suffix for the number of dimensions. This should make it very clear that the 2 stands for dimensionality and not shape. cdef numpy.ndarray[numpy.int64, 2D] arr Variations: a) numpy.ndarray[numpy.int64, 2D] b) numpy.ndarray[numpy.int64 2D] c) numpy.ndarray[2D numpy.int64] (Note that I want to leave the way open for adding mode="c" or mode="fortran", so using a "," seems good for that reason). - Require an ndim keyword: cdef numpy.ndarray[numpy.int64, ndim=2] - Other type of brackets. This boils down to (for various reasons) the <> brackets: cdef numpy.ndarray<numpy.int64, 2> (However, one should remember that Cython may want to support using C++ templates at some point). -- Dag Sverre

7 10

f2py: undefined symbol: zgemm_
by Harry Mangalam Aug. 6, 2008

Aug. 6, 2008

Hi All, Using: Python 2.5.2, f2py ver 2_4422, gfortran --version GNU Fortran (GCC) 4.2.3 (Ubuntu 4.2.3-2ubuntu7) on kubuntu 8.04 (Hardy) up-to-date on a Thinkpad T60. After building and using an f2py-generated lib for a while with this command: f2py --opt="-O3" -c -m fd_rrt1d --fcompiler=gnu95 \ --link-lapack_opt *.f I had to regenerate it after some code changes (just commenting out some debugging info). Altho the shared lib gets built correctly, when the calling python is run, I now get this: ./1d.py Traceback (most recent call last): File "./1d.py", line 27, in <module> from fd_rrt1d import * ImportError: /home/hjm/shaka/1D-Mangalam-py/fd_rrt1d.so: undefined symbol: zgemm_ sure enough, nm reports it as undefined: nm fd_rrt1d.so |tail 000065b0 t string_from_pyobj U strlen@@GLIBC_2.0 U strncpy@@GLIBC_2.0 0001df40 b u.1294 00013bee T umatrix1d_dms0_ 000128c3 T umatrix1d_dms1_ 00015458 T umatrix1d_ms0_ U zgemm_ U zgemv_ U zgesv_ but why the recent change? Google and the numpy list seem not to have heard about this one.. -- Harry Mangalam - Research Computing, NACS, E2148, Engineering Gateway, UC Irvine 92697 949 824-0084(o), 949 285-4487(c) -- ..Kick at the darkness til it bleeds daylight. (Lovers in a Dangerous Time) - Bruce Cockburn

1 1

dot/inner with kwargs?
by Hans Meine Aug. 6, 2008

Aug. 6, 2008

Hi, I am trying to do something like the following as efficiently as possible: result[...,0] = (dat * cos(arange(100))).sum(axis = -1) result[...,1] = (dat * sin(arange(100))).sum(axis = -1) Thus, I was looking for dot / inner with an 'output' arg. Q1: What is the difference between dot and inner (apart from the fact that dot does not seem to have docs)? In http://www.scipy.org/Numpy_Example_List_With_Doc#inner, I read that: > Like the generic NumPy equivalent the product sum is over > the last dimension of a and b. (BTW: "innerproduct" should be "inner") Q2: What is "the generic NumPy equivalent" here? -- Ciao, / / /--/ / / ANS

1 0

cubic spline function in numarray or numpy
by Gong, Shawn (Contractor) Aug. 6, 2008

Aug. 6, 2008

hi list, I am trying to find 1-D cubic spline function. Google search yields that there is a function called spline_filter1d in numarray. But I am not able to call it. Would someone point me to 1-D cubic spline function either in numarray or numpy? I can find source codes and enter in Python. But I prefer the built-in function. thanks, Shawn

1 0

Segfault in PyArray_Item_XDECREF when using recarray object references titles
by Michael Droettboom Aug. 6, 2008

Aug. 6, 2008

I've run into a segfault that occurs in the array destructor with arrays containing object references with both names and titles. When a field contains both and name and a title, the fields dictionary contains two entries for that field. This means that the array item destructor (which iterates through the fields dictionary) will decref the pointed-to object twice. If the first decref causes the object to be deleted, the second decref has the potential to segfault. It seems the simplest patch is to set the object pointer to NULL after decref'ing, so the second decref will do nothing. However, perhaps there is a way to avoid decref'ing twice in the first place. I've attached a script that exercises the segfault, a gdb backtrace, and a patch. You may need to adjust the number of rows until it is high enough to create a segfault on your system. This is on: RHEL4 Python 2.5.2 Numpy SVN r5497 Cheers, Mike > gdb python GNU gdb Red Hat Linux (6.3.0.0-1.153.el4_6.2rh) Copyright 2004 Free Software Foundation, Inc. GDB is free software, covered by the GNU General Public License, and you are welcome to change it and/or distribute copies of it under certain conditions. Type "show copying" to see the conditions. There is absolutely no warranty for GDB. Type "show warranty" for details. This GDB was configured as "i386-redhat-linux-gnu"...Using host libthread_db library "/lib/tls/libthread_db.so.1". (gdb) run segfault.py Starting program: /wonkabar/data1/usr/bin/python segfault.py [Thread debugging using libthread_db enabled] [New Thread -1208489312 (LWP 30028)] len(dtype) = 1, len(dtype.fields) = 2 {'name': (dtype('object'), 0, 'title'), 'title': (dtype('object'), 0, 'title')} Program received signal SIGSEGV, Segmentation fault. [Switching to Thread -1208489312 (LWP 30028)] 0x0097285e in PyArray_Item_XDECREF ( data=0xb7a3e780 "\uffff_\224\uffff `\214\uffff(`\214\uffff0`\214\uffff8`\214\uffff@`\214\uffffH`\214\uffffP`\214\uffffX`\214\uffff``\214\uffffh`\214\uffffp`\214\uffffx`\214\uffff\200`\214\uffff\210`\214\uffff\220`\214\uffff\230`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff", descr=0x9d4680) at numpy/core/src/arrayobject.c:198 198 Py_XDECREF(*temp); (gdb) bt #0 0x0097285e in PyArray_Item_XDECREF ( data=0xb7a3e780 "\uffff_\224\uffff `\214\uffff(`\214\uffff0`\214\uffff8`\214\uffff@`\214\uffffH`\214\uffffP`\214\uffffX`\214\uffff``\214\uffffh`\214\uffffp`\214\uffffx`\214\uffff\200`\214\uffff\210`\214\uffff\220`\214\uffff\230`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff\uffff`\214\uffff", descr=0x9d4680) at numpy/core/src/arrayobject.c:198 #1 0x00991bc7 in PyArray_XDECREF (mp=0xb7ae4f0c) at numpy/core/src/arrayobject.c:211 #2 0x009a579b in array_dealloc (self=0xb7ae4f0c) at numpy/core/src/arrayobject.c:2089 #3 0x0809781f in subtype_dealloc (self=0xb7ae4f0c) at Objects/typeobject.c:709 #4 0x08082a02 in PyDict_SetItem (op=0xb7f56acc, key=0xb7ea7d80, value=0x81379c0) at Objects/dictobject.c:416 #5 0x08085a1e in _PyModule_Clear (m=0xb7f3e0ec) at Objects/moduleobject.c:136 #6 0x080d7138 in PyImport_Cleanup () at Python/import.c:439 #7 0x080e4343 in Py_Finalize () at Python/pythonrun.c:399 #8 0x08056633 in Py_Main (argc=1, argv=0xbff1ca24) at Modules/main.c:545 #9 0x08056323 in main (argc=2, argv=0xbff1ca24) at ./Modules/python.c:23 -- Michael Droettboom Science Software Branch Operations and Engineering Division Space Telescope Science Institute Operated by AURA for NASA Index: numpy/core/src/arrayobject.c =================================================================== --- numpy/core/src/arrayobject.c (revision 5497) +++ numpy/core/src/arrayobject.c (working copy) @@ -196,6 +196,7 @@ if (descr->type_num == PyArray_OBJECT) { temp = (PyObject **)data; Py_XDECREF(*temp); + *temp = NULL; } else if PyDescr_HASFIELDS(descr) { PyObject *key, *value, *title=NULL; import numpy as np ROWS = 200000 ra = np.recarray(ROWS, dtype=[(('title', 'name'), 'O')]) print "len(dtype) = %d, len(dtype.fields) = %d" % (len(ra.dtype), len(ra.dtype.fields)) print ra.dtype.fields for i in range(ROWS): ra[i]["name"] = object()

1 2

Slicing without a priori knowledge of dimension
by Matthew Czesarski Aug. 6, 2008

Aug. 6, 2008

Dear list. I've got a feeling that what I'm trying to do *should* be easy but at the moment I can't find a non-brute-force method. I'm working with quite a high-rank matrix; 7 dimensions filled with chi**2 values. It's form is something like this: chi2 = numpy.ones((3,4,5,6,7,8,9)) What I need to do is slice out certain 2 dimensional grids that I then want to plot for confidence estimation; make a nice graph. This is fine: as easy as it gets if the 2 dimensions are known. E.g. for the 2nd and 5th axes, one could hardcode something like this: subChi2 = chi2[ ia, ib, :, id, ie, :, ig ] The difficulty is that the user is to state which plane (s)he wants to slice out and we can't code the above. The function itself has to convert 2 rank numbers into an expression for a slice and I can't currently figure out how to do that. There are a huge number of options (well, there are 7*6=42). If one could just manually write a colon into a tuple like (2,2,:,2,2,:,2) or something even like 2,2,0:len(2ndDim),2,2,0:len(5thDim),2, things would be fine. But that doesn't seem to be an option. An alternative would be to set up 2 loops and incrementally read out the individual elements of the slice to a new numpy.zeros ( ( len(2ndDim), len(5thDim) ) ), but again we seem to be diverging from elegance. There must be something that I'm missing. Could somebody have a pointer as to what it is? Thanks in advance, Matt

2 2

member1d and unique elements
by Greg Novak Aug. 6, 2008

Aug. 6, 2008

I have two arrays of integers, and would like to know _where_ they have elements in common, not just _which_ elements are in common. This is because the entries in the integer array are aligned with other arrays. This seems very close to what member1d advertises as its function. However, member1d says that it expects arrays with only unique elements. First of all, my desired operation is well-posed: I'd like f(ar1, ar2) to return something in the shape of ar1 with True if the value at that position appears anywhere in ar2 (regardless of duplication) and False otherwise. So I looked at the code and have two questions: 1) What is this code trying to achieve? aux = perm[ii+1] perm[ii+1] = perm[ii] perm[ii] = aux Here perm is the stable argsort of the two concatenated arguments: perm = concatenate((ar1, ar2)).argsort(kind='mergesort'). arr is the array of combined inputs in sorted order: arr = concatenate((ar1, ar2))[perm] and ii is a list of indices into arr where the value of arr is equal to the next value in the array (arr[ii] == arr[ii+1]) _and_ arr[ii] came from the _second_ input (ar2). Now, this last bit (looking for elements of arr that are equal and both came from the second array) is clearly trying to deal with duplication, which is why I'm interested... So, the code snippet is trying to swap perm[ii+1] with perm[ii], but I don't see why. Furthermore, there are funny results if a value is duplicated three times, not just twice -- perm is no longer a permutation vector. Eg, member1d([1], [2,2,2]) results perm=[0,1,2,3] and ii=[1,2] before the above snippet, and the above snippet makes perm into [0,2,3,2] I've commented those three lines, and I've never seen any changes to the output of member1d. The new value of perm is used to compute the expression: perm.argsort(kind='mergesort')[:len( ar1 )], but the changes to that expression as a result of the above three lines are always at the high end of the array, which is sliced off by the last [:len(ar1)]. Finally, my second question is: 2) Does anyone have a test case where member1d fails as a result of duplicates in the input? So far I haven't found any, with the above lines commented or not. Upon reflection and review of the changelog, another theory occurs to me: member1d did not originally use a stable sort. What I've written above for interpretation of the value ii (indicates duplication within ar2) is true for a stable sort, but for an unstable sort the same condition has the interpretation that ii holds the values where the sorting algorithm swapped the order of equal values unstably. Then the code snippet in question 1) looks like an attempt to swap those values in the permutation array to make the sort stable again. The attempt would fail if there was duplication in either array. So, I would propose deleting those three lines (since they seem to be a non-functional relic) and declaring in the docstring that member1d doesn't require unique elements. Also, if this is correct, then the function simplifies considerably since several values don't need to be computed anymore: def setmember1d( ar1, ar2 ): ar = nm.concatenate( (ar1, ar2 ) ) perm = ar.argsort(kind='mergesort') aux = ar[perm] flag = nm.concatenate( (aux[1:] == aux[:-1], [False] ) ) indx = perm.argsort(kind='mergesort')[:len( ar1 )] return flag[indx] Corrections to the above are welcome since I'm going to start using member1d without regard for uniqueness, and I'd like to know if I'm making a big mistake... Thanks, Greg

4 5

Numpy random: independent streams
by Ludwig Aug. 6, 2008

Aug. 6, 2008

The Python standard random API allows me to define multiple independent random streams, which I can control with their own seed, e.g. import random generator_1 = random.Random() generator_2 = random.Random() generator_1.seed(100) generator_2.seed(100) So now generator_1 and 2 will produce the same sequence. Is there an equivalent for this in numpy? (Apart from, as it seems, the oldnumeric package http://numpy.scipy.org/numpydoc/numpy-20.html)? Ludwig

2 1

Bilteral filter
by Nadav Horesh Aug. 5, 2008

Aug. 5, 2008

Hi all, Is there any efficient implementation of bilateral filter that works smoothly with numpy? Nadav.

4 5

Changeset 5609 of numpy breaks scipy build
by Fernando Perez Aug. 5, 2008

Aug. 5, 2008

Howdy, I just noticed that scipy wasn't building via my little script that does a numpy/mpl/scipy svn update+clean build. Chris Burn's excellent hunch pinpointed the culprit: http://projects.scipy.org/scipy/numpy/changeset/5609 Innocent looking as it may be, with numpy after that revision, building current scipy SVN gives: creating build/temp.linux-x86_64-2.5/scipy/sparse/sparsetools compile options: '-I/home/fperez/usr/opt/lib/python2.5/site-packages/numpy/core/include -I/usr/include/python2.5 -c' g++: scipy/sparse/sparsetools/csr_wrap.cxx scipy/sparse/sparsetools/csr_wrap.cxx: In function 'int require_size(PyArrayObject*, npy_intp*, int)': scipy/sparse/sparsetools/csr_wrap.cxx:2902: error: expected `)' before 'PRIdPTR' scipy/sparse/sparsetools/csr_wrap.cxx:2902: warning: spurious trailing '%' in format scipy/sparse/sparsetools/csr_wrap.cxx:2902: warning: too many arguments for format scipy/sparse/sparsetools/csr_wrap.cxx:2909: error: expected `)' before 'PRIdPTR' scipy/sparse/sparsetools/csr_wrap.cxx:2909: warning: spurious trailing '%' in format scipy/sparse/sparsetools/csr_wrap.cxx:2909: warning: too many arguments for format scipy/sparse/sparsetools/csr_wrap.cxx: In function 'int require_size(PyArrayObject*, npy_intp*, int)': scipy/sparse/sparsetools/csr_wrap.cxx:2902: error: expected `)' before 'PRIdPTR' scipy/sparse/sparsetools/csr_wrap.cxx:2902: warning: spurious trailing '%' in format scipy/sparse/sparsetools/csr_wrap.cxx:2902: warning: too many arguments for format scipy/sparse/sparsetools/csr_wrap.cxx:2909: error: expected `)' before 'PRIdPTR' scipy/sparse/sparsetools/csr_wrap.cxx:2909: warning: spurious trailing '%' in format scipy/sparse/sparsetools/csr_wrap.cxx:2909: warning: too many arguments for format error: Command "g++ -pthread -fno-strict-aliasing -DNDEBUG -O2 -g -pipe -Wall -Wp,-D_FORTIFY_SOURCE=2 -fexceptions -fstack-protector --param=ssp-buffer-size=4 -m64 -mtune=generic -D_GNU_SOURCE -fPIC -fPIC -I/home/fperez/usr/opt/lib/python2.5/site-packages/numpy/core/include -I/usr/include/python2.5 -c scipy/sparse/sparsetools/csr_wrap.cxx -o build/temp.linux-x86_64-2.5/scipy/sparse/sparsetools/csr_wrap.o" failed with exit status 1 ### I don't know what the right solution needs to be, but I figured we'd at least report it quickly so it can be looked into. Cheers, f

2 2

Error creating an array
by Sameer DCosta Aug. 5, 2008

Aug. 5, 2008

Im having a little trouble creating a numpy array with a specific dtype. I can create the array b with dtype=int, but not with the dtype I want. Any idea what I am doing wrong here? In [450]: import numpy as np In [451]: print np.__version__ 1.2.0.dev5243 In [452]: dtype=np.dtype([('spam', '<i4',1),('ham', '<i4',1)] ) In [453]: data = [ [1,2], [3,4], [5,6] ] In [454]: b = np.array( data, dtype=dtype) --------------------------------------------------------------------------- TypeError Traceback (most recent call last) /home/titan/sameer/tmp/test.py ----> 1 2 3 4 5 TypeError: expected a readable buffer object In [455]: b = np.array( data, dtype=int) In [456]: print b [[1 2] [3 4] [5 6]] Thanks. Sameer

2 3

64bit issue?
by Charles Doutriaux Aug. 5, 2008

Aug. 5, 2008

Hello I'm running into some strange error on a 64bit machine, I tracked it down to this line returning a NULL pointer, any idea why is this? I tried both numpy1.1.1 and what in trunk, numpy.test() passes for both Ok here's the uname of the machine and the offending line: Linux quartic.txcorp.com 2.6.20-1.2320.fc5 #1 SMP Tue Jun 12 18:50:49 EDT 2007 x86_64 x86_64 x86_64 GNU/Linux array = (PyArrayObject *)PyArray_SimpleNew(d, dims, self->type); where d is 3 and dims: 120,46,72 self->type is 11 it seems to pass with d=1, but i'm not 100% positive. Any idea on why it could be failing? Thanks, C. PS: just in case here are numpy.test results: >>> import numpy >>> numpy.test() Numpy is installed in /home/facets/doutriaux1/cdat/latest/lib/python2.5/site-packages/numpy Numpy version 1.1.1 Python version 2.5.2 (r252:60911, Aug 4 2008, 13:47:12) [GCC 4.1.1 20070105 (Red Hat 4.1.1-51)] Found 3/3 tests for numpy.core.tests.test_memmap Found 145/145 tests for numpy.core.tests.test_regression Found 12/12 tests for numpy.core.tests.test_records Found 3/3 tests for numpy.core.tests.test_errstate Found 72/72 tests for numpy.core.tests.test_numeric Found 36/36 tests for numpy.core.tests.test_numerictypes Found 290/290 tests for numpy.core.tests.test_multiarray Found 18/18 tests for numpy.core.tests.test_defmatrix Found 63/63 tests for numpy.core.tests.test_unicode Found 16/16 tests for numpy.core.tests.test_umath Found 7/7 tests for numpy.core.tests.test_scalarmath Found 2/2 tests for numpy.core.tests.test_ufunc Found 5/5 tests for numpy.distutils.tests.test_misc_util Found 4/4 tests for numpy.distutils.tests.test_fcompiler_gnu Found 2/2 tests for numpy.fft.tests.test_fftpack Found 3/3 tests for numpy.fft.tests.test_helper Found 15/15 tests for numpy.lib.tests.test_twodim_base Found 1/1 tests for numpy.lib.tests.test_machar Found 1/1 tests for numpy.lib.tests.test_regression Found 43/43 tests for numpy.lib.tests.test_type_check Found 1/1 tests for numpy.lib.tests.test_ufunclike Found 15/15 tests for numpy.lib.tests.test_io Found 25/25 tests for numpy.lib.tests.test__datasource Found 10/10 tests for numpy.lib.tests.test_arraysetops Found 1/1 tests for numpy.lib.tests.test_financial Found 4/4 tests for numpy.lib.tests.test_polynomial Found 6/6 tests for numpy.lib.tests.test_index_tricks Found 49/49 tests for numpy.lib.tests.test_shape_base Found 55/55 tests for numpy.lib.tests.test_function_base Found 5/5 tests for numpy.lib.tests.test_getlimits Found 3/3 tests for numpy.linalg.tests.test_regression Found 89/89 tests for numpy.linalg.tests.test_linalg Found 96/96 tests for numpy.ma.tests.test_core Found 19/19 tests for numpy.ma.tests.test_mrecords Found 15/15 tests for numpy.ma.tests.test_extras Found 4/4 tests for numpy.ma.tests.test_subclassing Found 36/36 tests for numpy.ma.tests.test_old_ma Found 1/1 tests for numpy.oldnumeric.tests.test_oldnumeric Found 7/7 tests for numpy.tests.test_random Found 16/16 tests for numpy.testing.tests.test_utils Found 6/6 tests for numpy.tests.test_ctypeslib ............................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................ ---------------------------------------------------------------------- Ran 1292 tests in 1.237s OK <unittest._TextTestResult run=1292 errors=0 failures=0> and: >>> import numpy >>> numpy.test() Running unit tests for numpy NumPy version 1.2.0.dev5611 NumPy is installed in /home/facets/doutriaux1/cdat/latest/lib/python2.5/site-packages/numpy Python version 2.5.2 (r252:60911, Aug 4 2008, 13:47:12) [GCC 4.1.1 20070105 (Red Hat 4.1.1-51)] nose version 0.10.3 ................................................................................................................................................................................................................................................................................................................................................................................................................................S.SSS.SSSSSS.SS..S.SSSSSS.SS.S.SSSS.SS.SSS.SSS.SS.S..........SS.SSS.S..S..................S......................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................................... ---------------------------------------------------------------------- Ran 1924 tests in 7.046s OK (SKIP=45) <nose.result.TextTestResult run=1924 errors=0 failures=0>

2 5

1.2 tasks
by Jarrod Millman Aug. 5, 2008

Aug. 5, 2008

Here are the remaining tasks that I am aware of that need to be done before tagging 1.2.0b1 on the 8th. Median ====== The call signature for median needs to change from def median(a, axis=0, out=None, overwrite_input=False): to def median(a, axis=None, out=None, overwrite_input=False): in both numpy/ma/extras.py and numpy/lib/function_base.py Histogram ======== The call signature for histogram needs to change from def histogram(a, bins=10, range=None, normed=False, weights=None, new=False): to def histogram(a, bins=10, range=None, normed=False, weights=None, new=True): in numpy/lib/function_base.py Documentation ============ The documentation project needs to merge in its changes. Stefan will take care of this on the 5th. Please let me know ASAP if there is anything I am missing. Thanks, -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

6 14

"import numpy" is slow
by Andrew Dalke Aug. 5, 2008

Aug. 5, 2008

(Trying again now that I'm subscribed. BTW, there's no link to the subscription page from numpy.scipy.org .) The initial 'import numpy' loads a huge number of modules, even when I don't need them. Python 2.5 (r25:51918, Sep 19 2006, 08:49:13) [GCC 4.0.1 (Apple Computer, Inc. build 5341)] on darwin Type "help", "copyright", "credits" or "license" for more information. >>> import sys >>> len(sys.modules) 28 >>> import numpy >>> len(sys.modules) 256 >>> len([s for s in sorted(sys.modules) if 'numpy' in s]) 127 >>> numpy.__version__ '1.1.0' >>> As a result, I assume that's the reason my program's startup cost is quite high. [josiah:~/src/fp] dalke% time python -c 'a=4' 0.014u 0.038s 0:00.05 80.0% 0+0k 0+1io 0pf+0w [josiah:~/src/fp] dalke% time python -c 'import numpy' 0.161u 0.279s 0:00.44 97.7% 0+0k 0+9io 0pf+0w My total runtime is something like 1.4 seconds, and the only thing I'm using NumPy for is to make an array of doubles that I can pass to a C extension. (I could use the array module or ctypes, but figured numpy is more useful for downstream code.) Why does numpy/__init__.py need to import all of these other modules and submodules? Any chance of cutting down on the number, in order to improve startup costs? Andrew dalke(a)dalkescientific.com

17 83

Toward a numpy build with warning: handling unused variable
by David Cournapeau Aug. 5, 2008

Aug. 5, 2008

Hi, I have recently tried various warning flags with numpy, and already discover several potential errors in the C code. Now, because of the huge amount of warning we get with e.g. -W -Wall, it would be nice to clean the code for this kind of warning 'level' because we have a huge amount of garbage right now at this level. One big proportion of warning with -W is the unused variable, which happens with the dummy/self variable in the typical C functions. I know of two ways to handle this: - the simple: use the variable in some way, e.g. Add a (void)dummy; in the function code. Problem: it is not really nice, and if later you want to use the variable, it may hide other useful warnings, although I am not sure about the later. - more sophisticated way: using compiler specific handling (gcc attribute) + mangling: http://sourcefrog.net/weblog/software/languages/C/unused.html It is more complicated, is compiler specific, but it prevents accidental use of the variable (thanks to the mangling). The compiler-specific is not that much of a problem IMHO, since this is for warnings (and most of use use gcc ?). What do people think ? If nobody has anything against it, I will prepare a patch to support this. cheers, David

7 15

boolean defaults for svd
by Keith Goodman Aug. 4, 2008

Aug. 4, 2008

svd uses 1 for its defaults: svd(a, full_matrices=1, compute_uv=1) Anyone interested in changing 1 to True? It shouldn't break any code, right?

1 0

arrays : c aligmen to fortran and back
by Thomas Hrabe Aug. 4, 2008

Aug. 4, 2008

Hi all, I need to convert a 3d array from c alingment to fortran and was wandering for the simplest code available. Its all about python and an interfaced C++ program, which, however, processes fortran aligned arrays. Is there a simple code from converting an array arr = numpy.array((30,20,10)) to an array arr2 with the same shape (arr.shape == arr2.shape) = true where the data is aligned in fortran style. Furthermore, has anybody ever seen a C library which converts arrays from C alignment to Fortran without swapping the dimensions ([30,20,10] to [10,20,30]). Thank you in advance for your help, Thomas

2 3

WTFM!
by Joe Harrington Aug. 4, 2008

Aug. 4, 2008

SciPy Documentation Marathon 2008 Status Report We are now nearing the end of the summer. We have a ton of great docstrings, a nice PDF and HTML reference guide, a new package with pages on general topics like slicing, and a glossary. We had hoped to have all the numpy docstrings in first-draft form in time for the pre-fall (1.2) release. The actual number of pages was more than double our quick goal-setting assessment, so we won't make it. As of this moment, we have: status % pages Needs editing 52 430 Being written / Changed 27 226 Needs review 18 152 Needs review (revised) 0 1 Needs work (reviewed) 0 3 Reviewed (needs proof) 2 19 Proofed 0 0 Unimportant 1531 Our current status can always be seen at: http://sd-2116.dedibox.fr/pydocweb/stats/ Definitions of the categories are also on the wiki, but "being written" is our first-draft category. So, we're just shy of halfway there, and since the goal more than doubled, we can say we have not failed our expectations. But, we haven't succeeded, either, and we certainly haven't finished. So, this being a marathon, we're not going to stop! Please join us if you haven't already for the... PRE-CONFERENCE DOC BLITZ! We can, quite realistically, get up to 60% for numpy 1.2. We've had several 8% weeks this summer and we've got several weeks to go. Stefan will merge docstrings into the beta for 1.2 on 5 August, and will continue merging from the wiki for the release candidates and final cut. However, Writing has slowed to a crawl in recent weeks. Please pitch in to help those who are still writing so we can get to 60% by release 1.2. Looking further ahead, I hope all the volunteers will continue writing for the rest of the summer and fall, so that we can put 100% decent drafts into 1.3, and a 100% reviewed set of docstrings into 1.4. Then we can turn our attention to scipy. Enough stick, here's some carrot: This is the design for this year's documentation prize, a T-shirt in Robert Kern black, designed by Teresa Jeffcott: http://physics.ucf.edu/~jh/scipyshirt-2008-2.png We'll hand these out at SciPy '08 to anyone who has written 1000 words or more (according to the stats page) or who has made an equivalent contribution in other ways (reviewing, wiki creation, etc.; Stefan and I will judge). So far, 11 contributors qualify, but several more could easily reach that goal in time. In fact, several of our volunteer writers have produced 1000 words in one week. The offer remains good through the first-draft phase, though you'll need to act quickly to get your docs into 1.2 and be recognized at the conference! If you won't be at SciPy '08 or if you qualify later, we'll mail one to you. As always, further discussion belongs on the scipy-dev mailing list, as does your request to be added to the doc wiki editors group, so head over to the doc wiki main page, established an account there, and give us a shout! WTFM, --jh--

1 0

What/when is the next release ?
by Pierre GM Aug. 4, 2008

Aug. 4, 2008

All, * What will be the next release: 1.1.x (w/ support of Python 2.3) or 1.2.0 (w/o support to Python 2.3) ? I need to correct an elusive bug in numpy.ma.MaskedArray and wondered whether I needed to backport it or not (1.1.x would likely be affected as well). Thanks a lot in advance, P.

4 8

concatenate an array with a number
by Nicolas Chopin Aug. 4, 2008

Aug. 4, 2008

Hi list, I want to do this: y = concatenate(x[0],x) where x is a 1d array. However, the only way I managed to make this work is: y = concatenate( array((x[0],),ndmin=1), x) which is kind of cryptic. (If I remove ndmin, numpy still complains.) 1. Is there a better way? 2. Would it be possible to allow numbers as arguments for concatenate? Many thanks NC

4 3

Renaming MaskedArray._basedict
by Pierre GM Aug. 4, 2008

Aug. 4, 2008

All, I intend to rename a little-known attribute of MaskedArray, _basedict, to something more understandable such as _optinfo. This attribute is a dictionary that comes handy when subclassing MaskedArray or to save information. _basedict would still be available for a while as an alias to _optinfo, but may disappear in the long run. I intend to put a DeprecationWarning at the beginning of numpy.ma.core. A side-effect is that this warning is issued each time numpy.ma is imported, which clutters the output of the tests. Am I missing something ? P.

3 4

ANN: NumPy 1.1.1
by Jarrod Millman Aug. 2, 2008

Aug. 2, 2008

I'm pleased to announce the release of NumPy 1.1.1. NumPy is the fundamental package needed for scientific computing with Python. It contains: * a powerful N-dimensional array object * sophisticated (broadcasting) functions * basic linear algebra functions * basic Fourier transforms * sophisticated random number capabilities * tools for integrating Fortran code. Besides it's obvious scientific uses, NumPy can also be used as an efficient multi-dimensional container of generic data. Arbitrary data-types can be defined. This allows NumPy to seamlessly and speedily integrate with a wide-variety of databases. Numpy 1.1.1 is a bug fix release featuring major improvements in Python 2.3.x compatibility and masked arrays. For information, please see the release notes: http://sourceforge.net/project/shownotes.php?group_id=1369&release_id=617279 Thank you to everybody who contributed to this release. Enjoy, -- Jarrod Millman Computational Infrastructure for Research Labs 10 Giannini Hall, UC Berkeley phone: 510.643.4014 http://cirl.berkeley.edu/

2 2

User-defined data-types
by Francesc Alted Aug. 1, 2008

Aug. 1, 2008

Hi, We need to make a dtype that needs to be defined in the user's application space. In the NumPy book there is a section devoted to this subject, but also warns: """ Adding data-types is one of the less well-tested areas for NumPy 1.0, so there may be bugs remaining in the approach. """ Anybody out there has tried this successfully? Thanks, -- Francesc Alted

3 2

[ANNOUNCE] ETS 3.0.0b1 released!
by Dave Peterson July 31, 2008

July 31, 2008

Hello, I am very pleased to announce that ETS 3.0.0b1 has just been tagged and released! All ETS sub-projects, including Traits, Chaco, Mayavi, and Envisage, have been tagged and released in one form or another (alpha, beta, final) as part of this effort. All ETS projects have now been registered with PyPi and source tarballs have been uploaded! As a result, for the first time ever, you can now easy_install all of ETS via a simple command like: easy_install ETS Or, alternatively, you can get only specific projects with commands like: easy_install Mayavi Or easy_install Chaco We anticipate the final release of ETS 3.0, and all of its sub-projects at their current version number, in a matter of weeks. Please join us in helping make this final release the most robust release it can be by participating in this beta shakeout period. Discussions about the status of ETS and its sub-projects happens on the enthought-dev mailing list, who's home page is at: https://mail.enthought.com/mailman/listinfo/enthought-dev Thank you to all who have worked so hard on making this release possible! -- Dave Enthought Tool Suite --------------------------- The Enthought Tool Suite (ETS) is a collection of components developed by Enthought and open source participants, which we use every day to construct custom scientific applications. It includes a wide variety of components, including: * an extensible application framework * application building blocks * 2-D and 3-D graphics libraries * scientific and math libraries * developer tools The cornerstone on which these tools rest is the Traits package, which provides explicit type declarations in Python; its features include initialization, validation, delegation, notification, and visualization of typed attributes. More information is available for all the packages within ETS from the Enthought Tool Suite development home page at http://code.enthought.com/projects/tool-suite.php. Testimonials ---------------- "I set out to rebuild an application in one week that had been developed over the last seven years (in C by generations of post-docs). Pyface and Traits were my cornerstones and I knew nothing about Pyface or Wx. It has been a hectic week. But here ... sits in front of me a nice application that does most of what it should. I think this has been a huge success. ... Thanks to the tools Enthought built, and thanks to the friendly support from people on the [enthought-dev] list, I have been able to build what I think is the best application so far. I have built similar applications (controlling cameras for imaging Bose-Einstein condensate) in C+MFC, Matlab, and C+labWindows, each time it has taken me at least four times longer to get to a result I regard as inferior. So I just wanted to say a big "thank you". Thank you to Enthought for providing this great software open-source. Thank you for everybody on the list for your replies." — Gaël Varoquaux, Laboratoire Charles Fabry, Institut d’Optique, Palaiseau, France "I'm currently writing a realtime data acquisition/display application … I'm using Enthought Tool Suite and Traits, and Chaco for display. IMHO, I think that in five years ETS/Traits will be the most comonly used framework for scientific applications." — Gary Pajer, Department of Chemistry, Biochemistry and Physics, Rider University, Lawrenceville NJ

2 1