From b3d3fe670ed133d09e7e8a340326a0919442ae0e Mon Sep 17 00:00:00 2001 From: HenryTux <h.schoeller@fu-berlin.de> Date: Wed, 13 Dec 2023 17:27:11 +0100 Subject: [PATCH] updates --- calc_E.py | 9 +++------ calc_bound.py | 6 ++---- calc_dist.py | 3 --- main.py | 34 +++++++++++++++++++--------------- 4 files changed, 24 insertions(+), 28 deletions(-) diff --git a/calc_E.py b/calc_E.py index 142fccc..6846bdd 100644 --- a/calc_E.py +++ b/calc_E.py @@ -16,12 +16,9 @@ chdir(dirname(abspath(__file__))) from sys import argv, path import numpy as np path.append("../wp21/pyscripts/LIB") -import plot as pp import calc as cc import data as dd -import matplotlib.pyplot as plt from datetime import datetime -import itertools # %% Get the data @@ -37,7 +34,7 @@ fname = "traj" + date_0.strftime('_%Y%m%d_%H') + ".npy" trajs = np.load(trapath + fname) -t_sel = np.arange(0,trajs.shape[1]) +t_sel = np.arange(0,int(trajs.shape[1]/2)+1) z_coord_type = "p" @@ -61,13 +58,13 @@ T = trajs['T'][::n] H = cc.ptoh(P, 1013.25, 8.435) # standard atmosphere -k_p = cc.calc_k(U/1000, V/1000, Omg/100) +# k_p = cc.calc_k(U/1000, V/1000, Omg/100) k_p=15 # k_h = cc.calc_k(U, V, cc.omg2w(Omg/100, T, P)) D = list((list(), list(), list())) -epsilon_opt = np.array([25000, 50000, 100000, 250000, 500000]) +epsilon_opt = np.array([20000, 50000, 100000, 200000, 500000]) bound_meth = "$\\alpha$" alpha = 1e-3 diff --git a/calc_bound.py b/calc_bound.py index c69d0a0..21ba9ef 100644 --- a/calc_bound.py +++ b/calc_bound.py @@ -16,10 +16,8 @@ chdir(dirname(abspath(__file__))) from sys import argv, path, stdout import numpy as np path.append("../wp21/pyscripts/LIB") -import plot as pp import calc as cc import data as dd -import matplotlib.pyplot as plt from datetime import datetime import itertools @@ -81,8 +79,8 @@ elif coord == "x, y, p (stereo)": x, y, z = cc.coord_trans(Lon, Lat, P, "None", 1, proj="stereo") d_path = path + "stereop" if (not norm): - z = z * cc.calc_k(trajs['U']/1000, trajs['V']/1000, - trajs['OMEGA']/100) + z = z * 15 #cc.calc_k(trajs['U']/1000, trajs['V']/1000, + # trajs['OMEGA']/100) elif coord == "x, y, h (stereo)": x, y, z = cc.coord_trans(Lon, Lat, P, "std_atm", 1, proj="stereo") d_path = path + "stereoh" diff --git a/calc_dist.py b/calc_dist.py index 8a7e3a8..0e0131a 100644 --- a/calc_dist.py +++ b/calc_dist.py @@ -16,12 +16,9 @@ chdir(dirname(abspath(__file__))) from sys import argv, path import numpy as np path.append("../wp21/pyscripts/LIB") -import plot as pp import calc as cc import data as dd -import matplotlib.pyplot as plt from datetime import datetime -import itertools # %% Get the data diff --git a/main.py b/main.py index e4bebf7..1d40f41 100644 --- a/main.py +++ b/main.py @@ -32,11 +32,14 @@ n = 1 force_calc = False -date_0 = datetime(2016,4,30,12) +date_0 = datetime(2017,1,25,0) fname = "../wp21/era5/traj/" + date_0.strftime('%Y/traj_%Y%m%d_%H') + ".npy" B = dd.get_block_dat(date_0) -coasts = dd.get_coasts(0, 90, -180, 180) +if date_0.strftime("%Y") == "2016": + coasts = dd.get_coasts(10, 90, -60, -240) +elif date_0.strftime("%Y") == "2017": + coasts = dd.get_coasts(10, 90, 120, -60) #fname = "traj" + date_0.strftime('_%Y%m%d_%H') + ".npy" path = "./" + date_0.strftime("dat_%Y%m%d_%H") + "/" + str(n) + "/" @@ -46,7 +49,7 @@ if not exists(path): trajs = np.load(fname) if argv[0]=='': - direc = 'whole' + direc = 'backward' else: direc = argv[1] @@ -59,8 +62,8 @@ elif direc=="forward": plot_0 = 0 elif direc=="backward": - t_sel = np.arange(0,int(trajs.shape[1]/2)) - plot_0 = -1 + t_sel = np.arange(0,int(trajs.shape[1]/2)+1) + plot_0 = int(trajs.shape[1]/2) # t_sel = np.arange(73, 74) @@ -86,14 +89,14 @@ elev = 45 dist = 7 t_i = 0 alpha = 1e-3 -e = 100000 -epsilon_opt = np.array([25000, 50000, 100000, 250000, 500000]) +e = 20000 +epsilon_opt = np.array([20000, 50000, 100000, 200000, 500000]) alph_opt = np.array([1e-4, 2e-4, 5e-4, 0.001, 0.002, 0.005, 0.01, 0.05, 0.1, 0.5, 1, 5, 10]) -N_k = 3 +N_k = 6 # u and v is in m/s, omg in P/s; result is in km/hPa -k_p = cc.calc_k(trajs['U']/1000, trajs['V']/1000, trajs['OMEGA']/100) - +# k_p = cc.calc_k(trajs['U']/1000, trajs['V']/1000, trajs['OMEGA']/100) +k_p = 15 # Coordinates for visualization X = lon Y = lat @@ -103,6 +106,7 @@ k_h = cc.calc_k(trajs['U'], trajs['V'], cc.omg2w(trajs['OMEGA']/100, trajs['T'], trajs['p'])) # %% Boundary + # if False: if argv[0] == '' or argv[1] == "Boundary": @@ -307,11 +311,11 @@ if True: kclust = cc.kcluster_idx(E, N_k) # upper case letters for visualization - Xcs, Ycs, Zcs = cc.coord_trans(lon, lat, p, "None", k_p, proj="stereo") + Xcs, Ycs, Zcs = cc.coord_trans(lon, lat, p, "None", 1, proj="stereo") points = pp.plot_clustering(ax3dcs, Xcs, Ycs, Zcs, azim, elev, dist, t_i, c=kclust.labels_, plot_0=False, - coord="stereo", mean_traj=True, block=B, - coast=coasts) + coord="stereo", mean_traj=True, block=None, + coast=None) ax3dcs.invert_zaxis() # Add widget items @@ -338,7 +342,7 @@ if True: null_ax = fig_cs.add_axes([0.02, 0.43, 0.1, 0.05]) null_check = CheckButtons(null_ax, ["Plot start?"], [False]) geo_ax = fig_cs.add_axes([0.1, 0.43, 0.1, 0.05]) - geo_check = CheckButtons(geo_ax, ["Plot features?"], [True]) + geo_check = CheckButtons(geo_ax, ["Plot features?"], [False]) run_ax = fig_cs.add_axes([0.02, 0.35, 0.1, 0.05]) run_check2 = Button(run_ax, "Calculate") @@ -474,7 +478,7 @@ if True: ax_spag = fig_spag.add_subplot(111) fig_spag.subplots_adjust(left=0.2) var_ax = fig_spag.add_axes([0.025, 0.1, 0.1, 0.8]) - var_check = RadioButtons(var_ax, trajs.dtype.names, 3) + var_check = RadioButtons(var_ax, trajs.dtype.names, 6) var_ax.set_title("Variable") dat = trajs[var_check.value_selected][::n, t_sel] -- GitLab