Compare commits
No commits in common. "d4e9c3fd7840b9b627fbf9427d37a00145af8995" and "25e0f91bf0acce0ed5d16359bf50d6c864fc139a" have entirely different histories.
d4e9c3fd78
...
25e0f91bf0
19 changed files with 49 additions and 618 deletions
|
|
@ -1,29 +0,0 @@
|
|||
[inp]
|
||||
root=data
|
||||
base=Database_20220224.xyz
|
||||
hires=bathyhires.dat
|
||||
hstru=Hstru.dat
|
||||
poro=Poro.dat
|
||||
psize=Psize.dat
|
||||
raw_ts=cerema/raw/201702281700.raw
|
||||
raw_spec=cerema/spt/201702281715.spt
|
||||
hires_step=0.5
|
||||
cycle=14400
|
||||
|
||||
[out]
|
||||
margin=0.005
|
||||
#no_breakwater=True
|
||||
root=out2
|
||||
sub=bathy_sub.npy
|
||||
out=bathy.npy
|
||||
step=1
|
||||
left=-300
|
||||
right=150
|
||||
|
||||
[artha]
|
||||
lat=43.398450
|
||||
lon=-1.673097
|
||||
|
||||
[buoy]
|
||||
lat=43.408333
|
||||
lon=-1.681667
|
||||
|
|
@ -4,7 +4,6 @@ import logging
|
|||
import pathlib
|
||||
|
||||
import matplotlib.pyplot as plt
|
||||
from matplotlib.ticker import MultipleLocator
|
||||
import numpy as np
|
||||
|
||||
parser = argparse.ArgumentParser(description="Plot orbitals")
|
||||
|
|
@ -86,18 +85,16 @@ fig2dv, ax2dv = plt.subplots(figsize=(5/2.54, 2/3*10/2.54), dpi=200, constrained
|
|||
x0 = ts_flt["x"] * np.cos(theta) + ts_flt["y"] * np.sin(theta)
|
||||
#ax2dv.plot(x0, z0, c="#0066ff", lw=1)
|
||||
ax2dv.quiver(
|
||||
x0[:-1] * 1e-2,
|
||||
z0[:-1] * 1e-2,
|
||||
np.diff(x0)[:] * 1e-2,
|
||||
np.diff(z0)[:] * 1e-2,
|
||||
x0[:-1],
|
||||
z0[:-1],
|
||||
np.diff(x0)[:],
|
||||
np.diff(z0)[:],
|
||||
color="k",
|
||||
scale_units="xy",
|
||||
scale=1,
|
||||
)
|
||||
ax2dv.grid(c="k", alpha=.2)
|
||||
ax2dv.set(aspect="equal", xlabel="x (m)", ylabel="z (m)")
|
||||
ax2dv.set(ylim=(-10, 10))
|
||||
ax2dv.yaxis.set_minor_locator(MultipleLocator(1))
|
||||
ax2dv.set(aspect="equal", xlabel="x (cm)", ylabel="z (cm)")
|
||||
fig2dv.savefig("out_orbitals.pdf")
|
||||
fig2dv.savefig("out_orbitals.jpg")
|
||||
|
||||
|
|
|
|||
|
|
@ -131,12 +131,8 @@ np.savetxt(out_root.joinpath("hstru.dat"), hstru[::-1], newline=" ")
|
|||
np.savetxt(out_root.joinpath("poro.dat"), poro[::-1], newline=" ")
|
||||
np.savetxt(out_root.joinpath("psize.dat"), psize[::-1], newline=" ")
|
||||
|
||||
fig, ax = plt.subplots(figsize=(16 / 2.54, 2 / 3 * 10 / 2.54), constrained_layout=True)
|
||||
fig, ax = plt.subplots()
|
||||
ax.plot(-x, z, color="k")
|
||||
ax.fill_between(-x, z + hstru, z, color="k", alpha=0.2)
|
||||
#ax.set_title(f"N={z.size-1}, x=[{-x.max()};{-x.min()}]")
|
||||
ax.set(ylim=(-30, 15))
|
||||
ax.set(xlabel="x (m)", ylabel="z (m)")
|
||||
ax.autoscale(True, "x", True)
|
||||
ax.grid()
|
||||
ax.set_title(f"N={z.size-1}, x=[{-x.max()};{-x.min()}]")
|
||||
fig.savefig(out_root.joinpath("bathy.pdf"))
|
||||
|
|
|
|||
|
|
@ -4,7 +4,6 @@ import logging
|
|||
import pathlib
|
||||
|
||||
import matplotlib.pyplot as plt
|
||||
from matplotlib.ticker import MultipleLocator
|
||||
import numpy as np
|
||||
|
||||
parser = argparse.ArgumentParser(description="Pre-process time-series")
|
||||
|
|
@ -48,37 +47,7 @@ log.debug(f"{t=}")
|
|||
log.info(f"Saving timeseries to '{out_ts}'")
|
||||
np.savetxt(out_ts, np.stack((t, raw_ts["z"] / 100), axis=1))
|
||||
|
||||
fig, ax = plt.subplots(figsize=(8 / 2.54, 2 / 3 * 10 / 2.54), constrained_layout=True)
|
||||
tp = np.datetime64("2017-02-28T17:00:00") + t.astype(np.timedelta64)[-(t.size // 3) :]
|
||||
ax.plot(
|
||||
tp,
|
||||
raw_ts["z"][-(t.size // 3) :] * 1e-2,
|
||||
color="k",
|
||||
lw=1,
|
||||
)
|
||||
ax.axvline(
|
||||
np.datetime64("2017-02-28T17:00:00") + np.timedelta64(23 * 60 + 8),
|
||||
color="k",
|
||||
alpha=0.1,
|
||||
lw=20,
|
||||
)
|
||||
ax.autoscale(True, "x", True)
|
||||
ax.set(xlabel="t (s)", ylabel="z (m)")
|
||||
yabs_max = abs(max(ax.get_ylim(), key=abs))
|
||||
ax.set(ylim=(-10, 10))
|
||||
ax.set(
|
||||
xticks=(
|
||||
np.datetime64("2017-02-28T17:20:00"),
|
||||
np.datetime64("2017-02-28T17:25:00"),
|
||||
np.datetime64("2017-02-28T17:30:00"),
|
||||
),
|
||||
xticklabels=(
|
||||
"17:20",
|
||||
"17:25",
|
||||
"17:30",
|
||||
),
|
||||
)
|
||||
ax.yaxis.set_minor_locator(MultipleLocator(1))
|
||||
ax.grid(color="k", alpha=0.2)
|
||||
fig, ax = plt.subplots()
|
||||
ax.plot(t, raw_ts["z"])
|
||||
ax.set(xlabel="t (s)", ylabel="z (cm)")
|
||||
fig.savefig(out_root.joinpath("ts.pdf"))
|
||||
fig.savefig(out_root.joinpath("ts.jpg"), dpi=200)
|
||||
|
|
|
|||
|
|
@ -1,29 +0,0 @@
|
|||
/*--------------------------------*- C++ -*----------------------------------*\
|
||||
| ========= | |
|
||||
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
|
||||
| \\ / O peration | Version: 2.1.0 |
|
||||
| \\ / A nd | Web: www.OpenFOAM.org |
|
||||
| \\/ M anipulation | |
|
||||
\*---------------------------------------------------------------------------*/
|
||||
FoamFile
|
||||
{
|
||||
version 2.0;
|
||||
format ascii;
|
||||
class dictionary;
|
||||
location "constant";
|
||||
object porosityDict;
|
||||
}
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
// Materials: clear region, core, secondary armour layer, primary armour layer
|
||||
// a,b,c: tuning parameters
|
||||
a 2(0 50);
|
||||
b 2(0 1.2);
|
||||
c 2(0 0.34);
|
||||
|
||||
// D50: mean nominal diameter
|
||||
D50 2(1 4);
|
||||
// porosity (phi)
|
||||
porosity 2(1 0.4);
|
||||
|
||||
// ************************************************************************* //
|
||||
|
|
@ -1,29 +0,0 @@
|
|||
/*--------------------------------*- C++ -*----------------------------------*\
|
||||
| ========= | |
|
||||
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
|
||||
| \\ / O peration | Version: 2.1.0 |
|
||||
| \\ / A nd | Web: www.OpenFOAM.org |
|
||||
| \\/ M anipulation | |
|
||||
\*---------------------------------------------------------------------------*/
|
||||
FoamFile
|
||||
{
|
||||
version 2.0;
|
||||
format ascii;
|
||||
class dictionary;
|
||||
location "constant";
|
||||
object turbulenceProperties;
|
||||
}
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
simulationType RAS;
|
||||
|
||||
RAS
|
||||
{
|
||||
RASModel kOmegaSST;
|
||||
|
||||
turbulence on;
|
||||
|
||||
printCoeffs on;
|
||||
}
|
||||
|
||||
// ************************************************************************* //
|
||||
|
|
@ -1,81 +0,0 @@
|
|||
/*--------------------------------*- C++ -*----------------------------------*\ | ========= | |
|
||||
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
|
||||
| \\ / O peration | Version: 1.7.1 |
|
||||
| \\ / A nd | Web: www.OpenFOAM.com |
|
||||
| \\/ M anipulation | |
|
||||
\*---------------------------------------------------------------------------*/
|
||||
FoamFile
|
||||
{
|
||||
version 2.0;
|
||||
format ascii;
|
||||
class dictionary;
|
||||
object blockMeshDict;
|
||||
}
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
scale 1;
|
||||
|
||||
vertices
|
||||
(
|
||||
(-150 0 -30)
|
||||
(0 0 -30)
|
||||
(0 0 30)
|
||||
(-150 0 30)
|
||||
(-150 1 -30)
|
||||
(0 1 -30)
|
||||
(0 1 30)
|
||||
(-150 1 30)
|
||||
);
|
||||
|
||||
blocks
|
||||
(
|
||||
hex (0 1 5 4 3 2 6 7) (750 1 300) simpleGrading (1 1 1)
|
||||
);
|
||||
|
||||
edges
|
||||
(
|
||||
);
|
||||
|
||||
boundary
|
||||
(
|
||||
inlet
|
||||
{
|
||||
type patch;
|
||||
faces
|
||||
(
|
||||
(0 4 7 3)
|
||||
);
|
||||
}
|
||||
/*outlet
|
||||
{
|
||||
type patch;
|
||||
faces
|
||||
(
|
||||
(1 5 6 2)
|
||||
);
|
||||
}*/
|
||||
wall1
|
||||
{
|
||||
type wall;
|
||||
faces
|
||||
(
|
||||
(0 1 5 4)
|
||||
);
|
||||
}
|
||||
atmosphere
|
||||
{
|
||||
type patch;
|
||||
faces
|
||||
(
|
||||
(3 2 6 7)
|
||||
(1 5 6 2)
|
||||
);
|
||||
}
|
||||
);
|
||||
|
||||
mergePatchPairs
|
||||
(
|
||||
);
|
||||
|
||||
// ************************************************************************* //
|
||||
|
|
@ -1,138 +0,0 @@
|
|||
/*---------------------------------------------------------------------------*\
|
||||
| ========= | |
|
||||
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
|
||||
| \\ / O peration | Version: 1.3 |
|
||||
| \\ / A nd | Web: http://www.openfoam.org |
|
||||
| \\/ M anipulation | |
|
||||
\*---------------------------------------------------------------------------*/
|
||||
FoamFile
|
||||
{
|
||||
version 2.0;
|
||||
format ascii;
|
||||
location "system";
|
||||
class dictionary;
|
||||
object controlDict;
|
||||
}
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
application olaFlow;
|
||||
|
||||
startFrom latestTime;
|
||||
|
||||
startTime 0;
|
||||
|
||||
stopAt endTime;
|
||||
|
||||
endTime 400;
|
||||
|
||||
deltaT 0.1;
|
||||
|
||||
writeControl adjustableRunTime;
|
||||
|
||||
writeInterval 0.5;
|
||||
|
||||
purgeWrite 0;
|
||||
|
||||
writeFormat ascii;
|
||||
|
||||
writePrecision 6;
|
||||
|
||||
compression on;
|
||||
|
||||
timeFormat general;
|
||||
|
||||
timePrecision 6;
|
||||
|
||||
runTimeModifiable yes;
|
||||
|
||||
adjustTimeStep yes;
|
||||
|
||||
maxCo 0.45;
|
||||
maxAlphaCo 0.45;
|
||||
|
||||
maxDeltaT 0.5;
|
||||
/*
|
||||
functions
|
||||
{
|
||||
gaugesVOF
|
||||
{
|
||||
type sets;
|
||||
libs ("libsampling.so");
|
||||
writeControl outputTime;
|
||||
writeInterval 1;
|
||||
setFormat raw;
|
||||
surfaceFormat raw;
|
||||
interpolationScheme cell;
|
||||
fields ( alpha.water );
|
||||
sets
|
||||
(
|
||||
GaugeVOF01
|
||||
{
|
||||
type lineCellFace;
|
||||
axis xyz;
|
||||
start ( 0.5 0.001 0 );
|
||||
end ( 0.5 0.001 1.2 );
|
||||
}
|
||||
GaugeVOF02
|
||||
{
|
||||
type lineCellFace;
|
||||
axis xyz;
|
||||
start ( 9.25 0.001 0 );
|
||||
end ( 9.25 0.001 1.2 );
|
||||
}
|
||||
GaugeVOF03
|
||||
{
|
||||
type lineCellFace;
|
||||
axis xyz;
|
||||
start ( 15.75 0.001 0 );
|
||||
end ( 15.75 0.001 1.2 );
|
||||
}
|
||||
GaugeVOF04
|
||||
{
|
||||
type lineCellFace;
|
||||
axis xyz;
|
||||
start ( 17.75 0.001 0 );
|
||||
end ( 17.75 0.001 1.2 );
|
||||
}
|
||||
GaugeVOF05
|
||||
{
|
||||
type lineCellFace;
|
||||
axis xyz;
|
||||
start ( 21.1 0.001 0 );
|
||||
end ( 21.1 0.001 1.2 );
|
||||
}
|
||||
);
|
||||
}
|
||||
gaugesP
|
||||
{
|
||||
type sets;
|
||||
libs ("libsampling.so");
|
||||
writeControl outputTime;
|
||||
writeInterval 1;
|
||||
setFormat raw;
|
||||
surfaceFormat raw;
|
||||
interpolationScheme cellPointFace;
|
||||
fields ( p );
|
||||
sets
|
||||
(
|
||||
GaugesP
|
||||
{
|
||||
type boundaryPoints;
|
||||
axis xyz;
|
||||
patches 1(caisson);
|
||||
points ((18.0 0.01 0.75)
|
||||
(18.00 0.01 0.80)
|
||||
(18.00 0.01 0.85)
|
||||
(18.00 0.01 0.95)
|
||||
(18.01 0.01 0.70)
|
||||
(18.25 0.01 0.70)
|
||||
(18.50 0.01 0.70)
|
||||
(18.75 0.01 0.70));
|
||||
maxDistance 0.01;
|
||||
}
|
||||
);
|
||||
}
|
||||
}
|
||||
*/
|
||||
// ************************************************************************* //
|
||||
|
|
@ -1,25 +0,0 @@
|
|||
/*--------------------------------*- C++ -*----------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration | Website: https://openfoam.org
|
||||
\\ / A nd | Version: 9
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
Description
|
||||
Writes graph data for specified fields along a line, specified by start and
|
||||
end points. A specified number of graph points are used, distributed
|
||||
uniformly along the line.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
start (-50 0.5 -15);
|
||||
end (-50 0.5 15);
|
||||
nPoints 100;
|
||||
|
||||
fields (alpha.water U);
|
||||
|
||||
axis z;
|
||||
|
||||
#includeEtc "caseDicts/postProcessing/graphs/graphUniform.cfg"
|
||||
|
||||
// ************************************************************************* //
|
||||
|
|
@ -1,25 +0,0 @@
|
|||
/*--------------------------------*- C++ -*----------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration | Website: https://openfoam.org
|
||||
\\ / A nd | Version: 9
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
Description
|
||||
Writes graph data for specified fields along a line, specified by start and
|
||||
end points. A specified number of graph points are used, distributed
|
||||
uniformly along the line.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
start (-20 0.5 -15);
|
||||
end (-20 0.5 15);
|
||||
nPoints 100;
|
||||
|
||||
fields (alpha.water U);
|
||||
|
||||
axis z;
|
||||
|
||||
#includeEtc "caseDicts/postProcessing/graphs/graphUniform.cfg"
|
||||
|
||||
// ************************************************************************* //
|
||||
|
|
@ -8,7 +8,6 @@ import pickle
|
|||
import matplotlib.pyplot as plt
|
||||
import matplotlib.animation as animation
|
||||
from matplotlib.gridspec import GridSpec
|
||||
from matplotlib.ticker import MultipleLocator
|
||||
import numpy as np
|
||||
from scipy import interpolate
|
||||
|
||||
|
|
@ -29,12 +28,6 @@ parser.add_argument(
|
|||
help="Only compute maximum rather than animation",
|
||||
action="store_true",
|
||||
)
|
||||
parser.add_argument(
|
||||
"-i",
|
||||
"--initial",
|
||||
help="Only compute initial domain",
|
||||
action="store_true",
|
||||
)
|
||||
args = parser.parse_args()
|
||||
|
||||
logging.basicConfig(level=max((10, 20 - 10 * args.verbose)))
|
||||
|
|
@ -68,19 +61,11 @@ AW[:, iz0[idz0], ix0[idx0]] = model.fields["alpha.water"]
|
|||
U = np.full((model.t.size, *X.shape), np.nan)
|
||||
U[:, iz0[idz0], ix0[idx0]] = np.linalg.norm(model.fields["U"], axis=1)
|
||||
|
||||
fig = plt.figure(
|
||||
figsize=(15 / 2.54, 4 / 2.54), dpi=200, constrained_layout=True
|
||||
)
|
||||
gs = GridSpec(
|
||||
3 if not args.initial else 1,
|
||||
1,
|
||||
figure=fig,
|
||||
height_ratios=[1, 0.1, 0.1] if not args.initial else [1],
|
||||
)
|
||||
fig = plt.figure()
|
||||
gs = GridSpec(3, 1, figure=fig, height_ratios=[1, 0.05, 0.05])
|
||||
ax = fig.add_subplot(gs[0])
|
||||
if not args.initial:
|
||||
cax1 = fig.add_subplot(gs[1])
|
||||
cax2 = fig.add_subplot(gs[2])
|
||||
cax1 = fig.add_subplot(gs[1])
|
||||
cax2 = fig.add_subplot(gs[2])
|
||||
aw_m = ax.imshow(
|
||||
AW[0],
|
||||
vmin=0,
|
||||
|
|
@ -101,15 +86,10 @@ p_m = ax.imshow(
|
|||
ax.axhline(4.5, ls="-.", lw=1, c="k", alpha=0.2, zorder=1.2)
|
||||
|
||||
|
||||
if not args.initial:
|
||||
fig.colorbar(
|
||||
aw_m, label=r"$\alpha_w$", cax=cax1, shrink=0.6, orientation="horizontal"
|
||||
)
|
||||
fig.colorbar(p_m, label=r"Porosity", cax=cax2, shrink=0.6, orientation="horizontal")
|
||||
fig.colorbar(aw_m, label=r"$\alpha_w$", cax=cax1, shrink=0.6, orientation="horizontal")
|
||||
fig.colorbar(p_m, label=r"Porosity", cax=cax2, shrink=0.6, orientation="horizontal")
|
||||
ax.set(xlabel="x (m)", ylabel="z (m)", aspect="equal", facecolor="#000000")
|
||||
ax.grid(c="k", alpha=0.2)
|
||||
ax.xaxis.set_minor_locator(MultipleLocator(5))
|
||||
ax.yaxis.set_minor_locator(MultipleLocator(5))
|
||||
|
||||
|
||||
figU = plt.figure()
|
||||
|
|
@ -157,13 +137,6 @@ if args.max:
|
|||
|
||||
fig.savefig(out.joinpath("max_aw.pdf"))
|
||||
figU.savefig(out.joinpath("max_U.pdf"))
|
||||
elif args.initial:
|
||||
aw_m.set_array(AW[0])
|
||||
ax.vlines(-20, -15, 15, color="k", lw=1, ls="--", label="Measurements")
|
||||
ax.text(-20, 15, "Measurements", ha="right", va="bottom")
|
||||
|
||||
fig.savefig(out.joinpath("aw_t0.pdf"))
|
||||
fig.savefig(out.joinpath("aw_t0.jpg"), dpi=200)
|
||||
else:
|
||||
fig.set(figwidth=19.2, figheight=10.8, dpi=100)
|
||||
figU.set(figwidth=19.2, figheight=10.8, dpi=100)
|
||||
|
|
|
|||
|
|
@ -74,7 +74,7 @@ u_m = axU.quiver(
|
|||
UU[0],
|
||||
alpha=alp[0],
|
||||
cmap="inferno_r",
|
||||
clim=(0, 20),
|
||||
clim=(0, np.nanmax(UU)),
|
||||
)
|
||||
# (wat_p,) = axU.plot(x0, watl[0])
|
||||
|
||||
|
|
|
|||
|
|
@ -9,7 +9,7 @@ import sys
|
|||
|
||||
from cycler import cycler
|
||||
import matplotlib.pyplot as plt
|
||||
from matplotlib.ticker import MultipleLocator
|
||||
import matplotlib.gridspec as gridspec
|
||||
import numpy as np
|
||||
from scipy import interpolate
|
||||
|
||||
|
|
@ -67,12 +67,11 @@ def get_pickle(out):
|
|||
|
||||
models = list(map(get_pickle, args.output))
|
||||
|
||||
figsize = 15 / 2.54, 4 / 2.54 * len(models)
|
||||
|
||||
fig, ax_ = plt.subplots(
|
||||
len(models),
|
||||
figsize=figsize,
|
||||
dpi=200,
|
||||
figsize=(6, 1.5 * len(models)),
|
||||
dpi=100,
|
||||
constrained_layout=True,
|
||||
squeeze=False,
|
||||
)
|
||||
|
|
@ -91,71 +90,32 @@ if args.timestep is None:
|
|||
)
|
||||
case "U":
|
||||
for i, (_ax, _model) in enumerate(zip(ax, models)):
|
||||
v = np.nanmax(np.abs(np.where(
|
||||
_model.post_fields[args.func]["alpha.water"] > 0.5,
|
||||
#np.linalg.norm(_model.post_fields[args.func][args.field], axis=2),
|
||||
_model.post_fields[args.func][args.field][..., 0],
|
||||
np.nan,
|
||||
)))
|
||||
v150 = np.nanmax(np.abs(np.where(
|
||||
(_model.post_fields[args.func]["alpha.water"] > 0.5) & (_model.t[:, None] > 170) & (_model.t[:, None] < 200),
|
||||
#np.linalg.norm(_model.post_fields[args.func][args.field], axis=2),
|
||||
_model.post_fields[args.func][args.field][..., 0],
|
||||
np.nan,
|
||||
)))
|
||||
_data = _model.post_fields[args.func][args.field][..., 0].T
|
||||
#_c = _ax.contourf(
|
||||
# _model.t,
|
||||
# _model.post_fields[args.func][f"x_{args.field}"],
|
||||
# _data,
|
||||
# cmap="PiYG",
|
||||
# #levels=[-15, -10, -5, -2, -1, 0, 1, 2, 5, 10, 15],
|
||||
# vmin=-np.nanmax(np.abs(_data)),
|
||||
# vmax=np.nanmax(np.abs(_data)),
|
||||
# extend="both",
|
||||
#)
|
||||
_c = _ax.imshow(
|
||||
_data[::-1],
|
||||
cmap="PiYG",
|
||||
alpha=np.clip(_model.post_fields[args.func]["alpha.water"], 0, 1).T[::-1],
|
||||
np.where(_model.post_fields[args.func]["alpha.water"] > 0.5, np.linalg.norm(_model.post_fields[args.func][args.field], axis=2), np.nan).T[::-1],
|
||||
vmin=0,
|
||||
vmax=20,
|
||||
cmap="inferno_r",
|
||||
extent=(
|
||||
_model.t.min(),
|
||||
_model.t.max(),
|
||||
_model.post_fields[args.func][f"x_{args.field}"].min(),
|
||||
_model.post_fields[args.func][f"x_{args.field}"].max(),
|
||||
),
|
||||
vmin=-v150,
|
||||
vmax=v150,
|
||||
aspect="auto",
|
||||
)
|
||||
_ax.set(xlim=(100, 300))
|
||||
_ax.set(facecolor="k")
|
||||
_ax.xaxis.set_minor_locator(MultipleLocator(5))
|
||||
_ax.yaxis.set_minor_locator(MultipleLocator(1))
|
||||
fig.colorbar(_c, label=f"{args.field} (m/s)", ax=_ax)
|
||||
log.info(f"Vitesse max: {v}m/s")
|
||||
log.info(f"Vitesse max [170,200]: {v150}m/s")
|
||||
log.info(f"Écart: {abs(np.nanmax(_data)-17.7)/17.7:%}")
|
||||
case _:
|
||||
log.error(f"Cannot plot field {args.field} from {args.func}")
|
||||
sys.exit(1)
|
||||
|
||||
for i, (_ax, _model) in enumerate(zip(ax, models)):
|
||||
_ax.set(xlabel="t (s)", ylabel="z (m)")
|
||||
if len(models) > 1:
|
||||
_ax.set(title=f"Case {i}")
|
||||
#_ax.grid(color="#797979", alpha=0.5)
|
||||
_ax.set(xlabel="t (s)", ylabel="z (m)", title=f"Case {i}")
|
||||
_ax.grid(color="k", alpha=0.2)
|
||||
|
||||
fig.savefig(
|
||||
args.output[0].joinpath(
|
||||
f"diff_{args.func}_{args.field}_{'_'.join([o.name for o in args.output])}.pdf"
|
||||
)
|
||||
)
|
||||
fig.savefig(
|
||||
args.output[0].joinpath(
|
||||
f"diff_{args.func}_{args.field}_{'_'.join([o.name for o in args.output])}.jpg"
|
||||
)
|
||||
)
|
||||
else:
|
||||
match args.field:
|
||||
case "alpha.water":
|
||||
|
|
@ -163,9 +123,7 @@ else:
|
|||
_ax.tricontour(
|
||||
_model.x,
|
||||
_model.z,
|
||||
_model.fields[args.field][np.where(_model.t == args.timestep)[0]][
|
||||
0
|
||||
],
|
||||
_model.fields[args.field][np.where(_model.t == args.timestep)[0]][0],
|
||||
levels=(0.5,),
|
||||
colors="k",
|
||||
)
|
||||
|
|
|
|||
|
|
@ -10,5 +10,5 @@ mpi=8
|
|||
[post]
|
||||
inp=inp_post/ts_4lay_1h
|
||||
out=out_post/ts_4lay_1h
|
||||
x0=-50
|
||||
x0=-1250
|
||||
t0=180
|
||||
|
|
|
|||
|
|
@ -9,7 +9,7 @@ mpi=8
|
|||
|
||||
[post]
|
||||
inp=inp_post/real_spec_interp
|
||||
compare=inp_post/real_spec_interp_nb
|
||||
#compare=inp_post/real_spec_interp_nb
|
||||
out=out_post/real_spec_interp
|
||||
x0=-1250
|
||||
t0=180
|
||||
|
|
|
|||
|
|
@ -7,6 +7,8 @@ import matplotlib.pyplot as plt
|
|||
import numpy as np
|
||||
import scipy.signal as sgl
|
||||
|
||||
from .read_swash import *
|
||||
|
||||
parser = argparse.ArgumentParser(description="Post-process swash output")
|
||||
parser.add_argument("-v", "--verbose", action="count", default=0)
|
||||
parser.add_argument("-c", "--config", default="config.ini")
|
||||
|
|
@ -28,8 +30,7 @@ t = np.load(inp.joinpath("t.npy"))
|
|||
|
||||
botl = np.load(inp.joinpath("botl.npy"))
|
||||
watl = np.load(inp.joinpath("watl.npy"))
|
||||
#vel = np.load(inp.joinpath("vel_x.npy"))
|
||||
vel = np.load(inp.joinpath("vel.npy"))[0]
|
||||
vel = np.load(inp.joinpath("vel_x.npy"))
|
||||
|
||||
# Cospectral calculations
|
||||
x0 = config.getint("post", "x0")
|
||||
|
|
@ -64,25 +65,22 @@ R = np.sqrt(
|
|||
|
||||
if config.has_option("post", "compare"):
|
||||
inp_comp = pathlib.Path(config.get("post", "compare"))
|
||||
x_ = np.load(inp_comp.joinpath("x.npy"))
|
||||
t_ = np.load(inp_comp.joinpath("t.npy"))
|
||||
x_ = np.load(inp_comp.joinpath("xp.npy"))
|
||||
t_ = np.load(inp_comp.joinpath("tsec.npy"))
|
||||
|
||||
botl_ = np.load(inp_comp.joinpath("botl.npy"))
|
||||
watl_ = np.load(inp_comp.joinpath("watl.npy"))
|
||||
vel_ = np.load(inp_comp.joinpath("vel_x.npy"))
|
||||
|
||||
# Cospectral calculations
|
||||
arg_x0_ = np.abs(x_ - x0).argmin()
|
||||
arg_t0_ = np.abs(t_ - t0).argmin()
|
||||
dt_ = np.diff(t_).mean() * 1e-3
|
||||
f_ = 1 / dt_
|
||||
|
||||
eta_ = watl_[t_ > t0, arg_x0_]
|
||||
u_ = vel_[t_ > t0, arg_x0_]
|
||||
|
||||
phi_eta_ = sgl.welch(eta_, f_, nperseg=nperseg)
|
||||
phi_u_ = sgl.welch(u_, f_, nperseg=nperseg)
|
||||
phi_eta_u_ = sgl.csd(eta_, u_, f_, nperseg=nperseg)
|
||||
phi_eta_ = sgl.welch(eta_, f, nperseg=nperseg)
|
||||
phi_u_ = sgl.welch(u_, f, nperseg=nperseg)
|
||||
phi_eta_u_ = sgl.csd(eta_, u_, f, nperseg=nperseg)
|
||||
|
||||
H_ = np.sqrt(np.abs(phi_eta_[1]))
|
||||
U_ = np.sqrt(np.abs(phi_u_[1]))
|
||||
|
|
@ -93,6 +91,10 @@ if config.has_option("post", "compare"):
|
|||
(np.abs(phi_eta_[1]) + np.abs(phi_u_[1]) - 2 * phi_eta_u_[1].real)
|
||||
/ (np.abs(phi_eta_[1]) + np.abs(phi_u_[1]) + 2 * phi_eta_u_[1].real)
|
||||
)
|
||||
# R_ = np.sqrt(
|
||||
# (1 + G_**2 - 2 * G_ * np.cos(th_eta_u_))
|
||||
# / (1 + G_**2 + 2 * G_ * np.cos(th_eta_u_))
|
||||
# )
|
||||
|
||||
|
||||
# Plotting
|
||||
|
|
|
|||
|
|
@ -30,7 +30,7 @@ botl = np.load(inp.joinpath("botl.npy"))
|
|||
watl = np.load(inp.joinpath("watl.npy"))
|
||||
vel = np.load(inp.joinpath("vel.npy"))[0]
|
||||
|
||||
t0 = np.linspace(23 * 60 + 8, 23 * 60 + 8 + 100, 6)
|
||||
t0 = np.linspace(23 * 60 + 8, 23 * 60 + 8 + 100, 5)
|
||||
|
||||
# Plotting
|
||||
log.info("Plotting results")
|
||||
|
|
@ -38,10 +38,10 @@ log.info("Plotting results")
|
|||
vlim = np.nanmin(np.maximum(watl, -botl)), np.nanmax(np.maximum(watl, -botl))
|
||||
|
||||
fig_x, ax = plt.subplots(
|
||||
3, 2, figsize=(15 / 2.54, 2 / 3 * 10 / 2.54), constrained_layout=True
|
||||
len(t0), figsize=(10 / 2.54, 4 / 3 * 10 / 2.54), constrained_layout=True
|
||||
)
|
||||
i0 = np.argmin(np.abs(t * 1e-3 - t0[0]))
|
||||
for ax_x, t0_x in zip(ax.reshape(-1), t0):
|
||||
for ax_x, t0_x in zip(ax, t0):
|
||||
ax_x.plot(x, -botl, color="k")
|
||||
i = np.argmin(np.abs(t * 1e-3 - t0_x))
|
||||
ax_x.plot(
|
||||
|
|
@ -67,6 +67,5 @@ log.info(f"Saving plots in '{out}'")
|
|||
out.mkdir(parents=True, exist_ok=True)
|
||||
|
||||
fig_x.savefig(out.joinpath("x.pdf"))
|
||||
fig_x.savefig(out.joinpath("x.jpg"), dpi=200)
|
||||
|
||||
log.info("Finished post-processing")
|
||||
|
|
|
|||
|
|
@ -1,87 +0,0 @@
|
|||
import argparse
|
||||
import configparser
|
||||
import logging
|
||||
import pathlib
|
||||
|
||||
import matplotlib.pyplot as plt
|
||||
import numpy as np
|
||||
import scipy.signal as sgl
|
||||
|
||||
parser = argparse.ArgumentParser(description="Post-process swash output")
|
||||
parser.add_argument("-v", "--verbose", action="count", default=0)
|
||||
parser.add_argument("-c", "--config", default="config.ini")
|
||||
args = parser.parse_args()
|
||||
|
||||
logging.basicConfig(level=max((10, 20 - 10 * args.verbose)))
|
||||
log = logging.getLogger("post")
|
||||
|
||||
log.info("Starting post-processing")
|
||||
config = configparser.ConfigParser()
|
||||
config.read(args.config)
|
||||
|
||||
inp = pathlib.Path(config.get("post", "inp"))
|
||||
root = pathlib.Path(config.get("swash", "out"))
|
||||
|
||||
log.info(f"Reading data from '{inp}'")
|
||||
x = np.load(inp.joinpath("x.npy"))
|
||||
t = np.load(inp.joinpath("t.npy")) * 1e-3
|
||||
|
||||
botl = np.load(inp.joinpath("botl.npy"))
|
||||
watl = np.load(inp.joinpath("watl.npy"))
|
||||
vel = np.load(inp.joinpath("vel.npy"))[0]
|
||||
|
||||
# Plotting
|
||||
log.info("Plotting results")
|
||||
|
||||
vlim = np.nanmin(np.maximum(watl, -botl)), np.nanmax(np.maximum(watl, -botl))
|
||||
|
||||
x0 = np.linspace(-600, -200, 5)
|
||||
i0 = np.argmin(np.abs(x[:, None] - x0), axis=0)
|
||||
|
||||
fig_x, ax = plt.subplots(
|
||||
5, 1, figsize=(15 / 2.54, 15/ 2.54), constrained_layout=True
|
||||
)
|
||||
dt = np.mean(np.diff(t))
|
||||
N = t.size
|
||||
s0 = 2 * dt
|
||||
dj = 0.5
|
||||
J = 1 / dj * np.log2(N * dt / s0)
|
||||
j = np.arange(0, J)
|
||||
sj = s0 * 2 ** (j * dj)
|
||||
Mw = sj / dt
|
||||
sig = np.var(watl[:, i0])
|
||||
M = np.stack([(np.abs(sgl.cwt(watl[:, i], sgl.morlet2, Mw))/sig)**2 for i in i0])
|
||||
v = np.max(M)
|
||||
|
||||
for ax_x, M_, x_ in zip(ax.reshape(-1), M, x[i0]):
|
||||
c = ax_x.contourf(t, sj, M_, cmap="Greys", vmin=0, levels=[1, 2.5, 5, 10, 20, 40], extend="both")
|
||||
fig_x.colorbar(c, ax=ax_x, label="NWPS")
|
||||
ax_x.grid(color="k", alpha=0.2)
|
||||
ax_x.text(
|
||||
0.95,
|
||||
0.95,
|
||||
f"x={x_:.0f}m",
|
||||
horizontalalignment="right",
|
||||
verticalalignment="top",
|
||||
transform=ax_x.transAxes,
|
||||
#c="w",
|
||||
)
|
||||
ax_x.semilogy()
|
||||
ax_x.autoscale(True, "both", True)
|
||||
ax_x.set_rasterization_zorder(1.5)
|
||||
ax_x.set(ylabel="T (s)", ylim=(sj[0], sj[-1]))
|
||||
|
||||
if ax_x != ax.reshape(-1)[-1]:
|
||||
ax_x.axes.set_xticklabels([])
|
||||
else:
|
||||
ax_x.set(xlabel="t (s)")
|
||||
|
||||
out = pathlib.Path(config.get("post", "out")).joinpath(f"trans")
|
||||
log.info(f"Saving plots in '{out}'")
|
||||
out.mkdir(parents=True, exist_ok=True)
|
||||
|
||||
fig_x.savefig(out.joinpath("wavelet.pdf"), dpi=300)
|
||||
fig_x.savefig(out.joinpath("wavelet.jpg"), dpi=200)
|
||||
fig_x.savefig(out.joinpath("wavelet.png"), dpi=200)
|
||||
|
||||
log.info("Finished post-processing")
|
||||
|
|
@ -38,9 +38,6 @@ w0_comp = watl_comp[:, arg_x0]
|
|||
cr0 = np.where(np.diff(np.sign(w0)))[0]
|
||||
cr0_comp = np.where(np.diff(np.sign(w0_comp)))[0]
|
||||
|
||||
log.info(f"1: {cr0.size} waves")
|
||||
log.info(f"2: {cr0_comp.size} waves")
|
||||
|
||||
wave = np.fromiter(
|
||||
(
|
||||
np.abs(
|
||||
|
|
@ -75,22 +72,9 @@ ax.autoscale(True, "x", True)
|
|||
ax.grid()
|
||||
fig.savefig(out.joinpath("wsize.pdf"))
|
||||
|
||||
fig2, ax2 = plt.subplots(
|
||||
figsize=(10 / 2.54, 2 / 3 * 10 / 2.54), constrained_layout=True
|
||||
)
|
||||
ax2.plot(
|
||||
t[cr0[i0 - 5] : cr0[i0 + 7]] * 1e-3,
|
||||
w0[cr0[i0 - 5] : cr0[i0 + 7]],
|
||||
color="k",
|
||||
label="Cas 1",
|
||||
)
|
||||
ax2.plot(
|
||||
t[cr0[i0 - 5] : cr0[i0 + 7]] * 1e-3,
|
||||
w0_comp[cr0[i0 - 5] : cr0[i0 + 7]],
|
||||
color="k",
|
||||
ls="-.",
|
||||
label="Cas 2",
|
||||
)
|
||||
fig2, ax2 = plt.subplots(figsize=(10/2.54, 2/3*10/2.54), constrained_layout=True)
|
||||
ax2.plot(t[cr0[i0 - 5] : cr0[i0 + 7]] * 1e-3, w0[cr0[i0 - 5] : cr0[i0 + 7]], color="k", label="Case 1")
|
||||
ax2.plot(t[cr0[i0 - 5] : cr0[i0 + 7]] * 1e-3, w0_comp[cr0[i0 - 5] : cr0[i0 + 7]], color="k", ls="-.", label="Case 2")
|
||||
ax2.set(xlabel="t (s)", ylabel="z (m)")
|
||||
ax2.autoscale(True, "x", True)
|
||||
ax2.grid()
|
||||
|
|
@ -98,11 +82,7 @@ ax2.legend()
|
|||
fig2.savefig(out.joinpath("maxw.pdf"))
|
||||
fig2.savefig(out.joinpath("maxw.jpg"), dpi=200)
|
||||
|
||||
log.info(
|
||||
f"RMS difference: {np.sqrt(np.mean((w0_comp-w0)**2))}m ; {np.sqrt(np.mean((w0_comp-w0)**2))/(w0.max()-w0.min()):%}"
|
||||
)
|
||||
log.info(f"RMS difference: {np.sqrt(np.mean((w0_comp-w0)**2))}m ; {np.sqrt(np.mean((w0_comp-w0)**2))/(w0.max()-w0.min()):%}")
|
||||
log.info(f"Bias: {np.mean(w0_comp-w0)}m")
|
||||
log.info(f"Maximum wave size: {wave.max()}m ; {wave_comp.max()}m")
|
||||
log.info(
|
||||
f"Maximum wave size difference: {abs(wave_comp.max()-wave.max())/wave.max():%}"
|
||||
)
|
||||
log.info(f"Maximum wave size difference: {abs(wave_comp.max()-wave.max())/wave.max():%}")
|
||||
|
|
|
|||
Loading…
Reference in a new issue