Setup

library(tidyverse)

## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──

## ✓ ggplot2 3.3.3     ✓ purrr   0.3.4
## ✓ tibble  3.1.2     ✓ dplyr   1.0.6
## ✓ tidyr   1.1.3     ✓ stringr 1.4.0
## ✓ readr   1.4.0     ✓ forcats 0.5.1

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(sf)

## Linking to GEOS 3.8.1, GDAL 3.2.1, PROJ 7.2.1

library(raster)

## Loading required package: sp

## 
## Attaching package: 'raster'

## The following object is masked from 'package:dplyr':
## 
##     select

## The following object is masked from 'package:tidyr':
## 
##     extract

Shoreline

Read it and crop it.

shoreline=st_read("Shoreline.shp")

## Reading layer `Shoreline' from data source `/Users/haroldnelson/Dropbox/QGIS/Tobias/Shoreline.shp' using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 6 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -123.1086 ymin: 37.69325 xmax: -122.3277 ymax: 37.92975
## Geodetic CRS:  WGS84(DD)

shoreline = st_crop(shoreline,c(xmin = -122.6,xmax = 0,ymin = 0,ymax  = 37.9))

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

Streets

streets = st_read("StreetCenterlines.shp")

## Reading layer `StreetCenterlines' from data source `/Users/haroldnelson/Dropbox/QGIS/Tobias/StreetCenterlines.shp' using driver `ESRI Shapefile'
## Simple feature collection with 16241 features and 21 fields
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: -122.5136 ymin: 37.70702 xmax: -122.3583 ymax: 37.83213
## Geodetic CRS:  WGS84(DD)

streets = streets %>% dplyr::select(classcode)
streets = streets %>% 
  filter(classcode == "1")

Trees

trees = read_csv("Street_Tree_Map.csv")

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   TreeID = col_double(),
##   qLegalStatus = col_character(),
##   qSpecies = col_character(),
##   qAddress = col_character(),
##   SiteOrder = col_double(),
##   qSiteInfo = col_character(),
##   PlantType = col_character(),
##   qCaretaker = col_character(),
##   qCareAssistant = col_character(),
##   PlantDate = col_character(),
##   DBH = col_double(),
##   PlotSize = col_character(),
##   PermitNotes = col_character(),
##   XCoord = col_double(),
##   YCoord = col_double(),
##   Latitude = col_double(),
##   Longitude = col_double(),
##   Location = col_character()
## )

trees = trees %>% 
  dplyr::select(PlantType, qSpecies,Latitude,Longitude) %>% na.omit()
trees = trees %>% filter(Latitude < 37.9 & Latitude > 37.6)

sycamores = trees %>% 
  filter(str_detect(qSpecies,"Sycamore"))
sycamores_sf = sycamores %>% st_as_sf(coords=c("Longitude","Latitude"),crs = 4326)

DEM

DEM = raster("DEM_SF.tif")

DEM_tr = projectRaster(DEM,crs = crs(shoreline))

## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj
## = prefer_proj): Discarded datum Unknown based on WGS84 ellipsoid in Proj4
## definition

## Warning in showSRID(SRS_string, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded datum WGS84 in Proj4 definition

DEM_tr_df = as.data.frame(DEM_tr, xy = TRUE)

A Map

ggplot(data = shoreline) +
  geom_sf() +
  geom_sf(data = streets) +
  geom_sf(data = sycamores_sf,size=.01,color = "green") +
  geom_raster(data = DEM_tr_df,aes(x=x,y=y,fill = DEM_SF),
              alpha=.5,)

Make Themes Available

library(ggthemes)

Crop the DEM

DEM_crop = DEM_tr_df %>% 
  filter(x > -122.54 &
         x < -122.35 &
         y > 37.7 &
         y < 37.82)

Plot the Cropped Map

I also played with the size parameter of shoreline and the alpha parameter of the raster.

ggplot(data = shoreline) +
  geom_sf(size = .2) +
  geom_sf(data = streets) +
  geom_sf(data = sycamores_sf,size=.01,color = "green") +
  geom_raster(data = DEM_crop,aes(x=x,y=y,fill = DEM_SF),
              alpha=.8) +
  theme_bw()

tmap

library(tmap)

Shoreline with tmap

tm_shape(shoreline) +
  tm_polygons()

Add the DEM + streets + trees

tm_shape(shoreline) +
  tm_polygons() +  
tm_shape(DEM_tr) +
  tm_raster() +
tm_shape(streets) +
  tm_lines() +
tm_shape(sycamores_sf) +
  tm_dots()

## stars object downsampled to 1099 by 910 cells. See tm_shape manual (argument raster.downsample)

## Variable(s) "NA" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.