Numeric help!

Sheldon shejo284 at gmail.com
Fri Jun 30 11:46:14 CEST 2006


Carl Banks wrote:
> Sheldon wrote:
> >Carl Banks wrote:
> >> I'm not sufficiently sure this isn't a homework problem, so here's a
> >> partial answer.
> [snip]
> >
> > My days as a student is over for the most part. I am learning python on
> > my own and Numeric is not properly documented so I am learning by doing
> > and copying from others.
>
> Are you aware of this guide?
>
> http://numeric.scipy.org/numpydoc/numdoc.htm
>
> I've seen better documents but it's fairly complete--I wouldn't call it
> improper.
>
>
> > I thought about the problem and a solution
> > another problem given to me earlier using "putmask" is the solution but
> > there is a bug that I cannot figure out:
> > ************************
> >             index = 0
> >             for z in range_va:
> >                 wk     = msgva # working arrary
> >                 sattmp = sat_id # working array
>
> This appears to be doing something you don't expect.  In Python all
> names are references; assignments don't create new objects but rather
> new names for the same object.  In the above, wk and msvga are the same
> array: any changes in wk also appear in msvga.  So wk[1,1] = 4 also
> causes msgva[1,1] to be 4.
>
> Numeric arrays even share data when slicing.  If you were to take a
> slice of wk, for example, a = wk[0:10,0:10], then modifying a would
> also show up in wk (and msvga).  (However, slicing most other objects
> copies rather than shares the data.)
>
> It's completely different from Matlab, which always copies data.
>
> Here's what you should do:
>
>     wk = array(msgva) # this creates a new array for wk
>     sattmp = array(sat_id)
>
> >                 mask = where(equal(wk,z),1,0) # creating a mask of
> > valid data pixs
>
> Note that the "where" is unnecessary here.  equal creates an array of
> ones and zeros.
>
> >                 putmask(sattmp,mask==0,-999)  # since zero is valid
> > data, -999 is used instead
>
> This would overwrite sat_id unless you copy the array as I've shown
> above.
>
> >                 rdata =
> > compress(ravel(not_equal(sattmp,-999)),ravel(sattmp))
>
> I believe you could avoid the above putmask above step and just use
> ravel(mask) as the first argument of compress here, or even just
> equal(wk,z).
>
> "equal(wx,z)", "mask", and "not_equal(sattmp,-999)" are all equal
> arrays.
>
> >                 if sum(sum(rdata)) == 0:
> >                     av = 0
>
> rdata is one-dimensional, so you only need one sum call; but why do
> this?  average will still work if the sum is zero.
>
> >                 else:
> >                     av = average(rdata,axis=None)
>
> axis argument isn't necessary here since rdata is one-dimesional.
>
> >                 tmparray[0,index] = av
> >                 tmparray[1,index] = z
> >                 index += 1
> >             print tmparray
> > ***********************************
> > But the tmparray is returning zeros as averages. When I try just one
> > value for z everything works. I can't see where I going wrong. I am not
> > using the original arrays, only the copies and when a new z is chosen
> > then these are recreated.
> > Care to help out with this?
> > /Sheldon
>
> Other than the array sharing mistake, it looks like it should work.
>
> Also, note that this group disdains top-posting; replies should go
> below the quoted text because in these sorts of discussions it's best
> to keep things in conversational order, especially for the sake of
> other interested readers.
>
> Good luck.
>
>
> Carl Banks

Thanks for the tips about the array and how it is copied. I figured
this out late yesterday and wrote a dirty solution to this problem but
your is more elegant and concise.
Yes, where statement is not needed and I removed it. The replacement is
the following:

average(compress(ravel(equal(wk,z)),ravel(sattmp)),axis=None)
This is much more compact and elegant. Thanks for pointing this out.
I don't know why average() returned a divide by zero error and to avoid
this I inserted this if statement.   Now it works much better !

Much obliged,
Sheldon




More information about the Python-list mailing list