Python code for selecting and averaging raster layers within subfolders

DeeMcElligott1

Hi there,

Basically I have a folder with ten subfolders (one for each year of data...1992-2012), within each of these subfolders are 52 satellite images, one for each week of the corresponding year. I want to convert the weekly data into monthly data, thereby resulting in a single output for each month and twelve for each year. So I want to run a script that picks every four images consecutively, i.e tell it: the first four images = Jan, second four = Feb ............and so on, and get the average of these four images.


I can get the average by running the tool 'Cell Statistic' in Model builder and use 'Calculate Value' to allow me to enter the Python script to tell 'Cell Statistic' which four images to average.

I'm afraid my Python scripting capabilities are rather limited. Does anyone know a Python code to enable me to do this?

I hope it makes sense!

Thanks

Dee

ThomasEmge

for ArcGIS 10.0 the code should look somewhat like this

# we are using the python calendar module to help us out with the weeks
import calendar 

# define your current workspace as the arcpy listing function use this information
arcpy.env.workspace = r'd:\archive\1992'

# let's generate a dictionary with the names for the months
month_dict = dict((v,k) for k,v in enumerate(calendar.month_name))

# within the workspace let's get all the rasters
# the assumption here is that the naming of the individual rasters actually orders 
# them by weeks starting at the beginning of the year all the way to the end
# (this might or might not be true in your case)
listOfRastersByWeek = arcpy.ListRasters()

# now this getting a little funky and might be too much
# how do you distribute the 52 weeks over the year, or how many weeks are in each month?
year = []
weekCounter = 0
for index in range(1,13):
    # now month contains an arrangement of calendar days in weekly chunks
    month = calendar.monthcalendar(1992, index)
    # this array tells us how many days in a week for the specified month
    daysinweek = [len[1 for day in weeks if day] for weeks in month]
    # if a weekly chuck has more than 4 days we'll keep it and count it as a full week
    weeksInMonth = sum(1 for daycount in daysinweek if daycount > 3)
    # since the assumption of 4 days is somewhat random, let's doublecheck the number
    # of weeks
    weekCounter = weekCounter + weeksInMonth
    if (weekCounter > 52):
        year.append(52 - (weekCounter - weeksInMonth))
    else:
        year.append(weeksInMonth)    

# at this point the year array has the number of weeks per month
# with the understanding of how many weeks per month let's compute the average

startIndex = 0
for index in range(12):
    # initialize the monthly raster with the first week
    monthRaster = arcpy.Raster(listOfRastersByWeek[startIndex])
    for rasterIndex in range(1,year[index]):
        if startIndex < 52:
            # add each week into the monthly raster
            monthRaster = monthRaster + arcpy.Raster(listOfRastersByWeek[startIndex])
            # push the index to the next week
            startIndex = startIndex + 1

    # average the added weekly values
    monthRaster = monthRaster / year[index]
    # save the raster with the name of the month
    monthRaster.save(month_dict[index + 1])