[GRASS-user] MODIS HDF conversion in TIF tool.

[Nikos Alexandris]

The one-line answer is here: http://trac.osgeo.org/gdal/ticket/1239

[  $ gdalwarp -of HFA \

    -s_srs '+proj=sinu +R=6371007.181 +nadgrids=@null +wktext' \ -t_srs
TARGET-SRS \ 'HDF4_EOS:EOS_GRID:"fpar8.hdf":MOD_Grid_MOD15A2:Fpar_1km' \
fpar8.hfa

        Note the key `+nadgrids=@null +wktext' options. In TARGET-SRS,
`+ellps=sphere' becomes unnecessary (unless, of course, one wants the
dataset to be warped into a projection on a sphere.)  ]

Using Markus script and replacing the necessary parameters as mentioned in
the end of the abovementioned link I get for example MOD09GA converted in
tif and reprojected in HGRS87 in one step.


Yooohoo!

Thank you for the script Markus!
---

The script ( ...slighlty modified, using 460x460 pixel resolution for the
output):

#!/bin/sh

#MODIS HDF/Sinusoidal to GeoTiff/HGRS87
#Based on:
#2005, Markus Neteler
# warp MODIS HDF to Erdas/Img - LatLong/WGS84

PROG="$0"

if [ $# -lt 1 -o "$1" = "-h" -o "$1" = "--help" -o "$1" = "-help" ] ; then
 echo "Script to warp all layers of MODIS HDF to Erdas/Img - LatLong/WGS84"
 echo "Usage:"
 echo "      $PROG M[OD|MYD]blablabla.hdf"
 exit 1
fi

FILE=$1
HDFlayerNAMES="`gdalinfo $FILE | grep SUBDATASET_ | grep _NAME | cut -d'='
-f2`"

#aarg! NDVI/EVI MOD13 contains spaces:
HDFlayerNAMES_NOSPACES="`gdalinfo $FILE | grep SUBDATASET_ | grep _NAME |
cut -d'=' -f2 | tr -s ' ' '|'`"

#MODIS does not use the "Normal Sphere (r=6370997)"!!!
# MODIS sphere, radius of 6371007.181 --> Instead use: '+proj=sinu
+R=6371007.181 +nadgrids=@null +wktext'
# http://edcdaac.usgs.gov/landdaac/tools/mrtswath/info/ReleaseNotes.pdf
#
#GGRS87 (or HGRS87 or EGSA87)
#<2100> +proj=tmerc +lat_0=0 +lon_0=24 +k=0.9996 +x_0=500000 +y_0=0
+ellps=GRS80 +towgs84=-199.87,74.79,246.62,0,0,0,0 +units=m +no_defs <>

#loop over all layers:
for i in ${HDFlayerNAMES_NOSPACES} ; do
        #transform back:
        INPUT="`echo ${i} | tr -s '|' ' '`"
	NEWNAME="`echo ${i} | cut -d':' -f4,5 | tr ':' '_' | tr -s '|' '_'`"
	# apply GCPs
	gdalwarp -s_srs '+proj=sinu +R=6371007.181 +nadgrids=@null +wktext' -t_srs
'+init=epsg:2100' -tr 460 460 "${INPUT}" $NEWNAME.HGRS87.tif
done





Nikos Alexandris wrote:
> 
> Well,
> 
> after searching and reading I found a script written by Markus Neteler. I
> tried to change the parameters in it and get directly from sinusoidal to
> HGRS87 (below the script only with the parameters changed). Initially it
> looked like it works... but when comparing with some vectors data
> (projected in HGRS87) a big shift appears.
> 
> After bothering Markus again about "how-to reproject MODIS surface
> reflectance data into the HGRS87 (EPSG:2100)" and following his
> recommendation to make it in 2 steps (1. sinusoidal to LatLong_WGS84 and
> 2. (re-)reprojecting into HGRS87) it worked.
> 
> I don't understand why can't it be done in one step.
> 
> The script:
> 
> #!/bin/sh
> 
> #Warp MODIS HDF to GeoTiff - HGRS87
> #Based on "modis_hdf2erdas_ll_wgs84.sh" script of Markus Neteler, 2005:
> warp MODIS HDF to Erdas/Img - LatLong/WGS84
> 
> PROG="$0"
> 
> if [ $# -lt 1 -o "$1" = "-h" -o "$1" = "--help" -o "$1" = "-help" ] ; then
>  echo "Script to warp all layers of MODIS HDF to GeoTiff - HGRS87"
>  echo "Usage:"
>  echo "      $PROG M[OD|MYD]blabla.hdf"
>  exit 1
> fi
> 
> FILE=$1
> HDFlayerNAMES="`gdalinfo $FILE | grep SUBDATASET_ | grep _NAME | cut -d'='
> -f2`"
> 
> #aarg! NDVI/EVI MOD13 contains spaces:
> HDFlayerNAMES_NOSPACES="`gdalinfo $FILE | grep SUBDATASET_ | grep _NAME |
> cut -d'=' -f2 | tr -s ' ' '|'`"
> 
> #MODIS does not use the "Normal Sphere (r=6370997)"!!!
> # MODIS sphere, radius of 6371007.181
> # http://edcdaac.usgs.gov/landdaac/tools/mrtswath/info/ReleaseNotes.pdf
> #
> #GGRS87 (or HGRS87 or EGSA87)
> #<2100> +proj=tmerc +lat_0=0 +lon_0=24 +k=0.9996 +x_0=500000 +y_0=0
> +ellps=GRS80 +towgs84=-199.87,74.79,246.62,0,0,0,0 +units=m +no_defs <>
> 
> 
> # cs2cs -le
> SPHERE="+ellps=GRS80"
> 
> #loop over all layers:
> for i in ${HDFlayerNAMES_NOSPACES} ; do
>         #transform back:
>         INPUT="`echo ${i} | tr -s '|' ' '`"
> 	NEWNAME="`echo ${i} | cut -d':' -f4,5 | tr ':' '_' | tr -s '|' '_'`"
> 	# apply GCPs
> 	gdalwarp -t_srs "+proj=tmerc +lat_0=0 +lon_0=24 +k=0.9996 +x_0=500000
> +y_0=0 +towgs84=-199.87,74.79,246.62,0,0,0,0 +units=m +no_defs $SPHERE"
> "${INPUT}" $NEWNAME.tmp.tif
> 	# fix the LatLong reference cut to fit into world map (ulx uly lrx lry):
> ????
> 	# -projwin -180 90 180 -90 ???
> 	# not desired for NDVI etc ???
> 	gdal_translate -of HFA -a_srs "EPSG:2100" \
>             $NEWNAME.tmp.tif $NEWNAME.HGRS87.tif
> 	rm -f $NEWNAME.tmp.tif
> done
> 

-- 
View this message in context: http://www.nabble.com/MODIS-HDF-conversion-in-TIF-tool.-tf4799564.html#a13900375
Sent from the Grass - Users mailing list archive at Nabble.com.