CODEWITHSUNDEEP
rgisprojectionraster
modestprog
modestprog

Reputation: 73

Using R to change raster projection

I'm trying to change raster projection using R and the raster package. The input raster projection is Lambert Azimuthal; the parameters are here:

Coordinate System:
Lambert_Azimuthal_Equal_Area
False_Easting: 4321000,000000
False_Northing: 3210000,000000
Central_Meridian: 10,000000
Latitude_Of_Origin: 52,000000
GCS_ETRS_1989
Datum: D_ETRS_1989
Prime Meridian: 0


PROJCS
 ["ETRS_1989_LAEA",
   GEOGCS ["GCS_ETRS_1989",
           DATUM ["D_ETRS_1989",
                  SPHEROID ["GRS_1980",6378137.0,298.257222101]],
           PRIMEM["Greenwich",0.0],
          UNIT["Degree",0.0174532925199433]],
   PROJECTION["Lambert_Azimuthal_Equal_Area"],
   PARAMETER["False_Easting",4321000.0],
   PARAMETER["False_Northing",3210000.0],
   PARAMETER["Central_Meridian",10.0],
   PARAMETER["Latitude_Of_Origin",52.0],
   UNIT["Meter",1.0]]

I need to convert them to simple rasters in ESRI ASCII format, using longitude and latitude coordinates, with a Mercator-style projection, with a cell size of 0.1 degrees (I hope I explain myself well enough, because I don't have enough GIS skills, sorry). What I need are rasters in format .ASC where each value of the raster corresponds to a single cell of size N x N, where N is in degrees (e.g. 0.1 degrees), and the raster coordinates are in longitude/latitude.

I tried to use raster library in R, and follow examples found for the projectRaster function. But after many attempts using many several parameters, I wasn't be able to make it correctly. I think I'm not using the correct parameters for projection, datum, or something like this.

Here's what I tried. I load the raster in R, then I set its projection using:

>crs(r)<-"+proj=laea +lat_1=52 +lon_0=-10 +ellps=GRS80"

Then I define the output projection and I try the conversion and save:

>newproj <- "+proj=lonlat +lat_1=52 +lon_0=-10 +ellps=WGS84"
>pr2 <- projectRaster(r, crs=newproj, res=0.1)
>writeRaster(pr2, "newraster.asc", overwrite=TRUE)

No error messages, but resulting raster is not correctly projected (country borders don’t match, countries are slightly distorted).

Thanks for any help!

Upvotes: 2

Views: 20623

Answers (1)

Robert Hijmans
Robert Hijmans

Reputation: 47676

Given the description of the projection you provide, this seems wrong:

crs(r) <- "+proj=laea +lat_1=52 +lon_0=-10 +ellps=GRS80"

As you do not include the false northing and easting parameters; and lat_1 should probably be lat_0. This may be better:

crs(r) <- "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m"

This also looks odd:

newproj <- "+proj=lonlat +lat_1=52 +lon_0=-10 +ellps=WGS84"

How about

newproj <- "+proj=longlat +datum=WGS84"

Upvotes: 6

Related Questions