1.本发明涉及一种基于geotools的等值面生成方法。
背景技术：2.目前，使用geotools实现idw等值面绘制结果需要兼顾精度和效率，支持等值面实时生成，并结合不同插值方法、借助气象工具重新生成等值面，以及需深入研究等值面实现算法，比较耗时。同时，目前生成等值面方法过程中网格构建普遍基于三角形网格的算法，此算法适用于雨量控制点少，面积较小的区域，对雨量控制点较多、面积较大的区域并不适用。因此，有必要设计一种既能同时平衡系统效率与计算精确度，又能适用于雨量控制点多、面积大区域的等值面生成的方案。
技术实现要素：3.针对上述现有技术的不足，本发明提供了一种基于geotools的等值面生成方法。
4.为实现上述目的，本发明采用的技术方案如下：
5.一种基于geotools的等值面生成方法，包括以下步骤：
6.(1)获取从前端传入的固定入参格式，参数包括控制点经纬度以及数值、研究区域边界数据、是否裁剪、插值格网数量大小、数据间隔；
7.(2)根据裁剪范围构建正方形网格，并为每个网格分配索引；
8.(3)采用反距离加权方法，传入已知点数据，计算出未知格网的属性值；
9.(4)根据计算出的格网数据以及已知分级数据，通过网格索引，获取连线的两个顶点信息，通过在点信息和编号，保存连线的最终点信息，生成等值线；
10.(5)计算等值线平滑度，并根据等值线生成等值面；
11.(6)判断是否需要根据研究区范围进行裁剪，是，则筛选出研究区外的坐标点，并返回步骤(5)重新生成等值面；否，则执行步骤(7)；
12.(7)将等值面结果以同样的格式返回前端，包括组成等值面的点经纬度坐标以及等值面的属性值；
13.(8)将等值面结果以可视化效果表达。
14.具体地，所述步骤(2)包括以下步骤：
15.(21)利用二次误差函数生成顶点坐标，并计算极点值，所述二次误差函数表达式如下：
[0016][0017]
式中，pi为交点的位置，ni为交点的法向，x为偏移量；
[0018]
(22)生成网格面片，其包括以下流程：
[0019]
a)创建1个4*4的零矩阵用于存放qr矩阵分解的结果；
[0020]
b)对于体素单元的每条相交边，计算交点的位置pi和对应的法向ni；
[0021]
c)将向量[ni.x,ni.y,ni.z,dot(pi,ni)]添加到4*4的零矩阵底部；
[0022]
d)通过qr矩阵分解得到3*3的上三角矩阵a'和向量b'；
[0023]
e)求解线性方程组a'ta'x＝(a'tb'-a'tb'c),其中c是体素单元中所有交点的质心位置；
[0024]
f)将计算得到的偏移量x加上质心位置c即为体素单元中的顶点坐标；
[0025]
g)如果计算得到的顶点坐标位于体素单元之外，则顶点坐标用质心位置c来代替；
[0026]
h)对于每一条相交的体素边，将其周围4个体素单元内的顶点连接生成1个四边形面片。
[0027]
具体地，所述步骤(3)中，反距离加权的表达式如下：
[0028][0029]
式中，z(s0)为s0处预测值；n为预测计算过程中要使用的预测点周围的样点数量；λi为预测计算过程中使用的各样点权重，该值随着预测点与样点之间的距离的增加而减少；z(si)为si处获得的测量值。
[0030]
进一步地，所述样点权重的确定公式如下：
[0031][0032]
式中，d
i0
是s0与各已知样点si之间的距离。
[0033]
与现有技术相比，本发明具有以下有益效果：
[0034]
现有的等值面算法针对控制点密度较大的研究区域来说，计算效率不佳，主要原因是传统方法很难平衡系统效率与计算准确度的关系，而本发明将三角形网格改进为正方形网格，并以此为基础，结合反距离加权方法，生成等值线，最终通过计算等值线平滑度，并根据等值线即可生成等值面。本发明方案简单、高效，既能解决雨量控制点多、面积大区域的等值面生成问题，同时还兼顾了平衡系统效率与计算精确度，很好地适用于雨量控制点较多、面积较大的区域。
附图说明
[0035]
图1为本发明-实施例的流程示意图。
[0036]
图2为本发明-实施例中的正方形网格示意图。
[0037]
图3为根据图2中的正方形网格所列出的16种可能结果的示意图。
[0038]
图4为三角形格网的8种可能结果的示意图。
具体实施方式
[0039]
下面结合附图说明和实施例对本发明作进一步说明，本发明的实施包含但不限于以下实施例。
[0040]
实施例
[0041]
本实施例提供了一种基于geotools的等值面生成方法，适用于雨量控制点较多、面积较大区域的等值面生成。如图1所示，本实施例的主要流程包括：
[0042]
第一步：获取从前端传入的固定入参格式(json格式)，参数包括控制点经纬度以及数值、研究区域边界数据、是否裁剪、插值格网数量大小、数据间隔。
[0043]
第二步：根据裁剪范围构建正方形网格，并为每个网格分配索引，具体流程如下：
[0044]
(1)利用二次误差函数生成顶点坐标，并计算极点值，所述二次误差函数表达式如下：
[0045][0046]
式中，pi为交点的位置，ni为交点的法向，x为偏移量；
[0047]
误差函数可以写成矩阵形式：
[0048][0049]
其中矩阵a的行向量为交点的法向ni，向量b的每个元素为ni·
pi。
[0050]
极值点可以通过求解正则方程得到：
[0051][0052]
但是这种方式会存在数值不稳定，本实施例采用了另一种解法，即：基于qr矩阵分解计算正交矩阵q，使得q与[a b]相乘为如下上三角矩阵形式：
[0053][0054]
其中a'为3*3的上三角矩阵，b'为长度为3的向量，r为标量。
[0055]
那么误差函数可以变化为：
[0056]
(ax-b)
t
(ax-b)＝(ax-b)
tqt
q(ax-b)
[0057]
＝(qax-qb)
t
(qax-qb)
[0058]
＝(a
′
x-b
′
)
t
(a
′
x-b
′
)+r2[0059]
然后再根据上式计算极值点；
[0060]
(2)生成网格面片，其包括以下流程：
[0061]
a)创建1个4*4的零矩阵用于存放qr矩阵分解的结果；
[0062]
b)对于体素单元的每条相交边，计算交点的位置pi和对应的法向ni；
[0063]
c)将向量[ni.x,ni.y,ni.z,dot(pi,ni)]添加到4*4的零矩阵底部；
[0064]
d)通过qr矩阵分解得到3*3的上三角矩阵a'和向量b'；
[0065]
e)求解线性方程组a'ta'x＝(a'tb'-a'tb'c),其中c是体素单元中所有交点的质心位置；
[0066]
f)将计算得到的偏移量x加上质心位置c即为体素单元中的顶点坐标；
[0067]
g)如果计算得到的顶点坐标位于体素单元之外，则顶点坐标用质心位置c来代替；
[0068]
h)对于每一条相交的体素边，将其周围4个体素单元内的顶点连接生成1个四边形面片。
[0069]
正方形网格与三角形网格在本实施例方案中的对比：
[0070]
如图2所示，正方形网格上每个顶点代表一个属性值，在本实施例中为降雨量。格网顶点属性值为使用idw插值方法根据样本点属性计算而来。令格网顶点的值大于等值线阈值时为1，小于阈值时为0，这样即可把结果保存在一个四位二进制内。
[0071]
将所有可能结果都列出来，既可得到如图3所示的16种可能结果。而传统三角形格网的所有可能结果如图4所示。
[0072]
可以看出，三角形网格的可能结果只有8种，且一个图形内无两条线的情况。
[0073]
综上，可以看出，正方形格网对于三角形格网来说，等值线与格网的交点更多，对于研究区范围较大的数据来说，减少了格网密度，增加了计算效率。
[0074]
第三步：采用反距离加权方法，传入已知点数据，计算出未知格网的属性值。
[0075]
本实施例中，反距离加权的表达式如下：
[0076][0077]
式中，z(s0)为s0处预测值；n为预测计算过程中要使用的预测点周围的样点数量；λi为预测计算过程中使用的各样点权重，该值随着预测点与样点之间的距离的增加而减少；z(si)为si处获得的测量值。
[0078]
样点权重的确定公式如下：
[0079][0080]
式中，d
i0
是s0与各已知样点si之间的距离。
[0081]
样点在预测点值的计算过程中所占权重的大小受p的影响；也就是说随着采样点与预测值之间距离的增加，标准样点对预测点影响的权重按指数规律减少。在预测过程中，各样点值对预测点值作用的权重大小是成比例的，这些权重值的总和为1。
[0082]
第四步：根据计算出的格网数据以及已知分级数据，通过网格索引，获取连线的两个顶点信息，通过在点信息和编号，保存连线的最终点信息，生成等值线；
[0083]
第五步：计算等值线平滑度，并根据等值线生成等值面。所采用的平滑函数如下：
[0084]
prop＝0
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
当prop≦1
[0085]
prop＝(leve1(i-1)+leve1(i))/2
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
当leve1(i-1)≤prop≤leve1(i)
[0086]
prop＝leve1(max)*1.5
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
当prop≥leve1(max)
[0087]
其中prop为属性值，本实施例中指降雨量。level为降雨量级别阈值；i为第i个级别。
[0088]
第六步：判断是否需要根据研究区范围进行裁剪，是，则筛选出研究区外的坐标点，并重新生成等值面；否，则将等值面结果以json格式返回前端，包括组成等值面的点经
纬度坐标以及等值面的属性值。
[0089]
第七步：采用可视化手段，将等值面结果以可视化效果表达。
[0090]
示例：
[0091]
根据研究区大小，经过大量实验，本实施例中，自定义了一张研究区面积与划分格网数量的对照表(计算时间小于30s为合格)，自动划分格网密度，达到系统效率与准确度直接最大的平衡。研究区面积与格网数量平衡的表格如下：
[0092]
研究区面积/万平方公里格网数量/个计算时间/秒1-101000*10002610-30500*5002830-50300*30028》50100*10023
[0093]
根据区间阈值大小进行已知点极限值的取舍，进行已知点属性的平滑，若已知点属性值位于阈值区间边缘，则代入平滑函数，平滑已知点属性，使其落在阈值中间，减少了入参的误差，增加了最终生成等值面的准确度。
[0094]
上述实施例仅为本发明的优选实施方式，不应当用于限制本发明的保护范围，凡在本发明的主体设计思想和精神上作出的毫无实质意义的改动或润色，其所解决的技术问题仍然与本发明一致的，均应当包含在本发明的保护范围之内。
技术特征：1.一种基于geotools的等值面生成方法，其特征在于，包括以下步骤：(1)获取从前端传入的固定入参格式，参数包括控制点经纬度以及数值、研究区域边界数据、是否裁剪、插值格网数量大小、数据间隔；(2)根据裁剪范围构建正方形网格，并为每个网格分配索引；(3)采用反距离加权方法，传入已知点数据，计算出未知格网的属性值；(4)根据计算出的格网数据以及已知分级数据，通过网格索引，获取连线的两个顶点信息，通过在点信息和编号，保存连线的最终点信息，生成等值线；(5)计算等值线平滑度，并根据等值线生成等值面；(6)判断是否需要根据研究区范围进行裁剪，是，则筛选出研究区外的坐标点，并返回步骤(5)重新生成等值面；否，则执行步骤(7)；(7)将等值面结果以同样的格式返回前端，包括组成等值面的点经纬度坐标以及等值面的属性值；(8)将等子面结果以可视化效果表达。2.根据权利要求1所述的一种基于geotools的等值面生成方法，其特征在于，所述步骤(2)包括以下步骤：(21)利用二次误差函数生成顶点坐标，并计算极点值，所述二次误差函数表达式如下：式中，p
i
为交点的位置，n
i
为交点的法向，x为偏移量；(22)生成网格面片，其包括以下流程：a)创建1个4*4的零矩阵用于存放qr矩阵分解的结果；b)对于体素单元的每条相交边，计算交点的位置p
i
和对应的法向n
i
；c)将向量[n
i
.x,n
i
.y,n
i
.z,dot(p
i
,n
i
)]添加到4*4的零矩阵底部；d)通过qr矩阵分解得到3*3的上三角矩阵a'和向量b'；e)求解线性方程组a'ta'x＝(a'tb'-a'tb'c),其中c是体素单元中所有交点的质心位置；f)将计算得到的偏移量x加上质心位置c即为体素单元中的顶点坐标；g)如果计算得到的顶点坐标位于体素单元之外，则顶点坐标用质心位置c来代替；h)对于每一条相交的体素边，将其周围4个体素单元内的顶点连接生成1个四边形面片。3.根据权利要求1或2所述的一种基于geotools的等值面生成方法，其特征在于，所述步骤(3)中，反距离加权的表达式如下：式中，z(s0)为s0处预测值；n为预测计算过程中要使用的预测点周围的样点数量；λ
i
为预测计算过程中使用的各样点权重，该值随着预测点与样点之间的距离的增加而减少；z(s
i
)为s
i
处获得的测量值。
4.根据权利要求2或3所述的一种基于geotools的等值面生成方法，其特征在于，所述样点权重的确定公式如下：式中，d
i0
是s0与各已知样点s
i
之间的距离。
技术总结本发明公开了一种基于geotools的等值面生成方法，包括：(1)获取从前端传入的固定入参格式；(2)根据裁剪范围构建正方形网格，并为每个网格分配索引；(3)采用反距离加权方法，传入已知点数据，计算出未知格网的属性值；(4)根据计算出的格网数据以及已知分级数据，通过网格索引，获取连线的两个顶点信息，通过在点信息和编号，保存连线的最终点信息，生成等值线；(5)计算等值线平滑度，并根据等值线生成等值面；(6)将等值面结果以json格式返回前端，包括组成等值面的点经纬度坐标以及等值面的属性值；(7)将等值面结果以可视化效果表达。本发明既能解决雨量控制点多、面积大区域的等值面生成问题，还兼顾了平衡系统效率与计算精确度。还兼顾了平衡系统效率与计算精确度。还兼顾了平衡系统效率与计算精确度。
技术研发人员：徐永全 杨胜利 方世禄 周定春 李超 邱辛 魏慧 吴旭欣 付尔康 刘孔震 傅康 贺加贝 雷双舟 冯杰
受保护的技术使用者：四川华能宝兴河水电有限责任公司
技术研发日：2022.07.06
技术公布日：2022/11/1
