基于直方圖分布的航天相機(jī)動(dòng)態(tài)范圍調(diào)整方法

李春梅 郭元榮 董書莉 胡雨婷 張斐然

?

基于直方圖分布的航天相機(jī)動(dòng)態(tài)范圍調(diào)整方法

李春梅 郭元榮 董書莉 胡雨婷 張斐然

（北京空間機(jī)電研究所，北京 100094）

動(dòng)態(tài)范圍自動(dòng)調(diào)整對航天相機(jī)有重要意義，尤其是相機(jī)動(dòng)態(tài)范圍調(diào)整頻繁，實(shí)時(shí)性要求較高時(shí)，需要?jiǎng)討B(tài)范圍實(shí)時(shí)自動(dòng)調(diào)整才能滿足要求。以灰度均值、最值作為動(dòng)態(tài)范圍調(diào)整的依據(jù)，都缺乏對圖像信息分布的整體理解；而以圖像熵作為圖像動(dòng)態(tài)范圍調(diào)整的依據(jù)，效率不足。為此，文章提出了基于直方圖分布特征的相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整方法。該方法在統(tǒng)計(jì)圖像直方圖的基礎(chǔ)上，進(jìn)一步實(shí)施直方圖信息提取，設(shè)計(jì)了直方圖區(qū)域面積函數(shù)作為動(dòng)態(tài)范圍評價(jià)依據(jù)，并給出動(dòng)態(tài)范圍調(diào)整策略，調(diào)整的目的是使圖像直方圖更接近于理想的分布。基于直方圖分布特征的方法相較于灰度均值、最值和圖像熵的方法，能得到更優(yōu)的調(diào)整效果。經(jīng)大量實(shí)驗(yàn)驗(yàn)證，基于直方圖分布特征的相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整方法，能夠獲得較好的相機(jī)動(dòng)態(tài)范圍調(diào)整效果，且效率較高，運(yùn)行穩(wěn)定。

直方圖 動(dòng)態(tài)范圍 自動(dòng)調(diào)整 評價(jià)依據(jù) 航天相機(jī)

0 引言

動(dòng)態(tài)范圍是航天相機(jī)的一項(xiàng)重要性能指標(biāo)，它是指相機(jī)能夠產(chǎn)生線性輸出電平的入瞳輻亮度范圍，以臨界飽和的入瞳輻亮度與形成噪聲（均方根）等效曝光量的入瞳輻亮度之比表示[1]。動(dòng)態(tài)范圍設(shè)置是否合理關(guān)系到衛(wèi)星圖像的層次、亮度和對比度，最終影響像質(zhì)。如果設(shè)置不合理，會(huì)在圖像中出現(xiàn)高端飽和、丟失信息或高端閑置、信息壓縮在低端的情況。動(dòng)態(tài)范圍自動(dòng)調(diào)整對航天相機(jī)有重要意義。某型號相機(jī)需要24小時(shí)連續(xù)開機(jī)成像，成像的時(shí)間跨度、軌道跨度都很大，甚至成像瞬間跨越晨昏線，相機(jī)動(dòng)態(tài)范圍調(diào)整頻繁，實(shí)時(shí)性要求較高，需要實(shí)時(shí)自動(dòng)調(diào)整動(dòng)態(tài)范圍才能滿足要求。

目前常用的相機(jī)動(dòng)態(tài)范圍調(diào)整方法有亮度均值法[2-11]（還有一般均值法、權(quán)重均值法、基于直方圖的加權(quán)亮度均值法）、亮度最值法、圖像熵法[12]。均值法無論是基于整幅圖像的均值，還是局部感興趣圖像區(qū)域的均值，亦或者感興趣灰度區(qū)域的均值，用于航天相機(jī)動(dòng)態(tài)范圍調(diào)整有明顯不足。首先，基于均值的調(diào)整，可能導(dǎo)致圖像中出現(xiàn)過多飽和像元；其次，遙感圖像各區(qū)域都需要適度曝光，沒有感興趣區(qū)域或不感興趣區(qū)域的劃分。亮度最值法，對圖像噪聲、探測器中的奇異點(diǎn)，或者圖像中的個(gè)別景物反光點(diǎn)異常敏感，容易出現(xiàn)誤判。而圖像熵法對于調(diào)整相機(jī)動(dòng)態(tài)范圍效率較低。因此，以上三種動(dòng)態(tài)范圍調(diào)整的方法，均不適用于航天相機(jī)成像應(yīng)用。

針對航天相機(jī)的成像特點(diǎn)，本文提出了基于直方圖分布特征的動(dòng)態(tài)范圍自動(dòng)調(diào)整方法。該方法以直方圖中幾個(gè)特定區(qū)域的面積為評價(jià)圖像亮暗的依據(jù)，調(diào)整的目的是使圖像直方圖更接近于理想的分布。經(jīng)大量實(shí)驗(yàn)驗(yàn)證，本文所提出的方法能夠獲得較好的相機(jī)動(dòng)態(tài)范圍調(diào)整效果，且效率較高，運(yùn)行穩(wěn)定。

1 動(dòng)態(tài)范圍調(diào)整方法

1.1 直方圖統(tǒng)計(jì)

灰度直方圖是對圖像亮度的一種統(tǒng)計(jì)表達(dá)，它反映了該圖中不同灰度級出現(xiàn)的統(tǒng)計(jì)概率。一幅數(shù)字圖像在范圍[0,]內(nèi)總共有個(gè)灰度級，其直方圖定義為離散函數(shù)[13]，見式（1）：

通常，我們會(huì)用到歸一化的直方圖，即直方圖中所有元素除以圖像中的像素總數(shù)所得到的數(shù)據(jù)，見式（2）：

式中(r)表示灰度級r出現(xiàn)的頻率；為圖像中的像素總數(shù)。

直方圖統(tǒng)計(jì)算法由于占用資源較大，在實(shí)時(shí)在軌圖像處理中應(yīng)用較少。改進(jìn)的直方圖統(tǒng)計(jì)方法，將圖像數(shù)據(jù)位寬轉(zhuǎn)化為5bit（截取高5位），設(shè)置的值為31（人眼能夠分辨的只有32個(gè)級別），同樣實(shí)現(xiàn)圖像灰度級統(tǒng)計(jì)的目的，又節(jié)省算法的資源消耗。

位寬為8bit的圖像數(shù)據(jù)，計(jì)算256級直方圖統(tǒng)計(jì)結(jié)果如圖1（a）與截取高5bit、計(jì)算32級直方圖統(tǒng)計(jì)結(jié)果如圖1（b），兩者對比如圖1所示。

1.2 直方圖信息提取

灰度直方圖是評價(jià)圖像質(zhì)量的重要客觀依據(jù)，對圖像的灰度分布給出了總體描述。從直方圖中可以得到圖像的灰度范圍，每個(gè)灰度級的出現(xiàn)概率、灰度級的分布、整體圖像的平均明暗和對比度等信息，為圖像評價(jià)提供重要依據(jù)。在暗色的圖像中，直方圖的組成成分集中在灰度級低的一側(cè)；反之，明亮的圖像的直方圖的組成成分集中在灰度級高的一側(cè)。圖像直方圖比較集中，則表明該圖像的動(dòng)態(tài)范圍小，如圖2所示；如果圖像的直方圖分布比較均勻地分布在較大的灰度范圍內(nèi)，則表明圖像具有較大的動(dòng)態(tài)范圍，如圖3所示。

（a）256級歸一化直方圖 （b）32級歸一化直方圖

（a）較小動(dòng)態(tài)范圍圖像 （b）較小動(dòng)態(tài)范圍圖像的直方圖

（a）較大動(dòng)態(tài)范圍圖像 （b）較大動(dòng)態(tài)范圍圖像的直方圖

為了從直方圖中獲得圖像灰度分布更直觀的信息，這里定義了直方圖區(qū)域面積函數(shù)，見式（3）。

式中(r)為直方圖中從灰度r到灰度31的面積；(31)代表了圖像中的飽和（包括已飽和和近飽和）信息的總量；(15)代表了圖像中中高輻亮度信息的總量；1-(6)代表了圖像中的低灰度信息的總量。從變量(r)中，能夠快速獲得直方圖分布的特征信息，可以作為相機(jī)動(dòng)態(tài)范圍調(diào)整的依據(jù)。

1.3 直方圖信息的使用

依據(jù)一行圖像（或一幀圖像）的直方圖信息調(diào)整相機(jī)動(dòng)態(tài)范圍容易引起算法的震蕩，不穩(wěn)定。本文提出將連續(xù)采集的多行（或多幀）圖像的直方圖信息作為參考，對多行圖像的直方圖信息分別設(shè)置不同的權(quán)值，統(tǒng)計(jì)后的結(jié)果再作為相機(jī)動(dòng)態(tài)范圍調(diào)整的依據(jù)，見式（4）：

2 動(dòng)態(tài)范圍調(diào)整

2.1 動(dòng)態(tài)范圍調(diào)整的策略

根據(jù)線陣掃描航天相機(jī)的成像特點(diǎn)，設(shè)定如下動(dòng)態(tài)范圍調(diào)整策略：1）由于成像地物的連續(xù)性和成像時(shí)間的連續(xù)性，相機(jī)成像圖像的直方圖信息不會(huì)發(fā)生突變，因此相機(jī)的動(dòng)態(tài)范圍調(diào)整（除開機(jī)時(shí)刻）不應(yīng)過于頻繁。2）具備“就高”和“就低”兩種工作模式，“就高”模式的工作原則是，優(yōu)先保證圖像中高輻亮度景物的成像性能，在盡量避免圖像出現(xiàn)飽和或允許少量飽和區(qū)域的情況下，調(diào)整相機(jī)動(dòng)態(tài)范圍，擴(kuò)展圖像灰度直方圖，增加景物的灰度細(xì)節(jié)；“就低”模式的工作原則是，優(yōu)先保證圖像中低輻亮度景物的成像性能，使其具有合適的灰度分層，不能保證或不完全保證高輻亮度景物成像，圖像可能出現(xiàn)一定面積的飽和。用戶根據(jù)目標(biāo)特性和成像需求選擇工作模式。3）動(dòng)態(tài)范圍調(diào)整應(yīng)快速高效。

“就高”模式調(diào)整策略如圖4所示，圖中W1、W2、W3、W4、W5、W6分別為就高模式策略參數(shù)。

“就低”模式調(diào)整策略如圖5所示，圖中V1、V2、V3、V4、V5、V6、V7分別為就低模式策略參數(shù)。

圖5 “就低”模式流程圖

2.2 動(dòng)態(tài)范圍調(diào)整的實(shí)現(xiàn)

航天相機(jī)獲取圖像，其圖像的動(dòng)態(tài)范圍是隨成像參數(shù)而變化[14-15]。某航天相機(jī)所使用的探測器利用高性能的低噪聲輸出放大技術(shù)和電子倍增增益調(diào)節(jié)技術(shù)，實(shí)現(xiàn)了高動(dòng)態(tài)范圍的成像。其中電子倍增可以實(shí)現(xiàn)高達(dá)500倍的放大[16-18]，為相機(jī)提供了寬動(dòng)態(tài)范圍微光成像的能力。微光相機(jī)結(jié)合微光探測器的積分時(shí)間和電子倍增調(diào)整[19]，設(shè)置了10個(gè)檔位（相鄰兩檔為2倍），實(shí)現(xiàn)了高達(dá)1 000倍的動(dòng)態(tài)范圍調(diào)整能力。由于成像地物的連續(xù)性和成像時(shí)間的連續(xù)性，相機(jī)成像圖像的直方圖信息不會(huì)發(fā)生突變，因此相機(jī)的動(dòng)態(tài)范圍調(diào)整步長固定為1檔，每次檔位增加或降低，只變化一個(gè)檔位。

基于直方圖分布特征的航天相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整方法是航天CCD相機(jī)、CMOS相機(jī)等可見光和多光譜相機(jī)能夠使用的通用方法，動(dòng)態(tài)范圍調(diào)整的手段依據(jù)相機(jī)所使用的探測器的特征而定。

3 實(shí)驗(yàn)結(jié)果

某航天相機(jī)利用該算法，通過調(diào)整積分時(shí)間和電子倍增，采取“就高”和“就低”模式，實(shí)現(xiàn)了相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整的目的。算法在當(dāng)前幀結(jié)束前輸出相機(jī)調(diào)整參數(shù)，滿足了實(shí)時(shí)性要求。

為驗(yàn)證所提出算法的有效性和可靠性，實(shí)驗(yàn)選取了多類遙感圖像完成動(dòng)態(tài)范圍調(diào)整的仿真。其中遙感云圖的動(dòng)態(tài)范圍調(diào)整仿真效果如圖6所示。

在圖6中，標(biāo)準(zhǔn)云圖（a）經(jīng)過人為動(dòng)態(tài)范圍壓縮為云圖（c）；對云圖（c）使用“就高”模式動(dòng)態(tài)范圍調(diào)整后的云圖為（e），在保證圖像不飽和或少量飽和的前提下實(shí)現(xiàn)了直方圖拉伸；對云圖6（c）使用“就低”模式動(dòng)態(tài)范圍調(diào)整后的云圖為（g），保證圖像中低灰度直方圖拉伸的前提下做到了盡可能少量像元飽和。結(jié)果顯示，經(jīng)算法調(diào)整后的圖像具有較好的主觀評價(jià)結(jié)果，其圖像直方圖也表現(xiàn)出較寬的動(dòng)態(tài)范圍，達(dá)到了動(dòng)態(tài)范圍調(diào)整的目的。

（a）標(biāo)準(zhǔn)云圖 （b）標(biāo)準(zhǔn)云圖的直方圖

（a）Standard cloud picture （b）Histogram of standard cloud picture

（c）動(dòng)態(tài)范圍設(shè)置偏小時(shí)的云圖 （d）動(dòng)態(tài)范圍偏小時(shí)的直方圖

（c）Cloud picture with narrow dynamic range （d）Histogram of cloud picture with narrow dynamic range

（e）“就高”模式動(dòng)態(tài)范圍調(diào)整后的云圖 （f）“就高”模式動(dòng)態(tài)范圍調(diào)整后的直方圖

（e）Cloud picture adjusted using“High”mode （f）Histogram of cloud picture adjusted using“High”mode

（g）“就低”模式動(dòng)態(tài)范圍調(diào)整后的云圖 （h）“就低”模式動(dòng)態(tài)范圍調(diào)整后的直方圖

4 結(jié)束語

本文根據(jù)航天相機(jī)的成像特點(diǎn)提出一種基于直方圖分布特征的相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整方法，完成了從圖像灰度數(shù)據(jù)中提取更可靠的數(shù)據(jù)信息作為相機(jī)動(dòng)態(tài)范圍調(diào)節(jié)的依據(jù)，并設(shè)計(jì)動(dòng)態(tài)范圍調(diào)整的策略，實(shí)現(xiàn)自動(dòng)調(diào)整的目的。實(shí)驗(yàn)結(jié)果表明，該算法能夠較好的完成相機(jī)動(dòng)態(tài)范圍自動(dòng)調(diào)整的目的，且執(zhí)行效率高?！熬透摺焙汀熬偷汀眱煞N工作模式，既滿足大動(dòng)態(tài)范圍寬覆蓋成像的需求，也滿足小動(dòng)態(tài)范圍低反射率景物灰度細(xì)分高質(zhì)量成像的要求。

[1] 陳世平,楊秉新,王懷義, 等.空間相機(jī)設(shè)計(jì)與試驗(yàn)[M]. 宇航出版社, 2003. CHEN Shiping, YANG Bingxin, WANG Huaiyi, et al. Space Camera Design and Experiment[M]. Aerospace press, 2003. (in Chinese)

[2] 隋龍, 趙健, 畢國玲. 變參數(shù)CCD相機(jī)自動(dòng)曝光調(diào)節(jié)算法[J]. 儀器儀表學(xué)報(bào), 2011, 32(6): 199-202. SUI Long, ZHAO Jian, BI Guoling. Varing Parameters Exposure Adjustment Algorithm for CCD Cameras[J]. Chinese Journal of Scientific Instrument, 2011, 32(6): 199-202. (in Chinese)

[3] 安然, 陳彥, 解靜. CMOS圖像傳感器自適應(yīng)動(dòng)態(tài)范圍調(diào)整曝光控制[J]. 紅外與激光工程, 2013, 42(S1): 88-92. AN Ran, CHEN Yan, XIE Jing. Exposure Algorithm for CMOS Image Sensor with Adaptive Dynamic Range[J]. Infrared and Laser Engineering, 2013, 42(S1): 88-92. (in Chinese)

[4] 仵宗欽, 劉鵬. 基于FPGA的數(shù)字相機(jī)自動(dòng)曝光技術(shù)研究[J]. 電視技術(shù), 2015, 39(11): 117-121. WU Zongqin, LIU Peng. Study on Technology of Digital Camera Automatic Exposure Based on FPGA[J]. Video Engineering, 2015, 39(11): 117-121. (in Chinese)

[5] 陳紅圓. 基于圖像處理的自動(dòng)曝光算法專利技術(shù)綜述[J]. 電視技術(shù), 2015, 39(S1): 5-7. CHEN Hongyuan. Patent Overview of Auto-exposure Algorithm Based on Image Processing Technique[J]. Video Engineering, 2015, 39(S1): 5-7. (in Chinese)

[6] 楊海濤, 常義林, 王靜, 等. 一種基于亮度直方圖的自動(dòng)曝光控制方法[J]. 光學(xué)學(xué)報(bào), 2007, 27(5): 841-847. YANG Haitao, CHANG Yilin, WANG Jing, et al. A New Automatic Exposure Algorithm for Video Cameras Using Luminance Histogram[J]. Acta Optica Sinica, 2007, 27(5): 841-847. (in Chinese)

[7] 黃輝, 任國強(qiáng), 孫健, 等. 基于圖像直方圖統(tǒng)計(jì)的CCD相機(jī)自動(dòng)調(diào)光算法[J]. 半導(dǎo)體光電, 2013, 34(4), 698-701. HUANG Hui, REN Guoqiang, SUN Jian, et al. Auto-exposure Algorithm for CCD Camera Using Image Histogram[J]. Semi-conductor Optoelectronics, 2013, 34(4), 698-701. (in Chinese)

[8] 梁佳毅. 高性能數(shù)碼相機(jī)自動(dòng)曝光算法研究與實(shí)現(xiàn)[D]. 上海: 復(fù)旦大學(xué), 2008. LIANG Jiayi. Study and Implementation of Exposure Algorithms in Advanced Digital Camera[D]. Shanghai: Fudan University, 2008. (in Chinese)

[9] 楊作廷, 阮萍, 翟波. 基于圖像熵的高動(dòng)態(tài)范圍場景的自動(dòng)曝光算法[J]. 光子學(xué)報(bào), 2013, 42(6): 742-746. YANG Zuoting, RUAN Ping, ZHAI Bo. Auto-exposure Algorithm for Scenes with High Dynamic Range Based on Image Entropy[J]. Acta Photonica Sinica, 2013, 42(6): 742-746. (in Chinese)

[10] AHMED S, ANTONIO R. Exposure Time Calculation for Spectral Cameras[J]. Journal of Electronic Imaging, 2015, 24(5): 053025.

[11] 匡艷, 張振環(huán), 鄭裕林. 一種基于圖像多區(qū)域的變步長自動(dòng)曝光算法[J]. 電視技術(shù), 2015, 39(19), 28-31. KUANG Yan, ZHANG Zhenhuan, ZHENG Yulin. Auto Exposure Algorithm Using Image Subregion and Variable Step Size[J]. Video Engineering, 2015, 39(19), 28-31. (in Chinese)

[12] 孟夢. 高清攝像機(jī)自動(dòng)白平衡和自動(dòng)曝光算法研究及實(shí)現(xiàn)[D]. 西安: 西安電子科技大學(xué), 2014. MENG Meng. Study and Implementation of Automatic White Balance and Exposure Algorithms in HD Camera[D]. Xian: Xidian University, 2014. (in Chinese)

[13] RAFAEL C G, RICHARD E W, STEVEN L E. 數(shù)字圖像處理[M]. 電子工業(yè)出版社, 2005. RAFAEL C G, RICHARD E W, STEVEN L E. Digital Image Processing[M]. Publishing House of Electronics Industry, 2005. (in Chinese)

[14] 何紅艷, 王小勇, 付興科. 遙感衛(wèi)星CCD相機(jī)的動(dòng)態(tài)范圍設(shè)計(jì)思考[J]. 航天返回與遙感, 2008, 29(1): 39-49. HE Hongyan, WANG Xiaoyong, FU Xingke. Study on the Designing of the Dynamic Range of Remote Satellite’s CCD Camera[J]. Spacecraft Rocovery & Remote Sensing, 2008, 29(1): 39-49. (in Chinese)

[15] 李智勇, 楊校軍. 關(guān)于遙感衛(wèi)星TDICCD相機(jī)動(dòng)態(tài)范圍設(shè)計(jì)的思考[J]. 航天返回與遙感, 2011, 32(1): 24-27. LI Zhiyong, YANG Xiaojun. Consideration on Designing of Dynamic Range for TDICCD Camera on Remote Sensing Satellites[J]. Spacecraft Recovery & Remote Sensing, 2011, 32(1): 24-27. (in Chinese)

[16] 何家維, 何昕, 魏仲慧, 等. EMCCD相機(jī)在全天時(shí)星敏感器中的設(shè)計(jì)與應(yīng)用[J]. 紅外技術(shù), 2013, 35(11), 718-722. HE Jiawei, HE Xin, WEI Zhonghui, et al. Design and Application of EMCCD Camera for Full-time Star Sensor[J]. Infrared Technology, 2013,35(11), 718-722. (in Chinese)

[17] 盧家莉, 李彬華, 胡泊. 一種電子倍增CCD電子倍增增益的確定方法[J]. 兵工學(xué)報(bào), 2015, 36(4), 710-715. LU Jiali, LI Binhua, HU Po. A Method of Determining EMCCD Electron Multiplication Gain[J]. Acta Armamentarii, 2015, 36(4), 710-715. (in Chinese)

[18] 胡泊, 李彬華. 低溫下EMCCD電子倍增模型[J]. 電子學(xué)報(bào), 2013, 41(9), 1826-1830.HU Po, LI Binhua. Electron Multiplication Model of EMCCD in Low Temperature[J]. Acta Electronica Sinica, 2013, 41(9), 1826-1830. (in Chinese)

[19] 宋秋冬, 林秋月, 劉會(huì)通, 等. 基于FPGA的CCD相機(jī)自動(dòng)增益控制方法及實(shí)現(xiàn)[J]. 紅外與激光工程, 2006, 35(增刊)，366-370 .SONG Dongqiu, LIN Qiuyue, LIU Huitong, et al. Automatic Gain Control Method and Achievement of CCD Digital Camera Based on FPGA[J]. Infrared and Laser Engineering, 2006, 35(Supplement): 366-370. (in Chinese)

（編輯：劉穎）

Dynamic Range Adjustment Method of the Aerospace Camera Based on Histogram Distribution

LI Chunmei GUO Yuanrong DONG Shuli HU Yuting ZHANG Feiran

（Beijing Institute of Space Mechanics & Electricity, Beijing 100094, China）

It is very important for the aerospace camera to adjust the dynamic range automatically. Especially when the camera dynamic range requires high-frequency and real-time adjustments, the automatic adjustment of dynamic range will be indispensable. Methods using image gray scale average or maximum value as the basis of dynamic range adjustment, lack understanding the overall image information distribution, and methods using image entropy as the basis are inefficient in computing. Therefore, a new automatic adjustment method of camera dynamic range was proposed based on histogram distribution characteristics in this paper. On the basis of gathering image histogram, the information of histogram was further extracted in this method. A histogram regional area function was designed as the evaluation criterion of dynamic range, and then the dynamic range adjustment strategy was given according to this function. The purpose of adjustment was to make the image histogram closer to ideal distribution. Better performance can be achieved by the proposed method in this paper than those three ones mentioned previously. Lots of experiments showed that the automatic adjustment method of camera dynamic range based on histogram distribution characteristics can obtain better dynamic range adjustment results with high efficiency and stability.

histogram; dynamic range; automatic adjustment; evaluation criterion; Aerospace camera

TN401

A

1009-8518(2017)05-0036-08

10.3969/j.issn.1009-8518.2017.05.005

李春梅，女，1980年生，2005年獲哈爾濱工業(yè)大學(xué)儀器科學(xué)與技術(shù)專業(yè)碩士學(xué)位，高級工程師。研究方向?yàn)檫b感圖像處理與應(yīng)用。E-mail：9757290@qq.com。

2017-03-08

國家重大科技專項(xiàng)工程