[Numpy-discussion] Help making better use of numpy array functions
mdekauwe
mdekauwe@gmail....
Tue Nov 24 15:33:55 CST 2009
Hi I have written some code and I would appreciate any suggestions to make
better use of the numpy arrays functions to make it a bit more efficient and
less of a port from C. Any tricks are thoughts would be much appreciated.
The code reads in a series of images, collects a running total if the value
is valid and averages everything every 8 images.
Code below...
Thanks in advance...
#!/usr/bin/env python
"""
Average the daily LST values to make an 8 day product...
Martin De Kauwe
18th November 2009
"""
import sys, os, glob
import numpy as np
import matplotlib.pyplot as plt
def averageEightDays(files, lst, count, numrows, numcols):
day_count = 0
# starting day - tag for output file
doy = 122
# loop over all the images and average all the information
# every 8 days...
for fname in glob.glob(os.path.join(path, files)):
current_file = fname.split('/')[8]
year = current_file.split('_')[2][0:4]
try:
f = open(fname, 'rb')
except IOError:
print "Cannot open outfile for read", fname
sys.exit(1)
data = np.fromfile(f, dtype=np.float32).reshape(numrows, numcols)
f.close()
# if it is day 1 then we need to fill up the holding array...
# copy the first file into the array...
if day_count == 0:
lst = data
for i in xrange(numrows):
for j in xrange(numcols):
# increase the pixel count if we got vaild data
(undefined = -999.0
if lst[i,j] > -900.:
count[i,j] = 1
day_count += 1
# keep a running total of valid lst values in an 8 day sequence
elif day_count > 0 and day_count <= 7:
for i in xrange(numrows):
for j in xrange(numcols):
# lst valid pixel?
if data[i,j] > -900.:
# was the existing pixel valid?
if lst[i,j] < -900.:
lst[i,j] = data[i,j]
count[i,j] = 1
else:
lst[i,j] += data[i,j]
count[i,j] += 1
day_count += 1
# need to average previous 8 days and write to a file...
if day_count == 8:
for i in xrange(numrows):
for j in xrange(numcols):
if count[i,j] > 0:
lst[i,j] = lst[i,j] / count[i,j]
else:
lst[i,j] = -999.0
day_count = 0
doy += 8
lst, count = write_outputs(lst, count, year, doy)
# it is possible we didn't have enough slices to do the last 8day avg...
# but if we have more than one day we shall use it
# need to average previous 8 days and write to a file...
if day_count > 1:
for i in xrange(numrows):
for j in xrange(numcols):
if count[i,j] > 0:
lst[i,j] = lst[i,j] / count[i,j]
else:
lst[i,j] = -999.0
day_count = 0
doy += 8
lst, count = write_outputs(lst, count, year, doy)
def write_outputs(lst, count, year, doy):
path = "/users/eow/mgdk/research/HOFF_plots/LST/8dayLST"
outfile = "lst_8day1030am_" + str(year) + str(doy) + ".gra"
try:
of = open(os.path.join(path, outfile), 'wb')
except IOError:
print "Cannot open outfile for write", outfile
sys.exit(1)
outfile2 = "pixelcount_8day1030am_" + str(year) + str(doy) + ".gra"
try:
of2 = open(os.path.join(path, outfile2), 'wb')
except IOError:
print "Cannot open outfile for write", outfile2
sys.exit(1)
# empty stuff and zero 8day counts
lst.tofile(of)
count.tofile(of2)
of.close()
of2.close()
lst = 0.0
count = 0.0
return lst, count
if __name__ == "__main__":
numrows = 332
numcols = 667
path = "/users/eow/mgdk/research/HOFF_plots/LST/gridded_03/"
lst = np.zeros((numrows, numcols), dtype=np.float32)
count = np.zeros((numrows, numcols), dtype=np.int)
averageEightDays('lst_scr_2006*.gra', lst, count, numrows, numcols)
lst = 0.0
count = 0.0
averageEightDays('lst_scr_2007*.gra', lst, count, numrows, numcols)
--
View this message in context: http://old.nabble.com/Help-making-better-use-of-numpy-array-functions-tp26503657p26503657.html
Sent from the Numpy-discussion mailing list archive at Nabble.com.
More information about the NumPy-Discussion
mailing list