From 4052380b5fb8b154b2a7c4420c9c2d2f08449f46 Mon Sep 17 00:00:00 2001 From: Scott Stewart Date: Wed, 6 Nov 2024 12:52:16 -0700 Subject: [PATCH 1/4] do not use melt_onset QA flag mask for south --- seaice_ecdr/intermediate_monthly.py | 94 +++++++++++++++-------------- 1 file changed, 49 insertions(+), 45 deletions(-) diff --git a/seaice_ecdr/intermediate_monthly.py b/seaice_ecdr/intermediate_monthly.py index f28dde23..1eb52008 100644 --- a/seaice_ecdr/intermediate_monthly.py +++ b/seaice_ecdr/intermediate_monthly.py @@ -204,7 +204,7 @@ def get_daily_ds_for_month( ) # TODO: rename. This is actually a bit mask (except the fill_value) -CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS = OrderedDict( +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH = OrderedDict( { "average_concentration_exceeds_0.15": 1, "average_concentration_exceeds_0.30": 2, @@ -217,6 +217,18 @@ def get_daily_ds_for_month( } ) +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH = OrderedDict( + { + "average_concentration_exceeds_0.15": 1, + "average_concentration_exceeds_0.30": 2, + "at_least_half_the_days_have_sea_ice_conc_exceeds_0.15": 4, + "at_least_half_the_days_have_sea_ice_conc_exceeds_0.30": 8, + "invalid_ice_mask_applied": 16, + "at_least_one_day_during_month_has_spatial_interpolation": 32, + "at_least_one_day_during_month_has_temporal_interpolation": 64, + } +) + def _qa_field_has_flag(*, qa_field: xr.DataArray, flag_value: int) -> xr.DataArray: """Returns a boolean DataArray indicating where a flag value occurs in the given `qa_field`.""" @@ -232,8 +244,14 @@ def calc_cdr_seaice_conc_monthly_qa_flag( *, daily_ds_for_month: xr.Dataset, cdr_seaice_conc_monthly: xr.DataArray, + hemisphere: Hemisphere, ) -> xr.DataArray: """Create `cdr_seaice_conc_monthly_qa_flag`.""" + if hemisphere == "north": + qa_flag_bitmasks = CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH + else: + qa_flag_bitmasks = CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH + # initialize the variable cdr_seaice_conc_monthly_qa_flag = xr.full_like( cdr_seaice_conc_monthly, @@ -246,18 +264,14 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~average_exceeds_15, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "average_concentration_exceeds_0.15" - ], + + qa_flag_bitmasks["average_concentration_exceeds_0.15"], ) average_exceeds_30 = cdr_seaice_conc_monthly > 0.30 cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~average_exceeds_30, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "average_concentration_exceeds_0.30" - ], + + qa_flag_bitmasks["average_concentration_exceeds_0.30"], ) days_in_ds = len(daily_ds_for_month.time) @@ -269,9 +283,7 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~at_least_half_have_sic_gt_15, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "at_least_half_the_days_have_sea_ice_conc_exceeds_0.15" - ], + + qa_flag_bitmasks["at_least_half_the_days_have_sea_ice_conc_exceeds_0.15"], ) at_least_half_have_sic_gt_30 = (daily_ds_for_month.cdr_seaice_conc > 0.30).sum( @@ -280,9 +292,7 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~at_least_half_have_sic_gt_30, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "at_least_half_the_days_have_sea_ice_conc_exceeds_0.30" - ], + + qa_flag_bitmasks["at_least_half_the_days_have_sea_ice_conc_exceeds_0.30"], ) # Use "invalid_ice_mask_applied", which is actually the invalid ice mask. @@ -292,8 +302,7 @@ def calc_cdr_seaice_conc_monthly_qa_flag( ).any(dim="time") cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~invalid_ice_mask_applied, - cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS["invalid_ice_mask_applied"], + cdr_seaice_conc_monthly_qa_flag + qa_flag_bitmasks["invalid_ice_mask_applied"], ) at_least_one_day_during_month_has_spatial_interpolation = _qa_field_has_flag( @@ -305,9 +314,7 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~at_least_one_day_during_month_has_spatial_interpolation, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "at_least_one_day_during_month_has_spatial_interpolation" - ], + + qa_flag_bitmasks["at_least_one_day_during_month_has_spatial_interpolation"], ) at_least_one_day_during_month_has_temporal_interpolation = _qa_field_has_flag( @@ -319,34 +326,30 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( ~at_least_one_day_during_month_has_temporal_interpolation, cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "at_least_one_day_during_month_has_temporal_interpolation" - ], - ) - - at_least_one_day_during_month_has_melt_detected = _qa_field_has_flag( - qa_field=daily_ds_for_month.cdr_seaice_conc_qa_flag, - flag_value=CDR_SEAICE_CONC_QA_FLAG_DAILY_BITMASKS["start_of_melt_detected"], - ).any(dim="time") - cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( - ~at_least_one_day_during_month_has_melt_detected, - cdr_seaice_conc_monthly_qa_flag - + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ - "at_least_one_day_during_month_has_melt_detected" - ], - ) + + qa_flag_bitmasks["at_least_one_day_during_month_has_temporal_interpolation"], + ) + + if hemisphere == "north": + at_least_one_day_during_month_has_melt_detected = _qa_field_has_flag( + qa_field=daily_ds_for_month.cdr_seaice_conc_qa_flag, + flag_value=CDR_SEAICE_CONC_QA_FLAG_DAILY_BITMASKS["start_of_melt_detected"], + ).any(dim="time") + cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.where( + ~at_least_one_day_during_month_has_melt_detected, + cdr_seaice_conc_monthly_qa_flag + + qa_flag_bitmasks["at_least_one_day_during_month_has_melt_detected"], + ) + monthly_qa_valid_range = (np.uint8(0), np.uint8(255)) + else: + monthly_qa_valid_range = (np.uint8(0), np.uint8(127)) cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.assign_attrs( long_name="NOAA/NSIDC CDR of Passive Microwave Monthly Northern Hemisphere Sea Ice Concentration QA flags", standard_name="status_flag", - flag_meanings=" ".join( - k for k in CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS.keys() - ), - flag_masks=[ - np.uint8(v) for v in CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS.values() - ], + flag_meanings=" ".join(k for k in qa_flag_bitmasks.keys()), + flag_masks=[np.uint8(v) for v in qa_flag_bitmasks.values()], grid_mapping="crs", - valid_range=(np.uint8(0), np.uint8(255)), + valid_range=monthly_qa_valid_range, ) cdr_seaice_conc_monthly_qa_flag.encoding = dict( @@ -612,6 +615,7 @@ def make_intermediate_monthly_ds( cdr_seaice_conc_monthly_qa_flag = calc_cdr_seaice_conc_monthly_qa_flag( daily_ds_for_month=daily_ds_for_month, cdr_seaice_conc_monthly=cdr_seaice_conc_monthly, + hemisphere=hemisphere, ) # Set stdev to invalid and QA to all-missing if this is an empty month @@ -663,10 +667,10 @@ def make_intermediate_monthly_ds( temporality="monthly", aggregate=False, source=", ".join([fp.item().name for fp in daily_ds_for_month.filepaths]), - # TODO: consider providing all platforms that went into month? This would be - # consistent with how we handle the aggregate filenames. Is it - # misleading to indicate that a month is a single platform when it may not - # really be? + # TODO: consider providing all platforms that went into month? + # This would be consistent with how we handle the aggregate filenames. + # Is it misleading to indicate that a month is a single platform + # when it may not really be? platform_ids=[platform_id], resolution=resolution, hemisphere=hemisphere, From 97539a589b227d21f9cfa57ade8ab508fa328ea2 Mon Sep 17 00:00:00 2001 From: Scott Stewart Date: Wed, 6 Nov 2024 12:56:55 -0700 Subject: [PATCH 2/4] fix test_monthly unit test --- seaice_ecdr/tests/unit/test_monthly.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/seaice_ecdr/tests/unit/test_monthly.py b/seaice_ecdr/tests/unit/test_monthly.py index 64400192..901f3b46 100644 --- a/seaice_ecdr/tests/unit/test_monthly.py +++ b/seaice_ecdr/tests/unit/test_monthly.py @@ -12,7 +12,7 @@ from seaice_ecdr.constants import ECDR_PRODUCT_VERSION from seaice_ecdr.intermediate_daily import get_ecdr_dir from seaice_ecdr.intermediate_monthly import ( - CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS, + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH, CDR_SEAICE_CONC_QA_FLAG_DAILY_BITMASKS, _get_daily_complete_filepaths_for_month, _platform_id_for_month, @@ -301,6 +301,9 @@ def _mock_daily_ds_for_month(): def test_calc_cdr_seaice_conc_monthly_qa_flag(): _mock_daily_ds = _mock_daily_ds_for_month() + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS = ( + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH + ) expected_flags = np.array( [ CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS[ @@ -357,6 +360,7 @@ def test_calc_cdr_seaice_conc_monthly_qa_flag(): actual = calc_cdr_seaice_conc_monthly_qa_flag( daily_ds_for_month=_mock_daily_ds, cdr_seaice_conc_monthly=_mean_daily_conc, + hemisphere="north", ) nptesting.assert_array_equal(expected_flags, actual.values) From e984265f720611d4da6b1beb0d3588490db4ef47 Mon Sep 17 00:00:00 2001 From: Trey Stafford Date: Thu, 7 Nov 2024 09:24:18 -0700 Subject: [PATCH 3/4] De-duplicate bitmask definition for SH --- seaice_ecdr/intermediate_monthly.py | 13 ++----------- 1 file changed, 2 insertions(+), 11 deletions(-) diff --git a/seaice_ecdr/intermediate_monthly.py b/seaice_ecdr/intermediate_monthly.py index 1eb52008..9a1a3ee8 100644 --- a/seaice_ecdr/intermediate_monthly.py +++ b/seaice_ecdr/intermediate_monthly.py @@ -217,17 +217,8 @@ def get_daily_ds_for_month( } ) -CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH = OrderedDict( - { - "average_concentration_exceeds_0.15": 1, - "average_concentration_exceeds_0.30": 2, - "at_least_half_the_days_have_sea_ice_conc_exceeds_0.15": 4, - "at_least_half_the_days_have_sea_ice_conc_exceeds_0.30": 8, - "invalid_ice_mask_applied": 16, - "at_least_one_day_during_month_has_spatial_interpolation": 32, - "at_least_one_day_during_month_has_temporal_interpolation": 64, - } -) +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH = CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH.copy() +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH.pop("at_least_one_day_during_month_has_melt_detected") def _qa_field_has_flag(*, qa_field: xr.DataArray, flag_value: int) -> xr.DataArray: From 53accf2b04ab9f435701bf52b61604d3d9969b63 Mon Sep 17 00:00:00 2001 From: Trey Stafford Date: Thu, 7 Nov 2024 09:26:55 -0700 Subject: [PATCH 4/4] Avoid hard-coding monthly qa range. Calculate max from bitmasks --- seaice_ecdr/intermediate_monthly.py | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/seaice_ecdr/intermediate_monthly.py b/seaice_ecdr/intermediate_monthly.py index 9a1a3ee8..69927cc6 100644 --- a/seaice_ecdr/intermediate_monthly.py +++ b/seaice_ecdr/intermediate_monthly.py @@ -217,8 +217,12 @@ def get_daily_ds_for_month( } ) -CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH = CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH.copy() -CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH.pop("at_least_one_day_during_month_has_melt_detected") +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH = ( + CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_NORTH.copy() +) +CDR_SEAICE_CONC_MONTHLY_QA_FLAG_BITMASKS_SOUTH.pop( + "at_least_one_day_during_month_has_melt_detected" +) def _qa_field_has_flag(*, qa_field: xr.DataArray, flag_value: int) -> xr.DataArray: @@ -330,17 +334,16 @@ def calc_cdr_seaice_conc_monthly_qa_flag( cdr_seaice_conc_monthly_qa_flag + qa_flag_bitmasks["at_least_one_day_during_month_has_melt_detected"], ) - monthly_qa_valid_range = (np.uint8(0), np.uint8(255)) - else: - monthly_qa_valid_range = (np.uint8(0), np.uint8(127)) + + monthly_qa_values = [np.uint8(value) for value in qa_flag_bitmasks.values()] cdr_seaice_conc_monthly_qa_flag = cdr_seaice_conc_monthly_qa_flag.assign_attrs( long_name="NOAA/NSIDC CDR of Passive Microwave Monthly Northern Hemisphere Sea Ice Concentration QA flags", standard_name="status_flag", flag_meanings=" ".join(k for k in qa_flag_bitmasks.keys()), - flag_masks=[np.uint8(v) for v in qa_flag_bitmasks.values()], + flag_masks=monthly_qa_values, grid_mapping="crs", - valid_range=monthly_qa_valid_range, + valid_range=(np.uint8(0), np.uint8(sum(monthly_qa_values))), ) cdr_seaice_conc_monthly_qa_flag.encoding = dict(