Numeric help!

[Sheldon]

Hi,

I have the following loop that I think can be written to run faster in
Numeric. I am currently using Numeric.
range_va = [60,61,62,63,64, 65,66,67,68,69, 70,71,72]
main.xsize= 600
main.ysize= 600
#msgva is an (600x600) Numeric array with mutiple occurrences of the
values in range_va
#sat_id is an (600x600) Numeric array with values ranging from -2 to 2
for z in range_va:
count = 0
mbias = 0
for i in range(main.xsiz e):
for j in range(main.ysiz e):
if msgva[i,j] == z:
mbias += sat_id[i,j] # take the sum of the
values
count += 1 # count the occurrences
tmp_array[0,index] = round(mbias/count,1) # store the mean
tmp_array[1,index] = z
index += 1

Any help would be greatly appreciated!

Sincerely,
Sheldon

[Carl Banks]

Re: Numeric help!

Sheldon wrote:[color=blue]
> Hi,
>
> I have the following loop that I think can be written to run faster in
> Numeric. I am currently using Numeric.
> range_va = [60,61,62,63,64, 65,66,67,68,69, 70,71,72]
> main.xsize= 600
> main.ysize= 600
> #msgva is an (600x600) Numeric array with mutiple occurrences of the
> values in range_va
> #sat_id is an (600x600) Numeric array with values ranging from -2 to 2
> for z in range_va:
> count = 0
> mbias = 0
> for i in range(main.xsiz e):
> for j in range(main.ysiz e):
> if msgva[i,j] == z:
> mbias += sat_id[i,j] # take the sum of the
> values
> count += 1 # count the occurrences
> tmp_array[0,index] = round(mbias/count,1) # store the mean
> tmp_array[1,index] = z
> index += 1
>
> Any help would be greatly appreciated![/color]

I'm not sufficiently sure this isn't a homework problem, so here's a
partial answer.

Your intuition is correct--there's almost always a faster way to do it
when you're cycling through Numeric array indices in a Python for loop.
Numeric and successors exist mostly to get rid of them.

Given z, you can calculate mbias in one line:

mbias = Numeric.sum(
Numeric.ravel(N umeric.where(ms gva==z,sat_id,0 )))

Here, the Numeric.where function allows you to filter out entries: it
returns an array with the elements of sat_id, but only where the
corresponding element of msvga==z; otherwise the entry is zero.
Numeric.ravel flattens an multidimensiona l array into one dimension,
and, of course, Numeric.sum adds up all the elements in the array.

How to get count and to work this into your loop are left as an
exercise.


Carl Banks

[Sheldon]

Re: Numeric help!

Hi Carl,

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. 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
mask = where(equal(wk, z),1,0) # creating a mask of
valid data pixs
putmask(sattmp, mask==0,-999) # since zero is valid
data, -999 is used instead
rdata =
compress(ravel( not_equal(sattm p,-999)),ravel(sat tmp))
if sum(sum(rdata)) == 0:
av = 0
else:
av = average(rdata,a xis=None)
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

Carl Banks wrote:[color=blue]
> Sheldon wrote:[color=green]
> > Hi,
> >
> > I have the following loop that I think can be written to run faster in
> > Numeric. I am currently using Numeric.
> > range_va = [60,61,62,63,64, 65,66,67,68,69, 70,71,72]
> > main.xsize= 600
> > main.ysize= 600
> > #msgva is an (600x600) Numeric array with mutiple occurrences of the
> > values in range_va
> > #sat_id is an (600x600) Numeric array with values ranging from -2 to 2
> > for z in range_va:
> > count = 0
> > mbias = 0
> > for i in range(main.xsiz e):
> > for j in range(main.ysiz e):
> > if msgva[i,j] == z:
> > mbias += sat_id[i,j] # take the sum of the
> > values
> > count += 1 # count the occurrences
> > tmp_array[0,index] = round(mbias/count,1) # store the mean
> > tmp_array[1,index] = z
> > index += 1
> >
> > Any help would be greatly appreciated![/color]
>
> I'm not sufficiently sure this isn't a homework problem, so here's a
> partial answer.
>
> Your intuition is correct--there's almost always a faster way to do it
> when you're cycling through Numeric array indices in a Python for loop.
> Numeric and successors exist mostly to get rid of them.
>
> Given z, you can calculate mbias in one line:
>
> mbias = Numeric.sum(
> Numeric.ravel(N umeric.where(ms gva==z,sat_id,0 )))
>
> Here, the Numeric.where function allows you to filter out entries: it
> returns an array with the elements of sat_id, but only where the
> corresponding element of msvga==z; otherwise the entry is zero.
> Numeric.ravel flattens an multidimensiona l array into one dimension,
> and, of course, Numeric.sum adds up all the elements in the array.
>
> How to get count and to work this into your loop are left as an
> exercise.
>
>
> Carl Banks[/color]

[Carl Banks]

Re: Numeric help!

Sheldon wrote:[color=blue]
>Carl Banks wrote:[color=green]
>> I'm not sufficiently sure this isn't a homework problem, so here's a
>> partial answer.[/color][/color]
[snip][color=blue]
>
> 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.[/color]

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.

[color=blue]
> 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[/color]

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)
[color=blue]
> mask = where(equal(wk, z),1,0) # creating a mask of
> valid data pixs[/color]

Note that the "where" is unnecessary here. equal creates an array of
ones and zeros.
[color=blue]
> putmask(sattmp, mask==0,-999) # since zero is valid
> data, -999 is used instead[/color]

This would overwrite sat_id unless you copy the array as I've shown
above.
[color=blue]
> rdata =
> compress(ravel( not_equal(sattm p,-999)),ravel(sat tmp))[/color]

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(satt mp,-999)" are all equal
arrays.
[color=blue]
> if sum(sum(rdata)) == 0:
> av = 0[/color]

rdata is one-dimensional, so you only need one sum call; but why do
this? average will still work if the sum is zero.
[color=blue]
> else:
> av = average(rdata,a xis=None)[/color]

axis argument isn't necessary here since rdata is one-dimesional.
[color=blue]
> 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[/color]

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

[Alan Isaac]

Re: Numeric help!

On Thu, 29 Jun 2006 05:25:25 -0400, Sheldon <shejo284@gmail .com> wrote:[color=blue]
> I am learning python on
> my own and Numeric is not properly documented[/color]

1. Use NumPy (Numeric's successor): http://www.numpy.org/
2. Documentation is excellent: http://www.tramy.us/
(Also see http://www.scipy.org/Cookbook )

Cheers,
Alan Isaac

[Sheldon]

Re: Numeric help!

Carl Banks wrote:[color=blue]
> Sheldon wrote:[color=green]
> >Carl Banks wrote:[color=darkred]
> >> I'm not sufficiently sure this isn't a homework problem, so here's a
> >> partial answer.[/color][/color]
> [snip][color=green]
> >
> > 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.[/color]
>
> 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.
>
>[color=green]
> > 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[/color]
>
> 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)
>[color=green]
> > mask = where(equal(wk, z),1,0) # creating a mask of
> > valid data pixs[/color]
>
> Note that the "where" is unnecessary here. equal creates an array of
> ones and zeros.
>[color=green]
> > putmask(sattmp, mask==0,-999) # since zero is valid
> > data, -999 is used instead[/color]
>
> This would overwrite sat_id unless you copy the array as I've shown
> above.
>[color=green]
> > rdata =
> > compress(ravel( not_equal(sattm p,-999)),ravel(sat tmp))[/color]
>
> 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(satt mp,-999)" are all equal
> arrays.
>[color=green]
> > if sum(sum(rdata)) == 0:
> > av = 0[/color]
>
> rdata is one-dimensional, so you only need one sum call; but why do
> this? average will still work if the sum is zero.
>[color=green]
> > else:
> > av = average(rdata,a xis=None)[/color]
>
> axis argument isn't necessary here since rdata is one-dimesional.
>[color=green]
> > 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[/color]
>
> 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[/color]

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(compres s(ravel(equal(w k,z)),ravel(sat tmp)),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

[Carl Banks]

Re: Numeric help!

Sheldon wrote:[color=blue]
> average(compres s(ravel(equal(w k,z)),ravel(sat tmp)),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 ![/color]

Probably you had a case where the array length was zero, but it
wouldn't happen in the present case unless your input arrays are zero
by zero.


Carl Banks

[Sheldon]

Re: Numeric help!


Carl Banks wrote:

Sheldon wrote:

average(compres s(ravel(equal(w k,z)),ravel(sat tmp)),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 !

>
Probably you had a case where the array length was zero, but it
wouldn't happen in the present case unless your input arrays are zero
by zero.
>
>
Carl Banks

Thanks for you help kind Sir!

/Sheldon