How do I change the projection of raster data in MATLAB?
65 views (last 30 days)
Is it possible to change the projection of georeferenced raster data in Matlab? I am interested in both reprojecting the my spatial data (geotiffs) between different projected coordinate systems and "unprojecting" it to geographic coordinates. Basically, I'm looking to do the equivalent "project raster" in ArcGIS. I can project and un-project shapefiles without issue using minvtran and mfwdtran, but I can't find an equivalent function for rasters.
Chad Greene on 28 Feb 2016
The mfwdtran and minvtran functions will work. If you have some x and y vectors associated with your geotiff, grid the coordinates like this:
[X,Y] = meshgrid(x,y);
then use minvtran to get geocoordinates of each pixel by
[lat,lon] = minvtran(X,Y);
Then you'll be able to use pcolorm(lat,lon,Z) or transform lat and lon into a different projection using mfwdtran.
Matthew Cooper on 17 May 2021
I recommend using gdal, which can be called from Matlab using the system command or the shell escape character. Here is an example that pretends you want to reproject your data from North America Lambert Conformal Conic projection to Polar Stereographic North projection:
ras_in = 'MyRaster_Lambert.tif';
ras_out = 'MyRaster_Stereo.tif';
func = '/path/to/gdalwarp ';
opts = '-s_srs ESRI:102009 -t_srs EPSG:3413 -r bilinear -of GTiff -co "TFW=YES" ';
command = [func opts ras_in ' ' ras_out];
status = system(command);
You will have to replace the path to gdalwarp with the path where it exists on your system. Mine is /usr/local/bin/gdalwarp. You can customize this by defining path variables or adding the path to your matlab search path. In the 'opts' variable, I included some common gdalwarp options that might be helpful, but check out the gdal documention to understand what is available.
The same command using the shell escape character would look like this:
!/path/to/gdalwarp -s_srs ESRI:102009 -t_srs EPSG:3413 -r bilinear -of GTiff -co "TFW=YES" MyRaster_Lambert.tif MyRaster_Stereo.tif
You can run the command above directly in the matlab editor/terminal.
Performing the reprojection with Matlab is possible but quite difficult and involves careful thought, way beyond the scope of this answer. In its simplest form, you need to interpolate the original data onto the new coordinates, but you need to consider what method of interpolation is appropriate and how the projection affects the mapping between the original data and the 'projected' data.
One last point. Matlab will not allow you to change the R.RasterInterpretation property between 'Cells' and 'Postings'. This has important implications for all of the other functions in the mapping toolbox. Make sure you understand which interpretation is correct. If you need to modify the interpretation before reading the data into matlab, you can use something like this:
!/path/to/gdal_translate -mo AREA_OR_POINT=AREA MyRaster_Stereo.tif MyRaster_Stereo_asCells.tif