04-projection-section.Rmd

---
title: "Geospatial tutorial - section 4"
author: "JR Ferrer-Paris"
---
## Explore lines with `sf`

### Hiking in Franconia

This data set has a selection of hiking trails in Franconia, the geometry here is a _multistring_, which are collections of lines. Notice that in this case we are using a different EPSG code, it refers to a specific zone of the Universal Transverse Mercator projection for this region of the world.

```{r linesex1, exercise=TRUE, paged.print=FALSE}
trails
```

Again, if you see a warning about old-style _crs_, then reassign the same _crs_.

```{r}
st_crs(trails) <- "EPSG:32632"
```

---
The [Universal Transverse Mercator (UTM)](https://en.wikipedia.org/wiki/Universal_Transverse_Mercator_coordinate_system) is a map projection system for assigning coordinates to locations on the surface of the Earth. The UTM system divides the Earth into 60 zones, each 6° of longitude in width, and uses a projection that can map a region of large north-south extent with low distortion.
---


### Length and area of spatial objects

We are often interested in the spatial dimensions of object, for example in this case we want to know the length of each trail, we use the function `st_length`:

```{r linesex2, exercise=TRUE}
trails %>% st_length() %>% hist(main="Length of trails")
```

###
Often you want to add information to your spatial object, so that you can use that information latter in your map or in other analysis, here we will use the function `mutate` to add the length as a variable of the trails object, and then find the longest trail for a hiking excursion.

```{r linesex3, exercise=TRUE}
trails %>% mutate(length=st_length(geometry)) -> trails
trails %>% slice(which.max(length))
```

###
We can also filter by a desired length and create objects for short and long trails:

```{r}
trails %>% mutate(length=st_length(geometry)) %>% filter(length>set_units(3000,'m'),length<set_units(3400,'m')) -> short_trails

trails %>% mutate(length=st_length(geometry)) %>% filter(length>set_units(40000,'m')) -> long_trails
```


### Projections and plots (maps)

When plotting object, `geom_sf` will handle the projections automatically, in this case the trails object is projected to geographical coordinates (some people would call this "inverse projection"). Notice the [graticule](https://www.thefreedictionary.com/graticule) of geographical coordinates in the background is rectilinear (straight lines in the background)?

```{r}
ggplot() +
  geom_sf(data = franconia, aes(fill=district)) +
  geom_sf(data = UBT, cex=2, col="yellow") +
  geom_sf(data = long_trails, col="whitesmoke")
```

###
If we want to have all our data in the UTM projection, we could invert the order in which we call the objects:

```{r linesex4, exercise=TRUE}
ggplot() +
  geom_sf(data=long_trails, col="whitesmoke") +
  geom_sf(data = franconia , aes(fill=district)) +
  geom_sf(data=UBT , cex=2, col="yellow")
```

How does this look? Do you notice the graticule behind our object has changed?

Unfortunately the lines are overlaid with polygons and we can not see them. So this is NOT the best way to do this.

###
If we want to have all our data projected on the fly, we can add a call to `coord_sf` and specify our choice of _crs_:

```{r linesex5, exercise=TRUE}
ggplot() +
  geom_sf(data = franconia, aes(fill=district)) +
  geom_sf(data = UBT, cex=2, col="yellow") +
  geom_sf(data = long_trails, col="whitesmoke") +
  coord_sf(crs=st_crs(long_trails))
```


###

However, if we want to work with the projected data more than once, we can save time by doing the transformation once and then using the projected objects instead of the originals.

Let's try a projection, just hit the run button to see what this means:

```{r projectionex1, exercise = TRUE, paged.print=FALSE}
(UBT_utm = UBT %>% st_transform(st_crs(trails)))
```

```{r projectionex1-check}
"UBT is now projected in Universal Transverse Mercator!"
```

Can you do the same with the franconia dataset?

```{r projectionex2, exercise = TRUE, paged.print=FALSE}
# franconia_utm = ____ %>% ____(____)
```

```{r projectionex2-solution}
franconia_utm = franconia %>% st_transform(st_crs(trails))
```

```{r projectionex2-check}
"the municipalities polygon should now be projected in Universal Transverse Mercator!"
```

###

Check how this look like

```{r}

ggplot() +
  geom_sf(data = franconia_utm, aes(fill=district)) +
  geom_sf(data = UBT_utm, cex=2, col="yellow") +
  geom_sf(data =long_trails, col="whitesmoke")

```