Generar archivos HDF4 con GDAL

El formato HDF es utilizado para la distribución de datos de sensores remotos, como por ejemplo, imágenes de satélite de varias misiones de la NASA. Aunque la mayoría de los datos de años recientes se distribuye en formato HDF5 o netCDF, todavía se pueden encontrar datos en la versión HDF4 que ya está cayendo en desuso. Un ejemplo son los datos de la misión MODIS al menos hasta la versión 6.1.

Dado que HDF4 es un formato viejo, no es fácil encontrar información al respecto sobre como manipularlos. Por ello cuando me encontré en la situación de generar este tipo de archivos buscar en la web no dio muy buenos resultados.

La situación era convertir un archivo shapefile a un archivo HDF, el shapefile contenía un campo con valores enteros que quería convertir a raster en formato HDF4.

GDAL es una librería frecuentemente utilizada en la manipulación de datos geoespaciales en formato raster. En mi caso utilicé dos herramientas, primero gdal_rasterize para convertir el shapefile a GeoTIFF, y después en un segundo paso con gdal_translate convertirlo finalmente a HDF4.

Convertir de SHP a GeoTIFF

El comando que se utiliza es el siguiente (cambiar las variables que empiezan con $ por tus valores):

gdal_rasterize -l $Tile -a $SiteNum -tr 30.0 30.0 -a_nodata 0.0 -te $West $South $East $North -ot Byte -of GTiff $FileName $FileOut

Las opciones que utilicé son:

• -l: para especificar el nombre del layer
• -a: el nombre del campo de atributo o field, estos son los valores que van a contener los pixeles del raster resultante
• -tr <xres> <yres>: para especificar el tamaño del pixel a generar, en este caso como se desea una resolución compatible con Landsat, se quiere un tamaño de pixel de 30×30 metros
• -a_nodata <valor>: especifica el valor para NoData, en este caso 0.0
• -te <xmin> <ymin> <xmax> <ymax>: la extensión, es decir las coordenadas extremas, del raster en orden Oeste, Sur, Este, y Norte; o bien punto inferior izquierdo (x, y) y punto superior derecho (x, y)
• -ot: especifica el tipo de datos, en este caso enteros tipo byte
• -of: el formato de archivo de salida, GeoTIFF
• Al final se colocan los nombres de archivo con ruta completa o relativa, primero el archivo de entrada (el shapefile) y después el de salida (el raster).

Dado que yo utilicé un script de shell en GNU/Linux para ejecutar el comando utilicé variables con el símbolo $ las cuales se tienen que definir previamente. Para ejecutar el comando directamente se omite esto y colocar directamente las opciones.

Convertir de GeoTIFF a HDF4

Finalmente a lo que venimos. Para convertir de GeoTIFF a HDF4 se utiliza el comando gdal_translate, sin embargo, al usar las opciones por default se generará un archivo de tres dimensiones, lo cual es ridículo para una imagen que es 2D. Esto además creará un archivo raster que no será visualizado correctamente en en QGIS or ArcGIS.

Para forzar a gdal_translate a utilizar solo dos dimensiones se pueden configurar las opciones, sin embargo, la forma de hacerlo se encuentra muy escondida en la documentación (https://gdal.org/en/stable/drivers/raster/hdf4.html, de nada 😉 )

gdal_translate -co rank=2 -of HDF4Image $FileIn $FileOutHDF

• -co: esta opción se utiliza para especificar opciones durante la creación del archivo. En este caso rank=2 especifica que se cree una archivo 2D en lugar de un 3D que por alguna razón es la opción por defecto.
• -of: el formato especificado es HDF4
• Finalmente los archivos de entrada y de salida.

Acerca de la proyección

Los dos comandos de GDAL que utilizamos no manipulan la proyección o el sistema de coordenadas por ello es necesario que las coordenadas extremas del raster se encuentren en el mismo sistema de referencia geoespacial que el archivo shapefile de entrada, así como también asegurarse que el tamaño de pixel sea adecuado.

Espero que la información se a de utilidad para todos aquellos que tengan un problema similar.