diff --git a/data/wrappers/Makefile.am b/data/wrappers/Makefile.am index d2f83ff125..deb919438e 100644 --- a/data/wrappers/Makefile.am +++ b/data/wrappers/Makefile.am @@ -24,8 +24,10 @@ wrappers_DATA = \ set_python_env.py \ read_tmp_dataplane.py \ read_tmp_ascii.py \ + read_tmp_point_nc.py \ write_tmp_dataplane.py \ write_tmp_point.py \ + write_tmp_point_nc.py \ write_tmp_mpr.py EXTRA_DIST = ${wrappers_DATA} diff --git a/data/wrappers/Makefile.in b/data/wrappers/Makefile.in index 9ab052cca5..da04b2b2a0 100644 --- a/data/wrappers/Makefile.in +++ b/data/wrappers/Makefile.in @@ -279,6 +279,7 @@ MET_HDFLIB = @MET_HDFLIB@ MET_NETCDF = @MET_NETCDF@ MET_NETCDFINC = @MET_NETCDFINC@ MET_NETCDFLIB = @MET_NETCDFLIB@ +MET_PYTHON_BIN_EXE = @MET_PYTHON_BIN_EXE@ MET_PYTHON_CC = @MET_PYTHON_CC@ MET_PYTHON_LD = @MET_PYTHON_LD@ MKDIR_P = @MKDIR_P@ @@ -359,8 +360,10 @@ wrappers_DATA = \ set_python_env.py \ read_tmp_dataplane.py \ read_tmp_ascii.py \ + read_tmp_point_nc.py \ write_tmp_dataplane.py \ write_tmp_point.py \ + write_tmp_point_nc.py \ write_tmp_mpr.py EXTRA_DIST = ${wrappers_DATA} diff --git a/data/wrappers/read_tmp_point_nc.py b/data/wrappers/read_tmp_point_nc.py new file mode 100644 index 0000000000..8afaf16684 --- /dev/null +++ b/data/wrappers/read_tmp_point_nc.py @@ -0,0 +1,30 @@ +######################################################################## +# +# Reads temporary file into memory. +# +# usage: /path/to/python read_tmp_dataplane.py dataplane.tmp +# +######################################################################## + +import os +import sys + +met_base_dir = os.getenv('MET_BASE',None) +if met_base_dir is not None: + sys.path.append(os.path.join(met_base_dir, 'python')) + +# add share/met/python directory to system path to find met_point_obs +sys.path.append(os.path.abspath(os.path.join(os.path.dirname(__file__), + os.pardir, 'python'))) +from met_point_obs import met_point_obs +from met_point_obs_nc import nc_point_obs + +netcdf_filename = sys.argv[1] + +# read NetCDF file +print('{p} reading{f}'.format(p=met_point_obs.get_prompt(), f=netcdf_filename)) +point_obs_data = nc_point_obs() +point_obs_data.read_data(netcdf_filename) + +met_point_data = point_obs_data.get_point_data() +met_point_data['met_point_data'] = point_obs_data diff --git a/data/wrappers/write_tmp_point_nc.py b/data/wrappers/write_tmp_point_nc.py new file mode 100644 index 0000000000..1dbda19807 --- /dev/null +++ b/data/wrappers/write_tmp_point_nc.py @@ -0,0 +1,59 @@ +######################################################################## +# +# Adapted from a script provided by George McCabe +# Adapted by Randy Bullock +# +# usage: /path/to/python write_tmp_point.py \ +# tmp_output_filename .py +# +######################################################################## + +import os +import sys +import importlib.util + +met_base_dir = os.getenv('MET_BASE',None) +if met_base_dir is not None: + sys.path.append(os.path.join(met_base_dir, 'python')) + +# add share/met/python directory to system path to find met_point_obs +sys.path.append(os.path.abspath(os.path.join(os.path.dirname(__file__), + os.pardir, 'python'))) + +from met_point_obs import met_point_obs +from met_point_obs_nc import nc_point_obs + +PROMPT = met_point_obs.get_prompt() +print("{p} Python Script:\t".format(p=PROMPT) + repr(sys.argv[0])) +print("{p} User Command:\t".format(p=PROMPT) + repr(' '.join(sys.argv[2:]))) +print("{p} Temporary File:\t".format(p=PROMPT) + repr(sys.argv[1])) + +tmp_filename = sys.argv[1] +pyembed_module_name = sys.argv[2] +sys.argv = sys.argv[2:] + +# append user script dir to system path +pyembed_dir, pyembed_file = os.path.split(pyembed_module_name) +if pyembed_dir: + sys.path.insert(0, pyembed_dir) + +if not pyembed_module_name.endswith('.py'): + pyembed_module_name += '.py' + +user_base = os.path.basename(pyembed_module_name).replace('.py','') + +spec = importlib.util.spec_from_file_location(user_base, pyembed_module_name) +met_in = importlib.util.module_from_spec(spec) +spec.loader.exec_module(met_in) + +if hasattr(met_in, 'point_obs_data'): + met_in.point_obs_data.save_ncfile(tmp_filename) +else: + if hasattr(met_in.met_point_data, 'point_obs_data'): + met_in.met_point_data['point_obs_data'].save_ncfile(tmp_filename) + else: + tmp_point_obs = nc_point_obs() + tmp_point_obs.put_data(met_in.met_point_data) + tmp_point_obs.save_ncfile(tmp_filename) + +#print('{p} writing {f}'.format(p=PROMPT, f=tmp_filename)) diff --git a/internal/test_unit/xml/unit_python.xml b/internal/test_unit/xml/unit_python.xml index 05f9d2650d..5a519d9212 100644 --- a/internal/test_unit/xml/unit_python.xml +++ b/internal/test_unit/xml/unit_python.xml @@ -495,6 +495,25 @@ + + + &MET_BIN;/point2grid + + MET_PYTHON_EXE &MET_PYTHON_EXE; + + \ + 'PYTHON_NUMPY=&MET_BASE;/python/read_met_point_obs.py &OUTPUT_DIR;/pb2nc/ndas.20120409.t12z.prepbufr.tm00.nc' \ + G212 \ + &OUTPUT_DIR;/python/pb2nc_TMP_user_python.nc \ + -field 'name="TMP"; level="*"; valid_time="20120409_120000"; censor_thresh=[ <0 ]; censor_val=[0];' \ + -name TEMP \ + -v 1 + + + &OUTPUT_DIR;/python/pb2nc_TMP_user_python.nc + + + &MET_BIN;/plot_point_obs diff --git a/scripts/python/Makefile.am b/scripts/python/Makefile.am index 76aa04cdf1..689708e4c3 100644 --- a/scripts/python/Makefile.am +++ b/scripts/python/Makefile.am @@ -27,6 +27,7 @@ pythonscriptsdir = $(pkgdatadir)/python pythonscripts_DATA = \ met_point_obs.py \ + met_point_obs_nc.py \ read_ascii_numpy.py \ read_ascii_numpy_grid.py \ read_ascii_xarray.py \ diff --git a/scripts/python/Makefile.in b/scripts/python/Makefile.in index 80da2e8f84..6d85ed81f9 100644 --- a/scripts/python/Makefile.in +++ b/scripts/python/Makefile.in @@ -221,6 +221,7 @@ MET_HDFLIB = @MET_HDFLIB@ MET_NETCDF = @MET_NETCDF@ MET_NETCDFINC = @MET_NETCDFINC@ MET_NETCDFLIB = @MET_NETCDFLIB@ +MET_PYTHON_BIN_EXE = @MET_PYTHON_BIN_EXE@ MET_PYTHON_CC = @MET_PYTHON_CC@ MET_PYTHON_LD = @MET_PYTHON_LD@ MKDIR_P = @MKDIR_P@ @@ -298,6 +299,7 @@ top_srcdir = @top_srcdir@ pythonscriptsdir = $(pkgdatadir)/python pythonscripts_DATA = \ met_point_obs.py \ + met_point_obs_nc.py \ read_ascii_numpy.py \ read_ascii_numpy_grid.py \ read_ascii_xarray.py \ diff --git a/scripts/python/met_point_obs.py b/scripts/python/met_point_obs.py index 5cfc5f7141..41bc6462f6 100755 --- a/scripts/python/met_point_obs.py +++ b/scripts/python/met_point_obs.py @@ -1,21 +1,25 @@ - +#!/usr/bin/env python3 ''' Created on Nov 10, 2021 @author: hsoh - This is the base class and the customized script should extend the met_point_obs. -- The customized script (for example "custom_reader") must implement "def read_data(self, args)" - which fills the array (python list or numpy array) variables at __init__(). -- The args can be 1) single string argument, 2) the list of arguments, and 3) the dictionary of arguments. -- The variable "met_point_data" must be set for MET tools +- The customized script (for example "custom_reader") must implement + "def read_data(self, args)" which fills the array variables at __init__(). +- The args can be 1) single string argument, 2) the list of arguments, + or 3) the dictionary of arguments. +- A python objects, met_point_data, must set: + + "point_obs_data" is an optional to use custom python EXE. + It's a python instance which processes the point observation data - The customized script is expected to include following codes: # prepare arguments for the customized script args = {'input', sys.argv[1]} # or args = [] point_obs_data = custom_reader() - point_obs_data.read_data() + point_obs_data.read_data(args) met_point_data = point_obs_data.get_point_data() + ''' import os @@ -24,22 +28,29 @@ COUNT_SHOW = 30 -MET_PYTHON_OBS_ARGS = "MET_POINT_PYTHON_ARGS" +def get_prompt(): + return " python:" + +def met_is_python_prefix(user_cmd): + return user_cmd.startswith(base_met_point_obs.python_prefix) class base_met_point_obs(object): ''' classdocs ''' - ERROR_P = " ==ERROR_PYTHON==" - INFO_P = " ==INFO_PYTHON==" + ERROR_P = " ==PYTHON_ERROR==" + INFO_P = " ==PYTHON_INFO==" + + python_prefix = 'PYTHON_POINT_USER' - python_prefix = 'PYTHON_POINT_RAW' + FILL_VALUE = -9999. def __init__(self, use_var_id=True): ''' Constructor ''' + self.count_info = "" self.input_name = None self.ignore_input_file = False self.use_var_id = use_var_id # True if variable index, False if GRIB code @@ -74,6 +85,8 @@ def __init__(self, use_var_id=True): self.obs_val = [] # (nobs) float self.obs_qty_table = [] # (nobs_qty, mxstr) string self.obs_var_table = [] # (nobs_var, mxstr2) string, required if self.use_var_id is True + self.obs_var_unit = [] # (nobs_var, mxstr2) string, optional if self.use_var_id is True + self.obs_var_desc = [] # (nobs_var, mxstr3) string, optional if self.use_var_id is True # Optional variables for PREPBUFR, not supported yet self.hdr_prpt_typ = [] # optional @@ -98,17 +111,12 @@ def check_data_member_float(self, local_var, var_name): self.add_error_msg("{v} is empty (float)".format(v=var_name)) elif isinstance(local_var, list): if isinstance(local_var[0], str) and not self.is_number(local_var[0]): - self.add_error_msg("Not supported data type ({v} string) for {n}[0] (int or float only)".format( + self.add_error_msg("Not supported data type: {n}[0]={v}, string type, not a number (int or float only)".format( n=var_name, v=local_var[0])) - elif 0 <= str(type(local_var[0])).find('numpy'): - self.log_info("Recommend using numpy instead of python list for {v} ({t}) to avoid side effect".format( - v=var_name, t=type(local_var[0]))) - elif not isinstance(local_var[0], (int, float)): - # self.add_error_msg("Not supported data type ({t}) for {v}[0] (int or float only)".format( - # v=var_name, t=local_var[0].type)) + elif 0 > str(type(local_var[0])).find('numpy') and not isinstance(local_var[0], (int, float)): self.add_error_msg("Not supported data type ({t}) for {v}[0] (int or float only)".format( v=var_name, t=type(local_var[0]))) - elif not isinstance(local_var, np.ndarray): + elif not self.is_numpy_array(local_var): self.add_error_msg("Not supported data type ({t}) for {v} (list and numpy.ndarray)".format( v=var_name, t=type(local_var))) @@ -117,15 +125,12 @@ def check_data_member_int(self, local_var, var_name): self.add_error_msg("{v} is empty (int)".format(v=var_name)) elif isinstance(local_var, list): if isinstance(local_var[0], str) and not self.is_number(local_var[0]): - self.add_error_msg("Not supported data type ({v} string) for {n}[0] (int or float only)".format( + self.add_error_msg("Not supported data type: {n}[0]={v}, string type, not a number (int only)".format( n=var_name, v=local_var[0])) - elif 0 <= str(type(local_var[0])).find('numpy'): - self.log_info("Recommend using numpy instead of python list for {v} ({t}) to avoid side effect".format( - v=var_name, t=type(local_var[0]))) - elif not isinstance(local_var[0], int): + elif 0 > str(type(local_var[0])).find('numpy') and not isinstance(local_var[0], int): self.add_error_msg("Not supported data type ({t}) for {v}[0] (int only)".format( v=var_name, t=type(local_var[0]))) - elif not isinstance(local_var, np.ndarray): + elif not self.is_numpy_array(local_var): self.add_error_msg("Not supported data type ({t}) for {v} (list and numpy.ndarray)".format( v=var_name, t=type(local_var))) @@ -160,6 +165,15 @@ def check_point_data(self): if self.use_var_id: self.check_data_member_string(self.obs_var_table,'obs_var_table') + def convert_to_numpy(self, value_list): + return np.array(value_list) + + def dump(self): + base_met_point_obs.print_point_data(self.get_point_data()) + + def get_count_string(self): + return f' nobs={self.nobs} nhdr={self.nhdr} ntyp={self.nhdr_typ} nsid={self.nhdr_sid} nvld={self.nhdr_vld} nqty={self.nobs_qty} nvar={self.nobs_var}' + def get_point_data(self): if self.nhdr <= 0: self.nhdr = len(self.hdr_lat) @@ -178,34 +192,108 @@ def get_point_data(self): if self.nobs_var <= 0: self.nobs_var = len(self.obs_var_table) self.check_point_data() - return self.__dict__ - def get_type(self, value): - return 'string' if isinstance('str') else type(value) + if not self.is_numpy_array(self.hdr_typ): + self.hdr_typ = self.convert_to_numpy(self.hdr_typ) + if not self.is_numpy_array(self.hdr_sid): + self.hdr_sid = self.convert_to_numpy(self.hdr_sid) + if not self.is_numpy_array(self.hdr_vld): + self.hdr_vld = self.convert_to_numpy(self.hdr_vld) + if not self.is_numpy_array(self.hdr_lat): + self.hdr_lat = self.convert_to_numpy(self.hdr_lat) + if not self.is_numpy_array(self.hdr_lon): + self.hdr_lon = self.convert_to_numpy(self.hdr_lon) + if not self.is_numpy_array(self.hdr_elv): + self.hdr_elv = self.convert_to_numpy(self.hdr_elv) + + if not self.is_numpy_array(self.obs_qty): + self.obs_qty = self.convert_to_numpy(self.obs_qty) + if not self.is_numpy_array(self.obs_hid): + self.obs_hid = self.convert_to_numpy(self.obs_hid) + if not self.is_numpy_array(self.obs_vid): + self.obs_vid = self.convert_to_numpy(self.obs_vid) + if not self.is_numpy_array(self.obs_lvl): + self.obs_lvl = self.convert_to_numpy(self.obs_lvl) + if not self.is_numpy_array(self.obs_hgt): + self.obs_hgt = self.convert_to_numpy(self.obs_hgt) + if not self.is_numpy_array(self.obs_val): + self.obs_val = self.convert_to_numpy(self.obs_val) + + self.count_info = self.get_count_string() + self.met_point_data = self + return self.__dict__ def is_number(self, num_str): return num_str.replace('-','1').replace('+','2').replace('.','3').isdigit() + def is_numpy_array(self, var): + return isinstance(var, np.ndarray) + def log_error_msg(self, err_msg): - print('{p} {m}'.format(p=self.ERROR_P, m=err_msg)) + base_met_point_obs.error_msg(err_msg) def log_error(self, err_msgs): print(self.ERROR_P) for err_line in err_msgs.split('\n'): - print('{p} {m}'.format(p=self.ERROR_P, m=err_line)) + self.log_error_msg(err_line) print(self.ERROR_P) def log_info(self, info_msg): - print('{p} {m}'.format(p=self.INFO_P, m=info_msg)) + base_met_point_obs.info_msg(info_msg) + + def put_data(self, point_obs_dict): + self.hdr_typ = point_obs_dict['hdr_typ'] + self.hdr_sid = point_obs_dict['hdr_sid'] + self.hdr_vld = point_obs_dict['hdr_vld'] + self.hdr_lat = point_obs_dict['hdr_lat'] + self.hdr_lon = point_obs_dict['hdr_lon'] + self.hdr_elv = point_obs_dict['hdr_elv'] + self.hdr_typ_table = point_obs_dict['hdr_typ_table'] + self.hdr_sid_table = point_obs_dict['hdr_sid_table'] + self.hdr_vld_table = point_obs_dict['hdr_vld_table'] + + #Observation data + self.obs_qty = point_obs_dict['obs_qty'] + self.obs_hid = point_obs_dict['obs_hid'] + self.obs_lvl = point_obs_dict['obs_lvl'] + self.obs_hgt = point_obs_dict['obs_hgt'] + self.obs_val = point_obs_dict['obs_val'] + self.obs_vid = point_obs_dict['obs_vid'] + self.obs_var_table = point_obs_dict['obs_var_table'] + self.obs_qty_table = point_obs_dict['obs_qty_table'] + po_array = point_obs_dict.get('obs_unit', None) + if po_array is not None: + self.obs_var_unit = po_array + po_array = point_obs_dict.get('obs_desc', None) + if po_array is not None: + self.obs_var_desc = po_array + + po_array = point_obs_dict.get('hdr_prpt_typ', None) + if po_array is not None: + self.hdr_prpt_typ = po_array + po_array = point_obs_dict.get('hdr_irpt_typ', None) + if po_array is not None: + self.hdr_irpt_typ = po_array + po_array = point_obs_dict.get('hdr_inst_typ', None) + if po_array is not None: + self.hdr_inst_typ = po_array @staticmethod - def is_python_script(arg_value): - return arg_value.startswith(met_point_obs.python_prefix) + def error_msg(msg): + print(f'{get_prompt()} {base_met_point_obs.ERROR_P} {msg}') + + @staticmethod + def info_msg(msg): + print(f'{get_prompt()} {base_met_point_obs.INFO_P} {msg}') @staticmethod def get_python_script(arg_value): return arg_value[len(met_point_obs.python_prefix)+1:] + @staticmethod + def is_python_script(arg_value): + return arg_value.startswith(met_point_obs.python_prefix) + @staticmethod def print_data(key, data_array, show_count=COUNT_SHOW): if isinstance(data_array, list): @@ -274,22 +362,18 @@ def print_point_data(met_point_data, print_subset=True): class csv_point_obs(ABC, base_met_point_obs): def __init__(self, point_data): - super(csv_point_obs, self).__init__() - - hdr_cnt = obs_cnt = len(point_data) self.point_data = point_data - self.nhdr = self.nobs = obs_cnt - self.nhdr_typ = self.nhdr_sid = self.nhdr_vld = hdr_cnt - self.nobs_qty = self.nobs_var = obs_cnt - self.input_file = None - self.ignore_input_file = True + super(csv_point_obs, self).__init__() + self.obs_cnt = obs_cnt = len(point_data) self.obs_qty = [ 0 for _ in range(0, obs_cnt) ] # (nobs_qty) integer, index of self.obs_qty_table self.obs_hid = [ 0 for _ in range(0, obs_cnt) ] # (nobs) integer self.obs_vid = [ 0 for _ in range(0, obs_cnt) ] # (nobs) integer, veriable index from self.obs_var_table or GRIB code - self.obs_lvl = [ nc_tools.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float - self.obs_hgt = [ nc_tools.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float - self.obs_val = [ nc_tools.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float + self.obs_lvl = [ self.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float + self.obs_hgt = [ self.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float + self.obs_val = [ self.FILL_VALUE for _ in range(0, obs_cnt) ] # (nobs) float + + self.convert_point_data() def check_csv_record(self, csv_point_data, index): error_msgs = [] @@ -364,33 +448,33 @@ def convert_point_data(self): obs_qty_map = {} index = 0 - # Build map - # names=['typ', 'sid', 'vld', 'lat', 'lon', 'elv', 'var', 'lvl', 'hgt', 'qc', 'obs'] - for csv_point_data in self.point_data: - hdr_typ_str = csv_point_data[0] + #names=['typ', 'sid', 'vld', 'lat', 'lon', 'elv', 'var', 'lvl', 'hgt', 'qc', 'obs'] + for csv_point_record in self.point_data: + # Build header map. + hdr_typ_str = csv_point_record[0] hdr_typ_idx = hdr_typ_map.get(hdr_typ_str,-1) if hdr_typ_idx < 0: hdr_typ_idx = hdr_typ_cnt hdr_typ_map[hdr_typ_str] = hdr_typ_idx hdr_typ_cnt += 1 - hdr_sid_str = csv_point_data[1] + hdr_sid_str = csv_point_record[1] hdr_sid_idx = hdr_sid_map.get(hdr_sid_str,-1) if hdr_sid_idx < 0: hdr_sid_idx = hdr_sid_cnt hdr_sid_map[hdr_sid_str] = hdr_sid_idx hdr_sid_cnt += 1 - hdr_vld_str = csv_point_data[2] + hdr_vld_str = csv_point_record[2] hdr_vld_idx = hdr_vld_map.get(hdr_vld_str,-1) if hdr_vld_idx < 0: hdr_vld_idx = hdr_vld_cnt hdr_vld_map[hdr_vld_str] = hdr_vld_idx hdr_vld_cnt += 1 - lat = csv_point_data[3] - lon = csv_point_data[4] - elv = self.get_num_value(csv_point_data[5] ) + lat = csv_point_record[3] + lon = csv_point_record[4] + elv = self.get_num_value(csv_point_record[5] ) hdr_key = (hdr_typ_idx,hdr_sid_idx,hdr_vld_idx,lat,lon,elv) hdr_idx = hdr_map.get(hdr_key,-1) if hdr_idx < 0: @@ -398,7 +482,7 @@ def convert_point_data(self): hdr_map[hdr_key] = hdr_idx hdr_cnt += 1 - var_id_str = csv_point_data[6] + var_id_str = csv_point_record[6] if self.use_var_id: var_id = obs_var_map.get(var_id_str,-1) if var_id < 0: @@ -408,7 +492,7 @@ def convert_point_data(self): else: var_id = int(var_id_str) - qc_str = csv_point_data[9] + qc_str = csv_point_record[9] qc_id = obs_qty_map.get(qc_str,-1) if qc_id < 0: qc_id = qc_cnt @@ -418,9 +502,9 @@ def convert_point_data(self): # names=['typ', 'sid', 'vld', 'lat', 'lon', 'elv', 'var', 'lvl', 'hgt', 'qc', 'obs'] self.obs_vid[index] = var_id self.obs_hid[index] = hdr_idx - self.obs_lvl[index] = self.get_num_value(csv_point_data[7]) - self.obs_hgt[index] = self.get_num_value(csv_point_data[8]) - self.obs_val[index] = self.get_num_value(csv_point_data[10]) + self.obs_lvl[index] = self.get_num_value(csv_point_record[7]) + self.obs_hgt[index] = self.get_num_value(csv_point_record[8]) + self.obs_val[index] = self.get_num_value(csv_point_record[10]) self.obs_qty[index] = qc_id index += 1 @@ -436,9 +520,9 @@ def convert_point_data(self): self.hdr_typ = [ 0 for _ in range(0, hdr_cnt) ] self.hdr_sid = [ 0 for _ in range(0, hdr_cnt) ] self.hdr_vld = [ 0 for _ in range(0, hdr_cnt) ] - self.hdr_lat = [ nc_tools.FILL_VALUE for _ in range(0, hdr_cnt) ] - self.hdr_lon = [ nc_tools.FILL_VALUE for _ in range(0, hdr_cnt) ] - self.hdr_elv = [ nc_tools.FILL_VALUE for _ in range(0, hdr_cnt) ] + self.hdr_lat = [ self.FILL_VALUE for _ in range(0, hdr_cnt) ] + self.hdr_lon = [ self.FILL_VALUE for _ in range(0, hdr_cnt) ] + self.hdr_elv = [ self.FILL_VALUE for _ in range(0, hdr_cnt) ] for key, idx in hdr_map.items(): self.hdr_typ[idx] = key[0] self.hdr_sid[idx] = key[1] @@ -463,16 +547,15 @@ def convert_point_data(self): for key, idx in obs_var_map.items(): self.obs_var_table[idx] = key - return self.get_point_data() - def get_num_value(self, column_value): num_value = column_value if isinstance(column_value, str): - if column_value.lower() == 'na': - num_value = nc_tools.FILL_VALUE - elif self.is_number(column_value): + if self.is_number(column_value): num_value = float(column_value) - num_value = nc_tools.FILL_VALUE + else: + num_value = self.FILL_VALUE + if column_value.lower() != 'na' and column_value.lower() != 'n/a': + self.log_info(f'{column_value} is not a number, converted to the missing value') return num_value def is_grib_code(self): @@ -488,7 +571,7 @@ def is_grib_code(self): def is_num_string(self, column_value): is_string = isinstance(column_value, str) if is_string: - is_num = True if self.is_number(column_value) or column_value.lower() == 'na' else False + is_num = True if self.is_number(column_value) or column_value.lower() == 'na' or column_value.lower() == 'n/a' else False else: is_num = True return is_string, is_num @@ -496,6 +579,8 @@ def is_num_string(self, column_value): class met_point_obs(ABC, base_met_point_obs): + MET_ENV_RUN = 'MET_FORCE_TO_RUN' + @abstractmethod def read_data(self, args): # args can be input_file_name, list, or dictionary @@ -514,31 +599,13 @@ def read_data(self, args): ''' pass - -# Note: caller should import netCDF4 -# the argements nc_group(dataset) and nc_var should not be None -class nc_tools(): - - FILL_VALUE = -9999. - met_missing = -99999999. - - @staticmethod - def get_num_array(nc_group, var_name): - nc_var = nc_group.variables.get(var_name, None) - #return [] if nc_var is None else nc_var[:].filled(nc_var._FillValue if '_FillValue' in nc_var.ncattrs() else nc_tools.met_missing) - return [] if nc_var is None else nc_var[:] - @staticmethod - def get_ncbyte_array_to_str(nc_var): - nc_str_data = nc_var[:] - if nc_var.datatype.name == 'bytes8': - nc_str_data = [ str(s.compressed(),"utf-8") for s in nc_var[:] ] - return nc_str_data + def get_prompt(): + return get_prompt() @staticmethod - def get_string_array(nc_group, var_name): - nc_var = nc_group.variables.get(var_name, None) - return [] if nc_var is None else nc_tools.get_ncbyte_array_to_str(nc_var) + def is_python_prefix(user_cmd): + return user_cmd.startswith(base_met_point_obs.python_prefix) # This is a sample drived class @@ -568,12 +635,15 @@ def read_data(self, arg_map={}): self.obs_var_table = [ "TMP", "RH" ] self.obs_qty_table = [ "NA" ] -def convert_point_data(point_data, check_all_records=False): - _csv_point_data = csv_point_obs(point_data) - if _csv_point_data.is_grib_code(): - _csv_point_data.use_var_id = False - _csv_point_data.check_csv_point_data(check_all_records) - return _csv_point_data.convert_point_data() +def convert_point_data(point_data, check_all_records=False, input_type='csv'): + tmp_point_data = {} + if 'csv' == input_type: + csv_point_data = csv_point_obs(point_data) + csv_point_data.check_csv_point_data(check_all_records) + tmp_point_data = csv_point_data.get_point_data() + else: + base_met_point_obs.error_msg('Not supported input type: {input_type}') + return tmp_point_data def main(): args = {} # or args = [] diff --git a/scripts/python/met_point_obs_nc.py b/scripts/python/met_point_obs_nc.py new file mode 100644 index 0000000000..cb86be5011 --- /dev/null +++ b/scripts/python/met_point_obs_nc.py @@ -0,0 +1,276 @@ +#!/usr/bin/env python3 + +''' +Separated from read_met_point_obs on Feb 09, 2023 + +@author: hsoh + +This script reads the MET point observation NetCDF file like MET tools do. +''' + +import os +import sys +from datetime import datetime +import numpy as np +import netCDF4 as nc + +from met_point_obs import met_point_obs, base_met_point_obs, get_prompt + +DO_PRINT_DATA = False +ARG_PRINT_DATA = 'show_data' + +# Note: caller should import netCDF4 +# the argements nc_group(dataset) and nc_var should not be None +class nc_tools(): + + met_missing = -99999999. + + @staticmethod + def get_num_array(nc_group, var_name): + nc_var = nc_group.variables.get(var_name, None) + return [] if nc_var is None else nc_var[:] + + @staticmethod + def get_ncbyte_array_to_str(nc_var): + nc_str_data = nc_var[:] + if nc_var.datatype.name == 'bytes8': + nc_str_data = [ str(s.compressed(),"utf-8") for s in nc_var[:] ] + return nc_str_data + + @staticmethod + def get_string_array(nc_group, var_name): + nc_var = nc_group.variables.get(var_name, None) + return [] if nc_var is None else nc_tools.get_ncbyte_array_to_str(nc_var) + + +class nc_point_obs(met_point_obs): + + # args should be string, list, or dictionary + def get_nc_filename(self, args): + nc_filename = None + if isinstance(args, dict): + nc_filename = args.get('nc_name',None) + elif isinstance(args, list): + nc_filename = args[0] + elif args != ARG_PRINT_DATA: + nc_filename = args + + return nc_filename + + def read_data(self, nc_filename): + if nc_filename is None: + self.log_error_msg("The input NetCDF filename is missing") + elif not os.path.exists(nc_filename): + self.log_error_msg(f"input NetCDF file ({nc_filename}) does not exist") + else: + dataset = nc.Dataset(nc_filename, 'r') + # Header + self.hdr_typ = dataset['hdr_typ'][:] + self.hdr_sid = dataset['hdr_sid'][:] + self.hdr_vld = dataset['hdr_vld'][:] + self.hdr_lat = dataset['hdr_lat'][:] + self.hdr_lon = dataset['hdr_lon'][:] + self.hdr_elv = dataset['hdr_elv'][:] + self.hdr_typ_table = nc_tools.get_string_array(dataset, 'hdr_typ_table') + self.hdr_sid_table = nc_tools.get_string_array(dataset, 'hdr_sid_table') + self.hdr_vld_table = nc_tools.get_string_array(dataset, 'hdr_vld_table') + + nc_var = dataset.variables.get('obs_unit', None) + if nc_var: + self.obs_var_unit = nc_var[:] + nc_var = dataset.variables.get('obs_desc', None) + if nc_var: + self.obs_var_desc = nc_var[:] + + nc_var = dataset.variables.get('hdr_prpt_typ', None) + if nc_var: + self.hdr_prpt_typ = nc_var[:] + nc_var = dataset.variables.get('hdr_irpt_typ', None) + if nc_var: + self.hdr_irpt_typ = nc_var[:] + nc_var = dataset.variables.get('hdr_inst_typ', None) + if nc_var: + self.hdr_inst_typ =nc_var[:] + + #Observation data + self.hdr_sid = dataset['hdr_sid'][:] + self.obs_qty = np.array(dataset['obs_qty'][:]) + self.obs_hid = np.array(dataset['obs_hid'][:]) + self.obs_lvl = np.array(dataset['obs_lvl'][:]) + self.obs_hgt = np.array(dataset['obs_hgt'][:]) + self.obs_val = np.array(dataset['obs_val'][:]) + nc_var = dataset.variables.get('obs_vid', None) + if nc_var is None: + self.use_var_id = False + nc_var = dataset.variables.get('obs_gc', None) + else: + self.obs_var_table = nc_tools.get_string_array(dataset, 'obs_var') + if nc_var: + self.obs_vid = np.array(nc_var[:]) + + self.obs_qty_table = nc_tools.get_string_array(dataset, 'obs_qty_table') + + def save_ncfile(self, nc_filename): + met_data = self.get_point_data() + with nc.Dataset(nc_filename, 'w') as nc_dataset: + self.set_nc_data(nc_dataset) + return met_data + + def set_nc_data(self, nc_dataset): + return nc_point_obs.write_nc_data(nc_dataset, self) + + @staticmethod + def write_nc_file(nc_filename, point_obs): + with nc.Dataset(nc_filename, 'w') as nc_dataset: + nc_point_obs.set_nc_data(nc_dataset, point_obs) + + @staticmethod + def write_nc_data(nc_dataset, point_obs): + do_nothing = False + if 0 == point_obs.nhdr: + do_nothing = True + base_met_point_obs.info_msg("the header is empty") + if 0 == point_obs.nobs: + do_nothing = True + base_met_point_obs.info_msg("the observation data is empty") + if do_nothing: + print() + return + + # Set global attributes + nc_dataset.MET_Obs_version = "1.02" ; + nc_dataset.use_var_id = "true" if point_obs.use_var_id else "false" + + # Create dimensions + nc_dataset.createDimension('mxstr', 16) + nc_dataset.createDimension('mxstr2', 40) + nc_dataset.createDimension('mxstr3', 80) + nc_dataset.createDimension('nhdr', point_obs.nhdr) + nc_dataset.createDimension('nobs', point_obs.nobs) + #npbhdr = len(point_obs.hdr_prpt_typ) + if 0 < point_obs.npbhdr: + nc_dataset.createDimension('npbhdr', point_obs.npbhdr) + nc_dataset.createDimension('nhdr_typ', point_obs.nhdr_typ) + nc_dataset.createDimension('nhdr_sid', point_obs.nhdr_sid) + nc_dataset.createDimension('nhdr_vld', point_obs.nhdr_vld) + nc_dataset.createDimension('nobs_qty', point_obs.nobs_qty) + nc_dataset.createDimension('obs_var_num', point_obs.nobs_var) + + type_for_string = 'S1' # np.byte + dims_hdr = ('nhdr',) + dims_obs = ('nobs',) + + # Create header and observation variables + var_hdr_typ = nc_dataset.createVariable('hdr_typ', np.int32, dims_hdr, fill_value=-9999) + var_hdr_sid = nc_dataset.createVariable('hdr_sid', np.int32, dims_hdr, fill_value=-9999) + var_hdr_vld = nc_dataset.createVariable('hdr_vld', np.int32, dims_hdr, fill_value=-9999) + var_hdr_lat = nc_dataset.createVariable('hdr_lat', np.float32, dims_hdr, fill_value=-9999.) + var_hdr_lon = nc_dataset.createVariable('hdr_lon', np.float32, dims_hdr, fill_value=-9999.) + var_hdr_elv = nc_dataset.createVariable('hdr_elv', np.float32, dims_hdr, fill_value=-9999.) + + var_obs_qty = nc_dataset.createVariable('obs_qty', np.int32, dims_obs, fill_value=-9999) + var_obs_hid = nc_dataset.createVariable('obs_hid', np.int32, dims_obs, fill_value=-9999) + var_obs_vid = nc_dataset.createVariable('obs_vid', np.int32, dims_obs, fill_value=-9999) + var_obs_lvl = nc_dataset.createVariable('obs_lvl', np.float32, dims_obs, fill_value=-9999.) + var_obs_hgt = nc_dataset.createVariable('obs_hgt', np.float32, dims_obs, fill_value=-9999.) + var_obs_val = nc_dataset.createVariable('obs_val', np.float32, dims_obs, fill_value=-9999.) + + if 0 == point_obs.npbhdr: + var_hdr_prpt_typ = None + var_hdr_irpt_typ = None + var_hdr_inst_typ = None + else: + dims_npbhdr = ('npbhdr',) + var_hdr_prpt_typ = nc_dataset.createVariable('hdr_prpt_typ', np.int32, dims_npbhdr, fill_value=-9999.) + var_hdr_irpt_typ = nc_dataset.createVariable('hdr_irpt_typ', np.int32, dims_npbhdr, fill_value=-9999.) + var_hdr_inst_typ = nc_dataset.createVariable('hdr_inst_typ', np.int32, dims_npbhdr, fill_value=-9999.) + + var_hdr_typ_table = nc_dataset.createVariable('hdr_typ_table', type_for_string, ('nhdr_typ','mxstr2')) + var_hdr_sid_table = nc_dataset.createVariable('hdr_sid_table', type_for_string, ('nhdr_sid','mxstr2')) + var_hdr_vld_table = nc_dataset.createVariable('hdr_vld_table', type_for_string, ('nhdr_vld','mxstr')) + var_obs_qty_table = nc_dataset.createVariable('obs_qty_table', type_for_string, ('nobs_qty','mxstr')) + var_obs_var_table = nc_dataset.createVariable('obs_var', type_for_string, ('obs_var_num','mxstr2')) + var_obs_var_unit = nc_dataset.createVariable('obs_unit', type_for_string, ('obs_var_num','mxstr2')) + var_obs_var_desc = nc_dataset.createVariable('obs_desc', type_for_string, ('obs_var_num','mxstr3')) + + # Set variables + var_hdr_typ[:] = point_obs.hdr_typ[:] + var_hdr_sid[:] = point_obs.hdr_sid[:] + var_hdr_vld[:] = point_obs.hdr_vld[:] + var_hdr_lat[:] = point_obs.hdr_lat[:] + var_hdr_lon[:] = point_obs.hdr_lon[:] + var_hdr_elv[:] = point_obs.hdr_elv[:] + for i in range(0, point_obs.nhdr_typ): + for j in range(0, len(point_obs.hdr_typ_table[i])): + var_hdr_typ_table[i,j] = point_obs.hdr_typ_table[i][j] + for i in range(0, point_obs.nhdr_sid): + for j in range(0, len(point_obs.hdr_sid_table[i])): + var_hdr_sid_table[i,j] = point_obs.hdr_sid_table[i][j] + for i in range(0, point_obs.nhdr_vld): + for j in range(0, len(point_obs.hdr_vld_table[i])): + var_hdr_vld_table[i,j] = point_obs.hdr_vld_table[i][j] + if 0 < point_obs.npbhdr: + var_hdr_prpt_typ[:] = point_obs.hdr_prpt_typ[:] + var_hdr_irpt_typ[:] = point_obs.hdr_irpt_typ[:] + var_hdr_inst_typ[:] = point_obs.hdr_inst_typ[:] + + var_obs_qty[:] = point_obs.obs_qty[:] + var_obs_hid[:] = point_obs.obs_hid[:] + var_obs_vid[:] = point_obs.obs_vid[:] + var_obs_lvl[:] = point_obs.obs_lvl[:] + var_obs_hgt[:] = point_obs.obs_hgt[:] + var_obs_val[:] = point_obs.obs_val[:] + for i in range(0, point_obs.nobs_var): + for j in range(0, len(point_obs.obs_var_table[i])): + var_obs_var_table[i,j] = point_obs.obs_var_table[i][j] + var_obs_var_unit[i] = "" if i >= len(point_obs.obs_var_unit) else point_obs.obs_var_unit[i] + var_obs_var_desc[i] = "" if i >= len(point_obs.obs_var_desc) else point_obs.obs_var_desc[i] + for i in range(0, point_obs.nobs_qty): + for j in range(0, len(point_obs.obs_qty_table[i])): + var_obs_qty_table[i,j] = point_obs.obs_qty_table[i][j] + + # Set variable attributes + var_hdr_typ.long_name = "index of message type" + var_hdr_sid.long_name = "index of station identification" + var_hdr_vld.long_name = "index of valid time" + var_hdr_lat.long_name = "latitude" + var_hdr_lat.units = "degrees_north" + var_hdr_lon.long_name = "longitude" + var_hdr_lon.units = "degrees_east" + var_hdr_elv.long_name = "elevation" + var_hdr_elv.units = "meters above sea level (msl)" + + var_obs_qty.long_name = "index of quality flag" + var_obs_hid.long_name = "index of matching header data" + var_obs_vid.long_name = "index of BUFR variable corresponding to the observation type" + var_obs_lvl.long_name = "pressure level (hPa) or accumulation interval (sec)" + var_obs_hgt.long_name = "height in meters above sea level (msl)" + var_obs_val.long_name = "observation value" + var_hdr_typ_table.long_name = "message type" + var_hdr_sid_table.long_name = "station identification" + var_hdr_vld_table.long_name = "valid time" + var_hdr_vld_table.units = "YYYYMMDD_HHMMSS UTC" + var_obs_qty_table.long_name = "quality flag" + var_obs_var_table.long_name = "variable names" + var_obs_var_unit.long_name = "variable units" + var_obs_var_desc.long_name = "variable descriptions" + + +def main(argv): + if len(argv) != 1 and argv[1] != ARG_PRINT_DATA: + netcdf_filename = argv[1] + tmp_nc_name = 'tmp_met_point.nc' + point_obs_data = nc_point_obs() + point_obs_data.read_data(point_obs_data.get_nc_filename(netcdf_filename)) + met_point_data = point_obs_data.save_ncfile(tmp_nc_name) + print(f'{get_prompt()} saved met_point_data to {tmp_nc_name}') + met_point_data['met_point_data'] = point_obs_data + + if DO_PRINT_DATA or ARG_PRINT_DATA == argv[-1]: + met_point_obs.print_point_data(met_point_data) + +if __name__ == '__main__': + start_time = datetime.now() + main(sys.argv) + run_time = datetime.now() - start_time + print(f'{get_prompt()} Done python script {sys.argv[0]} took {run_time}') diff --git a/scripts/python/read_met_point_obs.py b/scripts/python/read_met_point_obs.py index 60f2c62065..57ccd22e7a 100755 --- a/scripts/python/read_met_point_obs.py +++ b/scripts/python/read_met_point_obs.py @@ -1,121 +1,91 @@ +#!/usr/bin/env python3 ''' Created on Nov 10, 2021 @author: hsoh This script reads the MET point observation NetCDF file like MET tools do. + +Usage: + + python3 read_met_point_obs.py + python3 read_met_point_obs.py + : 11 columns + 'typ', 'sid', 'vld', 'lat', 'lon', 'elv', 'var', 'lvl', 'hgt', 'qc', 'obs' + string columns: 'typ', 'sid', 'vld', 'var', , 'qc' + numeric columns 'lat', 'lon', 'elv', 'lvl', 'hgt', 'qc', 'obs' + python3 read_met_point_obs.py + ''' import os import sys from datetime import datetime -import numpy as np -import netCDF4 as nc + +met_base_dir = os.getenv('MET_BASE',None) +if met_base_dir is not None: + sys.path.append(os.path.join(met_base_dir, 'python')) from met_point_obs import met_point_obs, sample_met_point_obs +from met_point_obs_nc import nc_point_obs DO_PRINT_DATA = False ARG_PRINT_DATA = 'show_data' -class nc_point_obs(met_point_obs): - - def read_data(self, args): - # args should be list or dictionaryu - if isinstance(args, dict): - nc_filename = args.get('nc_name',None) - elif isinstance(args, list): - nc_filename = args[0] - else: - nc_filename = args - if nc_filename is None: - print("==ERROR== The input NetCDF filename is missing") - elif not os.path.exists(nc_filename): - print("==ERROR== input NetCDF file ({f}) does not exist".format(f=nc_filename)) - else: - dataset = nc.Dataset(nc_filename, 'r') - # Header - self.hdr_typ = dataset['hdr_typ'][:] - self.hdr_sid = dataset['hdr_sid'][:] - self.hdr_vld = dataset['hdr_vld'][:] - self.hdr_lat = dataset['hdr_lat'][:] - self.hdr_lon = dataset['hdr_lon'][:] - self.hdr_elv = dataset['hdr_elv'][:] - self.hdr_typ_table = nc_tools.get_string_array(dataset, 'hdr_typ_table') - self.hdr_sid_table = nc_tools.get_string_array(dataset, 'hdr_sid_table') - self.hdr_vld_table = nc_tools.get_string_array(dataset, 'hdr_vld_table') - - nc_var = dataset.variables.get('hdr_prpt_typ', None) - if nc_var: - self.hdr_prpt_typ = nc_var[:] - nc_var = dataset.variables.get('hdr_irpt_typ', None) - if nc_var: - self.hdr_irpt_typ = nc_var[:] - nc_var = dataset.variables.get('hdr_inst_typ', None) - if nc_var: - self.hdr_inst_typ =nc_var[:] - - #Observation data - self.hdr_sid = dataset['hdr_sid'][:] - self.obs_qty = np.array(dataset['obs_qty'][:]) - self.obs_hid = np.array(dataset['obs_hid'][:]) - self.obs_lvl = np.array(dataset['obs_lvl'][:]) - self.obs_hgt = np.array(dataset['obs_hgt'][:]) - self.obs_val = np.array(dataset['obs_val'][:]) - nc_var = dataset.variables.get('obs_vid', None) - if nc_var is None: - self.use_var_id = False - nc_var = dataset.variables.get('obs_gc', None) - else: - self.obs_var_table = nc_tools.get_string_array(dataset, 'obs_var') - if nc_var: - self.obs_vid = np.array(nc_var[:]) - - self.obs_qty_table = nc_tools.get_string_array(dataset, 'obs_qty_table') - - -class nc_tools(): - - @staticmethod - def get_num_array(dataset, var_name): - nc_var = dataset.variables.get(var_name, None) - return nc_var[:].filled(nc_var._FillValue if '_FillValue' in nc_var.ncattrs() else -9999) if nc_var else [] - - @staticmethod - def get_ncbyte_array_to_str(nc_var): - nc_str_data = nc_var[:] - if nc_var.datatype.name == 'bytes8': - nc_str_data = [ str(s.compressed(),"utf-8") for s in nc_var[:]] - return nc_str_data - - @staticmethod - def get_string_array(dataset, var_name): - nc_var = dataset.variables.get(var_name, None) - return nc_tools.get_ncbyte_array_to_str(nc_var) if nc_var else [] - - start_time = datetime.now() +prompt = met_point_obs.get_prompt() point_obs_data = None -if len(sys.argv) == 1: +if len(sys.argv) == 1 or ARG_PRINT_DATA == sys.argv[1]: point_obs_data = sample_met_point_obs() point_obs_data.read_data([]) +elif met_point_obs.is_python_prefix(sys.argv[1]): + import importlib.util + + print("{p} Python Script:\t".format(p=prompt) + repr(sys.argv[0])) + print("{p} User Command:\t".format(p=prompt) + repr(' '.join(sys.argv[2:]))) + + pyembed_module_name = sys.argv[2] + sys.argv = sys.argv[1:] + + # append user script dir to system path + pyembed_dir, pyembed_file = os.path.split(pyembed_module_name) + if pyembed_dir: + sys.path.insert(0, pyembed_dir) + + if not pyembed_module_name.endswith('.py'): + pyembed_module_name += '.py' + os.environ[met_point_obs.MET_ENV_RUN] = "TRUE" + + user_base = os.path.basename(pyembed_module_name).replace('.py','') + + spec = importlib.util.spec_from_file_location(user_base, pyembed_module_name) + met_in = importlib.util.module_from_spec(spec) + spec.loader.exec_module(met_in) + + met_point_obs = met_in.met_point_obs + print("met_point_obs: ", met_point_obs) + met_point_data = met_in.met_point_data + print("met_point_data: ", met_point_data) + #print(hasattr("met_in: ", dir(met_in))) + #met_point_data = met_point_obs.get_point_data() + #met_point_data = None if met_in.get('met_point_data', None) else met_in.met_point_data + #met_data = None if met_in.get('met_data', None) else met_in.met_data + print(met_point_data) else: netcdf_filename = sys.argv[1] args = [ netcdf_filename ] #args = { 'nc_name': netcdf_filename } point_obs_data = nc_point_obs() - point_obs_data.read_data(args) - - if ARG_PRINT_DATA == sys.argv[-1]: - DO_PRINT_DATA = True + point_obs_data.read_data(point_obs_data.get_nc_filename(args)) if point_obs_data is not None: met_point_data = point_obs_data.get_point_data() met_point_data['met_point_data'] = point_obs_data - if DO_PRINT_DATA: - met_point_obs.print_point_data(met_point_data) + if DO_PRINT_DATA or ARG_PRINT_DATA == sys.argv[-1]: + point_obs_data.dump() run_time = datetime.now() - start_time -print('Done python script {s} took {t}'.format(s=sys.argv[0], t=run_time)) +print('{p} Done python script {s} took {t}'.format(p=prompt, s=sys.argv[0], t=run_time)) diff --git a/src/libcode/vx_data2d_python/python_dataplane.cc b/src/libcode/vx_data2d_python/python_dataplane.cc index f588db9c03..af886d7efe 100644 --- a/src/libcode/vx_data2d_python/python_dataplane.cc +++ b/src/libcode/vx_data2d_python/python_dataplane.cc @@ -23,7 +23,7 @@ //////////////////////////////////////////////////////////////////////// -GlobalPython GP; // this needs external linkage +extern GlobalPython GP; // this needs external linkage //////////////////////////////////////////////////////////////////////// @@ -35,12 +35,8 @@ static const char write_tmp_nc [] = "MET_BASE/wrappers/write_tmp_datapla static const char read_tmp_nc [] = "read_tmp_dataplane"; // NO ".py" suffix -static const char tmp_nc_base_name [] = "tmp_met_nc"; - static const char tmp_nc_var_name [] = "met_info"; -static const char tmp_nc_file_var_name [] = "tmp_nc_filename"; - //////////////////////////////////////////////////////////////////////// diff --git a/src/libcode/vx_pointdata_python/python_pointdata.cc b/src/libcode/vx_pointdata_python/python_pointdata.cc index e99cc49596..e8bd3731d2 100644 --- a/src/libcode/vx_pointdata_python/python_pointdata.cc +++ b/src/libcode/vx_pointdata_python/python_pointdata.cc @@ -13,15 +13,30 @@ #include "python_pointdata.h" #include "pointdata_from_array.h" #include "vx_util.h" +#include "vx_python3_utils.h" #include "global_python.h" #include "wchar_argv.h" //////////////////////////////////////////////////////////////////////// +extern GlobalPython GP; // this needs external linkage + //////////////////////////////////////////////////////////////////////// +static const char * user_ppath = 0; + +static const char write_tmp_nc [] = "MET_BASE/wrappers/write_tmp_point_nc.py"; + +static const char read_tmp_nc [] = "read_tmp_point_nc"; // NO ".py" suffix + +//////////////////////////////////////////////////////////////////////// + + +static bool tmp_nc_point_obs(const char * script_name, int user_script_argc, + char ** user_script_argv, MetPointDataPython &met_pd_out); + static bool straight_python_point_data(const char * script_name, int script_argc, char ** script_argv, const bool use_xarray, MetPointDataPython &met_pd_out); @@ -52,13 +67,172 @@ bool python_point_data(const char * script_name, int script_argc, char ** script { -bool status = straight_python_point_data(script_name, script_argc, script_argv, - use_xarray, met_pd_out); +bool status = false; + +if ( user_ppath == 0 ) user_ppath = getenv(user_python_path_env); + +if ( user_ppath != 0 ) { // do_tmp_nc = true; + + status = tmp_nc_point_obs(script_name, script_argc, script_argv, + met_pd_out); +} +else { + status = straight_python_point_data(script_name, script_argc, script_argv, + use_xarray, met_pd_out); +} return ( status ); } +//////////////////////////////////////////////////////////////////////// + +bool process_python_point_data(PyObject *module_obj, MetPointDataPython &met_pd_out) +{ + +int int_value; +PyObject *module_dict_obj = 0; +PyObject *python_value = 0; +PyObject *python_met_point_data = 0; +ConcatString cs, user_dir, user_base; +const char *method_name = "process_python_point_data -> "; +const char *method_name_s = "process_python_point_data()"; + + // + // get the namespace for the module (as a dictionary) + // + +module_dict_obj = PyModule_GetDict (module_obj); + + // + // get handles to the objects of interest from the module_dict + // + +python_met_point_data = PyDict_GetItemString (module_dict_obj, python_key_point_data); + +python_value = PyDict_GetItemString (python_met_point_data, python_use_var_id); + +bool use_var_id = pyobject_as_bool(python_value); +met_pd_out.set_use_var_id(use_var_id); +mlog << Debug(9) << method_name << "use_var_id: \"" << use_var_id + << "\" from python. is_using_var_id(): " << met_pd_out.is_using_var_id() << "\n"; + +python_value = PyDict_GetItemString (python_met_point_data, python_key_nhdr); + +int_value = pyobject_as_int(python_value); +if (int_value == 0) { + mlog << Error << "\n" << method_name + << "The header is empty. Please check if python input exists\n\n"; + exit (1); +} +met_pd_out.set_hdr_cnt(int_value); + +python_value = PyDict_GetItemString (python_met_point_data, python_key_nobs); +int_value = pyobject_as_int(python_value); +if (int_value == 0) { + mlog << Error << "\n" << method_name + << "The point observation data is empty. Please check if python input is processed properly\n\n"; + exit (1); +} + +met_pd_out.allocate(int_value); + +MetPointObsData *obs_data = met_pd_out.get_point_obs_data(); +MetPointHeader *header_data = met_pd_out.get_header_data(); + + + // look up the data array variable name from the dictionary + + set_array_from_python(python_met_point_data, numpy_array_hdr_typ, &header_data->typ_idx_array); + set_array_from_python(python_met_point_data, numpy_array_hdr_sid, &header_data->sid_idx_array); + set_array_from_python(python_met_point_data, numpy_array_hdr_vld, &header_data->vld_idx_array); + set_array_from_python(python_met_point_data, numpy_array_hdr_lat, &header_data->lat_array); + set_array_from_python(python_met_point_data, numpy_array_hdr_lon, &header_data->lon_array); + set_array_from_python(python_met_point_data, numpy_array_hdr_elv, &header_data->elv_array); + if (header_data->typ_idx_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_typ is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->sid_idx_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_sid is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->vld_idx_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_vld is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->lat_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_lat is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->lon_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_lon is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->elv_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_elv is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + + set_str_array_from_python(python_met_point_data, numpy_array_hdr_typ_table, &header_data->typ_array); + set_str_array_from_python(python_met_point_data, numpy_array_hdr_sid_table, &header_data->sid_array); + set_str_array_from_python(python_met_point_data, numpy_array_hdr_vld_table, &header_data->vld_array); + if (header_data->typ_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_typ_table is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->sid_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_sid_table is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (header_data->vld_array.n() == 0) { + mlog << Error << "\n" << method_name + << "The hdr_vld_table is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + set_array_from_python(python_met_point_data, numpy_array_prpt_typ_table, &header_data->prpt_typ_array, false); + set_array_from_python(python_met_point_data, numpy_array_irpt_typ_table, &header_data->irpt_typ_array, false); + set_array_from_python(python_met_point_data, numpy_array_inst_typ_table, &header_data->inst_typ_array, false); + + set_array_from_python(python_met_point_data, numpy_array_obs_qty, obs_data->obs_qids); + set_array_from_python(python_met_point_data, numpy_array_obs_hid, obs_data->obs_hids); + set_array_from_python(python_met_point_data, numpy_array_obs_vid, obs_data->obs_ids); + set_array_from_python(python_met_point_data, numpy_array_obs_lvl, obs_data->obs_lvls); + set_array_from_python(python_met_point_data, numpy_array_obs_hgt, obs_data->obs_hgts); + set_array_from_python(python_met_point_data, numpy_array_obs_val, obs_data->obs_vals); + + set_str_array_from_python(python_met_point_data, numpy_array_obs_qty_table, &obs_data->qty_names); + set_str_array_from_python(python_met_point_data, numpy_array_obs_var_table, &obs_data->var_names); + if (obs_data->qty_names.n() == 0) { + mlog << Error << "\n" << method_name + << "The obs_qty_table is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + if (use_var_id && obs_data->var_names.n() == 0) { + mlog << Error << "\n" << method_name + << "The obs_var_table is empty. Please check if python input is processed properly\n\n"; + exit (1); + } + + if(mlog.verbosity_level()>=point_data_debug_level) print_met_data(obs_data, header_data, method_name_s); + + // + // done + // + +return ( true ); + +} + + //////////////////////////////////////////////////////////////////////// @@ -170,153 +344,175 @@ if ( ! module_obj ) { } + +return process_python_point_data(module_obj, met_pd_out); + +} + + +//////////////////////////////////////////////////////////////////////// + + +bool tmp_nc_point_obs(const char * user_script_name, int user_script_argc, + char ** user_script_argv, MetPointDataPython &met_pd_out) + +{ + +int j; +int status; +ConcatString command; +ConcatString path; +ConcatString tmp_nc_path; +const char * tmp_dir = 0; +Wchar_Argv wa; +const char *method_name = "tmp_nc_point_obs() -> "; + // - // get the namespace for the module (as a dictionary) + // if the global python object has already been initialized, + // we need to reload the module // -module_dict_obj = PyModule_GetDict (module_obj); +bool do_reload = GP.is_initialized; + +GP.initialize(); // - // get handles to the objects of interest from the module_dict + // start up the python interpreter // -python_met_point_data = PyDict_GetItemString (module_dict_obj, python_key_point_data); +if ( PyErr_Occurred() ) { -python_value = PyDict_GetItemString (python_met_point_data, python_use_var_id); + PyErr_Print(); -bool use_var_id = pyobject_as_bool(python_value); -met_pd_out.set_use_var_id(use_var_id); -mlog << Debug(9) << method_name << "use_var_id: \"" << use_var_id << "\" from python. is_using_var_id(): " << met_pd_out.is_using_var_id() << "\n"; + mlog << Warning << "\n" << method_name + << "an error occurred initializing python\n\n"; -python_value = PyDict_GetItemString (python_met_point_data, python_key_nhdr); + return ( false ); -int_value = pyobject_as_int(python_value); -if (int_value == 0) { - mlog << Error << "\n" << method_name - << "The header is empty. Please check if python input exists\n\n"; - exit (1); } -met_pd_out.set_hdr_cnt(int_value); -python_value = PyDict_GetItemString (python_met_point_data, python_key_nobs); -int_value = pyobject_as_int(python_value); -if (int_value == 0) { +run_python_string("import sys"); +command << cs_erase + << "sys.path.append(\"" + << replace_path(python_dir) + << "\")"; +run_python_string(command.text()); + +mlog << Debug(3) << "Calling " << user_ppath + << " to run user's python script (" << user_script_name + << ").\n"; + +tmp_dir = getenv ("MET_TMP_DIR"); + +if ( ! tmp_dir ) tmp_dir = default_tmp_dir; + +path << cs_erase + << tmp_dir << '/' + << tmp_nc_base_name; + +tmp_nc_path = make_temp_file_name(path.text(), 0); + +command << cs_erase + << user_ppath << ' ' // user's path to python + << replace_path(write_tmp_nc) << ' ' // write_tmp_nc.py + << tmp_nc_path << ' ' // tmp_nc output filename + << user_script_name; // user's script name + +for (j=1; jtyp_idx_array); - set_array_from_python(python_met_point_data, numpy_array_hdr_sid, &header_data->sid_idx_array); - set_array_from_python(python_met_point_data, numpy_array_hdr_vld, &header_data->vld_idx_array); - set_array_from_python(python_met_point_data, numpy_array_hdr_lat, &header_data->lat_array); - set_array_from_python(python_met_point_data, numpy_array_hdr_lon, &header_data->lon_array); - set_array_from_python(python_met_point_data, numpy_array_hdr_elv, &header_data->elv_array); - if (header_data->typ_idx_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_typ is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->sid_idx_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_sid is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->vld_idx_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_vld is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->lat_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_lat is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->lon_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_lon is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->elv_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_elv is empty. Please check if python input is processed properly\n\n"; - exit (1); - } + module_obj = PyImport_ReloadModule (module_obj); - set_str_array_from_python(python_met_point_data, numpy_array_hdr_typ_table, &header_data->typ_array); - set_str_array_from_python(python_met_point_data, numpy_array_hdr_sid_table, &header_data->sid_array); - set_str_array_from_python(python_met_point_data, numpy_array_hdr_vld_table, &header_data->vld_array); - if (header_data->typ_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_typ_table is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->sid_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_sid_table is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (header_data->vld_array.n() == 0) { - mlog << Error << "\n" << method_name - << "The hdr_vld_table is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - set_array_from_python(python_met_point_data, numpy_array_prpt_typ_table, &header_data->prpt_typ_array, false); - set_array_from_python(python_met_point_data, numpy_array_irpt_typ_table, &header_data->irpt_typ_array, false); - set_array_from_python(python_met_point_data, numpy_array_inst_typ_table, &header_data->inst_typ_array, false); +} - set_array_from_python(python_met_point_data, numpy_array_obs_qty, obs_data->obs_qids); - set_array_from_python(python_met_point_data, numpy_array_obs_hid, obs_data->obs_hids); - set_array_from_python(python_met_point_data, numpy_array_obs_vid, obs_data->obs_ids); - set_array_from_python(python_met_point_data, numpy_array_obs_lvl, obs_data->obs_lvls); - set_array_from_python(python_met_point_data, numpy_array_obs_hgt, obs_data->obs_hgts); - set_array_from_python(python_met_point_data, numpy_array_obs_val, obs_data->obs_vals); +if ( PyErr_Occurred() ) { - set_str_array_from_python(python_met_point_data, numpy_array_obs_qty_table, &obs_data->qty_names); - set_str_array_from_python(python_met_point_data, numpy_array_obs_var_table, &obs_data->var_names); - if (obs_data->qty_names.n() == 0) { - mlog << Error << "\n" << method_name - << "The obs_qty_table is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - if (use_var_id && obs_data->var_names.n() == 0) { - mlog << Error << "\n" << method_name - << "The obs_var_table is empty. Please check if python input is processed properly\n\n"; - exit (1); - } - - if(mlog.verbosity_level()>=point_data_debug_level) print_met_data(obs_data, header_data, method_name_s); + PyErr_Print(); + + mlog << Warning << "\n" << method_name + << "an error occurred importing module " + << '\"' << path << "\"\n\n"; + + return ( false ); } +if ( ! module_obj ) { + + mlog << Warning << "\n" << method_name + << "error running python script\n\n"; + + return ( false ); + +} + + // + // read the tmp_nc file + // + + // + // get the namespace for the module (as a dictionary) + // + + +process_python_point_data(module_obj, met_pd_out); + + + // + // cleanup + // + +remove_temp_file(tmp_nc_path); + // // done // diff --git a/src/libcode/vx_python3_utils/python3_util.cc b/src/libcode/vx_python3_utils/python3_util.cc index ab17acd7b8..f454c837b1 100644 --- a/src/libcode/vx_python3_utils/python3_util.cc +++ b/src/libcode/vx_python3_utils/python3_util.cc @@ -15,8 +15,13 @@ using namespace std; #include "vx_math.h" #include "python3_util.h" +#include "global_python.h" +//////////////////////////////////////////////////////////////////////// + +GlobalPython GP; // this needs external linkage + //////////////////////////////////////////////////////////////////////// diff --git a/src/libcode/vx_python3_utils/python3_util.h b/src/libcode/vx_python3_utils/python3_util.h index b1681463c0..e81dfe260d 100644 --- a/src/libcode/vx_python3_utils/python3_util.h +++ b/src/libcode/vx_python3_utils/python3_util.h @@ -30,6 +30,13 @@ static const char user_python_path_env [] = "MET_PYTHON_EXE"; static const char wrappers_dir [] = "MET_BASE/wrappers"; +static const char python_dir [] = "MET_BASE/python"; + +static const char tmp_nc_base_name [] = "tmp_met_nc"; + +static const char tmp_nc_file_var_name [] = "tmp_nc_filename"; + +static const char tmp_nc_point_var_name[] = "met_point_data"; //////////////////////////////////////////////////////////////////////// @@ -56,7 +63,7 @@ extern PyObject * get_attribute(PyObject *, const char * attribute_name); extern int pyobject_as_int (PyObject *); -extern bool pyobject_as_bool (PyObject *); +extern bool pyobject_as_bool (PyObject *); extern double pyobject_as_double (PyObject *); extern std::string pyobject_as_string (PyObject *); extern ConcatString pyobject_as_concat_string (PyObject *);