# Hengistbury

```{r setup_hh, include=FALSE}
knitr::opts_chunk$set(echo = TRUE,message = FALSE)
```


```{r}
library(raster)
library(mapview)
library(tmap)
library(giscourse)
library(sf)
library(tidyverse)

conn<-connect()


mapviewOptions(basemaps =c("Esri.WorldImagery","OpenStreetMap","Esri.WorldStreetMap", "CartoDB.Positron", "OpenTopoMap", "OpenStreetMap.BlackAndWhite", "OpenStreetMap.HOT","Stamen.TonerLite", "CartoDB.DarkMatter", "Esri.WorldShadedRelief", "Stamen.Terrain"))
```



```{r}


query<-"select geom,lnr_name from local_reserves where lnr_name like 'Hengist%'"
lnr<-st_read(conn,query=query)
bound<-as(lnr, "Spatial")
mapview(lnr,  alpha.regions = 0, legend=FALSE, lwd =3) %>% extras()

```


```{r}
hh_habitat<-sprintf("select main_habit mhabitat, 
s1habtype,s2habtype,s3habtype habitat,
st_intersection(h.geom,s.geom) geom from 
(%s) s,
ph_v2_1 h 
where st_intersects(h.geom,s.geom)",query)

habitats<-st_read(conn, query=hh_habitat)

habitats$area <- st_area(habitats)

habitats%>% st_buffer(0.001)%>% group_by(mhabitat) %>% summarise(area=sum(area))-> habitats



qtm(habitats, "mhabitat") +tm_style("classic")
```

```{r}
mapview(habitats, zcol="mhabitat", burst=TRUE,legend=FALSE,  hide = TRUE) %>% extras()
```



```{r}
dsm<- pgGetRast(conn, "dsm2m",boundary=bound)
dtm<- pgGetRast(conn, "dtm2m",boundary=bound)

dsm[dsm<0.1]<-NA
dtm[dtm<0.1]<-NA
```

```{r}
slope<-terrain(dtm, opt='slope', unit='degrees')
aspect<-terrain(dtm, opt='aspect',unit='degrees')
sloper<-terrain(dtm, opt='slope', unit='radians')
aspectr<-terrain(dtm, opt='aspect',unit='radians')
hillshade<-hillShade(sloper,aspectr,angle=5)


plot(hillshade, col = gray(0:100 / 100), legend = FALSE)
plot(dtm, col = terrain.colors(25), alpha = 0.3, legend = TRUE, add = TRUE)
```



