From ccc90cf276321af30c189ed82740da590f705afd Mon Sep 17 00:00:00 2001 From: pertft Date: Tue, 1 Oct 2024 11:38:18 +0200 Subject: [PATCH] calc normfactor with ifm_model --- src/ifm.jl | 70 +++++++++++++++++++++++++++++++++++++++++++------ src/io.jl | 2 +- src/local_db.jl | 7 +++-- 3 files changed, 68 insertions(+), 11 deletions(-) diff --git a/src/ifm.jl b/src/ifm.jl index 943bcb5..a68c065 100644 --- a/src/ifm.jl +++ b/src/ifm.jl @@ -61,8 +61,11 @@ mutable struct TwoStateIfmHandler{P <: AbstractTwoStateIfmPredictor, prev_t::Union{TuLiPa.ProbTime, Nothing} ndays_forecast_used::Int + scen_start::Any + scen_stop::Any + function TwoStateIfmHandler(predictor, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) @assert ndays_forecast >= 0 @assert ndays_pred >= ndays_forecast @assert ndays_obs > 0 @@ -78,7 +81,7 @@ mutable struct TwoStateIfmHandler{P <: AbstractTwoStateIfmPredictor, T3 = typeof(hist_Lday) return new{P, U, T1, T2, T3}(predictor, updater, basin_area, m3s_per_mm, hist_P, hist_T, hist_Lday, ndays_pred, ndays_obs, ndays_forecast, data_pred, - data_obs, data_forecast, nothing, 0) + data_obs, data_forecast, nothing, 0, scen_start, scen_stop) end end @@ -259,10 +262,10 @@ struct TwoStateBucketIfm{H} <: AbstractInflowModel handler::H function TwoStateBucketIfm(id, model_params, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) predictor = TwoStateBucketIfmPredictor(model_params) handler = TwoStateIfmHandler(predictor, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) return new{typeof(handler)}(id, handler) end end @@ -312,10 +315,10 @@ struct TwoStateNeuralODEIfm{H} <: AbstractInflowModel handler::H function TwoStateNeuralODEIfm(id, model_params, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) predictor = TwoStateNeuralODEIfmPredictor(model_params) handler = TwoStateIfmHandler(predictor, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) return new{typeof(handler)}(id, handler) end end @@ -400,9 +403,17 @@ function common_includeTwoStateIfm!(Constructor, toplevel::Dict, lowlevel::Dict, data_forecast = nothing ndays_forecast = 0 + + + periodkey = TuLiPa.Id(TuLiPa.TIMEPERIOD_CONCEPT, "ScenarioTimePeriod") + period = lowlevel[periodkey] + scen_start = period["Start"] + scen_stop = period["Stop"] + + id = TuLiPa.getobjkey(elkey) toplevel[id] = Constructor(id, model_params, updater, basin_area, hist_P, hist_T, hist_Lday, - ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast) + ndays_pred, ndays_obs, ndays_forecast, data_obs, data_forecast, scen_start, scen_stop) return (true, deps) end @@ -459,6 +470,42 @@ function save_ifm_Q(db, inflow_name, stepnr, Q) end end +function calculate_normalize_factor(ifm_model) + + S0 = 0 + G0 = 0 + itp_Lday = ifm_model.handler.hist_Lday + itp_P = ifm_model.handler.hist_P + itp_T = ifm_model.handler.hist_T + + + days = Dates.value( Day(ifm_model.handler.scen_stop - ifm_model.handler.scen_start)) + + timepoints = collect((1: days)) + + P = Vector{Float64}([i for i in timepoints]) + T = Vector{Float64}([i for i in timepoints]) + Lday =Vector{Float64}([i for i in timepoints]) + for i in timepoints + start = ifm_model.handler.scen_start + Day(i - 1) + P[i] = TuLiPa.getweightedaverage(itp_P, start, JulES.ONEDAY_MS_TIMEDELTA) + T[i] = TuLiPa.getweightedaverage(itp_T, start, JulES.ONEDAY_MS_TIMEDELTA) + Lday[i] = TuLiPa.getweightedaverage(itp_Lday, start, JulES.ONEDAY_MS_TIMEDELTA) + end + + itp_method = JulES.SteffenMonotonicInterpolation() + itp_P = JulES.interpolate(timepoints, P, itp_method) + itp_T = JulES.interpolate(timepoints, T, itp_method) + itp_Lday = JulES.interpolate(timepoints, Lday, itp_method) + + res = JulES.predict(ifm_model.handler.predictor, S0, G0, itp_Lday, itp_P, itp_T, timepoints); + (Q, _) = res + Q = Float64.(Q) + Q .= Q .* ifm_model.handler.m3s_per_mm + + return 1 / sum(Q) +end + """ Sequentially solve inflow models stored locally. Each inflow model is solved for each scenario. @@ -485,7 +532,14 @@ function solve_ifm(t, stepnr) db.ifm_output[inflow_name] = Dict() end inflow_model = db.ifm[inflow_name] - normalize_factor = normfactors[inflow_name] + + + #if haskey(normfactors, inflow_name) == false + normalize_factor = calculate_normalize_factor(inflow_model) + # TODO: bruk info om hvor lang hist perioden er, start ett år før det, regn ut snittet av døgnverdiene, 1/snitet = normalizefactor + #else + #normalize_factor = normfactors[inflow_name] + #end u0 = estimate_u0(inflow_model, t) diff --git a/src/io.jl b/src/io.jl index e2a18a8..2420099 100644 --- a/src/io.jl +++ b/src/io.jl @@ -148,7 +148,7 @@ end get_iprogtype(input::DefaultJulESInput) = get(input.dataset, "iprogtype", "direct") has_ifm_results(input::DefaultJulESInput) = get_iprogtype(input) != "direct" -get_ifm_normfactors(input::DefaultJulESInput) = get(input.dataset, "ifm_normfactors", Dict{String, Float64}()) +#get_ifm_normfactors(input::DefaultJulESInput) = get(input.dataset, "ifm_normfactors", Dict{String, Float64}()) get_ifm_elements(input::DefaultJulESInput) = get(input.dataset, "ifm_elements", JulES.TuLiPa.DataElement[]) function get_ifm_names(input::DefaultJulESInput) diff --git a/src/local_db.jl b/src/local_db.jl index efba80a..431d8fb 100644 --- a/src/local_db.jl +++ b/src/local_db.jl @@ -76,6 +76,9 @@ mutable struct LocalDB ifm_output::Dict{String, Dict{ScenarioIx, Tuple{Vector{DateTime}, Vector{Float64}}}} ifm_derived::Dict{String, Dict{ScenarioIx, Tuple{Vector{DateTime}, Vector{Float64}}}} + ifm_normfactors::Dict{String, Float64} + + function LocalDB() return new( -1, @@ -115,7 +118,7 @@ mutable struct LocalDB Dict(), # ifm_output Dict(), # ifm_weighted - + Dict(), ) end end @@ -150,7 +153,7 @@ get_ifm_output(db::LocalDB) = db.ifm_output get_ifm_derived(db::LocalDB) = db.ifm_derived get_ifm_weights(db::LocalDB) = get_ifm_weights(get_input(db)) -get_ifm_normfactors(db::LocalDB) = get_ifm_normfactors(get_input(db)) +get_ifm_normfactors(db::LocalDB) = db.ifm_normfactors #get_ifm_normfactors(get_input(db)) get_ifm_elements(db::LocalDB) = get_ifm_elements(get_input(db)) get_cores(db::LocalDB) = get_cores(get_input(db))