EMCCD真彩色成像方法與驗證

白 玉，陳遠(yuǎn)金,2，張猛蛟,2，王嶺雪，蔡 毅,3

?

EMCCD真彩色成像方法與驗證

白 玉1，陳遠(yuǎn)金1,2，張猛蛟1,2，王嶺雪1，蔡 毅1,3

（1. 北京理工大學(xué)光電學(xué)院 光電成像技術(shù)與系統(tǒng)教育部重點實驗室，北京 100081；2. 華東光電集成器件研究所，江蘇 蘇州 215163；3. 中國兵器科學(xué)研究院，北京 100089）

針對夜天光包含的大量近紅外輻射，本文從理論和實驗上研究了基于近紅外波段全透的紅、綠、藍(lán)濾光片的EMCCD真彩色成像方法。夜天光在近紅外波段的能量均被EMCCD吸收，再結(jié)合EMCCD的電子倍增功能獲得較高亮度和信噪比的紅、綠、藍(lán)分量圖像，最終通過彩色融合等圖像處理方法獲得真彩色EMCCD圖像。本文進(jìn)行了全月、半月、1/4月和晴朗星光照明條件下的真彩色成像色彩仿真，在此基礎(chǔ)上制作了彩色濾光片，測試了典型樣本的光譜反射率，搭建了真彩色EMCCD實驗裝置并開展外場試驗。工作結(jié)果表明：該方法獲得了色彩豐富、空間深度感較好、信噪比較高的真彩色EMCCD圖像。

真彩色EMCCD；夜天光；真彩色成像；近紅外；彩色濾光片

0 前言

基于電子倍增技術(shù)的低照度電荷耦合器件（Electron-Multiplying CCD, EMCCD）具有靈敏度高、信號增益高、噪聲低、體積小、分辨力高、性能穩(wěn)定的特點，而且還有片上防飽和功能，可在全日光到星光的大范圍照度條件下連續(xù)工作，真正實現(xiàn)全天時使用。色彩是圖像的重要描繪子之一，使得場景分割、圖像識別等任務(wù)更加容易，Cavanillas等人[1]的研究表明：色彩能提高場景理解能力、縮短反應(yīng)時間、幫助人們更準(zhǔn)確和快速地識別目標(biāo)。因此，色彩被引入夜視成像技術(shù)并相繼發(fā)展了偽彩色和假彩色技術(shù)路線的彩色夜視技術(shù)。近年來，隨著夜視成像器件靈敏度的提高和微光學(xué)加工技術(shù)的發(fā)展，真彩色夜視成像逐漸成為研究熱點。真彩色圖像與人類的視覺認(rèn)知一致，有助于使用者利用大腦中存儲的色彩經(jīng)驗來識別物體。Tenebraex公司[2]使用像增強器和彩色濾光片旋轉(zhuǎn)掃描的方式在2012年開發(fā)了ColorPath真彩色夜視成像儀。SPI公司[3]宣稱使用固體成像器件和寬光譜薄膜陣列（Broad Spectrum Thin Film Array, BSTFA）技術(shù)開發(fā)了低照度彩色攝像機(jī)。

EMCCD真彩色成像技術(shù)也在持續(xù)發(fā)展。2006年，Ball航空公司[4]提出將彩色攝像機(jī)使用的經(jīng)典Bayer彩色濾光陣列（Color Filter Array, CFA）進(jìn)行改裝，通過低于4%稀疏分布的微型彩色濾光片和色彩渲染、空間濾波等圖像處理方法開發(fā)了彩色EMCCD攝像機(jī)。Ball的方法具有超高的光子利用率，圖像色調(diào)正確，能在晴朗的無月夜天空下（約1.0×10－3lx照度）成彩色圖像，但色彩的飽和度有損失，圖像顏色具有空間突變性。2010年，Jason M等人[5]研制了非制冷紅外熱像儀和真彩色EMCCD的彩色融合裝置，其中真彩色EMCCD采用青-黃-品-綠（Cyan-Yellow- Magenta-Green, CYMG）模式的彩色濾光陣列，在約1.0×10－2照度下能獲得真彩色圖像，但該彩色模式下由于C、Y、M、G通道的光譜響應(yīng)曲線存在較多的重疊區(qū)域，會帶來色域小、色彩豐富程度較低的問題，同時，在處理近紅外成分對色彩保真度的影響時，需要“額外的濾光片以及色彩優(yōu)化算法”。2015年，韓國Son等人[6]的實驗表明：采用RGB模式Bayer濾光陣列的普通彩色CCD在照度降到5 lx時就無法成像，而采用同樣濾光陣列的EMCCD能在0.3 lx照度下進(jìn)行彩色成像。在前人研究的基礎(chǔ)上，本文針對夜天光包含的大量近紅外輻射，研究基于近紅外波段全透的紅、綠、藍(lán)濾光片的EMCCD真彩色成像方法。夜天光在近紅外波段的能量均被EMCCD吸收，再結(jié)合EMCCD的電子倍增功能獲得較高亮度和信噪比的紅、綠、藍(lán)分量圖像，最終通過彩色融合等圖像處理方法獲得真彩色EMCCD圖像。

1 夜天光的亮度和光譜

在24h的一個地球自轉(zhuǎn)周期中，地面景物環(huán)境光照的強度和光譜成分都在發(fā)生巨大變化。太陽光的照度高達(dá)1.0×105lx，無月星光的照度則低至1.0×10－3lx。月光（Moonlight）、明亮行星光（Bright Planets）、黃道光（Zodiacal Light）、恒星光（Integrated Starlight）、大氣輝光（Airglow）、散射的銀河系光（Diffuse Galactic Light）和宇宙光（Cosmic Light）等共同組成夜天光。由于月光的亮度是星光和大氣輝光的100多倍，因此夜天空亮度受月相的影響很大[7]，表1是不同月相時地面景物的照度值[8]。

表1 不同月相時的地面景物的照度值

夜天光光譜的測量研究首先集中在可見光波段，例如1922年Meinel[9]在《自然》上發(fā)表論文，稱在夜天光中測量到黃綠波段的極光譜線以及H和K夫朗和費譜線；之后隨著探測器技術(shù)的發(fā)展，夜天光光譜的測量延伸至紅外波段，例如1953年A. Vallance Jones[10]使用硫化鉛紅外探測器測量了夜天光在1200～2000 nm波段的光譜分布情況。相比太陽光光譜，夜天光光譜更為復(fù)雜，至今仍沒有統(tǒng)一的表征方法。圖1給出了0.89月、0.50月、0.25月和無月星光下夜天光的光子輻射通量分布圖[11]，可以看出：①夜天光包含大量近紅外輻射（780～1100nm）；②月相從0.89月變?yōu)?.50月、0.25月、無月時，夜天光中的可見光（380～780nm）迅速減弱，而1000～1100 nm波段的近紅外輻射變化不明顯。

圖1 不同月相時夜天光的光子輻射通量分布圖

2 EMCCD真彩色成像方法

2.1 EMCCD真彩色成像原理

為獲得豐富的色彩，傳統(tǒng)真彩色成像系統(tǒng)要求彩色濾光片對近紅外截止，但是這種系統(tǒng)用于夜天光條件時，圖像信噪比低，甚至?xí)捎诠饽芰可俣鵁o法成像，另外沒有充分利用夜天光中的近紅外成分。針對這些情況，本文一是使用帶有電子倍增功能的EMCCD（光譜響應(yīng)范圍：380～1100nm），二是采用近紅外波段全透的紅、綠、藍(lán)濾光片，即R＋NIR、G＋NIR、B＋NIR濾光片，真彩色成像原理如圖2所示，透過3個濾光片的光能量到達(dá)EMCCD產(chǎn)生光電轉(zhuǎn)換后獲得R￠（R＋NIR1）、G￠（G＋NIR2）、B￠（B＋NIR3）通道的數(shù)據(jù)，這3個通道的數(shù)據(jù)再映射到RGB顏色空間作為紅、綠、藍(lán)分量進(jìn)行融合。經(jīng)反復(fù)分析，R￠、G￠和B￠通道的近紅外成分滿足NIR1＝NIR2＝NIR3時，EMCCD在3個通道的近紅外響應(yīng)不會影響景物的基本色調(diào)。

圖2 EMCCD真彩色成像原理框圖

2.2 EMCCD真彩色成像的濾光片設(shè)計

在設(shè)計EMCCD真彩色成像濾光片的透過率曲線時主要考慮以下兩點：

1）可見光波段的透過率符合傳統(tǒng)真彩色成像的要求。

為獲得較好的顏色復(fù)現(xiàn)（color reproduction）效果，傳統(tǒng)真彩色成像系統(tǒng)的濾光片透過率曲線是根據(jù)人眼3種彩色視覺細(xì)胞的光譜靈敏度曲線來設(shè)計的[12]，常用的紅、綠、藍(lán)濾光片的峰值波長分別為550nm、540nm、450nm，峰值半寬高分別為50nm、50nm、40nm[13]。

2）近紅外波段的透過率盡量高，以便充分利用夜天光的近紅外能量，而且近紅外波段的起始波長盡量重合。

基于以上考慮，本文制作了可見光波段具有帶通濾光性質(zhì)、近紅外波段具有長波通性質(zhì)的濾光片，在K9玻璃表面鍍制多種組合膜系達(dá)到相應(yīng)的光譜透過率，并在表面鍍制增透膜提高濾光片的透光能力。濾光片的物理參數(shù)如表2所示，光譜透過率曲線如圖3(a)所示，實物圖如圖3(b)所示。忽略大氣傳輸、光學(xué)鏡頭等因素的影響，EMCCD對紅、綠、藍(lán)通道的光譜靈敏度函數(shù)由濾光片的光譜透過率（圖3(a)）與EMCCD的相對光譜靈敏度曲線（圖3(c)）相乘得到，如圖3(d)所示。

3 EMCCD真彩色效果的色彩仿真

為比較本文方法在不同照度夜天光下的色彩還原能力，暫不考慮光源的能量、僅考慮光源的相對光譜功率分布，并假設(shè)光學(xué)鏡頭透過率為1，以標(biāo)準(zhǔn)24色麥克白色卡（Macbeth Color Checker chart）與6種典型物體為樣本，對本文方法的EMCCD真彩色效果進(jìn)行仿真。

3.1 樣品的光譜反射率的測量

本文使用標(biāo)準(zhǔn)24色麥克白彩色卡作為測試目標(biāo)樣本，該色卡采用特殊配置的涂料，能模擬大范圍的顏色和色調(diào)。除此之外，本文還使用“混凝土”、“干草”、“綠色涂料”、“紅色玫瑰”、“黃色玫瑰”、“綠葉（新鮮的）”6種實際物體進(jìn)行研究。新鮮的植物葉片由于對近紅外輻射有較高的光譜反射率，因此重點加以研究。標(biāo)準(zhǔn)24色麥克白色卡的24個色塊、“紅色玫瑰”、“黃色玫瑰”、“綠葉”總共27個樣品在380～1100nm波段的光譜反射率使用美國Perkin Elmer公司的Lambda 19分光光度計測量，測量前首先測量標(biāo)準(zhǔn)白板并用其數(shù)據(jù)校正其他顏色的數(shù)據(jù)。27種樣品的光譜反射率曲線如圖4(a)所示。“混凝土”、“干草”、“綠色涂料”3種物品的光譜反射率從加利福尼亞技術(shù)研究所的ASTER光譜庫[13]查到。6種實際物品的光譜反射率如圖4(b)所示。圖4(c)為24色麥克白色卡和6種實際物品在CIE D65光源、CIE1931標(biāo)準(zhǔn)觀察者條件下的色彩效果。

表2 本文彩色濾光片的物理參數(shù)

圖3 本文彩色濾光片透過率光譜和EMCCD響應(yīng)光譜

3.2 色彩仿真過程

令表示波長，()表示夜天光的相對光譜功率分布，()表示樣本的光譜反射率，()表示濾光片的光譜透過率，()表示EMCCD的歸一化光譜靈敏度，則反射率為()的樣本在()照明下經(jīng)過濾光片()濾光后，EMCCD輸出的R￠G￠B￠通道的信號由下式表示：

式中：1、2分別為積分的上限和下限波長；式中變量的下標(biāo)R￠、G￠、B￠分別表示R￠、G￠、B￠通道，常數(shù)R′、G′、B′分別是紅、綠、藍(lán)通道相互獨立的增益系數(shù)。

色彩仿真時所取的數(shù)據(jù)為：

1）()取自圖1，并對數(shù)據(jù)做關(guān)于550nm波長的歸一化處理；

2）()取自圖4(a)和圖4(b)樣品的光譜反射率；

3）R′、G′、B′分別取圖3(a)濾光片R＋NIR、G＋NIR、B＋NIR的光譜透過率；

圖4 樣品的顏色及光譜反射率曲線圖

Fig.4 Color and spectral reflectance curve of samples

4）1、2選（380nm，780nm）和（380nm，1100nm）兩組；

5）R′、G′、B′是R￠、G￠B￠通道相互獨立的增益系數(shù)，設(shè)置為固定值，以便系統(tǒng)在拍攝中性物體時R￠、G￠和B￠通道可以得到相等的信號值，即R′＝G′＝B′。本文對增益系數(shù)的求解方法是：針對麥克白色卡的第19個色塊（純白色），R￠、G￠和B￠通道的信號值取R′＝G′＝B′＝255，利用式(1)解出R′、G′、B′。

3.3 色彩仿真結(jié)果

圖5給出了0.89月、0.50月、0.25月和無月星光4種夜天光條件下真彩色EMCCD的色彩仿真效果,其中圖5(a)、(c)、(e)、(g)為380～780nm波段的仿真結(jié)果，代表了真彩色EMCCD所能達(dá)到的最理想的色彩效果（相當(dāng)于系統(tǒng)增加了近紅外截止濾光片）；圖5(b)、(d)、(f)、(h)為380～1100nm波段的仿真結(jié)果，代表真彩色EMCCD真實的色彩效果，反映了圖像色彩受近紅外的影響。

由圖5可以看出，380～780nm的仿真結(jié)果與標(biāo)準(zhǔn)24色卡（圖4(c)）對比，樣品的顏色接近真實，說明濾光片在可見光波段的參數(shù)設(shè)計合理；從380～1100nm的仿真結(jié)果中發(fā)現(xiàn)，EMCCD增加的近紅外響應(yīng)對顏色有“沖淡”作用，色彩的飽和度降低，這種影響隨著夜天光中月光照明的減弱（近紅外在整個響應(yīng)波段中占的比重增加）而增強，但樣品基本的色調(diào)（紅、綠、藍(lán)）是正確的。

綠色植物對近紅外輻射的反射率較高，其在夜天光條件下的R￠、G￠B￠通道信號主要是近紅外輻射的響應(yīng)，因此3個通道的信號值基本相等，所以綠色植物的顏色偏向于灰色。圖5(b)、5(d)、5(f)、5(h)中最后1行代表“新鮮綠葉”的最后1個色塊的顏色始終是灰色的，而在同一行中代表“綠色涂料”的第3個色塊在0.89月和0.50月時還能呈現(xiàn)綠色，這說明本文彩色成像方法在一定照度條件下還具有識別新鮮綠葉和偽裝目標(biāo)綠色涂料的能力。

圖5 本文真彩色EMCCD實驗裝置的色彩仿真結(jié)果圖

4 EMCCD真彩色成像方法的實驗驗證

4.1 EMCCD真彩色成像實驗裝置

本文設(shè)計的R＋NIR、G＋NIR、和B＋NIR濾光片加裝在濾光輪上，實驗時針對特定場景旋轉(zhuǎn)濾光輪一周，EMCCD相機(jī)分別獲得R￠、G￠、B￠通道的模擬視頻圖像，再由數(shù)字?jǐn)z像機(jī)（DV）顯示并存儲，實驗裝置示意圖見圖6。使用MATLAB軟件對R￠、G￠、B￠通道的圖像數(shù)據(jù)進(jìn)行融合得到真彩色圖像。

4.2 實驗結(jié)果與分析

圖7是實驗裝置全月夜晚在北京鳳凰嶺地區(qū)拍攝的圖像效果（圖7(a)）以及相同照度條件下的色彩仿真結(jié)果（圖7(c)），圖7(b)是圖7(a)中色卡的放大圖?？梢钥闯觯弘m然真實拍攝的圖像噪聲較大，但真彩色效果與仿真結(jié)果一致，這表明本文的真彩色EMCCD成像方法具有較好的色彩復(fù)現(xiàn)能力。

圖8～圖11給出了更多的室外拍攝效果。圖8拍攝于2016年7月6日21:30，地點在北京理工大學(xué)6號教學(xué)樓附近，照度：～10－1lx，共進(jìn)行了3組實驗，每一組實驗的圖(a)～(c)分別對應(yīng)R￠通道、G￠通道和B￠通道的圖像，圖(d)為本文EMCCD真彩色成像方法結(jié)果，圖(e)為白天的彩色圖像（使用普通彩色CCD在白天拍攝得到）。圖8中EMCCD拍攝R￠、G￠、B￠通道圖像時的增益相同。

圖9拍攝于2016年3月20日23:30，地點為北京理工大學(xué)中心花園，照度：～10－1lx，EMCCD拍攝R￠、G￠、B￠通道圖像時的增益有變化。圖10拍攝于2016年3月18日21:00，地點為北京理工大學(xué)6號教學(xué)樓附近，照度：～10–1lx，EMCCD拍攝R￠、G￠、B￠通道圖像時的增益有變化。圖11拍攝于2016年3月21日21:30，地點為北京郊外的鳳凰嶺，照度：～10–1lx，EMCCD拍攝R￠、G￠、B￠通道圖像時的增益有變化。

圖6 實驗裝置的結(jié)構(gòu)示意圖

Fig.6 Schematic diagram of experimental device

圖7 麥克白24色卡的色彩仿真結(jié)果以及實驗裝置拍攝效果圖

（時間：2016-07-06 21:30；地點：北京理工大學(xué)6號教學(xué)樓附近；照度：～10－1lx）

Fig.8 True-color imaging results in this paper and reference images taken by a color CCD in daytime

(Time: 2016-07-06 21:30; Location: No.6 Teaching Building of BIT;Illuminance: ～10－1lx)

從圖8～圖11的實驗效果可以看出：圖(c)的B￠通道圖像最暗、噪聲最大，說明夜晚藍(lán)色波段的能量最低；相比較于圖(a)～圖(c)的灰度圖像，圖(d)的本文真彩色結(jié)果具有較豐富的紋理細(xì)節(jié)、較高的圖像亮度和信噪比，色彩也較為豐富?？偟膩碚f，本文EMCCD真彩色成像方法獲得了較好的夜視彩色成像效果，豐富的色彩使得圖像更具空間深度感，使復(fù)雜的場景更具層次感，不同物體間的區(qū)分度增加，從而降低復(fù)雜場景中目標(biāo)識別的難度，而且真彩色效果有利于人眼長時間觀察避免視覺疲勞。

圖9 本文真彩色成像方法的其他更多結(jié)果（時間：2016-03-20 23:30；地點：北京理工大學(xué)中心花園；照度：～10－1 lx）

(Time: 2016-03-20 23:30; Location: Center Garden of BIT; Illuminance: ～10－1lx)

圖10 本文真彩色成像方法的其他更多結(jié)果（時間：2016-03-18 21:00；地點：北京理工大學(xué)6號教學(xué)樓附近；照度：～10－1 lx）

(Time: 2016-03-18 21:00; Location: No.6 Teaching Building of BIT; Illuminance: ～10－1lx)

圖11 本文真彩色成像方法的其他更多結(jié)果（時間：2016-03-21 21:30；地點：北京鳳凰嶺（郊區(qū)）；照度：～10－1 lx）

(Time: 2016-03-18 21:00; Location: Beijing Phoenix Valley (suburb) ; Illuminance: ～10－1lx)

5 結(jié)論

本文EMCCD真彩色成像方法在紅、綠、藍(lán)3個通道中增加近紅外響應(yīng)，雖然損失了一定的色彩真實性，但系統(tǒng)的光子利用率和圖像信噪比提高，方法簡便易實施，并能獲得較好的夜視彩色效果，對提高夜視系統(tǒng)性能具有一定的理論意義和實用價值。下一步，將完善真彩色EMCCD的色彩仿真工作，參考文獻(xiàn)[14]中Farrell等人的研究成果，綜合考慮照明場景的光譜能量信息、光學(xué)系統(tǒng)的參數(shù)及探測器的參數(shù)等許多重要因素，以得到更為完備的仿真結(jié)果。此外，后續(xù)還將對成像效果的色彩優(yōu)化和真彩色評價[15-16]進(jìn)行研究。

[1] Cavanillas J A A. The Role of Color and False Color in Object Recognition with Degraded and non-degraded Images[D]. Monterey: Naval Postgraduate School, 1999.

[2] The Official Magazine of the United States Marine Corps. See color at night[EB/OL]. http://marinesmagazine. dodlive.mil/ 2012/03/12/ color- nvgs/.

[3] SPI Corp. Color Low Light Night Vision Camera [EB/OL]. http:// www. x20.org/color-night-vision/.

[4] Heim G B, Burkepile J. Low-light-level EMCCD color camera[C]//, 2006, 6209: 62090F.

[5] Kriesel J M, Gat N. True-color night vision (TCNV) fusion system using a VNIR EMCCD and a LWIR microbolometer camera[C]//, 2010, 7697: 76970Z.

[6] Son J, Kang M, Min D, et al. EMCCD color correction based on spectral sensitivity analysis[J]., 2016, 75(13): 7589-7604.

[7] 李蔚, 常本康. 夜天光下景物反射光譜特性的研究[J]. 兵工學(xué)報, 2000, 21(2): 177-179.

LI Wei, CHANG Benkang. A study on the characteristics of reflection spectrum of field of view under night-sky radiation[J]., 2000, 21(2): 177-179.

[8] 金偉其. 輻射度 光度與色度及其測量[M]. 北京: 北京理工大學(xué)出版社, 2006.

JIN Weiqi.,[M]. Beijing: Beijing Institute of Technology Press, 2006.

[9] Meinel A B. Spectrum of the night sky[J]., 1958, 182(4632): 361.

[10] Jones A V, Gush H. Spectrum of the night sky in the range 1.2-2mm[J]., 1953, 172(4376): 496-496.

[11] Carapezza E M. Multi-capability color night vision HD camera for defense, surveillance, and security[C]//, 2015, 9451:94511Y.

[12] Luther R, Aus Dem Gebiet der Farbreizmetrik, Technol Z.[M]. Michigan: University of Michigan library, 1927, 8: 540-558.

[13] Quan S. Evaluation and Optimal Design of Spectral Sensitivities for Digital Color Imaging[D]. Rochester: Rochester Institute of Technology, 2004.

[14] Farrell J E, XIAO F, Catrysse P B, et al. A simulation tool for evaluating digital camera image quality[C]//, 2003, 5294: 124-131.

[15] Hawley C T. Technical assessment of low light color camera technology[C]//, 2010, 7687: 76870H.

[16] Richards A, Richardson M. EMCCD and CMOS imaging systems confront color night vision[J]., 2013, 49(2): 38-42.

Method and Demonstration of True-Color Imaging by EMCCD

BAI Yu1，CHEN Yuanjin1,2，ZHANG Mengjiao1,2，WANG Lingxue1，CAI Yi1,3

（1. Key Laboratory of Photoelectronic Imaging Technology and System, Ministry of Education of China, School of Optoelectronics, Beijing Institute of Technology, Beijing 100081, China; 2. East China Institute of Optoelectronic Integrated Device, Suzhou 215163, China; 3. China Research and Development Academy of Machinery Equipment, Beijing 100089, China）

In the light of a great deal of near-infrared radiationin the night sky, the true-color imaging method of EMCCD based on the red, green and blue optical filters with all fully transmitting NIR is studied through experiment and theory in this paper. EMCCD with electron multiplier function, can take full use of near infrared radiation in the night sky, and obtain the images of red, green, and blue channel with high brightness and signal-to-noise ratio. Ultimately, true-color EMCCD image can be obtained through comprehensive image processing methods including the color image fusion. This article makes color simulation under the conditions of the full-moon, half-moon, quarter-moon and bright-star night sky and carries out a series of studies on the production of color filters, the spectral reflectance measurement of samples, the construction of true-color EMCCD experimental device and the field tests. The experimental results show that based on the method in this paper, true-color EMCCD images with rich color, good spatial depth and high signal-noise ratio can be obtained successfully.

true-color EMCCD，night sky，true-color imaging，NIR，color filters

TN223

A

1001-8891(2017)04-0329-12

2016-08-22；

2016-09-14.

白玉（1991-），女，云南玉溪人，碩士研究生，主要從事彩色夜視成像和圖像處理的研究。E-mail：yuxiaobai927@163.com。

王嶺雪（1973-），女，云南石屏人，副教授，博士，主要從事紅外成像、圖像處理和紅外光譜等方面的研究。E-mail：neobull@bit.edu.cn。

國家自然科學(xué)基金項目（61471044）。