diff --git a/parm/atm/obs/config/seviri_m08.yaml.j2 b/parm/atm/obs/config/seviri_m08.yaml.j2 new file mode 100644 index 000000000..e92277d66 --- /dev/null +++ b/parm/atm/obs/config/seviri_m08.yaml.j2 @@ -0,0 +1,488 @@ +- obs space: + name: seviri_m08 + obsdatain: + engine: + type: H5File + obsfile: '{{DATA}}/obs/{{OPREFIX}}seviri_m08.tm00.nc' + obsdataout: + engine: + type: H5File + obsfile: '{{DATA}}/diags/diag_seviri_m08_{{ current_cycle | to_YMDH }}.nc' + simulated variables: [brightnessTemperature] + channels: &seviri_m08_channels 4-11 + + obs operator: + name: CRTM + Absorbers: [H2O,O3,CO2] + obs options: + Sensor_ID: seviri_m08 + EndianType: little_endian + CoefficientPath: {{DATA}}/crtm/ + linear obs operator: + Absorbers: [H2O, O3] + + obs bias: + input file: '{{DATA}}/obs/{{GPREFIX}}seviri_m08.satbias.nc' + output file: '{{DATA}}/bc/{{APREFIX}}seviri_m08.satbias.nc' + variational bc: + predictors: + - name: constant + - name: lapse_rate + order: 2 + tlapse: &seviri_m08_tlapse '{{DATA}}/obs/{{GPREFIX}}seviri_m08.tlapse.txt' + - name: lapse_rate + tlapse: *seviri_m08_tlapse + - name: emissivity + - name: scan_angle + var_name: sensorScanPosition + order: 4 + - name: scan_angle + var_name: sensorScanPosition + order: 3 + - name: scan_angle + var_name: sensorScanPosition + order: 2 + - name: scan_angle + var_name: sensorScanPosition + + covariance: + minimal required obs number: 20 + variance range: [1.0e-6, 10.0] + step size: 1.0e-4 + largest analysis variance: 10000.0 + prior: + input file: '{{DATA}}/obs/{{GPREFIX}}seviri_m08.satbias_cov.nc' + inflation: + ratio: 1.1 + ratio for small dataset: 2.0 + output file: '{{DATA}}/bc/{{APREFIX}}seviri_m08.satbias_cov.nc' + + obs prior filters: + # Step 0: Assign obs error for each channel + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + action: + name: assign error + error parameter vector: [1.80, 2.50, 2.25, 1.25, 1.25, 1.25, 1.45, 1.25] + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/surfaceFlag + type: int + function: + name: IntObsFunction/Conditional + options: + defaultvalue: 4 + firstmatchingcase: true + cases: + - where: + - variable: + name: GeoVaLs/water_area_fraction + minvalue: 0.99 + value: 0 + - where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + value: 1 + - where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + value: 2 + - where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + value: 3 + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/surfaceParam + type: int + function: + name: IntObsFunction/Conditional + options: + defaultvalue: 0 + firstmatchingcase: true + cases: + - where: + - variable: + name: GeoVaLs/water_area_fraction + minvalue: 0.99 + value: 30 + - where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + value: 15 + - where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + value: 20 + - where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + value: 15 + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/thinningCriteria + type: int + function: + name: IntObsFunction/Arithmetic + options: + variables: + - name: DerivedMetaData/surfaceParam + # - name: MetaData/fractionOfClearPixelsInFOV + coefs: [1] + # intercept: 300 + + obs post filters: + # Step 1: Satellite Zenith Angle Check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + test variables: + - name: MetaData/sensorZenithAngle + maxvalue: 65. + action: + name: reject + + # Step 2: Clear-sky fraction check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + test variables: + - name: MetaData/cloudAmount + maxvalue: 30. + action: + name: reject + + # Step 3: Scene homogeneous check + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 5 + is_defined: + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 6 + is_defined: + where operator: and + + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 5 + maxvalue: 1.3 + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 6 + maxvalue: 1.3 + where operator: and + + # Step 4: Thinning + - filter: Gaussian Thinning + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + horizontal_mesh: 145 + use_reduced_horizontal_grid: true + distance_norm: maximum + time_mesh: 'PT06H' + time_min: '2021-07-31T21:00:00Z' + time_max: '2021-08-01T03:00:00Z' + priority_variable: DerivedMetaData/thinningCriteria + action: + name: reject + + + # Step 5: Surface type check + # Reject channels 5-6 over land-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: 4, 7-11 + where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + + # Reject all channels over snow-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + + # Reject all channels over ice-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + + # Reject all channelsover mixed surface + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: GeoVaLs/land_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/water_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/ice_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/surface_snow_area_fraction + maxvalue: 0.99 + max_exclusive: true + + # Step 6: Terrain Check: Do not use when height > 1km + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + where: + - variable: + name: GeoVaLs/surface_geopotential_height + maxvalue: 1000.0 + + # Step 7: Observation Range Sanity Check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + minvalue: 0.0000 + maxvalue: 1000.0 + action: + name: reject + + # Step 8: Error Inflation based on topography + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorFactorTopo + channels: *seviri_m08_channels + type: float + function: + name: ObsFunction/ObsErrorFactorTopoRad + channels: *seviri_m08_channels + options: + sensor: seviri_m08 + channels: *seviri_m08_channels + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + action: + name: inflate error + inflation variable: + name: DerivedMetaData/ObsErrorFactorTopo + channels: *seviri_m08_channels + + # Step 9: Error Inflation based on TOA transmittance + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorTransmitTopRad + channels: *seviri_m08_channels + options: + channels: *seviri_m08_channels + + # Step 10: Cloud detection check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + test variables: + - name: ObsFunction/CloudDetectMinResidualIR + channels: *seviri_m08_channels + options: + channels: *seviri_m08_channels + use_flag: [ -1, 1, 1, -1, -1, -1, -1, -1 ] + use_flag_clddet: [ -2, -2, -2, -2, -2, -2, -2, -2 ] + obserr_dtempf: [0.50, 2.00, 4.00, 2.00, 4.00] + error parameter vector: [1.80, 2.50, 2.25, 1.25, 1.25, 1.25, 1.45, 1.25] + maxvalue: 1.0e-12 + action: + name: reject + + # Step 11: Scene consistency check using channel 9 + # Reject channels 4, 6-11 if channel 9 if scene consistency is greated than 0.5 + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 9 + maxvalue: 0.5 + max_exclusive: true + + # Step 12: Gross check + # Reject channels 4, 6-11 if omf > 2 + - filter: Background Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + absolute threshold: 2.0 + action: + name: reject + + # Step 13: Error inflation for channels 3-4 based on scene consistency from channel 5 + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.5 + minvalue: 0.4 + min_exclusive: true + action: + name: inflate error + inflation factor: 1.14891 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.6 + max_exclusive: true + minvalue: 0.5 + min_exclusive: true + action: + name: inflate error + inflation factor: 1.29228 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.7 + max_exclusive: true + minvalue: 0.6 + action: + name: inflate error + inflation factor: 1.49666 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + minvalue: 0.7 + action: + name: inflate error + inflation factor: 1.51987 + + # Step 14: + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + test variables: + - name: ObsFunction/NearSSTRetCheckIR + channels: *seviri_m08_channels + options: + channels: *seviri_m08_channels + use_flag: [ -1, 1, 1, -1, -1, -1, -1, -1 ] + maxvalue: 1.0e-12 + action: + name: reject + + # Step 15: Error inflation based on surface jacobian check + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorSurfJacobianRad + channels: *seviri_m08_channels + options: + channels: *seviri_m08_channels + sensor: seviri_m08 + obserr_demisf: [0.01, 0.02, 0.03, 0.02, 0.03] + obserr_dtempf: [0.50, 2.00, 4.00, 2.00, 4.00] + + # Step 16: Cloud fraction cehck + # Reject channels 4, 6-11 if Cloud fraction (percent) > 2 + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + where: + - variable: + name: MetaData/cloudAmount + maxvalue: 2 + + # Step 17: Final gross check + - filter: Background Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m08_channels + function absolute threshold: + - name: ObsFunction/ObsErrorBoundIR + channels: *seviri_m08_channels + options: + sensor: seviri_m08 + channels: *seviri_m08_channels + obserr_bound_latitude: + name: ObsFunction/ObsErrorFactorLatRad + options: + latitude_parameters: [0.0, 0.0, 0.0, 0.0] + obserr_bound_transmittop: + name: ObsFunction/ObsErrorFactorTransmitTopRad + channels: *seviri_m08_channels + options: + channels: *seviri_m08_channels + obserr_bound_max: [ 2.0, 4.0, 3.5, 2.0, 2.0, 2.0, 2.0, 3.0 ] + error parameter vector: [1.80, 2.50, 2.25, 1.25, 1.25, 1.25, 1.45, 1.25] + action: + name: reject diff --git a/parm/atm/obs/config/seviri_m11.yaml.j2 b/parm/atm/obs/config/seviri_m11.yaml.j2 new file mode 100644 index 000000000..3b215613e --- /dev/null +++ b/parm/atm/obs/config/seviri_m11.yaml.j2 @@ -0,0 +1,486 @@ +- obs space: + name: seviri_m11 + obsdatain: + engine: + type: H5File + obsfile: '{{DATA}}/obs/{{OPREFIX}}seviri_m11.tm00.nc' + obsdataout: + engine: + type: H5File + obsfile: '{{DATA}}/diags/diag_seviri_m11_{{ current_cycle | to_YMDH }}.nc' + simulated variables: [brightnessTemperature] + channels: &seviri_m11_channels 4-11 + + obs operator: + name: CRTM + Absorbers: [H2O,O3,CO2] + obs options: + Sensor_ID: seviri_m11 + EndianType: little_endian + CoefficientPath: {{DATA}}/crtm/ + linear obs operator: + Absorbers: [H2O, O3] + + obs bias: + input file: '{{DATA}}/obs/{{GPREFIX}}seviri_m11.satbias.nc' + output file: '{{DATA}}/bc/{{APREFIX}}seviri_m11.satbias.nc' + variational bc: + predictors: + - name: constant + - name: lapse_rate + order: 2 + tlapse: &seviri_m11_tlapse '{{DATA}}/obs/{{GPREFIX}}seviri_m11.tlapse.txt' + - name: lapse_rate + tlapse: *seviri_m11_tlapse + - name: emissivity + - name: scan_angle + var_name: sensorScanPosition + order: 4 + - name: scan_angle + var_name: sensorScanPosition + order: 3 + - name: scan_angle + var_name: sensorScanPosition + order: 2 + - name: scan_angle + var_name: sensorScanPosition + + covariance: + minimal required obs number: 20 + variance range: [1.0e-6, 10.0] + step size: 1.0e-4 + largest analysis variance: 10000.0 + prior: + input file: '{{DATA}}/obs/{{GPREFIX}}seviri_m11.satbias_cov.nc' + inflation: + ratio: 1.1 + ratio for small dataset: 2.0 + output file: '{{DATA}}/bc/{{APREFIX}}seviri_m11.satbias_cov.nc' + + obs prior filters: + # Step 0: Assign obs error for each channel + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + action: + name: assign error + error parameter vector: [0.75, 2.50, 2.25, 1.25, 1.25, 0.75, 0.80, 1.25] + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/surfaceFlag + type: int + function: + name: IntObsFunction/Conditional + options: + defaultvalue: 4 + firstmatchingcase: true + cases: + - where: + - variable: + name: GeoVaLs/water_area_fraction + minvalue: 0.99 + value: 0 + - where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + value: 1 + - where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + value: 2 + - where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + value: 3 + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/surfaceParam + type: int + function: + name: IntObsFunction/Conditional + options: + defaultvalue: 0 + firstmatchingcase: true + cases: + - where: + - variable: + name: GeoVaLs/water_area_fraction + minvalue: 0.99 + value: 30 + - where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + value: 15 + - where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + value: 20 + - where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + value: 15 + + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/thinningCriteria + type: int + function: + name: IntObsFunction/Arithmetic + options: + variables: + - name: DerivedMetaData/surfaceParam + # - name: MetaData/fractionOfClearPixelsInFOV + coefs: [1] + # intercept: 300 + + obs post filters: + # Step 1: Satellite Zenith Angle Check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + test variables: + - name: MetaData/sensorZenithAngle + maxvalue: 65. + action: + name: reject + + # Step 2: Clear-sky fraction check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + test variables: + - name: MetaData/cloudAmount + maxvalue: 30. + action: + name: reject + + # Step 3: Scene homogeneous check + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 5 + is_defined: + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 6 + is_defined: + where operator: and + + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 5 + maxvalue: 1.3 + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 6 + maxvalue: 1.3 + where operator: and + + ## Step 4: Thinning + - filter: Gaussian Thinning + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + horizontal_mesh: 145 + use_reduced_horizontal_grid: true + distance_norm: maximum + time_mesh: 'PT06H' + time_min: '2021-07-31T21:00:00Z' + time_max: '2021-08-01T03:00:00Z' + action: + name: reject + + # Step 5: Surface type check + # Reject channels 5-6 over land-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: 4, 7-11 + where: + - variable: + name: GeoVaLs/land_area_fraction + minvalue: 0.99 + + # Reject all channels over snow-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: GeoVaLs/surface_snow_area_fraction + minvalue: 0.99 + + # Reject all channels over ice-dominant area + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: GeoVaLs/ice_area_fraction + minvalue: 0.99 + + # Reject all channelsover mixed surface + - filter: RejectList + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: GeoVaLs/land_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/water_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/ice_area_fraction + maxvalue: 0.99 + max_exclusive: true + - variable: + name: GeoVaLs/surface_snow_area_fraction + maxvalue: 0.99 + max_exclusive: true + + # Step 6: Terrain Check: Do not use when height > 1km + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + where: + - variable: + name: GeoVaLs/surface_geopotential_height + maxvalue: 1000.0 + + # Step 7: Observation Range Sanity Check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + minvalue: 0.0000 + maxvalue: 1000.0 + action: + name: reject + + # Step 8: Error Inflation based on topography + - filter: Variable Assignment + assignments: + - name: DerivedMetaData/ObsErrorFactorTopo + channels: *seviri_m11_channels + type: float + function: + name: ObsFunction/ObsErrorFactorTopoRad + channels: *seviri_m11_channels + options: + sensor: seviri_m08 + channels: *seviri_m11_channels + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + action: + name: inflate error + inflation variable: + name: DerivedMetaData/ObsErrorFactorTopo + channels: *seviri_m11_channels + + # Step 9: Error Inflation based on TOA transmittance + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorTransmitTopRad + channels: *seviri_m11_channels + options: + channels: *seviri_m11_channels + + # Step 10: Cloud detection check + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + test variables: + - name: ObsFunction/CloudDetectMinResidualIR + channels: *seviri_m11_channels + options: + channels: *seviri_m11_channels + use_flag: [ -1, 1, 1, -1, -1, -1, -1, -1 ] + use_flag_clddet: [ -2, -2, -2, -2, -2, -2, -2, -2 ] + obserr_dtempf: [0.50, 2.00, 4.00, 2.00, 4.00] + error parameter vector: [0.75, 2.50, 2.25, 1.25, 1.25, 0.75, 0.80, 1.25] + maxvalue: 1.0e-12 + action: + name: reject + + # Step 11: Scene consistency check using channel 9 + # Reject channels 4, 6-11 if channel 9 if scene consistency is greated than 0.5 + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature + channels: 9 + maxvalue: 0.5 + max_exclusive: true + + # Step 12: Gross check + # Reject channels 4, 6-11 if omf > 2 + - filter: Background Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + absolute threshold: 2.0 + action: + name: reject + + # Step 13: Error inflation for channels 3-4 based on scene consistency from channel 5 + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.5 + minvalue: 0.4 + min_exclusive: true + action: + name: inflate error + inflation factor: 1.14891 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.6 + max_exclusive: true + minvalue: 0.5 + min_exclusive: true + action: + name: inflate error + inflation factor: 1.29228 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + maxvalue: 0.7 + max_exclusive: true + minvalue: 0.6 + action: + name: inflate error + inflation factor: 1.49666 + + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: 6-7 + where: + - variable: + name: ClearSkyStdDev/brightnessTemperature_5 + minvalue: 0.7 + action: + name: inflate error + inflation factor: 1.51987 + + # Step 14: + - filter: Bounds Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + test variables: + - name: ObsFunction/NearSSTRetCheckIR + channels: *seviri_m11_channels + options: + channels: *seviri_m11_channels + use_flag: [ -1, 1, 1, -1, -1, -1, -1, -1 ] + maxvalue: 1.0e-12 + action: + name: reject + + # Step 15: Error inflation based on surface jacobian check + - filter: Perform Action + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + action: + name: inflate error + inflation variable: + name: ObsFunction/ObsErrorFactorSurfJacobianRad + channels: *seviri_m11_channels + options: + channels: *seviri_m11_channels + sensor: seviri_m11 + obserr_demisf: [0.01, 0.02, 0.03, 0.02, 0.03] + obserr_dtempf: [0.50, 2.00, 4.00, 2.00, 4.00] + + # Step 16: Cloud fraction cehck + # Reject channels 4, 6-11 if Cloud fraction (percent) > 2 + - filter: Domain Check + filter variables: + - name: brightnessTemperature + channels: 4, 6-11 + where: + - variable: + name: MetaData/cloudAmount + maxvalue: 2 + + # Step 17: Final gross check + - filter: Background Check + filter variables: + - name: brightnessTemperature + channels: *seviri_m11_channels + function absolute threshold: + - name: ObsFunction/ObsErrorBoundIR + channels: *seviri_m11_channels + options: + sensor: seviri_m08 + channels: *seviri_m11_channels + obserr_bound_latitude: + name: ObsFunction/ObsErrorFactorLatRad + options: + latitude_parameters: [0.0, 0.0, 0.0, 0.0] + obserr_bound_transmittop: + name: ObsFunction/ObsErrorFactorTransmitTopRad + channels: *seviri_m11_channels + options: + channels: *seviri_m11_channels + obserr_bound_max: [ 2.0, 4.0, 3.5, 2.0, 2.0, 2.0, 2.0, 3.0 ] + error parameter vector: [0.75, 2.50, 2.25, 1.25, 1.25, 0.75, 0.80, 1.25] + action: + name: reject diff --git a/parm/atm/obs/lists/gdas_prototype_3d.yaml.j2 b/parm/atm/obs/lists/gdas_prototype_3d.yaml.j2 index f576fc25e..11d43e2d0 100644 --- a/parm/atm/obs/lists/gdas_prototype_3d.yaml.j2 +++ b/parm/atm/obs/lists/gdas_prototype_3d.yaml.j2 @@ -7,4 +7,6 @@ observers: {% include 'atm/obs/config/ompsnp_npp.yaml.j2' %} {% include 'atm/obs/config/ompstc_npp.yaml.j2' %} {% include 'atm/obs/config/omi_aura.yaml.j2' %} +{% include 'atm/obs/config/seviri_m08.yaml.j2' %} +{% include 'atm/obs/config/seviri_m11.yaml.j2' %} {% endfilter %} diff --git a/parm/ioda/bufr2ioda/bufr2ioda_sevcsr.json b/parm/ioda/bufr2ioda/bufr2ioda_sevcsr.json new file mode 100644 index 000000000..15cc0ee63 --- /dev/null +++ b/parm/ioda/bufr2ioda/bufr2ioda_sevcsr.json @@ -0,0 +1,16 @@ +{ + "data_format" : "bufr_d", + "data_type" : "sevcsr", + "cycle_type" : "{{ RUN }}", + "cycle_datetime" : "{{ current_cycle | to_YMDH }}", + "dump_directory" : "{{ DMPDIR }}", + "ioda_directory" : "{{ COM_OBS }}", + "subsets" : [ ], + "data_description" : "NC021043 SEVIRI, Meteosat-11; NC021043 SEVIRI, Meteosat-10, NC021043 SEVIRI, Meteosat-08", + "data_provider" : "U.S. NOAA/NESDIS", + "sensor_info" : { "sensor_name": "SEVIRI", "sensor_full_name": "Spinning Enhanced Visible Infra-Red Imager", "sensor_id": 207 }, + "satellite_info" : [ + { "satellite_name": "m11", "satellite_full_name": "Meteosat-11", "satellite_id": 70, "launch time": "20150715" }, + { "satellite_name": "m10", "satellite_full_name": "Meteosat-10", "satellite_id": 57, "launch time": "20120705" }, + { "satellite_name": "m08", "satellite_full_name": "Meteosat-08", "satellite_id": 55, "launch time": "20020828" } ] +} diff --git a/ush/ioda/bufr2ioda/bufr2ioda_sevcsr.py b/ush/ioda/bufr2ioda/bufr2ioda_sevcsr.py new file mode 100755 index 000000000..b1a4e3274 --- /dev/null +++ b/ush/ioda/bufr2ioda/bufr2ioda_sevcsr.py @@ -0,0 +1,560 @@ +#!/usr/bin/env python3 +import argparse +import calendar +import datetime +import json +import math +import os +import time +from datetime import datetime + +import numpy as np +import numpy.ma as ma +from wxflow import Logger + +from pyioda import ioda_obs_space as ioda_ospace +from pyiodaconv import bufr + +# Define and initialize global variables +global float32_fill_value +global int32_fill_value +global int64_fill_value + +float32_fill_value = np.float32(0) +int32_fill_value = np.int32(0) +int64_fill_value = np.int64(0) + + +def bufr_to_ioda(config, logger): + subsets = config["subsets"] + logger.debug(f"Checking subsets = {subsets}") + + # Get parameters from configuration + subsets = config["subsets"] + data_format = config["data_format"] + data_type = config["data_type"] + data_description = config["data_description"] + data_provider = config["data_provider"] + cycle_type = config["cycle_type"] + dump_dir = config["dump_directory"] + ioda_dir = config["ioda_directory"] + cycle = config["cycle_datetime"] + yyyymmdd = cycle[0:8] + hh = cycle[8:10] + + satellite_info_array = config["satellite_info"] + sensor_name = config["sensor_info"]["sensor_name"] + sensor_full_name = config["sensor_info"]["sensor_full_name"] + sensor_id = config["sensor_info"]["sensor_id"] + + # 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") + + # General informaton + converter = "BUFR to IODA Converter" + process_level = "Level-2" + platform_description = "NOAA Series of Geostationary Operational Environmental Satellites - 3rd generation since 2016" + sensor_description = "Spinning Enhanced Visible and InfraRed Imager;12 channels, 1 narrow-bandwidth, 1 high-resolution broad-bandwidth VIS" + + logger.info(f"sensor_name = {sensor_name}") + logger.info(f"sensor_full_name = {sensor_full_name}") + logger.info(f"sensor_id = {sensor_id}") + logger.info(f"reference_time = {reference_time}") + + bufrfile = f"{cycle_type}.t{hh}z.{data_type}.tm00.{data_format}" + DATA_PATH = os.path.join( + dump_dir, f"{cycle_type}.{yyyymmdd}", str(hh), "atmos", bufrfile + ) + if not os.path.isfile(DATA_PATH): + logger.info(f"The DATA_PATH is: {DATA_PATH}") + + # ============================================ + # Make the QuerySet for all the data we want + # ============================================ + start_time = time.time() + + logger.info("Making QuerySet") + q = bufr.QuerySet(subsets) + + # MetaData + q.add("latitude", "*/CLATH") + q.add("longitude", "*/CLONH") + q.add("satelliteId", "*/SAID") + q.add("year", "*/YEAR") + q.add("month", "*/MNTH") + q.add("day", "*/DAYS") + q.add("hour", "*/HOUR") + q.add("minute", "*/MINU") + q.add("second", "*/SECO") + q.add("sensorId", "*/SIID[1]") + q.add("sensorZenithAngle", "*/SAZA") + q.add("sensorCentralFrequency", "*/RPSEQ7/SCCF") + q.add("solarZenithAngle", "*/SOZA") + q.add("cloudFree", "*/RPSEQ7{5}/NCLDMNT") + q.add("brightnessTemperature", "*/RPSEQ7/TMBRST") + q.add("ClearSkyStdDev", "*/RPSEQ7/SDTB") + + end_time = time.time() + running_time = end_time - start_time + logger.debug(f"Processing 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.info("Executing QuerySet to get ResultSet") + with bufr.File(DATA_PATH) as f: + try: + r = f.execute(q) + except Exception as err: + logger.info(f'Return with {err}') + return + # MetaData + satid = r.get("satelliteId") + instid = r.get("sensorId") + year = r.get("year") + month = r.get("month") + day = r.get("day") + hour = r.get("hour") + minute = r.get("minute") + second = r.get("second") + lat = r.get("latitude") + lon = r.get("longitude") + satzenang = r.get("sensorZenithAngle") + chanfreq = r.get("sensorCentralFrequency", type="float") + BT = r.get("brightnessTemperature") + clrStdDev = r.get("ClearSkyStdDev") + cldFree = r.get("cloudFree", type="float") + solzenang = r.get("solarZenithAngle") + # DateTime: seconds since Epoch time + # IODA has no support for numpy datetime arrays dtype=datetime64[s] + timestamp = r.get_datetime( + "year", "month", "day", "hour", "minute", "second" + ).astype(np.int64) + + # Global variables declaration + # Set global fill values + float32_fill_value = satzenang.fill_value + int32_fill_value = satid.fill_value + int64_fill_value = timestamp.fill_value.astype(np.int64) + + end_time = time.time() + running_time = end_time - start_time + logger.info( + f"Processing time for executing QuerySet to get ResultSet : {running_time} seconds" + ) + + # ========================= + # Create derived variables + # ========================= + start_time = time.time() + + nfov = satzenang.shape[0] + scanpos = np.zeros(nfov, dtype=np.int32) + scanpos = satzenang.astype(np.int32) + 1 + + cloudAmount = 100. - cldFree + + sataziang = np.full_like(solzenang, float32_fill_value, dtype=np.float32) + solaziang = np.full_like(solzenang, float32_fill_value, dtype=np.float32) + viewang = np.full_like(solzenang, float32_fill_value, dtype=np.float32) + # Define Channel dimension for channels 4 to 11 since the other channel values are missing + channel_start = 4 + channel_end = 11 + channum = np.arange(channel_start, channel_end + 1) + + # Define wavenumbers for each satellite ID + wavenum_values_dict = { + 70: np.array( + [ + 257301, + 159467.6, + 136116.2, + 114714.5, + 103474.5, + 92997.55, + 83834.94, + 74750.95, + ], + dtype=np.float32, + ), + 55: np.array( + [ + 256550.4, + 159490.2, + 136128.6, + 114870.3, + 103420.4, + 92938.72, + 83886.68, + 75122.19, + ], + dtype=np.float32, + ), + } + wavenum_fill_value = float32_fill_value + + logger.info("Creating derived variables") + + end_time = time.time() + running_time = end_time - start_time + logger.info( + f"Processing time for creating derived variables : {running_time} seconds" + ) + + # ===================================== + # Split output based on satellite id + # Create IODA ObsSpace + # Write IODA output + # ===================================== + logger.info("Create IODA ObsSpace and Write IODA output based on satellite ID") + + # Find nique satellite identifiers in data to process + unique_satids = np.unique(satid) + logger.info(f"Number of Unique satellite identifiers: {len(unique_satids)}") + logger.info(f"Unique satellite identifiers: {unique_satids}") + logger.debug(f"Loop through unique satellite identifier {unique_satids}") + total_ob_processed = 0 + for sat in unique_satids.tolist(): + start_time = time.time() + + matched = False + for satellite_info in satellite_info_array: + if satellite_info["satellite_id"] == sat: + matched = True + satellite_id = satellite_info["satellite_id"] + satellite_name = satellite_info["satellite_name"] + satinst = sensor_name.lower() + "_" + satellite_name.lower() + logger.debug(f"Split data for {satinst} satid = {sat}") + + if matched: + if satellite_id in wavenum_values_dict: + # Extract the wavenum values for the current satellite ID + Wavenum = wavenum_values_dict[satellite_id] + else: + # If the satellite ID is not in the dictionary + logger.debug(f"satellite ID is not in the dictionary {satellite_id}") + + # Define a boolean mask to subset data from the original data object + satelite_mask = satid == sat + + # Define a boolean mask based on the condition 0 < satzenang2 < 80 + satzenang_mask = np.logical_and(0 < satzenang, satzenang < 80) + + combined_mask = satzenang_mask * satelite_mask + + # MetaData + lon2 = lon[combined_mask] + lat2 = lat[combined_mask] + timestamp2 = timestamp[combined_mask] + satid2 = satid[combined_mask] + instid2 = instid[combined_mask] + satzenang2 = satzenang[combined_mask] + chanfreq2 = chanfreq[3:11] + + # Convert scanpos to np.int32 before applying the mask + scanpos2 = scanpos.astype(np.int32)[combined_mask] + # Replace masked values with the fill value before writing to IODA variable + scanpos2 = np.where(scanpos2.mask, int32_fill_value, scanpos2) + solzenang2 = solzenang[combined_mask] + cldFree2 = cldFree[combined_mask] + cloudAmount2 = cloudAmount[combined_mask] + BT2 = BT[combined_mask] + + # Extract only channels 4 to 11 + BT2 = BT2[:, 3:11] + clrStdDev2 = clrStdDev[combined_mask] + viewang2 = viewang.flatten()[combined_mask] + sataziang2 = sataziang.flatten()[combined_mask] + solaziang2 = solaziang.flatten()[combined_mask] + + # Timestamp Range + timestamp2_min = datetime.fromtimestamp(timestamp2.min()) + timestamp2_max = datetime.fromtimestamp(timestamp2.max()) + + # Check unique observation time + unique_timestamp2 = np.unique(timestamp2) + logger.debug(f"Processing output for satid {sat}") + logger.info(f"number of unique_timestamp2 {len(unique_timestamp2)}") + logger.info(f"unique_timestamp2 {unique_timestamp2}") + + # Create the dimensions + dims = { + "Location": np.arange(0, BT2.shape[0]), + "Channel": np.arange(channel_start, channel_end + 1), + } + + # Create IODA ObsSpace + iodafile = f"{cycle_type}.t{hh}z.{satinst}.tm00.nc" + OUTPUT_PATH = os.path.join(ioda_dir, iodafile) + logger.info(f"Create output file : {OUTPUT_PATH}") + obsspace = ioda_ospace.ObsSpace(OUTPUT_PATH, mode="w", dim_dict=dims) + + # Create Global attributes + logger.debug("Write global attributes") + obsspace.write_attr("Converter", converter) + obsspace.write_attr("sourceFiles", bufrfile) + obsspace.write_attr("description", data_description) + obsspace.write_attr("datetimeReference", reference_time) + obsspace.write_attr( + "datetimeRange", [str(timestamp2_min), str(timestamp2_max)] + ) + obsspace.write_attr("sensor", sensor_id) + obsspace.write_attr("platform", satellite_id) + obsspace.write_attr("platformCommonName", satellite_name) + obsspace.write_attr("sensorCommonName", sensor_name) + obsspace.write_attr("processingLevel", process_level) + obsspace.write_attr("platformLongDescription", platform_description) + obsspace.write_attr("sensorLongDescription", sensor_description) + + # Create IODA variables + logger.debug("Write variables: name, type, units, and attributes") + + # Sensor Channel Number + obsspace.create_var( + "MetaData/sensorChannelNumber", + dim_list=["Channel"], + dtype=np.int32, + fillval=int32_fill_value, + ).write_attr("long_name", "Sensor Channel Number").write_data(channum) + + # Sensor Central Frequency + obsspace.create_var( + "MetaData/sensorCentralFrequency", + dim_list=["Channel"], + dtype=chanfreq2.dtype, + fillval=chanfreq2.fill_value, + ).write_attr("units", "Hz").write_attr( + "long_name", "Satellite Channel Center Frequency" + ).write_data( + chanfreq2 + ) + + # Sensor Central Wavenumber + obsspace.create_var( + "MetaData/sensorCentralWavenumber", + dim_list=["Channel"], + dtype=Wavenum.dtype, + fillval=wavenum_fill_value, + ).write_attr("units", "m-1").write_attr( + "long_name", "Sensor Central Wavenumber" + ).write_data( + Wavenum + ) + + if np.any(combined_mask): + # Longitude + obsspace.create_var( + "MetaData/longitude", dtype=lon2.dtype, fillval=lon2.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( + lon2 + ) + + # Latitude + obsspace.create_var( + "MetaData/latitude", dtype=lat2.dtype, fillval=lat2.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( + lat2 + ) + + # Datetime + obsspace.create_var( + "MetaData/dateTime", dtype=np.int64, fillval=int64_fill_value + ).write_attr("units", "seconds since 1970-01-01T00:00:00Z").write_attr( + "long_name", "Datetime" + ).write_data( + timestamp2 + ) + + # Satellite Identifier + obsspace.create_var( + "MetaData/satelliteIdentifier", + dtype=satid2.dtype, + fillval=satid2.fill_value, + ).write_attr("long_name", "Satellite Identifier").write_data(satid2) + + # Instrument Identifier + obsspace.create_var( + "MetaData/instrumentIdentifier", + dtype=instid2.dtype, + fillval=instid2.fill_value, + ).write_attr("long_name", "Satellite Instrument Identifier").write_data( + instid2 + ) + + # Scan Position (derived variable, need to specified fill value explicitly) + obsspace.create_var( + "MetaData/sensorScanPosition", + dtype=scanpos2.astype(np.int32).dtype, + fillval=int32_fill_value, + ).write_attr("long_name", "Sensor Scan Position").write_data(scanpos2) + + # Sensor Zenith Angle + obsspace.create_var( + "MetaData/sensorZenithAngle", + dtype=satzenang2.dtype, + fillval=satzenang2.fill_value, + ).write_attr("units", "degree").write_attr( + "valid_range", np.array([0, 90], dtype=np.float32) + ).write_attr( + "long_name", "Sensor Zenith Angle" + ).write_data( + satzenang2 + ) + + # Sensor Azimuth Angle + obsspace.create_var( + "MetaData/sensorAzimuthAngle", + dtype=np.float32, + fillval=sataziang2.fill_value, + ).write_attr("units", "degree").write_attr( + "valid_range", np.array([0, 360], dtype=np.float32) + ).write_attr( + "long_name", "Sensor Azimuth Angle" + ).write_data( + sataziang2 + ) + + # Solar Azimuth Angle + obsspace.create_var( + "MetaData/solarAzimuthAngle", + dtype=np.float32, + fillval=solaziang2.fill_value, + ).write_attr("units", "degree").write_attr( + "valid_range", np.array([0, 360], dtype=np.float32) + ).write_attr( + "long_name", "Solar Azimuth Angle" + ).write_data( + solaziang2 + ) + + # Sensor View Angle + obsspace.create_var( + "MetaData/sensorViewAngle", + dtype=np.float32, + fillval=viewang2.fill_value, + ).write_attr("units", "degree").write_attr( + "long_name", "Sensor View Angle" + ).write_data( + viewang2 + ) + + # Solar Zenith Angle + obsspace.create_var( + "MetaData/solarZenithAngle", + dtype=solzenang2.dtype, + fillval=solzenang2.fill_value, + ).write_attr("units", "degree").write_attr( + "valid_range", np.array([0, 180], dtype=np.float32) + ).write_attr( + "long_name", "Solar Zenith Angle" + ).write_data( + solzenang2 + ) + + # Cloud free + obsspace.create_var( + "MetaData/cloudFree", + dtype=cldFree2.dtype, fillval=int32_fill_value + ).write_attr("units", "1").write_attr( + "valid_range", np.array([0, 100], dtype=np.int32) + ).write_attr( + "long_name", "Amount Segment Cloud Free" + ).write_data( + cldFree2 + ) + + # Cloud amount based on computation + obsspace.create_var( + "MetaData/cloudAmount", + dtype=cloudAmount2.dtype, + fillval=cloudAmount2.fill_value, + ).write_attr("units", "1").write_attr( + "valid_range", np.array([0, 100], dtype=np.float32) + ).write_attr( + "long_name", "Amount of cloud coverage in layer" + ).write_data( + cloudAmount2 + ) + + # ObsType based on computation method/spectral band + obsspace.create_var( + "ObsValue/brightnessTemperature", + dim_list=["Location", "Channel"], + dtype=np.float32, + fillval=BT2.fill_value, + ).write_attr("units", "k").write_attr( + "long_name", "Brightness Temperature" + ).write_data( + BT2 + ) + + obsspace.create_var( + "ClearSkyStdDev/brightnessTemperature", + dim_list=["Location", "Channel"], + dtype=np.float32, + fillval=clrStdDev2.fill_value, + ).write_attr( + "long_name", "Standard Deviation Brightness Temperature" + ).write_data( + clrStdDev2[:, 3:11] + ) + + else: + logger.debug( + "No valid values (0