[Grass-de] Import von NOAA Daten nach GRASS

Do Aug 20 23:17:19 CEST 2009

Hallo Thomas,

2009/8/20 Becker, Thomas <thob at dmu.dk>:
> Hallo,
>
> ich habe hier einige Daten von NOAA-15 die ich gern in eine GRASS Location importieren würde. Leider bekomme ich das auch über Umwege nicht hin. Die Daten liegen im HDF-EOS Format vor und wurden von dieser Seite heruntergeladen (ftp://ghrc.nsstc.nasa.gov/pub2/data/amsu-a/noaa-15/2000/). Der HDF Explorer zeigt mir die Daten auch wunderbar an und so wie sie dort dargestellt werden bekomme ich sie auch mittels 'gdal_translate' in ein GTiff, aber es fehlt die räumliche Referenz.

Die Dateien sind in HDF4:

file amsua15_2000.104_09964_0701_0856_WI.eos
amsua15_2000.104_09964_0701_0856_WI.eos: Hierarchical Data Format
(version 4) data

Mit etwas Muehe habe ich mir HDF4.2 kompiliert und sehe nun:

gdalinfo amsua15_2000.104_09964_0701_0856_WI.eos
Driver: HDF4/Hierarchical Data Format Release 4
Files: amsua15_2000.104_09964_0701_0856_WI.eos
Size is 512, 512
Coordinate System is `'
Metadata:
  HDFEOSVersion=HDFEOS_V2.4
Subdatasets:
  SUBDATASET_1_NAME=HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:23800.37 MHz
  SUBDATASET_1_DESC=[857x30] 23800.37 MHz Orbit 9964 (32-bit
floating-point)
  SUBDATASET_2_NAME=HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:31400.42 MHz
  SUBDATASET_2_DESC=[857x30] 31400.42 MHz Orbit 9964 (32-bit
floating-point)
  SUBDATASET_3_NAME=HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:50299.91 MHz
...
  SUBDATASET_15_NAME=HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:88997.00 MHz
  SUBDATASET_15_DESC=[857x30] 88997.00 MHz Orbit 9964 (32-bit floating-point)
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0,  512.0)
Upper Right (  512.0,    0.0)
Lower Right (  512.0,  512.0)
Center      (  256.0,  256.0)

Einen einzelne Karte kann man so finden:
gdalinfo 'HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:57290.33 +- 322.2 +- 4.5 MHz'
Driver: HDF4Image/HDF4 Dataset
Files: amsua15_2000.104_09964_0701_0856_WI.eos
Size is 30, 857
Coordinate System is `'
Metadata:
  HDFEOSVersion=HDFEOS_V2.4
  _FillValue=-999
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0,  857.0)
Upper Right (   30.0,    0.0)
Lower Right (   30.0,  857.0)
Center      (   15.0,  428.5)
Band 1 Block=30x1 Type=Float32, ColorInterp=Gray
  NoData Value=-999

Aber, wie man sieht, nur Pixelkoordinaten :(

Die Doku sagt
ftp://ghrc.nsstc.nasa.gov/pub2/data/amsu-a/doc/amsu-a_dataset.html
-> " Geolocation fields"

Aber mir scheint, das GDAL sie nicht sieht (denn ansonsten wuerden sie
hier angewendet):

gdalwarp 'HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:57290.33 +- 322.2 +- 4.5 MHz' out.tif
ERROR 1: Unable to compute a transformation between pixel/line
and georeferenced coordinates for
HDF4_EOS:EOS_SWATH:"amsua15_2000.104_09964_0701_0856_WI.eos":Orbit
9964:57290.33 +- 322.2 +- 4.5 MHz.
There is no affine transformation and no GCPs.

> Ich weiß nicht was ich noch machen soll. Um eine räumliche Referenz für 'gdal_translate' zur verfügung zu stellen habe ich mir ein Python-Script geschrieben, welches für jedes 15. Pixel die Koordinaten aus den Latitude/Longitude Tabellen der HDF Daten ausliesst und in ein 'optfile' schreibt. Der Beginn des 'optfile' sieht dann so aus:
>
> -of GTiff -gcp 1 1 -91.868202 -6.065200 -gcp 15 1 -83.118896 -4.622600 ...
>
> Auf grund der Koordinaten geh ich davon aus, dass die Daten in WGS84 vorliegen, was ja EPSG:3395 entsprechen sollte.

http://www.spatialreference.org/ref/epsg/3395/
->     PROJECTION["Mercator_1SP"],

Es ist also nicht LatLong EPSG 4326.

> Selbst wenn ich das 'optfile' um '-a_srs EPSG:3395' erweitere bekomme ich die Daten nicht in die von mir gewünschte Referenz. Wenn ich das resultierende GTiff nach Grass importiere dann werden die Daten weiterhin in Form des Scan stripes dargestellt und die linke untere Ecke liegt bei 0 0, die rechte untere Ecke bei 30 0 usw. Mit anderen Worten, die Sensordaten sind nicht mit den Koordinaten verbunden, denn die Pixel werden einfach nur entsprechend ihrer Position im Scan stripe angeordnet. Der Scan stripe ist immer genau 30 Spalten breit.

Ich glaube, dass die '-a_srs EPSG:3395' Geschichte erst geht, wenn die Daten
als richtige Matrix gelesen werden.

> Der Versuch die Daten mittels HEG (http://newsroom.gsfc.nasa.gov/sdptoolkit/HEG/HEGHome.html) in GTiff zu konvertieren scheitert, da die Software nicht mit Gruppennamen mit Leerzeichen umgehen kann. Ich habe mit den Leuten Kontakt aufgenommen und man wird versuchen das zu ändern, aber das wird mir rein zeitlich gesehen für mein Projekt nicht mehr helfen.
>
> Kann mir einer von Euch in der Geschichte weiterhelfen?

Ich habs versucht, leider ohne Erfolg. Vorschlag: An gdal-dev schreiben...
Vielleicht waeren sie bei einem Hinweis auf
ftp://ghrc.nsstc.nasa.gov/pub/doc/amsu_a/amsuareader.c
bereit, einen Treiber zu schreiben?

Viele Gruesse aus Trento,
Markus