#* Functions for plotting hurricane impacts modelled in Impact.jl
#*------------------------------------------------------------------------------

#=
Plots histograms of the number of line failures for given values of γ 
=#
function plothist_tl_failures(
        γ_values::Array{Float64,1}, 
        N_runs::Int64, 
        rstep = 20; # number of realizations that were saved together
        datapath = "", # path to folder containing different γ scans
        plotpath = "" # where to save the plot
    )
    
    N_γ = length(γ_values) # number of given γ values
    mean_Nf_arr = zeros(N_γ) # array for mean number of line failures
    σ_arr = zeros(N_γ) # array for standard deviations

    ### Go through γ values and plot the respective histogram
    for i in 1:N_γ
        mean_Nf, σ = plothist_tl_failures(
            γ_values[i], N_runs, rstep, datapath=datapath, plotpath=plotpath
        )
        mean_Nf_arr[i] = mean_Nf # save mean number of failures
        σ_arr[i] = σ # save standard deviation
    end

    return mean_Nf_arr, σ_arr
end

#=
Plots a histogram of the number of line failures for a given γ value
=#
function plothist_tl_failures(
        γ::Float64,
        N_runs::Int64,
        rstep = 20; # number of realizations that were saved together
        datapath::String, # path to folder containing data files
        plotpath::String # where to save the plot
    )

    Nf = zeros(Int64, N_runs) # array for number of failures
    
    ### Read data and save the number of failures that occured in each run
    for i in rstep:rstep:N_runs
        res = load(datapath * "gamma$γ/" * "results_$γ" * "_$i.jld2", "result")
        for j in 1:rstep
            Nf[i-rstep+j] = length(res[j])
        end
    end

    ### Plot a historgram of the number of failures
    N_bins = maximum(Nf) - minimum(Nf) # number of bins
    n, bins, patches = plt[:hist](Nf, bins=N_bins, edgecolor="k")
    
    ### Calculate mean number of failures and standard deviation
    mean_Nf = round(mean!([1.], Nf)[1], digits=1)
    σ = round(std(Nf, mean=mean_Nf), digits=1)

    ### Legend with γ, mean number of failures and standard deviation
    mlines = pyimport("matplotlib.lines")
    no_marker = mlines.Line2D([], [], ls="None")
    plt.legend(
        [no_marker for i in 1:4], 
        [L"γ = " * "$γ", 
        L"N_r = " * "$N_runs", 
        L"\bar{N}_\mathrm{f} \approx" * "$mean_Nf", 
        L"\sigma_\mathrm{f} \approx" * "$σ"], 
        handletextpad = 0, handlelength = 0
    )

    ### Axes labels
    plt.xlabel("Number of line failures " * L"N_\mathrm{f}")
    plt.ylabel("Frequency")

    ### Save and close figure
    figpath = plotpath * "Texas_hist_tl_failures_gamma$γ.png"
    plt.savefig(figpath, bbox_inches="tight")
    plt.clf()

    return mean_Nf, σ
end

#*------------------------------------------------------------------------------

#=
Plots all frames contained in the provided wind data together with the current power grid. A hurricane may damage the grid according to the data in tl_failures and the corresponding lines are removed from variable_ndd.
=#
function plot_pg_impact(
        network_data::Dict{String,<:Any},
        tl_failures::Array{Tuple{String,Int64}}, # time series of tl failures
        settings = Dict{String,Any}(); # plot settings
        windfile = "", # .nc file containing wind data
        figpath = ""
    )
    
    variable_ndd = deepcopy(network_data) # NDD that is changed during storm
    
    ### Get frames in which transmission lines failed
    failure_frames = unique([tl_failures[i][2] for i in 1:length(tl_failures)])
    N_frames, max_ws = get_winddata(windfile) # get number of frames
    for f in 1:N_frames # plot every time step
        if f in failure_frames # line failure occured
            ### Disable all lines that failed in this frame
            destroyed_tl = [tlf[1] for tlf in tl_failures if tlf[2] == f]
            destroy_tl!(variable_ndd, destroyed_tl)
        end

        figpath_f = figpath * "_frame$f.png"
        plot_pg_map(
            variable_ndd, settings, wind=(windfile, f), figpath=figpath_f
        )
    end

    return variable_ndd
end