張建偉,王根,張華,黃靜,陳靖,吳玲玲
(1.南京信息工程大學數理學院,江蘇南京210044;2.國家氣象中心,北京100081;3.天津市氣象科學研究所,天津300074)
星載大氣垂直探測器遙感資料在數值預報和氣候研究中起著重要的作用。以前由于技術的限制,一般的衛(wèi)星探測儀器只能提供數量有限的通道,在這種情況下通道選擇的原則是通過利用通道的光譜響應特性(如通道的中心頻率、帶寬等),控制通道的權重函數,使得觀測結果能最有效地對大氣參數進行反演,此時求解的問題是欠定的。
近年來,隨著探測技術的發(fā)展,越來越多的高光譜探測器被搭載在氣象衛(wèi)星上。此時求解的問題是超定的。例如,地球觀測系統(tǒng)(EOS)第二顆衛(wèi)星Aqua上攜帶的大氣紅外探測儀A IRS,采用紅外光柵陣分光術按NASA設計要求,2 378個通道(實際在使用2 378個通道的過程中,已經經過了通道選擇,選出了324個通道,以下只考慮這324個通道)覆蓋650~2 700 cm-1(3.7~15.4μm)紅外譜區(qū)域,其光譜分辨率v/Δv高于1 200,掃描寬度約1 650km,星下點分辨率約為13km,垂直分辨率為1km,輻射絕對精度優(yōu)于0.2K,溫度垂直探測精度在l km垂直分辨率下可達1K,而濕度垂直探測精度在2km垂直分辨率下可達10%。真正實現了高精度探測和高光譜分辨率,被用來探測精細的大氣溫度、濕度廓線等。
與A TOVS僅有20個通道相比,A IRS在實際使用324個通道仍顯太多,由此引起新的問題,一個是計算量太大,另一個是通道之間存在相關性。計算量太大不利于業(yè)務應用。通道之間的相關性會造成反演或同化的不適定性。如何盡可能多地提取有用信息,去除通道間的相關性,考慮到實際業(yè)務中需要從324個通道中選出一些對特定目標(例如,溫度和濕度)起主要作用的通道子集。
對高光譜紅外探測器的通道選擇國內外學者做了大量富有成效的工作。Rodgers(1996)提出用基于信息熵的分步迭代法(以下簡記為分步迭代法)進行通道選擇。Florence et al.(2002)實現了Rodgers提出的分步迭代法進行通道選擇的方法,并且提出了雅可比矩陣法,對比了這兩種方法得出雅可比方法比迭代方法省時,但迭代方法的效果較好。最后綜合了兩者的優(yōu)點,提出了常量迭代法,常量迭代法的目的是同一通道子集可用于不同的參考廓線,從而其實用性較強。Fourrie et al.(2002)提出手動選取通道的思想,主要采用實時數據評估了Florence et al.(2002)提出的常量迭代法和自己提出的手動選取通道的方法,得出手動選取通道更符合實際情況。Andrew et al.(2003)和Jam es et al.(2005)從同化A IRS觀測亮溫的角度出發(fā),采用分步迭代法,給出了所選的通道組合,白天48,夜晚62個通道(他們文章中給出的是白天45,夜晚60,后來有所改變,與Jam es Cameron博士私下交流)。
杜華棟等(2008)和張水平(2009)提出的通道選擇方法的核心思想是借用Florence et al.(2002)的基于信息熵的分步迭代法,劉輝(2008)采用的是基于權重函數所對應的峰值層進行通道選擇。
國內外通行的方法是基于信息熵的分步迭代法,此方法每次都是選信息量較大的通道,且前面選出的通道影響到后面的通道選擇,但其存在的主要缺陷是在用迭代法進行通道選擇時,依賴于背景誤差和觀測誤差協方差矩陣的給定。然而,背景誤差協方差的具體信息缺乏;分步迭代法每次是從未被選中的通道中選出一個,往往比較耗時。Florence et al.(2002)提出的雅可比矩陣法,Filipe et al.(2002)給出了具體的實施過程,首先,用背景誤差和觀測誤差協方差矩陣對其雅可比矩陣進行標準化,然后,類似于權重函數的思想選每層峰值最大時所對應的通道。雅可比方法也依賴于背景誤差和觀測誤差協方差矩陣的給定。在具體選取通道的過程中,只考慮了通道峰值最大貢獻層,沒有考慮通道之間的相關性。
針對上面方法的缺陷,本文提出基于主成分累計影響系數的通道選擇方法。主成分累計影響系數法,既能保住最主要信息,又能得到影響及敏感性較大的通道,且去除了通道之間的相關性;直接對Jacobi矩陣的結構進行分析,不需要過多的先驗知識;可以同時考慮多個通道。
變分同化方法可以有效克服單純反演計算不適定的困難,實現用正演方法求解反演問題,從方法論上避開了反演問題的復雜性。
變分同化方法的基本思想(張華等,2004)是將資料同化歸結為一個表征分析場與觀測場、分析場與背景場偏差的二次泛函極小值問題。該泛函(目標函數)一般定義為
其中:xb為背景場或初始猜測值;B是背景誤差協方差矩陣;yobs為實際的觀測亮溫;Robs為觀測誤差協方差矩陣;H為觀測算子及文中下面介紹的輻射傳輸模式。
變分同化問題可以歸結為一個目標函數的極小化問題,根據最優(yōu)化理論,通過調用某種下降算法來獲得最優(yōu)解xa,即由變分同化后得到的分析場,xa的求解實現了反演和同化的一體化。是Jacobi矩陣,表示各通道觀測或模擬亮溫相對待反演大氣參數溫度、濕度等的敏感性。又分別稱為溫度Jacobi矩陣、濕度Jacobi矩陣等其他。
對于上述所提及的輻射傳輸模式H,文中采用RTTOV模式(Saunders,2002)進行計算。此模式在垂直方向上,從0.1到1 013.3hPa,共分為43層。RTTOV可考慮大氣溫度、濕度廓線、云廓線、雨廓線參數的影響。除正模式外,同時還具有切線性模式、伴隨模式以及Jacobi矩陣模式。Jacobi矩陣結構如下:
其中:m在文中為通道個數;n為模式層數;T為溫度;q為濕度。
主成分分析是通過構造原變量的線性組合,產生一系列互不相關的新變量,從中選出少數幾個新變量使他們含有盡可能多的原變量信息,從而以較少新變量代替原來較多變量,消除信息冗余,實現模型的簡化。也可以利用主成分分析計算原變量的重要性系數進行排序,達到篩選變量降低模型維數的目的。本文用后一種思想,采用累計影響系數法進行通道選擇。
PCA數學模型:主成分分析把原始數據(施能,2002)X=Xj=(xij),i=1,2,…,n;j=1,2,…,m映射到新的變量空間Y=Y=(yij),i=1,2,…,n;j=1,2,…,m。其中
U=(ukj),k=1,2,…,m;j=1,2,…,m,U的行向量是原始數據X的協方差矩陣的特征向量。n在文中為模式層數,m為通道個數。
如果系數ukj滿足2,…,m;而且系數ukj的確定要使得Yi與Yj(i≠j)相互無關,并使Y1是X1,X2,…,Xm的一切線性組合中方差最大者,Y2是與Y1不相關的X1,X2,…,Xm的所有線性組合中方差次大者,……,Ym是與Y1,Y2,…,Ym-1都不相關的X1,X2,…,Xm的所有線性組合中方差第m大者,則稱Y1,Y2,…,Ym為原變量的第1,第2,…,第m主成分。
方差貢獻率定義為各主成分相應的特征值占總特征值之和的比例即為該主成分的方差貢獻率。
影響系數s定義為主成分的特征向量與各主成分的方差貢獻率的乘積,則累計影響系數表示為
sjm為第j個變量對第m個主成分的影響系數。Sj為第j個變量在m個主成分中的影響系數的累計,稱為累計影響系數。
剔除噪聲較大和衡量氣體的強敏感通道,考慮非局部熱動平衡(non-L TE)效應的影響,波數在2 220cm-1和2 287cm-1之間的通道全部剔除,考慮到使用模式頂層的不確定性和同化系統(tǒng)的時效性,剔除長波區(qū)域CO2和水汽的高層通道,參考ECMW F在使用A IRS時給出的黑名單。其他限制,短波通道(λ≤5μm)由于在白天受到太陽光的影響,所以文中分別考慮白天和夜晚的通道組合。
下面給出在溫度、濕度變分同化反演過程中,通道選擇的具體實施方案。首先進行通道的預處理,然后采用主成分累計影響系數法選取通道。
首先剔除黑名單中的通道,下面的步驟只針對剩余通道進行操作。
第1步:先考慮白天的通道組合。
第2步:初始只進行溫度分析,為了確保溫度信息主要來自CO2通道,屏蔽掉受水汽影響和黑名單中的通道。這樣確保溫度的最小信息都是來自CO2通道,而不是其他通道。受太陽光影響的通道被排除。
第3步:執(zhí)行第2步,對CO2通道所對應的溫度雅可比進行主成分累計影響系數分析。
主成分累計影響系數分析具體操作如下:
1)取(2)式中CO2通道所對應的溫度Jacobi矩陣,把H(x)各通道變量按下式
標準化得到標準化矩陣,記為:RPCA=(rji),i=1,2,…,n;j=1,2,…,m。其中:rji為標準化矩陣的元素;n在文中為模式層數;m為通道個數。
2)計算RPCA的協方差矩陣記為Cov(RPCA)。
3)計算Cov(RPCA)的特征值λj和特征矩陣,特征值λj按從大到小排列。
4)對Cov(RPCA)進行主成分分析,取前p個主成分。
挑選前p個主成分的原則:給定一個確定的百分比ε(本文取ε等于90%),選擇p使
代表前p個主成分占總方差的百分比,也就是p個主成分的總的貢獻率。其中λj表示第j個最大特征值。
5)由(3)式得到新的矢量,根據(6)式,取前面幾個影響較大的主要成分,進入累計影響系數的分析。
6)如果有相同的通道在不同的主成分中都被選入,則采用式(4)計算累計影響系數。否則,按照計算影響系數s的定義計算影響系數,最后進行影響系數的排序,選出對溫度影響較大的通道,當所選的通道數目達到初始給定的通道數目,或某個通道累計影響系數小于給定的閾值時,則停止此時的通道選擇,轉到第4步。
第4步:考慮水汽吸收帶附近的通道,對H2O通道所對應的濕度雅可比矩陣進行主成分累計影響系數分析,具體過程和第3步溫度雅可比矩陣分析類似。由于水汽通道實質對溫度的垂直結構也敏感,所以此步驟包括了溫度和濕度信息。
第5步:考慮夜晚的通道組合。
第6步:加入受太陽光影響的通道,此時則可以考慮從受到非局部熱動平衡(non-L TE)效應影響的通道中選一部分(原始黑名單被放寬)。
第7步:如果要進一步考慮地表信息,則可以加入窗區(qū)通道。
第8步:當所選的通道總數達到初始給定的總數或通道的累計影響系數較小時,則停止整個通道選擇過程。最終得到白天和夜晚的通道組合。
在進行通道選擇試驗時,為了說明本文方法可行性,將本文方法和目前公認比較好的基于信息熵分步迭代法的通道選擇進行比較。分步迭代法的具體細節(jié),可參閱Florence et al.(2002)所做的工作。
本文和分步迭代法的方法都是采用M et O ffice 1D-V ar(一維變分同化反演)溫度反演范圍為RTTOV模式分層的整層(43層,從0.1hPa到1 013.3 hPa)。由于在大氣層頂部濕度為零,因此本文在進行大氣濕度廓線反演時,反演范圍為RTTOV模式分層中的第18層到第43層(從122.04hPa到1 013.25hPa),而對于122.04hPa以上的大氣濕度不進行反演。對于本文研究的A IRS紅外觀測,當同一通道從不同觀測角度對大氣進行探測時,受臨邊效應的影響,Jacobi矩陣H的值會發(fā)生變化。因此,同一儀器針對不同的觀測角度進行通道選擇時將會產生不同的結果,但本文著重方法的研究,所有的觀測選為星下點,方位角為零。對于參考大氣廓線,本文采用RTTOV模式自帶的中緯度冬季廓線。
對通道個數的選取,不同文獻給的通道數目不同,本文是為后續(xù)進行A IRS觀測亮溫同化做探索試驗的,所以參考A ndrew et al.(2003)和Jam es et al.(2005)給的通道數目,白天48,夜晚62。
首先,考慮到地表信息的不確定性,本文先剔除所有峰值在地表的通道(進一步規(guī)范黑名單),選出來的通道組合如圖1所示,深黃色星號和品紅色點號分別表示分步迭代法和累計影響系數法所選出的通道分布。其他不同的顏色代表不同的氣體吸收帶或窗區(qū)通道。
對所選出的通道進行分析,本文的方法和分步迭代法選出的通道不同,在長波CO2區(qū)域,兩種方法都有通道選入;在窗區(qū),分步迭代法有通道選入,而本文的方法沒有通道選入的原因是考慮地表信息的不確定性,文中一開始就剔除了權重峰值在地表的所有通道。本文方法在水汽的強吸收帶附近有通道選入,原因相比較于分步迭代法,累計影響系數法更好的得到了通道的信息量(水汽強吸收帶附近,濕度Jacobi矩陣的幅度值更大);在短波區(qū)域,文中剔除了N2O吸收帶處的通道。
圖1 分步迭代法和累計影響系數法所選的通道分布(剔除所有峰值在地表的通道)Fig.1 Distribution of the channels selected through the iterative method and cumulativeeffect coefficient method(without channels whose peak values appear on the surface)
其次,考慮到地表信息在后期同化過程中起一定的作用,保留峰值在地表的通道,選出來的通道組合如圖2所示。
圖2 分步迭代法和累計影響系數法所選的通道分布(保留峰值在地表的通道)Fig.2 Distribution of the channels selected through the iterative method and cumulativeeffect coefficient method(with channels w hose peak values appear on the surface)
從圖2可以看出,本文方法和分步迭代法選出的通道大體分布相同,但在CO2和H2O的強吸收帶處入選的通道很少或者沒有通道入選,原因是在先前的預處理過程中剔除了部分吸收帶處的通道(剔除原因是發(fā)現有些通道的雅可比敏感線在對流層和平流層都出現了峰值從而導致平流層的信息影響到對流層,造成反演或同化的不適定性),剩余的通道的信息量相對較少,從而未被選入。
采用分步迭代法和主成分累計影響系數法,對中緯度冬季廓線進行反演,用分析誤差標準差進行比較。
首先考慮剔除所有峰值在地表的通道得到的通道組合,反演溫度、濕度的分析誤差標準差如圖3所示。圖3a、c為白天(48個通道),圖3b、d為夜晚(62個通道)不同方法下的溫度和濕度分析誤差標準差對比。
從圖3a中可以看出,對于溫度的反演,白天通道組合累計影響系數法在700hPa以上的效果勝于分步迭代法。而在700hPa以下效果比分步迭代法差的原因是,先前的預處理過程中,剔除了權重函數峰值在地面的所有通道。圖3b夜晚通道組合,效果和圖3a白天的相似。圖3c對于濕度反演,白天通道組合累計影響系數法在800hPa以上優(yōu)于分步迭代法。
圖3d對于濕度反演,夜晚通道組合累計影響系數法整體優(yōu)于分步迭代法。
其次,考慮保留峰值在地表的通道得到的通道組合,反演溫度、濕度的分析誤差標準差。如圖4所示。
從圖4a可見,對于溫度的反演,白天通道組合分步迭代法在200~300hPa的誤差小于累計影響系數法,但在600~900hPa累計影響系數法優(yōu)于分步迭代法。圖4b夜晚通道組合,分步迭代法在200~300hPa的誤差小于累計影響系數法,但在500~600hPa累計影響系數法優(yōu)于分步迭代法。從圖4c白天通道組合和圖4d夜晚通道組合中可以看出,對于濕度反演,累計影響系數法整體優(yōu)于分步迭代法。
相比較于分步迭代法依賴背景場誤差、觀測誤差協方差和由傳輸模式得到的雅可比矩陣,累計影響系數法只依賴于雅可比矩陣;分步迭代法每次只考慮一個通道,而累計影響系數法是多通道同時考慮。整體而言累計系數法優(yōu)于分步迭代法。
實際業(yè)務中只考慮A IRS的324個通道,但考慮到通道之間的相關性、同化的時效性(通道數目)以及同化或反演效果,必須進行通道選擇,選出主要的通道信息,本文提出了主成分累計影響系數法,主成分分析是對雅可比矩陣進行操作的,不依賴于過多的先驗知識,對累計影響系數進行排序,可達到通道選擇的目的。通道預處理可以首先剔除性能較差的通道,分別考慮白天和夜晚的通道組合,更能較好地利用受太陽光影響的短波通道。從理論分析和試驗效果圖可以看出此方法用于通道選擇是可行的。本文只是著重進行了A IRS通道選擇的理論研究,在后期的工作中考慮A IRS實際觀測亮溫資料同化對模式效果的影響。
圖3 白天(a,c)和夜晚(b,d)時溫度(a,b;單位:K)及濕度(c,d;單位:g/kg)的分析誤差標準差(剔除所有峰值在地表的通道)Fig.3 Analysis of(a,b)the temperature(K)and(c,d)hum idity(g/kg)error standard deviation during the(a,c)daytim e and(b,d)nighttim e(w ithout channels whose peak values appear on the surface)
致謝:Jam es Cam eron博士提供了資料和一些理論解釋,謹致謝忱!
圖4 白天(a,c)和夜晚(b,d)時溫度(a,b;單位:K)及濕度(c,d;單位:g/kg)的分析誤差標準差(保留峰值在地表的通道)Fig.4 Analysis of(a,b)the temperature(K)and(c,d)humidity(g/kg)error standard deviation during the(a,c)daytime and(b,d)nighttime(with channels whose peak values appear on the surface)
杜華棟,黃思訓,石漢青.2008.高光譜分辨率遙感資料通道最優(yōu)選擇方法及試驗[J].物理學報,57(12):7685-7692.
劉輝.2008.A IRS晴空大氣溫度廓線反演試驗[J].氣象學報,66(4):513-519.
施能.2002.氣象科研與預報中的多元分析方法[M].2版.北京:氣象出版社.
張華,薛紀善,莊世宇,等.2004.GRAPeS三維變分同化系統(tǒng)的理想試驗[J].氣象學報,62(1):31-41.
張水平.2009.A IRS資料反演大氣溫度廓線的通道選擇研究[J].氣象科學,29(4):475-481.
Andrew C,Saunders R,Cameron J,et al.2003.Assim ilation of data from A IRS for improved numerical weather prediction[C]∥Thirteenth Interational TOVS Study Conference.Adele,Canada.
Filipe A,Chedin A,Scott N A,et al.2002.A regularized neural net approach for retrieval of atmospheric and surface temperatures with the IASI instrument[J].J Appl Meteor,41(2):144-159.
Florence R,Fourrie N,Ourrie D,et al.2002.Channel selection methods for infrared atmospheric sounding interferometer radiances[J].Quart J Roy Meteor Soc,128(581):1011-1027.
Fourrie N,Jean-Noel thepaut research department.2002.Validation of the NESD IS near real time A IRS channel selection[C]∥ECMW F Tech.M emo.European Center for Medium Range.
James C,Andrew C,English S.2005.Operational use of A IRS observations at the met office[C]∥14th Internat Tovs Study Conference.Beijing,China.
Rodgers C D.1996.Information content and optimisation of high spec-tralresolution measurements[R]∥Paul B H,Wang Jinxue.Optical spectroscopic techniques and Instrumentation for atmospheric and space research II:136-147.
Saunders R.2002.RTTOV-7Users Guide[M].http:∥www.metoffice.com/research/interproj/nwpsaf/rtm.