########## 
##########
# Pulling in all data from simulations, getting network metrics, and writing
# the results to a file for analysis
##########
##########
# AUTHOR: Cole B. Brookson
# DATE OF CREATION: 2021-11-20
##########
##########

# set up =======================================================================
import Pkg
Pkg.instantiate()
using DelimitedFiles

include("./05_network_metrics.jl")
include("./04_eqn_component_functions.jl")
include("./01_species_list_creation.jl")

# set all positional values 
organism_id_pos = 1
organism_type_pos = 2
body_size_pos = 3
lower_prey_limit_pos = 4
upper_prey_limit_pos = 5
habitat_midpoint_pos = 6
habitat_lower_limit_pos = 7
habitat_upper_limit_pos = 8
nutritional_value_pos = 9
biomass_pos = 10
network_id_pos = 11
alive_pos = 12 

# read in data =================================================================

# read in data for all the species lists
all_species_lists = readdlm("./data/intermediate-files/species_lists_for_simulations.csv", ',', Any)
network_data = readdlm("./data/intermediate-files/all_network_data_for_simulations.csv", ',', Any)

"""
    add_metrics(sims_data, network_num, sim_num, perturb_type)

This function takes in the information about a given simulation, and calculates 
the network metrics at every time point 
"""
function add_metrics(sims_data, 
                    network_num,
                    sim_num, 
                    perturb_type,
                    species_list)
    # for each row in the sims_data, calculate the network metrics and put them
    # in the data object
    data = Array(hcat(Array{Int64}(undef, size(sims_data, 1)),
                        Array{String}(undef, size(sims_data, 1)),
                        Array{Int64}(undef, size(sims_data, 1)),
                        Array{String}(undef, size(sims_data, 1)),
                        Array{Float64}(undef, (size(sims_data, 1), 6))))
    # get the network type 
    network_type = network_data[(network_data[:,1] .== network_num),:][1,3]

    # loop through each time point and get the network metrics for that time
    for i in 1:size(sims_data, 1)
        # set biomasses to be the ones at that time step 
        species_list[:,biomass_pos] = sims_data[i, 2:end]
        # get the network
        network = biomass_interaction_matrix(species_list)
        # set the row of the data to be the 
        data[i,:] = hcat(network_num, 
                        network_type,
                        sim_num, 
                        perturb_type, 
                        network_metrics(species_list, network)[1]')
    end 
    
    return data 
end 

# run through every file in each folder 
for perturb in ["pulse", "press", "period"]
    # set the path of the files to process
    path = string("./data/simulation-outputs/",perturb,"-sims/")
    # set empty data structure
    data_pre = Array(hcat(Array{Int64, 1},
                    Array{String, 1},
                    Array{Int64, 1},
                    Array{String, 1},
                    Array{Float64, 1}, 
                    Array{Float64, 1}, 
                    Array{Float64, 1}, 
                    Array{Float64, 1}, 
                    Array{Float64, 1}, 
                    Array{Float64, 1} 
                    ))
    for id in append!([2],collect(range(1000,9000, step = 500)))
        species_list = all_species_lists[all_species_lists[:,network_id_pos] .== id,:]
        for random in collect(1:10)
            # get the file names and read in the data one file at a time 
            file_name = string(perturb,"-",id,"-",random,"-sims.csv")
            path_file = string(path, file_name)
            sims_data = readdlm(path_file,',',Any)
            # fix data
            for row in 1:size(sims_data,1)
                sims_data[row,2] = parse(Float64, sims_data[row,2][2:end])
                sims_data[row,end] = parse(Float64, sims_data[row,end][2:end-1])
            end 
            # use function to get data and append it to the previous one 
            temp_data = add_metrics(sims_data, id, random, perturb, species_list)
            data_pre = cat(data_pre, temp_data, dims = 1)
        end 
    end 
    data_pre = data_pre[2:end,:]
    out_name = string(perturb,"-metrics.csv")
    out_path = string("./data/metric-outputs/")
    writedlm(string(out_path,out_name), data_pre, ",")
end
sims_data_float = convert(Matrix{Float64}, sims_data[:,2:end])
histogram(sims_data_float, bins = :scott)
#find(sims_data_float)

# testing
perturb = "pulse"
network_id = 501
id = 501
random = 1
path = string("./data/simulation-outputs/",perturb,"-sims/")
file_name = string(perturb,"-",id,"-",random,"-sims.csv")
path_file = string(path, file_name)
sims_data = readdlm(path_file,',',Any)
# fix data
for row in 1:size(sims_data,1)
    sims_data[row,2] = parse(Float64, sims_data[row,2][2:end])
    sims_data[row,end] = parse(Float64, sims_data[row,end][2:end-1])
end 
species_list = all_species_lists[all_species_lists[:,network_id_pos] .== network_id,:]
start_species_list = species_list
start_species_list[:,biomass_pos] = sims_data[1, 2:end]
end_species_list = species_list
end_species_list[:,biomass_pos] = sims_data[end, 2:end]
start_network = biomass_interaction_matrix(start_species_list)
end_network = biomass_interaction_matrix(end_species_list)
start_network == end_network

j = 1
i = 3

predator_id = j
prey_id = i
end_species_list[3,biomass_pos] = 0
   # extract the values for habitat, size, and predation limits for ease 
   predator = get_species(start_species_list, predator_id)
   prey = get_species(start_species_list, prey_id)

   pred_hab = predator[habitat_midpoint_pos:habitat_upper_limit_pos]
   prey_hab = prey[habitat_midpoint_pos:habitat_upper_limit_pos]
   prey_size = prey[body_size_pos]
   pred_min = predator[lower_prey_limit_pos]
   pred_max = predator[upper_prey_limit_pos]

   if (predator_id == prey_id) || (predator[organism_type_pos] == "Producer")
       value = 0 # assuming no canibalism and that producers don't predate
   elseif predator[biomass_pos] == 0.0
       value = 0
   elseif prey[biomass_pos] == 0.0
       value = 0
   elseif (((prey_hab[2] >= pred_hab[2] && prey_hab[2] <= pred_hab[3]) ||
       (prey_hab[3] >= pred_hab[2] && prey_hab[3] <= pred_hab[3])) && 
       (prey_size[1] <= pred_max && prey_size[1] >= pred_min))
       value = 1 # assign 1 if within body size and habitat range
   else
       value = 0
   end 

   return value 
start_network = fill(0, (size(start_species_list, 1), size(start_species_list, 1)))
for i in 1:size(start_network)[1] #rows
    for j in 1:size(start_network)[2] #cols
        start_network[i,j] = biomass_interaction(j, i, start_species_list)
    end
end
showall(start_network[:,22])


end_network = fill(0, (size(end_species_list, 1), size(end_species_list, 1)))
for i in 1:size(end_network)[1] #rows
    for j in 1:size(end_network)[2] #cols
        end_network[i,j] = biomass_interaction(j, i, end_species_list)
    end
end
histogram(end_species_list[:,biomass_pos])
end_species_list[end_species_list[:,biomass_pos] .== 0.0, :]

function myshowall(io, x, limit = false) 
    println(io, summary(x), ":")
    Base.print_matrix(IOContext(io, :limit => limit), x)
end

show(stdout, "text/plain", start_network[:,22])