PointStat Error: "percentile threshold "==FBIAS1" used before it was set"

[malloryprow]

I have been running PointStat for the EVS global_det wave portion. I have noticed an error and it seems to be sporadic. Last week the same job failed with an error, but this morning I ran the same job and no error. I attached the logs. I have seen it too where some runs for other users have no errors, but do for others.

I don't have an example job, config, and data files that do result in an error. Dogwood is down today but hopefully by the end of the week I can get something. I just wanted to get this started to see if maybe anything was seen between the good vs. bad log files.

DEBUG 2: Computing Categorical Statistics.
ERROR  :
ERROR  : Simple_Node::check(double, double, double) const -> percentile threshold "==FBIAS1" used before it was set.
ERROR  :

BAD_PointStat_fcstGFS_obsNDBC_climoERA5_Wave_Multifield_valid2023080212_f138.log
GOOD_PointStat_fcstGFS_obsNDBC_climoERA5_Wave_Multifield_valid2023080212_f138.log

[malloryprow]

I will note that the bad config was ran on WCOSS2-Dogwood and the good one was run on WCOSS2-Cactus. I would hope this wouldn't make a difference.

[JohnHalleyGotway]

@malloryprow, thanks for sending these log files. I've been staring at them and stepping through the code. There's definitely no obvious explanation for the difference in behavior.

I do notice one difference in the environment:
In the "GOOD" run on Cactus:

08/08 12:14:15.803 metplus INFO: OMP_NUM_THREADS=1

Whereas in the "BAD" run on Dogwood:

08/03 16:47:05.426 metplus (config_util.py:115) WARNING: Config variable OMP_NUM_THREADS (1) will be overridden by the environment variable OMP_NUM_THREADS (25)
08/03 16:47:05.449 metplus INFO: OMP_NUM_THREADS=25

So the $OMP_NUM_THREADS environment variable setting trumps the config file one. The trouble is that the Point-Stat code isn't parallelized right now. So that setting really should have no effect. I just noticed the difference.

I would like to see additional log output from MET for both the GOOD and BAD runs. This DEBUG(3) log message would provide more information about the percentile computation. In the METplus config file, you'd set the LOG_MET_VERBOSITY option... might as well use 4... LOG_MET_VERBOSITY = 4.

I'm mostly worried here about memory issues, such as an uninitialized variable or incomplete clear functions. Clearly this warning message that appears in the BAD run and not the GOOD run is a symptom of the problem:

WARNING: Simple_Node::set_perc() -> the requested percentile (100) for threshold "==FBIAS1" differs from the actual percentile (0) by 100.

Are you able to rerun both using:
export OMP_NUM_THREADS=1
and
LOG_MET_VERBOSITY = 4

[malloryprow]

I can try. I'll have to wait until Dogwood is back tomorrow though.

[JohnHalleyGotway]

@malloryprow, I thought of two more requests:

1. Please turn on the MPR output line type for both runs:

POINT_STAT_OUTPUT_FLAG_MPR = STAT

2. And please send the output .stat file from the successful run. I'd like to see the actual fcst and obs values over which these percentiles are computed.

It's very interesting that the ==FBIAS1 thresholds seems to work fine for the wave height evaluation over 5 masking regions but fail on the 6th and final region, for HAWAII which consists of only 22 matched pairs. I'll note that computing the 95-th percentile or trying to de-bias a set of only 22 values likely isn't going to produce very meaningful results.

Regardless, the software should work as expected.

[malloryprow]

Attached are BAD2a, BAD2b, and GOOD2 logs that were done on WCOSS2-Dogwood this morning. These were run to produce MPR stat file output and with LOG_MET_VERBOSITY=4. Attached is the output from GOOD2 (had to change it from .stat to .txt to get it to attach).

BAD2a - run in the queue on Dogwood; runs in parallel with other PointStat jobs
BAD2b - same as BAD2a but with export OMP_NUM_THREADS=1
GOOD2 - run on the command line, attached point_stat file is from this run (end in .txt so I could attach it)

It seems maybe the problem is coming from when the job is run in parallel with other jobs??

BAD2a_PointStat_fcstGFS_obsNDBC_climoERA5_Wave_Multifield_valid2023080212_f138.log
BAD2b_PointStat_fcstGFS_obsNDBC_climoERA5_Wave_Multifield_valid2023080212_f138.log
GOOD2_PointStat_fcstGFS_obsNDBC_climoERA5_Wave_Multifield_valid2023080212_f138.log
point_stat_fcstGFS_obsNDBC_climoERA5_1380000L_20230802_120000V.txt

[JohnHalleyGotway]

Mallory, thanks for sending this data.

I don't have a good answer for you yet, but I did grep out the 22 lines of MPR data for Wave Height in the Hawaii masking region.
See WVHT_HAWAII_mpr.txt. You'll notice that there are 22 pairs but only 7 unique buoy id's in the OBS_SID column:
51000, 51002, 51003, 51004, 51101, 51202, 51213

This is caused by buoys reporting multiple times in the observation matching time window.

I'd recommend coordinating with Deanna to modify her configuration to only use the observation that is closest in time to the forecast valid time. This is controlled by the obs_summary configuration option for Point-Stat:
obs_summary = NEAREST;.

Checking in the wrapped METplus config file for Point-Stat, I see on this line that obs_summary is not currently configurable via METplus. So for now, you'd need to use:

POINT_STAT_MET_CONFIG_OVERRIDES = obs_summary = NEAREST;

Computing percentile thresholds for a list of 22 values doesn't make much sense... and doing so with only 7 makes even less. So I'd recommend removing percentile thresholds from the operational runs.

But I'll continue investigating the source of runtime diffs.

[malloryprow]

Thank you!!

[JohnHalleyGotway]

@malloryprow I'm feeling a bit suspicious of this line of code.

This is a function to compute the value of the requested percentile from an array of numbers. If the requested percentile is 100 (i.e. the max), as in the case with this sample HAWAII set of 22 points being threshold <= 3.0, then the 0-based index = 21 and index + 1 = 22. So ordered_array[index+1] is reading beyond the range of valid data. In this case, delta should be 0 and 0 times anything is 0. So it shouldn't affect the results. But this just isn't good practice. It's never good to read beyond the range of valid data.

One option would be me revising this function with more range checking, but this is a somewhat wild guess. I really wish I could replicate the problem here so that I could actually test fixes.

[malloryprow]

Would you like me to send you any data? The input files I used?

[JohnHalleyGotway]

Mallory, yes please, you can post them to our anonymous ftp site. I have MET compiled on hera using the -ftrapuv intel compiler option to initialize memory to weird values which is usually for turning up memory issue. I could try running multiple instances at the same time to see if I can trigger the error.

[malloryprow]

Done!

[malloryprow]

Let me know if you have any questions!

[JohnHalleyGotway]

@malloryprow I tried running this data through hera, running 5 instances of the same command concurrently, but I'm not able to trigger the error.

I am missing the climo files though:

/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/climos/wave/era5/wave_climo_mean_0802.nc
/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/climos/wave/era5/wave_climo_stdev_0802.nc

And the mask files:

08/08 12:14:15.801 metplus INFO: METPLUS_MASK_POLY=poly = ["/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/latlon_0p25_GLOBAL.nc"];
08/08 12:14:15.801 metplus INFO: METPLUS_MASK_SID=sid = ["/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/SEUS_CARB.ndbc.buoys", "/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/GOM.ndbc.buoys", "/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/NEUS_CAN.ndbc.buoys", "/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/WCOAST_AK.ndbc.buoys", "/lfs/h2/emc/vpppg/noscrub/emc.vpppg/verification/EVS_fix/masks/HAWAII.ndbc.buoys"]; 

If you send me those, I could try adding them to the test, but I'm not very hopeful I'll be able to replicate.

[malloryprow]

Uploaded the climo files!

[JohnHalleyGotway]

Thanks Mallory. While I wasn't able to actually replicate the error on Hera using the data you sent, @jprestop and I were able to produce the exact same error messages you're getting with one small change to the code.

DEBUG 3: The observation threshold value ">20.0" represents the 100-th percentile.
DEBUG 3: For forecast threshold "==FBIAS1" with type ">" update the requested percentile from 100 to 100.
WARNING: 
WARNING: Simple_Node::set_perc() -> the requested percentile (100) for threshold "==FBIAS1" differs from the actual percentile (0) by 100.
WARNING: This is common for small samples or data that contains ties.
WARNING: 
DEBUG 3: The requested percentile (95) for threshold ">=SOP95(2.2)" includes 95.6522% of the data.
DEBUG 2: Computing Categorical Statistics.
ERROR  : 
ERROR  : Simple_Node::check(double, double, double) const -> percentile threshold "==FBIAS1" used before it was set.
ERROR  : 

If this percentile() function returns a value of bad data or nan, then we see the bad behavior. And it's very likely that the problem lies in the line of code that I referenced previously. If the data stored in ordered_array[index+1] is nan, then the whole expression will evaluate to nan. This function just needs more bounds checking.

I'm going to write up a MET bug issue and fix this in the main_v11.0, main_v11.1, and develop branches. Unfortunately, that means you won't be able to use this feature in the 11.0.2 release.

[JohnHalleyGotway]

FYI, here's the MET GitHub bugfix issue: dtcenter/MET#2644

[malloryprow]

Great, thanks John! Let me know when you need any testing from me.

[jprestop]

I am going to go ahead and close this Discussion since there is a related issue. If anything else needs to be addressed here, please comment on the issue: dtcenter/MET#2644