##############################################################################################
# Clean NWIS Suspended Sediment Data
#
# This script cleans Discrete Suspended Sediment data compiled in the USGS Sediment Data Portal
# <http://cida.usgs.gov/sediment/> based on a list of pre-determined sites.
# In some cases, instantaneous or daily average flow values are missing. I also identify sites
# with < 15 measurements of suspended sediment in which no sand was observed.
#
# Joel Sholtes, PhD Candidate, Colorado State University, March 2014
##############################################################################################

rm(list=ls(all=TRUE))
library(MASS)

setwd("C:/Users/jsholtes/Dropbox/Research/Qeff/Analysis/SuspLoad")

# import list of sites to use
sites <- read.table("Susp_Sed_Sites_v10.csv",
                               header = TRUE, sep = ",")
# load and clean data  - probably only need to do this once
data <- read.table("C:/Users/jsholtes/Dropbox/Research/Qeff/Data/Suspended Load/discrete_data.csv", 
                   header = T, stringsAsFactors = F, sep = ",", comment.char = "#")
# subset from original USGS database only the sites that match my culled list
# based on studied sites
select_site <- subset(data, SITE_NO %in% sites$site.no)

# Check if all sites in 'sites' are in 'select_site'
sum( sort(as.numeric(unique(select_site$SITE_NO))) != sort(sites$site.no) ) == 0

# # remove any rows with pct < 0.063 mm = 100 (i.e., no sand)
#select_site <- subset(select_site, P63 != 100)

# assign any blank values with NA
select_site[select_site == ""]  <- NA

# selects only certain columns needed
select_site <- select_site[c("SITE_NO", "STATION_NM", "DATETIME", "SSC", "DAILYFLOW",
                            "INSTFLOW", "P63")]

#Convert select columns to numeric
select_site$SSC <- as.numeric(select_site$SSC)
select_site$DAILYFLOW <- as.numeric(select_site$DAILYFLOW)
select_site$INSTFLOW <- as.numeric(select_site$INSTFLOW)
select_site$P63 <- as.numeric(select_site$P63)

# remove any NA vals from SSC column or P63 column
select_site <- subset(select_site, !is.na(select_site$SSC))
select_site <- subset(select_site, !is.na(select_site$P63))

# remove rows with both INST and Daily flows as NA
select_site <- subset(select_site, !is.na(select_site$INSTFLOW) & 
                        !is.na(select_site$DAILYFLOW))

# remove rows with zero in one flow column and NA in the other
select_site <- subset(select_site, !is.na(select_site$INSTFLOW) & 
                        !select_site$DAILYFLOW == 0)
select_site <- subset(select_site, !select_site$INSTFLOW == 0 & 
                        !is.na(select_site$DAILYFLOW))

# remove rows where P63 > 100 or less than 0
select_site <- subset(select_site, !select_site$P63 > 100)

# remove rows with negative values in DAILYFLOW or INSTFLOW column
select_site <- subset(select_site, !select_site$DAILYFLOW < 0)
select_site <- subset(select_site, !select_site$INSTFLOW < 0)

# Save cleaned data
write.csv(select_site, "susp_data_clean2.csv", row.names = F)

sites2 <- unique(select_site$SITE_NO)