Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feature 1954 semilatlon #2262

Merged
merged 40 commits into from Sep 18, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
40 commits
Select commit Hold shift + click to select a range
a445afb
Per #1954, adding initial unstructured grid library code.
JohnHalleyGotway Jul 7, 2022
cd6e32f
Per #1954, since I'm adding new files to the vx_grid library, I had t…
JohnHalleyGotway Jul 7, 2022
775ffed
Per #1954, adding operator== and serialize() functions to the NumArra…
JohnHalleyGotway Jul 7, 2022
88afd03
Merge branch 'feature_1954_ugrid' of https://github.com/dtcenter/MET …
JohnHalleyGotway Jul 7, 2022
d32c648
Per #1954, switch from storing times in TimeArray to NumArray to simp…
JohnHalleyGotway Jul 7, 2022
8966e39
Per #1954, when the compilation fails add the full make_install.log f…
JohnHalleyGotway Jul 7, 2022
70361a9
Per #1954, update python embedding to parse the unstructured grid type.
JohnHalleyGotway Jul 7, 2022
b63e821
Per #1954, instantiate unstructured grid from data.
JohnHalleyGotway Jul 7, 2022
e447a19
Per #1954, fix up unstructured grid code.
JohnHalleyGotway Jul 8, 2022
d3ec1b6
Merge remote-tracking branch 'origin/develop' into feature_1954_ugrid
JohnHalleyGotway Jul 26, 2022
424a6f1
Merge remote-tracking branch 'origin/develop' into feature_1954_ugrid
JohnHalleyGotway Aug 11, 2022
083487d
Per #1954, rename Dim1/Dim2 to xDim/yDim for greater clarity.
JohnHalleyGotway Aug 11, 2022
4a1a2d4
Per #1954, UnstructuredData structs require a deep copy rather than j…
JohnHalleyGotway Aug 11, 2022
1914219
Per #1954, fix the logic to skip the plotting of map data for unstruc…
JohnHalleyGotway Aug 11, 2022
eb17db0
Per #1954, move the logic for skipping the plotting of map data from …
JohnHalleyGotway Aug 11, 2022
d8c3b84
Per #1954, link -lvx_cal into the compilation of test_grid_area utili…
JohnHalleyGotway Aug 11, 2022
36611f6
Merge branch 'develop' into feature_1954_semilatlon
JohnHalleyGotway Aug 17, 2022
76dafdb
Per #1954, rename from unstructured grid to semilatlon grid to avoid …
JohnHalleyGotway Aug 17, 2022
95d9bdf
Per #1954, add NumArray::summarize() function and update the SemiLatL…
JohnHalleyGotway Aug 19, 2022
ed85c8a
Per #1954, switch flag from Is2Dim to IsLatLon to avoid confusion.
JohnHalleyGotway Aug 19, 2022
7d3efaf
Merge remote-tracking branch 'origin/develop' into feature_1954_semil…
JohnHalleyGotway Aug 23, 2022
4b22445
Merge remote-tracking branch 'origin/develop' into feature_1954_semil…
JohnHalleyGotway Aug 30, 2022
803d01b
Per #1954, define new semilatlon projection type string.
JohnHalleyGotway Aug 31, 2022
b7db013
Per #1954, define logic to read semilatlon data via python embedding …
JohnHalleyGotway Aug 31, 2022
f544c1e
Per #1954, update the library code to have the write_netcdf_proj() fu…
JohnHalleyGotway Aug 31, 2022
aec6c45
Per #1954, update all the application code to remove the creation of …
JohnHalleyGotway Aug 31, 2022
d2019aa
Per #1954, fix login in met_file.cc to define x_dim_name and y_dim_na…
JohnHalleyGotway Aug 31, 2022
d4605af
Merge remote-tracking branch 'origin/develop' into feature_1954_semil…
JohnHalleyGotway Sep 2, 2022
82ec2f7
Per #1954, trying to debug difference in behavior on my Mac vs Seneca.
JohnHalleyGotway Sep 6, 2022
eac02b7
Per #1954, unrelated to development, but correcting an apparent typo …
JohnHalleyGotway Sep 6, 2022
9528d82
Per #1954, removing development print statement.
JohnHalleyGotway Sep 6, 2022
95fea19
Per #1954, update the grid_from_python_dict.cc logic to parse the lat…
JohnHalleyGotway Sep 6, 2022
092bd34
Merge remote-tracking branch 'origin/develop' into feature_1954_semil…
JohnHalleyGotway Sep 15, 2022
a3bc30c
Per #1954, ci-run-unit add unit tests to demonstrate the newly added …
JohnHalleyGotway Sep 15, 2022
c60dd6b
ci-run-unit Fix expected file type
JohnHalleyGotway Sep 16, 2022
c164c34
Per #1954, fix quote in log message.
JohnHalleyGotway Sep 16, 2022
fe6b03d
Per #1954, ci-run-unit fix failure from unit_point2grid.xml. Had to a…
JohnHalleyGotway Sep 16, 2022
2f92e78
Per #1954, GHA doesn't allow for colons in output file names. Switchi…
JohnHalleyGotway Sep 16, 2022
0631729
Per #1954, add documentation about the semilatlon grid.
JohnHalleyGotway Sep 16, 2022
9c9858b
Small grammar update
j-opatz Sep 16, 2022
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
3 changes: 3 additions & 0 deletions .github/jobs/build_docker_image.sh
Expand Up @@ -15,5 +15,8 @@ time_command docker build -t ${DOCKERHUB_TAG} \
-f $DOCKERFILE_PATH ${GITHUB_WORKSPACE}
if [ $? != 0 ]; then
cat ${GITHUB_WORKSPACE}/docker_build.log
# Append the full make_install.log file
echo "Appending make_install.log to docker_build.log. See the logs artifact for details."
cat ${GITHUB_WORKSPACE}/make_install.log >> ${GITHUB_WORKSPACE}/docker_build.log
exit 1
fi
6 changes: 6 additions & 0 deletions docs/Users_Guide/appendixB.rst
Expand Up @@ -23,6 +23,8 @@ The following map projections are currently supported in MET:

* Gaussian Projection

* Semi Lat/Lon

Grid Specification Strings
==========================

Expand Down Expand Up @@ -82,6 +84,10 @@ For a Gaussian grid, the syntax is

The parameters **Nx** and **Ny** are as before, while **lon_zero** defines the first longitude.

For a Semi Lat/Lon grid, no grid specification string is supported. This grid type is only supported via Python embedding or when reading NetCDF files generated by another MET tool. A Semi Lat/Lon grid defines the information about 2D field of data whose dimension are defined by arrays of latitude (**lats**), longitude (**lons**), level (**levels**), and time (**times**). Times are defined as unixtime, the number of seconds since January 1, 1970. Typically, the lats or lons array and the levels or times array has non-zero length. For example, a zonal mean field is defined using the lats and levels array. A meridional mean field is defined using the lons and levels array. A Hovmoeller field is defined using lats or lons versus times. An arbitrary cross-section is defined by specifying both the lats and lons array with exactly the same length versus levels or times.

Statistics can be computed from data on Semi Lat/Lon grids but only when all data resides on the same Semi Lat/Lon grid. Two Semi Lat/Lon grids are equal when their lats, lons, levels, and times arrays match. No functionality is provided to regrid Semi Lat/Lon data. The MET tools can plot Semi Lat/Lon data, however no map data is overlaid since these grids lack two spatial dimensions.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think I know the answer, but just clarifying: when it's stated "Two Semi Lat/Lon grids are equal when their lats, lons, levels, and times arrays match.", that means the values of the arrays match, correct? Didn't want the confusion of "equal" meaning lengths of the arrays match.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's correct, the code checks that the arrays have the same length and all the array values match.


Grids
=====

Expand Down
9 changes: 9 additions & 0 deletions docs/Users_Guide/appendixF.rst
Expand Up @@ -173,6 +173,15 @@ When specified as a dictionary, the contents of the **grid** dictionary vary bas
• lon_zero (double)
• nx, ny (int)

• **SemiLatLon** grid dictionary entries:

• type ("SemiLatLon")
• name (string)
• lats (list of doubles)
• lons (list of doubles)
• levels (list of doubles)
• times (list of doubles)

Additional information about supported grids can be found in :ref:`appendixB`.

**Using Xarray DataArrays**
Expand Down
41 changes: 41 additions & 0 deletions internal/test_unit/xml/unit_python.xml
Expand Up @@ -587,4 +587,45 @@
</output>
</test>

<test name="python_plot_data_plane_SEMILATLON_ZONAL_MEAN">
<exec>&MET_BIN;/plot_data_plane</exec>
<param> \
PYTHON_NUMPY \
&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_ZONAL_MEAN.ps \
'name="&MET_BASE;/python/derive_WRF_semilatlon.py &DATA_DIR_MODEL;/p_interp/wrfout_d01_2008-08-08_12:00:00_PLEV TT lat";' \
-title "WRF Zonal Mean" \
-v 1
</param>
<output>
<ps>&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_ZONAL_MEAN.ps</ps>
</output>
</test>

<test name="python_pcp_combine_SEMILATLON_MERIDIONAL_MEAN">
<exec>&MET_BIN;/pcp_combine</exec>
<param> \
-add PYTHON_NUMPY \
'name="&MET_BASE;/python/derive_WRF_semilatlon.py &DATA_DIR_MODEL;/p_interp/wrfout_d01_2008-08-08_12:00:00_PLEV TT lon";' \
&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_MERIDIONAL_MEAN.nc \
-name "TT_MERIDIONAL_MEAN" -v 1
</param>
<output>
<grid_nc>&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_MERIDIONAL_MEAN.nc</grid_nc>
</output>
</test>

<test name="python_plot_data_plane_SEMILATLON_MERIDIONAL_MEAN">
<exec>&MET_BIN;/plot_data_plane</exec>
<param> \
&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_MERIDIONAL_MEAN.nc \
&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_MERIDIONAL_MEAN.ps \
'name="TT_MERIDIONAL_MEAN"; level="(*,*)";' \
-title "WRF Meridional Mean" \
-v 1
</param>
<output>
<ps>&OUTPUT_DIR;/python/wrfout_d01_2008-08-08_12_00_00_PLEV_MERIDIONAL_MEAN.ps</ps>
</output>
</test>

</met_test>
1 change: 1 addition & 0 deletions internal/test_util/libcode/vx_grid/Makefile.am
Expand Up @@ -20,6 +20,7 @@ test_grid_area_LDADD = -lvx_grid \
-lvx_util_math \
-lvx_util \
-lvx_math \
-lvx_cal \
-lvx_log \
-lz \
-lm
1 change: 1 addition & 0 deletions internal/test_util/libcode/vx_grid/Makefile.in
Expand Up @@ -320,6 +320,7 @@ test_grid_area_LDADD = -lvx_grid \
-lvx_util_math \
-lvx_util \
-lvx_math \
-lvx_cal \
-lvx_log \
-lz \
-lm
Expand Down
3 changes: 2 additions & 1 deletion scripts/python/Makefile.am
Expand Up @@ -32,7 +32,8 @@ pythonscripts_DATA = \
read_ascii_xarray.py \
read_ascii_point.py \
read_ascii_mpr.py \
read_met_point_obs.py
read_met_point_obs.py \
derive_WRF_semilatlon.py

EXTRA_DIST = ${pythonscripts_DATA}

Expand Down
3 changes: 2 additions & 1 deletion scripts/python/Makefile.in
Expand Up @@ -303,7 +303,8 @@ pythonscripts_DATA = \
read_ascii_xarray.py \
read_ascii_point.py \
read_ascii_mpr.py \
read_met_point_obs.py
read_met_point_obs.py \
derive_WRF_semilatlon.py

EXTRA_DIST = ${pythonscripts_DATA}
MAINTAINERCLEANFILES = Makefile.in
Expand Down
103 changes: 103 additions & 0 deletions scripts/python/derive_WRF_semilatlon.py
@@ -0,0 +1,103 @@
from __future__ import print_function

import os
import sys
import numpy as np
import datetime as dt
from netCDF4 import Dataset,chartostring

###########################################

##
## input file specified on the command line
## load the data into the numpy array
##

if len(sys.argv) == 4:
# Read the input file as the first argument
input_file = os.path.expandvars(sys.argv[1])
var_name = sys.argv[2]
axis = sys.argv[3]

try:
# Print some output to verify that this script ran
print("Input File: " + repr(input_file))
print("Variable: " + repr(var_name))
print("Axis: " + repr(axis))

# Read input file
f = Dataset(input_file, 'r')

data = np.float64(f.variables[var_name][0,:,:,:])
data[data > 1.0e30] = np.nan

pvals = list(np.float64(f.variables["pressure"][:]))

if axis == "lon":
met_data = np.nanmean(data[::-1], axis=1).copy()
elif axis == "lat":
met_data = np.nanmean(data[::-1], axis=2).transpose().copy()
elif axis == "latlon":
met_data = np.nanmean(data[::-1], axis=1).copy()
else:
print("ERROR: Unsupported axis type: " + axis)
sys.exit(1)

print("Data Shape: " + repr(met_data.shape))
print("Data Type: " + repr(met_data.dtype))
except NameError:
print("Trouble reading data from input file")
else:
print("Must specify exactly one input file, variable name, and summary axis (lat, lon, latlon).")
sys.exit(1)

###########################################

##
## create the metadata dictionary
##

init = dt.datetime.strptime(getattr(f, "START_DATE"), "%Y-%m-%d_%H:%M:%S")
valid_ref = dt.datetime.strptime(getattr(f.variables["Time"], "units"), "hours since %Y-%m-%d %H:%M:%S")
add_hours = float(f.variables["Time"][:])
valid = valid_ref + dt.timedelta(hours=add_hours)
lead, rem = divmod((valid-init).total_seconds(), 3600)
accum = "00"

# Use the first column of lats

if axis == "lon":
lats = list()
lons = list(np.float64(f.variables["XLONG"][0,0,:]))
elif axis == "lat":
lats = list(np.float64(f.variables["XLAT"][0,:,0]))
lons = list()
elif axis == "latlon":
lats = list(np.float64(f.variables["XLONG"][0,0,:]))
lons = list(np.float64(f.variables["XLAT"][0,0,:]))

levels = list(pvals)
times = list()

attrs = {
'valid': valid.strftime("%Y%m%d_%H%M%S"),
'init': init.strftime("%Y%m%d_%H%M%S"),
'lead': str(int(lead)),
'accum': accum,

'name': var_name,
'long_name': str(getattr(f.variables[var_name], "description")),
'level': axis + "_mean",
'units': str(getattr(f.variables[var_name], "units")),

'grid': {
'type' : "SemiLatLon",
'name' : axis + "_mean",
'lats' : lats,
'lons' : lons,
'levels' : levels,
'times' : times
}
}

print("Attributes: " + repr(attrs))
45 changes: 45 additions & 0 deletions src/basic/vx_cal/time_array.cc
Expand Up @@ -94,6 +94,30 @@ return ( * this );
////////////////////////////////////////////////////////////////////////


bool TimeArray::operator==(const TimeArray & a) const

{

if ( Nelements != a.Nelements ) return ( false );

bool status = true;

for (int j=0; j<(Nelements); ++j) {

if ( e[j] != a.e[j] ) {
status = false;
break;
}
}

return ( status );

}


////////////////////////////////////////////////////////////////////////


void TimeArray::init_from_scratch()

{
Expand Down Expand Up @@ -460,6 +484,27 @@ return(u);
////////////////////////////////////////////////////////////////////////


ConcatString TimeArray::serialize() const

{

ConcatString s;

if(n_elements() == 0) return(s);

int j;

s << e[0];
for(j=1; j<n_elements(); j++) s << " " << unix_to_yyyymmdd_hhmmss(e[j]);

return(s);

}


////////////////////////////////////////////////////////////////////////


void TimeArray::sort_array()

{
Expand Down
5 changes: 4 additions & 1 deletion src/basic/vx_cal/time_array.h
Expand Up @@ -52,6 +52,7 @@ class TimeArray {
~TimeArray();
TimeArray(const TimeArray &);
TimeArray & operator=(const TimeArray &);
bool operator==(const TimeArray &) const;

void clear();

Expand Down Expand Up @@ -80,6 +81,8 @@ class TimeArray {
unixtime min() const;
unixtime max() const;

ConcatString serialize() const;

int n_elements() const;
int n() const;

Expand All @@ -96,7 +99,7 @@ inline int TimeArray::n() const { return ( Nelements ); }
////////////////////////////////////////////////////////////////////////


extern ConcatString write_css (const TimeArray &);
extern ConcatString write_css(const TimeArray &);


////////////////////////////////////////////////////////////////////////
Expand Down
55 changes: 55 additions & 0 deletions src/basic/vx_util/num_array.cc
Expand Up @@ -93,6 +93,31 @@ NumArray & NumArray::operator=(const NumArray & a)
////////////////////////////////////////////////////////////////////////


bool NumArray::operator==(const NumArray & a) const

{

if ( e.size() != a.e.size() ) return ( false );

bool status = true;
int n = e.size();

for (int j=0; j<n; ++j) {

if ( e[j] != a.e[j] ) {
status = false;
break;
}
}

return ( status );

}


////////////////////////////////////////////////////////////////////////


void NumArray::init_from_scratch()

{
Expand Down Expand Up @@ -984,6 +1009,36 @@ ConcatString NumArray::serialize() const
////////////////////////////////////////////////////////////////////////


ConcatString NumArray::summarize() const

{

ConcatString s;

s << "n = " << n_elements();

if(n_elements() > 0) {

double min_v, max_v;
min_v = max_v = e[0];

for(int j=0; j<n_elements(); j++) {
if(is_bad_data(e[j])) continue;
if(e[j] < min_v) min_v = e[j];
if(e[j] > max_v) max_v = e[j];
}

s << ", min = " << min_v << ", max = " << max_v;
}

return(s);

}


////////////////////////////////////////////////////////////////////////


NumArray NumArray::subset(int beg, int end) const

{
Expand Down