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-tp2650... Sent from the Numpy-discussion mailing list archive at Nabble.com.