2

I would like to extract values from multiple rasters with a single lat / long or Shapefile. I have tried every tutorial and stackoverflow post and am getting nowhere. Maybe someone could please help. Here is what I have tried so far:

# Imports
library(raster)
library(sp)

# Raster List
tmp <- lapply(list.files("D:/rasters", pattern = ".tif$", full.names = TRUE), raster)

# Lat and Lon
coords <- data.frame("lat" = c(29.940668), "lon" = c(48.363455))

# Create spatial points
pts <- SpatialPoints(coords = coords, 
              proj4string = CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0  "))

# Extract the points
as.data.frame(cbind(coords, do.call("rbind", lapply(tmp, extract, pts))))

The code runs but the return values are all NA and I checked and I should be getting back some values.

I have also tried to read the values in as a shapefile using readOGR as in Extracting values from raster according to Lat and long of values?, but I get an error

Error in ogrInfo(dsn = dsn, layer = layer, encoding = encoding, use_iconv = use_iconv, : Cannot open layer

I'm at a loss why nothing is working for me.

  • 1
    What spatial reference does your raster(s) have? You may need to project your points. – Michael Stimson Oct 29 '18 at 23:18
  • Thanks, I am guessing you are asking regarding the readOGR issue. Both are in Geographic Coordinate systems - GCS WGS1984. Should both file types be projected? – Rachel Rotz Oct 29 '18 at 23:20
  • 1
    No, projected or geographic shouldn't matter so long as they are both the same. – Michael Stimson Oct 29 '18 at 23:25
  • 1
    Can you confirm that Iraq is your area of interest (or, at least, this point should be in Iraq)? – phloem Oct 29 '18 at 23:28
  • 1
    According to this question, lat/lon should not be in quotes: stackoverflow.com/questions/49181715/… – phloem Oct 30 '18 at 15:47
1

Thanks all! I finally got it working per the code below.

# Imports
library(raster)
library(sp)

# Raster List
tmp <- lapply(list.files("D:/rasters", pattern = ".tif$", full.names = TRUE), raster)

# Lon and lat without quotes per phloem's comment above  
coords <- data.frame(lon = 47.667, lat = 29.317)

# Create spatial points
pts <- SpatialPoints(coords = coords, 
              proj4string = CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0  "))

# Extract the points
result <- as.data.frame(cbind(coords, do.call("rbind", lapply(tmp, extract, pts))))

# Write to csv file
write.csv(result, file="filename.csv")

I was also able to finally read in the shapefile too but abandoned that method since the first one finally worked. Nonetheless, here is the code that finally worked to read in the shapefile:

point_location <- readOGR("foldername","shapefile name without extension")

Thanks again!

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.