[Numpy-discussion] Help making better use of numpy array functions

josef.pktd@gmai... josef.pktd@gmai...
Wed Nov 25 15:22:34 CST 2009


On Wed, Nov 25, 2009 at 4:13 PM, mdekauwe <mdekauwe@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
>


More information about the NumPy-Discussion mailing list