徐曉文,聶蕓
(華北計(jì)算技術(shù)研究所北京100083)
隨著時(shí)代的發(fā)展,計(jì)算機(jī)的技術(shù)越發(fā)的先進(jìn),同時(shí)數(shù)據(jù)觀測(cè)和采集儀器的進(jìn)步,可以生成大量的氣象數(shù)據(jù)。氣象數(shù)據(jù)具有很強(qiáng)的時(shí)空性,表現(xiàn)在一個(gè)氣象信息處理系統(tǒng)中的數(shù)據(jù)源既有同時(shí)間不同空間的數(shù)據(jù),也有同空間不同時(shí)間的數(shù)據(jù)[1]。同時(shí),通過不同尺度的觀測(cè)也會(huì)有不同的比例和精度。
氣象可視化科學(xué)研究已經(jīng)比較成熟,能夠?qū)庀笥^測(cè)的要素?cái)?shù)據(jù)、分析預(yù)報(bào)的產(chǎn)品進(jìn)行按空間高度層次、按時(shí)次的分析和顯示,如等值線分析顯示、填色分析與顯示等,在氣象水文預(yù)報(bào)業(yè)務(wù)和氣象水文信息應(yīng)用系統(tǒng)中廣泛使用[2]。隨著氣象水文信息對(duì)日常生活、行業(yè)工作的影響日益被重視,在應(yīng)用系統(tǒng)中提供氣象水文環(huán)境信息的感知能力成為迫切的需求,同時(shí)也為了更好地進(jìn)行氣象環(huán)境應(yīng)用分析,更加充分、直觀地展現(xiàn)氣象水文信息,對(duì)氣象信息可視化提出了更高的要求。以往單層次、單時(shí)次的靜止數(shù)據(jù)展現(xiàn)不能充分體現(xiàn)氣象數(shù)據(jù)時(shí)空連續(xù)的特征;天氣圖等面向?qū)I(yè)人員的方式也難以被一般用戶所接受。
本論文研究基于已有的氣象可視化技術(shù)基礎(chǔ),在不改變可視化的具體實(shí)現(xiàn)的過程情況下,利用現(xiàn)有的可視化數(shù)據(jù),通過多種方法實(shí)現(xiàn)中間插值,從而獲得氣象時(shí)空上連續(xù)的可視化效果。
在實(shí)現(xiàn)氣象時(shí)空數(shù)據(jù)動(dòng)態(tài)可視化過程中,需要每一個(gè)中間圖像生成的數(shù)據(jù),通過這些圖像的連續(xù)播放,從而實(shí)現(xiàn)氣象時(shí)空數(shù)據(jù)的動(dòng)態(tài)可視化。在很多情況下,氣象數(shù)據(jù)都是單時(shí)次,不連續(xù)的,所以需要將這些不連續(xù)的單時(shí)次的數(shù)據(jù)組合起來,如果直接生成可視化圖像進(jìn)行連續(xù)播放,那么大部分的氣象數(shù)據(jù)在感官上都存在明顯的撕裂感。所以需要利用這些原本的單時(shí)次,不連續(xù)的數(shù)據(jù)生成一連串的數(shù)據(jù)來補(bǔ)充這個(gè)畫面減少撕裂感的產(chǎn)生。這時(shí)候就需要應(yīng)用到插值技術(shù)。
由于是做氣象時(shí)空數(shù)據(jù)動(dòng)態(tài)可視化,所以需要進(jìn)行插值做出中間動(dòng)畫來實(shí)現(xiàn)動(dòng)態(tài)展示。考慮到分段低次插值函數(shù)都有一致收斂性,但光滑性較差,對(duì)于氣象時(shí)空數(shù)據(jù)動(dòng)態(tài)可視化要求有一定的光滑度,即最好保證存在二階連續(xù)導(dǎo)數(shù)。所以本文采用最常用的三次樣條函數(shù)[3]。三次樣條插值(Cubic Spline Interpolation)簡(jiǎn)稱Spline插值,我們選擇的是三次樣條插值[3]。這種插值技術(shù)比較容易生成平滑的曲線,并且在多個(gè)關(guān)鍵點(diǎn)上數(shù)據(jù)相關(guān),在保證關(guān)鍵的準(zhǔn)確的獲取到的數(shù)據(jù)的精度上,可以合理的生成出相關(guān)的數(shù)據(jù)[4]。同時(shí)由于數(shù)據(jù)是矩陣的形式,所以利用近似最鄰近特征匹配方法進(jìn)行采樣,獲取出一連串的點(diǎn),從而構(gòu)建出合理的三次樣條函數(shù)[5]。
之后通過線性條插值補(bǔ)全數(shù)據(jù)。通過這種方法形成的新的矩陣,不僅在數(shù)據(jù)上連貫性良好,還在初始?xì)庀髷?shù)據(jù)矩陣十分精確。
圖1表示了實(shí)驗(yàn)的具體流程。根據(jù)流程圖示意,首先數(shù)據(jù)格式要求為矩陣形式。那么初始數(shù)據(jù)為給定的一連串的時(shí)間上的矩陣數(shù)據(jù)。
假定初始存在5個(gè)數(shù)據(jù)場(chǎng),數(shù)據(jù)矩陣分別為M1,…,M5。這是一個(gè)固定時(shí)間間隔的序列矩陣,分別依照時(shí)間順序編號(hào)1,2,…,5。如圖2所示。
圖1 系統(tǒng)流程圖
圖2 初始矩陣序列圖
根據(jù)圖2所示,存在一個(gè)矩陣序列,包含依照時(shí)間排序的5個(gè)數(shù)據(jù)矩陣(其中2和4沒有畫出來),矩陣大小全部都為m×n,它們依照時(shí)間順序編號(hào)為1,2,3,4,5,本文主要實(shí)現(xiàn)的就是在這一個(gè)連續(xù)的數(shù)據(jù)矩陣序列中插入中間數(shù)據(jù)并利用這些數(shù)據(jù)和原始的數(shù)據(jù)生成一個(gè)動(dòng)態(tài)的可視化效果。根據(jù)生成動(dòng)畫的原理,為了獲得流暢的動(dòng)畫效果,需要在這些數(shù)據(jù)矩陣插入中間數(shù)據(jù)。
文中研究中假定初始數(shù)據(jù)為矩陣型網(wǎng)格數(shù)據(jù)。這種類型的數(shù)據(jù)是實(shí)際業(yè)務(wù)中最常見的格式。實(shí)際使用中,通過氣象水文資料的數(shù)據(jù)處理,就能夠直接得到這類數(shù)據(jù),通常按照行列順序存儲(chǔ)。對(duì)于氣象實(shí)況等基于觀測(cè)站點(diǎn)的離散數(shù)據(jù),可以通過網(wǎng)格化處理形成矩陣型數(shù)據(jù)。
氣象數(shù)據(jù)的網(wǎng)格節(jié)點(diǎn)數(shù)(行列數(shù))一般是按照在一定地理范圍內(nèi)固定間隔分布來計(jì)算的,比如全球預(yù)報(bào)數(shù)據(jù),按照4度、1度劃分,可形成90*45,360*180的數(shù)據(jù)矩陣[6]。
本章以4度劃分的矩陣為例進(jìn)行算法說明。
同時(shí),為了獲取到比較合理的數(shù)據(jù),本文在實(shí)驗(yàn)中所采用的數(shù)據(jù)為NOAA(美國國家海洋和大氣管理局)的數(shù)據(jù)(https://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.pressure.html),下 載 Air Temperature數(shù)據(jù)并提取文件里的數(shù)據(jù),縮小矩陣規(guī)模進(jìn)行實(shí)驗(yàn)。
文中存在5個(gè)初始數(shù)據(jù)矩陣M1,…,M5,這是一個(gè)依照時(shí)間順序排序的序列矩陣,分別賦予編號(hào)1,2,3,4,5。現(xiàn)在本文需要在這些矩陣中建立一定的映射關(guān)系用來方便插值。取矩陣M1和M2作為樣本,都是同等大小的m×n矩陣,要求矩陣M1的M1(x,y){(0<x<90,0<y<45,且x,y為整數(shù),其中x,y對(duì)應(yīng)的是矩陣中的坐標(biāo)位置)}求對(duì)應(yīng)到矩陣M2的最近臨近相似點(diǎn)[7]。
考慮到數(shù)據(jù)矩陣空間上的對(duì)應(yīng)與時(shí)間上的連續(xù),所以需要考慮空間對(duì)應(yīng)關(guān)系以及時(shí)間上的變化關(guān)系[8]。
本文所采用的方法就是求取
這個(gè)由M2矩陣的一部分元素組成的新的矩陣與M1(x,y)最相近的,此處最相近不僅僅代表空間上的相近,也有時(shí)間上數(shù)值的相近。由于采用的是近似最臨近特征匹配,所以本文需要考慮兩個(gè)矩陣中的距離因素,即為空間要素,假設(shè)這個(gè)距離為D12。為了確定距離的具體數(shù)值,本文做了以下的定義:
在M2中求取M1(x,y)中的近似最臨近特征,依照本文設(shè)定距離D12的數(shù)值為M1(x,y)與M2(x,y)對(duì)應(yīng),那么距離D12的數(shù)值為0,然后與M2(x-1,y),M2(x,y-1),M2(x+1,y)與M2(x,y+1)對(duì)應(yīng),那么距離D12的數(shù)值為1,其他距離D12的數(shù)值為2。
之后應(yīng)該計(jì)算所有條件下對(duì)每個(gè)對(duì)應(yīng)匹配占有的影響數(shù)值I,所利用的公式I=(M2(x,y)-M1(x,y))+D12。其中M1(x,y)可以分別變換成
矩陣中任何的元素。之后選取I最小的點(diǎn)作為匹配。實(shí)現(xiàn)過程如下所示[9]:
1)首先將原始數(shù)據(jù)矩陣M1與M2進(jìn)行分割,通過多區(qū)域采樣來提高數(shù)據(jù)的精確度。如圖3所示,兩個(gè)矩陣按照紅線所示分成等分的9個(gè)部分,分別是取橫軸和縱軸的三分之一作為一個(gè)分劃位置。這些分給出來的矩陣依照先從左到右再從上到下的規(guī)則 分 別 命 名 為M1p1,M1p2,…,M1p9與M2p1,M2p2,…,M2p9。之后各個(gè)相應(yīng)的區(qū)域找對(duì)應(yīng)關(guān)系,比如M1p1與M2p1對(duì)應(yīng);
圖3 原始數(shù)據(jù)矩陣區(qū)域劃分示意圖
2)以M1p1與M2p1對(duì)應(yīng)為例,本文采用的對(duì)應(yīng)關(guān)系如圖4所示。為圓點(diǎn)到星點(diǎn)的匹配過程。分割后M1p1與M2p1均為m/3×n/3大小的數(shù)據(jù)矩陣,其中圓點(diǎn)的取值為這個(gè)區(qū)域的中心點(diǎn)M1p1(m/6,n/6),星點(diǎn)的位置對(duì)應(yīng)矩陣如下:
之后的對(duì)應(yīng)關(guān)系重新依照這個(gè)方法對(duì)應(yīng)。同理M2與M3之間的關(guān)系重新匹配。每個(gè)關(guān)鍵點(diǎn)取由前一個(gè)中心點(diǎn)確定的作為樣本數(shù)據(jù)。即第一個(gè)樣本點(diǎn)的數(shù)值是Nm1p1,根據(jù)這個(gè)關(guān)系找到的是Nm2p1,但是Nm2p1可能已經(jīng)不是矩陣M2p1(m/6,n/6)的數(shù)據(jù)了。但是尋找Nm3p1的時(shí)候仍然以矩陣M2p1的M2p1(m/6,n/6)元素開始參考運(yùn)算。Nm2p1并不會(huì)參與到M2p1與M3p1的運(yùn)算中;
圖4 數(shù)據(jù)關(guān)鍵點(diǎn)的查找示意圖
3)后運(yùn)算影響數(shù)值數(shù)據(jù)I,當(dāng)存在多個(gè)相同的I時(shí),并且剛好I的取值為最小,那么根據(jù)氣象時(shí)空數(shù)據(jù)時(shí)間和空間都具有變換的關(guān)系,采用匹配點(diǎn)的順時(shí)針方向順序選取第一個(gè)I值最小的點(diǎn)作為匹配點(diǎn)。
通過上述方法,可以在數(shù)據(jù)矩陣M1和M2之間建立出一個(gè)對(duì)應(yīng)的關(guān)系。本文開始分割矩陣M1與M2,得到了多個(gè)數(shù)據(jù)矩陣,分別為M1p1,M1p2,…,M1p9與M2p1,M2p2,…,M2p9。之后以M1p1與M2p1為例,通過M1p1的中心點(diǎn)找到對(duì)應(yīng)的矩陣M2p1的所需求的數(shù)據(jù)和位置,之后再利用矩陣M2p1的中心點(diǎn)找到矩陣M3p1對(duì)應(yīng)的數(shù)據(jù)和位置,之后就可以得到一連串的5個(gè)數(shù)值,分別為Nm1p1,Nm2p1,… ,Nm5p1。N就是全部 5個(gè)序列矩陣演變關(guān)系的一個(gè)基準(zhǔn)向量,建立了5個(gè)矩陣的對(duì)應(yīng)關(guān)系。同理可以獲得其他的分塊矩陣的基準(zhǔn)向量。本文以之前的矩陣獲取到的基準(zhǔn)向量為例。
根據(jù)2.2可知,通過近似最臨近特征匹配本文在初始數(shù)據(jù)上獲取到了一連串的5個(gè)數(shù)值,分別為Nm1p1,Nm2p1,… ,Nm5p1。那么接下來為了構(gòu)建出這些關(guān)鍵點(diǎn)的中間插值數(shù)據(jù),本文考慮到氣象數(shù)據(jù)的流體性質(zhì),采用三次樣條插值[10]。首先根據(jù)初始數(shù)據(jù)是一個(gè)存在時(shí)間順序關(guān)系的數(shù)據(jù)矩陣,那么根據(jù)這個(gè)信息,賦予一個(gè)變量T。那么假設(shè)T=[1,2,3,4,5]。同理存在隨著這個(gè)序數(shù)改變的變量,即為上一小節(jié)獲取的 5 個(gè)數(shù)值Nm=[Nm1p1,Nm2p1,…,Nm5p1],那么就可以根據(jù)這個(gè)條件得到三次樣條函數(shù)。主要步驟如下所示[11]:
1)三次樣條插值多項(xiàng)式Sn(x)在每個(gè)小區(qū)間[xi-1,xi]上是3次多項(xiàng)式,其在此區(qū)間上的表達(dá)式如下:
其中,hi=xi+1-xi。因此,只要確定了Mi的值,就確定了整個(gè)表達(dá)式;
2)在此處x的取值為T,S(x)的取值為Nm。令:
則Mi滿足如下n-1個(gè)方程:
本文采用第一種邊界方法,所以得到以下的方程:
3)根據(jù)上文信息,存在n+1個(gè)方程,有著n+1個(gè)Mi未知數(shù),那么根據(jù)得到的方程解出Mi,就可以得到方程的解。
根據(jù)2.3得到的三次樣條插值的方程,假設(shè)本文需要求取M1p1和M2p1之間的插值數(shù)據(jù)。假設(shè)兩個(gè)矩陣之間插值數(shù)據(jù)為20個(gè)。這些數(shù)據(jù)和M1p1和M2p1是等時(shí)間間隔的,所以假設(shè)序號(hào)也是等間隔的。那么假設(shè)這些數(shù)據(jù)為N1l2p1n1,N1l2p1n2,…,N1l2p1n20。那根據(jù)時(shí)間次序的序列的劃分T=[1,1.05,1.10,…,2]。根據(jù)上一小節(jié)的方程,可以得到M1p1和M2p1之間的三次樣條函數(shù),其中自變量為T,通過T的數(shù)值來確定對(duì)應(yīng)的插值點(diǎn)的數(shù)值。
如圖5所示,在兩個(gè)矩陣之間插入數(shù)據(jù),不過現(xiàn)在插入的只是一些關(guān)鍵點(diǎn)的數(shù)據(jù),所以插入矩陣的數(shù)據(jù)不完全,它插入的是星號(hào)的數(shù)據(jù),所以需要對(duì)剩余部分進(jìn)行填充。
圖5 中間數(shù)據(jù)的插入過程
以N1l2p1n1為例子,假設(shè)已經(jīng)取到了N1l2p1n1,現(xiàn)在利用這個(gè)數(shù)據(jù)進(jìn)行補(bǔ)全成插值矩陣M1l2p1n1,這個(gè)插值矩陣對(duì)應(yīng)的序號(hào)為1.05。
圖6 補(bǔ)全插入數(shù)據(jù)的過程
如上所示,M1l2n1已經(jīng)獲取到了,同理可以得到更多的插值數(shù)據(jù)。本文采用的是現(xiàn)有的可視化手段,插值矩陣和原始數(shù)據(jù)構(gòu)造相同,所以可以直接進(jìn)行可視化展示[14]。
本文采用最常用的等值線方法進(jìn)行可視化??梢暬Ч鐖D7所示。
為了方便觀察,圖8為圖7通過提取的左上角部分區(qū)域的可視化顯示圖片,可以觀測(cè)到等值線有變化過程,但是大致的區(qū)域的變動(dòng)不大。說明這個(gè)數(shù)據(jù)插值進(jìn)行可視化是有效果的。同時(shí)在運(yùn)行過程中并沒有感受到明顯的卡頓,而且動(dòng)畫的過度也比較流暢,證明了優(yōu)化的效果。
本文研究需要大量的成熟數(shù)值方法、圖形圖像處理方法支撐,為了能夠快速實(shí)驗(yàn)多種方法,使用了科學(xué)計(jì)算工具matlab作為實(shí)驗(yàn)平臺(tái)。
實(shí)驗(yàn)中采取5個(gè)連續(xù)45×90的初始數(shù)據(jù)矩陣,每個(gè)數(shù)據(jù)之間插入300個(gè)數(shù)據(jù),去除可視化步驟。在實(shí)驗(yàn)中內(nèi)存占用如下:
圖7 可視化顯示序列圖
圖8 可視化顯示部分序列圖
物理內(nèi)存:4135 MB。
交換頁面:5036 MB。
虛擬內(nèi)存:6341 MB。
在使用e3處理器,16 g內(nèi)存的環(huán)境下,所需要的運(yùn)行時(shí)間為15.051630 s;即為每個(gè)數(shù)據(jù)插入平均耗時(shí)不到0.02 s,滿足可視化要求。
本論文研究基于已有的氣象可視化技術(shù)基礎(chǔ),在不改變可視化的具體實(shí)現(xiàn)的過程情況下,利用現(xiàn)有的可視化數(shù)據(jù),通過多種方法實(shí)現(xiàn)中間插值,利用插值的數(shù)據(jù)進(jìn)行可視化展現(xiàn),然后連續(xù)播放這些可視化效果,從而獲得氣象時(shí)空上連續(xù)的可視化效果。論文設(shè)計(jì)的方法,分離了可視化方法和數(shù)據(jù)處理方法,使得動(dòng)態(tài)展現(xiàn)的效果能夠針對(duì)多種可視化方法實(shí)施,便于擴(kuò)展應(yīng)用;數(shù)據(jù)處理方法簡(jiǎn)單,可并行執(zhí)行,便于應(yīng)用實(shí)現(xiàn)。論文有些工作還有待進(jìn)一步展開[15],包括:
1)可視化效果可擴(kuò)展??紤]到與現(xiàn)有系統(tǒng)的結(jié)合使用,本文在說明氣象時(shí)空數(shù)據(jù)可視化方法時(shí),僅選取了等值線場(chǎng)和色彩場(chǎng)方法。對(duì)于時(shí)空數(shù)據(jù)可視化的呈現(xiàn)方式本身,也是值得研究的問題。
2)方法驗(yàn)證和改進(jìn)。論文采用的實(shí)驗(yàn)數(shù)據(jù)類型有限,還需要采用不同類型不同規(guī)模的真實(shí)數(shù)據(jù),對(duì)方法進(jìn)行測(cè)試,更好地確定方法的適應(yīng)性、有效性。