基于三次樣條插值的采樣數(shù)據(jù)光滑曲線形成法-

1、基于三次樣條插值的采樣數(shù)據(jù)光滑曲線形成法宋又廉(上海鐵道大學(xué)機(jī)械系上海,200331摘要計(jì)算機(jī)輔助測(cè)試所得記錄為一系列經(jīng)濾波的離散采樣數(shù)據(jù),要把它恢復(fù)成不失真的連續(xù)光滑曲線打印輸出,應(yīng)對(duì)數(shù)據(jù)進(jìn)行處理。對(duì)變化平緩的穩(wěn)態(tài)過程可采用擬合法;而對(duì)波動(dòng)較大,采樣點(diǎn)較多的動(dòng)態(tài)過程必須用插值法才能真實(shí)反映這種波動(dòng)。介紹了三次樣條插值函數(shù)的形成法。指出常規(guī)方法僅能指定過程始、終點(diǎn)被測(cè)物理量變化的速度或加速度,不能完整描述物理量的變化狀況。提出在始端和末端規(guī)定兩個(gè)不獨(dú)立的附加虛擬采樣節(jié)點(diǎn)可以增加方程組的兩個(gè)自由度,使其能同時(shí)指定過程始終點(diǎn)被測(cè)物理量的速度和加速度,從而得到符合實(shí)際情況的完美結(jié)果。關(guān)鍵詞:采樣數(shù)據(jù)

2、;樣條插值;光滑曲線中圖分類號(hào):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 系統(tǒng)對(duì)被測(cè)物理量變化過程的記錄是以離散點(diǎn)的形式實(shí)現(xiàn)的,因此實(shí)驗(yàn)記錄的輸出有一個(gè)把離散點(diǎn)重新以合理的方法轉(zhuǎn)化成連續(xù)曲線的問題。試驗(yàn)中的物理量變化過程有穩(wěn)態(tài)和瞬態(tài)兩種形態(tài),由于

6、物理量在穩(wěn)態(tài)過程中變化緩慢,因此設(shè)定的采樣周期較長(zhǎng),實(shí)測(cè)數(shù)據(jù)點(diǎn)少,在CA T 系統(tǒng)中可以采用擬合法生成穩(wěn)態(tài)曲線。對(duì)瞬態(tài)過程因其變化迅第14卷第2期1999年6月數(shù)據(jù)采集與處理Jou rnal of D ata A cqu isiti on &P rocessing V o l .14N o.2Jun .1999速且伴有波動(dòng),因此采樣周期必須很短,采樣數(shù)據(jù)一般都在幾十點(diǎn)以上。即使如此,若用直線連接采樣數(shù)據(jù)其結(jié)果仍是折線而不是平滑的曲線。然而采用擬合法是不合適的,因?yàn)橐獢M合出曲線的波動(dòng)勢(shì)必要提高曲線的方次,擬合方次過高不但求解困難而且曲線會(huì)產(chǎn)生與實(shí)際過程不符的波動(dòng)。合理的做法是以插值代擬合

7、,常用的插值法有代數(shù)插值,三角多項(xiàng)式插值及樣條函數(shù)插值等。其中三次樣條函數(shù)以其簡(jiǎn)單易求而獲廣泛采用。樣條插值的基本做法是給出了一系列離散點(diǎn)的坐標(biāo),用三次樣條函數(shù)聯(lián)系相鄰點(diǎn),相鄰樣條函數(shù)在銜接點(diǎn)不僅有相同的函數(shù)值而且有相同的切線和曲率,因此非常光順。整個(gè)插值函數(shù)就由這些分段三次樣條鏈接而成1。三次樣條函數(shù)的求解需要指定始點(diǎn)和終點(diǎn)的邊界條件,常見的方法有兩種:(1指定始、終點(diǎn)的加速度值a 1和a n ;(2指定始、終點(diǎn)的速度值v 1和v n ,然后求解24。這兩種方法都有片面性,因?yàn)樵囼?yàn)中被測(cè)物理量的瞬態(tài)變化過程實(shí)際就是它受激后的運(yùn)動(dòng)過程。在開始記錄時(shí)該物理量應(yīng)處于一種初始穩(wěn)定狀態(tài),即此時(shí)的速度與

8、加速度均為零;在記錄的終了則變?yōu)榱硪环N穩(wěn)定狀態(tài),其速度與加速度恢復(fù)為零。因此僅指定端點(diǎn)速度或僅指定端點(diǎn)加速度都是片面的,正確的做法是同時(shí)指定端點(diǎn)速度和加速度,即v 1,a 1和v n ,a n5,然而籍以求解分段三次樣條函數(shù)的方程組沒有足夠的自由度同時(shí)滿足這些邊界條件。本文將簡(jiǎn)敘常規(guī)三次樣條曲線的形成法及其局限性,并提出同時(shí)指定端點(diǎn)速度與加速度時(shí)動(dòng)態(tài)過程的采樣數(shù)據(jù)三次樣條插值函數(shù)的求解方法。1指定端點(diǎn)速度或加速度時(shí)的三次樣條函數(shù)求解法令三次樣條函數(shù)Y (t ,定義在t 1,t n 區(qū)間上,采樣時(shí)間節(jié)點(diǎn)序列為t 1,t 2,t n -1,t n ,相應(yīng)物理量幅值為y 1,y 2,y n -1,y

9、 n ,若采樣時(shí)間間隔為h i =t i +1-t i ,則時(shí)間間隔序列為h 1,h 2,h n -1。Y (t 是由定義在t i ,t i +1區(qū)間上的分段三次多項(xiàng)式Y(jié) i (t (i =1,2,n -1鏈接而成的。因三次多項(xiàng)式的二階導(dǎo)數(shù)呈線性,故有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積分兩次,積分常數(shù)可據(jù)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由此可見,只要求得各節(jié)點(diǎn)上的加速度a i 就能解出各段三次多項(xiàng)式,然后鏈接成整個(gè)區(qū)間t ,t

11、n 上的樣條函數(shù)Y (t 。式(2是各種邊界條件下求解樣條函數(shù)的基本式。實(shí)際上計(jì)算機(jī)采樣時(shí)間間隔是固定的,此時(shí)可令h =h i (i =1,2,n -1。111指定端點(diǎn)加速度a 1,a n 時(shí)的求解法分段的三次樣條在節(jié)點(diǎn)處一階導(dǎo)數(shù)應(yīng)連續(xù),即Y(t i -0=Y (t i +0。在子區(qū)間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在子區(qū)間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應(yīng)用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數(shù)據(jù)采集與處理第14卷很顯然,n -2個(gè)方程不能解n 個(gè)未知量a i ,必須指定兩項(xiàng)邊界條件。最通常的做法是指定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 數(shù)學(xué)處理最容易,但并不合理,因?yàn)樵谒矐B(tài)過程中被測(cè)物理量的初速度和末速度應(yīng)為零。因此比較合理的邊界條件是指定始、終點(diǎn)速度v 1和v n 。112指定端點(diǎn)速度v 1,v n 時(shí)的求解法式(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聯(lián)立得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同時(shí)指定兩端點(diǎn)速度與加速度的三次樣條函數(shù)求解法在瞬態(tài)過程的測(cè)試中,采樣開始在過程啟動(dòng)之前,采樣結(jié)束應(yīng)在過程穩(wěn)定之后,此時(shí)被測(cè)物理量處于不變化狀態(tài)下,也即其速度和加速度均為零。因此模擬曲線的始點(diǎn)不僅斜率應(yīng)為零(一階導(dǎo)數(shù)為零,曲率也應(yīng)為零(二階導(dǎo)數(shù)為零,終點(diǎn)的情況也是一樣。在求取三次樣條函數(shù)以復(fù)現(xiàn)試驗(yàn)曲線時(shí),設(shè)定始點(diǎn)的邊界條件為v 1=a 1=0和v n =a n =0是合理的。由前述討論可知以式(2,5構(gòu)成的三次樣條

16、求取法只能容納兩項(xiàng)邊界條件,如同時(shí)指定四項(xiàng)邊界條件則方程(5沒有足夠的求解自由度,必須另設(shè)法為方程增加兩個(gè)自由度。方法是人為地在t 1和t 2之間以及t n -1和t n 之間各增加一個(gè)虛擬采樣點(diǎn),該兩點(diǎn)的時(shí)間坐標(biāo)可以任定(例如取在t 1,t 2及t n -1,t n 的正中間,但縱坐標(biāo)值不獨(dú)立,分別落在Y 1(t 和Y n -1(t 這兩段三次樣條上。為了敘述問題方便清楚起見,對(duì)測(cè)試采樣序列重新編號(hào),嵌入虛擬采樣點(diǎn)作為第二點(diǎn)和第n -1點(diǎn)。即令時(shí)間節(jié)點(diǎn)序列為t 1=0,t 2=h2(虛擬采樣點(diǎn),t 3=h ,t 4=2h ,t n -2=(n -4h ,t n -1=(n -4+1 2h (

17、虛擬采樣點(diǎn),t n =(n -3h 。問題求解虛擬點(diǎn)縱坐標(biāo)y 2,y n -1y 2=y 1+h 2v (t d t由于Y (t 是三次樣條,在從y 1到y(tǒng) 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期宋又廉:基于三次樣條插值的采樣數(shù)據(jù)光滑曲線形成法類似可得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該式適用于采樣點(diǎn)大于7的情況。事實(shí)上動(dòng)態(tài)過程測(cè)試的采樣點(diǎn)均在數(shù)十至數(shù)百點(diǎn),因此該式的實(shí)際應(yīng)用并不受影響。下圖為某機(jī)電速度控制系統(tǒng)階躍輸入的動(dòng)態(tài)響應(yīng)測(cè)試實(shí)例。.點(diǎn)為采樣數(shù)據(jù),被測(cè)物理量在采樣開始和結(jié)束時(shí)均處于穩(wěn)定狀態(tài),即其一階和二階導(dǎo)數(shù)為零,因此速度和加速度邊界條件為零,應(yīng)用上述算法得光滑插值曲線。圖1同時(shí)指定速度與加速度邊界條件時(shí)的三次樣條插值實(shí)例4結(jié)論(1CA T所得記錄為被測(cè)物理量的一系列離散采樣數(shù)據(jù),要恢復(fù)成連續(xù)光滑的時(shí)域曲線輸出,必須對(duì)數(shù)據(jù)進(jìn)行數(shù)學(xué)處理。對(duì)象的穩(wěn)態(tài)性