葉晨，駱冬根，李怡心，姚萍萍，王羿，向光峰，李揚(yáng)，李雙，洪津

（1 中國科學(xué)院合肥物質(zhì)科學(xué)研究院 安徽光學(xué)精密機(jī)械研究所，合肥 230031）

（2 中國科學(xué)技術(shù)大學(xué)，合肥 230026）

（3 中國科學(xué)院通用光學(xué)定標(biāo)與表征技術(shù)重點(diǎn)實(shí)驗(yàn)室，合肥 230031）

（4 桂林電子科技大學(xué) 信息與通信學(xué)院，廣西 桂林 541004）

0 引言

中國科學(xué)院安徽光學(xué)精密機(jī)械研究所研制的星載偏振傳感器多角度偏振成像儀（Directional Polarimetric Camera，DPC）具有低畸變、大視場(chǎng)、畫幅式成像的特性，使其可獲得地氣系統(tǒng)的多角度、多光譜、偏振輻射信息［1-4］，將獲得的偏振信息和大氣特性反演模型相結(jié)合，可得到大氣氣溶膠和云的分布狀況及相關(guān)應(yīng)用參數(shù)［5-6］。DPC 載荷在出廠前，需要在實(shí)驗(yàn)室對(duì)DPC 響應(yīng)做穩(wěn)定性能測(cè)試［7-8］，實(shí)驗(yàn)要求光源的輻射能盡可能穩(wěn)定。輻射源通常采用鹵鎢燈積分球，周小麗等研制的近紅外光源在穩(wěn)定后10 min 以內(nèi)光強(qiáng)波動(dòng)約0.1%［9］；袁銀麟等研制的大孔徑可調(diào)光譜積分球光源2 h 內(nèi)的全譜段不穩(wěn)定性約為0.2%［10］；D'AMATO D 闡明了同一鹵鎢燈，其不同波段的相對(duì)光譜輻射率不同，且不穩(wěn)定性也不一致［11］。光源使用年限增加導(dǎo)致其性能衰變，結(jié)合環(huán)境溫濕度變化和操作規(guī)范性差異等因素影響，光源非穩(wěn)定性［12］參數(shù)會(huì)進(jìn)一步變大。因此有必要對(duì)光源進(jìn)行能量監(jiān)測(cè)，以降低多角度偏振成像儀的測(cè)量誤差，校正穩(wěn)定性測(cè)量參數(shù)［13］。本文通過對(duì)多角度偏振成像儀12 個(gè)通道的灰度值序列進(jìn)行小波分解，提取出各通道的光源能量變化量并在數(shù)據(jù)處理中進(jìn)行校正，以有效提升載荷穩(wěn)定性參數(shù)測(cè)量精度。

1 基本原理

1.1 多角度偏振成像儀和陷阱探測(cè)器觀測(cè)成像模型

搭載在某型號(hào)衛(wèi)星上的DPC 在一個(gè)探測(cè)周期內(nèi)可獲取3 個(gè)非偏波段（763 nm、765 nm、910 nm）和3 個(gè)偏振波段（490 nm、670 nm、865 nm）共12 個(gè)通道的圖像，每個(gè)偏振波段設(shè)置0°、60°、120°三個(gè)檢偏方向。剩余一個(gè)本底通道用于測(cè)量DPC 的暗電流數(shù)據(jù)，以消除由于探測(cè)器自身、電路盒溫度等因素造成的影響。

根據(jù)DPC 探測(cè)原理，其圖像數(shù)字量（Digital Number，DN）可表示為

式中，DNDPC為像元的DN 值，GDPC為DPC 的電子增益系數(shù)，λ1、λ2為DPC 探測(cè)波段的上下限，Ls(λ)為積分球光源的光譜輻亮度，Apixel為像元面積，F(xiàn)為等效光圈數(shù)，Moptic為光學(xué)系統(tǒng)放大倍率，Tatm(λ)為環(huán)境光譜透過率，Toptic-DPC(λ)為DPC 光學(xué)系統(tǒng)的光譜透過率，η(λ)為探測(cè)器量子效率，tint為積分時(shí)間，Ndark為暗電流電子數(shù)總和，Nnoise1為DPC 各種噪聲的總和。因此從式（1）可以看到主要影響DPC 長期穩(wěn)定性的是其中4 個(gè)參數(shù)：受積分球輻亮度穩(wěn)定性影響的Ls(λ)、受環(huán)境影響的Τatm(λ)、受探測(cè)器溫度影響的暗電流Ndark和探測(cè)器的其它噪聲Nnoise1。

陷阱探測(cè)器的探測(cè)原理可描述為

式中，s為陷阱探測(cè)器靈敏度，λ3、λ4為其探測(cè)帶寬上下限，對(duì)于同一款陷阱探測(cè)器來說，λ3、λ4是固定值，MTRAP、ATRAP分別為陷阱探測(cè)器光學(xué)系統(tǒng)參數(shù)、探測(cè)器面積，Noffset、Nnoise2分別為陷阱探測(cè)器的固定偏置、噪聲總和。因此陷阱探測(cè)器可直接測(cè)量獲得光源的穩(wěn)定性。

1.2 載荷非穩(wěn)定性評(píng)價(jià)

載荷的一項(xiàng)重要指標(biāo)是像元響應(yīng)不確定性，即載荷非穩(wěn)定性（Instability，NS），其計(jì)算公式為

式中，k為載荷在穩(wěn)定性實(shí)驗(yàn)期間響應(yīng)的總次數(shù)，y(i)為載荷第i次響應(yīng)時(shí)的輸出信號(hào)，為全時(shí)間段內(nèi)輸出信號(hào)的算術(shù)平均值。

1.3 小波分解原理

小波分解可由粗到精地對(duì)信號(hào)進(jìn)行多尺度分解，所以能夠很好地提取出信號(hào)中的各層低頻成分［14-16］，待分解信號(hào)通過低通濾波器得到低頻系數(shù)A，通過高通濾波器得到高頻系數(shù)D，如圖1 所示，信號(hào)經(jīng)第一層小波分解后得到A1和D1，然后再將第一層的低頻系數(shù)A1作為第二層的輸入，又得到一組低頻系數(shù)A2和高頻系數(shù)D2。以此類推分解，直到滿足想要的頻率要求。

圖1 信號(hào)小波分解過程Fig.1 Signal wavelet decomposition process

待分解的原始信號(hào)可看做0 級(jí)低頻系數(shù)a0=(f0，f1，...，fn)，那么

式中，G為低通濾波器矩陣，H為高通濾波器矩陣，am、dm分別為當(dāng)前分解m層時(shí)的低頻系數(shù)、高頻系數(shù)，因此

式中，G*、H*為G、H的共軛矩陣。

小波基函數(shù)為

式中，φ(t)為基小波或母小波函數(shù)，經(jīng)過尺度因子a和平移因子b變換后的φa，b(t)統(tǒng)稱小波。

對(duì)于離散信號(hào)，有

2 各觀測(cè)波段光源變化量提取

2.1 實(shí)驗(yàn)數(shù)據(jù)的獲取

以搭載在某型號(hào)衛(wèi)星上的多角度偏振儀DPC 穩(wěn)定性測(cè)試實(shí)驗(yàn)為例。其穩(wěn)定性實(shí)驗(yàn)測(cè)試現(xiàn)場(chǎng)示意圖如圖2 所示，將陷阱探測(cè)器和DPC 平行放置于積分球的出光口，控制柜控制積分球的點(diǎn)燈數(shù)和電流，采用直徑為2.5 m 的鹵鎢燈IS2500-1000 型積分球，工作波段是350～2 500 nm，標(biāo)稱全譜非穩(wěn)定性為0.2%/h，所采用的SI19008 型陷阱探測(cè)器的輻射不確定度低于0.65%。先將DPC 和陷阱探測(cè)器分別正對(duì)積分球輻射源，打開積分球光源及DPC 進(jìn)行預(yù)熱至穩(wěn)定狀態(tài)，設(shè)置好相應(yīng)的儀器參數(shù)，各自記錄響應(yīng)值，根據(jù)陷阱探測(cè)器監(jiān)測(cè)的光源穩(wěn)定段所對(duì)應(yīng)的DPC 成像情況，經(jīng)數(shù)據(jù)處理判斷DPC 響應(yīng)是否穩(wěn)定。

圖2 DPC 穩(wěn)定性實(shí)驗(yàn)示意圖Fig.2 Schematic diagram of DPC stability experiment

DPC 的同一通道在6 595 s 里采了1 024 次，采樣頻率為0.155 Hz，陷阱探測(cè)器監(jiān)測(cè)頻率是2 Hz，由Nyquist 采樣定律知：DPC 各通道的DN 幀序列能檢測(cè)到的最高頻率是0.077 5 Hz，陷阱探測(cè)器能識(shí)別光源變化最高的頻率是1 Hz。

在穩(wěn)定性實(shí)驗(yàn)中，重點(diǎn)關(guān)注的是鹵鎢燈積分球的低頻輻射變化信號(hào)［17］，而陷阱探測(cè)器在監(jiān)測(cè)時(shí)引入的噪聲及儀器測(cè)量誤差主要是相對(duì)高頻噪聲，因此利用小波分解可以從含有高頻噪聲的曲線中提取出低頻的光源緩變信號(hào)。同時(shí)，在DPC 的CCD 記錄數(shù)字量灰度值DN 值過程中，CCD 的光電響應(yīng)的噪聲較光源能量變化屬于高頻成分，陷阱探測(cè)器的監(jiān)測(cè)靈敏度變化等噪聲較光源本身的能量變化也屬于高頻成分，因此光源能量的低頻變化曲線是暗藏在DN 值幀序列和陷阱探測(cè)器監(jiān)測(cè)曲線中的。又由于陷阱探測(cè)器的采樣頻率比CCD 同一通道的采樣頻率高，因此陷阱探測(cè)器中蘊(yùn)含的光源變化信息更多，但是穩(wěn)定性實(shí)驗(yàn)只要求光源能量穩(wěn)定在一定范圍內(nèi)，因此理論上可從DN 值幀序列中提取出光源變化趨勢(shì)。

在穩(wěn)定性實(shí)驗(yàn)中，積分球輻射源在DPC 的CCD 上所成的某一幀圖像如圖3 所示，由于觀測(cè)位置和條件均保持不變，因此每幀的被輻射區(qū)域位置是固定的，其內(nèi)接正方形尺寸約是160 pixels×160 pixels?？紤]用多像元平均法提高數(shù)據(jù)信噪比，因此本次選取的像元監(jiān)測(cè)區(qū)域?yàn)橹醒雲(yún)^(qū)域140 pixels×140 pixels。DPC 采用的是幀轉(zhuǎn)移型面陣CCD，因此圖像的預(yù)處理包括扣除當(dāng)前幀及前后各3 幀本底平均和幀轉(zhuǎn)移效應(yīng)校正。

圖3 DPC 某一通道接收到的一幀圖像Fig.3 A frame of image received by one channel of DPC

2.2 光源變化曲線的提取

基于像元的光子信號(hào)和轉(zhuǎn)換過程中產(chǎn)生的噪聲變換特性，本次小波基選取的是db2 小波，分解層數(shù)為8，從陷阱探測(cè)器監(jiān)測(cè)曲線中提取出光源的低頻變化趨勢(shì)A1～A8，也從各通道的DN 值幀曲線中提取出A1～A8，將兩者之間的低頻曲線做Pearson 相關(guān)，取Pearson 相關(guān)系數(shù)最大所對(duì)應(yīng)的低頻曲線，認(rèn)為是各通道的光能量變化趨勢(shì)，提取光源變化曲線流程如圖4 所示。Pearson 相關(guān)系數(shù)是用來反映兩變量之間相似程度的統(tǒng)計(jì)量，衡量兩者之間的線性關(guān)聯(lián)性，其對(duì)兩組數(shù)據(jù)的異常值較為敏感，即數(shù)據(jù)的異常點(diǎn)會(huì)導(dǎo)致兩組數(shù)據(jù)間的Pearson 相關(guān)系數(shù)降低。

圖4 提取光源變化曲線流程Fig.4 Process of extracting light source change curve

依次對(duì)12 個(gè)通道的DN 值幀曲線做8 層小波分解，將每個(gè)通道每層的低頻曲線An和陷阱探測(cè)器的An計(jì)算Pearson 相關(guān)系數(shù)。以670 nm（P1）通道為例，此通道的DN 值幀曲線的A1～A8和對(duì)應(yīng)陷阱探測(cè)器的A1～A8之間的Pearson 相關(guān)系數(shù)如圖5 所示，DPC 的670（P1）、670（P2）、670（P3）之間的Pearson 相關(guān)系數(shù)如圖6 所示。

圖5 670 nm（P1）通道的DN 值幀序列的A1～A8和陷阱探測(cè)器的A1～A8的相關(guān)系數(shù)Fig.5 Correlation coefficient between A1～A8 of trap detector and A1～A8 of DPC 670 nm（P1）

圖6 670 nm（P1）、670 nm（P2）、670 nm（P3）的DN 值幀序列的A1～A8之間的Pearson 相關(guān)系數(shù)Fig.6 Pearson correlation coefficient between A1～A8 of DN value frame sequence of 670 nm（P1）、670 nm（P2）、670 nm（P3）channel

由圖5 可看出，對(duì)于670 nm（P1）而言，在陷阱探測(cè)器的A1～A8曲線中，Pearson 相關(guān)系數(shù)最大值均是在DN 值幀序列的A5處取得，且均在0.95 以上，說明670 nm（P1）和光源的670 nm 波段輻射能的變化極度相關(guān)。由圖6 可看出670 nm 的三個(gè)偏振通道的DN 值幀序列在相同的小波分解層數(shù)下（即對(duì)角線方向），其Pearson 相關(guān)系數(shù)取得最大值，且在A5時(shí)的兩兩之間的Pearson 相關(guān)系數(shù)接近0.998，說明三個(gè)偏振通道的DN 值幀序列變化趨勢(shì)一致，是由光源670 nm 波段的輻射能變化造成的。因此對(duì)670 nm 波段而言，其對(duì)應(yīng)的DN 值幀序列的A5最能體現(xiàn)光源此波段的變化趨勢(shì)。同理計(jì)算得到其余8 個(gè)成像通道的光源變化曲線是各自通道DN 值幀序列的A4或A5，說明在整個(gè)穩(wěn)定性觀測(cè)實(shí)驗(yàn)期間，光源的光輻射能主要以低于0.004 8 Hz的頻率緩慢變化。因此，每個(gè)通道的DN 值幀序列小波分解的A4或A5即為實(shí)驗(yàn)期間光源變化量。

以測(cè)量的初始時(shí)刻的光能量為基準(zhǔn)，各通道光源能量變化趨勢(shì)和變化量分別如圖7 和表1 所示。光源各波段的能量穩(wěn)定性不一致，其中490 nm、910 nm 和865 nm 波段變化超過0.5%，因此有必要對(duì)DPC 的觀測(cè)數(shù)據(jù)進(jìn)行去光源變化處理，以提高穩(wěn)定性參數(shù)的真實(shí)性和有效性。

圖7 各通道的光源能量變化Fig.7 Variation of light source energy in all channels

表1 各波段的光能量變化Table 1 Variation of light energy in each wavelength band

利用所提取的各波段光源能量變化情況進(jìn)行穩(wěn)定性校正，得到利用陷阱探測(cè)器監(jiān)測(cè)方法和小波分解方法前后的穩(wěn)定性參數(shù)計(jì)算結(jié)果如表2 所示，可以看到兩種方法獲得的穩(wěn)定性參數(shù)計(jì)算結(jié)果均有一定程度的提升，而小波分解方法獲得的穩(wěn)定性參數(shù)的提升更為明顯，較好地去除了光源自身的波動(dòng)、環(huán)境光譜透過率不穩(wěn)定性、DPC 自身溫度變化導(dǎo)致的暗電流不穩(wěn)定性，計(jì)算出的不穩(wěn)定性更接近DPC 實(shí)際工作時(shí)的不穩(wěn)定性。

表2 兩種校正方法的非穩(wěn)定性參數(shù)計(jì)算結(jié)果對(duì)比Table 2 Comparison of calculation results of the instability parameter by two correction methods

3 基于小波分解方法提取穩(wěn)定性參數(shù)的驗(yàn)證

3.1 實(shí)驗(yàn)室信噪比計(jì)算方法

1）多幀圖像數(shù)據(jù)信噪比計(jì)算法

CCD 同一通道接收到了k幀圖像，在被輻照區(qū)域內(nèi)選取計(jì)算區(qū)域m×n，對(duì)于計(jì)算區(qū)域內(nèi)的某一個(gè)像元(i，j)，在k幀圖像中，該像素點(diǎn)對(duì)應(yīng)k個(gè)DN 值，求其平均值μ(i，j)及標(biāo)準(zhǔn)差std(i，j)，定義該像元的信噪比為

以此類推，可求得計(jì)算區(qū)域內(nèi)其它像元的信噪比。將m×n個(gè)像元的信噪比SNR(i，j)取平均，作為CCD此通道接收到的圖像的信噪比SNR，即

2）相鄰兩幀信噪比計(jì)算法

歐洲機(jī)器視覺協(xié)會(huì)制定的圖像傳感器與相機(jī)性能測(cè)試標(biāo)準(zhǔn)EMVA Standard 1288［18］中提出用相鄰兩幀計(jì)算信噪比，取CCD 相鄰兩幀圖像A、B，假設(shè)噪聲源靜止且具有均勻性，則在被輻照區(qū)域內(nèi)選取計(jì)算區(qū)域m×n，定義平均灰度值即信號(hào)為

式中，yA(i，j)表示圖像A中坐標(biāo)(i，j)處的像元灰度值，yB(i，j)表示圖像B中坐標(biāo)(i，j)處的像元灰度值。

兩圖像A、B相減，得到新圖像C，其m×n個(gè)像元灰度值的標(biāo)準(zhǔn)差記為stdC，則原圖像A或B的噪聲大小視為

3.2 陷阱探測(cè)器修正光源計(jì)算信噪比

在穩(wěn)定性實(shí)驗(yàn)中，常用陷阱探測(cè)器監(jiān)測(cè)光源變化曲線來對(duì)DPC 所成的圖像進(jìn)行去光源變化，將當(dāng)前圖像的建立時(shí)間與所對(duì)應(yīng)的光源曲線的點(diǎn)對(duì)應(yīng)起來，然后再去光源變化，即可計(jì)算信噪比。設(shè)某一通道的第1幀、第2 幀、...、第k幀圖像建立時(shí)對(duì)應(yīng)的陷阱探測(cè)器監(jiān)測(cè)的光輻射能依次為Q1、Q2、...、Qk，將其光能量歸一化處理，則此通道第i幀的光能量修正系數(shù)為

第i幀圖像修正后的DNTRAP_corrected_i值為

將每幀圖像根據(jù)陷阱探測(cè)器的光能量進(jìn)行DN 值修正后，計(jì)算出的信噪比記為陷阱探測(cè)器修正后信噪比。

3.3 基于小波分解法修正光源計(jì)算信噪比

為證明提出的用小波分解法提取光源的有效性，用基于小波分解法提取出12 個(gè)成像通道各自的光源變化曲線，依據(jù)光源變化修正圖像DN 值，得到修正后的幀序列信噪比以及信噪比提升比例分別見圖8 和9，具體提升的信噪比數(shù)值比較見表3。由圖8 可看出，小波分解法修正光源后計(jì)算出的信噪比接近相鄰兩幀信噪比，證明提取光源方法的正確性。

表3 陷阱探測(cè)器去光源變化法和小波分解法去光源變化法提升的信噪比數(shù)值比較Table 3 Numerical comparison of the improved SNR between trap detector de-illumination variation method and wavelet decomposition de-illumination variation method

圖8 不同的方法計(jì)算的信噪比Fig.8 SNR calculated by different methods

陷阱探測(cè)器去光源修正后較未修正前信噪比平均提升了0.72%，而小波分解法修正后信噪比平均提升了2.52%。在865 nm 和910 nm 波段，小波分解去光源法將信噪比分別提升了4.34%和12.9%，遠(yuǎn)高于陷阱探測(cè)器監(jiān)測(cè)數(shù)據(jù)修正的結(jié)果，反映了實(shí)驗(yàn)期間光源在865 nm 波段和910 nm 波段能量波動(dòng)較大，也說明用寬譜的陷阱探測(cè)器監(jiān)測(cè)數(shù)據(jù)來修正865 nm 和910 nm 波段的變化存在局限。在490 nm 波段，小波分解法將信噪比提升了0.54%，但陷阱探測(cè)器修正法僅將信噪比提升了0.09%，主要是因?yàn)?90 nm 波段的光源能量較低，獲取的圖像DN 值低。雖然光源在490 nm 波段的變化最大，但信噪比提升并沒有910 nm 波段顯著。兩種方法在670 nm 波段信噪比提升幅度小則是因?yàn)榇瞬ǘ蔚墓庠醋兓∮?.205%。

圖9 陷阱探測(cè)器去光源變化法和小波分解法去光源變化法提升的信噪比Fig.9 The improved SNR of the trap detector de-illumination variation method and the wavelet decomposition deillumination variation method

4 結(jié)論

積分球光源不同光譜波段的變化趨勢(shì)不同，且與光源的使用年限、環(huán)境控制、操作規(guī)范等因素有一定的相關(guān)性，但總的來說，大口徑積分球光源的變化趨勢(shì)較緩慢。本文提出用小波分解法從DPC 各成像通道的DN 值幀曲線中提取出光源各波段的能量變化趨勢(shì)，進(jìn)而扣除光源變化趨勢(shì)，使穩(wěn)定性測(cè)試的結(jié)果更為客觀和真實(shí)。光源變化校正的效果優(yōu)于傳統(tǒng)陷阱探測(cè)器監(jiān)測(cè)數(shù)據(jù)校正，DPC 的信噪比較未扣除光源波動(dòng)前平均提升了2.52%，非穩(wěn)定性降低到0.031%，驗(yàn)證了此方法提取光源變化量的有效性。