基于三次样条插值的采样数据光滑曲线形成法-

1、基于三次样条插值的采样数据光滑曲线形成法宋又廉(上海铁道大学机械系上海,200331摘要计算机辅助测试所得记录为一系列经滤波的离散采样数据,要把它恢复成不失真的连续光滑曲线打印输出,应对数据进行处理。对变化平缓的稳态过程可采用拟合法;而对波动较大,采样点较多的动态过程必须用插值法才能真实反映这种波动。介绍了三次样条插值函数的形成法。指出常规方法仅能指定过程始、终点被测物理量变化的速度或加速度,不能完整描述物理量的变化状况。提出在始端和末端规定两个不独立的附加虚拟采样节点可以增加方程组的两个自由度,使其能同时指定过程始终点被测物理量的速度和加速度,从而得到符合实际情况的完美结果。关键词:采样数据

2、;样条插值;光滑曲线中图分类号:T P 27412;O 24113S m ooth Curve Output for D ynam ic Process Sam pl i ng Da taBa sed on Cub ic Polynom i a l I n terpola tionS ong Y ou lian(D epartm en t of M echan ical Engineering ,Shanghai T iedao U n iversity Shanghai ,200331Abstract T he reco rd of CA T is a series of sam p lin

3、g data .In o rder to p rin t ou t the reco rd as a s m oo th cu rve ,the data shou ld be fitted o r in terpo lated .T he fo rm er is fo r steady p rocess and the latter is esp ecially fo r dynam ic p rocess that m ay have b ig fluctuati on .A pp roaches fo r so lving cub ic po lynom ial sp line in t

4、erpo lati on are in tro 2duced ,that can on ly define the accelerati on s o r velocities at the beginn ing and the end .T h is is no t com p lete to describe the m oving statu s of a p hysical variab le .B y defin ing tw o ex tra kno ts the freedom of the equati on s is increased and bo th veloci 2t

5、ies and accelerati on s at the tw o ends be defined .T h is resu lts in a com p lete ou tp u t cu rve .Key words :sam p ling data ;sp line in terpo lati on ;s m oo th cu rve收稿日期:1998211211;修改稿收到日期:1999203220引言CA T 系统对被测物理量变化过程的记录是以离散点的形式实现的,因此实验记录的输出有一个把离散点重新以合理的方法转化成连续曲线的问题。试验中的物理量变化过程有稳态和瞬态两种形态,由于

6、物理量在稳态过程中变化缓慢,因此设定的采样周期较长,实测数据点少,在CA T 系统中可以采用拟合法生成稳态曲线。对瞬态过程因其变化迅第14卷第2期1999年6月数据采集与处理Jou rnal of D ata A cqu isiti on &P rocessing V o l .14N o.2Jun .1999速且伴有波动,因此采样周期必须很短,采样数据一般都在几十点以上。即使如此,若用直线连接采样数据其结果仍是折线而不是平滑的曲线。然而采用拟合法是不合适的,因为要拟合出曲线的波动势必要提高曲线的方次,拟合方次过高不但求解困难而且曲线会产生与实际过程不符的波动。合理的做法是以插值代拟合

7、,常用的插值法有代数插值,三角多项式插值及样条函数插值等。其中三次样条函数以其简单易求而获广泛采用。样条插值的基本做法是给出了一系列离散点的坐标,用三次样条函数联系相邻点,相邻样条函数在衔接点不仅有相同的函数值而且有相同的切线和曲率,因此非常光顺。整个插值函数就由这些分段三次样条链接而成1。三次样条函数的求解需要指定始点和终点的边界条件,常见的方法有两种:(1指定始、终点的加速度值a 1和a n ;(2指定始、终点的速度值v 1和v n ,然后求解24。这两种方法都有片面性,因为试验中被测物理量的瞬态变化过程实际就是它受激后的运动过程。在开始记录时该物理量应处于一种初始稳定状态,即此时的速度与

8、加速度均为零;在记录的终了则变为另一种稳定状态,其速度与加速度恢复为零。因此仅指定端点速度或仅指定端点加速度都是片面的,正确的做法是同时指定端点速度和加速度,即v 1,a 1和v n ,a n5,然而籍以求解分段三次样条函数的方程组没有足够的自由度同时满足这些边界条件。本文将简叙常规三次样条曲线的形成法及其局限性,并提出同时指定端点速度与加速度时动态过程的采样数据三次样条插值函数的求解方法。1指定端点速度或加速度时的三次样条函数求解法令三次样条函数Y (t ,定义在t 1,t n 区间上,采样时间节点序列为t 1,t 2,t n -1,t n ,相应物理量幅值为y 1,y 2,y n -1,y

9、 n ,若采样时间间隔为h i =t i +1-t i ,则时间间隔序列为h 1,h 2,h n -1。Y (t 是由定义在t i ,t i +1区间上的分段三次多项式Y i (t (i =1,2,n -1链接而成的。因三次多项式的二阶导数呈线性,故有Y i (t =t i +1-t h iY i (t i +t -t i h iY i (t i +1i =1,2,n -1(1积分两次,积分常数可据Y i (t i =y i ,Y i (t i +1=y i +1求得,记Y i (t i +1=Y i +1(t i +1=a i +1于是Y i =(t i +1-t 36h i a i +(t

10、 -t i 36h i a i +1+c 1(t i +1-t h i+c 2(t -t i h i由Y i (t i =y i ,Y i (t i +1=y i +1可得c 1=y i -a ih i 26c 2=y i +1-a i +1h i 26上式可写为Y i (t =(t i +1-t 36h i a i +(t -t i 36h ia i +1+y i -a ih i 26(t i +1-t h i +y i +1-a i +1h i 26(t -t i h ii =1,2,n -1(2由此可见,只要求得各节点上的加速度a i 就能解出各段三次多项式,然后链接成整个区间t ,t

11、n 上的样条函数Y (t 。式(2是各种边界条件下求解样条函数的基本式。实际上计算机采样时间间隔是固定的,此时可令h =h i (i =1,2,n -1。111指定端点加速度a 1,a n 时的求解法分段的三次样条在节点处一阶导数应连续,即Y(t i -0=Y (t i +0。在子区间t i ,t i +1上Y i (t =-a i(t i +1-t 22h +a i +1(t -t i 22h +y i +1-y ih-a i +1-a i6h (3在子区间t i -1,t i 上Yi -1(t =-a i -1(t i -t 22h +a i (t -t i -122h+y i -y i

12、-1h-a i -a i -16h(4应用Yi (t i =Y i -1(t i ,于是ha i -1+4ha i +ha i +1=6y i +1-2y i +y i -1hi =2,3,n -1(5该方程组可写成矩阵形式A(n -2×n X n ×1=B(n -2×1其中Xn ×1=a 1,a 2,a n T491数据采集与处理第14卷很显然,n -2个方程不能解n 个未知量a i ,必须指定两项边界条件。最通常的做法是指定a 1和a n (而且通常是a 1=a n =0,则方程可写为A (n -2×(n -2X (n -2×1=

13、B(n -2×1X(n -2×1=a 2,a 3,a n -1T其中A =4h h0h4h h00h4hhh4h hh4h(n -2×(n -2B =6( y 3-2y 2+y 1 h -ha 16(y 4-2y 3-y 2 h6(y n -2y n -1+y n -2 h -ha n(n -2×1以这样的边界条件求解Y (t 数学处理最容易,但并不合理,因为在瞬态过程中被测物理量的初速度和末速度应为零。因此比较合理的边界条件是指定始、终点速度v 1和v n 。112指定端点速度v 1,v n 时的求解法式(5中不含有速度变量,但从式(3,4中可以获得速

14、度信息。由式(3,令t =t i ,i =1得Y 1(t 1=v 1=-a 1h2+y 2-y 1h-a 2-a 16h整理得2ha 1+ha 2=6y 2-y 1h-v 1(6由式(4,令t =t i ,i =n 得y n -1(t n =v n =a nh2+y n -y n -1h-a n -a n -16h整理得ha n -1+2ha n =6v n -y n -y n -1h(7将式 (6,7与式(5联立得An ×nXn ×1=B n ×1式中A =2h h0h4h h00h4h hh4h h0h2hn ×nB =6(y 2-y 1 h -v

15、16(y 3-2y 2+y 1 h6(y n -2y n -1+y n -2 h6(v n -(y n -y n -1 h n ×1由于v 1,v n 已知,方程组可解。2同时指定两端点速度与加速度的三次样条函数求解法在瞬态过程的测试中,采样开始在过程启动之前,采样结束应在过程稳定之后,此时被测物理量处于不变化状态下,也即其速度和加速度均为零。因此模拟曲线的始点不仅斜率应为零(一阶导数为零,曲率也应为零(二阶导数为零,终点的情况也是一样。在求取三次样条函数以复现试验曲线时,设定始点的边界条件为v 1=a 1=0和v n =a n =0是合理的。由前述讨论可知以式(2,5构成的三次样条

16、求取法只能容纳两项边界条件,如同时指定四项边界条件则方程(5没有足够的求解自由度,必须另设法为方程增加两个自由度。方法是人为地在t 1和t 2之间以及t n -1和t n 之间各增加一个虚拟采样点,该两点的时间坐标可以任定(例如取在t 1,t 2及t n -1,t n 的正中间,但纵坐标值不独立,分别落在Y 1(t 和Y n -1(t 这两段三次样条上。为了叙述问题方便清楚起见,对测试采样序列重新编号,嵌入虚拟采样点作为第二点和第n -1点。即令时间节点序列为t 1=0,t 2=h2(虚拟采样点,t 3=h ,t 4=2h ,t n -2=(n -4h ,t n -1=(n -4+1 2h (

17、虚拟采样点,t n =(n -3h 。问题求解虚拟点纵坐标y 2,y n -1y 2=y 1+h 2v (t d t由于Y (t 是三次样条,在从y 1到y 2的变化过程中加速度作线性变化,即a (t =a 1+a 2-a 1(h 2t故有v (t =v 1+t 0a (t d t =v 1+a 1t +a 2-a 1(h 2t 22y 2=y 1+h 2v1+a 1t +a 2-a 1(h 2t 22d t=y 1+h2v 1+h 212a 1+h 224a 2(8591第2期宋又廉:基于三次样条插值的采样数据光滑曲线形成法类似可得y n-1=y n-h n-12v n+h212a n+h2

18、24a n-1(9将式(8,9代入式(5并令h1=h2=h n-1=h n-2=h 2得3ha2+h2a3=12h(y3-y1-12v1-52ha1(103ha3+ha4=6h(2y1-3y3+y4+6v1+ha1(11ha i-1+4ha i+ha i+1=6h(y i-1-2y i+y i+1(12i=4,5,n-3ha n-3+3ha n-2=6h(y n-3-3y n-2+2y n-6v n+ha n(13h2a n-2+3ha n-1=12h(y n-2-y n+12v n-5h2a n(14式(10,11包含了始端边界条件v1,a1,式(13,14包含了终端边界条件v n,a n,

19、方程组获得了足够的自由度。式(1014可写成矩阵形式A(n-2×(n-2X(n-2×1=B(n-2×1其中X(n-2 ×1=a2,a3,a n-1TA=3h h 2003h h00h4h h00h4h h00h4h h00h3h00h 23h(n-2×(n-2B=12h(y3-y1-12v1-52ha16h(2y1-3y3+y4+6v 1+ha16h(y i-1-2y i+y i+1i=4,n-36h(y n-3-3y n-2+2y n-6v n+ha n12h(y n-2-y n+12v n-5h2a n(n-2×1该式适用于采样点大于7的情况。事实上动态过程测试的采样点均在数十至数百点,因此该式的实际应用并不受影响。下图为某机电速度控制系统阶跃输入的动态响应测试实例。.点为采样数据,被测物理量在采样开始和结束时均处于稳定状态,即其一阶和二阶导数为零,因此速度和加速度边界条件为零,应用上述算法得光滑插值曲线。图1同时指定速度与加速度边界条件时的三次样条插值实例4结论(1CA T所得记录为被测物理量的一系列离散采样数据,要恢复成连续光滑的时域曲线输出,必须对数据进行数学处理。对象的稳态性