1

Goal- Building python script to compute monthly average from 8 days MODIS data

I have MODIS data which is contain 8 days product. Some of month has 3 images and some of month has four imgaes. So how i can build a Arc GIS python scrip to calculate monthly average from this data

I have listed below number of images in each month.

Jan - 4, Feb- 4, Mar- 4, Apr- 3, May- 3, June-4, July- 4, Aug - 4, Sept - 4, Oct - 3, Nov - 4, Dec - 4,

Total 45 images are in same foldar with 8 days interval name

e.g. MOD15A2.MRTWEB.A2005337.005.Fpar_1km, MOD15A2.MRTWEB.A2005345.005.Fpar_1km

After calculation i just i want to save it in month wise name.

I have Build one code but i can take only defined interval e.g 3 days or 4 days. not 3 days and 4 days at a same time.

Below i have shown my python code

import arcpy, os, sys  
from arcpy import env  
from arcpy.sa import *  
arcpy.env.overwriteOutput = True  
arcpy.CheckOutExtension("Spatial")  
arcpy.env.extent = "MAXOF"    
env.workspace = 'D:\MODIS-NDVI\FPAR-2005\MASKED-FPAR-05B'  
rasters = arcpy.ListRasters()  
out_ws = 'D:\MODIS-NDVI\FPAR-2005\AVG-FPAR-05'  

rastersCount = len(rasters)  
counter = 0  
dek = 1  
while counter < rastersCount:  
    dekad = []  
    for i in range(counter,(counter+2)):  
        print i  
        dekad.append(rasters[i])  
        outCellStatistics = CellStatistics(dekad, "MEAN", "NODATA")  
        outRasterName = out_ws + 'tmax_dek_{:03d}.tif'.format (dek)  
        #outRasterName = "00" + str(dek) + ".tif"  
        outCellStatistics.save(outRasterName)  
    counter += 2  
    dek += 1  

Above can take only defined interval e.g 3 days or 4 days. not 3 days and 4 days at a same time.

-1

I'm not familiar with the MODIS file naming conventions, but maybe you could try the regular expression module (re) to add rasters matching a pattern in your root folder, add them to "dekad" based on the pattern, then perform CellStatistics. Below, I've assumed the month is the 3 digits before the "Fpar_1km" (??) You could adjust it to suit your needs.

import arcpy, os, glob, re
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")

rootdir = 'D:/MODIS-NDVI/FPAR-2005/MASKED-FPAR-05B/'
env.workspace = rootdir
dekad = []
mons = ['005', '006']  # <-- adjust as necessary

for mon in mons:

    patt = r"MOD15A2\.MRTWEB\.A2005\d{3}\." + mon + r"\.Fpar_1km\.tif"
    prog = re.compile(patt)

    for raster in glob.glob(rootdir + "*.tif"):
        rastername = os.path.basename(raster)
        if prog.match(rastername):
            work_file = os.path.join(rootdir, rastername)
            print 'Work file: ' + work_file
            dekad.append(rastername)

    print dekad
    print "Creating a MEAN raster!"
    final_raster = os.path.join(rootdir, mon) + '_mean.tif'
    outfile = CellStatistics(dekad, "MEAN", "NODATA")
    outfile.save(final_raster)
    dekad = []    
  • Thank you Sir, which you assumed as month of images "3 digits before the "Fpar_1km" 005" its not month. it is common in every images. name format of images for different time is like A200500.005, A2005009.005, A2005017.005, A2005025.005, A2005033.005. name of different time period is stared after the year "2005" and its 8 days interval its increase . similarly i have the raster till A2005361 . First massage i have mentioned the full image name. Image name is changes 3 digit after the year (e.g 2005) – Shouvik Jha Jul 14 '16 at 5:21
1

It sounds like you just need to iterate over each month's folder(?), add each month's tif files to the list "dekad", perform the CellStatistics command and save result as the monthly average raster, correct?

Instead of the "for i in range... dekad.append[]..." method, you could do something like below using "glob"...

import arcpy, os, glob
from arcpy import env
from arcpy.sa import *
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")

rootdir = 'D:/MODIS-NDVI/FPAR-2005/MASKED-FPAR-05B'
dekad = []


for month in os.listdir(rootdir):
    if os.path.isdir(os.path.join(rootdir, item)):
        print month  
        work_dir = os.path.join(rootdir, month)
        final_raster = os.path.join(work_dir, month) + '_mean.tif'
        for raster in glob.glob(work_dir + "*.tif"):
            dekad.append(raster)
        print dekad
        print "Creating a MEAN raster!"
        outfile = CellStatistics(dekad, "MEAN", "NODATA")
        outfile.save(final_raster)
    dekad = []    
  • Thank you. I have call Raster in same folder only. Just for better understanding purpose i mentioned month wise number. But all images in same folder only. I have shared image details with folder – Shouvik Jha Jul 13 '16 at 19:19
  • You should use os.path.join() to concatenate your paths. – Paul Jul 13 '16 at 19:20
  • Thank you Sir. We have all data in same folder only. I just mentioned no of images in different month for better understanding purpose. Below i have mention my all input file names MOD15A2.MRTWEB.A2005001.005.Fpar_1km MOD15A2.MRTWEB.A2005009.005.Fpar_1km MOD15A2.MRTWEB.A2005017.005.Fpar_1km MOD15A2.MRTWEB.A2005025.005.Fpar_1km MOD15A2.MRTWEB.A2005033.005.Fpar_1km MOD15A2.MRTWEB.A2005041.005.Fpar_1km So on....... – Shouvik Jha Jul 13 '16 at 19:27
-1

I don't know of any other way to do this unless you could access the image's EXIF data. Test "datetime" or "datetimeoriginal" and check the month the image was taken. The Python image library (PIL) module can do this like below.

import arcpy, os, glob, re
from arcpy import env
from arcpy.sa import *
from PIL import Image

arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")

rootdir = 'D:/MODIS-NDVI/FPAR-2005/MASKED-FPAR-05B/'
env.workspace = rootdir
dekad = []
mons = ['05', '06']  # <-- adjust as necessary

def get_date_taken(path):
    # returns datetimeoriginal like "2014:02:18 14:25:24"
    return Image.open(path)._getexif()[36867]

for mon in mons:

    # match datetimeoriginal on month 
    patt = r"\d{4}\:" + mon + r"\:\d{2}\s\d{2}\:\d{2}\:\d{2}"
    prog = re.compile(patt)

    for raster in glob.glob(rootdir + "*.jpg"):

        image_date = get_date_taken(raster)
        print image_date
        rastername = os.path.basename(raster)
        print rastername
        normrast = os.path.normpath(raster)

        if prog.match(image_date):
            print 'it matched!'
            work_file = os.path.join(rootdir, rastername)
            dekad.append(rastername)

    print dekad
  • Glad to hear it. Regards - cm1 – cm1 Jul 16 '16 at 4:40
  • Happy to help, and welcome to Stack Overflow. If this answer or any other one solved your issue, please mark it as accepted. – cm1 Jul 16 '16 at 15:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.