### https://cran.r-project.org/web/packages/sf/index.html
library(sf)
## Linking to GEOS 3.7.0, GDAL 2.4.0, PROJ 5.2.0
nc <- st_read(system.file("shape/nc.shp", package="sf"))
## Reading layer `nc' from data source `C:\Users\subas\Documents\R\win-library\3.6\sf\shape\nc.shp' using driver `ESRI Shapefile'
## Simple feature collection with 100 features and 14 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## epsg (SRID): 4267
## proj4string: +proj=longlat +datum=NAD27 +no_defs
## Simple feature collection with 100 features and 6 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## epsg (SRID): 4267
## proj4string: +proj=longlat +datum=NAD27 +no_defs
## First 3 features:
## BIR74 SID74 NWBIR74 BIR79 SID79 NWBIR79 geometry
## 1 1091 1 10 1364 0 19 MULTIPOLYGON (((-81.47276 3...
## 2 487 0 10 542 3 12 MULTIPOLYGON (((-81.23989 3...
## 3 3188 5 208 3616 6 260 MULTIPOLYGON (((-80.45634 3...
## [1] $<- [ [[<-
## [4] aggregate as.data.frame cbind
## [7] coerce dbDataType dbWriteTable
## [10] identify initialize merge
## [13] plot print rbind
## [16] show slotsFromS3 st_agr
## [19] st_agr<- st_area st_as_sf
## [22] st_bbox st_boundary st_buffer
## [25] st_cast st_centroid st_collection_extract
## [28] st_convex_hull st_coordinates st_crop
## [31] st_crs st_crs<- st_difference
## [34] st_geometry st_geometry<- st_interpolate_aw
## [37] st_intersection st_intersects st_is
## [40] st_line_merge st_nearest_points st_node
## [43] st_normalize st_point_on_surface st_polygonize
## [46] st_precision st_segmentize st_set_precision
## [49] st_simplify st_snap st_sym_difference
## [52] st_transform st_triangulate st_union
## [55] st_voronoi st_wrap_dateline st_write
## [58] st_zm
## see '?methods' for accessing help and source code
(nc_geom <- st_geometry(nc))
## Geometry set for 100 features
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## epsg (SRID): 4267
## proj4string: +proj=longlat +datum=NAD27 +no_defs
## First 5 geometries:
par(mar = c(0,0,1,0))
plot(nc[1], reset = FALSE) # reset = FALSE: we want to add to a plot with a legend
plot(nc[1,1], col = 'grey', add = TRUE)

(p = st_point(c(0,2)))
## POINT (0 2)
p + 1
## POINT (1 3)
p + c(1,2)
## POINT (1 4)
p + p
## POINT (0 4)
p * p
## POINT (0 4)
rot = function(a) matrix(c(cos(a), sin(a), -sin(a), cos(a)), 2, 2)
p * rot(pi/4)
## POINT (1.414214 1.414214)
p * rot(pi/2)
## POINT (2 1.224647e-16)
p * rot(pi)
## POINT (2.449294e-16 -2)
nc = st_read(system.file("shape/nc.shp", package="sf"), quiet = TRUE)
ncg = st_geometry(nc)
plot(ncg, border = 'grey')
cntrd = st_centroid(ncg)
## Warning in st_centroid.sfc(ncg): st_centroid does not give correct
## centroids for longitude/latitude data
ncg2 = (ncg - cntrd) * rot(pi/2) * .75 + cntrd
plot(ncg2, add = TRUE)
plot(cntrd, col = 'red', add = TRUE, cex = .5)
library(sf)
## Linking to GEOS 3.5.0, GDAL 2.2.2, PROJ 4.8.0
demo(nc, ask = FALSE, echo = FALSE)
## Reading layer `nc.gpkg' from data source `C:\Users\subas\Documents\R\win-library\3.6\sf\gpkg\nc.gpkg' using driver `GPKG'
## Simple feature collection with 100 features and 14 fields
## Attribute-geometry relationship: 0 constant, 8 aggregate, 6 identity
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## epsg (SRID): 4267
## proj4string: +proj=longlat +datum=NAD27 +no_defs
## Reading layer `nc.gpkg' from data source `/home/travis/R/Library/sf/gpkg/nc.gpkg' using driver `GPKG'
## Simple feature collection with 100 features and 14 fields
## Attribute-geometry relationship: 0 constant, 8 aggregate, 6 identity
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## epsg (SRID): 4267
## proj4string: +proj=longlat +datum=NAD27 +no_defs
plot(st_geometry(nc))

plot(st_geometry(nc), col = sf.colors(12, categorical = TRUE), border = 'grey',
axes = TRUE)
plot(st_geometry(st_centroid(nc)), pch = 3, col = 'red', add = TRUE)
## Warning in st_centroid.sf(nc): st_centroid assumes attributes are constant
## over geometries of x
## Warning in st_centroid.sfc(st_geometry(x), of_largest_polygon =
## of_largest_polygon): st_centroid does not give correct centroids for
## longitude/latitude data
plot(nc)



options(sf_max.plot=1)
plot(nc)

library(mapview)
mapview(nc["BIR74"], col.regions = sf.colors(10))
mm= st_read("C:/Users/subas/Syncplicity/RuralSpeedSafety (L.D. White)/app/www/OH_shp/OH_Fin2.shp")
## Reading layer `OH_Fin2' from data source `C:\Users\subas\Syncplicity\RuralSpeedSafety (L.D. White)\app\www\OH_shp\OH_Fin2.shp' using driver `ESRI Shapefile'
## Simple feature collection with 8466 features and 37 fields
## geometry type: MULTILINESTRING
## dimension: XY
## bbox: xmin: 400586.2 ymin: 60008.65 xmax: 766904.5 ymax: 434522.5
## epsg (SRID): NA
## proj4string: +proj=lcc +lat_1=38.73333333333333 +lat_2=40.03333333333333 +lat_0=38 +lon_0=-82.5 +x_0=600000 +y_0=0 +datum=NAD83 +units=m +no_defs

plot(mm["Total"], breaks = c(0,5,10,15,20,100))

plot(mm["FI"], breaks = c(0,5,10,15,20,25))

library(mapview)
mapview(mm["Total"], col.regions = sf.colors(10))