From 9a3c226237c79b6a74312c8b1c3b99cf33ac4b28 Mon Sep 17 00:00:00 2001 From: emilyhcliu <36091766+emilyhcliu@users.noreply.github.com> Date: Tue, 31 Oct 2023 16:12:13 -0400 Subject: [PATCH 1/9] Update scatwind yamls for testing and end-to-end (#697) * Update scatwind YAML files based on UFO update accordingly. * Remove extra comments * Add YAML files for scatwind under config * Move filters two spacings left to line up with the filter header. --- parm/atm/obs/config/scatwind_metop-a.yaml | 280 +++++++++++++ parm/atm/obs/config/scatwind_metop-b.yaml | 280 +++++++++++++ parm/atm/obs/testing/scatwind.yaml | 481 ++++++++++++---------- parm/atm/obs/testing/scatwind_noqc.yaml | 39 +- 4 files changed, 862 insertions(+), 218 deletions(-) create mode 100644 parm/atm/obs/config/scatwind_metop-a.yaml create mode 100644 parm/atm/obs/config/scatwind_metop-b.yaml diff --git a/parm/atm/obs/config/scatwind_metop-a.yaml b/parm/atm/obs/config/scatwind_metop-a.yaml new file mode 100644 index 000000000..4735c31fe --- /dev/null +++ b/parm/atm/obs/config/scatwind_metop-a.yaml @@ -0,0 +1,280 @@ +obs space: + name: ascatw_ascat_metop-a + obsdatain: + engine: + type: H5File + obsfile: $(DATA)/obs/$(OPREFIX)ascatw.ascat_metop-a.{{ current_cycle | to_YMDH }}.nc4 + obsdataout: + engine: + type: H5File + obsfile: $(DATA)/diags/diag_ascatw_ascat_metop-a_{{ current_cycle | to_YMDH }}.nc4 + io pool: + max pool size: 1 + simulated variables: [windEastward, windNorthward] + +obs operator: + name: VertInterp + # Use height vertical coordinate first +# vertical coordinate: geometric_height + vertical coordinate: geopotential_height + observation vertical coordinate group: DerivedVariables + observation vertical coordinate: adjustedHeight + interpolation method: linear + hofx scaling field: SurfaceWindScalingHeight + hofx scaling field group: DerivedVariables + +obs prior filters: +# Apply variable changes needed for rescaled height coordinate +- filter: Variable Transforms + Transform: AdjustedHeightCoordinate + SkipWhenNoObs: False + +# Apply variable changes needed for wind scaling +- filter: Variable Transforms + Transform: SurfaceWindScalingHeight + SkipWhenNoObs: False + +# Assign the initial observation error (constant value, 1.5 m/s right now). +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + action: + name: assign error + error parameter: 1.5 + +# Calculate error inflation factor for duplicate observations +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windEastward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windEastward + +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windNorthward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windNorthward + +# Reject all obs with PreQC mark already set above 3 +# NOTE: All scatwinds have an automatic PreQC mark of 2 (hard-wired default from GSI) +# - filter: PreQC +# maxvalue: 3 +# action: +# name: reject + +obs post filters: +# Reject all ASCAT (Type 290) winds with tsavg <= 273.0 (surface temperature) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 +# - variable: GeoVaLs/surface_temperature + - variable: GeoVaLs/surface_temperature_where_land + maxvalue: 273. + action: + name: reject + +# Reject all ASCAT (Type 290) winds with isflg /= 0 (non-water surface) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + - variable: GeoVaLs/water_area_fraction + maxvalue: 0.99 + action: + name: reject + +# Reject ASCAT (Type 290) when observed component deviates from background by more than 5.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windEastward + - name: HofX/windEastward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 + +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windNorthward + - name: HofX/windNorthward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 + +# Reject OSCAT (Type 291) when observed component deviates from background by more than 6.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Background Check + filter variables: + - name: windEastward + - name: windNorthward + threshold: 6. + absolute threshold: 6. + where: + - variable: ObsType/windEastward + is_in: 291 + action: + name: reject + +# Reject ASCAT (Type 290) when ambiguity check fails (returned value is negative) +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + test variables: + - name: ObsFunction/ScatWindsAmbiguityCheck + options: + minimum_uv: 0.0001 # hard-coding a minimum-uv for transparancy, want this to basically be zero + maxvalue: 0. + action: + name: reject + +# All scatwinds must adjust errors based on ObsErrorFactorPressureCheck. +# This check will inflate errors for obs that are too close to either +# the model top or bottom. +- filter: Perform Action + filter variables: + - name: windEastward + where: + - variable: + name: ObsType/windEastward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck + options: + variable: windEastward + inflation factor: 4.0 + +- filter: Perform Action + filter variables: + - name: windNorthward + where: + - variable: + name: ObsType/windNorthward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck + options: + variable: windNorthward + inflation factor: 4.0 + +# All scatwinds subject to a gross error check. This is contained within +# the WindsSPDBCheck, although it is not exclusive to satwinds. +- filter: Background Check + filter variables: + - name: windEastward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windEastward + action: + name: reject + +- filter: Background Check + filter variables: + - name: windNorthward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windNorthward + action: + name: reject + +# The last error inflation check is for duplicate observations. This one needs +# to come last, because we don't want to inflate errors for duplication if one +# of the duplicates should be rejected. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windEastward + +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windNorthward + +# There is no across-the-board inflation for nvqc=.true. for scatwinds, presumably because for +# this inflation to take place both nvqc must be .true. AND ibeta must be >0, see: +# https://github.com/NOAA-EMC/GSI/blob/14ae595af1b03471287d322596d35c0665336e95/src/gsi/setupw.f90#L1229 +# GSI settings must have ibeta>0 for satwinds, but not for scatwinds. +# +# If the ibeta settings for scatwinds were to change while nvqc remained .true., we would extend YAML to +# an additional filter that inflates final ob-errors across-the-board by 1/0.8 = 1.25. NOTE: the nvqc setting +# is defaulted to .false. in GSI code, but is overridden in global operational configuration. See: +# configuration, see: https://github.com/NOAA-EMC/global-workflow/blob/d5ae3328fa4041b177357b1133f6b92e81c859d7/scripts/exglobal_atmos_analysis.sh#L750 +# This setting activates Line 1229 of setupw.f90 to scale ratio_errors by 0.8, which is applied in +# the denominator of the final ob-error, so 1/0.8 = 1.25 factor of ob-error inflation. +# +# If this functionality were to be activated for scatwinds, you would want to include this last inflation filter. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# where: +# - variable: ObsType/windEastward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# where: +# - variable: ObsType/windNorthward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 + +# END OF FILTERS# diff --git a/parm/atm/obs/config/scatwind_metop-b.yaml b/parm/atm/obs/config/scatwind_metop-b.yaml new file mode 100644 index 000000000..5ed7c38bc --- /dev/null +++ b/parm/atm/obs/config/scatwind_metop-b.yaml @@ -0,0 +1,280 @@ +obs space: + name: ascatw_ascat_metop-b + obsdatain: + engine: + type: H5File + obsfile: $(DATA)/obs/$(OPREFIX)ascatw.ascat_metop-b.{{ current_cycle | to_YMDH }}.nc4 + obsdataout: + engine: + type: H5File + obsfile: $(DATA)/diags/diag_ascatw_ascat_metop-b_{{ current_cycle | to_YMDH }}.nc4 + io pool: + max pool size: 1 + simulated variables: [windEastward, windNorthward] + +obs operator: + name: VertInterp + # Use height vertical coordinate first +# vertical coordinate: geometric_height + vertical coordinate: geopotential_height + observation vertical coordinate group: DerivedVariables + observation vertical coordinate: adjustedHeight + interpolation method: linear + hofx scaling field: SurfaceWindScalingHeight + hofx scaling field group: DerivedVariables + +obs prior filters: +# Apply variable changes needed for rescaled height coordinate +- filter: Variable Transforms + Transform: AdjustedHeightCoordinate + SkipWhenNoObs: False + +# Apply variable changes needed for wind scaling +- filter: Variable Transforms + Transform: SurfaceWindScalingHeight + SkipWhenNoObs: False + +# Assign the initial observation error (constant value, 1.5 m/s right now). +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + action: + name: assign error + error parameter: 1.5 + +# Calculate error inflation factor for duplicate observations +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windEastward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windEastward + +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windNorthward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windNorthward + +# Reject all obs with PreQC mark already set above 3 +# NOTE: All scatwinds have an automatic PreQC mark of 2 (hard-wired default from GSI) +# - filter: PreQC +# maxvalue: 3 +# action: +# name: reject + +obs post filters: +# Reject all ASCAT (Type 290) winds with tsavg <= 273.0 (surface temperature) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 +# - variable: GeoVaLs/surface_temperature + - variable: GeoVaLs/surface_temperature_where_land + maxvalue: 273. + action: + name: reject + +# Reject all ASCAT (Type 290) winds with isflg /= 0 (non-water surface) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + - variable: GeoVaLs/water_area_fraction + maxvalue: 0.99 + action: + name: reject + +# Reject ASCAT (Type 290) when observed component deviates from background by more than 5.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windEastward + - name: HofX/windEastward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 + +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windNorthward + - name: HofX/windNorthward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 + +# Reject OSCAT (Type 291) when observed component deviates from background by more than 6.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Background Check + filter variables: + - name: windEastward + - name: windNorthward + threshold: 6. + absolute threshold: 6. + where: + - variable: ObsType/windEastward + is_in: 291 + action: + name: reject + +# Reject ASCAT (Type 290) when ambiguity check fails (returned value is negative) +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + test variables: + - name: ObsFunction/ScatWindsAmbiguityCheck + options: + minimum_uv: 0.0001 # hard-coding a minimum-uv for transparancy, want this to basically be zero + maxvalue: 0. + action: + name: reject + +# All scatwinds must adjust errors based on ObsErrorFactorPressureCheck. +# This check will inflate errors for obs that are too close to either +# the model top or bottom. +- filter: Perform Action + filter variables: + - name: windEastward + where: + - variable: + name: ObsType/windEastward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck + options: + variable: windEastward + inflation factor: 4.0 + +- filter: Perform Action + filter variables: + - name: windNorthward + where: + - variable: + name: ObsType/windNorthward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck + options: + variable: windNorthward + inflation factor: 4.0 + +# All scatwinds subject to a gross error check. This is contained within +# the WindsSPDBCheck, although it is not exclusive to satwinds. +- filter: Background Check + filter variables: + - name: windEastward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windEastward + action: + name: reject + +- filter: Background Check + filter variables: + - name: windNorthward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windNorthward + action: + name: reject + +# The last error inflation check is for duplicate observations. This one needs +# to come last, because we don't want to inflate errors for duplication if one +# of the duplicates should be rejected. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windEastward + +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windNorthward + +# There is no across-the-board inflation for nvqc=.true. for scatwinds, presumably because for +# this inflation to take place both nvqc must be .true. AND ibeta must be >0, see: +# https://github.com/NOAA-EMC/GSI/blob/14ae595af1b03471287d322596d35c0665336e95/src/gsi/setupw.f90#L1229 +# GSI settings must have ibeta>0 for satwinds, but not for scatwinds. +# +# If the ibeta settings for scatwinds were to change while nvqc remained .true., we would extend YAML to +# an additional filter that inflates final ob-errors across-the-board by 1/0.8 = 1.25. NOTE: the nvqc setting +# is defaulted to .false. in GSI code, but is overridden in global operational configuration. See: +# configuration, see: https://github.com/NOAA-EMC/global-workflow/blob/d5ae3328fa4041b177357b1133f6b92e81c859d7/scripts/exglobal_atmos_analysis.sh#L750 +# This setting activates Line 1229 of setupw.f90 to scale ratio_errors by 0.8, which is applied in +# the denominator of the final ob-error, so 1/0.8 = 1.25 factor of ob-error inflation. +# +# If this functionality were to be activated for scatwinds, you would want to include this last inflation filter. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# where: +# - variable: ObsType/windEastward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# where: +# - variable: ObsType/windNorthward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 + +# END OF FILTERS# diff --git a/parm/atm/obs/testing/scatwind.yaml b/parm/atm/obs/testing/scatwind.yaml index 30a1af6f7..9befdda81 100644 --- a/parm/atm/obs/testing/scatwind.yaml +++ b/parm/atm/obs/testing/scatwind.yaml @@ -9,228 +9,283 @@ obs space: type: H5File obsfile: !ENV scatwind_diag_${CDATE}.nc4 simulated variables: [windEastward, windNorthward] + geovals: filename: !ENV scatwind_geoval_${CDATE}.nc4 -vector ref: GsiHofXBc -tolerance: 0.01 + obs operator: name: VertInterp - apply near surface wind scaling: true -# -obs pre filters: - # Assign the initial observation error (constant value, 1.5 m/s right now). - - filter: Perform Action - filter variables: - - name: windEastward - - name: windNorthward - action: - name: assign error - error parameter: 1.5 -#obs prior filters: - # - # Reject all obs with PreQC mark already set above 3 - # NOTE: All scatwinds have an automatic PreQC mark of 2 (hard-wired default from GSI) - # - filter: PreQC - # maxvalue: 3 - # action: - # name: reject - # + # Use height vertical coordinate first +# vertical coordinate: geometric_height + vertical coordinate: geopotential_height + observation vertical coordinate group: DerivedVariables + observation vertical coordinate: adjustedHeight + interpolation method: linear + hofx scaling field: SurfaceWindScalingHeight + hofx scaling field group: DerivedVariables + +obs prior filters: +# Apply variable changes needed for rescaled height coordinate +- filter: Variable Transforms + Transform: AdjustedHeightCoordinate + SkipWhenNoObs: False + +# Apply variable changes needed for wind scaling +- filter: Variable Transforms + Transform: SurfaceWindScalingHeight + SkipWhenNoObs: False + +# Assign the initial observation error (constant value, 1.5 m/s right now). +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + action: + name: assign error + error parameter: 1.5 + +# Calculate error inflation factor for duplicate observations +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windEastward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windEastward + +#- filter: Variable Assignment +# assignments: +# - name: ObsErrorFactorDuplicateCheck/windNorthward +# type: float +# function: +# name: ObsFunction/ObsErrorFactorDuplicateCheck +# options: +# use_air_pressure: true +# variable: windNorthward + +# Reject all obs with PreQC mark already set above 3 +# NOTE: All scatwinds have an automatic PreQC mark of 2 (hard-wired default from GSI) +# - filter: PreQC +# maxvalue: 3 +# action: +# name: reject + obs post filters: - # - # Reject all ASCAT (Type 290) winds with tsavg <= 273.0 (surface temperature) - # - - filter: Perform Action - filter variables: - - name: windEastward - - name: windNorthward - where: - - variable: ObsType/windEastward - is_in: 290 - - variable: GeoVaLs/surface_temperature - maxvalue: 273. - action: - name: reject - # - # Reject all ASCAT (Type 290) winds with isflg /= 0 (non-water surface) - # - - filter: Perform Action - filter variables: - - name: windEastward - - name: windNorthward - where: - - variable: ObsType/windEastward - is_in: 290 - - variable: GeoVaLs/water_area_fraction - maxvalue: 0.001 - action: - name: reject - # Reject ASCAT (Type 290) when observed component deviates from background by more than 5.0 m/s - # NOTE: This check can reject a u- or v-component of the same observation independently, which - # is fundamentally different from how GSI rejects obs (both components are rejected if - # either component fails a check). - - filter: Background Check - filter variables: - - name: windEastward - - name: windNorthward - threshold: 5. - absolute threshold: 5. - where: - - variable: ObsType/windEastward - is_in: 290 - action: - name: reject - # Reject OSCAT (Type 291) when observed component deviates from background by more than 6.0 m/s - # NOTE: This check can reject a u- or v-component of the same observation independently, which - # is fundamentally different from how GSI rejects obs (both components are rejected if - # either component fails a check). - - filter: Background Check - filter variables: - - name: windEastward - - name: windNorthward - threshold: 6. - absolute threshold: 6. - where: - - variable: ObsType/windEastward - is_in: 291 - action: - name: reject - # Reject ASCAT (Type 290) when ambiguity check fails (returned value is negative) - - filter: Bounds Check - filter variables: - - name: windEastward - - name: windNorthward - where: - - variable: ObsType/windEastward - is_in: 290 - test variables: - - name: ObsFunction/ScatWindsAmbiguityCheck - options: - test_hofx: GsiHofX - minimum_uv: 0.0001 # hard-coding a minimum-uv for transparancy, want this to basically be zero - maxvalue: 0. - action: - name: reject -# - # All scatwinds must adjust errors based on ObsErrorFactorPressureCheck. - # This check will inflate errors for obs that are too close to either - # the model top or bottom. - - filter: Perform Action - filter variables: - - name: windEastward - where: - - variable: - name: ObsType/windEastward - is_in: 290-291 - action: - name: inflate error - inflation variable: - name: ObsFunction/ObsErrorFactorPressureCheck - options: - variable: windEastward - inflation factor: 4.0 -# - - filter: Perform Action - filter variables: - - name: windNorthward - where: - - variable: - name: ObsType/windNorthward - is_in: 290-291 - action: - name: inflate error - inflation variable: - name: ObsFunction/ObsErrorFactorPressureCheck - options: - variable: windNorthward - inflation factor: 4.0 - # All scatwinds subject to a gross error check. This is contained within - # the WindsSPDBCheck, although it is not exclusive to satwinds. - - filter: Background Check - filter variables: - - name: windEastward - function absolute threshold: - - name: ObsFunction/WindsSPDBCheck +# Reject all ASCAT (Type 290) winds with tsavg <= 273.0 (surface temperature) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 +# - variable: GeoVaLs/surface_temperature + - variable: GeoVaLs/surface_temperature_where_land + maxvalue: 273. + action: + name: reject + +# Reject all ASCAT (Type 290) winds with isflg /= 0 (non-water surface) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + - variable: GeoVaLs/water_area_fraction + maxvalue: 0.99 +# - variable: GeoVaLs/surface_type +# minvalue: 0.001 + action: + name: reject +#passedBenchmark: 53558 # 2 variables (u,v), u=26779 passing, v=26779 passing + +# Reject ASCAT (Type 290) when observed component deviates from background by more than 5.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windEastward + - name: HofX/windEastward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 + +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + test variables: + - name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/windNorthward + - name: HofX/windNorthward + coefs: [1.0, -1.0] + minvalue: -5.0 + maxvalue: 5.0 +#passedBenchmark: 52362 # 2 variables (u,v), u=26181 passing, v=26181 passing + +# Reject OSCAT (Type 291) when observed component deviates from background by more than 6.0 m/s +# NOTE: This check can reject a u- or v-component of the same observation independently, which +# is fundamentally different from how GSI rejects obs (both components are rejected if +# either component fails a check). +- filter: Background Check + filter variables: + - name: windEastward + - name: windNorthward + threshold: 6. + absolute threshold: 6. + where: + - variable: ObsType/windEastward + is_in: 291 + action: + name: reject + +# Reject ASCAT (Type 290) when ambiguity check fails (returned value is negative) +- filter: Bounds Check + filter variables: + - name: windEastward + - name: windNorthward + where: + - variable: ObsType/windEastward + is_in: 290 + test variables: + - name: ObsFunction/ScatWindsAmbiguityCheck + options: +# test_hofx: GsiHofX + minimum_uv: 0.0001 # hard-coding a minimum-uv for transparancy, want this to basically be zero + maxvalue: 0. + action: + name: reject +#passedBenchmark: 51776 # 2 variables (u,v), u=25888 passing, v=25888 passing + +# All scatwinds must adjust errors based on ObsErrorFactorPressureCheck. +# This check will inflate errors for obs that are too close to either +# the model top or bottom. +- filter: Perform Action + filter variables: + - name: windEastward + where: + - variable: + name: ObsType/windEastward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck options: - wndtype: [ 290, 291] - cgross: [ 5.0, 5.0] - error_min: [1.4, 1.4] - error_max: [6.1, 6.1] variable: windEastward - action: - name: reject - - - filter: Background Check - filter variables: - - name: windNorthward - function absolute threshold: - - name: ObsFunction/WindsSPDBCheck + inflation factor: 4.0 + +- filter: Perform Action + filter variables: + - name: windNorthward + where: + - variable: + name: ObsType/windNorthward + is_in: 290-291 + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorPressureCheck options: - wndtype: [ 290, 291] - cgross: [ 5.0, 5.0] - error_min: [1.4, 1.4] - error_max: [6.1, 6.1] variable: windNorthward - action: - name: reject - # The last error inflation check is for duplicate observations. This one needs - # to come last, because we don't want to inflate errors for duplication if one - # of the duplicates should be rejected. - - filter: RejectList - filter variables: - - name: windEastward - action: - name: inflate error - inflation variable: - name: ObsFunction/ObsErrorFactorDuplicateCheck - options: - use_air_pressure: true - variable: windEastward - - filter: RejectList - filter variables: - - name: windNorthward - action: - name: inflate error - inflation variable: - name: ObsFunction/ObsErrorFactorDuplicateCheck - options: - use_air_pressure: true - variable: windNorthward - # There is no across-the-board inflation for nvqc=.true. for scatwinds, presumably because for - # this inflation to take place both nvqc must be .true. AND ibeta must be >0, see: - # https://github.com/NOAA-EMC/GSI/blob/14ae595af1b03471287d322596d35c0665336e95/src/gsi/setupw.f90#L1229 - # GSI settings must have ibeta>0 for satwinds, but not for scatwinds. - # - # - # If the ibeta settings for scatwinds were to change while nvqc remained .true., we would extend YAML to - # an additional filter that inflates final ob-errors across-the-board by 1/0.8 = 1.25. NOTE: the nvqc setting - # is defaulted to .false. in GSI code, but is overridden in global operational configuration. See: - # configuration, see: https://github.com/NOAA-EMC/global-workflow/blob/d5ae3328fa4041b177357b1133f6b92e81c859d7/scripts/exglobal_atmos_analysis.sh#L750 - # This setting activates Line 1229 of setupw.f90 to scale ratio_errors by 0.8, which is applied in - # the denominator of the final ob-error, so 1/0.8 = 1.25 factor of ob-error inflation. - # - # If this functionality were to be activated for scatwinds, you would want to include this last inflation filter. -# - filter: Perform Action -# filter variables: -# - name: windEastward -# where: -# - variable: ObsType/windEastward -# is_in: 290-291 -# action: -# name: inflate error -# inflation factor: 1.25 -# - filter: Perform Action -# filter variables: -# - name: windNorthward -# where: -# - variable: ObsType/windNorthward -# is_in: 290-291 -# action: -# name: inflate error -# inflation factor: 1.25 -linear obs operator: - name: Identity -#passedBenchmark: 52097 # 2 variables (u,v), u=26087 passing, v=26010 passing -passedBenchmark: 52109 # 2 variables (u,v), u=26093 passing, v=26016 passing + inflation factor: 4.0 +#passedBenchmark: 51776 # 2 variables (u,v), u=25888 passing, v=25888 passing + +# All scatwinds subject to a gross error check. This is contained within +# the WindsSPDBCheck, although it is not exclusive to satwinds. +- filter: Background Check + filter variables: + - name: windEastward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windEastward + action: + name: reject + +- filter: Background Check + filter variables: + - name: windNorthward + function absolute threshold: + - name: ObsFunction/WindsSPDBCheck + options: + wndtype: [ 290, 291] + cgross: [ 5.0, 5.0] + error_min: [1.4, 1.4] + error_max: [6.1, 6.1] + variable: windNorthward + action: + name: reject +#passedBenchmark: 51776 # 2 variables (u,v), u=25888 passing, v=25888 passing + +# The last error inflation check is for duplicate observations. This one needs +# to come last, because we don't want to inflate errors for duplication if one +# of the duplicates should be rejected. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windEastward + +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# action: +# name: inflate error +# inflation variable: +# name: ObsErrorFactorDuplicateCheck/windNorthward + +# There is no across-the-board inflation for nvqc=.true. for scatwinds, presumably because for +# this inflation to take place both nvqc must be .true. AND ibeta must be >0, see: +# https://github.com/NOAA-EMC/GSI/blob/14ae595af1b03471287d322596d35c0665336e95/src/gsi/setupw.f90#L1229 +# GSI settings must have ibeta>0 for satwinds, but not for scatwinds. +# +# If the ibeta settings for scatwinds were to change while nvqc remained .true., we would extend YAML to +# an additional filter that inflates final ob-errors across-the-board by 1/0.8 = 1.25. NOTE: the nvqc setting +# is defaulted to .false. in GSI code, but is overridden in global operational configuration. See: +# configuration, see: https://github.com/NOAA-EMC/global-workflow/blob/d5ae3328fa4041b177357b1133f6b92e81c859d7/scripts/exglobal_atmos_analysis.sh#L750 +# This setting activates Line 1229 of setupw.f90 to scale ratio_errors by 0.8, which is applied in +# the denominator of the final ob-error, so 1/0.8 = 1.25 factor of ob-error inflation. +# +# If this functionality were to be activated for scatwinds, you would want to include this last inflation filter. +#- filter: Perform Action +# filter variables: +# - name: windEastward +# where: +# - variable: ObsType/windEastward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 +#- filter: Perform Action +# filter variables: +# - name: windNorthward +# where: +# - variable: ObsType/windNorthward +# is_in: 290-291 +# action: +# name: inflate error +# inflation factor: 1.25 +passedBenchmark: 51764 # 2 variables (u,v), u=25882 passing, v=25882 passing # GSI rejects both u- and v-component in first-guess check, UFO does not, but when UFO rej is reconciled btwn u and v these match GSI rej # u: 207 obs pass with corresponding v being rejected, 25880 obs in UFO/GSI agreement # v: 130 obs pass with corresponding u being rejected, 25880 obs in UFO/GSI agreement diff --git a/parm/atm/obs/testing/scatwind_noqc.yaml b/parm/atm/obs/testing/scatwind_noqc.yaml index f8f36efeb..ec4aef452 100644 --- a/parm/atm/obs/testing/scatwind_noqc.yaml +++ b/parm/atm/obs/testing/scatwind_noqc.yaml @@ -9,13 +9,42 @@ obs space: type: H5File obsfile: !ENV scatwind_diag_${CDATE}.nc4 simulated variables: [windEastward, windNorthward] + geovals: filename: !ENV scatwind_geoval_${CDATE}.nc4 -vector ref: GsiHofXBc -tolerance: 0.01 + obs operator: name: VertInterp - apply near surface wind scaling: true + # Use height vertical coordinate first +# vertical coordinate: geometric_height + vertical coordinate: geopotential_height + observation vertical coordinate group: DerivedVariables + observation vertical coordinate: adjustedHeight + interpolation method: linear + hofx scaling field: SurfaceWindScalingHeight + hofx scaling field group: DerivedVariables + +obs prior filters: +# Apply variable changes needed for rescaled height coordinate +- filter: Variable Transforms + Transform: AdjustedHeightCoordinate + SkipWhenNoObs: False + +# Apply variable changes needed for wind scaling +- filter: Variable Transforms + Transform: SurfaceWindScalingHeight + SkipWhenNoObs: False + +# Assign the initial observation error (constant value, 1.5 m/s) +- filter: Perform Action + filter variables: + - name: windEastward + - name: windNorthward + action: + name: assign error + error parameter: 1.5 -vector ref: GsiHofXBc -tolerance: 1.e-5 +passedBenchmark: 53546 +#passedBenchmark: 53558 +#vector ref: GsiHofXBc +#tolerance: 1.e-5 From e609c1102db6450c6db03d5de0bc97c79d6530aa Mon Sep 17 00:00:00 2001 From: Guillaume Vernieres Date: Wed, 1 Nov 2023 08:02:46 -0400 Subject: [PATCH 2/9] Add functionality to locally change the weights used for the hybrid B (#700) * patch ens. B * saving ens. mean * adjusted weights locations * tidy * uups --- parm/soca/berror/soca_ensb.yaml | 6 +++ parm/soca/berror/soca_ensweights.yaml | 10 +++++ test/soca/testinput/socahybridweights.yaml | 15 ++++++- utils/soca/gdas_ens_handler.h | 11 +++-- utils/soca/gdas_socahybridweights.h | 51 +++++++++++++++++++++- 5 files changed, 86 insertions(+), 7 deletions(-) diff --git a/parm/soca/berror/soca_ensb.yaml b/parm/soca/berror/soca_ensb.yaml index 814dfc7ac..447d5b5fd 100644 --- a/parm/soca/berror/soca_ensb.yaml +++ b/parm/soca/berror/soca_ensb.yaml @@ -30,6 +30,12 @@ steric height: linear variable changes: - linear variable change name: BalanceSOCA # Only the steric balance is applied +ensemble mean output: + datadir: ./static_ens + date: '{{ATM_WINDOW_BEGIN}}' + exp: ens_mean + type: incr + ssh output: unbalanced: datadir: ./static_ens diff --git a/parm/soca/berror/soca_ensweights.yaml b/parm/soca/berror/soca_ensweights.yaml index c7ea40a98..094377f95 100644 --- a/parm/soca/berror/soca_ensweights.yaml +++ b/parm/soca/berror/soca_ensweights.yaml @@ -19,6 +19,16 @@ weights: # Need to provide weights^2 when reading from file ice: 0.0025 # 5% of original variance ocean: 0.0625 # 25% " " + # Apply localized weights to the ocean ens. B + ocean local weights: + - lon: -172.0 + lat: 11.0 + amplitude: -1.0 + length scale: 700.0 + - lon: -160.0 + lat: 12.0 + amplitude: -1.0 + length scale: 700.0 output: datadir: ./ diff --git a/test/soca/testinput/socahybridweights.yaml b/test/soca/testinput/socahybridweights.yaml index cdad1ab27..855dd6d60 100644 --- a/test/soca/testinput/socahybridweights.yaml +++ b/test/soca/testinput/socahybridweights.yaml @@ -16,8 +16,19 @@ background: read_from_file: 1 weights: - ice: 0.1 - ocean: 0.5 + # Need to provide weights^2 when reading from file + ice: 0.0025 # 5% of original variance + ocean: 0.0625 # 25% " " + # Apply localized weights to the ocean ens. B + ocean local weights: + - lon: -172.0 + lat: 11.0 + amplitude: -1.0 + length scale: 700.0 + - lon: -160.0 + lat: 12.0 + amplitude: -1.0 + length scale: 700.0 output: datadir: ./ diff --git a/utils/soca/gdas_ens_handler.h b/utils/soca/gdas_ens_handler.h index 44e994a0b..5cf862651 100644 --- a/utils/soca/gdas_ens_handler.h +++ b/utils/soca/gdas_ens_handler.h @@ -104,6 +104,10 @@ namespace gdasapp { gdasapp_ens_utils::ensMoments(ensMembers, ensMean, ensStd, ensVariance); oops::Log::info() << "mean: " << ensMean << std::endl; oops::Log::info() << "std: " << ensStd << std::endl; + if ( fullConfig.has("ensemble mean output") ) { + const eckit::LocalConfiguration ensMeanOutputConfig(fullConfig, "ensemble mean output"); + ensMean.write(ensMeanOutputConfig); + } // Remove mean from ensemble members for (size_t i = 0; i < postProcIncr.ensSize_; ++i) { @@ -115,6 +119,7 @@ namespace gdasapp { eckit::LocalConfiguration stericVarChangeConfig; fullConfig.get("steric height", stericVarChangeConfig); oops::Log::info() << "steric config 0000: " << stericVarChangeConfig << std::endl; + // Initialize trajectories const eckit::LocalConfiguration trajConfig(fullConfig, "trajectory"); soca::State cycleTraj(geom, trajConfig); // trajectory of the cycle @@ -148,7 +153,7 @@ namespace gdasapp { soca::Increment incr = postProcIncr.appendLayer(ensMembers[i]); // Save total ssh - oops::Log::info() << "ssh ensemble memnber " << i << std::endl; + oops::Log::info() << "ssh ensemble member " << i << std::endl; soca::Increment ssh_tmp(geom, socaSshVar, postProcIncr.dt_); ssh_tmp = ensMembers[i]; sshTotal.push_back(ssh_tmp); @@ -186,14 +191,14 @@ namespace gdasapp { // Add the unbalanced ssh to the recentered perturbation // this assumes ssh_u is independent of the trajectory - oops::Log::info() << "&&&&& before adding ssh_u " << incr << std::endl; + oops::Log::debug() << "&&&&& before adding ssh_u " << incr << std::endl; atlas::FieldSet incrFs; incr.toFieldSet(incrFs); atlas::FieldSet sshNonStericFs; sshNonSteric[i].toFieldSet(sshNonStericFs); util::addFieldSets(incrFs, sshNonStericFs); incr.fromFieldSet(incrFs); - oops::Log::info() << "&&&&& after adding ssh_u " << incr << std::endl; + oops::Log::debug() << "&&&&& after adding ssh_u " << incr << std::endl; // Save final perturbation, used in the offline EnVAR result = postProcIncr.save(incr, i+1); diff --git a/utils/soca/gdas_socahybridweights.h b/utils/soca/gdas_socahybridweights.h index adf20188a..219450a8e 100644 --- a/utils/soca/gdas_socahybridweights.h +++ b/utils/soca/gdas_socahybridweights.h @@ -3,12 +3,17 @@ #include #include #include +#include #include "eckit/config/LocalConfiguration.h" #include "atlas/field.h" +#include "atlas/util/Earth.h" +#include "atlas/util/Geometry.h" +#include "atlas/util/Point.h" #include "oops/base/PostProcessor.h" +#include "oops/generic/gc99.h" #include "oops/mpi/mpi.h" #include "oops/runs/Application.h" #include "oops/util/DateTime.h" @@ -20,6 +25,38 @@ #include "soca/State/State.h" namespace gdasapp { + // Create a simple mask based on a Gaussian function + void gaussianMask(const soca::Geometry & geom, soca::Increment & gaussIncr, + eckit::LocalConfiguration conf) { + // Get the 2D grid + std::vector lats; + std::vector lons; + bool halo = true; + geom.latlon(lats, lons, halo); + + // Prepare fieldset from increment + atlas::FieldSet gaussIncrFs; + gaussIncr.toFieldSet(gaussIncrFs); + + // Get the GC99 parameters from config + double amp = conf.getDouble("amplitude"); + double scale = conf.getDouble("length scale"); + const atlas::PointLonLat p0(conf.getDouble("lon"), conf.getDouble("lat")); + + // Recompute weights + for (auto & field : gaussIncrFs) { + oops::Log::info() << "---------- Field name: " << field.name() << std::endl; + auto view = atlas::array::make_view(field); + for (int jnode = 0; jnode < field.shape(0); ++jnode) { + atlas::PointLonLat p1(lons[jnode], lats[jnode]); + double d = atlas::util::Earth::distance(p0, p1)/1000.0; + for (int jlevel = 0; jlevel < field.shape(1); ++jlevel) { + view(jnode, jlevel) += amp * oops::gc99(d/scale); + } + } + } + gaussIncr.fromFieldSet(gaussIncrFs); + } class SocaHybridWeights : public oops::Application { public: @@ -45,8 +82,7 @@ namespace gdasapp { socaVars += socaOcnVars; /// Read the background - // TODO(guillaume): Use the ice extent to set the weights ... no clue if this is - // possible at this level + // TODO(guillaume): Use the ice extent to set the weights soca::State socaBkg(geom, socaVars, dt); const eckit::LocalConfiguration socaBkgConfig(fullConfig, "background"); socaBkg.read(socaBkgConfig); @@ -70,6 +106,17 @@ namespace gdasapp { /// Create fields of weights for the ocean soca::Increment socaOcnHW(geom, socaOcnVars, dt); socaOcnHW.ones(); + + /// Apply localized gaussians to the weights + eckit::LocalConfiguration localWeightsConfigs(fullConfig, "weights.ocean local weights"); + std::vector localWeightsList = + localWeightsConfigs.getSubConfigurations(); + for (auto & conf : localWeightsList) { + gaussianMask(geom, socaOcnHW, conf); + oops::Log::info() << "Local weights for socaOcnHW: " << std::endl << conf << std::endl; + oops::Log::info() << socaOcnHW << std::endl; + } + socaOcnHW *= wOcean; oops::Log::info() << "socaOcnHW: " << std::endl << socaOcnHW << std::endl; socaOcnHW.write(socaHWOutConfig); From 00201eff3e7491f79ae764a36c9a368f5e1e9885 Mon Sep 17 00:00:00 2001 From: NicholasEsposito-NOAA <62616739+nicholasesposito@users.noreply.github.com> Date: Thu, 2 Nov 2023 12:45:05 -0400 Subject: [PATCH 3/9] Feature/nick e adpsfc prepbufr (#702) * working sfcshp * initial changes. no py add yet * add adpsfc py * some updates * update README * get rid of sys.appendpath * finally works again with json * coding norms hopefully pass now * more coding norms * more coding norms * one line wont work * logger * few logger updates * ok now should be 100 * change logger filename line * print -> loggers * Update bufr2ioda_adpsfc_prepbufr.py heightOfStation, stationElevation, var order, obstype * zobqm added * update zob elv * zob elv confusion * Oba -> Obs * readme change * rm readme, add to run_bufr2ioda.py * rm whitespace * rm more whitespace * rm run_bufr2ioda_adpsfc.sh - not needed * NESDIS removed. logger -> debug --------- Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito --- .../bufr2ioda/bufr2ioda_adpsfc_prepbufr.json | 12 + .../bufr2ioda/bufr2ioda_adpsfc_prepbufr.py | 347 ++++++++++++++++++ ush/ioda/bufr2ioda/gen_bufr2ioda_json.py | 8 +- ush/ioda/bufr2ioda/run_bufr2ioda.py | 2 +- 4 files changed, 367 insertions(+), 2 deletions(-) create mode 100644 parm/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.json create mode 100755 ush/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.py diff --git a/parm/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.json b/parm/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.json new file mode 100644 index 000000000..bdde19135 --- /dev/null +++ b/parm/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.json @@ -0,0 +1,12 @@ +{ + "data_format" : "prepbufr", + "subsets" : [ "ADPSFC" ], + "source" : "prepBUFR", + "data_type" : "ADPSFC", + "cycle_type" : "{{ RUN }}", + "cycle_datetime" : "{{ current_cycle | to_YMDH }}", + "dump_directory" : "{{ DMPDIR }}", + "ioda_directory" : "{{ COM_OBS }}", + "data_description" : "ADPSFC_prepbufr", + "data_provider" : "U.S. NOAA" +} diff --git a/ush/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.py b/ush/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.py new file mode 100755 index 000000000..a5700574a --- /dev/null +++ b/ush/ioda/bufr2ioda/bufr2ioda_adpsfc_prepbufr.py @@ -0,0 +1,347 @@ +#!/usr/bin/env python3 +# (C) Copyright 2023 NOAA/NWS/NCEP/EMC +# +# This software is licensed under the terms of the Apache Licence Version 2.0 +# which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. + +import sys +import numpy as np +import os +import argparse +import math +import calendar +import time +from datetime import datetime +import json +from pyiodaconv import bufr +from collections import namedtuple +from pyioda import ioda_obs_space as ioda_ospace +from wxflow import Logger + + +def Compute_dateTime(cycleTimeSinceEpoch, dhr): + + dhr = np.int64(dhr*3600) + dateTime = dhr + cycleTimeSinceEpoch + + return dateTime + + +def bufr_to_ioda(config, logger): + + subsets = config["subsets"] + logger.debug(f"Checking subsets = {subsets}") + + # Get parameters from configuration + data_format = config["data_format"] + source = config["source"] + data_type = config["data_type"] + data_description = config["data_description"] + data_provider = config["data_provider"] + cycle_type = config["cycle_type"] + cycle_datetime = config["cycle_datetime"] + dump_dir = config["dump_directory"] + ioda_dir = config["ioda_directory"] + cycle = config["cycle_datetime"] + + # Get derived parameters + yyyymmdd = cycle[0:8] + hh = cycle[8:10] + + reference_time = datetime.strptime(cycle, "%Y%m%d%H") + reference_time = reference_time.strftime("%Y-%m-%dT%H:%M:%SZ") + reference_time_full = f"{yyyymmdd}{hh}00" + + logger.debug(f"reference_time = {reference_time}") + + # General informaton + converter = 'BUFR to IODA Converter' + platform_description = 'SFCSHP data from prepBUFR format' + + bufrfile = f"{cycle_type}.t{hh}z.{data_format}" + DATA_PATH = os.path.join(dump_dir, f"{cycle_type}.{yyyymmdd}", + str(hh), bufrfile) + + logger.debug(f"The DATA_PATH is: {DATA_PATH}") + + # ============================================ + # Make the QuerySet for all the data we want + # ============================================ + start_time = time.time() + + logger.debug('Making QuerySet ...') + q = bufr.QuerySet(subsets) + + # ObsType + q.add('observationType', '*/TYP') + + # MetaData + q.add('stationIdentification', '*/SID') + q.add('latitude', '*/YOB') + q.add('longitude', '*/XOB') + q.add('obsTimeMinusCycleTime', '*/DHR') + q.add('heightOfStation', '*/Z___INFO/Z__EVENT{1}/ZOB') + q.add('pressure', '*/P___INFO/P__EVENT{1}/POB') + +# # Quality Infomation (Quality Indicator) + q.add('qualityMarkerStationPressure', '*/P___INFO/P__EVENT{1}/PQM') + q.add('qualityMarkerStationElevation', '*/Z___INFO/Z__EVENT{1}/ZQM') + + # ObsValue + q.add('stationPressure', '*/P___INFO/P__EVENT{1}/POB') + q.add('stationElevation', '*/ELV') + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for making QuerySet: {running_time} seconds") + + # ============================================================== + # Open the BUFR file and execute the QuerySet to get ResultSet + # Use the ResultSet returned to get numpy arrays of the data + # ============================================================== + start_time = time.time() + + logger.debug(f"Executing QuerySet to get ResultSet ...") + with bufr.File(DATA_PATH) as f: + r = f.execute(q) + + logger.debug(" ... Executing QuerySet: get ObsType ...") + # ObsType + typ = r.get('observationType') + + logger.debug(" ... Executing QuerySet: get MetaData ...") + # MetaData + sid = r.get('stationIdentification') + lat = r.get('latitude') + lon = r.get('longitude') + lon[lon > 180] -= 360 + zob = r.get('heightOfStation', type='float') + pressure = r.get('pressure') + pressure *= 100 + + logger.debug(f" ... Executing QuerySet: get QualityMarker information ...") + # Quality Information + pobqm = r.get('qualityMarkerStationPressure') + zobqm = r.get('qualityMarkerStationElevation') + + logger.debug(f" ... Executing QuerySet: get ObsValue ...") + # ObsValue + elv = r.get('stationElevation', type='float') + pob = r.get('stationPressure') + pob *= 100 + + logger.debug(f" ... Executing QuerySet: get dateTime ...") + # DateTime: seconds since Epoch time + # IODA has no support for numpy datetime arrays dtype=datetime64[s] + dhr = r.get('obsTimeMinusCycleTime', type='int64') + + logger.debug(f" ... Executing QuerySet: Done!") + + logger.debug(f" ... Executing QuerySet: Check BUFR variable generic \ + dimension and type ...") + # Check BUFR variable generic dimension and type + logger.debug(f" typ shape = {typ.shape}") + logger.debug(f" sid shape = {sid.shape}") + logger.debug(f" dhr shape = {dhr.shape}") + logger.debug(f" lat shape = {lat.shape}") + logger.debug(f" lon shape = {lon.shape}") + logger.debug(f" zob shape = {zob.shape}") + logger.debug(f" pressure shape = {pressure.shape}") + + logger.debug(f" pobqm shape = {pobqm.shape}") + logger.debug(f" zobqm shape = {zobqm.shape}") + + logger.debug(f" elv shape = {elv.shape}") + logger.debug(f" pob shape = {pob.shape}") + + logger.debug(f" dhr type = {dhr.shape}") + + logger.debug(f" sid type = {sid.dtype}") + logger.debug(f" dhr type = {dhr.dtype}") + logger.debug(f" lat type = {lat.dtype}") + logger.debug(f" lon type = {lon.dtype}") + logger.debug(f" zob type = {zob.dtype}") + logger.debug(f" typ type = {typ.dtype}") + logger.debug(f" pressure type = {pressure.dtype}") + + logger.debug(f" pobqm type = {pobqm.dtype}") + logger.debug(f" zobqm type = {zobqm.dtype}") + + logger.debug(f" elv type = {elv.dtype}") + logger.debug(f" pob type = {pob.dtype}") + + logger.debug(f" dhr type = {dhr.dtype}") + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for executing QuerySet to get ResultSet: \ + {running_time} seconds") + + # ========================= + # Create derived variables + # ========================= + start_time = time.time() + + logger.debug(f"Creating derived variables - dateTime ...") + + cycleTimeSinceEpoch = np.int64(calendar.timegm(time.strptime( + reference_time_full, '%Y%m%d%H%M'))) + dateTime = Compute_dateTime(cycleTimeSinceEpoch, dhr) + + logger.debug(f" Check derived variables type ... ") + logger.debug(f" dateTime shape = {dateTime.shape}") + logger.debug(f" dateTime type = {dateTime.dtype}") + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for creating derived variables: \ + {running_time} seconds") + + # ===================================== + # Create IODA ObsSpace + # Write IODA output + # ===================================== + + # Create the dimensions + dims = {'Location': np.arange(0, lat.shape[0])} + + iodafile = f"{cycle_type}.t{hh}z.{data_type}.{data_format}.nc" + OUTPUT_PATH = os.path.join(ioda_dir, iodafile) + logger.debug(f" ... ... Create OUTPUT file: {OUTPUT_PATH}") + + path, fname = os.path.split(OUTPUT_PATH) + if path and not os.path.exists(path): + os.makedirs(path) + + obsspace = ioda_ospace.ObsSpace(OUTPUT_PATH, mode='w', dim_dict=dims) + + # Create Global attributes + logger.debug(f" ... ... Create global attributes") + + obsspace.write_attr('Converter', converter) + obsspace.write_attr('source', source) + obsspace.write_attr('sourceFiles', bufrfile) + obsspace.write_attr('dataProviderOrigin', data_provider) + obsspace.write_attr('description', data_description) + obsspace.write_attr('datetimeReference', reference_time) + obsspace.write_attr('datetimeRange', + [str(min(dateTime)), str(max(dateTime))]) + obsspace.write_attr('platformLongDescription', platform_description) + + # Create IODA variables + logger.debug(f" ... ... Create variables: name, type, units, & attributes") + + # Observation Type - Station Elevation + obsspace.create_var('ObsType/stationElevation', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Station Elevation Observation Type') \ + .write_data(typ) + + # Observation Type - Station Pressure + obsspace.create_var('ObsType/stationPressure', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Station Pressure Observation Type') \ + .write_data(typ) + + # Longitude + obsspace.create_var('MetaData/longitude', dtype=lon.dtype, + fillval=lon.fill_value) \ + .write_attr('units', 'degrees_east') \ + .write_attr('valid_range', np.array([-180, 180], dtype=np.float32)) \ + .write_attr('long_name', 'Longitude') \ + .write_data(lon) + + # Latitude + obsspace.create_var('MetaData/latitude', dtype=lat.dtype, + fillval=lat.fill_value) \ + .write_attr('units', 'degrees_north') \ + .write_attr('valid_range', np.array([-90, 90], dtype=np.float32)) \ + .write_attr('long_name', 'Latitude') \ + .write_data(lat) + + # Datetime + obsspace.create_var('MetaData/dateTime', dtype=dateTime.dtype, + fillval=dateTime.fill_value) \ + .write_attr('units', 'seconds since 1970-01-01T00:00:00Z') \ + .write_attr('long_name', 'Datetime') \ + .write_data(dateTime) + + # Station Identification + obsspace.create_var('MetaData/stationIdentification', dtype=sid.dtype, + fillval=sid.fill_value) \ + .write_attr('long_name', 'Station Identification') \ + .write_data(sid) + + # Height Of Station + obsspace.create_var('MetaData/heightOfStation', dtype=zob.dtype, + fillval=zob.fill_value) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Height Of Station') \ + .write_data(zob) + + # Pressure + obsspace.create_var('MetaData/pressure', dtype=pressure.dtype, + fillval=pressure.fill_value) \ + .write_attr('units', 'Pa') \ + .write_attr('long_name', 'Pressure') \ + .write_data(pressure) + + # QualityMarker - Station Elevation + obsspace.create_var('QualityMarker/stationElevation', dtype=zobqm.dtype, + fillval=zobqm.fill_value) \ + .write_attr('long_name', 'Station Elevation Quality Marker') \ + .write_data(zobqm) + + # QualityMarker - Station Pressure + obsspace.create_var('QualityMarker/stationPressure', dtype=pobqm.dtype, + fillval=pobqm.fill_value) \ + .write_attr('long_name', 'Station Pressure Quality Marker') \ + .write_data(pobqm) + + # Station Elevation + obsspace.create_var('ObsValue/stationElevation', dtype=elv.dtype, + fillval=elv.fill_value) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Station Elevation') \ + .write_data(elv) + + # Station Pressure + obsspace.create_var('ObsValue/stationPressure', dtype=pob.dtype, + fillval=pob.fill_value) \ + .write_attr('units', 'Pa') \ + .write_attr('long_name', 'Station Pressure') \ + .write_data(pob) + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for splitting and output IODA: \ + {running_time} seconds") + + logger.debug("All Done!") + + +if __name__ == '__main__': + + start_time = time.time() + + parser = argparse.ArgumentParser() + parser.add_argument('-c', '--config', type=str, + help='Input JSON configuration', + required=True) + parser.add_argument('-v', '--verbose', + help='print debug logging information', + action='store_true') + args = parser.parse_args() + + log_level = 'DEBUG' if args.verbose else 'INFO' + logger = Logger('bufr2ioda_adpsfc_prepbufr.py', level=log_level, + colored_log=True) + + with open(args.config, "r") as json_file: + config = json.load(json_file) + + bufr_to_ioda(config, logger) + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Total running time: {running_time} seconds") diff --git a/ush/ioda/bufr2ioda/gen_bufr2ioda_json.py b/ush/ioda/bufr2ioda/gen_bufr2ioda_json.py index 617a8b024..0748e20f6 100755 --- a/ush/ioda/bufr2ioda/gen_bufr2ioda_json.py +++ b/ush/ioda/bufr2ioda/gen_bufr2ioda_json.py @@ -6,7 +6,8 @@ import argparse import json import os -from wxflow import Logger, parse_j2yaml +from wxflow import Logger, parse_j2yaml, cast_strdict_as_dtypedict +from wxflow import add_to_datetime, to_timedelta # Initialize root logger logger = Logger('gen_bufr2ioda_json.py', level='INFO', colored_log=True) @@ -28,4 +29,9 @@ def gen_bufr_json(config, template, output): parser.add_argument('-t', '--template', type=str, help='Input JSON template', required=True) parser.add_argument('-o', '--output', type=str, help='Output JSON file', required=True) args = parser.parse_args() + # get the config from your environment + config = cast_strdict_as_dtypedict(os.environ) + # we need to add in current cycle from PDYcyc + config['current_cycle'] = add_to_datetime(config['PDY'], to_timedelta(f"{config['cyc']}H")) + # call the parsing function gen_bufr_json(config, args.template, args.output) diff --git a/ush/ioda/bufr2ioda/run_bufr2ioda.py b/ush/ioda/bufr2ioda/run_bufr2ioda.py index d7a7bb130..63071d838 100755 --- a/ush/ioda/bufr2ioda/run_bufr2ioda.py +++ b/ush/ioda/bufr2ioda/run_bufr2ioda.py @@ -30,7 +30,7 @@ def bufr2ioda(current_cycle, RUN, DMPDIR, config_template_dir, COM_OBS): } # Specify observation types to be processed by a script - BUFR_py = ["satwind_amv_goes", "satwind_scat", "adpupa_prepbufr"] + BUFR_py = ["satwind_amv_goes", "satwind_scat", "adpupa_prepbufr", "adpsfc_prepbufr"] for obtype in BUFR_py: logger.info(f"Convert {obtype}...") From 19701ac72b804e66f5eb66f7c916b0091ae61ddd Mon Sep 17 00:00:00 2001 From: emilyhcliu <36091766+emilyhcliu@users.noreply.github.com> Date: Thu, 2 Nov 2023 15:58:26 -0400 Subject: [PATCH 4/9] Update conventional station pressure (conv_ps) for testing and end-to-end (#698) * Add a temporary hack of setting data window length for conv_ps. Due to error inflation for duplicate observations and the JEDI window length definition, the window length in the script has to be relaxed so that all observations from the 6-hour obs data will be accounted for. * Remove comments * Add assignment based on GSI error table, and update error inflation for pressure check. * Move filters two spacings left to line up with the filter header. * Add conv_ps.yaml under config for end-to-end test * Fix typo in comments * Fix typo * modify comments * Modify comments * Add and set window shift to true for the conventional difference in data window between GSI and JEDI. * Add logics to set window shift parameter * Add handling for tests with and without QC * Remove unwated comments * Remove debugging lines --------- Co-authored-by: Cory Martin --- parm/atm/obs/config/conv_ps.yaml | 295 ++++++++++++++++ parm/atm/obs/testing/conv_ps.yaml | 455 +++++++++++++++---------- parm/atm/obs/testing/conv_ps_noqc.yaml | 93 ++++- ush/ufoeval/run_ufo_hofx_test.sh | 4 +- ush/ufoeval/test_yamls.sh | 19 +- 5 files changed, 670 insertions(+), 196 deletions(-) create mode 100644 parm/atm/obs/config/conv_ps.yaml diff --git a/parm/atm/obs/config/conv_ps.yaml b/parm/atm/obs/config/conv_ps.yaml new file mode 100644 index 000000000..8fa4ab185 --- /dev/null +++ b/parm/atm/obs/config/conv_ps.yaml @@ -0,0 +1,295 @@ +obs space: + name: surface_ps + obsdatain: + engine: + type: H5File + obsfile: $(DATA)/obs/$(OPREFIX)sondes.{{ current_cycle | to_YMDH }}.nc4 + obsdataout: + engine: + type: H5File + obsfile: $(DATA)/diags/diag_sondes_{{ current_cycle | to_YMDH }}.nc4 + io pool: + max pool size: 1 + simulated variables: [stationPressure] + +obs operator: + name: SfcPCorrected + variables: + - name: stationPressure + da_psfc_scheme: GSI + station_altitude: height + geovar_sfc_geomz: surface_altitude + geovar_geomz: geopotential_height + +obs prior filters: +# Initial Error Assignments for SFC Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [181] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [110, 120, 120, 120, 120, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [187] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [85000, 80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 120, 140, 140, 140, 140, 140, 1.0e+11] + +# Initial Error Assignments for SFCSHIP Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [180] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [60000, 55000 ] + errors: [ 130, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [183] + action: + name: assign error + error parameter: 1.0e+11 + +# Initial Error Assignments for Radiosonde +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [120] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 110, 120, 120, 120, 120, 1.0e+11] + +obs post filters: +# Observation range sanity check +- filter: Bounds Check + filter variables: + - name: stationPressure + minvalue: 37499.0 + maxvalue: 106999.0 + action: + name: reject + +# Reject all ObsType 183 +- filter: RejectList + where: + - variable: + name: ObsType/stationPressure + is_in: 183 + +# Reject surface pressure below 500 hPa +- filter: Bounds Check + filter variables: + - name: stationPressure + minvalue: 50000.00 + action: + name: reject + +- filter: RejectList + where: + - variable: + name: PreQC/stationPressure + is_in: 4-15 + +# Inflate obs error based on obs type +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_in: 3, 7 + action: + name: inflate error + inflation factor: 1.2 + +# Calculate obs error inflation factors for duplicated observations at the same location +- filter: Variable Assignment + assignments: + - name: ObsErrorFactorDuplicateCheck/stationPressure + type: float + function: + name: ObsFunction/ObsErrorFactorDuplicateCheck + options: + use_air_pressure: false + variable: stationPressure + +# Reduce effective observation error based on obs type and subtype +# In this case: reduce effective obs error for buoy +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: + name: ObsType/stationPressure + is_in: 180 + - variable: + name: ObsSubType/stationPressure + is_in: 0 + action: + name: inflate error + inflation factor: 0.7 + +# Reduce original observation error based on obs type and subtype +# In this case: reduce original obs error for buoy +- filter: Variable Assignment + where: + - variable: + name: ObsType/stationPressure + is_in: 180 + - variable: + name: ObsSubType/stationPressure + is_in: 0 + assignments: + - name: ObsError/stationPressure + type: float + function: + name: ObsFunction/Arithmetic + options: + variables: + - name: ObsError/stationPressure + coefs: [0.7] + +# Calculate obs error inflation factors for large discrepancies between model and observations +- filter: Variable Assignment + assignments: + - name: ObsErrorFactorSfcPressure/stationPressure + type: float + function: + name: ObsFunction/ObsErrorFactorSfcPressure + options: + geovar_sfc_geomz: surface_altitude + geovar_geomz: geopotential_height + station_altitude: height + +# Inflate surface pressure observation based on discrepancies between +# model and observations due to terrian +- filter: Perform Action + filter variables: + - name: stationPressure + action: + name: inflate error + inflation variable: + name: ObsErrorFactorSfcPressure/stationPressure + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/Innovation + type: float + function: + name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/stationPressure + - name: HofX/stationPressure + coefs: [1, -1] + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorBoundSfcPressure1 + type: float + function: + name: ObsFunction/ObsErrorBoundConventional + options: + obsvar: stationPressure + obserr_bound_min: 100 + obserr_bound_max: 300 + obserr_bound_factor: 5.0 + +- filter: Background Check + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_not_in: 3 + function absolute threshold: + - name: DerivedMetaData/ObsErrorBoundSfcPressure1 + action: + name: reject + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorBoundSfcPressure2 + type: float + function: + name: ObsFunction/ObsErrorBoundConventional + options: + obsvar: stationPressure + obserr_bound_min: 100 + obserr_bound_max: 300 + obserr_bound_factor: 3.5 + +- filter: Background Check + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_in: 3 + function absolute threshold: + - name: DerivedMetaData/ObsErrorBoundSfcPressure2 + action: + name: reject + +# Inflate obs error based on duplicate check +- filter: Perform Action + filter variables: + - name: stationPressure + action: + name: inflate error + inflation variable: + name: ObsErrorFactorDuplicateCheck/stationPressure + +# Reject data based on PreUseFlag (usage in GSI) +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: PreUseFlag/stationPressure + is_not_in: 0, 1 + action: + name: reject + +# End of Filters + diff --git a/parm/atm/obs/testing/conv_ps.yaml b/parm/atm/obs/testing/conv_ps.yaml index f43bc9f41..20a449a8a 100644 --- a/parm/atm/obs/testing/conv_ps.yaml +++ b/parm/atm/obs/testing/conv_ps.yaml @@ -8,8 +8,8 @@ obs space: engine: type: H5File obsfile: !ENV conv_ps_diag_${CDATE}.nc4 - overwrite: true simulated variables: [stationPressure] + geovals: filename: !ENV conv_ps_geoval_${CDATE}.nc4 @@ -22,191 +22,276 @@ obs operator: geovar_sfc_geomz: surface_altitude geovar_geomz: geopotential_height +obs prior filters: +# Initial Error Assignments for SFC Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [181] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [110, 120, 120, 120, 120, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [187] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [85000, 80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 120, 140, 140, 140, 140, 140, 1.0e+11] + +# Initial Error Assignments for SFCSHIP Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [180] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [60000, 55000 ] + errors: [ 130, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [183] + action: + name: assign error + error parameter: 1.0e+11 + +# Initial Error Assignments for Radiosonde +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [120] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 110, 120, 120, 120, 120, 1.0e+11] + obs post filters: - # Observation range sanity check - - filter: Bounds Check - filter variables: - - name: stationPressure - minvalue: 37499.0 - maxvalue: 106999.0 - action: - name: reject - - # Reject all ObsType 183 - - filter: RejectList - where: - - variable: - name: ObsType/stationPressure - is_in: 183 - - # Reject surface pressure below 500 hPa - - filter: Bounds Check - filter variables: - - name: stationPressure - minvalue: 50000.00 - action: - name: reject - - - filter: RejectList - where: - - variable: - name: PreQC/stationPressure - is_in: 4-15 - - # Inflate obs error based on obs type - - filter: Perform Action - filter variables: - - name: stationPressure - where: - - variable: PreQC/stationPressure - is_in: 3, 7 - action: - name: inflate error - inflation factor: 1.2 - - # Calculate obs error inflation factors for duplicated observations at the same location - - filter: Variable Assignment - assignments: - - name: ObsErrorFactorDuplicateCheck/stationPressure - type: float - function: - name: ObsFunction/ObsErrorFactorDuplicateCheck - options: - use_air_pressure: false - variable: stationPressure - - # Reduce effective observation error based on obs type and subtype - # In this case: reduce effective obs error for buoy - - filter: Perform Action - filter variables: - - name: stationPressure - where: - - variable: - name: ObsType/stationPressure - is_in: 180 - - variable: - name: ObsSubType/stationPressure - is_in: 0 - action: - name: inflate error - inflation factor: 0.7 - - # Reduce original observation error based on obs type and subtype - # In this case: reduce original obs error for buoy - - filter: Variable Assignment - where: - - variable: - name: ObsType/stationPressure - is_in: 180 - - variable: - name: ObsSubType/stationPressure - is_in: 0 - assignments: - - name: ObsError/stationPressure - type: float - function: - name: ObsFunction/Arithmetic - options: - variables: - - name: ObsError/stationPressure - coefs: [0.7] - - # Calculate obs error inflation factors for large discrepancies between model and observations - - filter: Variable Assignment - assignments: - - name: ObsErrorFactorSfcPressure/stationPressure - type: float - function: - name: ObsFunction/ObsErrorFactorSfcPressure - options: - geovar_sfc_geomz: surface_altitude - - # Inflate surface pressure observation based on discrepancies between - # model and observations due to terrian - - filter: Perform Action - filter variables: - - name: stationPressure - action: - name: inflate error - inflation variable: - name: ObsErrorFactorSfcPressure/stationPressure - - - filter: Variable Assignment - assignments: - - name: DerivedMetaData/Innovation - type: float - function: - name: ObsFunction/Arithmetic - options: - variables: - - name: ObsValue/stationPressure - - name: HofX/stationPressure - coefs: [1, -1] - - - filter: Variable Assignment - assignments: - - name: DerivedMetaData/ObsErrorBoundSfcPressure1 - type: float - function: - name: ObsFunction/ObsErrorBoundConventional - options: - obsvar: stationPressure - obserr_bound_min: 100 - obserr_bound_max: 300 - obserr_bound_factor: 5.0 - - - filter: Background Check - filter variables: - - name: stationPressure - where: - - variable: PreQC/stationPressure - is_not_in: 3 - function absolute threshold: - - name: DerivedMetaData/ObsErrorBoundSfcPressure1 - action: - name: reject - - - filter: Variable Assignment - assignments: - - name: DerivedMetaData/ObsErrorBoundSfcPressure2 - type: float - function: - name: ObsFunction/ObsErrorBoundConventional - options: - obsvar: stationPressure - obserr_bound_min: 100 - obserr_bound_max: 300 - obserr_bound_factor: 3.5 - - - filter: Background Check - filter variables: - - name: stationPressure - where: - - variable: PreQC/stationPressure - is_in: 3 - function absolute threshold: - - name: DerivedMetaData/ObsErrorBoundSfcPressure2 - action: - name: reject - - # Inflate obs error based on duplicate check - - filter: Perform Action - filter variables: - - name: stationPressure - action: - name: inflate error - inflation variable: - name: ObsErrorFactorDuplicateCheck/stationPressure - - # Reject data based on PreUseFlag (usage in GSI) - - filter: Perform Action - filter variables: - - name: stationPressure - where: - - variable: PreUseFlag/stationPressure - is_not_in: 0, 1 - action: - name: reject +# Observation range sanity check +- filter: Bounds Check + filter variables: + - name: stationPressure + minvalue: 37499.0 + maxvalue: 106999.0 + action: + name: reject + +# Reject all ObsType 183 +- filter: RejectList + where: + - variable: + name: ObsType/stationPressure + is_in: 183 + +# Reject surface pressure below 500 hPa +- filter: Bounds Check + filter variables: + - name: stationPressure + minvalue: 50000.00 + action: + name: reject + +- filter: RejectList + where: + - variable: + name: PreQC/stationPressure + is_in: 4-15 + +# Inflate obs error based on obs type +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_in: 3, 7 + action: + name: inflate error + inflation factor: 1.2 + +# Calculate obs error inflation factors for duplicated observations at the same location +- filter: Variable Assignment + assignments: + - name: ObsErrorFactorDuplicateCheck/stationPressure + type: float + function: + name: ObsFunction/ObsErrorFactorDuplicateCheck + options: + use_air_pressure: false + variable: stationPressure + +# Reduce effective observation error based on obs type and subtype +# In this case: reduce effective obs error for buoy +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: + name: ObsType/stationPressure + is_in: 180 + - variable: + name: ObsSubType/stationPressure + is_in: 0 + action: + name: inflate error + inflation factor: 0.7 + +# Reduce original observation error based on obs type and subtype +# In this case: reduce original obs error for buoy +- filter: Variable Assignment + where: + - variable: + name: ObsType/stationPressure + is_in: 180 + - variable: + name: ObsSubType/stationPressure + is_in: 0 + assignments: + - name: ObsError/stationPressure + type: float + function: + name: ObsFunction/Arithmetic + options: + variables: + - name: ObsError/stationPressure + coefs: [0.7] + +# Calculate obs error inflation factors for large discrepancies between model and observations +- filter: Variable Assignment + assignments: + - name: ObsErrorFactorSfcPressure/stationPressure + type: float + function: + name: ObsFunction/ObsErrorFactorSfcPressure + options: + geovar_sfc_geomz: surface_altitude + geovar_geomz: geopotential_height + station_altitude: height + +# Inflate surface pressure observation based on discrepancies between +# model and observations due to terrian +- filter: Perform Action + filter variables: + - name: stationPressure + action: + name: inflate error + inflation variable: + name: ObsErrorFactorSfcPressure/stationPressure + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/Innovation + type: float + function: + name: ObsFunction/Arithmetic + options: + variables: + - name: ObsValue/stationPressure + - name: HofX/stationPressure + coefs: [1, -1] + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorBoundSfcPressure1 + type: float + function: + name: ObsFunction/ObsErrorBoundConventional + options: + obsvar: stationPressure + obserr_bound_min: 100 + obserr_bound_max: 300 + obserr_bound_factor: 5.0 + +- filter: Background Check + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_not_in: 3 + function absolute threshold: + - name: DerivedMetaData/ObsErrorBoundSfcPressure1 + action: + name: reject + +- filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorBoundSfcPressure2 + type: float + function: + name: ObsFunction/ObsErrorBoundConventional + options: + obsvar: stationPressure + obserr_bound_min: 100 + obserr_bound_max: 300 + obserr_bound_factor: 3.5 + +- filter: Background Check + filter variables: + - name: stationPressure + where: + - variable: PreQC/stationPressure + is_in: 3 + function absolute threshold: + - name: DerivedMetaData/ObsErrorBoundSfcPressure2 + action: + name: reject + +# Inflate obs error based on duplicate check +- filter: Perform Action + filter variables: + - name: stationPressure + action: + name: inflate error + inflation variable: + name: ObsErrorFactorDuplicateCheck/stationPressure + +# Reject data based on PreUseFlag (usage in GSI) +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: PreUseFlag/stationPressure + is_not_in: 0, 1 + action: + name: reject +# End of Filters passedBenchmark: 85378 diff --git a/parm/atm/obs/testing/conv_ps_noqc.yaml b/parm/atm/obs/testing/conv_ps_noqc.yaml index 44a6ef4c8..8d469729f 100644 --- a/parm/atm/obs/testing/conv_ps_noqc.yaml +++ b/parm/atm/obs/testing/conv_ps_noqc.yaml @@ -8,8 +8,8 @@ obs space: engine: type: H5File obsfile: !ENV conv_ps_diag_${CDATE}.nc4 - overwrite: true simulated variables: [stationPressure] + geovals: filename: !ENV conv_ps_geoval_${CDATE}.nc4 @@ -22,10 +22,89 @@ obs operator: geovar_sfc_geomz: surface_altitude geovar_geomz: geopotential_height -vector ref: GsiHofXBc -tolerance: 1.e-4 -#linear obs operator test: -# coef TL: 0.1 -# tolerance TL: 1.0e-13 -# tolerance AD: 1.0e-11 +obs prior filters: +# Initial Error Assignments for SFC Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [181] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [110, 120, 120, 120, 120, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [187] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [85000, 80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 120, 140, 140, 140, 140, 140, 1.0e+11] + +# Initial Error Assignments for SFCSHIP Observations +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [180] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [60000, 55000 ] + errors: [ 130, 1.0e+11] + +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [183] + action: + name: assign error + error parameter: 1.0e+11 + + # Initial Error Assignments for Radiosonde +- filter: Perform Action + filter variables: + - name: stationPressure + where: + - variable: ObsType/stationPressure + is_in: [120] + action: + name: assign error + error function: + name: ObsFunction/ObsErrorModelStepwiseLinear + options: + round_to_the_nearest_integer: true + xvar: + name: ObsValue/stationPressure + xvals: [80000, 75000, 70000, 65000, 60000, 55000 ] + errors: [ 110, 120, 120, 120, 120, 1.0e+11] + +passedBenchmark: 92824 # total: 92842; missing: 18 +#vector ref: GsiHofXBc +#tolerance: 1.e-4 diff --git a/ush/ufoeval/run_ufo_hofx_test.sh b/ush/ufoeval/run_ufo_hofx_test.sh index ccc043bd4..230966284 100755 --- a/ush/ufoeval/run_ufo_hofx_test.sh +++ b/ush/ufoeval/run_ufo_hofx_test.sh @@ -33,6 +33,7 @@ run_filtering=YES run_eva=YES eva_stats_only=NO keep_output=NO +window_shift=False while getopts "c:hsxq" opt; do case $opt in @@ -99,7 +100,6 @@ if [ $run_filtering == NO ]; then else yamlpath=$GDASApp/parm/atm/obs/testing/${obtype}.yaml fi - exename=test_ObsFilters.x #-------------- Do not modify below this line ---------------- @@ -178,11 +178,13 @@ export OPREFIX=gdas.t${cyc}z export APREFIX=gdas.t${cyc}z export GPREFIX=gdas.t${gcyc}z +if [ $obtype == "conv_ps" ]; then window_shift=True; fi cat > $workdir/temp.yaml << EOF window begin: '{{ WINDOW_BEGIN | to_isotime }}' window end: '{{ WINDOW_END | to_isotime }}' observations: - !INC $yamlpath +window shift: ${window_shift} EOF $GDASApp/ush/genYAML --input $workdir/temp.yaml --output $workdir/${obtype}_${cycle}.yaml diff --git a/ush/ufoeval/test_yamls.sh b/ush/ufoeval/test_yamls.sh index 7e98d80e2..d46056bd8 100755 --- a/ush/ufoeval/test_yamls.sh +++ b/ush/ufoeval/test_yamls.sh @@ -13,10 +13,23 @@ if [ $? -ne 0 ]; then exit 1 fi -for file in `ls ../../parm/atm/obs/testing/*yaml`; do +# Process tests wiht QC +for file in `find ../../parm/atm/obs/testing/*.yaml -type f -not -name "*noqc*"`; do basefile=${file##*/} - inst="${basefile%.*}" - ./run_ufo_hofx_test.sh -x $inst > $WORKDIR/$inst.log 2> $WORKDIR/$inst.err + obtype="${basefile%.*}" + ./run_ufo_hofx_test.sh -x $obtype > $WORKDIR/$obtype.log 2> $WORKDIR/$obtype.err + if [ $? == 0 ]; then + echo $basefile Passes \(yay!\) + else + echo $basefile Fails \(boo!\) + fi +done + +# Process tests without QC (HofX + Observation error assignment) +for file in `ls ../../parm/atm/obs/testing/*_noqc.yaml`; do + basefile=${file##*/} + obtype="${basefile%_noqc.*}" + ./run_ufo_hofx_test.sh -x -q $obtype > $WORKDIR/${obtype}_noqc.log 2> $WORKDIR/${obtype}_noqc.err if [ $? == 0 ]; then echo $basefile Passes \(yay!\) else From a0616fdac531d778ef6d90bbe538bd43b3b36b0c Mon Sep 17 00:00:00 2001 From: NicholasEsposito-NOAA <62616739+nicholasesposito@users.noreply.github.com> Date: Thu, 2 Nov 2023 16:43:51 -0400 Subject: [PATCH 5/9] Add sfcship prepBUFR conversion (#703) * working sfcshp * rm some comments * small cleanup * update README * get rid of sys.path.append * now works again * logger and code norms * remove some spaces * 1 more thing * rm some remaining prints * reference time print * Update bufr2ioda_sfcshp_prepbufr.py height, var changes * some mods * some elv zob updates * aesthetic * reame chamge * changes * some small changes * revert and fix coding norms * rm unneeded files * json update * run fixes --------- Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito Co-authored-by: Nicholas Esposito Co-authored-by: Cory Martin --- .../bufr2ioda/bufr2ioda_sfcshp_prepbufr.json | 12 + .../bufr2ioda/bufr2ioda_sfcshp_prepbufr.py | 515 ++++++++++++++++++ ush/ioda/bufr2ioda/run_bufr2ioda.py | 3 +- 3 files changed, 529 insertions(+), 1 deletion(-) create mode 100644 parm/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.json create mode 100755 ush/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.py diff --git a/parm/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.json b/parm/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.json new file mode 100644 index 000000000..4bdc8f781 --- /dev/null +++ b/parm/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.json @@ -0,0 +1,12 @@ +{ + "data_format" : "prepbufr", + "subsets" : [ "SFCSHP" ], + "source" : "prepBUFR", + "data_type" : "SFCSHP", + "cycle_type" : "{{ RUN }}", + "cycle_datetime" : "{{ current_cycle | to_YMDH }}", + "dump_directory" : "{{ DMPDIR }}", + "ioda_directory" : "{{ COM_OBS }}", + "data_description" : "SFCSHP_prepbufr", + "data_provider" : "U.S. NOAA" +} diff --git a/ush/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.py b/ush/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.py new file mode 100755 index 000000000..47207a3e9 --- /dev/null +++ b/ush/ioda/bufr2ioda/bufr2ioda_sfcshp_prepbufr.py @@ -0,0 +1,515 @@ +#!/usr/bin/env python3 +# (C) Copyright 2023 NOAA/NWS/NCEP/EMC +# +# This software is licensed under the terms of the Apache Licence Version 2.0 +# which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. + +import sys +import numpy as np +import os +import argparse +import math +import calendar +import time +from datetime import datetime +import json +from pyiodaconv import bufr +from collections import namedtuple +from pyioda import ioda_obs_space as ioda_ospace +from wxflow import Logger + + +def Compute_dateTime(cycleTimeSinceEpoch, dhr): + + dhr = np.int64(dhr*3600) + dateTime = dhr + cycleTimeSinceEpoch + + return dateTime + + +def bufr_to_ioda(config, logger): + + subsets = config["subsets"] + logger.debug(f"Checking subsets = {subsets}") + + # Get parameters from configuration + data_format = config["data_format"] + source = config["source"] + data_type = config["data_type"] + data_description = config["data_description"] + data_provider = config["data_provider"] + cycle_type = config["cycle_type"] + cycle_datetime = config["cycle_datetime"] + dump_dir = config["dump_directory"] + ioda_dir = config["ioda_directory"] + cycle = config["cycle_datetime"] + + # Get derived parameters + yyyymmdd = cycle[0:8] + hh = cycle[8:10] + + reference_time = datetime.strptime(cycle, "%Y%m%d%H") + reference_time = reference_time.strftime("%Y-%m-%dT%H:%M:%SZ") + reference_time_full = f"{yyyymmdd}{hh}00" + + logger.debug(f"Reference time = {reference_time}") + + # General informaton + converter = 'BUFR to IODA Converter' + platform_description = 'SFCSHP data from prepBUFR format' + + bufrfile = f"{cycle_type}.t{hh}z.{data_format}" + DATA_PATH = os.path.join(dump_dir, f"{cycle_type}.{yyyymmdd}", + str(hh), bufrfile) + + logger.debug(f"The DATA_PATH is: {DATA_PATH}") + + # ============================================ + # Make the QuerySet for all the data we want + # ============================================ + start_time = time.time() + + logger.debug("Making QuerySet ...") + q = bufr.QuerySet(subsets) + + # ObsType + q.add('observationType', '*/TYP') + + # MetaData + q.add('stationIdentification', '*/SID') + q.add('latitude', '*/YOB') + q.add('longitude', '*/XOB') + q.add('obsTimeMinusCycleTime', '*/DHR') + q.add('heightOfStation', '*/Z___INFO/Z__EVENT{1}/ZOB') + q.add('pressure', '*/P___INFO/P__EVENT{1}/POB') + q.add('temperatureEventCode', '*/T___INFO/T__EVENT{1}/TPC') + +# # Quality Infomation (Quality Indicator) + q.add('qualityMarkerStationElevation', '*/Z___INFO/Z__EVENT{1}/ZQM') + q.add('qualityMarkerStationPressure', '*/P___INFO/P__EVENT{1}/PQM') + q.add('qualityMarkerAirTemperature', '*/T___INFO/T__EVENT{1}/TQM') + q.add('qualityMarkerSpecificHumidity', '*/Q___INFO/Q__EVENT{1}/QQM') + q.add('qualityMarkerWindNorthward', '*/W___INFO/W__EVENT{1}/WQM') + q.add('qualityMarkerSeaSurfaceTemperature', '*/SST_INFO/SSTEVENT{1}/SSTQM') + + # ObsValue + q.add('stationElevation', '*/ELV') + q.add('stationPressure', '*/P___INFO/P__EVENT{1}/POB') + q.add('airTemperature', '*/T___INFO/T__EVENT{1}/TOB') + q.add('specificHumidity', '*/Q___INFO/Q__EVENT{1}/QOB') + q.add('windNorthward', '*/W___INFO/W__EVENT{1}/VOB') + q.add('windEastward', '*/W___INFO/W__EVENT{1}/UOB') + q.add('seaSurfaceTemperature', '*/SST_INFO/SSTEVENT{1}/SST1') + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for making QuerySet: {running_time} seconds") + + # ============================================================== + # Open the BUFR file and execute the QuerySet to get ResultSet + # Use the ResultSet returned to get numpy arrays of the data + # ============================================================== + start_time = time.time() + + logger.debug(f"Executing QuerySet to get ResultSet ...") + with bufr.File(DATA_PATH) as f: + r = f.execute(q) + + logger.debug(f" ... Executing QuerySet: get metadata: basic ...") + # ObsType + logger.debug(" ... Executing QuerySet: get ObsType ...") + typ = r.get('observationType') + + # MetaData + logger.debug(" ... Executing QuerySet: MetaData ...") + sid = r.get('stationIdentification') + lat = r.get('latitude') + lon = r.get('longitude') + lon[lon > 180] -= 360 + zob = r.get('heightOfStation', type='float') + pressure = r.get('pressure') + pressure *= 100 + tpc = r.get('temperatureEventCode') + + # Quality Information + logger.debug(f" ... Executing QuerySet: QualityMarker ...") + zobqm = r.get('qualityMarkerStationElevation') + pobqm = r.get('qualityMarkerStationPressure') + tobqm = r.get('qualityMarkerAirTemperature') + tsenqm = np.full(tobqm.shape[0], tobqm.fill_value) + tsenqm = np.where(((tpc >= 1) & (tpc < 8)), tobqm, tsenqm) + tvoqm = np.full(tobqm.shape[0], tobqm.fill_value) + tvoqm = np.where((tpc == 8), tobqm, tvoqm) + qobqm = r.get('qualityMarkerSpecificHumidity') + wobqm = r.get('qualityMarkerWindNorthward') + sstqm = r.get('qualityMarkerSeaSurfaceTemperature') + + logger.debug(f" ... Executing QuerySet: ObsValue ...") + # ObsValue + elv = r.get('stationElevation', type='float') + pob = r.get('stationPressure') + pob *= 100 + tob = r.get('airTemperature') + tob += 273.15 + tsen = np.full(tob.shape[0], tob.fill_value) + tsen = np.where(((tpc >= 1) & (tpc < 8)), tob, tsen) + tvo = np.full(tob.shape[0], tob.fill_value) + tvo = np.where((tpc == 8), tob, tvo) + qob = r.get('specificHumidity', type='float') + qob *= 0.000001 + uob = r.get('windEastward') + vob = r.get('windNorthward') + sst1 = r.get('seaSurfaceTemperature') + + logger.debug(f" ... Executing QuerySet: get datatime: observation time ...") + # DateTime: seconds since Epoch time + # IODA has no support for numpy datetime arrays dtype=datetime64[s] + dhr = r.get('obsTimeMinusCycleTime', type='int64') + + logger.debug(f" ... Executing QuerySet: Done!") + + logger.debug(f" ... Executing QuerySet: Check BUFR variable generic \ + dimension and type ...") + # Check BUFR variable generic dimension and type + logger.debug(f" typ shape = {typ.shape}") + logger.debug(f" sid shape = {sid.shape}") + logger.debug(f" dhr shape = {dhr.shape}") + logger.debug(f" lat shape = {lat.shape}") + logger.debug(f" lon shape = {lon.shape}") + logger.debug(f" zob shape = {zob.shape}") + logger.debug(f" pressure shape = {pressure.shape}") + logger.debug(f" tpc shape = {tpc.shape}") + + logger.debug(f" zobqm shape = {zobqm.shape}") + logger.debug(f" pobqm shape = {pobqm.shape}") + logger.debug(f" tobqm shape = {tobqm.shape}") + logger.debug(f" tsenqm shape = {tsenqm.shape}") + logger.debug(f" tvoqm shape = {tvoqm.shape}") + logger.debug(f" qobqm shape = {qobqm.shape}") + logger.debug(f" wobqm shape = {wobqm.shape}") + logger.debug(f" sstqm shape = {sstqm.shape}") + + logger.debug(f" elv shape = {elv.shape}") + logger.debug(f" pob shape = {pob.shape}") + logger.debug(f" tob shape = {tob.shape}") + logger.debug(f" tsen shape = {tsen.shape}") + logger.debug(f" tvo shape = {tvo.shape}") + logger.debug(f" qob shape = {qob.shape}") + logger.debug(f" uob shape = {uob.shape}") + logger.debug(f" vob shape = {vob.shape}") + logger.debug(f" sst1 shape = {sst1.shape}") + + logger.debug(f" typ type = {typ.dtype}") + logger.debug(f" sid type = {sid.dtype}") + logger.debug(f" dhr type = {dhr.dtype}") + logger.debug(f" lat type = {lat.dtype}") + logger.debug(f" lon type = {lon.dtype}") + logger.debug(f" zob type = {zob.dtype}") + logger.debug(f" pressure type = {pressure.dtype}") + logger.debug(f" tpc type = {tpc.dtype}") + + logger.debug(f" pobqm type = {pobqm.dtype}") + logger.debug(f" tobqm type = {tobqm.dtype}") + logger.debug(f" tsenqm type = {tsenqm.dtype}") + logger.debug(f" tvoqm type = {tvoqm.dtype}") + logger.debug(f" qobqm type = {qobqm.dtype}") + logger.debug(f" wobqm type = {wobqm.dtype}") + logger.debug(f" sstqm type = {sstqm.dtype}") + + logger.debug(f" elv type = {elv.dtype}") + logger.debug(f" pob type = {pob.dtype}") + logger.debug(f" tob type = {tob.dtype}") + logger.debug(f" tsen type = {tsen.dtype}") + logger.debug(f" tvo type = {tvo.dtype}") + logger.debug(f" qob type = {qob.dtype}") + logger.debug(f" uob type = {uob.dtype}") + logger.debug(f" vob type = {vob.dtype}") + logger.debug(f" sst1 type = {sst1.dtype}") + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for executing QuerySet to get ResultSet: \ + {running_time} seconds") + + # ========================= + # Create derived variables + # ========================= + start_time = time.time() + + logger.debug(f"Creating derived variables - dateTime ...") + + cycleTimeSinceEpoch = np.int64(calendar.timegm(time.strptime( + reference_time_full, '%Y%m%d%H%M'))) + dateTime = Compute_dateTime(cycleTimeSinceEpoch, dhr) + + logger.debug(f" Check derived variables type ... ") + logger.debug(f" dateTime shape = {dateTime.shape}") + logger.debug(f" dateTime type = {dateTime.dtype}") + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for creating derived variables: {running_time} \ + seconds") + + # ===================================== + # Create IODA ObsSpace + # Write IODA output + # ===================================== + + # Create the dimensions + dims = {'Location': np.arange(0, lat.shape[0])} + + iodafile = f"{cycle_type}.t{hh}z.{data_type}.{data_format}.nc" + OUTPUT_PATH = os.path.join(ioda_dir, iodafile) + logger.debug(f" ... ... Create OUTPUT file: {OUTPUT_PATH}") + + path, fname = os.path.split(OUTPUT_PATH) + if path and not os.path.exists(path): + os.makedirs(path) + + obsspace = ioda_ospace.ObsSpace(OUTPUT_PATH, mode='w', dim_dict=dims) + + # Create Global attributes + logger.debug(f" ... ... Create global attributes") + + obsspace.write_attr('Converter', converter) + obsspace.write_attr('source', source) + obsspace.write_attr('sourceFiles', bufrfile) + obsspace.write_attr('dataProviderOrigin', data_provider) + obsspace.write_attr('description', data_description) + obsspace.write_attr('datetimeReference', reference_time) + obsspace.write_attr('datetimeRange', + [str(min(dateTime)), str(max(dateTime))]) + obsspace.write_attr('platformLongDescription', platform_description) + + # Create IODA variables + logger.debug(f" ... ... Create variables: name, type, units, & attributes") + + # ObsType: stationElevation + obsspace.create_var('ObsType/stationElevation', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Station Elevation Observation Type') \ + .write_data(typ) + + # ObsType: stationPressure + obsspace.create_var('ObsType/stationPressure', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Station Pressure Observation Type') \ + .write_data(typ) + + # ObsType: air Temperature + obsspace.create_var('ObsType/airTemperature', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Air Temperature Observation Type') \ + .write_data(typ) + + # ObsType: virtual Temperature + obsspace.create_var('ObsType/virtualTemperature', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Virtual Temperature Observation Type') \ + .write_data(typ) + + # ObsType: Specific Humidity + obsspace.create_var('ObsType/specificHumidity', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Specific Humidity Observation Type') \ + .write_data(typ) + + # ObsType: wind Eastward + obsspace.create_var('ObsType/windEastward', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Wind Eastward Observation Type') \ + .write_data(typ) + + # ObsType: wind Northward + obsspace.create_var('ObsType/windNorthward', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Wind Northward Observation Type') \ + .write_data(typ) + + # ObsType: Sea Surface Temperature + obsspace.create_var('ObsType/seaSurfaceTemperature', dtype=typ.dtype, + fillval=typ.fill_value) \ + .write_attr('long_name', 'Sea Surface Temperature Observation Type') \ + .write_data(typ) + + # Longitude + obsspace.create_var('MetaData/longitude', dtype=lon.dtype, + fillval=lon.fill_value) \ + .write_attr('units', 'degrees_east') \ + .write_attr('valid_range', np.array([-180, 180], dtype=np.float32)) \ + .write_attr('long_name', 'Longitude') \ + .write_data(lon) + + # Latitude + obsspace.create_var('MetaData/latitude', dtype=lat.dtype, + fillval=lat.fill_value) \ + .write_attr('units', 'degrees_north') \ + .write_attr('valid_range', np.array([-90, 90], dtype=np.float32)) \ + .write_attr('long_name', 'Latitude') \ + .write_data(lat) + + # Datetime + obsspace.create_var('MetaData/dateTime', dtype=dateTime.dtype, + fillval=dateTime.fill_value) \ + .write_attr('units', 'seconds since 1970-01-01T00:00:00Z') \ + .write_attr('long_name', 'Datetime') \ + .write_data(dateTime) + + # Station Identification + obsspace.create_var('MetaData/stationIdentification', dtype=sid.dtype, + fillval=sid.fill_value) \ + .write_attr('long_name', 'Station Identification') \ + .write_data(sid) + + # Station Elevation + obsspace.create_var('MetaData/heightOfStation', dtype=zob.dtype, + fillval=zob.fill_value) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Height Of Station') \ + .write_data(zob) + + # Pressure + obsspace.create_var('MetaData/pressure', dtype=pressure.dtype, + fillval=pressure.fill_value) \ + .write_attr('units', 'Pa') \ + .write_attr('long_name', 'Pressure') \ + .write_data(pressure) + + # Temperature Event Code + obsspace.create_var('MetaData/temperatureEventCode', dtype=tpc.dtype, + fillval=tpc.fill_value) \ + .write_attr('long_name', 'Temperature Event Code') \ + .write_data(tpc) + + # Quality Marker: Station Pressure + obsspace.create_var('QualityMarker/stationPressure', dtype=pobqm.dtype, + fillval=pobqm.fill_value) \ + .write_attr('long_name', 'Station Pressure Quality Marker') \ + .write_data(pobqm) + + # Quality Marker: Air Temperature + obsspace.create_var('QualityMarker/airTemperature', dtype=tobqm.dtype, + fillval=tobqm.fill_value) \ + .write_attr('long_name', 'Air Temperature Quality Marker') \ + .write_data(tsenqm) + + # Quality Marker: Virtual Temperature + obsspace.create_var('QualityMarker/virtualTemperature', dtype=tobqm.dtype, + fillval=tobqm.fill_value) \ + .write_attr('long_name', 'Virtual Temperature Quality Marker') \ + .write_data(tvoqm) + + # Quality Marker: Specific Humidity + obsspace.create_var('QualityMarker/specificHumidity', dtype=qobqm.dtype, + fillval=qobqm.fill_value) \ + .write_attr('long_name', 'Specific Humidity Quality Marker') \ + .write_data(qobqm) + + # Quality Marker: Eastward Wind + obsspace.create_var('QualityMarker/windEastward', dtype=wobqm.dtype, + fillval=wobqm.fill_value) \ + .write_attr('long_name', 'Eastward Wind Quality Marker') \ + .write_data(wobqm) + + # Quality Marker: Northward Wind + obsspace.create_var('QualityMarker/windNorthward', dtype=wobqm.dtype, + fillval=wobqm.fill_value) \ + .write_attr('long_name', 'Northward Wind Quality Marker') \ + .write_data(wobqm) + + # Quality Marker: Sea Surface Temperature + obsspace.create_var('QualityMarker/seaSurfaceTemperature', + dtype=sstqm.dtype, fillval=sstqm.fill_value) \ + .write_attr('long_name', 'Sea Surface Temperature Quality Marker') \ + .write_data(sstqm) + + # Station Elevation + obsspace.create_var('ObsValue/stationElevation', dtype=elv.dtype, + fillval=elv.fill_value) \ + .write_attr('units', 'm') \ + .write_attr('long_name', 'Station Elevation') \ + .write_data(elv) + + # Station Pressure + obsspace.create_var('ObsValue/stationPressure', dtype=pob.dtype, + fillval=pob.fill_value) \ + .write_attr('units', 'Pa') \ + .write_attr('long_name', 'Station Pressure') \ + .write_data(pob) + + # Air Temperature + obsspace.create_var('ObsValue/airTemperature', dtype=tob.dtype, + fillval=tob.fill_value) \ + .write_attr('units', 'K') \ + .write_attr('long_name', 'Air Temperature') \ + .write_data(tsen) + + # Virtual Temperature + obsspace.create_var('ObsValue/virtualTemperature', dtype=tob.dtype, + fillval=tob.fill_value) \ + .write_attr('units', 'K') \ + .write_attr('long_name', 'Virtual Temperature') \ + .write_data(tvo) + + # Specific Humidity + obsspace.create_var('ObsValue/specificHumidity', dtype=qob.dtype, + fillval=qob.fill_value) \ + .write_attr('units', 'kg kg-1') \ + .write_attr('long_name', 'Specific Humidity') \ + .write_data(qob) + + # Eastward Wind + obsspace.create_var('ObsValue/windEastward', dtype=uob.dtype, + fillval=uob.fill_value) \ + .write_attr('units', 'm s-1') \ + .write_attr('long_name', 'Eastward Wind') \ + .write_data(uob) + + # Northward Wind + obsspace.create_var('ObsValue/windNorthward', dtype=vob.dtype, + fillval=vob.fill_value) \ + .write_attr('units', 'm s-1') \ + .write_attr('long_name', 'Northward Wind') \ + .write_data(vob) + + # Sea Surface Temperature + obsspace.create_var('ObsValue/seaSurfaceTemperature', dtype=sst1.dtype, + fillval=sst1.fill_value) \ + .write_attr('units', 'K') \ + .write_attr('long_name', 'Sea Surface Temperature') \ + .write_data(sst1) + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Running time for splitting and output IODA: \ + {running_time} seconds") + + logger.debug("All Done!") + + +if __name__ == '__main__': + + start_time = time.time() + + parser = argparse.ArgumentParser() + parser.add_argument('-c', '--config', type=str, + help='Input JSON configuration', + required=True) + parser.add_argument('-v', '--verbose', + help='print debug logging information', + action='store_true') + args = parser.parse_args() + + log_level = 'DEBUG' if args.verbose else 'INFO' + logger = Logger('bufr2ioda_sfcshp_prepbufr.py', level=log_level, + colored_log=True) + + with open(args.config, "r") as json_file: + config = json.load(json_file) + + bufr_to_ioda(config, logger) + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Total running time: {running_time} seconds") diff --git a/ush/ioda/bufr2ioda/run_bufr2ioda.py b/ush/ioda/bufr2ioda/run_bufr2ioda.py index 63071d838..da01dd131 100755 --- a/ush/ioda/bufr2ioda/run_bufr2ioda.py +++ b/ush/ioda/bufr2ioda/run_bufr2ioda.py @@ -30,7 +30,8 @@ def bufr2ioda(current_cycle, RUN, DMPDIR, config_template_dir, COM_OBS): } # Specify observation types to be processed by a script - BUFR_py = ["satwind_amv_goes", "satwind_scat", "adpupa_prepbufr", "adpsfc_prepbufr"] + + BUFR_py = ["satwind_amv_goes", "satwind_scat", "adpupa_prepbufr", "adpsfc_prepbufr", "sfcshp_prepbufr"] for obtype in BUFR_py: logger.info(f"Convert {obtype}...") From a88e6e9ad1f3bdfea1aa8e23a7252f2248060f9c Mon Sep 17 00:00:00 2001 From: Shastri Paturi Date: Fri, 3 Nov 2023 15:47:59 +0000 Subject: [PATCH 6/9] insitu Obsyamls for monthly tank data (#704) * insitu yamls for subpfl:argo&glider;tesac & bathy (monthly tank) * updated obs_list.yaml with in situ yamls * removed extra comment lines in obs_list.yaml --- parm/soca/obs/config/salt_profile_argo.yaml | 31 +++++++++++++++++++ parm/soca/obs/config/salt_profile_glider.yaml | 31 +++++++++++++++++++ parm/soca/obs/config/salt_profile_tesac.yaml | 31 +++++++++++++++++++ parm/soca/obs/config/temp_profile_argo.yaml | 28 +++++++++++++++++ parm/soca/obs/config/temp_profile_bathy.yaml | 28 +++++++++++++++++ parm/soca/obs/config/temp_profile_glider.yaml | 28 +++++++++++++++++ parm/soca/obs/config/temp_profile_tesac.yaml | 28 +++++++++++++++++ parm/soca/obs/obs_list.yaml | 7 +++++ 8 files changed, 212 insertions(+) create mode 100644 parm/soca/obs/config/salt_profile_argo.yaml create mode 100644 parm/soca/obs/config/salt_profile_glider.yaml create mode 100644 parm/soca/obs/config/salt_profile_tesac.yaml create mode 100644 parm/soca/obs/config/temp_profile_argo.yaml create mode 100644 parm/soca/obs/config/temp_profile_bathy.yaml create mode 100644 parm/soca/obs/config/temp_profile_glider.yaml create mode 100644 parm/soca/obs/config/temp_profile_tesac.yaml diff --git a/parm/soca/obs/config/salt_profile_argo.yaml b/parm/soca/obs/config/salt_profile_argo.yaml new file mode 100644 index 000000000..9f3e2521a --- /dev/null +++ b/parm/soca/obs/config/salt_profile_argo.yaml @@ -0,0 +1,31 @@ +obs space: + name: salt_profile_argo + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}salt_profile_argo.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/salt_profile_argo.${PDY}${cyc}.nc4 + simulated variables: [salinity] + io pool: + max pool size: 1 +obs operator: + name: VertInterp + observation alias file: ./obsop_name_map.yaml + vertical coordinate: sea_water_depth + observation vertical coordinate: depth + interpolation method: linear +obs error: + covariance model: diagonal + +obs filter: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/salt_profile_glider.yaml b/parm/soca/obs/config/salt_profile_glider.yaml new file mode 100644 index 000000000..a1504886b --- /dev/null +++ b/parm/soca/obs/config/salt_profile_glider.yaml @@ -0,0 +1,31 @@ +obs space: + name: salt_profile_glider + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}salt_profile_glider.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/salt_profile_glider.${PDY}${cyc}.nc4 + simulated variables: [salinity] + io pool: + max pool size: 1 +obs operator: + name: VertInterp + observation alias file: ./obsop_name_map.yaml + vertical coordinate: sea_water_depth + observation vertical coordinate: depth + interpolation method: linear +obs error: + covariance model: diagonal + +obs filter: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/salt_profile_tesac.yaml b/parm/soca/obs/config/salt_profile_tesac.yaml new file mode 100644 index 000000000..53cc4496e --- /dev/null +++ b/parm/soca/obs/config/salt_profile_tesac.yaml @@ -0,0 +1,31 @@ +obs space: + name: salt_profile_tesac + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}salt_profile_tesac.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/salt_profile_tesac.${PDY}${cyc}.nc4 + simulated variables: [salinity] + io pool: + max pool size: 1 +obs operator: + name: VertInterp + observation alias file: ./obsop_name_map.yaml + vertical coordinate: sea_water_depth + observation vertical coordinate: depth + interpolation method: linear +obs error: + covariance model: diagonal + +obs filter: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/temp_profile_argo.yaml b/parm/soca/obs/config/temp_profile_argo.yaml new file mode 100644 index 000000000..ff1b8ce9d --- /dev/null +++ b/parm/soca/obs/config/temp_profile_argo.yaml @@ -0,0 +1,28 @@ +obs space: + name: temp_profile_argo + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}temp_profile_argo.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/temp_profile_argo.${PDY}${cyc}.nc4 + simulated variables: [waterTemperature] + io pool: + max pool size: 1 +obs operator: + name: InsituTemperature + +obs error: + covariance model: diagonal + +obs filters: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/temp_profile_bathy.yaml b/parm/soca/obs/config/temp_profile_bathy.yaml new file mode 100644 index 000000000..90c0569ea --- /dev/null +++ b/parm/soca/obs/config/temp_profile_bathy.yaml @@ -0,0 +1,28 @@ +obs space: + name: temp_profile_bathy + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}temp_profile_bathy.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/temp_profile_bathy.${PDY}${cyc}.nc4 + simulated variables: [waterTemperature] + io pool: + max pool size: 1 +obs operator: + name: InsituTemperature + +obs error: + covariance model: diagonal + +obs filters: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/temp_profile_glider.yaml b/parm/soca/obs/config/temp_profile_glider.yaml new file mode 100644 index 000000000..721fccafb --- /dev/null +++ b/parm/soca/obs/config/temp_profile_glider.yaml @@ -0,0 +1,28 @@ +obs space: + name: temp_profile_glider + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}temp_profile_glider.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/temp_profile_glider.${PDY}${cyc}.nc4 + simulated variables: [waterTemperature] + io pool: + max pool size: 1 +obs operator: + name: InsituTemperature + +obs error: + covariance model: diagonal + +obs filters: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/config/temp_profile_tesac.yaml b/parm/soca/obs/config/temp_profile_tesac.yaml new file mode 100644 index 000000000..8cbe84bd6 --- /dev/null +++ b/parm/soca/obs/config/temp_profile_tesac.yaml @@ -0,0 +1,28 @@ +obs space: + name: temp_profile_tesac + obsdatain: + engine: + type: H5File + obsfile: !ENV ${DATA}/obs/${OPREFIX}temp_profile_tesac.${PDY}${cyc}.nc4 + obsdataout: + engine: + type: H5File + obsfile: !ENV ${DATA}/diags/temp_profile_tesac.${PDY}${cyc}.nc4 + simulated variables: [waterTemperature] + io pool: + max pool size: 1 +obs operator: + name: InsituTemperature + +obs error: + covariance model: diagonal + +obs filters: + # Passivate obs where ocean fraction is > 90% + - filter: Domain Check + action: + name: passivate + where: + - variable: {name: GeoVaLs/sea_area_fraction} + maxvalue: 0.9 + diff --git a/parm/soca/obs/obs_list.yaml b/parm/soca/obs/obs_list.yaml index d710f3461..89b2659a1 100644 --- a/parm/soca/obs/obs_list.yaml +++ b/parm/soca/obs/obs_list.yaml @@ -13,3 +13,10 @@ observers: #- !INC ${OBS_YAML_DIR}/icec_ssmis_f17_south.yaml #- !INC ${OBS_YAML_DIR}/icec_ssmis_f18_north.yaml - !INC ${OBS_YAML_DIR}/icec_ssmis_f18_south.yaml +#- !INC ${OBS_YAML_DIR}/temp_profile_argo.yaml +#- !INC ${OBS_YAML_DIR}/temp_profile_bathy.yaml +#- !INC ${OBS_YAML_DIR}/temp_profile_glider.yaml +#- !INC ${OBS_YAML_DIR}/temp_profile_tesac.yaml +#- !INC ${OBS_YAML_DIR}/salt_profile_argo.yaml +#- !INC ${OBS_YAML_DIR}/salt_profile_glider.yaml +#- !INC ${OBS_YAML_DIR}/salt_profile_tesac.yaml \ No newline at end of file From 645e644e45be9b004a793ed98483d6487c159433 Mon Sep 17 00:00:00 2001 From: AndrewEichmann-NOAA <58948505+AndrewEichmann-NOAA@users.noreply.github.com> Date: Fri, 3 Nov 2023 15:41:09 -0400 Subject: [PATCH 7/9] tell JGLOBAL_PREP_OCEAN_OBS ctest location of ex-script in GDASApp (#706) --- test/soca/gw/run_jjobs.yaml.test | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/test/soca/gw/run_jjobs.yaml.test b/test/soca/gw/run_jjobs.yaml.test index f773a8c79..f908c841f 100644 --- a/test/soca/gw/run_jjobs.yaml.test +++ b/test/soca/gw/run_jjobs.yaml.test @@ -33,6 +33,9 @@ gw environement: OOPS_DEBUG: 1 OMP_NUM_THREADS: 1 + run scripts: + GDASPREPOCNOBSPY: @HOMEgfs@/sorc/gdas.cd/scripts/exglobal_prep_ocean_obs.py + setup_expt config: base: DO_JEDIATMVAR: "NO" @@ -51,6 +54,8 @@ setup_expt config: SABER_BLOCKS_YAML: @HOMEgfs@/sorc/gdas.cd/parm/soca/berror/saber_blocks.yaml NICAS_RESOL: 1 NICAS_GRID_SIZE: 150 + prepoceanobs: + SOCA_OBS_LIST: @HOMEgfs@/sorc/gdas.cd/parm/soca/obs/obs_list_small.yaml job options: account: da-cpu From c501fed4ec0607ddec36854705628ae5d0096bcc Mon Sep 17 00:00:00 2001 From: Jiarui Dong Date: Fri, 3 Nov 2023 16:06:19 -0400 Subject: [PATCH 8/9] Replace "MetaData/height" with "MetaData/stationElevation" (#707) Co-authored-by: Cory Martin --- test/land/letkfoi_land.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/land/letkfoi_land.yaml b/test/land/letkfoi_land.yaml index c9069f0e8..98cd3d171 100644 --- a/test/land/letkfoi_land.yaml +++ b/test/land/letkfoi_land.yaml @@ -64,7 +64,7 @@ observations: where: - minvalue: '-999.0' variable: - name: MetaData/height + name: MetaData/stationElevation - filter: Domain Check where: - maxvalue: '1.5' From 0d600e9c412db0f4c16849a8df89eb915f1caddc Mon Sep 17 00:00:00 2001 From: Andrew Collard <40322596+ADCollard@users.noreply.github.com> Date: Mon, 6 Nov 2023 08:35:17 -0500 Subject: [PATCH 9/9] Add error compare plot (#709) --- ush/eva/jedi_gsi_compare_conv.yaml | 55 ++++++++++++++++++++++++++++++ 1 file changed, 55 insertions(+) diff --git a/ush/eva/jedi_gsi_compare_conv.yaml b/ush/eva/jedi_gsi_compare_conv.yaml index 02293a7b0..c82a0dbd9 100644 --- a/ush/eva/jedi_gsi_compare_conv.yaml +++ b/ush/eva/jedi_gsi_compare_conv.yaml @@ -107,6 +107,24 @@ transforms: for: variable: *variables + # Obs Error that passed QC for JEDI + - transform: accept where + new name: experiment::EffectiveErrorQc::${variable} + starting field: experiment::EffectiveError::${variable} + where: + - experiment::EffectiveQC::${variable} == 0 + for: + variable: *variables + + # ObsError that passed QC for GSI + - transform: accept where + new name: experiment::GsiFinalObsErrorQc::${variable} + starting field: experiment::GsiFinalObsError::${variable} + where: + - experiment::GsiEffectiveQC::${variable} == 0 + for: + variable: *variables + graphics: # ---------- Scatter Plots ---------- @@ -276,6 +294,7 @@ graphics: markersize: 5 color: 'red' label: 'hofxdiff vs pressure' + do_linear_regression: False statistics: fields: - field_name: experiment::HofXDiff::${variable} @@ -354,6 +373,42 @@ graphics: markersize: 5 color: 'red' label: 'errordiff vs pressure' + do_linear_regression: False + # Error comparison as a function of pressure + - batch figure: + variables: *variables + @CHANNELSKEY@ + figure: + layout: [1,1] + title: 'JEDI-GSI Comparison (Passed QC) | @NAME@ @CYCLE@ | ${variable_title}' + output name: observation_scatter_plots/errorcomp_vs_pressure_@CYCLE@_@NAME@_${variable}@CHANNELVAR@.png + plots: + - add_xlabel: 'Observation Error' + add_ylabel: 'Observation Pressure' + add_grid: + add_legend: + loc: 'lower left' + layers: + - type: Scatter + x: + variable: experiment::GsiFinalObsErrorQc::${variable} + y: + variable: experiment::MetaData::pressure + @CHANNELKEY@ + markersize: 3 + color: 'red' + label: 'GSI Error' + do_linear_regression: False + - type: Scatter + x: + variable: experiment::EffectiveErrorQc::${variable} + y: + variable: experiment::MetaData::pressure + @CHANNELKEY@ + markersize: 5 + color: 'green' + label: 'JEDI Error' + do_linear_regression: False # ---------- Histograms ---------- # Histogram of h(x) difference