library(mapsapi)
library(sp)
library(sf)
library(tidyverse)
library(leaflet)

#Google API key
key='WRITE YOUR API KEY'

#load point layers

#layer of wastewater treatment plants
#with manually predefined waypoints of each kind of weekly route
wwtp <- st_read("wwtp.geojson")

#layer of lab1 & 2 (they are in the same building)
lab1_2 <- st_read("lab1_2.geojson")

#layer of lab3
lab3 <- st_read("lab3.geojson")

#reprojct all layers to WGS84
wwtp <- st_transform(wwtp, 4326)
lab1_2 <- st_transform(lab1_2, 4326)
lab3 <- st_transform(lab3, 4326)

#select the kind of weekly route: summer_1, summer_2, winter_1, winter_2
week <- "summer_1"

#select wwtp from that weekly route
wwtp_week <- wwtp %>% filter(!is.na(matches(week)))

#define origins of each route
origin1 <- filter(wwtp_week, wwtp_id == "DOLO")
origin2 <- filter(wwtp_week, wwtp_id == "DRSS")
origin3 <- filter(wwtp_week, wwtp_id == "DCPA")
origin4 <- filter(wwtp_week, wwtp_id == "DVEN")
origin5 <- filter(wwtp_week, wwtp_id == "DPUI")
origin51 <- filter(wwtp_week, wwtp_id == "DAMP")
origin52 <- filter(wwtp_week, wwtp_id == "DLLE")
origin53 <- filter(wwtp_week, wwtp_id == "DVIE")

#Run through each route
for (i in unique(wwtp_week[week])){
  print(i)
  #define parameters for mp_directions
  origin <- eval(parse(text=paste("origin",i,sep="")))
  waypoints <- filter(wwtp_week, matches(week) == i)

  #determine final destination using i value
  if (i < 50){
    destination <- lab1_2
  }
  else if (i > 50){
    destination <- lab3
  }else{print("i out of in range")}

  #define departure time as 8:00am
  departure_time <- as.POSIXct(strptime("2020-07-06 08:00:00", "%Y-%m-%d %H:%M:%OS"))

  #Calculate route i
  route <- mp_directions(
    origin = origin,
    waypoints = waypoints,
    destination = destination,
    departure_time = departure_time,
    region = "es",
    traffic_model = "best_guess",
    key = key,
    #quiet = TRUE
  )

  #convert to sf object
  st_route <- mp_get_routes(route)

  #create a field with route id
  st_route$route_id <- i

  #assign a name to the route
  assign(paste("route_",i,sep=""), st_route)

  #merge all routes
  if (i == 1){
    route_all <- route_1
  } else {
    route_all <- rbind(route_all, st_route)
  }
}
plot(st_geometry(route_all))

#export to geojson in ETRS89 (ESPG: 25831)
route_all <- st_transform(route_all, 25831)
st_write(route_all,
         dsn=paste(week, ".geojson", sep=""), delete_dsn = TRUE)

#create two fields with route's length in kn and hours
route_all$km <- route_all$distance_m / 1000
route_all$hours <- route_all$duration_s / 3600

#create data frame to summarize results
summary_route <- as.data.frame(table(route_all$route_id))

#summarize results for each route
for (i in summary_route$Var1){
  summary_route$km[summary_route$Var1==i] <- sum(route_all$km[route_all$route_id==i])
  summary_route$durada[summary_route$Var1==i] <- sum(route_all$hores[route_all$route_id==i])
}

#add 20 minutes for each stop in a route
summary_route$Freq <- summary_route$Freq-1
summary_route$dur_parades <- summary_route$durada+(0.33*summary_route$Freq)

#export summary_route to csv for excel
write.csv2(summary_route,
           file=paste(week, "summary.csv", sep="_"))

#create csv for excel with stops and coordinates
stops <- select(wwtp_week, matches(week), wwtp_id, county, province, x, y)
write.csv2(stops,
           file=paste(week, "stops.csv", sep="_"))