compare-load-reality

library(lubridate) # Include to avoid lubridate:: all the time
imsbasics::clc()
path <- paste0(getwd(), "/data/RData/SG/2014/")
t0_date = lubridate::ymd_hms("2014-01-01 00:00:00")


# ==============================================================================
# Files to load:
hist_scenario <- imsbasics::load_rdata("hist_scenario", path)
hist_scenario$missions <- hist_scenario$missions[!is.na(hist_scenario$missions$lat), ]
hist_scenario$missions <- hist_scenario$missions[!is.na(hist_scenario$missions$lng), ]
hist_scenario$vehicles <- hist_scenario$vehicles
hist_scenario$hist_events <- hist_scenario$events

# myvariables <- data911::get_ims_variables()$names[data911::get_ims_variables()$used_in_sim911]
select_col <- c("id", "vehicle_name", "vehicle_id", "event_id", "priority", "t_alarm_sec",
  "t_leave_poa_sec", "t_arrive_poa_sec", "dt_to_launch", "dt_to_poa", 
  "dt_to_destination", "dt_to_base", "dt_to_recondition",
  "dt_to_completion", "lat", "lng")
hist_scenario$missions <- hist_scenario$missions[, select_col]
missions <- hist_scenario$missions

# sim_missions
hist_vehicles <- hist_scenario$vehicles
hist_vehicles <- hist_vehicles[complete.cases(hist_vehicles$shift_weekday), ]
rm(path, hist_scenario)

# ==============================================================================
# Functions

# Use as:
# wdays_english <- weekdays_abbr()$english
# wdays_german <- weekdays_abbr()$german
weekdays_abbr <- function() { # ................................................ imsbasics
  german <- c("So", "Mo", "Di", "Mi", "Do", "Fr", "Sa")
  english <- c("Sun", "Mon", "Tues", "Wed", "Thurs", "Fri", "Sat") # lubridate!
  return(list(german = german, english = english))
}

g2e <- function(weekday_ger) {
  res <- character(length(weekday_ger))
  for (i in 1:length(weekday_ger)) {
    res[i] <- wdays_english[which(wdays_german == weekday_ger[i])]
  }
  return(res)
}

e2g <- function(weekday_eng) {
  res <- character(length(weekday_eng))
  for (i in 1:length(weekday_eng)) {
    res[i] <- wdays_german[which(wdays_english == weekday_eng[i])]
  }
  return(res)
}

is_busy <- function(t, t0, t1) {
  y <- fBasics::Heaviside(t,t0) * fBasics::Heaviside(-t,-t1)
  return(y)
}

hist_t_alarm_hour <- function(sim_missions, scale = 1) {
  t_alarm_hour <- plyr::count(trunc(simTimeR::simTime(sim_missions$t_alarm_sec)/60/60))
#   barplot(
#     t_alarm_hour$freq/scale,
#     main = "Einsatzverteilung über den Tag",
#     ylab = "Anzahl Einsätze",
#     xlab = "Stunde",
#     col = imsbasics::fhs(),
#     names.arg = c(8:18)
#     )
  plot(t_alarm_hour$x,
    t_alarm_hour$freq/scale,
    main = "Einsatzverteilung über den Tag",
    ylab = "Anzahl Einsätze",
    xlab = "Stunde",
    col = imsbasics::fhs(),
    # names.arg = c(0:23)
  )
}

compare_t_alarm_hour <- function(missions, available, scale = c(1,1)) {
  t_alarm_hour <- plyr::count(trunc(simTimeR::simTime(missions$t_alarm_sec)/60/60))
  bplt <- barplot(t_alarm_hour$freq/scale[1],
    main = "Einsatzverteilung über den Tag",
    ylab = "Anzahl Einsätze",
    xlab = "Stunde",
    col = imsbasics::fhs(),
    names.arg = c(0:23)
    )
  av_alarm_hour <- plyr::count(trunc(simTimeR::simTime(available$t_alarm_sec)/60/60))
  lines(x = bplt, y = av_alarm_hour$freq/scale[2], col = "red", lwd = 2)
  points(x = bplt, y = av_alarm_hour$freq/scale[2], col = "red", lwd = 2)
}

# ==============================================================================
# Schedule to availability matrix
wdays_english <- weekdays_abbr()$english
wdays_german <- weekdays_abbr()$german
# Genauigkeit maximal Minuten. Sonst: Stunden
t0_date = lubridate::ymd_hms("2014-01-01 00:00:00")
t1_date = lubridate::ymd_hms("2014-12-31 24:00:00")
interval <- lubridate::interval(t0_date, t1_date) # lubridate
range <- interval/lubridate::dminutes() # numeric dseconds(), dminutes(), dhours()
t_int <- seq(0, range, 1)

all_vehicle_availability <- c() # data.frame(t_alarm_sec = as.integer(0))
for (n_vehicle in 1:nrow(hist_vehicles)) { # nrow(hist_vehicles)
  print(n_vehicle)
  # Days of the year after schedule:
  days_of_year <- c((hist_vehicles$shift_from_simdate[n_vehicle]):(hist_vehicles$shift_to_simdate[n_vehicle] - 1)) # first days = 0!
  all_dates_of_year <- t0_date + lubridate::days(days_of_year)
  available_weekdays <- unlist(strsplit(hist_vehicles$shift_weekday[1], ", ", fixed = TRUE)) # See SCN is_on_duty
  dates_of_year <- all_dates_of_year[lubridate::wday(all_dates_of_year, label = T, abbr = T) %in% g2e(available_weekdays)]
  
  scale <- 3600 # 60=minutes, 3600=hours
  time_of_day_hours <- c((
    hist_vehicles$shift_from_simtime[n_vehicle]/scale):((hist_vehicles$shift_to_simtime[n_vehicle] - 1)/scale)) # lubridate::interval?
  
  all_times_of_year <- dates_of_year[1] + lubridate::hours(time_of_day_hours)
  for (i in 2:length(dates_of_year)) {
    all_times_of_year <- c(all_times_of_year, dates_of_year[i] + lubridate::hours(time_of_day_hours))
  }
  # vehicle_availability <- data.frame(t_alarm_sec = all_times_of_year)
  vehicle_availability <- as.integer(difftime(all_times_of_year,t0_date, units="secs"))
  all_vehicle_availability <- c(all_vehicle_availability, vehicle_availability)
}
df_all_vehicle_availability <- data.frame(t_alarm_sec = all_vehicle_availability)

# ==============================================================================
# Plots

hist_t_alarm_hour(missions[complete.cases(missions),])
# hist_t_alarm_day(missions[complete.cases(missions),])

hist_t_alarm_hour(df_all_vehicle_availability, scale = 365)

compare_t_alarm_hour(missions[complete.cases(missions),], df_all_vehicle_availability, c(365,nrow(hist_vehicles)*365))