diff --git a/internal/test_unit/xml/unit_point2grid.xml b/internal/test_unit/xml/unit_point2grid.xml index dd1792ccf0..5221ea4460 100644 --- a/internal/test_unit/xml/unit_point2grid.xml +++ b/internal/test_unit/xml/unit_point2grid.xml @@ -137,7 +137,7 @@ &OUTPUT_DIR;/point2grid/point2grid_GOES_16_AOD_TO_G212_compute.nc - + &MET_BIN;/point2grid @@ -155,7 +155,7 @@ &OUTPUT_DIR;/point2grid/point2grid_GOES_16_AOD_TO_G212_gaussian.nc - + &MET_BIN;/point2grid @@ -174,7 +174,7 @@ &OUTPUT_DIR;/point2grid/point2grid_GOES_16_ADP.nc - + &MET_BIN;/point2grid @@ -185,7 +185,7 @@ G212 \ &OUTPUT_DIR;/point2grid/point2grid_GOES_16_AOD_TO_G212_grid_map.nc \ -field 'name="AOD"; level="(*,*)";' \ - -qc 1,2,3 -method MAX \ + -qc 0,1,2 -method MAX \ -v 1 @@ -205,7 +205,7 @@ G212 \ &OUTPUT_DIR;/point2grid/point2grid_GOES_16_AOD_TO_G212.nc \ -field 'name="AOD"; level="(*,*)";' \ - -qc 1,2,3 -method MAX \ + -qc 0,1,2 -method MAX \ -v 1 @@ -213,6 +213,25 @@ + + &MET_BIN;/point2grid + + MET_TMP_DIR &OUTPUT_DIR;/point2grid + + \ + &DATA_DIR_MODEL;/goes_16/OR_ABI-L2-AODC-M6_G16_s20241100001171_e20241100003544_c20241100006242.nc \ + G212 \ + &OUTPUT_DIR;/point2grid/point2grid_GOES_16_ADP_Enterprise_high.nc \ + -field 'name="AOD_Smoke"; level="(*,*)";' \ + -adp &DATA_DIR_MODEL;/goes_16/OR_ABI-L2-ADPC-M6_G16_s20241100001171_e20241100003544_c20241100006361.nc \ + -qc 0,1 -method MAX \ + -v 1 + + + &OUTPUT_DIR;/point2grid/point2grid_GOES_16_ADP_Enterprise_high.nc + + + &MET_BIN;/point2grid diff --git a/src/libcode/vx_nc_util/nc_utils.cc b/src/libcode/vx_nc_util/nc_utils.cc index 6350b6576a..83a6f47950 100644 --- a/src/libcode/vx_nc_util/nc_utils.cc +++ b/src/libcode/vx_nc_util/nc_utils.cc @@ -537,6 +537,17 @@ bool get_nc_att_value(const NcVar *var, const ConcatString &att_name, //////////////////////////////////////////////////////////////////////// +bool get_nc_att_values(const NcVar *var, const ConcatString &att_name, + unsigned short *att_val, bool exit_on_error) { + static const char *method_name = "get_nc_att_value(NcVar,float) -> "; + bool status = get_nc_att_values_(var, att_name, att_val, exit_on_error, + (unsigned short)bad_data_int, method_name); + return status; +} + +//////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////// + bool get_nc_att_value(const NcVarAtt *att, ConcatString &att_val) { bool status = false; @@ -1729,18 +1740,21 @@ bool get_nc_data(NcVar *var, char **data) { //////////////////////////////////////////////////////////////////////// -bool get_nc_data(NcVar *var, uchar *data) { +bool get_nc_data(NcVar *var, uchar *data, bool allow_conversion) { bool return_status = false; int data_type = GET_NC_TYPE_ID_P(var); static const char *method_name = "get_nc_data(NcVar *, uchar *) -> "; if (NC_UBYTE == data_type) return_status = get_nc_data_t(var, data); - else if (NC_BYTE == data_type && has_unsigned_attribute(var)) { - int cell_count = get_data_size(var); - ncbyte *signed_data = new ncbyte[cell_count]; - return_status = get_nc_data_t(var, signed_data); - for (int idx=0; idx +bool get_nc_att_values_(const netCDF::NcVar *var, const ConcatString &att_name, + T *att_vals, bool exit_on_error, + T bad_data, const char *caller_name) { + // caller should initialize att_vals + + // + // Retrieve the NetCDF variable attribute. + // + netCDF::NcVarAtt *att = get_nc_att(var, att_name); + bool status = IS_VALID_NC_P(att); + if (status) att->getValues(att_vals); + else { + mlog << Error << "\n" << caller_name + << get_log_msg_for_att(att, GET_SAFE_NC_NAME_P(var), att_name); + } + if (att) delete att; + if (!status && exit_on_error) exit(1); + return status; +} + +//////////////////////////////////////////////////////////////////////// + template bool get_nc_att_value_(const netCDF::NcVarAtt *att, T &att_val, bool exit_on_error, T bad_data, const char *caller_name) { diff --git a/src/tools/other/point2grid/point2grid.cc b/src/tools/other/point2grid/point2grid.cc index 1861a1ca66..edf18accaa 100644 --- a/src/tools/other/point2grid/point2grid.cc +++ b/src/tools/other/point2grid/point2grid.cc @@ -124,6 +124,9 @@ static NcFile *nc_out = (NcFile *) 0; static NcDim lat_dim ; static NcDim lon_dim ; +static int adp_qc_high; /* 3 as baseline algorithm, 0 for enterpirse algorithm */ +static int adp_qc_medium; /* 1 as baseline algorithm, 1 for enterpirse algorithm */ +static int adp_qc_low; /* 0 as baseline algorithm, 2 for enterpirse algorithm */ //////////////////////////////////////////////////////////////////////// @@ -170,6 +173,7 @@ static void regrid_nc_variable(NcFile *nc_in, Met2dDataFile *fr_mtddf, static bool keep_message_type(const int mt_index); static bool has_lat_lon_vars(NcFile *nc_in); +static void set_adp_gc_values(NcVar var_adp_qc); //////////////////////////////////////////////////////////////////////// // for GOES 16 @@ -940,8 +944,6 @@ void process_point_met_data(MetPointData *met_point_obs, MetConfig &config, VarI int valid_count = 0; int absent_count = 0; int censored_count = 0; - int qc_filtered_count = 0; - int adp_qc_filtered_count = 0; float data_value; float from_min_value = 10e10; float from_max_value = -10e10; @@ -1185,7 +1187,6 @@ void process_point_file(NcFile *nc_in, MetConfig &config, VarInfo *vinfo, void process_point_python(string python_command, MetConfig &config, VarInfo *vinfo, const Grid to_grid, bool use_xarray) { - int idx, hdr_idx; ConcatString vname, vname_cnt, vname_mask; DataPlane fr_dp, to_dp; DataPlane cnt_dp, mask_dp; @@ -1193,8 +1194,6 @@ void process_point_python(string python_command, MetConfig &config, VarInfo *vin NcVar var_obs_gc, var_obs_var; clock_t start_clock = clock(); - bool has_prob_thresh = !prob_cat_thresh.check(bad_data_double); - unixtime requested_valid_time, valid_time; static const char *method_name = "process_point_python() -> "; static const char *method_name_s = "process_point_python()"; @@ -1863,6 +1862,24 @@ void check_lat_lon(int data_size, float *latitudes, float *longitudes) { << method_name << "LONG: " << min_lon << " to " << max_lon << "\n"; } +//////////////////////////////////////////////////////////////////////// +// QC flags: 0=high, 1=medium, 2=low +// Enterpise algorithm: 0=high, 1=medium, 2=low +// Baseline algorithm: 3=high, 1=medium, 0=low (high=12/48, medium=4/16) +// returns bad_data_int if it does not belong to high, mediuam, or low. + +int compute_adp_qc_flag(int adp_qc, int shift_bits) { + int particle_qc = ((adp_qc >> shift_bits) & 0x03); + int qc_for_flag; + + if (particle_qc == adp_qc_high) qc_for_flag = 0; + else if (particle_qc == adp_qc_medium) qc_for_flag = 1; + else if (particle_qc == adp_qc_low) qc_for_flag = 2; + else qc_for_flag = bad_data_int; + + return qc_for_flag; +} + //////////////////////////////////////////////////////////////////////// static unixtime compute_unixtime(NcVar *time_var, unixtime var_value) { @@ -2342,12 +2359,17 @@ int get_lon_count(NcFile *_nc) { static NcVar get_goes_nc_var(NcFile *nc, const ConcatString var_name, bool exit_if_error) { - NcVar var_data = get_nc_var(nc, var_name.c_str(), false); + NcVar var_data; + static const char *method_name = "get_goes_nc_var() -> "; + if (has_var(nc, var_name.c_str())) var_data = get_nc_var(nc, var_name.c_str(), false); if (IS_INVALID_NC(var_data)) { - var_data = get_nc_var(nc, var_name.split("_")[0].c_str()); + mlog << Debug(4) << method_name + << "The variable \"" << var_name << "\" does not exist. Find \"" + << var_name.split("_")[0] << "\" variable\n"; + var_data = get_nc_var(nc, var_name.split("_")[0].c_str()); } if (IS_INVALID_NC(var_data)) { - mlog << Error << "\nget_goes_nc_var() -> " + mlog << Error << "\n" << method_name << "The variable \"" << var_name << "\" does not exist\n\n"; if (exit_if_error) exit(1); } @@ -2424,6 +2446,7 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, bool has_qc_var = false; bool has_adp_qc_var = false; + const int log_debug_level = 4; clock_t start_clock = clock(); int to_lat_count = to_grid.ny(); int to_lon_count = to_grid.nx(); @@ -2443,6 +2466,10 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, // -99 is arbitrary number as invalid QC value memset(qc_data, -99, from_data_size*sizeof(uchar)); + adp_qc_high = 3; /* 3 as baseline algorithm, 0 for enterpirse algorithm */ + adp_qc_medium = 1; /* 1 as baseline algorithm, 1 for enterpirse algorithm */ + adp_qc_low = 0; /* 0 as baseline algorithm, 2 for enterpirse algorithm */ + NcVar var_qc; NcVar var_adp; NcVar var_adp_qc; @@ -2465,13 +2492,14 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, else if (is_smoke_only) var_adp = get_goes_nc_var(nc_adp, vname_smoke); if (IS_VALID_NC(var_adp)) { - get_nc_data(&var_adp, adp_data); + get_nc_data(&var_adp, adp_data, true); //Smoke:ancillary_variables = "DQF" ; ubyte DQF(y, x) ; if (get_att_value_string(&var_adp, (string)"ancillary_variables", qc_var_name)) { var_adp_qc = get_nc_var(nc_adp, qc_var_name.c_str()); if (IS_VALID_NC(var_adp_qc)) { get_nc_data(&var_adp_qc, adp_qc_data); + set_adp_gc_values(var_adp_qc); has_adp_qc_var = true; mlog << Debug(5) << method_name << "found QC var: " << qc_var_name << " for " << GET_NC_NAME(var_adp) << ".\n"; @@ -2521,7 +2549,6 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, int non_missing_count = 0; int qc_filtered_count = 0; int adp_qc_filtered_count = 0; - float data_value; float from_min_value = 10e10; float from_max_value = -10e10; float qc_min_value = 10e10; @@ -2534,6 +2561,18 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, missing_count = non_missing_count = 0; to_dp.set_constant(bad_data_double); + int shift_bits = 2; + if (is_dust_only) shift_bits += 2; + + int cnt_aod_qc_low = 0; + int cnt_aod_qc_high = 0; + int cnt_aod_qc_medium = 0; + int cnt_aod_qc_nr = 0; // no_retrieval_qf + int cnt_adp_qc_low = 0; + int cnt_adp_qc_high = 0; + int cnt_adp_qc_medium = 0; + int cnt_adp_qc_nr = 0; // no_retrieval_qf + int cnt_adp_qc_adjused = 0; for (int xIdx=0; xIdxcensor_thresh().n(); i++) { + // Break out after the first match. + if(vinfo->censor_thresh()[i].check(data_value)) { + data_value = vinfo->censor_val()[i]; + censored_count++; + break; + } + } + + // Check the data existance (always 1 if ADP variable does not exist) + if (0 == adp_data[from_index]) { + absent_count++; + continue; + } + // Filter by QC flag - qc_value = qc_data[from_index]; - if (!has_qc_var || !has_qc_flags || qc_flags.has(qc_value)) { - for(int i=0; icensor_thresh().n(); i++) { - // Break out after the first match. - if(vinfo->censor_thresh()[i].check(data_value)) { - data_value = vinfo->censor_val()[i]; - censored_count++; - break; + if (has_qc_var || has_adp_qc_var) { + qc_value = qc_data[from_index]; + if (mlog.verbosity_level() >= log_debug_level) { + if (qc_min_value > qc_value) qc_min_value = qc_value; + if (qc_max_value < qc_value) qc_max_value = qc_value; + switch (qc_value) { + case 0: cnt_aod_qc_high++; break; + case 1: cnt_aod_qc_medium++; break; + case 2: cnt_aod_qc_low++; break; + default: cnt_aod_qc_nr++; break; } } - if (0 == adp_data[from_index]) { - absent_count++; - continue; - } + if (has_adp_qc_var) { + int qc_for_flag = compute_adp_qc_flag(adp_qc_data[from_index], shift_bits); + bool filter_out = is_eq(qc_for_flag, bad_data_int); + + if (mlog.verbosity_level() >= log_debug_level) { + switch (qc_for_flag) { + case 0: cnt_adp_qc_high++; break; + case 1: cnt_adp_qc_medium++; break; + case 2: cnt_adp_qc_low++; break; + default: cnt_adp_qc_nr++; break; + } + } - if (has_adp_qc_var && has_qc_flags) { - int shift_bits = 2; - if (is_dust_only) shift_bits += 2; - particle_qc = ((adp_qc_data[from_index] >> shift_bits) & 0x03); - int qc_for_flag = 3 - particle_qc; // high = 3, qc_flag for high = 0 - if (!qc_flags.has(qc_for_flag)) { + if (!filter_out) { + /* Adjust the quality by AOD data QC */ + if (1 == qc_value && 0 == qc_for_flag) { + qc_for_flag = 1; /* high to medium quality */ + cnt_adp_qc_adjused++; + } + else if (2 == qc_value) { + qc_for_flag = 2; /* high/medium to low quality */ + cnt_adp_qc_adjused++; + } + if (has_qc_flags && !qc_flags.has(qc_for_flag)) filter_out = true; + } + if (filter_out) { adp_qc_filtered_count++; continue; } } - - dataArray.add(data_value); - if (mlog.verbosity_level() >= 4) { - if (qc_min_value > qc_value) qc_min_value = qc_value; - if (qc_max_value < qc_value) qc_max_value = qc_value; + else if (has_qc_var && has_qc_flags && !qc_flags.has(qc_value)) { + qc_filtered_count++; + continue; } } - else { - qc_filtered_count++; - } + dataArray.add(data_value); valid_count++; } @@ -2618,7 +2686,6 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, << data_count << " data values.\n"; } } - else {} } } @@ -2627,14 +2694,20 @@ void regrid_goes_variable(NcFile *nc_in, VarInfo *vinfo, delete [] from_data; delete [] adp_qc_data; - mlog << Debug(4) << method_name << "Count: actual: " << to_cell_count + mlog << Debug(log_debug_level) << method_name << "Count: actual: " << to_cell_count << ", missing: " << missing_count << ", non_missing: " << non_missing_count - << "\n\tFiltered: by QC: " << qc_filtered_count + << "\n Filtered: by QC: " << qc_filtered_count << ", by adp QC: " << adp_qc_filtered_count << ", by absent: " << absent_count << ", total: " << (qc_filtered_count + adp_qc_filtered_count + absent_count) - << "\n\tRange: data: [" << from_min_value << " - " << from_max_value - << "] QC: [" << qc_min_value << " - " << qc_max_value << "]\n"; + << "\n Range: data: [" << from_min_value << " - " << from_max_value + << "] QC: [" << qc_min_value << " - " << qc_max_value << "]" + << "\n AOD QC: high=" << cnt_aod_qc_high + << " medium=" << cnt_aod_qc_medium << ", low=" << cnt_aod_qc_low + << ", no_retrieval=" << cnt_aod_qc_nr + << "\n ADP QC: high=" << cnt_adp_qc_high << " medium=" << cnt_adp_qc_medium + << ", low=" << cnt_adp_qc_low << ", no_retrieval=" << cnt_adp_qc_nr + << ", adjusted=" << cnt_adp_qc_adjused << "\n"; if (to_cell_count == 0) { mlog << Warning << "\n" << method_name @@ -2835,6 +2908,38 @@ void usage() { exit(1); } +//////////////////////////////////////////////////////////////////////// +const static ConcatString att_name_values = "flag_values"; +const static ConcatString att_name_meanings = "flag_meanings"; + +void set_adp_gc_values(NcVar var_adp_qc) { + ConcatString att_flag_meanings; + + if (get_nc_att_value(&var_adp_qc, (ConcatString)"flag_meanings", att_flag_meanings)) { + StringArray flag_meanings = to_lower(att_flag_meanings).split(" "); + unsigned short flag_values[256]; + for (int i=0; i> 2) & 0x03; + } + if (flag_meanings.has("medium_confidence_smoke_detection_qf", idx)) { + adp_qc_medium = (flag_values[idx] >> 2) & 0x03; + } + if (flag_meanings.has("high_confidence_smoke_detection_qf", idx)) { + adp_qc_high = (flag_values[idx] >> 2) & 0x03; + } + } + mlog << Debug(4) << "set_adp_gc_values() " + << " high_confidence = " << adp_qc_high + << ", medium_confidence = " << adp_qc_medium + << ", low_confidence = " << adp_qc_low << "\n"; + } +} + //////////////////////////////////////////////////////////////////////// void set_field(const StringArray &a) {