README.Rmd

---
output: github_document
---

<!-- badges: start -->
[![CRAN](http://www.r-pkg.org/badges/version/rsgeo)](https://cran.r-project.org/package=rsgeo)
[![R-CMD-check](https://github.com/JosiahParry/rsgeo/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/JosiahParry/rsgeo/actions/workflows/R-CMD-check.yaml)
<!-- badges: end -->


```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/README-",
  out.width = "100%"
)
```

# rsgeo

`rsgeo` is an interface to the Rust libraries `geo-types` and `geo`. `geo-types` implements pure rust geometry primitives. The `geo` library adds additional algorithm functionalities on top of `geo-types`. This package lets you harness the speed, safety, and memory efficiency of these libraries. `geo-types` does not support Z or M dimensions. There is no support for CRS at this moment.

```{r}
# install.packages(
#   'rsgeo', 
#   repos = c('https://josiahparry.r-universe.dev', 'https://cloud.r-project.org')
# )
library(rsgeo)
```

rsgeo works with vectors of geometries. When we compare this to `sf` this is always the geometry column which is a class `sfc` object (simple feature column).

```{r}
# get geometry from sf
data(guerry, package = "sfdep")

polys <- guerry[["geometry"]] |>
  sf::st_cast("POLYGON")

# cast to rust geo-types
rs_polys <- as_rsgeo(polys)

head(rs_polys)
```

Cast geometries to sf

```{r}
sf::st_as_sfc(rs_polys)
```


Calculate the unsigned area of polygons. 

```{r}
bench::mark(
  rust = unsigned_area(rs_polys),
  sf = sf::st_area(polys),
  check = FALSE
)
```
Find centroids 

```{r}
bench::mark(
  centroids(rs_polys),
  sf::st_centroid(polys),
  check = FALSE
)
```
Extract points coordinates

```{r}
coords(rs_polys) |> 
  head()
```

Plot the polygons and their centroids

```{r}
plot(rs_polys)
plot(centroids(rs_polys), add = TRUE)
```

Calculate a distance matrix. Note that there is often floating point error differences so `check = FALSE` in this case.

```{r}
pnts <- centroids(rs_polys)
pnts_sf <- sf::st_as_sfc(pnts)

bench::mark(
  rust = distance_euclidean_matrix(pnts, pnts),
  sf = sf::st_distance(pnts_sf, pnts_sf),
  check = FALSE
)

```

Simplify geometries. 

```{r}
x <- rs_polys
x_simple <- simplify_geoms(x, 5000)

plot(x_simple)
```

```{r}
bench::mark(
  rust = simplify_geoms(rs_polys, 500),
  sf = sf::st_simplify(polys, FALSE, 500),
  check = FALSE
)
```

Union geometries with `union_geoms()`. Some things sf is better at! One of which is performing unary unions of complex geometries.


```{r}
plot(union_geoms(rs_polys))

bench::mark(
  union_geoms(rs_polys),
  sf::st_union(polys),
  check = FALSE
)
```

We can cast between geometries as well.

```{r}
lns <- cast_geoms(rs_polys, "linestring")
```

Some unions are faster when using rsgeo vectors like linestrings. 

```{r}
lns_sf <- sf::st_cast(polys, "LINESTRING")

bench::mark(
  union_geoms(lns),
  sf::st_union(lns_sf),
  check = FALSE
)
```


Find the closest point to a geometry

```{r}
close_pnt <- closest_point(
  rs_polys, 
  geom_point(800000, 2090000)
)

plot(rs_polys[1])
plot(close_pnt, pch = 15, add = TRUE)
```

Find the haversine destination of a point, bearing, and distance. Compare to the very fast geosphere destination point function. 

```{r}
bench::mark(
  rust = haversine_destination(geom_point(10, 10), 45, 10000),
  Cpp = geosphere::destPoint(c(10, 10), 45, 10000),
  check = FALSE
)
```


```{r}
origin <- geom_point(10, 10)

destination <- haversine_destination(origin, 45, 10000)

plot(c(origin, destination), col = c("red", "blue"))
```

Find intermediate point on a great circle.

```{r}
middle <- haversine_intermediate(origin, destination, 1/2)

plot(origin)
plot(destination, add = TRUE, col = "red")
plot(middle, add = TRUE, col = "blue")
```

<!-- Utilize the chaikin smoothing algorithm with 5 iterations. -->

<!-- ```{r} -->
<!-- region <- rs_polys[[2]] -->
<!-- plot(chaikin_smoothing(region, 5)) -->
<!-- ``` -->

Find extreme coordinates with `extreme_coords()`


```{r}
france <- union_geoms(rs_polys)

plot(france)
plot(extreme_coords(france)[[1]], add = TRUE, pch = 15)
```

Get bounding rectangles

```{r}
rects <- bounding_rect(rs_polys)
plot(rects)
```

Convex hulls

```{r}
convex_hull(rs_polys) |> 
  plot()
```

Expand into constituent geometries as a list of geometry vectors

```{r}
expand_geoms(rs_polys) |> 
  head()
```
We can flatten the resultant geometries into a single vector using `flatten_geoms()`

```{r}
expand_geoms(rs_polys) |> 
  flatten_geoms() |> 
  head()
```

Combine geometries into a single multi- geometry

```{r}
combine_geoms(lns)
```


Spatial predicates 

```{r}
x <- rs_polys[1:5]
intersects_sparse(x, rs_polys)
```


<!-- Convert to and from wkb and wkt -->

<!-- ```{r} -->
<!-- wkt <- wkt_from_geoms(x) -->
<!-- wkt_to_geoms(wkt) -->
<!-- ``` -->

<!-- ```{r} -->
<!-- wkb <- wkb_from_geoms(x) -->
<!-- head(wkb[[1]]) -->

<!-- wkb_to_geoms(wkb) -->
<!-- ``` -->


#### Notes

Right now plotting is done using `wk` by first casting the rsgeo into an sfc object.