[Numpy-discussion] MaskedArray problem?

[Álvaro Tejero Cantero]

Hello,

I'm trying to rewrite my simulation using masked arrays. For that I
define masks for triangular_inferior NxN and
diagonal_triangular_superior NxN arrays.

I'm unsuccesfully trying to know what happens with the program when I
operate over a MaskedArray: it always says

Traceback (most recent call last):
  File "<stdin>", line 31, in ?
  File "<stdin>", line 26, in run
  File "/usr/lib/python2.3/site-packages/numarray/ma/MA.py", line 990,
in __pow__
    return power(self, other, third)
  File "/usr/lib/python2.3/site-packages/numarray/ma/MA.py", line 1585,
in power
    return masked_array(Numeric.power(fa, fb), m)
  File "/usr/lib/python2.3/site-packages/numarray/ufunc.py", line 877,
in _cache_miss2
    mode, win1, win2, wout, cfunc, ufargs = \
  File "/usr/lib/python2.3/site-packages/numarray/ufunc.py", line 921,
in _setup
    intypes = (in1._type.name, in2._type.name)
AttributeError: 'ObjectArray' object has no attribute '_type'


For reference, I copy a simplified version of the program. The statement
that causes the AttributeError is the last from run(N,D,S):

---------------------------

import numarray.ma as mka
import numarray.random_array as rnd
import numarray as nmr

def relative(x, msk):
    (N,D) = x.shape
    return mka.array([nmr.subtract.outer(x[:,i],x[:,i]) for i in nmr.arange(D)],
                     mask = nmr.array([msk for i in nmr.arange(D)],typecode=mka.MaskType))

def scalar(x,y):
    return mka.add.reduce(x*y,axis=0)

def run(D,N,S):
    trig_inf = nmr.array(([[(0,1)[i<j] for i in nmr.arange(N)] for j in nmr.arange(N)]), typecode=mka.MaskType)
    diag_trig_sup = 1 - trig_inf #print mka.is_mask(diag_trig_sup) gives "1"

    # positions and velocities
    r = rnd.random([N,D]) #should be overlap-checked
    v = rnd.random([N,D])

    # relative positions and velocities
    r_rel = relative(r, diag_trig_sup)
    v_rel = relative(v, diag_trig_sup)

    RV, VV, RR = scalar(r_rel,v_rel), scalar(v_rel,v_rel), scalar(r_rel,r_rel)
    discr = RV**2 - VV*(RR-S**2)      #<<<<----------PROBLEM HERE----------------



if __name__ == '__main__':
    run(D=3,N=4,S=1e-8)


-------------------------

Any other advice that you may give me is very wellcome. My deadline is
approaching and I'm really reluctant to rewrite all this in C++ but I
don't see the light.

Thank you in advance,

 á.


-- 
Álvaro Tejero Cantero

http://alqua.org    --	documentos libres
			free documents