Extract  files from  MODIS HDF files: Processing in a loop

4664
1
Jump to solution
03-01-2013 01:48 PM
ZiaAhmed
New Contributor III
I am trying to extract a file (0=NDVI raster)  from  a MODIS-hdf file and save it as geo-tiff file and project it. I did this in ArcGIS/ModelBuilder with a limted number of files. Below is my phython code for this model. I want to run this code in a loop- that can help to  process (extraction, projection and save in a different directory) all HDF file in the working directory.  I have no experience in python  doing something in a loop. If someone help me to do this in  a loop, that will save  lot of my time.
Thanks
Zia


#-------------------------------------------------------------------------
# Extract_HDF_TO_TIF_PROJ.py
#
# Created on: 2013-02-28 11:40:59.00000
#   (generated by ArcGIS/ModelBuilder)
# Description:
# ---------------------------------------------------------------------------

# Import arcpy module
import arcpy


# Local variables:
MOD13Q1_A2012145_h09v04_005_2012166114230_hdf = "L:\\MODIS_DATA\\MODIS_RAW\\May_24\\MOD13Q1.A2012145.h09v04.005.2012166114230.hdf"
MOD13Q1_A2012145_h09v05_005_2012166111623_hdf = "L:\\MODIS_DATA\\MODIS_RAW\\May_24\\MOD13Q1.A2012145.h09v05.005.2012166111623.hdf"
MOD13Q1_A2012145_h10v04_005_2012166113110_hdf = "L:\\MODIS_DATA\\MODIS_RAW\\May_24\\MOD13Q1.A2012145.h10v04.005.2012166113110.hdf"
MOD13Q1_A2012145_h10v05_005_2012166112540_hdf = "L:\\MODIS_DATA\\MODIS_RAW\\May_24\\MOD13Q1.A2012145.h10v05.005.2012166112540.hdf"
MOD13Q1_A2012145_h11v04_005_2012166113905_hdf = "L:\\MODIS_DATA\\MODIS_RAW\\May_24\\MOD13Q1.A2012145.h11v04.005.2012166113905.hdf"
# Output  file and directory
ma24s02_tif = "L:\\MODIS_DATA\\MODIS_TIF\\May_24\\ma24s02.tif"
ma24s02_tif__2_ = "L:\\MODIS_DATA\\MODIS_TIF_PROJECTED\\May_24\\ma24s02.tif"
ma24s03_tif = "L:\\MODIS_DATA\\MODIS_TIF\\May_24\\ma24s03.tif"
ma24s03_tif__2_ = "L:\\MODIS_DATA\\MODIS_TIF_PROJECTED\\May_24\\ma24s03.tif"
ma24s04_tif = "L:\\MODIS_DATA\\MODIS_TIF\\May_24\\ma24s04.tif"
ma24s04_tif__2_ = "L:\\MODIS_DATA\\MODIS_TIF_PROJECTED\\May_24\\ma24s04.tif"
ma24s05_tif = "L:\\MODIS_DATA\\MODIS_TIF\\May_24\\ma24s05.tif"
ma24s05_tif__2_ = "L:\\MODIS_DATA\\MODIS_TIF_PROJECTED\\May_24\\ma24s05.tif"
ma24s06_tif = "L:\\MODIS_DATA\\MODIS_TIF\\May_24\\ma24s06.tif"
ma24s06_tif__2_ = "L:\\MODIS_DATA\\MODIS_TIF_PROJECTED\\May_24\\ma24s06.tif"


# Process: Extract Subdataset
arcpy.ExtractSubDataset_management(MOD13Q1_A2012145_h09v04_005_2012166114230_hdf, ma24s02_tif, "0")
# Process: Project Raster
arcpy.ProjectRaster_management(ma24s02_tif, ma24s02_tif__2_, "PROJCS['NAD_1983_Albers',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',23.0],UNIT['Meter',1.0]]", "NEAREST", "", "", "", "")

# Process: Extract Subdataset (2)
arcpy.ExtractSubDataset_management(MOD13Q1_A2012145_h09v05_005_2012166111623_hdf, ma24s03_tif, "0")

# Process: Project Raster (2)
arcpy.ProjectRaster_management(ma24s03_tif, ma24s03_tif__2_, "PROJCS['NAD_1983_Albers',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',23.0],UNIT['Meter',1.0]]", "NEAREST", "", "", "", "")

# Process: Extract Subdataset (3)
arcpy.ExtractSubDataset_management(MOD13Q1_A2012145_h10v04_005_2012166113110_hdf, ma24s04_tif, "0")

# Process: Project Raster (3)
arcpy.ProjectRaster_management(ma24s04_tif, ma24s04_tif__2_, "PROJCS['NAD_1983_Albers',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',23.0],UNIT['Meter',1.0]]", "NEAREST", "", "", "", "")

# Process: Extract Subdataset (4)
arcpy.ExtractSubDataset_management(MOD13Q1_A2012145_h10v05_005_2012166112540_hdf, ma24s05_tif, "0")

# Process: Project Raster (4)
arcpy.ProjectRaster_management(ma24s05_tif, ma24s05_tif__2_, "PROJCS['NAD_1983_Albers',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',23.0],UNIT['Meter',1.0]]", "NEAREST", "", "", "", "")

# Process: Extract Subdataset (5)
arcpy.ExtractSubDataset_management(MOD13Q1_A2012145_h11v04_005_2012166113905_hdf, ma24s06_tif, "0")

# Process: Project Raster (5)
arcpy.ProjectRaster_management(ma24s06_tif, ma24s06_tif__2_, "PROJCS['NAD_1983_Albers',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',23.0],UNIT['Meter',1.0]]", "NEAREST", "", "", "", "")
Tags (2)
1 Solution

Accepted Solutions
ZiaAhmed
New Contributor III
I solved the with following codes.
Zia

try:
    import arcpy
    arcpy.env.workspace = r"L:/Test_input"
    rasterList = arcpy.ListRasters()

    ##Extract NDVI from HDF
    for raster in rasterList:
        tifOut="L:\\Test_out\\" + raster[9:-22] # SET OUTPUT DIRECTORY WITH FILE NAMES("2000193.h08v05")
        arcpy.ExtractSubDataset_management(raster, tifOut + ".tif", "0")

except:
    print "Extract Subdataset example failed."
    print arcpy.GetMessages()

View solution in original post

0 Kudos
1 Reply
ZiaAhmed
New Contributor III
I solved the with following codes.
Zia

try:
    import arcpy
    arcpy.env.workspace = r"L:/Test_input"
    rasterList = arcpy.ListRasters()

    ##Extract NDVI from HDF
    for raster in rasterList:
        tifOut="L:\\Test_out\\" + raster[9:-22] # SET OUTPUT DIRECTORY WITH FILE NAMES("2000193.h08v05")
        arcpy.ExtractSubDataset_management(raster, tifOut + ".tif", "0")

except:
    print "Extract Subdataset example failed."
    print arcpy.GetMessages()
0 Kudos