General Numerical Python question

[Michael Ressler]

In article <%xVib.26249$gi2.16759 at fed1read01>, Tim Hochberg wrote:
> Michael Ressler wrote:
>> This is "pseudo-code", so don't try to run it, but see if the ideas
>> are useful. One way to approach "running" things, is extensive use of

How I defend myself :-)

> I agree, you probably can't get rid of all the loops in this case, and 
> if you could, the resulting code would probably be horrible. I have a 
> couple of minor quibles with the above code though. I think I'd write it as:
> 
> lenavg =  len(a) - n + 1
> avg  =np.zeros(lenavg, np.Float)
> for i in range(n) : # window size to smooth over
>      avg += a[i:lenavg+i]  # Using += reuses the same array every time
>                       # Instead of creating a new one each time
>                       # Through the loop.
> avg /= n  # Same here.
> 
> The important point being the use of += and /=. And, in order to make 
> that work, you need to set the type of avg appropriately, not let it 
> default to int.

You are right, of course, but I was trying to be clear, not fast. The
point is that you could replace that sum with the sum of the squares
and crossproducts, and thus get a running standard deviation, for
example. The trick of Numeric is thinking about problems in new ways
which are far different from what the original author is used to with
loops.

Another example of thinking things differently is suppose you have a
vector where the values are randomly positive or negative. Suppose for
reasons known only to you, you want to replace the negative values
with the sqrt of their absolute values. With Numeric, no loops are
involved.

from Numeric import *
a=array([1.,2.,-3.,4.,-5.,6.,-7.,-8.,9.])	# make up an array
idx=nonzero(a<0)			# indexes of the negative values
sqrs=sqrt(abs(take(a,idx)))		# get the sqrts of neg elements
put(a,idx,sqrs)				# put them back into a
print a					# works!  

You can make the whole thing a one-liner if you want to get carried
away with it. It's too bad "nonzero" isn't called "whereis" or
something like that - it would make the idx= line more obvious.

Mike

-- 
Dr. Michael Ressler
Research Scientist, Astrophysics Research Element, Jet Propulsion Laboratory
Email: ressler at cheetah.jpl.nasa.gov      Phone: (818)354-5576
"A bad night at the telescope is still better than the best day at the office."