Question

Forgive me for me being naive. I am having trouble with re-projecting SGDF.

I have an xyz (x=longitude, y=latitude, z=value) 3-column dataset. Longitude and Latitude values are coordinates in EPSG:3035 format. I converted the data frame into a grid with a spatial resolution of 5 km * 5km projected in EPSG:3035. I wish to re-project the grid to EPSG:4326 with a spatial resolution of 0.05*0.05. However i get the following warning messages:

Warning messages:
1: In spTransform(radon, CRS("+init=epsg:4326")) :   Grid warping not available, coercing to points
2: In spTransform(as(x, "SpatialPixelsDataFrame"), CRSobj, ...) :   Grid warping not available, coercing to points

Can anybody please tell me, how i can re-project the grid. Below is a small reproducible example:

library(sp)
library(rgdal)
library(raster)

x=c(5013500, 5018500, 4883500, 4888500, 4893500, 4898500, 4908500,4948500, 4953500, 4958500, 4963500, 4973500, 4978500, 4988500, 5008500, 5013500, 5028500, 4878500, 4883500, 4888500, 4893500,4898500, 4903500,4928500, 4963500, 4968500, 4973500, 4978500, 4983500, 4988500)

y=c(5395500, 5395500, 5390500, 5390500, 5390500, 5390500, 5390500,5390500, 5390500, 5390500, 5390500, 5390500, 5390500, 5390500, 5390500,5390500, 5390500, 5385500, 5385500, 5385500, 5385500, 5385500, 5385500,5385500, 5385500, 5385500, 5385500, 5385500, 5385500, 5385500)

z=c(1.74, 1.74, 1.82, 1.82, 1.82, 1.81, 1.81, 1.78, 1.77, 1.77, 1.76,1.76, 1.75, 1.74, 1.73, 1.73, 1.72, 1.82, 1.82, 1.81, 1.81, 1.80, 1.80, 1.78, 1.75, 1.75, 1.74, 1.74, 1.73, 1.73)

df1=data.frame(x,y,z)
coordinates(df1) <- c("x", "y")
proj4string(df1)=CRS("+init=epsg:3035")
gridded(df1)=TRUE
fullgrid(df1)=TRUE

getGridTopology(df1)

                     x       y
cellcentre.offset 4878500 5385500
cellsize             5000    5000
cells.dim              31       3

newgrid = spTransform(df1, CRS("+init=epsg:4326"))

Well, thanks to Paul, i was able to re-project the grid using gdalwarp. But still, the spatial resolution is different:

    Coordinates:
    min      max
x -31.34409 50.36650
y  34.07928 71.87206
Is projected: FALSE 
proj4string :
[+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0]
Grid attributes:
       cellcentre.offset   cellsize cells.dim
x         -31.31151 0.06515996      1254
y          34.11186 0.06515996       580
Data attributes:
Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
-0.11290  0.00000  0.00000  0.06136  0.00000  1.90300 

Any thoughts??

Was it helpful?

Solution

You are getting this message because after reprojection the points that where on a regular grid in lat-lon are no longer a regular grid in another projected system. The reprojection step does not preserve the regularness because for a lot of projections the change between lat-lon and a projected system is not constant. So the distortion, or change, is larger in one area than in another, leading to a non-regular grid. Imagine the squares of the grid in lat-long becoming a non-square shape in a projected system.

The solution that is most simple is just to project the points to EPSG:4326 first, and then perform the interpolation step. I'm not sure which interpolation routine you use, but e.g. gstat does not support interpolation in lat-lon. So, going to a projected system, not lat-lon, before interpolating the data is always a safe bet.

If you really need to reproject a grid, some kind of interpolation needs to take place. You do two things: using a tool like gdalwarp that can perform the coordinate transform and the subsequent warping of the resulting grid. Alternatively, you can interpolate the non-regular grid to a regular grid using e.g. nearest neighbor or any interpolation routine you like.

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top