mirror of https://gitlab.com/pamhyr/pamhyr2
1039 lines
36 KiB
Python
1039 lines
36 KiB
Python
# AdisTT.py -- Pamhyr
|
|
# Copyright (C) 2023-2025 INRAE
|
|
#
|
|
# This program is free software: you can redistribute it and/or modify
|
|
# it under the terms of the GNU General Public License as published by
|
|
# the Free Software Foundation, either version 3 of the License, or
|
|
# (at your option) any later version.
|
|
#
|
|
# This program is distributed in the hope that it will be useful,
|
|
# but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
# GNU General Public License for more details.
|
|
#
|
|
# You should have received a copy of the GNU General Public License
|
|
# along with this program. If not, see <https://www.gnu.org/licenses/>.
|
|
|
|
# -*- coding: utf-8 -*-
|
|
|
|
import os
|
|
import logging
|
|
|
|
import numpy as np
|
|
|
|
import shutil
|
|
|
|
from tools import (
|
|
trace, timer, logger_exception,
|
|
timestamp_to_old_pamhyr_date,
|
|
old_pamhyr_date_to_timestamp,
|
|
timestamp_to_old_pamhyr_date_adists
|
|
)
|
|
|
|
from Solver.CommandLine import CommandLineSolver
|
|
|
|
from Model.Results.ResultsAdisTS import Results
|
|
from Model.Results.River.River import River, Reach, Profile
|
|
|
|
from itertools import chain
|
|
|
|
logger = logging.getLogger()
|
|
|
|
|
|
def adistt_file_open(filepath, mode):
|
|
f = open(filepath, mode)
|
|
|
|
if "w" in mode:
|
|
# Write header
|
|
comment = "*"
|
|
if ".ST" in filepath:
|
|
comment = "#"
|
|
|
|
f.write(
|
|
f"{comment} " +
|
|
"This file is generated by PAMHYR, please don't modify\n"
|
|
)
|
|
|
|
return f
|
|
|
|
|
|
class AdisTT(CommandLineSolver):
|
|
_type = "adistt"
|
|
|
|
def __init__(self, name):
|
|
super(AdisTT, self).__init__(name)
|
|
|
|
self._type = "adistt"
|
|
|
|
self._cmd_input = ""
|
|
self._cmd_solver = "@path @input -o @output"
|
|
self._cmd_output = ""
|
|
|
|
@classmethod
|
|
def default_parameters(cls):
|
|
lst = super(AdisTT, cls).default_parameters()
|
|
|
|
lst += [
|
|
("adistt_implicitation_parameter", "0.5"),
|
|
("adistt_timestep_screen", "86400"),
|
|
("adistt_timestep_bin", "3600"),
|
|
("adistt_timestep_csv", "3600"),
|
|
("adistt_timestep_mage", "300"),
|
|
("adistt_initial_concentration", "0"),
|
|
]
|
|
|
|
return lst
|
|
|
|
@classmethod
|
|
def checkers(cls):
|
|
lst = [
|
|
]
|
|
|
|
return lst
|
|
|
|
##########
|
|
# Export #
|
|
##########
|
|
|
|
_alph = list(
|
|
map(
|
|
chr,
|
|
chain(
|
|
range(48, 58), # 0..9
|
|
range(65, 91), # A..Z
|
|
range(97, 123) # a..z
|
|
)
|
|
)
|
|
)
|
|
|
|
_l_alph = len(_alph)
|
|
|
|
_nodes_cnt = 0
|
|
_nodes_names = {}
|
|
_nodes_views = set()
|
|
|
|
def get_reach_name(self, reach):
|
|
index = self._study.river.get_edge_id(reach) + 1
|
|
return f"Reach_{index:03d}"
|
|
|
|
def get_node_name(self, node):
|
|
"""Generate a 3 char name for node
|
|
|
|
Args:
|
|
node: The node
|
|
|
|
Returns:
|
|
A 3 char name string
|
|
"""
|
|
n = node.pamhyr_id
|
|
|
|
if n in self._nodes_names:
|
|
return self._nodes_names[n]
|
|
|
|
name = ""
|
|
|
|
checked_new = False
|
|
while not checked_new:
|
|
self._nodes_cnt += 1
|
|
nc = self._nodes_cnt
|
|
|
|
name = "".join(
|
|
map(
|
|
lambda i: self._alph[i % self._l_alph],
|
|
[
|
|
int(nc / (self._l_alph * self._l_alph)),
|
|
int(nc / self._l_alph),
|
|
nc
|
|
]
|
|
)
|
|
)
|
|
|
|
checked_new = name not in self._nodes_views
|
|
|
|
self._nodes_views.add(name)
|
|
self._nodes_names[n] = name
|
|
|
|
return name
|
|
|
|
def cmd_args(self, study):
|
|
lst = super(AdisTT, self).cmd_args(study)
|
|
return lst
|
|
|
|
def input_param(self):
|
|
name = self._study.name.replace(" ", "_")
|
|
return f"{name}.REP"
|
|
|
|
def log_file(self):
|
|
name = self._study.name.replace(" ", "_")
|
|
return f"{name}.TRA"
|
|
|
|
def _export_ST(self, study, repertory, qlog, name="0"):
|
|
files = []
|
|
|
|
if qlog is not None:
|
|
qlog.put("Export ST file")
|
|
|
|
os.makedirs(os.path.join(repertory, "net"), exist_ok=True)
|
|
|
|
# Write header
|
|
edges = study.river.enable_edges()
|
|
for edge in edges:
|
|
name = self.get_reach_name(edge)
|
|
|
|
with adistt_file_open(
|
|
os.path.join(repertory, "net", f"{name}.ST"),
|
|
"w+"
|
|
) as f:
|
|
files.append(str(os.path.join("net", f"{name}.ST")))
|
|
|
|
cnt_num = 1
|
|
for profile in edge.reach.profiles:
|
|
self._export_ST_profile_header(
|
|
f, files, profile, cnt_num
|
|
)
|
|
cnt_num += 1
|
|
|
|
# Points
|
|
for point in profile.points:
|
|
self._export_ST_point_line(
|
|
f, files, point
|
|
)
|
|
|
|
# Profile last line
|
|
f.write(f" 999.9990 999.9990 999.9990\n")
|
|
|
|
def _export_ST_profile_header(self, wfile, files,
|
|
profile, cnt):
|
|
num = f"{cnt:>6}"
|
|
c1 = f"{profile.code1:>6}"
|
|
c2 = f"{profile.code2:>6}"
|
|
t = f"{len(profile.points):>6}"
|
|
rk = f"{profile.rk:>12f}"[0:12]
|
|
pname = profile.name
|
|
if profile.name == "":
|
|
# Generate name from profile id prefixed with
|
|
# 'p' (and replace space char with '0' char)
|
|
pname = f"p{profile.id:>3}".replace(" ", "0")
|
|
name = f"{pname:<19}"
|
|
|
|
# Profile header line
|
|
# AdisTT only reads the hydraulic geometry. Mage sediment-layer
|
|
# fields are deliberately omitted because AdisTT interprets them
|
|
# with a different grammar.
|
|
wfile.write(f"{num}{c1}{c2}{t} {rk} {pname}\n")
|
|
|
|
def _export_ST_point_line(self, wfile, files, point):
|
|
x = f"{point.x:<12.4f}"[0:12]
|
|
y = f"{point.y:<12.4f}"[0:12]
|
|
z = f"{point.z:<12.4f}"[0:12]
|
|
n = f"{point.name:<3}"
|
|
|
|
# Point line
|
|
wfile.write(f"{x} {y} {z} {n}\n")
|
|
|
|
def _export_NET(self, study, repertory, qlog=None, name="0"):
|
|
|
|
if qlog is not None:
|
|
qlog.put("Export NET file")
|
|
|
|
with adistt_file_open(
|
|
os.path.join(repertory, f"{name}.NET"), "w+") as f:
|
|
edges = study.river.enable_edges()
|
|
|
|
for e in edges:
|
|
name = self.get_reach_name(e)
|
|
id = name
|
|
|
|
n1 = f"{self.get_node_name(e.node1):3}".replace(" ", "x")
|
|
n2 = f"{self.get_node_name(e.node2):3}".replace(" ", "x")
|
|
file = os.path.join("net", name + ".ST")
|
|
|
|
f.write(f"{id} {n1} {n2} {file}\n")
|
|
|
|
def _export_fake_INI(self, study, repertory, qlog=None, name="0"):
|
|
if qlog is not None:
|
|
qlog.put("Export fake INI file")
|
|
|
|
with adistt_file_open(
|
|
os.path.join(
|
|
repertory, "Mage_fin.ini"
|
|
), "w+"
|
|
) as f:
|
|
edges = study.river.enable_edges()
|
|
for i, edge in enumerate(edges):
|
|
lst = list(filter(
|
|
lambda f: f.is_full_defined(),
|
|
edge.frictions.frictions
|
|
))
|
|
rk_min = 9999999.9
|
|
rk_max = -9999999.9
|
|
coeff_min = -1.0
|
|
coeff_max = -1.0
|
|
for s in lst: # TODO optimise ?
|
|
if s.begin_rk > rk_max:
|
|
rk_max = s.begin_rk
|
|
coeff_max = s.begin_strickler
|
|
if s.begin_rk < rk_min:
|
|
rk_min = s.begin_rk
|
|
coeff_min = s.begin_strickler
|
|
if s.end_rk > rk_max:
|
|
rk_max = s.end_rk
|
|
coeff_max = s.end_strickler
|
|
if s.end_rk < rk_min:
|
|
rk_min = s.end_rk
|
|
coeff_min = s.end_strickler
|
|
|
|
print("min max", rk_min, rk_max)
|
|
|
|
def get_stricklers_from_rk(rk, lst):
|
|
print("rk", rk)
|
|
|
|
coeff = None
|
|
if rk > rk_max:
|
|
coeff = coeff_max
|
|
elif rk < rk_min:
|
|
coeff = coeff_min
|
|
else:
|
|
for s in lst:
|
|
if (rk >= s.begin_rk and rk <= s.end_rk or
|
|
rk <= s.begin_rk and rk >= s.end_rk):
|
|
coeff = s.begin_strickler # TODO: inerpolate
|
|
break
|
|
|
|
# TODO interpolation if rk is not in frictons
|
|
|
|
if coeff is None:
|
|
logger.error(
|
|
"Study frictions are not fully defined"
|
|
)
|
|
return None
|
|
|
|
return coeff.minor, coeff.medium
|
|
|
|
for j, profile in enumerate(edge.reach.profiles):
|
|
coef_min, coef_moy = get_stricklers_from_rk(profile.rk,
|
|
lst)
|
|
f.write(
|
|
f" {i+1:3}{j+1:4}{' '*116}{coef_min:9}{coef_moy:9}\n")
|
|
|
|
def _export_REP_additional_lines(self, study, rep_file):
|
|
lines = filter(
|
|
lambda line: line.is_enabled(),
|
|
study.river.rep_lines.lines
|
|
)
|
|
|
|
for line in lines:
|
|
rep_file.write(line.line)
|
|
|
|
def _export_REP(self, study, repertory, mage_rep, files, qlog, name="0"):
|
|
|
|
if qlog is not None:
|
|
qlog.put("Export REP file")
|
|
|
|
mage_path = (mage_rep if os.path.isabs(mage_rep) else os.path.join(
|
|
os.path.abspath(os.path.join(repertory, os.pardir)),
|
|
mage_rep,
|
|
))
|
|
relative_mage_path = os.path.relpath(mage_path, repertory)
|
|
|
|
# Write header
|
|
with adistt_file_open(
|
|
os.path.join(
|
|
repertory, f"{name}.REP"
|
|
), "w+"
|
|
) as f:
|
|
path = os.path.join(relative_mage_path, name)
|
|
f.write(f"NET {path}.NET\n")
|
|
f.write(f"REP {path}.REP\n")
|
|
|
|
for file in files:
|
|
EXT = file.split('.')[1]
|
|
f.write(f"{EXT} {file}\n")
|
|
|
|
self._export_REP_additional_lines(study, f)
|
|
|
|
# fake mage_fin.ini:
|
|
|
|
self._export_fake_INI(study, mage_path,
|
|
qlog, name=name)
|
|
|
|
# path_mage_net = os.path.join(os.path.abspath(
|
|
# os.path.join(repertory, os.pardir)
|
|
# ), os.path.join(mage_rep, "net"))
|
|
# path_adistt_net = os.path.join(repertory, "net")
|
|
|
|
# if os.path.exists(path_mage_net):
|
|
# shutil.copytree(path_mage_net,
|
|
# path_adistt_net,
|
|
# dirs_exist_ok=True)
|
|
|
|
@timer
|
|
def export(self, study, repertory, qlog=None):
|
|
self._study = study
|
|
name = study.name.replace(" ", "_")
|
|
|
|
self.export_additional_files(study, repertory, qlog, name=name)
|
|
|
|
return True
|
|
|
|
###########
|
|
# RESULTS #
|
|
###########
|
|
|
|
def read_bin(self, study, repertory, results, qlog=None, name="0"):
|
|
return
|
|
|
|
@timer
|
|
def results(self, study, repertory, qlog=None, name="0", type_pol=None):
|
|
results = Results(
|
|
study=study,
|
|
solver=self,
|
|
repertory=repertory,
|
|
name=name,
|
|
type_pol=type_pol,
|
|
)
|
|
self.read_bin(study, repertory, results, qlog, name=name)
|
|
|
|
return results
|
|
|
|
def output_param(self):
|
|
name = ""
|
|
return f"{name}"
|
|
|
|
_alph = list(
|
|
map(
|
|
chr,
|
|
chain(
|
|
range(48, 58), # 0..9
|
|
range(65, 91), # A..Z
|
|
range(97, 123) # a..z
|
|
)
|
|
)
|
|
)
|
|
|
|
_l_alph = len(_alph)
|
|
|
|
_nodes_cnt = 0
|
|
_nodes_names = {}
|
|
_nodes_views = set()
|
|
|
|
#################################
|
|
# Adis-TT in weak coupling mode #
|
|
#################################
|
|
|
|
|
|
class AdisTTwc(AdisTT):
|
|
_type = "adisttwc"
|
|
|
|
def __init__(self, name):
|
|
super(AdisTTwc, self).__init__(name)
|
|
|
|
self._type = "adisttwc"
|
|
|
|
@classmethod
|
|
def default_parameters(cls):
|
|
lst = super(AdisTTwc, cls).default_parameters()
|
|
|
|
# Insert new parameters at specific position
|
|
names = list(map(lambda t: t[0], lst))
|
|
|
|
return lst
|
|
|
|
##########
|
|
# Export #
|
|
##########
|
|
|
|
def cmd_args(self, study):
|
|
lst = super(AdisTTwc, self).cmd_args(study)
|
|
|
|
return lst
|
|
|
|
# def _export_D90(self, study, repertory, qlog=None, name="0"):
|
|
|
|
# files = []
|
|
|
|
# if qlog is not None:
|
|
# qlog.put("Export D90 file")
|
|
|
|
# with adistt_file_open(
|
|
# os.path.join(repertory, f"{name}.D90"), "w+"
|
|
# ) as f:
|
|
# files.append(f"{name}.D90")
|
|
|
|
# f.write(f"*Diamètres caractéristiques du fond stable\n")
|
|
|
|
# d90AdisTT = study.river.d90_adistt.D90_AdisTT_List
|
|
|
|
# f.write(f"DEFAULT = {d90AdisTT[0].d90}\n")
|
|
|
|
# self._export_d90_spec(study, d90AdisTT[0]._data, f, qlog)
|
|
|
|
# return files
|
|
|
|
# def _export_d90_spec(self, study, d90_spec_data, f, qlog, name="0"):
|
|
|
|
# for d90_spec in d90_spec_data:
|
|
# if (d90_spec.name is None
|
|
# or d90_spec.reach is None
|
|
# or d90_spec.start_rk is None
|
|
# or d90_spec.end_rk is None
|
|
# or d90_spec.d90 is None):
|
|
# return
|
|
|
|
# edges = study.river.enable_edges()
|
|
|
|
# id_edges = list(map(lambda x: x.id, edges))
|
|
|
|
# id_reach = d90_spec.reach
|
|
# reach = next((x for x in edges if x.id == id_reach), None)
|
|
|
|
# if reach is None:
|
|
# return
|
|
|
|
# f.write(" ".join((f"{d90_spec.name}",
|
|
# "=",
|
|
# f"{study.river.get_edge_id(reach)+1}",
|
|
# f"{d90_spec.start_rk}",
|
|
# f"{d90_spec.end_rk}",
|
|
# f"{d90_spec.d90}\n")))
|
|
|
|
# def _export_DIF(self, study, repertory, qlog=None, name="0"):
|
|
|
|
# files = []
|
|
|
|
# if qlog is not None:
|
|
# qlog.put("Export DIF file")
|
|
|
|
# with adistt_file_open(
|
|
# os.path.join(repertory, f"{name}.DIF"), "w+"
|
|
# ) as f:
|
|
# files.append(f"{name}.DIF")
|
|
|
|
# f.write(f"*Définition des paramètres des fonctions de calcul du\n")
|
|
# f.write(f"*coefficient de diffusion\n")
|
|
|
|
# difAdisTT = study.river.dif_adistt.DIF_AdisTT_List
|
|
|
|
# if difAdisTT[0].method != "generique":
|
|
# f.write(" ".join((f"defaut = ",
|
|
# f"{difAdisTT[0].method}",
|
|
# f"{difAdisTT[0].dif}\n")))
|
|
# else:
|
|
# f.write(" ".join((f"defaut ="
|
|
# f"{difAdisTT[0].method}",
|
|
# f"{difAdisTT[0].dif}",
|
|
# f"{difAdisTT[0].b}",
|
|
# f"{difAdisTT[0].c}\n")))
|
|
|
|
# self._export_dif_spec(study, difAdisTT[0]._data, f, qlog)
|
|
|
|
# return files
|
|
|
|
# def _export_dif_spec(self, study, dif_spec_data, f, qlog, name="0"):
|
|
|
|
# for dif_spec in dif_spec_data:
|
|
# if (dif_spec.reach is None
|
|
# or dif_spec.start_rk is None
|
|
# or dif_spec.end_rk is None
|
|
# or dif_spec.dif is None
|
|
# or dif_spec.b is None
|
|
# or dif_spec.c is None):
|
|
# return
|
|
|
|
# edges = study.river.enable_edges()
|
|
|
|
# id_reach = dif_spec.reach
|
|
# reach = next((x for x in edges if x.id == id_reach), None)
|
|
|
|
# if reach is None:
|
|
# return
|
|
|
|
# if dif_spec.method != "generique":
|
|
# f.write(" ".join((f"{dif_spec.method}",
|
|
# "=",
|
|
# f"{study.river.get_edge_id(reach)+1}",
|
|
# f"{dif_spec.start_rk}",
|
|
# f"{dif_spec.end_rk}",
|
|
# f"{dif_spec.dif}\n")))
|
|
# else:
|
|
# f.write(" ".join((f"{dif_spec.method}",
|
|
# f"=" f"{study.river.get_edge_id(reach)+1}",
|
|
# f"{dif_spec.start_rk}",
|
|
# f"{dif_spec.end_rk}",
|
|
# f"{dif_spec.dif}",
|
|
# f"{dif_spec.b}",
|
|
# f"{dif_spec.c}\n")))
|
|
|
|
def _export_NUM(self, study, repertory, qlog=None, name="0"):
|
|
|
|
dict_names = {"init_time": "start_date",
|
|
"final_time": "end_date",
|
|
"timestep": "dt0",
|
|
"implicitation_parameter": "theta",
|
|
"timestep_screen": "dtscr",
|
|
"timestep_bin": "dtbin",
|
|
"timestep_csv": "dtcsv",
|
|
"timestep_mage": "dtMage",
|
|
"initial_concentration": "c_initiale"}
|
|
files = []
|
|
|
|
if qlog is not None:
|
|
qlog.put("Export NUM file")
|
|
|
|
with adistt_file_open(
|
|
os.path.join(repertory, f"{name}.NUM"), "w+"
|
|
) as f:
|
|
files.append(f"{name}.NUM")
|
|
|
|
params = study.river.get_params(self.type).parameters
|
|
for p in params:
|
|
name = p.name\
|
|
.replace("all_", "")\
|
|
.replace("adistt_", "")
|
|
value = p.value
|
|
|
|
logger.debug(
|
|
f"export: NUM: {name}: {value} ({p.value})"
|
|
)
|
|
|
|
if name != "command_line_arguments":
|
|
f.write(f"{dict_names[name]} = {value}\n")
|
|
|
|
return files
|
|
|
|
@timer
|
|
def read_bin(self, study, repertory, results, qlog=None, name="0"):
|
|
|
|
filelist = [f for f in os.listdir(repertory)
|
|
if os.path.isfile(os.path.join(repertory, f))
|
|
]
|
|
files_bin_names = [f for f in filelist if f[-4:] == ".bin"]
|
|
files_bin_names.insert(0, files_bin_names.pop(
|
|
files_bin_names.index("total_sediment.bin"))
|
|
)
|
|
|
|
ifilename = os.path.join(repertory, files_bin_names[0])
|
|
|
|
logger.info(f"read_bin: Start reading '{ifilename}' ...")
|
|
|
|
with open(ifilename, 'rb') as f:
|
|
# header
|
|
# first line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
data = np.fromfile(f, dtype=np.int32, count=3)
|
|
ibmax = data[0] # number of reaches
|
|
ismax = data[1] # total number of cross sections
|
|
kbl = data[2] * -1 # block size for .BIN header
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# second line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
ibu = np.fromfile(f, dtype=np.int32, count=ibmax)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# third line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
data = np.fromfile(f, dtype=np.int32, count=2 * ibmax)
|
|
is1 = np.zeros(ibmax, dtype=np.int32)
|
|
is2 = np.zeros(ibmax, dtype=np.int32)
|
|
|
|
logger.debug(f"read_bin: nb_reach = {ibmax}")
|
|
logger.debug(f"read_bin: nb_profile = {ismax}")
|
|
|
|
results.set("nb_reach", f"{ibmax}")
|
|
results.set("nb_profile", f"{ismax}")
|
|
|
|
reachs = []
|
|
iprofiles = {}
|
|
reach_offset = {}
|
|
|
|
for i in range(ibmax):
|
|
# Add results reach to reach list
|
|
r = results.river.add(i)
|
|
reachs.append(r)
|
|
|
|
is1[i] = data[2 * i] - 1 # first section of reach i
|
|
is2[i] = data[2 * i + 1] - 1 # last section of reach i
|
|
|
|
key = (is1[i], is2[i])
|
|
iprofiles[key] = r
|
|
|
|
reach_offset[r] = is1[i]
|
|
|
|
logger.debug(f"read_bin: iprofiles = {iprofiles}")
|
|
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# fourth line
|
|
pk = np.zeros(ismax, dtype=np.float32)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
pk[k:min(k + kbl, ismax)] = np.fromfile(f,
|
|
dtype=np.float32,
|
|
count=min(
|
|
k + kbl, ismax
|
|
) - k)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
|
|
# fifth line (useless)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
zmin_OLD = np.fromfile(f, dtype=np.float32, count=1)[0]
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# sixth line
|
|
zf = np.zeros(ismax, dtype=np.float32)
|
|
z = np.zeros(ismax * 3, dtype=np.float32)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
z[3 * k:3 * min(k + kbl, ismax)] = \
|
|
np.fromfile(f,
|
|
dtype=np.float32,
|
|
count=3 * (min(k + kbl, ismax) - k)
|
|
)
|
|
# z[i*3+1] and z[i*3+2] are useless
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
zf = [z[i * 3] for i in range(ismax)]
|
|
# seventh line (useless)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
zero = np.fromfile(f, dtype=np.int32, count=ismax)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# end header
|
|
|
|
def ip_to_r(i):
|
|
return iprofiles[
|
|
next(
|
|
filter(
|
|
lambda k: k[0] <= i <= k[1],
|
|
iprofiles
|
|
)
|
|
)
|
|
]
|
|
|
|
def ip_to_ri(r, i): return i - reach_offset[r]
|
|
|
|
path_files = map(lambda file: os.path.join(
|
|
repertory, file), files_bin_names)
|
|
|
|
data_tmp = {}
|
|
|
|
for file_bin in path_files:
|
|
key_pol = os.path.basename(file_bin)[0:-4]
|
|
data_tmp[key_pol] = {}
|
|
logger.info(f"read_bin: Start reading '{file_bin}' ...")
|
|
with open(file_bin, 'rb') as f:
|
|
# header
|
|
# first line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
data = np.fromfile(f, dtype=np.int32, count=3)
|
|
ibmax = data[0] # number of reaches
|
|
ismax = data[1] # total number of cross sections
|
|
kbl = data[2] * -1 # block size for .BIN header
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# second line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
ibu = np.fromfile(f, dtype=np.int32, count=ibmax)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# third line
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
data = np.fromfile(f, dtype=np.int32, count=2 * ibmax)
|
|
is1 = np.zeros(ibmax, dtype=np.int32)
|
|
is2 = np.zeros(ibmax, dtype=np.int32)
|
|
for i in range(ibmax):
|
|
# first section of reach i (FORTRAN numbering)
|
|
is1[i] = data[2 * i]
|
|
# last section of reach i (FORTRAN numbering)
|
|
is2[i] = data[2 * i + 1]
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# fourth line
|
|
pk = np.zeros(ismax, dtype=np.float32)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
pk[k:min(k + kbl, ismax)] = np.fromfile(
|
|
f, dtype=np.float32, count=min(k + kbl, ismax) - k
|
|
)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# fifth line (useless)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
zmin_OLD = np.fromfile(f, dtype=np.float32, count=1)[0]
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# sixth line
|
|
zf = np.zeros(ismax, dtype=np.float32)
|
|
z = np.zeros(ismax * 3, dtype=np.float32)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
z[3 * k:3 * min(k + kbl, ismax)] = np.fromfile(
|
|
f, dtype=np.float32,
|
|
count=3 * (min(k + kbl, ismax) - k)
|
|
)
|
|
# z[i*3+1] and z[i*3+2] are useless
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
zf = [z[i * 3] for i in range(ismax)]
|
|
# seventh line (useless)
|
|
for k in range(0, ismax, kbl):
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
zero = np.fromfile(
|
|
f, dtype=np.int32, count=min(k + kbl, ismax) - k
|
|
)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
# end header
|
|
# data
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
while data.size > 0:
|
|
ismax = np.fromfile(f, dtype=np.int32, count=1)[0]
|
|
t = np.fromfile(f, dtype=np.float64, count=1)[0]
|
|
if t not in data_tmp[key_pol]:
|
|
data_tmp[key_pol][t] = {}
|
|
c = np.fromfile(f, dtype=np.byte, count=1)
|
|
# possible values :
|
|
# sediment : C, G, M, D, L, N, R
|
|
# polutant : C
|
|
phys_var = bytearray(c).decode()
|
|
data_tmp[key_pol][t][phys_var] = {}
|
|
real_data = np.fromfile(f, dtype=np.float32, count=ismax)
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (end)
|
|
data_tmp[key_pol][t][phys_var] = real_data
|
|
data = np.fromfile(f, dtype=np.int32, count=1) # (start)
|
|
# end data
|
|
|
|
pollutants_keys = list(data_tmp.keys())
|
|
timestamps_keys = list(data_tmp[pollutants_keys[0]].keys())
|
|
phys_data_names = list(data_tmp[pollutants_keys[0]]
|
|
[timestamps_keys[0]].keys())
|
|
type_pol_index = len(phys_data_names)
|
|
|
|
for r, reach in enumerate(reachs):
|
|
for i in range(is1[r]-1, is2[r]):
|
|
p_i = ip_to_ri(reach, i)
|
|
for t_data in timestamps_keys:
|
|
pol_view = []
|
|
for pol in pollutants_keys:
|
|
pol_view.append(tuple(list(map(
|
|
lambda data_el: data_el[i],
|
|
list(data_tmp[pol][t_data].values())
|
|
))))
|
|
|
|
reach.set(p_i, t_data, "pols", pol_view)
|
|
|
|
results.set("timestamps", set(timestamps_keys))
|
|
|
|
@timer
|
|
def results(self, study, repertory, qlog=None, name=None):
|
|
self._study = study
|
|
if name is None:
|
|
name = study.name.replace(" ", "_")
|
|
|
|
results = super(AdisTTwc, self).results(study,
|
|
repertory,
|
|
qlog,
|
|
name=name)
|
|
|
|
return results
|
|
|
|
_weather_files = {
|
|
"AT": ("tda0", "TDA"),
|
|
"SH": ("hsa0", "HSA"),
|
|
"GR": ("ray0", "RAY"),
|
|
"RWS": ("ven0", "VEN"),
|
|
"GT": ("tna0", "TNA"),
|
|
"GFR": ("qna0", "QNA"),
|
|
"ALB": ("alb0", "ALB"),
|
|
"SC": ("SF0", "SF"),
|
|
"CCF": ("cld0", "CLD"),
|
|
}
|
|
|
|
def _export_temperature_initial_conditions(self, study, repertory):
|
|
initial_conditions = study.river.ic_temperature.Initial_Conditions_List
|
|
filename = "TEM.INI"
|
|
with adistt_file_open(os.path.join(repertory, filename), "w+") as f:
|
|
f.write("* Initial temperature\n")
|
|
initial = initial_conditions[0] if initial_conditions else None
|
|
temperature = initial.temperature if initial is not None else 0.0
|
|
if temperature is None:
|
|
temperature = 0.0
|
|
f.write(f"DEFAULT = {temperature} 0.0 0.0 0.0\n")
|
|
|
|
edges = study.river.enable_edges()
|
|
specifications = initial._data if initial is not None else []
|
|
for spec in specifications:
|
|
if spec.is_deleted():
|
|
continue
|
|
reach = next((edge for edge in edges
|
|
if edge.id == spec.reach), None)
|
|
if reach is None or spec.temperature is None:
|
|
continue
|
|
reach_id = study.river.get_edge_id(reach) + 1
|
|
f.write(
|
|
f"{spec.name} = {reach_id} {spec.start_rk} "
|
|
f"{spec.end_rk} {spec.temperature} 0.0 0.0 0.0 0.0\n"
|
|
)
|
|
return filename
|
|
|
|
def _export_temperature_boundary_conditions(self, study, repertory):
|
|
boundary_conditions = (
|
|
study.river.boundary_conditions_temperature
|
|
.BCs_Temperature_List
|
|
)
|
|
filename = "TEM.CDT"
|
|
with adistt_file_open(os.path.join(repertory, filename), "w+") as f:
|
|
for boundary_condition in boundary_conditions:
|
|
node = next((node for node in study.river.nodes()
|
|
if node.id == boundary_condition.node), None)
|
|
if node is None:
|
|
continue
|
|
f.write(f"${self.get_node_name(node)}\n")
|
|
f.write("*temps |temperature (degC)\n")
|
|
f.write("*---------++++++++++\n")
|
|
for value in boundary_condition.data:
|
|
date = timestamp_to_old_pamhyr_date_adists(int(value[0]))
|
|
f.write(f"{date} {value[1]}\n")
|
|
f.write("*\n")
|
|
return filename
|
|
|
|
def _export_weather_file(self, study, repertory, weather_parameter,
|
|
extension):
|
|
reach = weather_parameter.reach
|
|
if reach is None or len(weather_parameter.data) == 0:
|
|
return None
|
|
|
|
filename = f"TEM.{extension}"
|
|
path = os.path.join(repertory, filename)
|
|
mode = "a" if os.path.exists(path) else "w+"
|
|
with adistt_file_open(path, mode) as f:
|
|
reach_name = self.get_reach_name(reach)
|
|
f.write(
|
|
f"${reach_name} {weather_parameter.begin_rk} "
|
|
f"{weather_parameter.end_rk}\n"
|
|
)
|
|
f.write("*temps |value\n")
|
|
f.write("*---------++++++++++\n")
|
|
for value in weather_parameter.data:
|
|
date = timestamp_to_old_pamhyr_date_adists(int(value[0]))
|
|
f.write(f"{date} {value[1]}\n")
|
|
f.write("*\n")
|
|
return filename
|
|
|
|
def _export_TEM(self, study, repertory, qlog=None, name="0"):
|
|
if qlog is not None:
|
|
qlog.put("Export TEM files")
|
|
|
|
legacy_pol = os.path.join(repertory, "TEM.POL")
|
|
if os.path.exists(legacy_pol):
|
|
os.remove(legacy_pol)
|
|
|
|
# Weather files are appended range by range below. Remove files from
|
|
# a previous export first so rerunning an unchanged study is stable.
|
|
for _, extension in self._weather_files.values():
|
|
path = os.path.join(repertory, f"TEM.{extension}")
|
|
if os.path.exists(path):
|
|
os.remove(path)
|
|
|
|
initial_file = self._export_temperature_initial_conditions(
|
|
study, repertory
|
|
)
|
|
boundary_file = self._export_temperature_boundary_conditions(
|
|
study, repertory
|
|
)
|
|
|
|
weather_files = {}
|
|
for weather_parameter in study.river.weather_parameters.lst:
|
|
config = self._weather_files.get(weather_parameter.type)
|
|
if (config is None or weather_parameter.reach is None
|
|
or weather_parameter.reach.is_deleted()):
|
|
continue
|
|
_, extension = config
|
|
exported = self._export_weather_file(
|
|
study, repertory, weather_parameter, extension
|
|
)
|
|
if exported is not None:
|
|
weather_files[extension] = exported
|
|
|
|
filename = f"{name}.TEM"
|
|
with adistt_file_open(os.path.join(repertory, filename), "w+") as f:
|
|
defaults = study.river.weather_parameters
|
|
for type_, (parameter, _) in self._weather_files.items():
|
|
default = defaults.default_for_type(type_)
|
|
if default is not None and not default.is_deleted():
|
|
f.write(f"{parameter} = {default.value}\n")
|
|
|
|
f.write(f"file_ini = {initial_file}\n")
|
|
f.write(f"file_cl = {boundary_file}\n")
|
|
for _, extension in self._weather_files.values():
|
|
if extension in weather_files:
|
|
f.write(
|
|
f"file_{extension.lower()} = "
|
|
f"{weather_files[extension]}\n"
|
|
)
|
|
|
|
return [filename]
|
|
|
|
def export_func_dict(self):
|
|
return [
|
|
self._export_NUM,
|
|
self._export_TEM
|
|
]
|
|
|
|
def rm_previous_results(self, study, repertory, qlog):
|
|
|
|
if "resultats" in os.listdir(repertory):
|
|
repertory_results = os.path.join(repertory, "resultats")
|
|
else:
|
|
repertory_results = os.path.normpath(repertory)
|
|
|
|
if not os.path.isdir(repertory_results):
|
|
return
|
|
|
|
filelist = [f for f in os.listdir(repertory_results)
|
|
if os.path.isfile(os.path.join(repertory_results, f))
|
|
]
|
|
files_bin_names = [f for f in filelist if f[-4:] == ".bin"]
|
|
if len(files_bin_names) < 1:
|
|
return
|
|
|
|
for el in files_bin_names:
|
|
os.remove(os.path.join(repertory_results, el))
|
|
|
|
@timer
|
|
def export(self, study, repertory, mage_rep, qlog=None, name="0"):
|
|
logger.debug(f"cmd solver adisttwc : {self._cmd_solver}")
|
|
self._study = study
|
|
name = study.name.replace(" ", "_")
|
|
|
|
mage_path = os.path.join(
|
|
os.path.abspath(os.path.join(repertory, os.pardir)),
|
|
mage_rep
|
|
)
|
|
source_net = os.path.join(mage_path, "net")
|
|
destination_net = os.path.join(repertory, "net")
|
|
|
|
if qlog is not None:
|
|
qlog.put("Copy Mage net directory")
|
|
|
|
if os.path.exists(destination_net):
|
|
shutil.rmtree(destination_net)
|
|
shutil.copytree(source_net, destination_net)
|
|
|
|
# Node names must be deterministic and identical in NET and TEM.CDT.
|
|
self._nodes_cnt = 0
|
|
self._nodes_names = {}
|
|
self._nodes_views = set()
|
|
for edge in study.river.enable_edges():
|
|
self.get_node_name(edge.node1)
|
|
self.get_node_name(edge.node2)
|
|
|
|
# Generate files
|
|
files = []
|
|
|
|
self.rm_previous_results(study, repertory, qlog)
|
|
|
|
try:
|
|
for func in self.export_func_dict():
|
|
files = files + func(study, repertory, qlog, name=name)
|
|
|
|
self.export_additional_files(study, repertory, qlog, name=name)
|
|
self._export_REP(study, repertory, mage_rep,
|
|
files, qlog, name=name)
|
|
|
|
return True
|
|
except Exception as e:
|
|
logger.error(f"Exception occurred during export: {e}")
|
|
logger.error(f"Failed to export study to {self._type}")
|
|
logger_exception(e)
|
|
return False
|