李夢楚,胡春宏
(中國水利水電科學(xué)研究院 流域水循環(huán)模擬與調(diào)控國家重點(diǎn)實(shí)驗(yàn)室,北京 100038)
小波分析是20世紀(jì)80年代初發(fā)展起來的一種信號時(shí)、頻局部化分析新方法。目前在信號處理、圖像壓縮、語音編碼、模式識別、地震勘探、大氣科學(xué)以及許多非線性科學(xué)領(lǐng)域內(nèi)取得了大量的研究成果。隨著小波理論的形成和發(fā)展,它的優(yōu)勢逐漸引起許多水科學(xué)工作者的重視并且引入水文學(xué)科中[1-2]。水文時(shí)間序列是一個(gè)觀測到的樣本,是隨時(shí)間變化的時(shí)變序列,幾乎都是非平穩(wěn)的。對于這一類非平穩(wěn)序列,需要提取某一時(shí)段(或極小瞬間)的頻域信息,或某一頻段所對應(yīng)的時(shí)間信息。而小波變換是一種具有時(shí)頻多分辨率的分析方法,能把時(shí)間序列在時(shí)間和頻率域上同時(shí)展開,進(jìn)行多尺度細(xì)化分析,研究它們在不同周期隨時(shí)間變化的特性[1]。
從1993年Kumar和Foufoular-Gegious[3]將小波變換介紹到水文學(xué)中,到目前已取得了一定的研究成果,主要表現(xiàn)在水文系統(tǒng)多時(shí)間尺度分析、水文時(shí)間序列變化特性分析、水文系統(tǒng)預(yù)測預(yù)報(bào)和水文系統(tǒng)隨機(jī)模擬等方面。關(guān)于黃河水沙變化的小波分析多是關(guān)于徑流量的分析[4-5],也有一些關(guān)于上游輸沙量的分析[6],但對于黃河下游研究較少,也沒有人把小波分析引進(jìn)到河床演變研究中來。
本文選取黃河下游小浪底、花園口、高村和利津四站,利用小波理論研究此4站的水量、沙量序列多時(shí)間尺度變化,并以花園口、高村、利津三站平灘流量序列作為河道斷面形態(tài)特征參數(shù)與水沙序列變化進(jìn)行對比,對于黃河下游水沙變化及水沙變化對斷面形態(tài)的影響研究都有重要意義。
1.1 小波函數(shù) 小波函數(shù)是具有震蕩特性、能夠迅速衰減到零的一類函數(shù),即也稱基小波,其伸縮和平移構(gòu)成一簇函數(shù)系:
式中:Ψa,b(t)稱子小波;t為時(shí)間,a為尺度因子或頻率因子,反映了小波的周期長度;b為時(shí)間因子或平移因子,反映了在時(shí)間上的平移。
1.2 小波變換 對于能量有限信號f(t)∈L2(R),其連續(xù)小波變換定義為
當(dāng)研究離散信號f(iΔt)(i=1,2,…,N;N為樣本容量;Δt為取樣時(shí)間間隔)時(shí),式(2)的離散形式表達(dá)為
Wf(a,b)能同時(shí)反映時(shí)域參數(shù)b和頻域參數(shù)a的特性,它是時(shí)間序列f(t)或f(iΔt)通過單位脈沖響應(yīng)的濾波器的輸出。當(dāng)a較小時(shí),對頻域的分辨率低,對時(shí)域的分辨率高;當(dāng)a增大時(shí),對頻域的分辨率高,對時(shí)域的分辨率低。因此,小波變換實(shí)現(xiàn)了窗口的大小固定、形狀可變的時(shí)域局部化。
Wf(a,b)隨參數(shù)a和b變化,可作出以b為橫坐標(biāo),a為縱坐標(biāo)的關(guān)于Wf(a,b)的二維等值線圖,稱為小波變換系數(shù)圖。通過小波變換系數(shù)圖可得到關(guān)于時(shí)間序列變化的小波變化特征。在尺度a相同情況下,小波變換系數(shù)隨時(shí)間的變化過程反映了系統(tǒng)在該尺度下的變化特征;正的小波變換系數(shù)對應(yīng)于偏多期,負(fù)的小波變換系數(shù)對應(yīng)于偏少期,小波變換系數(shù)為零對應(yīng)于突變點(diǎn);小波變換系數(shù)絕對值越大,表明該時(shí)間尺度變化越顯著。
1.3 小波方差 將時(shí)間域上的關(guān)于a的所有小波變換系數(shù)的平方進(jìn)行積分,即為小波方差
小波方差隨尺度a的變化過程稱小波方差圖。它反映了波動(dòng)的能量隨尺度的分布。通過小波方差圖,可以確定一個(gè)水文序列中存在的主要時(shí)間尺度,即主周期。
水沙序列資料采用小浪底、花園口、高村、利津四站1952—2009年年水量、年沙量、汛期、非汛期水量、沙量資料。小浪底站年水沙量過程如圖1。
采用Morlet小波為母函數(shù)。Morlet小波為復(fù)數(shù)小波,定義為φ(t)=eiω0te-t22,式中ω0為常數(shù)且當(dāng)ω0≥5時(shí)Morlet小波就能近似滿足允許性條件,其小波變換系數(shù)的模和實(shí)部是兩個(gè)重要的變量。模的大小表示特征時(shí)間尺度信號的強(qiáng)弱,實(shí)部表示不同特征時(shí)間尺度信號在不同時(shí)間上的分布和位相兩方面的信息。將各序列和Morlet小波函數(shù)代入式(3),取不同的a和b,計(jì)算小波變換系數(shù)的模平方和實(shí)部,作它們的二維等值線圖,并求出小波方差,作圖。Morlet小波變換系數(shù)的模值表示能量密度,
由于運(yùn)用小波分析時(shí)間序列需較長資料序列,因此采用以端點(diǎn)為中心對稱延拓的方法將兩端數(shù)據(jù)外延,消除邊界效應(yīng)。同時(shí)為應(yīng)用方便將數(shù)據(jù)進(jìn)行距平處理。小波變換系數(shù)的模值越大,表明其所對應(yīng)的時(shí)段和尺度的周期性越明顯,模平方圖即能量圖;小波系數(shù)的實(shí)部包含著給定時(shí)間和尺度信號相對于其他時(shí)間和尺度信號的強(qiáng)度和位相兩方面的信息;小波方差的峰值則對應(yīng)主要周期,小波方差越大即表明在該尺度下周期變化明顯。
圖2(a)、(b)分別給出了小浪底站年水量距平序列Morlet小波變換系數(shù)的模平方和實(shí)部的時(shí)頻分布,(a)中虛線表示負(fù)位相,對應(yīng)偏少期,實(shí)線表示正位相,對應(yīng)偏多期,交替處即為突變點(diǎn)。由圖可以看出,18~30a時(shí)間尺度十分突出,周期中心在23a,在1952—2009年都表現(xiàn)明顯,正負(fù)位相交替出現(xiàn)。12~18a時(shí)間尺度也比較強(qiáng),周期中心在14a,在1952—2009年都表現(xiàn)明顯,正負(fù)位相交替出現(xiàn)。2~6a時(shí)間尺度在1952—1980年表現(xiàn)比較明顯,周期中心主要在3a左右,正負(fù)位相交替出現(xiàn)。6~11a時(shí)間尺度在1968—1995年表現(xiàn)比較明顯,其周期中心在9a左右,正負(fù)位相交替出現(xiàn)。
由圖3(a)、(b)可以看出年沙量周期變化、強(qiáng)弱、突變點(diǎn)分布及其位相結(jié)構(gòu)。8~16a時(shí)間尺度表現(xiàn)明顯在1952—1990年左右,周期中心在11a,1985—2009年在此時(shí)間尺度附近體現(xiàn)為11~19a及5~10a時(shí)間尺度,周期中心分別為14a和8a。17~30a時(shí)間尺度在1952—2009年全時(shí)段都體現(xiàn)比較明顯,周期中心在23a。3~7a周期在1952—1995年表現(xiàn)比較明顯,周期中心主要在5a左右。30a以上周期也有所體現(xiàn),但由于采用的時(shí)間序列較短,不很明顯。
用式(4)計(jì)算小波方差,并作出小浪底站年水沙量距平序列的小波方差圖如圖4、圖5。由圖可清楚看出,小浪底年水量具有3、9、14、23a的變化周期,其中23a是主要周期;小浪底年沙量具有5、11、23、36a的變化周期,其中11a是主要周期。1a自然周期能量也較大,但忽略之。
小浪底站年水量與年沙量序列總體來說變化規(guī)律是一致的,都體現(xiàn)出了4a左右、10a左右、23a左右的周期變化。同樣方法分析小浪底站汛期、非汛期水量、沙量以及花園口、高村、利津三站的年水沙量、汛期、非汛期水沙量,可得類似結(jié)果,將四站水沙變化周期列于表1(表中黑體表示為此序列主周期)??芍髡鞠嗤蛄凶兓疽恢拢晁颗c汛期、非汛期水量變化周期也基本一致,非汛期沙量與年沙量、汛期沙量的周期稍有差異。水沙總體變化周期都是一致的。
中國大陸華北地區(qū)年降水存在2、4、7、15、26a左右的主周期,15a左右為第一主周期,西北地區(qū)存在2、4、9、20a左右的主周期,20a左右為第一主周期[11]。而太陽黑子存在9a、41a的周期變化特征[10]。故可知黃河下游水沙變化周期與年降水變化相關(guān),可能也與太陽黑子活動(dòng)有關(guān)。
表1 四站水沙變化周期 (單位:a)
采用花園口、高村、利津三站1950—2008年平灘流量資料(來源于黃委會(huì))?;▓@口站平灘流量序列如圖6。
同樣的,對平灘流量序列也進(jìn)行延拓、距平,進(jìn)行如上分析。以花園口為例,分析如下:
圖7(a)、(b)分別給出了花園口站平灘流量距平序列mor?let小波變換系數(shù)的模平方和實(shí)部的時(shí)頻分布,由圖可以看出花園口站平灘流量周期變化、強(qiáng)弱、突變點(diǎn)分布及其位相結(jié)構(gòu)。18~30a時(shí)間尺度十分突出,周期中心在24a,在1950—2008年都表現(xiàn)明顯,正負(fù)位相交替出現(xiàn)。12~20a時(shí)間尺度在1950—1975年表現(xiàn)明顯,周期中心在16a,正負(fù)位相交替出現(xiàn)。7~12a時(shí)間尺度在1975—2009年表現(xiàn)比較明顯,周期中心主要在10a左右,正負(fù)位相交替出現(xiàn)。5~10a時(shí)間尺度在1950—1975年表現(xiàn)比較明顯,其周期中心在7a左右,正負(fù)位相交替出現(xiàn)。作小波方差圖如圖8??煽闯?4a時(shí)間尺度表現(xiàn)非常明顯,為花園口站平灘流量主要周期,其他周期表現(xiàn)較弱。
高村、利津兩站同理分析,將平灘流量變化周期列于表2(表中黑體表示為此序列主周期)。
由此可見,平灘流量的變化與水沙變化周期基本一致,平灘流量變化受到水沙變化周期影響。
為進(jìn)一步分析,選取小尺度3a、中尺度10a、大尺度23a的各序列小波變換系數(shù)實(shí)部變化作圖進(jìn)行分析(以小浪底年水沙量、花園口平灘流量為例)。
圖9給出了小浪底站年水沙量及花園口站平灘流量3a的小波變換實(shí)部變化過程,正相位對應(yīng)偏多時(shí)期,負(fù)相位對應(yīng)偏少時(shí)期,小尺度下變異點(diǎn)不具有代表性,有許多偽突變點(diǎn),應(yīng)在中尺度和大尺度下判斷,故變化點(diǎn)不一一列出。
表2 花園口、高村、利津平灘流量變化周期(單位:a)
由圖9可以看出,在20世紀(jì)70年代中期以前年水沙量序列與平灘流量序列變化基本一致;70年代中期到80年代中期這一時(shí)段內(nèi),3個(gè)序列變化不同步,影響因素比較復(fù)雜,由于多時(shí)間尺度是從整體序列分析出來的,而同一尺度下的豐枯變化又是相對于本身序列的數(shù)值變化得出的[6],初步分析認(rèn)為可能與1974年、1986年三門峽水庫調(diào)水調(diào)沙及1986年龍羊峽水庫建成運(yùn)用有關(guān);80年代后期以后3個(gè)序列變化逐漸趨于一致,90年代后期到2003年左右平灘流量序列與水沙序列變化稍有紊亂后歸于一致,可能與小浪底水庫建成運(yùn)用有關(guān)。
圖10給出了小浪底站年水沙量及花園口站平灘流量10a的小波變換實(shí)部變化過程,可以看出:小浪底年水量序列在1953年以前、1958—1963年、1968—1973年、1977—1982年、1986—1991年、1996—2002年、2007年以后各時(shí)段為負(fù)相位,表示年水量偏少;而1953—1958年、1963—1968年、1973—1977年、1982—1986年、1991—1996年、2002—2007年各時(shí)段為正相位,表示年水量偏豐。由此可推斷在此尺度下2009年以后兩三年內(nèi)仍然處于流量偏少時(shí)段,2012年以后會(huì)處于流量偏豐時(shí)段。小浪底年沙量序列在1954年以前、1960—1965年、1970—1976年、1981—1986年、1991—1994年、1998—2003年、2007年以后各時(shí)段為負(fù)相位,表示年沙量偏少;而1954-1960年、1965—1970年、1976—1981年、1986—1991年、1994—1998年、2003—2007年各時(shí)段為正相位,表示年沙量偏多。由此可推斷在此尺度下2009年以后兩三年內(nèi)仍然處于年沙量偏少時(shí)段,2012年會(huì)處于年沙量偏多時(shí)段?;▓@口平灘流量序列在1953年以前、1957—1963年、1969—1974年、1979—1984年、1989—1994年、2000—2006年各時(shí)段為負(fù)相位,表示平灘流量偏??;而1953—1957年、1963—1969年、1974—1979年、1984—1989年、1994—2000年、2006年以后各時(shí)段為正相位,表示平灘流量偏豐。由此可推斷在此尺度下2009年以后兩年內(nèi)仍然處于平灘流量偏豐時(shí)段,2011年以后會(huì)處于平灘流量偏豐時(shí)段。
在中尺度下,水沙變化及平灘流量變化是一致的。沙量變化及平灘流量變化始終跟隨水量變化,符合實(shí)際力學(xué)規(guī)律,且能夠體現(xiàn)出河床演變滯后于水沙變化的現(xiàn)象,沒有明顯的紊亂的時(shí)段。
圖11給出了小浪底站年水沙量及花園口站平灘流量23a的小波變換實(shí)部變化過程,可以看出:小浪底年水量序列在1957年以前、1968—1980年、1991—2003年各時(shí)段為負(fù)相位,表示年水量偏少;而1957—1968年、1980—1991年、2003年以后各時(shí)段為正相位,表示年水量偏豐。由此可推斷在此尺度上2009年以后5、6a內(nèi)仍然處于流量偏豐時(shí)段,2015年以后會(huì)處于流量偏少時(shí)段。小浪底年沙量序列在1957—1969年、1980—1992年、2003年以后各時(shí)段為負(fù)相位,表示年沙量偏少;而1957年以前、1969—1980年、1992—2003年各時(shí)段為正相位,表示年沙量偏多。由此可推斷此尺度上2009年以后5、6a內(nèi)仍然處于年沙量偏少時(shí)段,2015年以后會(huì)處于年沙量偏多時(shí)段?;▓@口平灘流量序列在1956年以前、1967—1979年、1991—2003年各時(shí)段為負(fù)相位,表示平灘流量偏小;而1956—1967年、1979—1991年、2003年以后各時(shí)段為正相位,表示平灘流量偏豐。由此可推斷2009年以后5、6a內(nèi)仍然處于平灘流量偏豐時(shí)段,2015年以后會(huì)處于平灘流量偏少時(shí)段。
在大尺度下可看到,水沙序列呈相反趨勢,水量偏多期對應(yīng)于沙量偏少期,反之亦然。而平灘流量的變化基本與年水量變化一致,略有所超前,在圖9上補(bǔ)充了小浪底站汛期水量的變化,發(fā)現(xiàn)平灘流量的變化與汛期水量的變化更為吻合,是跟隨汛期水量變化而變化的,所以在長時(shí)間尺度上汛期水量對于造床的貢獻(xiàn)可能更大。
綜上可知,在大、中、小不同時(shí)間尺度上,水沙序列、平灘流量序列的小波變換實(shí)部變化基本一致,在小尺度和大尺度上體現(xiàn)出的不同步可能有兩個(gè)主要原因:一是黃河流域水沙異源的特性,二是人類活動(dòng)特別是大型水庫建設(shè)的影響。
將小浪底年水沙序列及花園口平灘流量序列在中、大尺度上得變異點(diǎn)列于表3。同尺度下各序列變異點(diǎn)基本吻合,比較穩(wěn)定。序列的奇異點(diǎn)可能在不同尺度上是不一致的,與時(shí)間尺度有關(guān)。一般在較大尺度上尋找奇異點(diǎn),從而判斷奇異性的存在,但也要綜合各個(gè)尺度進(jìn)行系統(tǒng)分析。此結(jié)果基本吻合1960年三門峽水庫建成投運(yùn),1968年劉家峽水庫蓄水運(yùn)用,1969年開始的三門峽水庫第二期改建,2002年小浪底水庫建成完工等大的歷史事件,故人類活動(dòng)在很大程度上影響了水沙序列的變化,下游河道也因此產(chǎn)生很大影響。
表3 小浪底年水沙序列及花園口平灘流量序列在10a、23a時(shí)變異點(diǎn)
從1952年到2009年,下游來水來沙經(jīng)歷了3個(gè)階段:1960年以前的天然階段、1960—1999年三門峽水庫單獨(dú)運(yùn)用階段以及1999年以后小浪底水庫投入運(yùn)用。本文選取全時(shí)段進(jìn)行研究,先在多時(shí)間尺度層面對整個(gè)黃河下游水沙及平灘流量變化周期進(jìn)行了整體規(guī)律的研究,發(fā)現(xiàn)變化規(guī)律比較明顯,又在不同時(shí)間尺度上進(jìn)行了分析,從宏觀到微觀,總結(jié)了在不同自然條件及人類活動(dòng)影響下黃河下游的變化情況。通過對黃河下游小浪底、花園口、高村、利津4站水沙序列及花高利三站平灘流量序列進(jìn)行小波分析,發(fā)現(xiàn)了在多時(shí)間尺度下水沙與平灘流量的變換特征,結(jié)果表明:(1)黃河下游水量主要存在3a、9a、14a、23a變化周期,其中23a為最顯著周期;沙量主要存在5a、11a、23a、36a變化周期,其中11a為最顯著周期;平灘流量主要存在24a變化周期。水沙序列存在相似的多時(shí)間尺度變化特征,平灘流量序列受到水沙序列顯著周期的影響。水沙變化周期與華北、西北地區(qū)降水變化周期基本吻合,也可能與太陽黑子變化有關(guān)。(2)在小尺度(3a)下,年水量、年沙量、平灘流量曲線變化基本一致,在某些時(shí)段發(fā)生紊亂;在中尺度(10a)下,年水量、年沙量、平灘流量曲線變化有一定相關(guān)性及同步性,相位基本一致;而大尺度(23a)下,年水量、平灘流量曲線變化基本一致,平灘流量與汛期流量變化相關(guān)性更強(qiáng),與年沙量序列相位相反。各尺度下水量變化、沙量變化與平灘流量變化有所不一致可能是由于受到人類活動(dòng)及自然情況的影響,也可能是自然規(guī)律的反映,需進(jìn)一步研究。(3)各尺度下各序列的突變點(diǎn)基本一致,隨著時(shí)間尺度的增大,小時(shí)間尺度上的一些突變點(diǎn)逐漸退化為平常點(diǎn),相反,隨著時(shí)間尺度的降低,除原有的突變點(diǎn)外,還不斷增加突變點(diǎn),大時(shí)間尺度上的一些平常態(tài),在低層次看來則是突變態(tài)。因此,突變點(diǎn)與時(shí)間尺度有密切關(guān)系。黃河下游水沙及平灘流量序列的突變點(diǎn)主要為1957、1968、1980、1991、2003年左右。
[1]王文圣,丁晶,李躍清.水文小波分析[M].北京:化學(xué)工業(yè)出版社,2005.
[2]王文圣,丁晶,向紅蓮.水文時(shí)間序列多時(shí)間尺度分析的小波變換法[J].四川大學(xué)學(xué)報(bào):工程科學(xué)版,2002,34(6):14-17 .
[3]P-Kumar,E-Foufoula-Georgiou.A multi-component decomposition of spatial rainfall fields 1.Segregation of large and small scale features using wavelet transforms[J].Water Resources Research.1993,29(8):2515-2532.
[4]蔣曉輝,劉昌明,黃強(qiáng).黃河上中游天然徑流多時(shí)間尺度變化及動(dòng)因分析[J].自然資源學(xué)報(bào),2003,18(2):142-147.
[5]劉俊萍,田峰巍,黃強(qiáng).黃河上游河川徑流變化多時(shí)間尺度分析[J].應(yīng)用科學(xué)學(xué)報(bào),2003,21(2):117-121.
[6]王秀杰,練繼建.近43年黃河上游來水來沙變化特點(diǎn)[J].干旱區(qū)研究,2008,25(3):342-347.
[7]姜世中,梁川.龍川江年輸沙率時(shí)間序列的小波特征[J].哈爾濱工業(yè)大學(xué)學(xué)報(bào),2009,41(11):197-200.
[8]高鵬,穆興民,王煒.基于Morlet小波的太陽活動(dòng)與黃河輸沙量關(guān)系分析[J].水土保持研究,2010,17(3):181-183 .