add scripts to use Mascaret

setup.py
Theophile Terraz 2023-11-02 16:43:57 +01:00
parent 325bb9ed02
commit 32777456b9
2 changed files with 440 additions and 0 deletions

341
src/Scripts/mage_to_mascaret.py Executable file
View File

@ -0,0 +1,341 @@
#!/bin/env python3
import sys
import importlib.util as imp
import os
import datetime as dt
import lxml.etree as ET
def mage_to_mascaret(filename):
# rep file
with open(filename, encoding = "ISO-8859-1", mode ="r") as rep:
lines = rep.readlines()
for line in lines:
l = line.split()
if (l[0] == 'PAR'):
parfile = l[1]
elif (l[0] == 'NET'):
netfile = l[1]
elif (l[0] == 'HYD'):
hydfile = l[1]
elif (l[0] == 'RUG'):
rugfile = l[1]
elif (l[0] == 'LIM'):
limfile = l[1]
elif (l[0] == 'INI'):
inifile = l[1]
# net file
with open(netfile, encoding = "ISO-8859-1", mode ="r") as net:
lines = net.readlines()
for line in lines:
if (line[0] != '*' and line[0] != '$'):
l = line.split()
biefname = l[0]
stfile = l[3]
break # on n'utilise que le premier bief pour le moment
# geometry file
print("geometrie fie : "+stfile)
wb.init_bief_from_geo_file(stfile, 0, 0)
wb.set_bief_name(biefname)
nb_sect = wb.get_nb_sections()
print("nb sections: "+str(nb_sect))
# hydrogramme file
with open(hydfile, encoding = "ISO-8859-1", mode ="r") as hyd:
HYD = []
lines = hyd.readlines()
for line in lines:
if (line[0] != '*' and line[0] != '$'):
HYD.append([x for x in line.split()])
# friction file
with open(rugfile, encoding = "ISO-8859-1", mode ="r") as lim:
RUG = []
lines = lim.readlines()
for line in lines:
if (line[0] != '*' and line[0] != '$'):
RUG.append([x for x in line.split()])
# tarage file
with open(limfile, encoding = "ISO-8859-1", mode ="r") as tar:
LIM = []
lines = tar.readlines()
for line in lines:
if (line[0] != '*' and line[0] != '$'):
LIM.append([x for x in line.split()])
# cond init file
with open(inifile, encoding = "ISO-8859-1", mode ="r") as ini:
INI = []
lines = ini.readlines()
for line in lines:
if (line[0] != '*' and line[0] != '$'):
INI.append([x for x in line.split()])
# PAR file
with open(parfile, encoding = "ISO-8859-1", mode ="r") as par:
PAR = {}
lines = par.readlines()
for line in lines:
l = line.split()
if l[0] == 'init_time':
t = l[1].split(":")
if len(t)>1:
PAR['init_time'] = int(t[0])*24*3600+int(t[1])*3600+int(t[2])*60+int(t[3])
else:
PAR['init_time'] = int(t[0])*60
elif l[0] == 'final_time':
t = l[1].split(":")
if len(t)>1:
PAR['final_time'] = int(t[0])*24*3600+int(t[1])*3600+int(t[2])*60+int(t[3])
else:
PAR['final_time'] = int(t[0])*60
elif l[0] == 'timestep':
PAR['timestep'] = float(l[1])
elif l[0] == 'timestep_bin':
PAR['timestep_bin'] = float(l[1])
elif l[0] == 'timestep_tra':
PAR['timestep_tra'] = float(l[1])
elif l[0] == 'min_height':
PAR['min_height'] = float(l[1])
#os.chdir(current_dir)
with open("hydrogramme.loi", mode ="w") as hyd:
hyd.write(" S\n")
for line in HYD:
hyd.write(" "+line[0]+" "+line[1]+"\n")
if float(HYD[-1][0]) < PAR['final_time']:
hyd.write(" "+str(PAR['final_time'])+" "+str(float(HYD[-1][1]))+"\n")
with open("limnigramme.loi", mode ="w") as lim:
lim.write(" S\n")
for line in LIM:
lim.write(" "+line[0]+" "+line[1]+"\n")
if float(LIM[-1][0]) < PAR['final_time']:
lim.write(" "+str(PAR['final_time'])+" "+str(float(LIM[-1][1]))+"\n")
with open("init.lig", mode ="w") as lig:
lig.write("RESULTATS CALCUL,DATE : "+dt.datetime.now().strftime("%d/%m/%y %H:%M"+"\n"))
lig.write("FICHIER RESULTAT MASCARET\n")
lig.write("-"*71+"\n")
lig.write(" IMAX = "+f"{len(INI):4d}"+" NBBIEF= 1\n")
lig.write(" I1,I2 = 1 "+str(len(INI))+"\n")
lig.write(" X\n")
for i in range(int(len(INI)/5)+1):
lig.write(" "+"".join([f"{float(INI[j][4]):12.2f}" for j in range(i*5, min(i*5+5, len(INI)))])+"\n")
lig.write(" Z\n")
for i in range(int(len(INI)/5)+1):
lig.write(" "+"".join([f"{float(INI[j][3]):12.3f}" for j in range(i*5, min(i*5+5, len(INI)))])+"\n")
lig.write(" Q\n")
for i in range(int(len(INI)/5)+1):
lig.write(" "+"".join([f"{float(INI[j][2]):12.3f}" for j in range(i*5, min(i*5+5, len(INI)))])+"\n")
lig.write(" FIN\n")
#os.remove("geometrie")
wb.output_bief_mascaret("geometrie")
pk = [wb.get_pk_section(i+1) for i in range(nb_sect)]
print(pk)
print(max(pk))
print(min(pk))
out_name = 'FichierCas.xcas'
root = ET.Element('fichierCas')
pc = ET.SubElement(root, 'parametresCas')
pg = ET.SubElement(pc, 'parametresGeneraux')
ET.SubElement(pg, 'versionCode').text='3'
ET.SubElement(pg, 'code').text = '3'
ET.SubElement(pg, 'fichMotsCles').text = out_name
ET.SubElement(pg, 'dictionaire').text='dico.txt'
ET.SubElement(pg, 'progPrincipal').text='princi.f'
ET.SubElement(pg, 'sauveModele').text='false'
ET.SubElement(pg, 'fichSauvModele').text='mascaret_exp.tmp'
ET.SubElement(pg, 'validationCode').text='false'
ET.SubElement(pg, 'typeValidation').text='1'
ET.SubElement(pg, 'presenceCasiers').text='false'
b = ET.SubElement(pg, 'bibliotheques')
ET.SubElement(b, 'bibliotheque').text='mascaretV5P1.a damoV3P0.a'
pmp = ET.SubElement(pc, 'parametresModelePhysique')
ET.SubElement(pmp, 'perteChargeConf').text='false'
ET.SubElement(pmp, 'compositionLits').text='1'
ET.SubElement(pmp, 'conservFrotVertical').text='false'
ET.SubElement(pmp, 'elevCoteArrivFront').text='0.05'
ET.SubElement(pmp, 'interpolLinStrickler').text='false'
d = ET.SubElement(pmp, 'debordement')
ET.SubElement(d, 'litMajeur').text='false'
ET.SubElement(d, 'zoneStock').text='false'
pn = ET.SubElement(pc, 'parametresNumeriques')
ET.SubElement(pn, 'calcOndeSubmersion').text='false'
ET.SubElement(pn, 'decentrement').text='false'
ET.SubElement(pn, 'froudeLimCondLim').text='1000.0'
ET.SubElement(pn, 'traitImplicitFrot').text='false'
ET.SubElement(pn, 'hauteurEauMini').text=str(PAR['min_height'])
ET.SubElement(pn, 'implicitNoyauTrans').text='true'
ET.SubElement(pn, 'optimisNoyauTrans').text='false'
ET.SubElement(pn, 'perteChargeAutoElargissement').text='false'
ET.SubElement(pn, 'termesNonHydrostatiques').text='false'
ET.SubElement(pn, 'apportDebit').text='0'
ET.SubElement(pn, 'attenuationConvection').text='false'
pt = ET.SubElement(pc, 'parametresTemporels')
ET.SubElement(pt, 'pasTemps').text=str(PAR['timestep'])
ET.SubElement(pt, 'tempsInit').text=str(PAR['init_time'])
ET.SubElement(pt, 'critereArret').text='1'
ET.SubElement(pt, 'nbPasTemps').text='2000'
ET.SubElement(pt, 'tempsMax').text=str(PAR['final_time'])
ET.SubElement(pt, 'pasTempsVar').text='true'
ET.SubElement(pt, 'nbCourant').text='0.8'
ET.SubElement(pt, 'coteMax').text='0.0'
ET.SubElement(pt, 'abscisseControle').text='0.0'
ET.SubElement(pt, 'biefControle').text='1'
pgr = ET.SubElement(pc, 'parametresGeometrieReseau')
g = ET.SubElement(pgr, 'geometrie')
ET.SubElement(g, 'fichier').text='geometrie'
ET.SubElement(g, 'format').text='2'
ET.SubElement(g, 'profilsAbscAbsolu').text='true'
lb = ET.SubElement(pgr, 'listeBranches')
ET.SubElement(lb, 'nb').text='1'
ET.SubElement(lb, 'numeros').text='1'
ET.SubElement(lb, 'abscDebut').text=str(min(pk))
ET.SubElement(lb, 'abscFin').text=str(max(pk))
ET.SubElement(lb, 'numExtremDebut').text='1'
ET.SubElement(lb, 'numExtremFin').text='2'
ln = ET.SubElement(pgr, 'listeNoeuds')
ET.SubElement(ln, 'nb').text='0'
ET.SubElement(ln, 'noeuds')
el = ET.SubElement(pgr, 'extrLibres')
ET.SubElement(el, 'nb').text='2'
ET.SubElement(el, 'num').text='1 2'
ET.SubElement(el, 'numExtrem').text='1 2'
noms = ET.SubElement(el, 'noms')
ET.SubElement(noms, 'string').text='limite1'
ET.SubElement(noms, 'string').text='limite2'
ET.SubElement(el, 'typeCond').text='1 2'
ET.SubElement(el, 'numLoi').text='1 2'
pc2 = ET.SubElement(pc, 'parametresConfluents')
ET.SubElement(pc2, 'nbConfluents').text='0'
ET.SubElement(pc2, 'confluents')
pcm = ET.SubElement(pc, 'parametresPlanimetrageMaillage')
ET.SubElement(pcm, 'methodeMaillage').text='1'
planim = ET.SubElement(pcm, 'planim')
ET.SubElement(planim, 'nbPas').text='101'
ET.SubElement(planim, 'nbZones').text='1'
ET.SubElement(planim, 'valeursPas').text='0.1'
ET.SubElement(planim, 'num1erProf').text='1'
ET.SubElement(planim, 'numDerProf').text=str(nb_sect)
maillage = ET.SubElement(pcm, 'maillage')
ET.SubElement(maillage, 'modeSaisie').text='2'
ET.SubElement(maillage, 'sauvMaillage').text='false'
maillageClavier = ET.SubElement(maillage, 'maillageClavier')
ET.SubElement(maillageClavier, 'nbSections').text='0'
ET.SubElement(maillageClavier, 'nbPlages').text='1'
ET.SubElement(maillageClavier, 'num1erProfPlage').text='1'
ET.SubElement(maillageClavier, 'numDerProfPlage').text=str(nb_sect)
ET.SubElement(maillageClavier, 'pasEspacePlage').text='10.0'
ET.SubElement(maillageClavier, 'nbZones').text='0'
ps = ET.SubElement(pc, 'parametresSingularite')
ET.SubElement(ps, 'nbSeuils').text='0'
pad = ET.SubElement(pc, 'parametresApportDeversoirs')
pc3 = ET.SubElement(pc, 'parametresCalage')
frt = ET.SubElement(pc3, 'frottement')
ET.SubElement(frt, 'loi').text='1'
ET.SubElement(frt, 'nbZone').text=str(len(RUG))
ET.SubElement(frt, 'numBranche').text=' '.join([i[1] for i in RUG])
ET.SubElement(frt, 'absDebZone').text=' '.join([i[2] for i in RUG])
ET.SubElement(frt, 'absFinZone').text=' '.join([i[3] for i in RUG])
ET.SubElement(frt, 'coefLitMin').text=' '.join([i[4] for i in RUG])
ET.SubElement(frt, 'coefLitMaj').text=' '.join([i[5] for i in RUG])
zst = ET.SubElement(pc3, 'zoneStockage')
ET.SubElement(zst, 'nbProfils').text='0'
ET.SubElement(zst, 'numProfil').text='-0'
ET.SubElement(zst, 'limGauchLitMaj').text='-0'
ET.SubElement(zst, 'limDroitLitMaj').text='-0'
plh = ET.SubElement(pc, 'parametresLoisHydrauliques')
ET.SubElement(plh, 'nb').text='2'
lois = ET.SubElement(plh, 'lois')
spl = ET.SubElement(lois, 'structureParametresLoi')
ET.SubElement(spl, 'nom').text='loi_1_hydrogramme'
ET.SubElement(spl, 'type').text='1'
donnees = ET.SubElement(spl, 'donnees')
ET.SubElement(donnees, 'modeEntree').text='1'
ET.SubElement(donnees, 'fichier').text='hydrogramme.loi'
ET.SubElement(donnees, 'uniteTps').text='-0'
ET.SubElement(donnees, 'nbPoints').text='-0'
ET.SubElement(donnees, 'nbDebitsDifferents').text='-0'
spl2 = ET.SubElement(lois, 'structureParametresLoi')
ET.SubElement(spl2, 'nom').text='loi_2_limnigramme'
ET.SubElement(spl2, 'type').text='2'
donnees = ET.SubElement(spl2, 'donnees')
ET.SubElement(donnees, 'modeEntree').text='1'
ET.SubElement(donnees, 'fichier').text='limnigramme.loi'
ET.SubElement(donnees, 'uniteTps').text='-0'
ET.SubElement(donnees, 'nbPoints').text='-0'
ET.SubElement(donnees, 'nbDebitsDifferents').text='-0'
pci = ET.SubElement(pc, 'parametresConditionsInitiales')
re = ET.SubElement(pci, 'repriseEtude')
ET.SubElement(re, 'repriseCalcul').text='false'
le = ET.SubElement(pci, 'ligneEau')
ET.SubElement(le, 'LigEauInit').text='true'
ET.SubElement(le, 'modeEntree').text='1'
ET.SubElement(le, 'fichLigEau').text='init.lig'
ET.SubElement(le, 'formatFichLig').text='2'
ET.SubElement(le, 'nbPts').text='-0'
pir = ET.SubElement(pc, 'parametresImpressionResultats')
ET.SubElement(pir, 'titreCalcul').text='Etude hydraulique1d'
impression = ET.SubElement(pir, 'impression')
ET.SubElement(impression, 'impressionGeometrie').text='false'
ET.SubElement(impression, 'impressionPlanimetrage').text='false'
ET.SubElement(impression, 'impressionReseau').text='false'
ET.SubElement(impression, 'impressionLoiHydraulique').text='false'
ET.SubElement(impression, 'impressionligneEauInitiale').text='false'
ET.SubElement(impression, 'impressionCalcul').text='false'
ps = ET.SubElement(pir, 'pasStockage')
ET.SubElement(ps, 'premPasTpsStock').text='1'
ET.SubElement(ps, 'pasStock').text=str(int(PAR['timestep_bin']/PAR['timestep']))
ET.SubElement(ps, 'pasImpression').text=str(int(PAR['timestep_tra']/PAR['timestep']))
results = ET.SubElement(pir, 'resultats')
ET.SubElement(results, 'fichResultat').text='mascaret_exp.opt'
ET.SubElement(results, 'postProcesseur').text='2'
listing = ET.SubElement(pir, 'listing')
ET.SubElement(listing, 'fichListing').text='mascaret_exp.lis'
reprise = ET.SubElement(pir, 'fichReprise')
ET.SubElement(reprise, 'fichRepriseEcr').text='mascaret_exp.rep'
rub = ET.SubElement(pir, 'rubens')
ET.SubElement(rub, 'ecartInterBranch').text='1.0'
stockage = ET.SubElement(pir, 'stockage')
ET.SubElement(stockage, 'option').text='1'
ET.SubElement(stockage, 'nbSite').text='0'
pvc = ET.SubElement(pc, 'parametresVariablesCalculees')
ET.SubElement(pvc, 'variablesCalculees').text="false false false false false false false false false false false false false false false"
pvs = ET.SubElement(pc, 'parametresVariablesStockees')
ET.SubElement(pvs, 'variablesStockees').text="true false false false false true true true false false false false false false false false false false false true false false true false false false false false false false false false false false false false false false false false false false"
tree = ET.ElementTree(root)
xmlstr = ET.tostring(root, encoding='unicode', pretty_print=True)
tree.write(out_name)
print("fin du programme")
parser = ET.XMLParser(remove_blank_text=True)
tree = ET.parse(out_name, parser)
tree.write(out_name, pretty_print=True)
if __name__ == "__main__":
current_dir = os.getcwd()
os.chdir('/home/theophile.terraz/Codes/mage/src')
spec = imp.spec_from_file_location("wb", "./wrapper_bief.py")
wb = imp.module_from_spec(spec)
spec.loader.exec_module(wb)
os.chdir(current_dir)
print( 'Number of arguments:', len(sys.argv), 'arguments.')
print( 'Argument List:', str(sys.argv))
mage_to_mascaret(sys.argv[1])

99
src/Scripts/mascaret_to_mage.py Executable file
View File

@ -0,0 +1,99 @@
#!/bin/env python3
import sys
import importlib.util as imp
import os
import numpy as np
def lire_opt(filename, res):
with open(filename, encoding = "ISO-8859-1", mode ="r") as f:
variables = ['"t"','"Bief"','"Section"','"Pk"']
lines = f.readlines()
for start,line in enumerate(lines):
print("line n "+str(start)+" : ",line)
if line == "[variables]\n": continue
if line == "[resultats]\n": break
val = line.split(';')
variables.append(val[1])
if val[1] == '"ZREF"': izref = start+3
if val[1] == '"Q"': iq = start+3
if val[1] == '"Z"': iz = start+3
print ("variables : "+str(variables))
#resultats
res.mage_version = 83
# first pass to count sect and biefs
prev_timestep = float(lines[start+1].split(';')[0])
for line in lines[start+1:]:
ls = line.split(';')
t = float(ls[0])
if prev_timestep != t : break
ibief = int(ls[1].strip('"'))
isect = int(ls[2].strip('"'))
pk = float(ls[3])
if ibief > res.ibmax:
res.ibmax += 1
res.is1.append(isect)
res.is2.append(isect)
if isect > res.ismax:
res.ismax += 1
res.xgeo.append(pk)
res.ygeo.append(0.0)
res.zfd.append(float(ls[izref]))
res.ybas.append(0.0)
if isect > res.is2[ibief-1]:
res.is2[ibief-1] = isect
if isect < res.is1[ibief-1]:
res.is1[ibief-1] = isect
print(" nb sections : "+str(res.ismax))
# first pass to get the data
prev_bief = 0
prev_sect = 0
z = np.zeros(res.ismax, dtype=np.float32)
q = np.zeros(res.ismax, dtype=np.float32)
isect = 0
count = np.array([1],dtype=np.int32)
bZ = bytearray('Z'.encode())
bQ = bytearray('Q'.encode())
for line in lines[start+1:]:
ls = line.split(';')
t = float(ls[0])
ibief = int(ls[1].strip('"'))
isect = int(ls[2].strip('"'))
if prev_timestep != t :
res.values['Z'].append(prev_timestep, z)
res.values['Q'].append(prev_timestep, q)
count[0] = 4 + 8 + 1 + 8*len(z)
res.raw_data.append(count)
res.raw_data.append(np.array([len(z)],dtype=np.int32))
res.raw_data.append(np.array([prev_timestep],dtype=np.float64))
res.raw_data.append(np.array(bZ,dtype=np.byte))
res.raw_data.append(z)
res.raw_data.append(count)
res.raw_data.append(count)
res.raw_data.append(np.array([len(q)],dtype=np.int32))
res.raw_data.append(np.array([prev_timestep],dtype=np.float64))
res.raw_data.append(np.array(bQ,dtype=np.byte))
res.raw_data.append(q)
res.raw_data.append(count)
prev_timestep = t
z[isect-1] = float(ls[iz])
q[isect-1] = float(ls[iq])
if __name__ == "__main__":
current_dir = os.getcwd()
os.chdir('/home/theophile.terraz/Codes/mage/src')
spec = imp.spec_from_file_location("me", "./mage_extraire.py")
me = imp.module_from_spec(spec)
spec.loader.exec_module(me)
os.chdir(current_dir)
print( 'Number of arguments:', len(sys.argv), 'arguments.')
print( 'Argument List:', str(sys.argv))
res = me.data()
lire_opt(sys.argv[1], res)
res.write_bin_8(sys.argv[2][0:-4])