Thanks. Agreed it will break down under that scenario but that shouldn't be
encountered as I am simply checking the value is greater than what I have
set the undefined to be (-999.0). Here is the refjigged logic using numpy
functionality for anyone who it might help. Would appreciate any suggestions
how I can further improve this. The code is significantly more efficient now
(thankfully ;P)!
#!/usr/bin/env python
"""
Average the daily LST values to make an 8 day product...
Martin De Kauwe
26th November 2009
"""
import sys, os, glob
import numpy as np
def averageEightDays(files, lst, count, numrows, numcols):
day_count = 0
doy = 122
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()
# copy the first file into the array...
# increase the pixel count if we have some valid data at
timestep 1
if day_count == 0:
lst = data
count = np.where(lst > -900.0, 1, 0)
day_count += 1
# keep a running total of valid lst values in an 8 day sequence
# need to check whether the lst pixel is valid
elif day_count > 0 and day_count <= 7:
lst = np.where(data > -900.0, np.where(lst < -900.0,
data, lst + data),
lst)
count = np.where(data > -900.0, np.where(lst < -900.0,
1, count + 1),
count)
day_count += 1
# need to average previous 8 days and write to a file...
if day_count == 8:
print year, ':', doy
lst = np.where(count > 0, lst / np.float32(count),
-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:
lst = np.where(count > 0, lst / np.float32(count), -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
lst.tofile(of)
count.tofile(of2)
of.close()
of2.close()
lst = np.zeros((numrows, numcols), dtype=np.float32)
count = np.zeros((numrows, numcols), dtype=np.int)
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 = np.zeros((numrows, numcols), dtype=np.float32)
count = np.zeros((numrows, numcols), dtype=np.int)
averageEightDays('lst_scr_2007*.gra', lst, count, numrows, numcols)
Pierre GM-2 wrote:
>
> On Nov 25, 2009, at 4:13 PM, mdekauwe wrote:
>> I tried redoing the internal logic for example by using the where
>> function
>> but I can't seem to work out how to match up the logic. For example (note
>> slightly different from above). If I change the main loop to
>>
>> lst = np.where((data > -900.0) & (lst < -900.0), data, lst)
>> lst = np.where((data > -900.0) & (lst > -900.0), 5.0, lst)
>>
>> If I then had a data array value of 10.0 and lst array value of -999.0.
>> In
>> this scenario the first statement would set the LST value to 10.0.
>> However
>> when you run the second statement, data is still bigger than the
>> undefined
>> (-999.0) but now so is the LST, hence the lst is set to 5.0
>>
>> This doesn't match with the logic of
>>
>> if data[i,j] > -900.:
>> if lst[i,j] < -900.:
>> lst[i,j] = data[i,j]
>> elif lst[i,j] > -900.:
>> lst[i,j] = 5.0
>>
>> as it would never get to the second branch under this logic.
>>
>> Does anyone have any suggestions on how to match up the logic?
>>
>
>
> np.where(data>-900,np.where(lst<-900,data,5.),lst)
>
> Note that this should fail if lst=-900
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
--
View this message in context:
http://old.nabble.com/Help-making-better-use-of-numpy-array-functions-tp26503657p26525315.html
Sent from the Numpy-discussion mailing list archive at Nabble.com.
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
