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

Hi, I am in the process of converting some code from Numeric to numarray,
and it seems that numarray no longer has the sign() function -- is that
so?
ex:
sign(-30.0) = -1
sign(0) = 0
sign(1000) = 1
Is there a replacement?
Thanks,
Mike

Hi,
I am having trouble understanding how exactly "where" works in
numarray.
What I am trying to do:
I am preparing a two-level mask in an array and then assign values to
the array where both masks are true:
>>> from numarray import *
>>> a = arange(10)
>>> # First mask
>>> m1 = where(a > 5)
>>> a[m1]
array([6, 7, 8, 9])
>>> # Second mask
>>> m2 = where(a[m1] < 8)
>>> a[m1][m2]
array([6, 7])
>>> # So far so good
>>> # Now change some values
>>> a[m1][m2] = array([10, 20])
>>> a[m1][m2]
array([6, 7])
>>> a
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> # Didn't work
>>> # Let's try a temporary variable
>>> t = a[m1]
>>> t[m2]
array([6, 7])
>>> t[m2] = array([10, 20])
>>> t[m2], t
(array([10, 20]), array([10, 20, 8, 9]))
>>> a
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
So, my assignment to a[m1][m2] seems to work (no messages), but it
doesn't produce the effect I want it to.
I have read the documentation but I couldn't find something that would
explain this behavior.
So my questions:
- did I miss something important in the documentation,
- I am expecting something I shouldn't, or
- there is a bug in numarray?
Thanks,
Alok
--
Alok Singhal (as8ca(a)virginia.edu) __
Graduate Student, dept. of Astronomy / _
University of Virginia \_O \
http://www.astro.virginia.edu/~as8ca/ __/

I have a series of luminance images that I want to do some correlation
analyses on. Each image is an MxN frame of a movie. I want to
compute the correlation between a given pixel i,j, and every other
pixel in the image over each frame. That is, if I assume
xij is a numFrames length time series, and xkl is another numFrames
length time series (the pixel intensities at two points in the image
over time), I want to compute the
corrcoeff(xij, xkl) for every kl with ij fixed.
I know I could do this by looping over the pixels of the image, but
I'm hoping for something a bit faster.
Any suggestions?
John Hunter

I am wondering if there is any way to set an array as read-only,
so that any trial to modify values of elements in the array raises some
warning.
Is it a cool feature to prevent unexpected side effect?
Daehyok Shin (Peter)

Hi,
I'm willing to use a lot the searchsorted function in numarray, but I'm a
bit surprised about the poor performance of it. Look at that:
>>> from time import time
>>> import numarray
>>> import Numeric
>>> na=numarray.arange(1000*1000)
>>> nu=Numeric.arange(1000*1000)
>>> t1=time();numarray.searchsorted(na,200*1000);time()-t1
200000
0.00055098533630371094
>>> t1=time();Numeric.searchsorted(nu,200*1000);time()-t1
200000
7.7962875366210938e-05
It may seem that Numeric is better optimised, but the standard python module
bisect is even faster than numarray.searchsorted:
>>> import bisect
>>> t1=time();bisect.bisect_left(na,200*1000);time()-t1
200000
8.8930130004882812e-05
>>> t1=time();bisect.bisect_left(nu,200*1000);time()-t1
200000
8.6069107055664062e-05
So, bisect performance is similar to that of Numeric searchsorted and both
are almost an order of magnitude better than numarray searchsorted. This a
little bit surprising as the bisect_left module is written in plain python.
From the python 2.3.3 sources:
def bisect_left(a, x, lo=0, hi=None):
if hi is None:
hi = len(a)
while lo < hi:
mid = (lo+hi)//2
if a[mid] < x: lo = mid+1
else: hi = mid
return lo
I'm using python2.3.3, numarray 0.9.1 (latest CVS) and Debian Linux (sid).
Cheers,
--
Francesc Alted

Hi,
the random_array poisson functions returns negative values if mean=0:
>>>from numarray import random_array as ra
>>> ra.seed(x=1, y=1)
>>> ra.poisson(0)
5
>>> ra.poisson(0)
-2
My "math book" tells me that it should be always zero.
This seems to be a constructed case, but I'm using this to put
"quantum statistic" into a simulated image:
obj = na.array( something )
imageFromDetector = ra.poisson( obj ) + gaussianNoiseArray
The object array might have lots of zeros surrounding the "actual object".
Thinking of a fluorescent object sending out photons it makes sense to not get
any photons at all from 'empty' regions.
I'm using numarray 0.8;
Thanks for numarray,
Sebastian Haase

The simplest solution is to use "solve_linear_equations" functin from numarray.linear_algebra. For more sophisticate solution you may try scipy
(www.scipy.org)
Nadav
-----Original Message-----
From: Craig Rasmussen [mailto:crasmussen@lanl.gov]
Sent: Wed 19-May-04 17:59
To: numarray
Cc:
Subject: [Numpy-discussion] Linear solvers
I need to solve a matrix equation of the form
A u = q
for u, where u and q are N-vectors and A is a
symmetric-positive-definite matrix.
For now, at least, A is a 2-D matrix.
I would like to use Python and numarray. Are there Python modules that
I could use to obtain a solution for u?
Thanks,
Craig
-------------------------------------------------------
This SF.Net email is sponsored by: SourceForge.net Broadband
Sign-up now for SourceForge Broadband and get the fastest
6.0/768 connection for only $19.95/mo for the first 3 months!
http://ads.osdn.com/?ad_id=2562&alloc_id=6184&op=click
_______________________________________________
Numpy-discussion mailing list
Numpy-discussion(a)lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/numpy-discussion

I have a matrix of particle collision times: times[i,j] gives the time
for particle "i" to collide with particle "j".
If I do, in order to get the first expected collision time for each
particle, the following (random array for testing purposes) :
>>> N=30
>>> times = rnd.random([N,N])
>>> choose(argmin(times,transpose(times))
Segmentation fault (of the python interactive shell!!!)
With N=100 I get a more informative traceback, and rest within the
shell:
Traceback (most recent call last):
File "<stdin>", line 1, in ?
File "/usr/lib/python2.3/site-packages/numarray/ufunc.py", line 1670, in choose
return _choose(selector, population, outarr, clipmode)
File "/usr/lib/python2.3/site-packages/numarray/ufunc.py", line 1579, in __call__
result = self._doit(computation_mode, woutarr, cfunc, ufargs, 0)
File "/usr/lib/python2.3/site-packages/numarray/ufunc.py", line 1564, in _doit
blockingparameters)
ValueError: _operator_compute: too many inputs + outputs
For N=10,N=15 I get the expected output, but for N=20 I get again the
brutal segfault...
regards,
á.
--
Álvaro Tejero Cantero
http://alqua.org -- documentos libres
free documents