# Raster Tutorial
# Author: Billy Raseman based on "Raster_Tutorial.R" by Jenna Stewart

# Notes: Many types of Rasters (Tiff, GeoTiff, ascii, netcdf, etc.).
# Examples include (DEM, NAtional Land Cover, flow direction, soils data, etc.).

# Clear Environment
rm(list=ls())

# Load packages
library(raster)     
library(rasterVis)
library(rgl)
library(rgdal)
library(ggplot2)

# Set working directory to where data is
setwd("C:/Users/Billy/Google Drive/CU Boulder Courses/CVEN 5454 Statistics and Probability/RMarkdown_CVEN-5454/scripts/")

##############
# Tutorial 1 #
##############
# source: http://neondataskills.org/R/Raster-Data-In-R/

# open the raster files
mask = raster("../data/mask.asc")
dem = raster("../data/dem.asc")

# get information on the rasters
mask
dem

# get the coordinate system
mask@crs
dem@crs

# plot initial images
image(mask)
image(dem)

# change the colors
col1=topo.colors(5)
image(dem,col=col1)

# add more colors
col2=topo.colors(10)
image(dem,col=col2)

# change the color scheme
col3=terrain.colors(10)
image(dem,col=col3)

# calculate the histogram
hist(dem, main="Distribution of elevation values", 
     col= "purple", 
     maxpixels=21752940)

# crop the raster to the mask, then plot the new cropped raster
dem_crop = dem*mask
# plot the cropped extent
image(dem_crop,col=col3)

##################
# ggplot section #
##################

# convert raster to points
dem.p <- rasterToPoints(dem_crop)
# create a dataframe and rename the columns
df.p <- data.frame(dem.p)
colnames(df.p) <- c('Easting', 'Northing', 'Elev')

# plot!
ggplot(data=df.p, aes(y=Northing, x=Easting)) +
  geom_raster(aes(fill=Elev))