# R/costrasterGen.R In ipdw: Spatial Interpolation by Inverse Path Distance Weighting

```#'@name costrasterGen
#'@title Generate a cost Raster
#'@description Generate a cost raster from an object of class \code{SpatialPolygons}, \code{matrix}, or \code{SpatialPointsDataFrame}
#'@author Joseph Stachelek
#'@param xymat Matrix of coordinates or a \code{SpatialPointsDataFrame} object
#'@param pols SpatialPolygons object
#'@param extent Define extent based on extent of xymat/xyspdf (points) or pols (polys). Default is polys.
#'@param projstr proj4 string defining the inherent projection
#'@details Ensure that the projection of the xymat coordinates and pols match. This can be accomplished by running the \code{projection} command on both inputs. If they do not match use the \code{spTransform} command.
#'@return RasterLayer
#'@importFrom raster raster rasterize reclassify
#'@importFrom sp proj4string coordinates bbox
#'@export
#'@examples \dontrun{
#'Sr1 <- Polygon(cbind(c(0, 0, 1, 1, 0), c(0, 12, 12, 0, 0)))
#'Sr4 <- Polygon(cbind(c(9, 9, 10, 10, 9), c(0, 12, 12, 0, 0)))
#'Sr2 <- Polygon(cbind(c(1, 1, 9, 9, 1), c(11, 12, 12, 11, 11)))
#'Sr3 <- Polygon(cbind(c(1, 1, 9, 9, 1),c(0, 1, 1, 0, 0)))
#'Sr5 <- Polygon(cbind(c(4, 4, 5, 5, 4), c(4, 8, 8, 4, 4)))

#'Srs1 <- Polygons(list(Sr1), "s1")
#'Srs2 <- Polygons(list(Sr2), "s2")
#'Srs3 <- Polygons(list(Sr3), "s3")
#'Srs4 <- Polygons(list(Sr4), "s4")
#'Srs5 <- Polygons(list(Sr5), "s5")
#'
#'pols <- SpatialPolygons(list(Srs1, Srs2, Srs3, Srs4, Srs5), 1:5)
#'
#'#using a matrix object
#'xymat <- matrix(3, 3, nrow = 1, ncol = 2)
#'costras <- costrasterGen(xymat, pols, projstr = NULL)
#'
#'#plotting
#'plot(costras)
#'points(xymat)
#'}

'costrasterGen' <- function(xymat, pols, extent = "polys", projstr){
if(class(xymat) == "SpatialPointsDataFrame" |
class(xymat) == "SpatialPoints"){
xymat <- sp::coordinates(xymat)
}

#add check to see if projstr and projection(pols) match
if(!identical(projstr, sp::proj4string(pols))){
message("Warning, the projection of polygons does not
match projstr. See rgdal::spTransform")
}

#define spatial domain based on pnts or polys
if(extent == "polys"){
xmn <- min(sp::bbox(pols)[1,])
xmx <- max(sp::bbox(pols)[1,])
ymn <- min(sp::bbox(pols)[2,])
ymx <- max(sp::bbox(pols)[2,])
}

if(extent == "points" | extent == "pnts"){
ymn <- range(xymat[,2])[1]
ymx <- range(xymat[,2])[2]
xmn <- range(xymat[,1])[1]
xmx <- range(xymat[,1])[2]
}

nrow <- ymx - ymn
ncol <- xmx - xmn

#generate cost raster
r <- raster::raster(nrow = nrow, ncol = ncol, crs = projstr, xmn = xmn,
xmx = xmx, ymn = ymn, ymx = ymx)
costras <- raster::rasterize(pols, r, silent = TRUE)
m <- c(0, +Inf, 10000)
rclmat <- matrix(m, ncol = 3, byrow = TRUE)
costras <- raster::reclassify(costras, rclmat)
costras <- raster::reclassify(costras, cbind(NA, 1))

return(costras)
}
```

## Try the ipdw package in your browser

Any scripts or data that you put into this service are public.

ipdw documentation built on May 29, 2017, 1:53 p.m.