diff --git a/src/conventional/CMakeLists.txt b/src/conventional/CMakeLists.txt index 5f4efbe09..89d675b6a 100644 --- a/src/conventional/CMakeLists.txt +++ b/src/conventional/CMakeLists.txt @@ -10,6 +10,7 @@ list(APPEND programs sonde_bufr2ioda.py sonde_tac2ioda.py synop_bufr2ioda.py + tcv_txt2ioda.py ) file( RELATIVE_PATH SCRIPT_LIB_PATH ${CMAKE_BINARY_DIR}/bin ${PYIODACONV_BUILD_LIBDIR} ) diff --git a/src/conventional/tcv_txt2ioda.py b/src/conventional/tcv_txt2ioda.py new file mode 100755 index 000000000..2bcef541c --- /dev/null +++ b/src/conventional/tcv_txt2ioda.py @@ -0,0 +1,147 @@ +#!/usr/bin/env python3 + +# converter for tcvitals format (tropical storm central pressure) data + +import numpy as np +import netCDF4 as nc +import argparse +import pandas as pd +from datetime import datetime +import dateutil.parser +from builtins import str +from pyioda import ioda_obs_space as ioda_ospace + + +def read_file(file_name): + + names = ["center", "storm_id", "storm_name", "datestr", + "lat", "latns", "lon", "lonew", "stdir", "stspd", "pcen"] + +# ctr id name datestr lat latns + colspecs = [(0, 4), (5, 8), (9, 18), (19, 32), (33, 36), (36, 37), + (38, 42), (42, 43), (44, 47), (48, 51), (52, 56)] +# lon lonew dir spd pres + + data = pd.read_fwf(file_name, colspecs=colspecs, header=None, names=names) + + return data + + +def create_output(data, output_name): + + float_missing_value = nc.default_fillvals['f4'] + + lon = data.lon / 10.0 + lon[data.lonew == 'W'] = 360. - lon + lat = data.lat / 10.0 + lat[data.latns == 'S'] = - lat +# calculate 'obs error' as in GSI + tcp_refps = 1000.0 + tcp_width = 50.0 + tcp_ermin = 0.75 + tcp_ermax = 5.0 +# alpha=max(min(psdif/tcp_width,one),zero) + psdif = tcp_refps - data.pcen + alpha = psdif / tcp_width + alpha[alpha > 1.0] = 1.0 + alpha[alpha < 0.0] = 0.0 + oberr = tcp_ermin+(tcp_ermax-tcp_ermin)*alpha + nobs = len(lat) + height = np.zeros(nobs) + obstype = np.full(nobs, 112) + obssubtype = np.full(nobs, 0) + preqc = np.full(nobs, 2) + preuseflg = np.full(nobs, 1) + tempK = np.full(nobs, float_missing_value) + psminPa = data.pcen*100. + oberrPa = oberr*100. + epoch = datetime.fromisoformat("1970-01-01T00:00:00Z") + obdate = pd.to_datetime(data.datestr, format="%Y%m%d %H%M") + eparr = np.full_like(obdate, epoch) + obdiff = obdate - eparr + sec_since_epoch = obdiff.astype('int64') // 10**9 + station_id = data['center'].str.cat(data['storm_id'], sep='_') + dims = { + 'Location': np.arange(0, psminPa.shape[0]), + } + obsspace = ioda_ospace.ObsSpace(output_name, mode='w', dim_dict=dims) + obsspace.create_var('MetaData/dateTime', dtype='int64') \ + .write_attr('units', 'seconds since 1970-01-01T00:00:00Z') \ + .write_attr('long_name', 'dateTime') \ + .write_data(sec_since_epoch) + obsspace.create_var('MetaData/stationElevation', dtype=np.float32) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Height Of Station') \ + .write_data(height) + obsspace.create_var('MetaData/height', dtype=np.float32) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Height Of Station') \ + .write_data(height) + obsspace.create_var('MetaData/latitude', dtype=np.float32) \ + .write_attr('units', 'degrees_north') \ + .write_attr('valid_range', np.array([-90, 90], dtype=np.float32)) \ + .write_attr('long_name', 'Latitude') \ + .write_data(lat) + obsspace.create_var('MetaData/longitude', dtype=np.float32) \ + .write_attr('units', 'degrees_east') \ + .write_attr('valid_range', np.array([0, 360], dtype=np.float32)) \ + .write_attr('long_name', 'Longitude') \ + .write_data(lon) + obsspace.create_var('MetaData/pressure', dtype=np.float32) \ + .write_attr('units', 'Pa') \ + .write_attr('valid_range', np.array([105000, 80000], dtype=np.float32)) \ + .write_attr('long_name', 'Pressure') \ + .write_data(psminPa) + obsspace.create_var('MetaData/stationIdentification', dtype=station_id.dtype) \ + .write_attr('long_name', 'stationIdentification') \ + .write_data(station_id) + obsspace.create_var('ObsType/stationPressure', dtype=np.int32) \ + .write_attr('long_name', 'Station Pressure Observation Type') \ + .write_data(obstype) + + obsspace.create_var('ObsSubType/stationPressure', dtype=np.int32) \ + .write_attr('long_name', 'Station Pressure Observation subType') \ + .write_data(obssubtype) + +# ObsError: initial error values loaded from the input ioda file + obsspace.create_var('ObsError/stationPressure', dtype=np.float32) \ + .write_attr('units', 'Pa') \ + .write_attr('long_name', 'ObsError') \ + .write_attr('coordinates', 'longitude latitude') \ + .write_data(oberrPa) + + obsspace.create_var('ObsValue/stationPressure', dtype=np.float32) \ + .write_attr('units', 'Pa') \ + .write_attr('valid_range', np.array([105000, 80000], dtype=np.float32)) \ + .write_attr('long_name', 'ObsValue') \ + .write_attr('coordinates', 'longitude latitude') \ + .write_data(psminPa) + +# temporary fake array for temperature requested by SfcCorrected + obsspace.create_var('ObsValue/airTemperature', dtype=np.float32) \ + .write_attr('units', 'K') \ + .write_attr('valid_range', np.array([250, 380], dtype=np.float32)) \ + .write_attr('long_name', 'ObsValue') \ + .write_attr('coordinates', 'longitude latitude') \ + .write_data(tempK) + + +def main(): + + desc = "Reads a tcvitals (tropical storm center data) file and converts into IODA format" + parser = argparse.ArgumentParser(description=desc) + parser.add_argument( + '-i', '--input', help='Input tcvitals text file', + type=str, required=True, default=None) + parser.add_argument( + '-o', '--output', help='Output tcvitals ioda nc4 file', + type=str, required=True, default=None) + args = parser.parse_args() + + data = read_file(args.input) + + create_output(data, args.output) + + +if __name__ == '__main__': + main() diff --git a/src/gsi_ncdiag/gsi_ncdiag.py b/src/gsi_ncdiag/gsi_ncdiag.py index 18d17911c..8d41f5527 100755 --- a/src/gsi_ncdiag/gsi_ncdiag.py +++ b/src/gsi_ncdiag/gsi_ncdiag.py @@ -53,6 +53,9 @@ ], "conv_sst": [ 'sst', + ], + "conv_tcp": [ + 'tcp', ] } @@ -81,6 +84,7 @@ 722, 723, 740, 741, 742, 743, 744, 745, \ 750, 751, 752, 753, 754, 755, 786, 803, 804, 820, 821, 825], "sst": [181, 182, 183, 189, 190, 191, 192, 193, 194, 195, 196, 197, 198, 199, 200, 201, 202], + "tcp": [112], # 132 are dropsondes } @@ -195,6 +199,7 @@ "bend": ["bendingAngle"], "refract": ["atmosphericRefractivity"], "sst": ["seaSurfaceTemperature"], + "tcp": ["stationPressure"], } conv_gsivarnames = { @@ -206,6 +211,7 @@ "bend": ["Observation"], "refract": ["Observation"], "sst": ["Observation"], + "tcp": ["Observation"], } gsi_add_vars_allsky = { @@ -782,6 +788,9 @@ def toGeovals(self, OutDir, clobber=True): if (v == 'sst'): outname = OutDir + '/' + v + '_geoval_' + \ self.validtime.strftime("%Y%m%d%H") + '.nc4' + if (v == 'tcp'): + outname = OutDir + '/' + v + '_geoval_' + \ + self.validtime.strftime("%Y%m%d%H") + '.nc4' if (p == 'windprof' or p == 'satwind' or p == 'scatwind' or p == 'vadwind' or p == 'pibal'): outname = OutDir + '/' + p + '_geoval_' + \ self.validtime.strftime("%Y%m%d%H") + '.nc4' @@ -905,6 +914,9 @@ def toIODAobs(self, OutDir, TotalBias=False, clobber=True, platforms=None): if (v == 'sst'): outname = OutDir + '/' + v + '_obs_' + \ self.validtime.strftime("%Y%m%d%H") + '.nc4' + if (v == 'tcp'): + outname = OutDir + '/' + v + '_obs_' + \ + self.validtime.strftime("%Y%m%d%H") + '.nc4' if (p == 'windprof' or p == 'satwind' or p == 'scatwind' or p == 'vadwind' or p == 'pibal'): outname = OutDir + '/' + p + '_obs_' + \ self.validtime.strftime("%Y%m%d%H") + '.nc4' diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 1ec21ab5a..5b1e571d9 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -112,6 +112,7 @@ list( APPEND test_input testinput/balloonsonde_sample_20241001T00Z_PT15M.json testinput/tec_groundBased_TENET_20201001T00Z.tec testinput/COWVR_TSDR.020421.20240430T215845.20240430T230345.V1001.sample.h5 + testinput/gfs.231010.t00z.syndata.tcvitals.txt ) list( APPEND test_output @@ -215,6 +216,7 @@ list( APPEND test_output testoutput/obs.20241001T000000Z_PT15M_balloonsonde.nc4 testoutput/obs.20201001T000000Z_PT1M_tec_groundbased.nc4 testoutput/obs.20240430T220600Z_PT2M_cowvr_iss.nc4 + testoutput/tcp_obs_2023101000.nc4 ) if( iodaconv_gnssro_ENABLED ) @@ -873,15 +875,15 @@ ecbuild_add_test( TARGET test_${PROJECT_NAME}_pace_radiance_L1B -d 2019032112" pace_radiance_L1B.nc ${IODA_CONV_COMP_TOL}) -ecbuild_add_test( TARGET test_${PROJECT_NAME}_modis_water_leaving_radiance_L2 - TYPE SCRIPT - ENVIRONMENT "PYTHONPATH=${IODACONV_PYTHONPATH}" - COMMAND bash - ARGS ${CMAKE_BINARY_DIR}/bin/iodaconv_comp.sh - netcdf - "${Python3_EXECUTABLE} ${CMAKE_BINARY_DIR}/bin/oc_l2_radiance2ioda.py - -i testinput/modis_aqua_Rrs_OC_L2.nc - -o testrun/modis_water_leaving_radiance_L2.nc +ecbuild_add_test( TARGET test_${PROJECT_NAME}_modis_water_leaving_radiance_L2 + TYPE SCRIPT + ENVIRONMENT "PYTHONPATH=${IODACONV_PYTHONPATH}" + COMMAND bash + ARGS ${CMAKE_BINARY_DIR}/bin/iodaconv_comp.sh + netcdf + "${Python3_EXECUTABLE} ${CMAKE_BINARY_DIR}/bin/oc_l2_radiance2ioda.py + -i testinput/modis_aqua_Rrs_OC_L2.nc + -o testrun/modis_water_leaving_radiance_L2.nc -d 2021080112" modis_water_leaving_radiance_L2.nc ${IODA_CONV_COMP_TOL}) @@ -979,6 +981,21 @@ ecbuild_add_test( TARGET test_${PROJECT_NAME}_sondetac --date 2021-08-16T12:00:00Z --netcdf" 2021081612_sonde_small.nc ${IODA_CONV_COMP_TOL}) +#=============================================================================== +# Conventional data, tropical storm central pressure +#=============================================================================== + +ecbuild_add_test( TARGET test_${PROJECT_NAME}_tcvitals + TYPE SCRIPT + ENVIRONMENT "PYTHONPATH=${IODACONV_PYTHONPATH}" + COMMAND bash + ARGS ${CMAKE_BINARY_DIR}/bin/iodaconv_comp.sh + netcdf + "${Python3_EXECUTABLE} ${CMAKE_BINARY_DIR}/bin/tcv_txt2ioda.py + -i testinput/gfs.231010.t00z.syndata.tcvitals.txt + -o testrun/tcp_obs_2023101000.nc4" + tcp_obs_2023101000.nc4 ${IODA_CONV_COMP_TOL}) + #=============================================================================== # MRMS (multi-radar, multi-sensor) radar reflectivity CONUS converter #=============================================================================== @@ -1021,7 +1038,7 @@ ecbuild_add_test( TARGET test_iodaconv_tropics -o testrun/20220216T00Z_PT3H_tms_tropics-01.nc4 -d 2022021600" 20220216T00Z_PT3H_tms_tropics-01.nc4 ${IODA_CONV_COMP_TOL}) - + ecbuild_add_test( TARGET test_iodaconv_amsr2_gcom-w1 TYPE SCRIPT ENVIRONMENT "PYTHONPATH=${IODACONV_PYTHONPATH}" @@ -1033,7 +1050,7 @@ ecbuild_add_test( TARGET test_iodaconv_amsr2_gcom-w1 -o testrun/20220216T00Z_PT1H_amsr2_sample.nc4 -d 2022021600" 20220216T00Z_PT1H_amsr2_sample.nc4 ${IODA_CONV_COMP_TOL}) - + ecbuild_add_test( TARGET test_iodaconv_gmi_gpm TYPE SCRIPT ENVIRONMENT "PYTHONPATH=${IODACONV_PYTHONPATH}" @@ -1479,7 +1496,7 @@ ecbuild_add_test( TARGET test_${PROJECT_NAME}_modis_aod -i testinput/modis_aod.hdf --platform Terra --date_range 2021080100 2021080101 - -o testrun/modis_aod.nc" + -o testrun/modis_aod.nc" modis_aod.nc ${IODA_CONV_COMP_TOL}) ecbuild_add_test( TARGET test_${PROJECT_NAME}_aeronet_aod @@ -2307,8 +2324,8 @@ endif() COMMAND bash ARGS ${CMAKE_BINARY_DIR}/bin/iodaconv_comp.sh netcdf - "${Python3_EXECUTABLE} ${CMAKE_BINARY_DIR}/bin/cowvr_hdf5_2ioda.py - -i testinput/COWVR_TSDR.020421.20240430T215845.20240430T230345.V1001.sample.h5 + "${Python3_EXECUTABLE} ${CMAKE_BINARY_DIR}/bin/cowvr_hdf5_2ioda.py + -i testinput/COWVR_TSDR.020421.20240430T215845.20240430T230345.V1001.sample.h5 -o testrun/obs.20240430T220600Z_PT2M_cowvr_iss.nc4 -d 2024050100" obs.20240430T220600Z_PT2M_cowvr_iss.nc4 ${IODA_CONV_COMP_TOL}) diff --git a/test/testinput/gfs.231010.t00z.syndata.tcvitals.txt b/test/testinput/gfs.231010.t00z.syndata.tcvitals.txt new file mode 100644 index 000000000..c13a3f366 --- /dev/null +++ b/test/testinput/gfs.231010.t00z.syndata.tcvitals.txt @@ -0,0 +1,4 @@ +NHC 15E LIDIA 20231010 0000 182N 1106W 065 041 0983 1005 0408 33 037 0130 0130 0093 0093 D +NHC 16E MAX 20231010 0000 179N 1011W 020 026 1002 1005 0111 18 037 0056 0056 0056 -999 S +JTWC 15W BOLAVEN 20231010 0000 137N 1465E 310 067 0990 1006 0509 36 027 0213 0222 0158 0139 D +JTWC 14W KOINU 20231010 0000 208N 1103E 250 047 1011 1013 0111 13 018 -999 -999 -999 -999 S diff --git a/test/testoutput/tcp_obs_2023101000.nc4 b/test/testoutput/tcp_obs_2023101000.nc4 new file mode 100644 index 000000000..985880c3f --- /dev/null +++ b/test/testoutput/tcp_obs_2023101000.nc4 @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:31b5bbc90bcb608dcf10483e063982ef6fec1d7455b228087f0f4111e2ddcbfd +size 20131