I will have to study that (!) carefully. Thank you for your time. BTW: I tried to reproject MODIS surface reflectance data with Markus script first from sinusoidal to WGS84 and then to HGRS87. And it works... ! But I can't get it in one line command successfuly.
I posted also in gdal (as I unterstand -with my limited experiences- it is a gdal-issue). Link to gdal-dev: http://www.nabble.com/warping-MODIS-data-with-gdal-tf2688531.html Ivan Shmakov wrote: > >>>>>> Nikos Alexandris <[EMAIL PROTECTED]> writes: > > [...] > > >>> you can use gdal_translate to add (per scanline) GCPs to the metadata > >>> before running gdalwarp. gdal_translate -gcp pixel line easting > >>> northing [elevation] > > >> Indeed, it works, at least for the AIRS data. (Although I've > >> had to upgrade GDAL 1.3.2, as supplied with Debian 4.0 r1, to > >> 1.4.3.) > > >> For those who are interested, the script I've used looks roughly > >> as follows: > > > You think your script would be applicable to Surface Reflectance > > products (L2) ? > > Do you mean MODIS L2 data? Well, I'll assume it for now... > > > Sorry for bothering again, but since I've limited experience with > > scripts... I ask before I I try to modify it for my need. > > I'll try to comment on some of the points where the script is > tied closely to the AIRS L2 data format. > > >> --cut-- > >> #!/bin/bash > >> file=AIRS.2007.07.09.079.L2.RetStd.v4.0.9.102.D07190164053.hdf > >> file_gcp=airs-2007-07-09-gran-79-l2-retstd-tsurfstd-gcp+0.5.tiff > >> > file_out=airs-2007-07-09-gran-79-l2-retstd-tsurfstd-gcp+0.5-laea-10km.tiff > > >> ## FIXME: specify ``no projection'' instead? > >> source_srs='+proj=latlong' > > First of all, this seems to specify the projection system for > the GCPs. Thus, it's essential. No FIXME: needed. > > >> target_srs='+proj=laea +lat_0=55.0 +lon_0=90' > > >> ## FIXME: more clean error handling > >> set -x -e > > >> long_lat_to_gcp () { > >> gawk '$1 != -9999 && $2 != -9999 { > >> print "-gcp", (NR - 1) % 30 + .5, int ((NR - 1) / 30) + > .5, $1, $2; > > Here, the `-gcp MINOR MAJOR LONGITUDE LATITUDE' lines are made. > It's assumed that the data is stored as an 30x45 (30 pixels > along the minor dimension) array, i. e.: > > a [1, 1] ... a [1, 30] a [2, 1] ... a [2, 30] ... a [45, 30] > > These are to be tailored for the data to be processed (e. g., > 1354 x (10 * scans) for 1 km MODIS L2 swath.) > > >> }' > >> } > > >> gdal_translate \ > >> $(paste <(hdfdump --text "$file" {Long,Lat}itude) | > long_lat_to_gcp) \ > > Extraction of the geolocation arrays is done with `hdfdump' [1], > yet another ``HDF to plain-something'' conversion tool. It > could probably be done with `hdp' or GRASS (through the ``x, y'' > location) as well. > > Now, the `-gcp ...' lines are put into the command line for > `gdal_translate'. Command line arguments arrays are subject to > OS- (and configuration-)dependent limits, and it's assumed that > there's enough room for that many options. > > This assumption could easily become false for MODIS data on most > OS, with some 1000000 points (or more.) Either some points are > to be discarded, or the GCP setting is to be applied in multiple > passes. I haven't tried doing either. > > >> > "HDF4_EOS:EOS_SWATH:\"$file\":L2_Standard_atmospheric&surface_product:TSurfStd" > \ > > This line have to be tailored according to the way GDAL ``sees'' > the file. One could use `gdalinfo' on the file in question to > find the pattern one needs. > > >> "$file_gcp" > >> gdalwarp \ > >> -tps -s_srs "$source_srs" -t_srs "$target_srs" \ > >> -tr 1e4 1e4 \ > > Since MODIS data has finer resolution than AIRS, the destination > resolution should be changed as well. E. g., `-tr 1e3 1e3' > could probably be used for 1 km MODIS L2 data. > > >> "$file_gcp" "$file_out" > >> --cut-- > > >> Is there any chance that GDAL will obtain ``GCPs'' from the > >> HDF-EOS file directly? > > Of course, this would eliminate all the aforementioned issues, > and the script itself. > > >> [...] > > >>>> The one more problem with the MODIS L2 data is that it consists of > >>>> overlapping scans -- MODIS could ``see'' some parts of the Earth in > >>>> consequent scans. > > This is one more question I've had to check. May be the > conversion has to be done scan-wise -- each scan is to be > selected and georeferenced (with `gdal_translate' and, e. g., > the above script), then reprojected with `gdalwarp'. The > resulting reprojected strips are to be somehow merged after. > > [...] > > [1] http://theory.asu.ru/~ivan/devel/hdfdump/ > > PS. I feel that the gdal-dev mailing list is a better place to discuss > this kind of problems. Could we move there? > > _______________________________________________ > grass-user mailing list > [email protected] > http://lists.osgeo.org/mailman/listinfo/grass-user > > -- View this message in context: http://www.nabble.com/MODIS-HDF-conversion-in-TIF-tool.-tf4799564.html#a13877276 Sent from the Grass - Users mailing list archive at Nabble.com. _______________________________________________ grass-user mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/grass-user
