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.

Improve this question

3 Answers 3

Reset to default
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 = []
Improve this answer
3
  • 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 Commented Jul 13, 2016 at 19:19
  • You should use os.path.join() to concatenate your paths. Commented Jul 13, 2016 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....... Commented Jul 13, 2016 at 19:27
-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 = []
Improve this answer
1
  • 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) Commented Jul 14, 2016 at 5:21
-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
Improve this answer
2
  • Glad to hear it. Regards - cm1 Commented Jul 16, 2016 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. Commented Jul 16, 2016 at 15:02

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.