I have satellite swath data in the form of an HDF5 file containing 3 arrays: 'data_array', 'latitude' and 'longitude'. All three are 32-bit floats, are arrays of the same dimensions (X-by-Y) and the names are self-explanatory. All three datasets use a fill value of -32767

What I want to do it reproject data_array from native satellite swath onto a latitude/longitude grid with a pixel size of 0.009 degrees. I want to save the resulting data array as a geoTIFF containing project information (for use in QGis/ENVI/etc).

I have tried using gdal based upon this but gdalwarp segfaults (maybe because there are some fill value pixels in the lat/lon grid? Unsure how to deal with this).

I also tried using python and the pytroll libraries, but find the documentation somewhat impenetrable. All my attempts thus far have failed.

So, can anyone recommend a simple method to perform this reprojection?

I assumed it would be something quite common, but am having difficulty finding a solution.

I can plot project the data onto a grid using NASA Panoply, but this appears to be for visualisation only. I can't see how to save a GEOTiff output.

  • 1
    I'm trying to understand your question and something is making me wonder : if you have lat/long value arrays... why would you want to reproject the data ? From what I understand, the data is already in lat long. no ? To help us understand, maybe you could try and upload some data, some extract of data or at least some header file content ? – gisnside Sep 3 '17 at 10:02
  • it would help to know what data you are using. Is it MODIS data ? other NASA data ? – radouxju Sep 7 '17 at 9:57
  • @gisnside The data is in raw instrument swath projection, which is unsuited for merging with other data sources. Hence why I wish to reproject onto a lat/lon grid. – os1 Sep 8 '17 at 11:30
  • @AndreJ, yes - i saw both of those topics before posting my own question. I can't really figure out how they are of use, though. – os1 Sep 8 '17 at 11:31

You might be interested by the Orfeo Toolbox which might have some useful feature for you (i only know about the tool, not the detailed procedure).

Orfeo ToolBox (OTB) is an open-source project for state-of-the-art remote sensing. Built on the shoulders of the open-source geospatial community, it can process high resolution optical, multispectral and radar images at the terabyte scale. A wide variety of applications are available: from ortho-rectification or pansharpening, all the way to classification, SAR processing, and much more!

All of OTB’s algorithms are accessible from Monteverdi, QGIS, Python, the command line or C++. Monteverdi is an easy to use visualization tool with an emphasis on hardware accelerated rendering for high resolution imagery (optical and SAR). With it, end-users can visualize huge raw imagery products and access all of the applications in the toolbox. From resource limited laptops to high performance MPI clusters, OTB is available on Windows, Linux and Mac. It is community driven, extensible and heavily documented.


You say in one of your comment : "the data is in raw instrument swath projection". You might specifically be interested by the Orfeo Toolbox recipes structured in articles and with code examples :

From raw image to calibrated product

Pre-processing tasks

This section presents various pre-processing tasks that are presented in a classical order to obtain a calibrated, pan-sharpened image.

  • Optical radiometric calibration

  • Pan-sharpening

  • Digital Elevation Model management

  • Ortho-rectification and map projections

In the same area, one more page should help you out :

Residual registration

Image registration is a fundamental problem in image processing. The aim is to align two or more images of the same scene often taken at different times, from different viewpoints, or by different sensors. It is a basic step for orthorectification, image stitching, image fusion, change detection, and others. [...] Sensor model is generally not sufficient to provide image registrations. Indeed, several sources of geometric distortion can be contained in optical remote sensing images including earth rotation, platform movement, non linearity, etc.

They result in geometric errors on scene level, image level and pixel level. It is critical to rectify the errors before a thematic map is generated, especially when the remote sensing data need to be integrated together with other GIS data.

| improve this answer | |
  • Thanks. I will investigate this and see if it is helpful, if so I'll accept your answer. – os1 Sep 8 '17 at 13:40
  • No problem. I regret it's not more detailed but I have limited knowledge on image preprocessing like this. I have more knowledge on the image classification side, when the image has been pre-processed ;) – gisnside Sep 8 '17 at 14:49

Try this with GDAL:

gdal_translate hdf5_file -sds bands.tif

Or this in R:

system(paste0('gdal_translate ', get_subdatasets(file_path)[band],' myfile_', band','.tif'))
| improve this answer | |

You could try georeferencing through thinplate spline (tps).

In Python should it be something like this code line:

import gdal
img2=gdal.Open(img1.GetSubDatasets()[3][0]) #where 3 is the band of interest (is an example)
gdal.Warp(FOLD+'/FILE.tif',img2,tps=True) #the -tps algorithm georreferences your data
| improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.