4

I am trying to reproject a raster from Cylindrical Equal-Area (Lambert). Central Meridian: -160. Datum: WGS 1984 to latitude and longitude But I'm getting multiple error messages. The original raster looks like this when plotted with rasterVis: enter image description here

My code is as follows:

library(raster)
#read raster
raster(r)

# Define new Proj.4 spatial reference 
sr <- "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"

# Project Raster
projected_raster <- projectRaster(r, crs = sr)

#But I get the following errors: 
Error in if (maxy == miny) { : missing value where TRUE/FALSE needed
In addition: Warning message:
In rgdal::rawTransform(projfrom, projto, nrow(xy), xy[, 1], xy[,  :
53 projected point(s) not finite

I am not sure whether this is a problem with the raster or my code? Please help!

5
  • Try CRS in your code: sr <- CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0")
    – Guz
    Dec 7, 2016 at 14:43
  • Hi Guzman, I tried this and it made no difference. Dec 7, 2016 at 15:07
  • Use gdaltransform() from gdalUtils package
    – aldo_tapia
    Dec 7, 2016 at 15:15
  • @PhoebeStewart-Sinclair please try to provide a reproducible example, so we can understand better your question and give a better answer!
    – Guz
    Dec 7, 2016 at 15:40
  • The website this data comes from wri.org/resources/data-sets/reefs-risk-revisited under Global Threats:Data. I want to reproject the basic w001001.adf raster files from Equal-Area (Lambert) to lat/long ("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0") in R Dec 7, 2016 at 16:00

1 Answer 1

1

Error probably came from its extent.

library(raster)
r <- raster("C:/GIS/Global_Threats/Acidification/arag_500/w001001.adf")
crs(r) <- "+proj=cea +datum=WGS84 +lon_0=-160.0 +lat_ts=0.0 +x_0=0.0 +y_0=0.0"
r

r returns extent : -20037507, 19962493, -6363885, 6436115 (xmin, xmax, ymin, ymax) so just as a test I have cut y-value (latitude) at 90% to north and south.

extent(r) <- c(xmin= -20037507, xmax= 19962493, ymin= 0.9*(-6363885), ymax= 0.9*(6436115))
sr <- "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0" 
projected_raster <- projectRaster(r, crs= sr, method = 'bilinear')
plot(projected_raster)

enter image description here

How much you can cut northern/southern edge will depend on your study area.

[EDIT] This discussion will be better and complete answer to your question about what happened. Re-projecting raster in R: gives warning that projected point(s) not finite?

1
  • I'm glad it worked out.
    – Kazuhito
    Dec 8, 2016 at 14:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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