Reputation: 67
I am looking to calculate distances from NA cells to a spatial object. However, the distance() function is not recognizing NA cells that I have created. Here is reproducible code:
To make spatial object:
library(raster)
CAN <- getData("GADM", country = "Canada", level = 1)#level 1 just gives borders for provinces
library(sf)#to make sf object
sf_CAN <- as(CAN,"sf")
st_crop_CAN <- st_crop(sf_CAN, xmin = -60.17, xmax = -59.65, ymin = 43.92, ymax = 44.05) #bounding box gives sable island
Create raster file and use mask to then use distance():
library(marmap)
shelf <- getNOAA.bathy(lon1 = -75, lon2 = -45,
lat1 = 35, lat2 = 60, resolution = 4)
library(raster)
CAN <- getData("GADM", country = "Canada", level = 1)#level 1 just gives borders for provinces
library(sf)#to make sf object
sf_CAN <- as(CAN,"sf")
st_crop_CAN <- st_crop(sf_CAN, xmin = -60.17, xmax = -59.65, ymin = 43.92, ymax = 44.05) #bounding box gives sable island
shelf_dis <- marmap::as.raster(shelf)
values(shelf_dis) <- NA
newR <- mask(shelf_dis,as(st_crop_CAN$geometry, "Spatial"))
dist_sable <- distance(newR)
#error code after distance()
Error in .local(x, y, ...) :
RasterLayer has no NA cells (for which to compute a distance)
One issue that I see is that in the min and max values of newR object are NA as shown below:
> newR
class : RasterLayer
dimensions : 375, 450, 168750 (nrow, ncol, ncell)
resolution : 0.06666667, 0.06666667 (x, y)
extent : -75, -45, 35, 60 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs
source : memory
names : layer
values : NA, NA (min, max)
So it is unable to distinguish between NA cells and the spatial object because there are only NA cells to begin with?
I am not sure if my approach to this is the best way. Any advice or help appreciated. At the end of the day I am just trying to get distances from points (not provided) to the spatial object. That part seems easy enough with extract(). I just need to get this new raster file! This is my first go at spatial work in R.
Upvotes: 0
Views: 244
Reputation: 47411
Yes, there must be cells that are not NA
to compute distance to. The general approach is usually like this
library(raster)
# example polygons
filename <- system.file("external/lux.shp", package="raster")
p <- shapefile(filename)
pp <- crop(p, extent(5.698, 6.231, 49.442, 49.781))
# create or use existing raster
r <- raster(p, res=0.01)
r <- rasterize(pp, r)
d <- distance(r)
d
#class : RasterLayer
#dimensions : 73, 78, 5694 (nrow, ncol, ncell)
#resolution : 0.01, 0.01 (x, y)
#extent : 5.74414, 6.52414, 49.45162, 50.18162 (xmin, xmax, ymin, ymax)
#crs : +proj=longlat +datum=WGS84 +no_defs
#source : memory
#names : layer
#values : 0, 49114.11 (min, max)
Upvotes: 1