基于埃尔米特运动预测的地形实时绘制

COMPUTER ENGINEERING AND DESIGN Feb. 2017

Vol. 38 No. 2

基于埃尔米特运动预测的地形实时绘制

高保祿，高锐军，王倩，窦明亮，张冰

(太原理工大学计算机科学与技术学院，山西太原030024)

摘要：在大规模地形实时绘制时，地形规模越大，帧速率越小，演染速度越慢，影响了漫游者在虚拟空间中的沉浸感。 基于此，提出一种基于埃尔米特运动预测的地形实时绘制方法。采用埃尔米特插值算法进行视点运动预测，提前加栽下一 视点位置及视线方向的地形数据，减少实时绘制时的动态加栽量，有效提高帧速率，避免漫游过程中画面停滞的现象；该 算法计算量较小，稳定性相对其它算法较高，在取适当步长的情况下，能获得期望的精度，不会出现龙格现象。实验结果 表明，在大规模地形的实时绘制时，该方法帧速率变化小，在保证渲染画面平滑和地形绘制实时性的同时，提高了地形渲 染速度，避免了视觉上的跳跃感。

关键词：地形；埃尔米特算法；运动预测；龙格；实时绘制中图法分类号：TP391. 9

文献标识号：A

文章编号：1000-7024

(2017) 02-0483-05

doi： 10. 16208/j. issnl000-7024. 2017. 02. 038

Real time rendering of terrain based on Hermite prediction

(College of Computer Science and Technology, Taiyuan University of Technology, Taiyuan 030024, China)

GAO Bao-lu, GAO Rui-jun, WANG Qian, DOU Ming-liang, ZHANG Bing

Abstract： Small frame rates affect the immersion of the roamer in the virtual space, because the terrain data are too large in real­

time rendering. Based on this? the method of viewpoint prediction based on Hermite was proposed. Hermite interpolation algo­

rithm was used to predict viewpoint movement, so that the frame rate was improved effectively and the rendering was steady. Compared with other algorithms? this method is more stable, and its calculation is smaller. The Runge phenomenon is avoided when an appropriate step size is chosen. The experimental results show that the variation of frame rates is small when rendering the large-scale terrain, and the rendering is steady.Keywords： terrain； Hermite interpolation； motion prediction； Runge； real-time rendering

〇引言

随着大数据时代的到来，以及虚拟现实领域的快速发 展，对大规模地形绘制的速度、逼真度有了更高的要求。 因此，漫游时视点的运动预测成为了研究热点。

本文在地形数据的预处理阶段，首先对地形数据进行 分块，其次对其构造合理的包围球结构，之后依据视锥裁 剪、屏幕误差等技术对大规模地形进行处理。同时，结合 当前视点位置的情况，根据埃尔米特插值算法对下一视点 位置及视线方向进行精确预测，将预测的可见地形提前调 人内存，减少内外存读写速度巨大差异对实时性的影响，收稿日期：2015-12-30;修订日期：2016-09-13

基金项目：虚拟现实技术与系统国家重点实验室开放基金项目（BUAA-VR-15KF-17)

作者简介：高保禄（19几-），男，山西太原人，博士，讲师，CCF会员，研究方向为智能信息处理；高锐军（1990 -)，男，山西吕梁人， 硕士研究生，研究方向为智能信息处理；王倩（1990 -），女，山西襄汾人，硕士研究生，研究方向为虚拟现实；窦明亮（1993 -)，男， 山西太原人，硕士研究生，研究方向为虚拟现实；张冰（1993 -），女，山西运城人，硕士研究生，研究方向为虚拟现实。

E-mail： 8538930l@qq.com

提高了大规模地形的渲染速度，保证了地形绘制的实时性； 埃尔米特插值算法相比拉格朗日插值算法，计算量比较小， 稳定性却很高，在取适当步长的情况下，能获得较高 的精度。

1地形数据预处理

大规模地形在三维实时绘制时，首先将大量的地形数

据调人内存，由于地形数据量较大，将地形数据进行分块， 实现部分数据加载。同时，结合基于包围球结构的视锥裁 剪算法，判断当前视点位置的可见范围。选择建立基于四 叉树和包围球结构的屏幕误差的动态LOD模型进一

• 484 •

计算机工程与设计

1.3四叉树LOD生成算法

2017 年

步保证了大规模地形的实时渲染。1.1地形数据分块

在大规模地形进行四叉树层次结构分割前，对地形瓦 片进行规则网格（regular square grid，RSG)的标准化操 作，使其满足如图1所示的分块结果。分块后的瓦片作为 预处理阶段构造包围球的基本单位，对每一个数据块的信 息进行保存。

依据文献[3]中的视锥裁剪策略，利用投影技术判断 视点的可见范围，再利用之前构造的包围球结构，计算出 所有落在投影区域的地形块。而这些地形块的集就是视点 在某一时刻、某一位置时的地形可见区域。而后通过相应 的坐标计算，加载可见地形块集。简言之，就是依据视相 关技术与四叉树包围球结构，来构建动态的LOD地 形模型。

令所有地形瓦片的中心点为进一步细分时的剖分点， 1.2基于四叉树的包围球构造

在地形分块后，对每个地形块的信息参数进行计算， 构造基于四叉树的包围球结构，以供实时渲染时视锥裁剪 的判断以及屏幕误差的计算判断。构造一个如图2所示的 正方形包围球。

球心处的空间几何误差计算如下。

根据球心Q处的高程值厶及正方形4个顶点的高程值 乙、厶、乙、厶计算包围球的半径，得下式

et =max( \\ Za — Zt | , \\ Zb — Zt | , \\ Zc — Zt | , \\ Zd — Zt | ,

Za +

+ Zc + Zd

4

(1)

根据A点在水平面的坐标（兄，Ya)，B点在水平面 的坐标（X，^)，计算正方形ABCD对角线的一半得 下式

AOt

V(xa-xcy + (Ya-Ycy

(2)

最终求得叶子节点z的半径为

i?f

=ef+Aa

(3)

而式（1)到式（3)计算任一地形块的几何误差与屏幕误 差，再结合下式判断其重要度，即是否需要进行细分

pt

,R,e) = /(dt —Rt)

(4)li = qoI

(2tan -^)

(5)

只，d[4]表示屏幕误差，co是视域广角p的像素个

数，忒为视点位置到包围球球心的距离。t为屏幕误差阈 值，即为判断地形是否需要进行细分的临界值。相比于T， 如果屏幕误差较大，则表明地形的划分未达到视觉所要求

的层次，需要继续划分，直到其层次满足要求。反之则表 明其层次满足要求，不必进行进一步的细分。对于某一地 形块节点z的细分判断如下

divide(z

)4~^pi C'Si，Ri，c)

/Cdi Ri)>^r

s

^ ) > 2视点运动预测

在虚拟场景空间中漫游时，视点往往会沿着某个轨迹

运动，为了保证视觉上地形的流畅渲染，避免视觉跳跃感， 进行视点的运动预测必不可少。采用埃尔米特插值算法[6,7] 进行下一个视点位置和视线方向的预测，提前将可见地形 数据调入内存，提高地形渲染速度，保证地形渲染的 流畅性。

文献[8]利用在拉格朗日算法基础上推导出的外推算 法进行预测，采用一个变量作为外插的步长，从而根据视 点运动的快慢调节步长大小以预测下一视点运动情况。但 这种方法有很大程度上的局限性，预测出下一视点的位置 及方向范围较大。因此，本文采用埃尔米特算法进行预测， 在拉格朗日算法的基础之上，在节点处求一阶导，保证在 节点上有一阶光滑性，从而保证视点运动预测的准确性。 2.1埃尔米特算法

由于埃尔米特插值算法是在拉格朗日插值算法基础之 上再次进行了一阶导数运算，故计算量也比较小，并且稳 定性很高，在取适当步长的情况下，能获得期望的精度， 不会出现龙格现象。

埃尔米特算法的数学定义请参见文献[10]，求出^十

2个

插值基函数及

(^ = 0，1

，…，

72)，每一

第38卷第2期 高保禄，高锐军，王倩，等：基于埃尔米特运动预测的地形实时绘制• 485 •

个基函数都是2〃十1次多项式，且满足条件[1°]

aj ( j〇k) ~ ^jkI〇，j•尹左，\\l,j = k,V t\\ — t〇 ,V t〇 — t\\ ,v h—t〇(19)

i (xk ) = 0；(7)

Pj (j〇k) = 0,^j (j〇k) = djk(j = 0,1, • • •,n)

示的形式

7 = 〇

由插值基函数所满足的条件式（7)，有

将，3=&+3/z代入得

X3 = 28X〇 - 27Xi +

2. 2

步长选取

将插值多项式H2„+1(:t) =H(:r)写成用插值基函数表

!2hX\\ + ld>hX\\

(20)

Hz^iCx) =

2 Vyiaj (x) + nij^j (x)] (8)

在平均帧速率为200 fps左右的情况下，分别设置/z为 20帧和40帧，随机改变视点的位置和视线方向，进行20 次埃尔米特插值运算，即A为40帧时运行800帧，/z为20 帧时运行400帧，对比不同步长的运行帧数得到如图3、图 H2^i (xk ) = yk, H'lltVx ( j〇k ) =

= 0,1, • • •, 77) (9)

利用拉格朗日插值基函数(z)，令

a； (x) ~

(ax~\\~b) Ij (x)

(10)

另有

aj Cscj ) ~ CciXj ~\\~ b) /) Caxj ~Cxj ) — 1

a j C-Tj ) — lj Cxj ) \\_a lj Cxj ) ~\\~ 2\\~ b) l j Cxj ) ^|

— 0(ID

整理得

I axj Jrb=l

1 a 2/ j (. Xj) — 0

(12)

解出

=-2/〇 (x7 ),b=l + 2xJlfJ ixj)

(13)

由于lj ( jo)(x — JT〇 ) ••• (JT — JT；-1 ) (x — X；+1 ) • • • (X — Xn )

{ Xj — JC〇 ) * • * { Xj — Xj-i ) { Xj — Xj+\\ ) • • * { Xj — Xn )

(14)

两端取对数再求导得

k=Q工厂工k(15)

k^j

于是

aj(x) = (1 — 2(jc — Xj)k-鄭 x)

(16)

k =^ j

〇 工j

工k

同理，得

j3j (x) = (jc — Xj) l) (jc)

(17)

可以写出分段三次Hermite插值函数的分段表达式

HU

) = (1 V + 2 jcz-+i — x i ))

(r\\Xi _1^-J^i+i) / 2^+(v

1 + 2

Xi — J〇i+i% / \\J〇i+l — Xi , yi+i

X ——JC/+1

(n

)(f^

) ’ 十(m

V Xi — J〇i+i ,y z+ii G [

;汁i]G = 〇，1，…，”一 1)

(18)

假设每帧渲染的时间是近似相等的，取变量A作为外

插步长，用

(f =

0，1，2，3，…）表系 ^ (f = 0，1，2，

3,…）时刻漫游时视点的不同状态，现假设刚进入场景的 初始时刻为广。，初始状态为X。，另外，广1 =广。十^，广2 =

h

+2/z。根据上述计算得时刻的Hermite插值公式为X

3

-2X

h — t〇 \\ /13 — t\\

t〇

—

-2X

h — h

t\\ ) V tt\\ ,X〇 +〇 —

t\\ — t〇 ,4所示曲线。

_

<#RWX

图4 A为20帧时预测结果比较曲线

由实验结果得，采用埃尔米特运动预测方法时，步长

的选取将对预测结果有影响。若选取较大步长，根据算法 预测的结果将产生较大误差，选取较小步长将得到更接近 于真实值的结果。但是由于步长选择过小会导致计算量的 增加，影响帧速率，因此要根据实际情况选择合适的步长。

3基于四叉树孤立分割的裂缝处理

3.1瓦片间及瓦片内的裂缝

大规模地形进行分块后，对单一地形块进行孤立分割。 假设根据RSG将地形分割成mXn个瓦片，某一时刻，视

• 486 •

计算机工程与设计

2017 年

点可见的瓦片集合为E = {en，eu，ei3，…，e44}，对每 个瓦片进行四叉树孤立分割。对每t分割时，由于分割自上 而下、自左而右，ei2和egi先分割，因此会在€22和eif、egg 和e21之间出现裂缝，即瓦片间裂缝。瓦片内部进行孤立分 割，结点连接处也会出现裂缝。

地形瓦片的依次分割，如图5所示^

图5

地形瓦片的依次分割

3.2裂缝消除

针对出现的裂缝，采用结点依赖关系消除导致相邻区 域结点层次相互影响，算法复杂度较高，很难用于大规模 地形可视化。因此，本文采用高程平均值法和添加拆分点 对存在于相邻结点间以及瓦片内部结点的裂缝进行消除[2]。

假设当前瓦片内结点Ba的分割层次为1，从而可以判 断出同层次的上邻Bb和左邻氏是否达到了精度要求，是 否需要继续分割。而Bb和a可能出现3种分割情况，如 图6所示。

图6

四叉树相邻结点分割情况

基于此，对Ba的分割情况进行分析，通过添加拆分点 消除瓦片内的裂缝[2:s

Ta为当前瓦片，Tb和都已完成分割，此时在点c、

d、e、h、j处出现裂隙。若边界结点的分割层次Ta大于

I

，裂隙点如图7中点c、d、e所示；若边界结点分割层

次Ta小于Tb，裂隙点如图8中点h、j所示。

对于以上情况，需要通过改变当前瓦片的边界高程值 或在其上添加拆分点来消除裂隙。对分割层办：Ta大于T。 的情况，采用高程平均值法，令c、d、e点的高程为点a和

b

对应边的线性插值，则插值后的c、d、e点必定在边ab

上，从而消除了裂隙。对于分割层次1小于Tb的情况， 在Ta中添加拆分点h、j来消除裂隙。

图7分割层次Ta大于的裂缝消除

4实验结果分析

本文采用基于埃尔米特的运动预测方法，实现了大规 模地形的实时绘制。实验系统为Windows 7旗舰版，实验 软件平台为Microsoft Visual Studio 2008和三维图形标准库

OpenGL〇其中地形数据采用的是The Puget Sound Area，

来自于 Large Geometric Models Archive at Georgia Tech 0

与文献[8]拉格朗日插值算法、文献[9]卡尔曼滤 波法进行预测准确性比较，在平均帧速率为200 fps左右的

情况下，一分钟内平均误差的折线如图9所示。

实验将文献[幻基于拉格朗日的外推算法，以及文

o7

o

oo6o _o班 o5 软

t7# 4

3

2 160

一•一文献

[8]

步长顺

本文算法

-

文献

[9]

图9 一分钟内平均误差的折线

第38卷第2期 高保禄，高锐军，王倩，等：基于埃尔米特运动预测的地形实时绘制

量，保证了地形绘制的实时性和流畅性。

• 487 •

献[9]卡尔曼滤波法，同埃尔米特插值算法进行了比较， 当步长值取10、20、30、40、50，漫游时每帧中预测值 的平均误差如图9所示。由上图可见，拉格朗日插值算法 有很大程度上的局限性，预测结果误差较大；卡尔曼滤波 法通过系统的状态预测方程和当前状态的协方差预测方 程，计算滤波增益系数矩阵、状态向量的最优估计以及状 态向量的协方差，虽然更加准确一些，但是平均帧速相对 较低（见表1),导致实时渲染的性能比较差；而本文采 用的埃尔米特插值算法预测误差比较小，相对于拉格朗日 插值算法，埃尔米特插值算法预测的视点位置与真实值更 实验结果进一步表明，本文提出的埃尔米特插值算法在对帧速率影响不大的情况下，对视点运动的预测较为准确，预测效果更逼真，并能保证大规模地形绘制的实时性。在下一步的研究工作中，将对GPU单指令多数据强大的并行计算能力问题进行进一步的研究，以使大规模地形渲染 更将快速，画面更加平滑。

参考文献：

[1] Judnich J, Ling N. Symmetric cluster set level of detail for 加接近。

设定漫游者在漫游空间的视点运动情况为匀速直线运 动，采用埃尔米特方法对漫游者的视点位置及视线方向进 行插值计算，每间隔30帧计算一次，将计算结果与不进行 预测的情况、文献[8]基于拉格朗日的外推算法、文献[9]卡尔曼滤波法进行比较，得到表1结果。

表1预测算法与帧速率的关系

视点运动轨迹

预测方法平均帧速

/£ps帧速率 平均方差

不进行预测267. 41433. 24文献

[8]拉格朗曰插值预测265. 34425.36匀速直线

文献

[9]

卡尔曼

滤波预测230. 57392. 74本文埃尔米

特插值预测

261. 32

412. 53

实验结果表明，虽然有视点运动预测的帧速率略微低 于无视点运动预测，但是前者帧速率方差要比后者小，与 文献[8]中采用的拉格朗日插值算法的视点运动预测相 比，埃尔米特插值算法的帧速率方差要小于前者，帧速与 前者相差不多，故得到的画面更为平滑；与文献[9]中采 用的卡尔曼滤波算法相比，埃尔米特插值算法在保证帧速 率方差与其相差不多时，更大程度保证了帧速，画面流畅 度更好，符合人的正常视觉感受。

5结束语

本文在地形数据加载之前通过视锥裁剪预测可见地形

块，将埃尔米特插值算法作为一种视点运动预测算法，避 免了帧速的巨大波动。大型场景中漫游时，确定初始位置 及初始视线，据此利用埃尔米特算法预测下一视点具体情 况，将在预测范围内的地形数据调人内存，可以有效避免 因地形的巨大反差而使内存调人量产生巨大的变化，且保 证了帧速率的相对稳定，防止地形变化较大而产生画面卡 滞的现象。埃尔米特插值算法是在拉格朗日插值算法基础 上求一阶导，计算量比较小，稳定性却很高，在取适当步 长的情况下，能获得较高的精度，减少了实时渲染的计算

real-time terrain rendering [C] //IEEE International Confe­

rence on Multimedia and Expo. IEEE, 2012：

320-324.[2] ZHANG Junfeng, YAO Zhihong. LOD algorithm 〇£ terrain based on conservative screen error and isolated division of quad­tree [J], Journal of Southwest Jiaotong University, 2013, 48 (4): 666-671 (in Chinese).[张俊峰，姚志宏.基于四叉树孤 立分割和屏幕误差的地形L0D算法[J].西南交通大学学报，2013, 48 (4)： 666-671.]

[3] WANGChaoyi, DING Ying, ZHANG Xin, et al. Algorithm based on projection for view frustum culling [J], Journal of System Simulation，2014 (10): 2412-2417 (in Chinese).[王

超一，丁莹，张鑫，等.一种基于投影的视锥体裁剪算法

[J].系统仿真学报，2014 (10): 2412-2417.]

[4] YANG Ying, FENG Liyin, ZHAO Jing, et al. Terrain rendering method based on blocking and enclosing ball error metric [J]. Com­puter Engineering, 2010，36 (15) ： 199-201 (in Chinese). [杨

莹，冯立颖，赵静，等.基于分块和包围球误差函数的地形绘

制方法

[J].计算机工程，2010, 36 (15): 199-201.]

[5] Sun Y，Dong Y, Mao K, et al. View-dependent progressive transmission and rendering for lunar model based on bicubic sub­division-surface wavelet [J], Advances in Space Research,2014, 53 (12)： 1848-1857.

[6] ZENG Changxiong. 3n+2 orders hermite interpolation polyno­mials and interpolation error [J]. Journal of Shaoyang Univer­sity (Natural Science Edition), 2010，7 (2) ： 9-12 (in Chi- nese).[曾长雄.3n + 2次Hennite插值多项式及插值误差 [J].邵阳学院学报（自然科学版），2010, 7 (2): 9-12.][7] gate the WANG Xiuyou, properties and applications of difference quotient FAN Jianzhong, XU Dongqing. To investi­[J], Journal of Fuyang Normal College： Natural Science Edition, 2013, 30 (3): 15-18 (in Chinese).[王秀友，范建中，徐冬 青.差商的性质及其应用探讨[J].阜阳师范学院学报：自然 科学版，[8] WANG Xiang, 2013, LEI Xiaoyong, 30 (3): 15-18.]

DAi Shuling. Real time rendering of large scale terrain based on viewpoint prediction [J], Journal of System Simulation, 2013，25 (6) ： 1202-1206 (in Chinese).[王

响，雷小永，戴树岭.基于视点预测的大规模地形的实时渲染

[J].系统仿真学报，2013, 25 (6): 1202-1206.]

(下转第556页）

• 556 •

计算机工程与设计

2017 年

frastructure management information system with GIS and IOT[J]. Geomatics World, 2010, 8 (5)： 17-21 (in Chinese). [乔彦友，李广文，常原飞，等.基于GIS和物联网技术的基 础设施管理信息系统[J].地理信息世界，2010，8 (5): 17-21.]

[8] XU Zheng, ZHAO Yuexuan, LIU Chong, et al. GIS based architecture for WSNs and its application for railway track moni­toring [J], Electronic Measurement Technology, 2014，37(10) ： 100-103 (in Chinese).[徐征，赵悦镟，刘冲，等.基 于GIS的无线传感网架构及在轨道监测的应用[J].电子测量 技术，2014，37 (10): 100-103.]

[9] YIN Fang, FENG Min, LIU Lei, et al. WebGIS solutionbased on OGC and Open Source projects [J], Computer Engi­neering and Applications, 2010，46 (30): 68-70 (in Chi­nese). [尹芳，冯敏，刘磊，等.基于OGC规范和OpenSource项目的WebGIS开发与应用[J].计算机工程与应用， 2010, 46 (30)： 68-70.]

[10] ZHANG Liming, YAN Haowea Development and application of

WebGIS based on OGC specifications [J], Geomatics Spatial Information Technology, 2011，34 (6) ： 41-43 (in Chinese). [张黎明，闫浩文.基于OGC规范的WebGIS设计与应用[J]. 测绘与空间地理信息，2011，34 (6): 41-43.]

[11] SHI Wentao, Research on key technologies of Html5 Web-

Socket protocol and implementation of real-time Web conmmu- nication system based on WebSocket protocol [D]. Nanjing：Nanjing University 〇£ Posts and Telecommunications, 2014 (in Chinese).[石文涛.Html5中WebSocket协议关键技术 的研究及基于WebSocket协议的实时Web通信系统的实现 [D].南京：南京邮电大学，2014.]

(上接第487页）

[12] MA Hui, HUANG Qiubo, SU Houqin. Efficient and reliable

data distribution released by incorporating TCP socket with Web services [J], Computer Applications and Software, 2011，28 (9): 144-147 (in Chinese).[马辉，黄秋波，苏厚 勤.TCP Socket通信与Web服务结合实现高效可靠的数据分 发[J].计算机应用与软件，2011，28 (9): 144-147.][13] ZHOU Leqin, YAN Cairong, SU Houqin. An application

study on web-socket protocol-based push-sending data tech­niques in monitoring system [J]. Computer Applications and Software, 2013, 30 (5)： 229-232 (in Chinese).[周乐钦， 燕彩蓉，苏厚勤.基于Web-Socket协议的推送数据技术在监 控系统中的应用研究[J].计算机应用与软件，2013, 30 (5)： 229-232.]

[14] CAI Yi, Design of network management software based on

MINA framework [J], Communications Technology, 2013 (3): 115-117 (in Chinese).[蔡谊.基于 MINA框架的网络 管理软件设计[J].通信技术，2013 (3): 115-117.]

[15] GUO Lixia, LI Guoqing, XIE Jibo, et al. Research on dy­

namic WMS publishing technology of spatial data [J], Geo­matics World, 2011, 9 (2)： 58-62 (in Chinese).[郭丽霞， 李国庆，解吉波，等.空间数据对象的动态WMS发布技术 研究[J].地理信息世界，2011，9 (2): 58-62.]

[16] GUO Pengwei, FENG Wei. Design and implementation of

machine vision intelligent terminal based on ModBus/Tcp pro­tocol and custom protocol [J], Computer Digital Enginee-ring, 2010, 38 (12) ： 85-88 (in Chinese).[郭朋軎伟，冯伟.基于ModBus/Tcp协议和自定义协议的机器视觉智能终端的设计与实现[J].计算机与数字工程，2010，38 (12):85-88.]

[9] JIANG Yi. The key technology of terrain rendering based on [11] Hildebrandt D. Imag亡based styling [J], Visual Computer, LOD dynamic [D]. Wuhan： Wuhan University of Science 2015, 32 (4)： 445-463.

[12] Gunther T, Theisel H. Finite-time mass separation for com­and Technology, 2015 (in Chinese). [江益，基于动态

parative visualizations of inertial particles [C] //Computer LOD地形渲染的关键技术研究[D].武汉：武汉科技大学，

Graphics Forum, 2015 ： 471-480.2015.]

[10] ZHU Jingmei, Constructed formula of Hermite interpolate [13] Dong YF，Sun YK，Tang ZS, Interactive visualization of 3D

[J], Journal of Shandong University of Technology, 1998 lunar model with texture and labels, using Chang? E-l data (5)： 488-495 (in Chinese).[祝精美.Hermite 插值的构造 [J], Science China Physics Mechanics Astronomy, 2013, 型公式[J].山东工业大学学报，1998 (5): 488-495.]56 (10)： 2002-2008.