DeltaPull.Rmd

---
title: "ReflectancePull_Walkthrough"
author: "Simon Topp"
date: "5/8/2019"
output: html_document
editor_options: 
  chunk_output_type: console
---

```{r setup, include=FALSE}
library(tidyverse)
library(googledrive)
library(feather)
library(sf)
library(lubridate)

knitr::opts_chunk$set(echo = TRUE)
```

## Here we pull down the landsat stacks for all the Deltas. 

```{r}
## Shapefile from JM, correct the Dnieper which is mislabeled
geo <- st_read('data/in/DeltaSamplesJM.shp')

# Fix Dnieper
geo %>% filter(Delta == 'Dnieper')
geo$surface[geo$id %in% c('Dnieper_Water_Downstream','Dnieper_Water_Midstream','Dnieper_Water_Upstream')] <- 'Land'
geo$surface[geo$id %in% c('Dnieper_Land_Downstream','Dnieper_Land_Midstream','Dnieper_Land_Upstream')] <- 'Water'

geo <- geo %>% mutate(id = paste0(Delta,'_',surface,'_',location))

mapview(geo)

## Check for other goofs
unique(geo$location)
geo$location[geo$location == 'MIdstream'] <- 'Midstream'
unique(geo$surface)
geo$surface[geo$surface == "Land'"] <- 'Land'

## Re-write the shapefile for the pull
st_write(geo, 'data/in/DeltaSampleLocations_96.shp')
```

## Read in the data from EE, code for pull script can be found
## at https://code.earthengine.google.com/1d647a3d92568d22fe0b2f528b9e5ac2?noload=true

```{r}
deltas <- read_csv('data/in/DeltaEEExports/landExport.csv') %>%
  bind_rows(read_csv('data/in/DeltaEEExports/WaterExport.csv')) %>%
  rename_all(tolower) %>%
  select(ublue, blue, green, red, nir, swir1,swir2, cScore = cscore_clouds, date, id, pixelCount = pcount_prop, location, surface, Delta = delta)

deltas <- deltas %>%
  filter(pixelCount > 8) %>%
  pivot_longer(cols = c(blue, green, red, nir, swir1,swir2), names_to = 'band', values_to = 'value') %>%
  filter(value > 0,
         value < 10000,
         cScore < .1) %>%
  pivot_wider(names_from = band, values_from = value)

# Create clean dataset for Evan
dClean <- deltas %>%
  mutate(month = month(date),
         year = year(date),
         ndvi = (nir-red)/(nir+red),
         evi = 2.5*((nir-red)/(nir + 6*red - 7.5*blue + 1)),
         evi = ifelse(abs(evi) > 10000, NA, evi), 
         savi = ((nir-red)/(nir+red+0.5))*1.5,
         gr = green/red,
         ndssi = (blue - nir)/(nir + blue), #Reversed order of this from lit convention to make higher NDSSI = more sediment
         nsmi = (red + green - blue)/(red + green + blue)) %>% #Added NSMI to these indices
  group_by(Delta, location, surface, year, month) %>%
  summarize(ndvi = mean(ndvi, na.rm = T),
            red = mean(red, na.rm = T),
            evi = mean(evi, na.rm = T),
            savi = mean(savi, na.rm = T),
            gr = mean(gr, na.rm = T),
            ndssi = mean(ndssi, na.rm = T),
            nsmi = mean(nsmi, na.rm = T))

write_csv(dClean, 'data/out/deltas_clean_96.csv')

## Pull Centerpoints for all the deltas.
shp <- st_read('data/in/DeltaSampleLocations_96.shp') %>%
  st_cast('MULTIPOINT') %>%
  group_by(Delta) %>%
  summarise(dummy = 1) %>%
  st_centroid() %>%
  mutate(Lon = map_dbl(geometry,1),
         Lat = map_dbl(geometry,2)) %>%
  select(-dummy, Deltas = Delta) %>%
  st_set_geometry(NULL)

write_csv(shp, 'data/DeltaLocations_96.csv')




```