如何在gnuplot绘图循环中检查column(x)/bin_width的有效值？

背景与现有代码 #

数据文件a_test.dat内容：

# block_1
x y
1 1979.2
1 1598.6
1 1430.8
1 1544.7
1 1953.3
1 1288.8
1 1741.3
1 1791.1


# block_2
x y
2 2027
2 1591.4
2 1401.6
2 1923.7
2 1677
2 1311.1
2 1689.3
2 2015.8

已定义的参数与初始化代码：

x_min = 1250
x_max = 2250
n_bins = 20
bin_width = 1.*(x_max - x_min)/n_bins

array histo[n_bins]
do for [i=1:n_bins] {
    histo[i] = 0
}

leftInclude(x) = ( _i=int( floor( column(x)/bin_width ) ), histo[_i] = histo[_i] + 1 )

当前执行计算的命令：

data_file = "a_test.dat"
data_blocks = "block_1 block_2"

set table $Temp
    plot for [i=1:words(data_blocks)] data_file index word(data_blocks, i) using ( leftInclude("y") ) w table 
unset table

解决方法 #

要检查column(x)/bin_width的无效结果（如NaN），可以利用gnuplot内置的valid()函数，同时修正原代码中存在的数组越界风险：

1. 基础版：仅检查无效值 #

修改leftInclude函数，加入有效性判断：

leftInclude(x) = ( _val = column(x), 
                   valid(_val) ? ( 
                       _i = int(floor( (_val - x_min)/bin_width )) + 1,  # 修正索引计算，适配数组1起始索引
                       histo[_i] = histo[_i] + 1 
                   ) : 0 )

• valid(_val)会判断数值是否为非NaN、非无穷大的有效值
• 修正索引计算逻辑：通过(_val - x_min)/bin_width将数值映射到0_{19区间，加1后对应数组1}20的索引，避免越界

2. 增强版：同时过滤超出分箱范围的值 #

如果需要排除[x_min, x_max)范围外的有效值，可加入范围判断：

leftInclude(x) = ( _val = column(x), 
                   valid(_val) && _val >= x_min && _val < x_max ? ( 
                       _i = int(floor( (_val - x_min)/bin_width )) + 1, 
                       histo[_i] = histo[_i] + 1 
                   ) : 0 )

3. 可选：统计无效值数量 #

若需要记录无效值的总数，可额外定义计数变量：

invalid_count = 0  # 初始化无效值计数器

leftInclude(x) = ( _val = column(x), 
                   valid(_val) ? ( 
                       _val >= x_min && _val < x_max ? ( 
                           _i = int(floor( (_val - x_min)/bin_width )) + 1, 
                           histo[_i] = histo[_i] + 1 
                       ) : 0 
                   ) : ( 
                       invalid_count = invalid_count + 1, 
                       0 
                   ) )

执行完直方图计算后，invalid_count变量会存储所有无效值的数量。

内容的提问来源于stack exchange，提问作者magfan