Skip to content
21 changes: 18 additions & 3 deletions src/grid/zgr.py
Original file line number Diff line number Diff line change
Expand Up @@ -212,7 +212,6 @@ def fill_zgrid_vars(zgr_type, grid, hgr_type, e_dict, missing):
grid (dict) : vertical grid data dictionary
"""
# gdep

t_done = "gdept" not in missing
if t_done is False:
# Fill in the 3D gdept data from 1D gdept_0
Expand Down Expand Up @@ -248,10 +247,25 @@ def fill_zgrid_vars(zgr_type, grid, hgr_type, e_dict, missing):
if w_done is False:
if "e3w" not in missing:
grid["gdepw"] = np.cumsum(grid["e3w"], axis=1)
elif zgr_type == "zco":
# Don't waste time calculating gdepw everywhere if it all the same.
gdepw_0 = calc_gdepw(
grid["gdept_0"][np.newaxis, :, np.newaxis, np.newaxis]
).squeeze()
grid["gdepw"] = np.tile(
gdepw_0,
(
grid["mbathy"].shape[0],
grid["mbathy"].shape[1],
grid["mbathy"].shape[2],
1,
),
)
grid["gdepw"] = np.transpose(grid["gdepw"], axes=[0, 3, 1, 2])
else:
# Fill in the 3D gdepw data from gdept
grid["gdepw"] = calc_gdepw(grid["gdept"])
missing = sorted(list(set(missing) - set(["gdepw"])))
missing = sorted(list(set(missing) - set(["gdepw"])))

# Calculate other gdep values
gdep = horiz_interp_lev(grid["gdept"], grid["gdepw"], zgr_type, hgr_type)
Expand All @@ -264,7 +278,6 @@ def fill_zgrid_vars(zgr_type, grid, hgr_type, e_dict, missing):
grid[vi] = gdep[vi[4:]]

# e3

e3t_done = "e3t" not in missing
if e3t_done is False:
grid["e3t"] = vert_calc_e3(grid["gdept"], grid["gdepw"], "e3t")
Expand Down Expand Up @@ -301,8 +314,10 @@ def calc_gdepw(gdept):
dep_out = np.zeros((gdept.shape))
indxt = np.arange(gdept.shape[1])
indxw = indxt[:-1] + 0.5

for j in range(gdept.shape[2]):
for i in range(gdept.shape[3]):
# dep_out[0, 1:, j, i] = np.interp(indxw, indxt, gdept[0, :, j, i])
func = interp.interp1d(indxt, gdept[0, :, j, i], kind="cubic")
# top gdepw is zero
dep_out[0, 1:, j, i] = func(indxw)
Expand Down
4 changes: 2 additions & 2 deletions src/pybdy/nemo_bdy_chunk.py
Original file line number Diff line number Diff line change
Expand Up @@ -361,8 +361,8 @@ def chunk_large(ibdy, jbdy, chunk_number):
for i in range(len(all_chunk)):
i_chk = ibdy[chunk_number == all_chunk[i]]
j_chk = jbdy[chunk_number == all_chunk[i]]
i_range = np.max(i_chk) - np.min(i_chk)
j_range = np.max(j_chk) - np.min(j_chk)
i_range = (np.max(i_chk) - np.min(i_chk)) + 1
j_range = (np.max(j_chk) - np.min(j_chk)) + 1
ratio_box[i] = np.sum(chunk_number == all_chunk[i]) / (i_range * j_range)

# Calculate which boundaries are too big
Expand Down
8 changes: 4 additions & 4 deletions src/pybdy/nemo_bdy_setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -86,9 +86,9 @@ def _get_var_name_value(self, line):
index = "-1"
value = ""
if name_prefix == "ln":
if value_str.find("true") is not -1:
if value_str.find("true") != -1:
value = True
elif value_str.find("false") is not -1:
elif value_str.find("false") != -1:
value = False
else:
raise ValueError(
Expand Down Expand Up @@ -242,11 +242,11 @@ def _get_val(vars_dictionary, bool_vars_dictionary, line):
value = line[1].strip()

if name_prefix == "ln":
if value.find("true") is not -1:
if value.find("true") != -1:
if (name in vars_dictionary) is not True:
vars_dictionary[name] = True
bool_vars_dictionary[name] = True
elif value.find("false") is not -1:
elif value.find("false") != -1:
if (name in vars_dictionary) is not True:
vars_dictionary[name] = False
bool_vars_dictionary[name] = False
Expand Down
60 changes: 29 additions & 31 deletions src/pybdy/profiler.py
Original file line number Diff line number Diff line change
Expand Up @@ -90,7 +90,6 @@ def process_bdy(setup_filepath=0, mask_gui=False):

"""
# Start Logger

logger.info("Start NRCT Logging: " + time.asctime())
logger.info("============================================")

Expand All @@ -100,10 +99,11 @@ def process_bdy(setup_filepath=0, mask_gui=False):
Setup = setup.Setup(setup_filepath) # default settings file
settings = Setup.settings

logger.info("Reading grid completed")
logger.info("Reading setup completed")

bdy_msk = _get_mask(Setup, mask_gui)
DstCoord.bdy_msk = bdy_msk == 1
logger.info("Reading mask completed")

DstCoord.hgr = hgr.H_Grid(settings["dst_hgr"], settings["nme_map"], logger, dst=1)
DstCoord.zgr = zgr.Z_Grid(
Expand All @@ -114,8 +114,7 @@ def process_bdy(setup_filepath=0, mask_gui=False):
logger,
dst=1,
)

logger.info("Reading mask completed")
logger.info("Reading dst grid completed")

bdy_ind = {} # define a dictionary to hold the grid information

Expand All @@ -136,38 +135,11 @@ def process_bdy(setup_filepath=0, mask_gui=False):
logger.info("File: coordinates.bdy.nc generated and populated")

# Idenitify number of boundary points

nbdy = {}

for grd in ["t", "u", "v"]:
nbdy[grd] = len(bdy_ind[grd].bdy_i[:, 0])

# Gather grid information

logger.info("Gathering grid information")
SourceCoord.hgr = hgr.H_Grid(
settings["src_hgr"], settings["nme_map"], logger, dst=0
)
SourceCoord.zgr = zgr.Z_Grid(
settings["src_zgr"],
settings["nme_map"],
SourceCoord.hgr.grid_type,
SourceCoord.hgr.grid,
logger,
dst=0,
)

# Fill horizontal grid information

try: # if they are masked array convert them to normal arrays
SourceCoord.hgr.grid["glamt"] = SourceCoord.hgr.grid["glamt"].filled() # lon
except Exception:
logger.debug("Not a masked array.")
try:
SourceCoord.hgr.grid["gphit"] = SourceCoord.hgr.grid["gphit"].filled() # lat
except Exception:
logger.debug("Not a masked array.")

# Assign horizontal grid data

DstCoord.lonlat = {"t": {}, "u": {}, "v": {}}
Expand Down Expand Up @@ -196,6 +168,32 @@ def process_bdy(setup_filepath=0, mask_gui=False):

logger.info("BDY lons/lats identified from %s", settings["dst_hgr"])

# Gather grid information

logger.info("Gathering src grid information")
SourceCoord.hgr = hgr.H_Grid(
settings["src_hgr"], settings["nme_map"], logger, dst=0
)
SourceCoord.zgr = zgr.Z_Grid(
settings["src_zgr"],
settings["nme_map"],
SourceCoord.hgr.grid_type,
SourceCoord.hgr.grid,
logger,
dst=0,
)

# Fill horizontal grid information

try: # if they are masked array convert them to normal arrays
SourceCoord.hgr.grid["glamt"] = SourceCoord.hgr.grid["glamt"].filled() # lon
except Exception:
logger.debug("Not a masked array.")
try:
SourceCoord.hgr.grid["gphit"] = SourceCoord.hgr.grid["gphit"].filled() # lat
except Exception:
logger.debug("Not a masked array.")

# Define z at t/u/v points

DstCoord.depths = {"t": {}, "u": {}, "v": {}}
Expand Down