On Wed, Nov 25, 2009 at 4:13 PM, mdekauwe <mdeka...@gmail.com> 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)
does an intermediate array help? lst2 = np.where((data > -900.0) & (lst < -900.0), data, lst) lst = np.where((data > -900.0) & (lst > -900.0), 5.0, lst2) I didn't read through all the loops to see what would be a better way to do it. Josef > > 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? > > > > mdekauwe wrote: >> >> 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-tp26503657p26520383.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 > _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion