export CLI file with extract_contour
Summary
export CLI files with extract_contour.py
Why is this feature useful?
The current version of extract_contour does a great job for calculating boundaries on a selafin file.
A small implementation could greatly improve the python API and/or the inline capacities by also exporting the CLI file.
How to implement it?
Here is the current implementation of that feature in pyposeidon where we just add the index of the BND points, along with their location:
def export_cli(ds: xr.Dataset, tel_path: str, outCli: str, tel_module: str = "telemac2d"):
"""
(This function is a modification of the existing extract_contour() function
in scripts/python3/pretel/extract_contour.py of the TELEMAC scripts)
Generic function for extraction of contour from a mesh (with our without
boundary file)
@param ds (xr.Dataset) xarray Dataset of the mesh file (used to extract the boundary types)
@param geo (str) path to the Telemac geo file
@param outCli (str) Path to the output contour file
@returns (list) List of polygons
"""
# normalise
ds = ds.rename(
{
"SCHISM_hgrid_node_x": "node_x",
"SCHISM_hgrid_node_y": "node_y",
"SCHISM_hgrid_face_nodes": "face_nodes",
}
)
#
node_to_type = dict(zip(ds.bnodes.values, ds.type.values))
domains_bnd = []
lines = []
bnd_node = 0
contours, contours_idx = extract_contour(ds)
for bnd in contours_idx:
poly_bnd = []
coord_bnd = []
for i, glo_node in enumerate(bnd[:-1]): # not taking the last node (not repeating)
x, y = ds.node_x[glo_node], ds.node_y[glo_node]
coord_bnd.append((x, y))
# Determine boundary type for the current node
boundary_type = node_to_type.get(glo_node, "Unknown")
if boundary_type == "open":
# Get previous and next node indices in a circular manner
prev_node = bnd[i - 1] if i > 0 else bnd[-2] # -2 to skip the repeated last node
next_node = bnd[i + 1] if i < len(bnd) - 2 else bnd[0] # Wrapping around to the first node
# Get boundary types for previous and next nodes
prev_boundary_type = node_to_type.get(prev_node, "Unknown")
next_boundary_type = node_to_type.get(next_node, "Unknown")
# If both adjacent nodes are not 'open', then bnd is closed
if prev_boundary_type != "open" and next_boundary_type != "open":
boundary_type = "Unknown"
boundary_settings = get_boundary_settings(boundary_type, glo_node, bnd_node, tel_module)
keys_order = [
"lihbor",
"liubor",
"livbor",
"hbor",
"ubor",
"vbor",
"aubor",
"litbor",
"tbor",
"atbor",
"btbor",
"nbor",
"k",
]
line = " ".join(str(boundary_settings[key]) for key in keys_order)
lines.append(line)
bnd_node += 1
poly_bnd.append((coord_bnd, bnd))
domains_bnd.append(poly_bnd)
# Writing to file
with open(outCli, "w") as f:
for line in lines:
formatted_line = " ".join(
format_value(value, 3, is_float=isinstance(value, float)) for value in line.split()
)
f.write(f"{formatted_line}\n")
return domains_bndAdditional resources
Note the get_boundary_settings function that can be adapted for each TELEMAC module (currently using default T2D settings)
def get_boundary_settings(boundary_type, glo_node, bnd_node, module):
if module in ["telemac2d", "telemac3d", "tomawac"]:
settings = {
"lihbor": 5 if boundary_type == "open" else 2,
"liubor": 6 if boundary_type == "open" else 2,
"livbor": 6 if boundary_type == "open" else 2,
"hbor": 0.0,
"ubor": 0.0,
"vbor": 0.0,
"aubor": 0.0,
"litbor": 5 if boundary_type == "open" else 2,
"tbor": 0.0,
"atbor": 0.0,
"btbor": 0.0,
"nbor": glo_node + 1,
"k": bnd_node + 1,
}
return settings