# Copyright (c) 2018 Kelvin Say
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in all
# copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
# SOFTWARE.

library(data.table)

# ---------------------------------------------------------------------------------------------------------------------------
# Determine the overall operation of a given household across each yearly time series and also prepare the technical
# yearly summary of the system's operation
# ---------------------------------------------------------------------------------------------------------------------------
process_MeterLoad_Years <- function(simulation_info, timeline, baseload, solar, battery) {

  # Set up calculation constants
  time_interval <- timeline$Info$'Time Step (s)'  #as.numeric(difftime(timeline$Data$Time[2], timeline$Data$Time[1], units="secs"))
  convert_W_to_Wh <- time_interval / 3600
  convert_Wh_to_W <- 1 / convert_W_to_Wh
  num_years <- timeline$Info$'Operational Years'
  n_samples <- length(timeline$Data$Time)
  n_samples_per_day <- n_samples / 365

  # Configure the annual meter summary elements
  annual_meter_summary <- data.table('Load (kWh)' = numeric(num_years),
                                     'Imported (kWh)' = numeric(num_years),
                                     'Exported (kWh)' = numeric(num_years),
                                     'Generated (kWh)' = numeric(num_years),
                                     'PVUsed (kWh)' = numeric(num_years),
                                     'B.Charged (kWh)' = numeric(num_years),
                                     'B.Discharged (kWh)' = numeric(num_years),
                                     'B.InvLosses (kWh)' = numeric(num_years),
                                     'B.CapLosses (Wh)' = numeric(num_years),
                                     'B.SoC.Initial (Wh)' = numeric(num_years),
                                     'B.SoC.Final (Wh)' = numeric(num_years))

  # Configure the daily operating parameters to store
  daily_operation <- data.table('B.Soc.Initial (Wh)' = numeric(num_years * 365))

  # Create dummy meter_operation if not required
  if (!g.save_detailed) {
    meter_operation <- data.table()
  }

  # Process each year of operation to produce a summary table
  for (i in 1:timeline$Info$'Operational Years') {

    # Remember the last state of the battery from year to year
    if (i == 1) {
      battery$Initial$'SoC (Wh)' <- 0
    } else {
      battery$Initial$'SoC (Wh)' <- last(pv_battery_operation$'B.SoC (Wh)')
    }

    # Attenuate the solar insolation with the expected pannel performance degradation for each nominal capacity
    if (solar$Info$'Nominal Capacity (Wp)' > 0 ) {
      solar_degradation_start <- (1-(i-1)*(1-solar$Info$'25 Year Capacity (%)')/25)
      solar_degradation_end <- (1-i*(1-solar$Info$'25 Year Capacity (%)')/25)
      solar_degradation <- solar_degradation_start + (solar_degradation_end - solar_degradation_start)/n_samples*0:(n_samples-1)
      solar$Data$'Generation (W)' <- solar$Info$'Nominal Capacity (Wp)' * solar$Data$'Insolation Reference (W/Wp)' * solar_degradation # xxx - Could technically precalculate this earlier and pass as an input
    } else {
      solar$Data$'Generation (W)' <- NULL
    }

    # Obtain the battery capacity decay curve for this year
    if (battery$Info$'Rated Capacity (Wh)' > 0) {
      capacity_start_Wh <- battery$Info$'Rated Capacity (Wh)'*(1-(i-1)*(1-battery$Info$'10 Year Capacity (%)')/10)
      capacity_end_Wh <- battery$Info$'Rated Capacity (Wh)'*(1-i*(1-battery$Info$'10 Year Capacity (%)')/10)
      battery$Data$'Capacity (Wh)' <- capacity_start_Wh + (capacity_end_Wh - capacity_start_Wh)/n_samples*0:(n_samples-1) # xxx - Could technically precalculate this earlier and pass as an input
    } else {
      battery$Data$'Capacity (Wh)' <- NULL
    }

    # Take the baseload, solar and battery data, and operate the customer loads to determine the meter operation
    pv_battery_operation <- process_PV_Battery_Load(timeline, baseload, solar, battery)

    # Only store the detailed data if required
    if (g.save_detailed) {
      meter_operation <- data.table('Time' = timeline$Data$Time,
                                    'Load (W)' = baseload$Data$'Load (W)',
                                    'PV (W)' = solar$Data$'Generation (W)',
                                    'NetLoad (W)' = pv_battery_operation$'NetLoad (W)',
                                    'PVUsed (W)' = pv_battery_operation$'PVUsed (W)',
                                    'B.Charging (W)' = pv_battery_operation$'B.Charging (W)',
                                    'B.Discharging (W)' = pv_battery_operation$'B.Discharging (W)',
                                    'B.InvLosses (W)' = pv_battery_operation$'B.InvLosses (W)',
                                    'B.CapLosses (Wh)' = pv_battery_operation$'B.CapLosses (Wh)',
                                    'B.SoC (Wh)' = pv_battery_operation$'B.SoC (Wh)',
                                    'Imported (W)' = pv_battery_operation$'Imported (W)',
                                    'Exported (W)' = pv_battery_operation$'Exported (W)')
    }

    # Store away the daily condition of the system
    daily_operation$'B.Soc.Initial (Wh)'[1+(i-1)*365] <- battery$Initial$'SoC (Wh)'
    daily_operation$'B.Soc.Initial (Wh)'[(2:365)+(i-1)*365] <- pv_battery_operation$'B.SoC (Wh)'[n_samples_per_day*(1:364)]

    # Summarise the annual performance figures for this house
    annual_meter_summary$'Load (kWh)'[i] <- sum(baseload$Data$'Load (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'Imported (kWh)'[i] <- sum(pv_battery_operation$'Imported (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'Exported (kWh)'[i] <- sum(pv_battery_operation$'Exported (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'Generated (kWh)'[i] <- sum(solar$Data$'Generation (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'PVUsed (kWh)'[i] <- sum(pv_battery_operation$'PVUsed (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'B.Charged (kWh)'[i] <- sum(pv_battery_operation$'B.Charging (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'B.Discharged (kWh)'[i] <- sum(pv_battery_operation$'B.Discharging (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'B.InvLosses (kWh)'[i] <- sum(pv_battery_operation$'B.InvLosses (W)') * convert_W_to_Wh / 1000
    annual_meter_summary$'B.CapLosses (Wh)'[i] <- sum(pv_battery_operation$'B.CapLosses (Wh)')
    annual_meter_summary$'B.SoC.Initial (Wh)'[i] <- battery$Initial$'SoC (Wh)'
    annual_meter_summary$'B.SoC.Final (Wh)'[i] <- last(pv_battery_operation$'B.SoC (Wh)')
  }
  return(list('Operation' = meter_operation,
              'Daily' = daily_operation,
              'Summary' = annual_meter_summary))
}

# ---------------------------------------------------------------------------------------------------------------------------
# Operate the battery using the input PV and baseload profiles to determine the netload of the household
# ---------------------------------------------------------------------------------------------------------------------------
process_PV_Battery_Load <- function(timeline, baseload, solar, battery) {

  # setup the battery's configuration
  time_interval <- timeline$Info$'Time Step (s)' # as.numeric(difftime(timeline$Data$Time[2], timeline$Data$Time[1], units="secs"))
  convert_W_to_Wh <- time_interval / 3600
  convert_Wh_to_W <- 1 / convert_W_to_Wh
  charge_limit_Wh <- battery$Info$'Input Limit (W)' * convert_W_to_Wh
  discharge_limit_Wh <- battery$Info$'Output Limit (W)' * convert_W_to_Wh

  num_steps <- length(timeline$Data$Time)
  if (solar$Info$'Nominal Capacity (Wp)' > 0) {
    v_netload_W <- baseload$Data$'Load (W)' - solar$Data$'Generation (W)'
  } else {
    v_netload_W <- baseload$Data$'Load (W)'
  }
  v_state_of_charge_Wh <- rep(0, num_steps)
  v_capacity_losses_Wh <- rep(0, num_steps)
  v_losses_W <- rep(0, num_steps)
  v_charging_W <- rep(0, num_steps)
  v_discharging_W <- rep(0, num_steps)
  v_PV_used_W <- rep(0, num_steps)

  # Battery actually has capacity to work with
  if (battery$Info$'Rated Capacity (Wh)' > 0) {
    for (i in 1:num_steps) {
      netload_W <- v_netload_W[i]
      if (i == 1) {
        prev_state_of_charge_Wh <- battery$Initial$'SoC (Wh)'
      } else {
        prev_state_of_charge_Wh <- v_state_of_charge_Wh[i-1]
      }
      curr_state_of_charge_Wh <- prev_state_of_charge_Wh

      # ---------------------- capacity loss cycle ---------------------
      # Skim off any 'state of charge' greater than the reduced battery capacity, and put into the 'capacity losses'
      max_state_of_charge_Wh <- battery$Data$'Capacity (Wh)'[i]
      if (prev_state_of_charge_Wh > max_state_of_charge_Wh) {
        v_capacity_losses_Wh[i] <- prev_state_of_charge_Wh - max_state_of_charge_Wh
        prev_state_of_charge_Wh <- max_state_of_charge_Wh
        curr_state_of_charge_Wh <- prev_state_of_charge_Wh
      }

      # ------------------------- charge cycle -------------------------
      if (netload_W < 0) {

        # Battery is not full
        available_capacity_Wh <- ( (max_state_of_charge_Wh - prev_state_of_charge_Wh) / (1 - battery$Info$'Input Loss (%)') ) # factor input losses
        if (available_capacity_Wh > 0) {
          committed_energy_Wh <- 0
          excess_energy_Wh <- -1 * netload_W * convert_W_to_Wh

          # Determine the maximum amount of energy the battery could accept
          max_charge_energy_Wh <- min(charge_limit_Wh, available_capacity_Wh)

          # Self-generation exceeds available battery capacity - store what you can, export the rest
          if (excess_energy_Wh > max_charge_energy_Wh) {
            committed_energy_Wh <- max_charge_energy_Wh * (1 - battery$Info$'Input Loss (%)')
            lost_energy_Wh <- max_charge_energy_Wh - committed_energy_Wh
            curr_state_of_charge_Wh <- prev_state_of_charge_Wh + committed_energy_Wh

            excess_energy_Wh <- excess_energy_Wh - max_charge_energy_Wh
            v_netload_W[i] <- -1 * excess_energy_Wh * convert_Wh_to_W
            v_charging_W[i] <- max_charge_energy_Wh * convert_Wh_to_W
            v_losses_W[i] <- lost_energy_Wh * convert_Wh_to_W
          }
          # Self-generation fits into the battery - store all of it
          else {
            committed_energy_Wh <- excess_energy_Wh * (1 - battery$Info$'Input Loss (%)')
            lost_energy_Wh <- excess_energy_Wh - committed_energy_Wh
            curr_state_of_charge_Wh <- prev_state_of_charge_Wh + committed_energy_Wh

            v_netload_W[i] <- 0
            v_charging_W[i] <- excess_energy_Wh * convert_Wh_to_W
            v_losses_W[i] <- lost_energy_Wh * convert_Wh_to_W
          }
        }
      }
      # ----------------------- discharge cycle ------------------------
      else {

        # Battery is not empty
        available_capacity_Wh <- prev_state_of_charge_Wh
        if (available_capacity_Wh > 0) {
          extracted_energy_Wh <- 0 # Amount of energy that will be extracted from the battery
          required_energy_Wh <- (netload_W * convert_W_to_Wh) / (1 - battery$Info$'Output Loss (%)')

          # Determine the maximum amount of energy the battery could send out
          max_discharge_energy_Wh <- min(discharge_limit_Wh, available_capacity_Wh)

          # Battery is able to meet the energy requirements - discharge the battery, import nothing
          if (max_discharge_energy_Wh > required_energy_Wh) {
            extracted_energy_Wh <- required_energy_Wh * (1 - battery$Info$'Output Loss (%)') # factor output losses
            lost_energy_Wh <- required_energy_Wh - extracted_energy_Wh
            curr_state_of_charge_Wh <- prev_state_of_charge_Wh - required_energy_Wh

            v_netload_W[i] <- 0
            v_discharging_W[i] <- extracted_energy_Wh * convert_Wh_to_W
            v_losses_W[i] <- lost_energy_Wh * convert_Wh_to_W
          }
          # Battery is unable to meet the energy requirements - discharge the battery, import the rest
          else {
            extracted_energy_Wh <- max_discharge_energy_Wh * (1 - battery$Info$'Output Loss (%)')
            lost_energy_Wh <- max_discharge_energy_Wh - extracted_energy_Wh
            curr_state_of_charge_Wh <- prev_state_of_charge_Wh - max_discharge_energy_Wh

            required_energy_Wh <- (netload_W * convert_W_to_Wh) - extracted_energy_Wh
            v_netload_W[i] <- required_energy_Wh * convert_Wh_to_W
            v_discharging_W[i] <- extracted_energy_Wh * convert_Wh_to_W
            v_losses_W[i] <- lost_energy_Wh * convert_Wh_to_W
          }
        }
      }
      v_state_of_charge_Wh[i] <- curr_state_of_charge_Wh
    }
  } # END IF: Rated Capacity > 0

  # Determine how much of the PV energy is directly used
  if (solar$Info$'Nominal Capacity (Wp)' > 0) {
    temp_mask <- sapply(baseload$Data$'Load (W)' - solar$Data$'Generation (W)', function(x) if(x > 0) 1 else 0)
    v_PV_used_W <- (temp_mask * solar$Data$'Generation (W)') + ((1-temp_mask) * baseload$Data$'Load (W)')
  }

  # Determine the amount of imported and exported energy
  v_imported_W <- sapply(v_netload_W, function(x) if(x > 0) x else 0)
  v_exported_W <- sapply(v_netload_W, function(x) if(x < 0) -x else 0)

  return(data.table('NetLoad (W)' = v_netload_W,
                    'PVUsed (W)' = v_PV_used_W,
                    'B.Charging (W)' = v_charging_W,
                    'B.Discharging (W)' = v_discharging_W,
                    'B.InvLosses (W)' = v_losses_W,
                    'B.CapLosses (Wh)' = v_capacity_losses_Wh,
                    'B.SoC (Wh)' = v_state_of_charge_Wh,
                    'Imported (W)' = v_imported_W,
                    'Exported (W)' = v_exported_W))
}

# ---------------------------------------------------------------------------------------------------------------------------
#
# ---------------------------------------------------------------------------------------------------------------------------
# technical_summary <- process_Summarise_Technical(technical_data)
process_Summarise_Technical <- function(technical_data) {

  # Setup constants
  max_elements_pv <- length(technical_data$'Operational'[,1])
  max_elements_battery <- length(technical_data$'Operational'[1,])
  op_years <- nrow(technical_data$Operational[1,1][[1]]$'Summary')

  # Setup the Energy summary tables
  pv_sizes <- rep(0, max_elements_pv)
  battery_sizes <- rep(-1, max_elements_battery)

  yearly.load <- rep(-1, op_years)
  yearly.imported <- array(rep(0, max_elements_pv*max_elements_battery*op_years), dim = c(max_elements_pv, max_elements_battery, op_years))
  yearly.exported <- array(rep(0, max_elements_pv*max_elements_battery*op_years), dim = c(max_elements_pv, max_elements_battery, op_years))
  yearly.generated <- array(rep(0, max_elements_pv*op_years), dim=c(max_elements_pv, op_years))
  yearly.pv_used <- array(rep(0, max_elements_pv*op_years), dim=c(max_elements_pv, op_years))
  yearly.batt_discharged <- array(rep(0, max_elements_pv*max_elements_battery*op_years), dim = c(max_elements_pv, max_elements_battery, op_years))
  yearly.batt_loss <- array(rep(0, max_elements_pv*max_elements_battery*op_years), dim = c(max_elements_pv, max_elements_battery, op_years))
  yearly.batt_soc_offset <- array(rep(0, max_elements_pv*max_elements_battery*op_years), dim = c(max_elements_pv, max_elements_battery, op_years))

  total.load <- -1
  total.imported <- array(rep(0, max_elements_pv*max_elements_battery), dim = c(max_elements_pv, max_elements_battery))
  total.exported <- array(rep(0, max_elements_pv*max_elements_battery), dim = c(max_elements_pv, max_elements_battery))
  total.generated <- rep(0, max_elements_pv)
  total.pv_used <- rep(0, max_elements_pv)
  total.batt_discharged <- array(rep(0, max_elements_pv*max_elements_battery), dim = c(max_elements_pv, max_elements_battery))
  total.batt_loss <- array(rep(0, max_elements_pv*max_elements_battery), dim = c(max_elements_pv, max_elements_battery))
  total.batt_soc_offset <- array(rep(0, max_elements_pv*max_elements_battery), dim = c(max_elements_pv, max_elements_battery))

  # Compute the summary fields for both the yearly and total data elements
  yearly.load <- technical_data$Operational[1,1][[1]]$Summary$'Load'
  total.load <- sum(yearly.load)
  for (j in 1:max_elements_battery) {
    battery_sizes[j] <- technical_data$Operational[1,j][[1]]$'Battery Capacity (Wh)'
  }

  for (i in 1:max_elements_pv) {
    pv_sizes[i] <- technical_data$Operational[i,1][[1]]$'PV Capacity (Wp)'
    for (j in 1:max_elements_battery) {
      yearly.imported[i,j,] <- technical_data$Operational[i,j][[1]]$Summary$'Imported (kWh)'
      yearly.exported[i,j,] <- technical_data$Operational[i,j][[1]]$Summary$'Exported (kWh)'
      yearly.batt_discharged[i,j,] <- technical_data$Operational[i,j][[1]]$Summary$'B.Discharged (kWh)'
      yearly.batt_loss[i,j,] <- technical_data$Operational[i,j][[1]]$Summary$'B.InvLosses (kWh)' + (technical_data$Operational[i,j][[1]]$Summary$'B.CapLosses (Wh)' + technical_data$Operational[i,j][[1]]$Summary$'B.SoC.Final (Wh)' - technical_data$Operational[i,j][[1]]$Summary$'B.SoC.Initial (Wh)') / 1000

      total.imported[i,j] <- sum(yearly.imported[i,j,])
      total.exported[i,j] <- sum(yearly.exported[i,j,])
      total.batt_discharged[i,j] <- sum(yearly.batt_discharged[i,j,])
      total.batt_loss[i,j] <- sum(yearly.batt_loss[i,j,])
    }
    yearly.generated[i,] <- technical_data$Operational[i,1][[1]]$Summary$'Generated'
    yearly.pv_used[i,] <- technical_data$Operational[i,1][[1]]$Summary$'PVUsed'

    total.generated[i] <- sum(yearly.generated[i,])
    total.pv_used[i] <- sum(yearly.pv_used[i,])
  }

  yearly <- list('Load (kWh)' = yearly.load,
                 'Imported (kWh)' = yearly.imported,
                 'Exported (kWh)' = yearly.exported,
                 'Generated (kWh)' = yearly.generated,
                 'Directly Used (kWh)' = yearly.pv_used,
                 'Battery Discharged (kWh)' = yearly.batt_discharged,
                 'Battery Losses (kWh)' = yearly.batt_loss,
                 'Battery Charged (kWh)' = yearly.batt_discharged + yearly.batt_loss,
                 'Self-Generation (%)' = (crossfill_2D_3D_matrix(yearly.pv_used, max_elements_battery) + yearly.batt_discharged + yearly.batt_loss) / crossfill_2D_3D_matrix(yearly.generated, max_elements_battery),
                 'Self-Sufficiency (%)' = (crossfill_2D_3D_matrix(yearly.pv_used, max_elements_battery) + yearly.batt_discharged) / crossfill_1D_3D_matrix(yearly.load, max_elements_pv, max_elements_battery),
                 'Dependence (%)' = yearly.imported / crossfill_1D_3D_matrix(yearly.load, max_elements_pv, max_elements_battery))

  total <- list('Load (kWh)' = total.load,
                'Imported (kWh)' = total.imported,
                'Exported (kWh)' = total.exported,
                'Generated (kWh)' = total.generated,
                'Directly Used (kWh)' = total.pv_used,
                'Battery Discharged (kWh)' = total.batt_discharged,
                'Battery Losses (kWh)' = total.batt_loss,
                'Battery Charged (kWh)' = total.batt_discharged + total.batt_loss,
                'Self-Generation (%)' = (total.pv_used + total.batt_discharged + total.batt_loss) / total.generated,
                'Self-Sufficiency (%)' = (total.pv_used + total.batt_discharged) / total.load,
                'Dependence (%)' = total.imported / total.load)

  return(list('Inputs' = technical_data$Inputs,
              'PV Capacity (Wp)' = pv_sizes,
              'Battery Capacity (Wh)' = battery_sizes,
              'Total' = total,
              'Yearly' = yearly))
}