陳青青1,2,，湯志亞1，楊玲1,2
?。?.成都信息工程大學,，四川 成都 610225；2.中國氣象局 大氣探測重點開放實驗室,，四川 成都 610225）
摘要：采用地基鏡頭旋轉掃描拍攝不同俯仰角的天空圖像,，通過全景圖拼接獲得全天空圖像，避免了魚眼透視畸變,，保留了低仰角信息,。目前全景圖拼接技術主要局限于單層柱面，通過空間坐標轉換及以漸入漸出加權平均融合算法和插值算法,，借鑒地圖學投影原理,，將28張圖像在垂直投影面和等積柱面投影面進行無縫拼接得到360°×180°的全天空全景圖，為識別云狀,、等積計算云量提供了連續(xù)完整的全天空圖像,。
0引言
圖1原始圖像隨著大氣探測技術的發(fā)展，云的觀測逐漸自動化,，國內外研發(fā)了大量不同波段不同視角的測云儀器［1］,，全景拍攝方式簡單快捷但低仰角畸變大,；另外還有全景圖拼接方式［2］。本文采用監(jiān)控攝像機,，無盲區(qū)掃描拍攝可見光波段天空圖像,，得到36張具有重疊部分的子圖，通過還原物體在天球中的相對位置,，并借鑒地圖投影轉換為平面圖像,，對重疊部分進行融合，可得到高清連續(xù)的全天空云圖,。
1觀測平臺
云臺型號為HY900B,，相機為索尼FCBEX1010P，采集器為HB18204T,，圖像傳感器CCD為EXview HAD 1/4英寸,，成像面積為3.2 mm×2.4 mm，固定焦距為3.4 mm,，視角為57.8°,，圖像尺寸為704×576，長邊視角為44.70°,，寬邊視角為36.58°,，圖像上下左右都存在重疊，無盲區(qū),。實驗用機配置Intel Core i3 CPU,，2 GB內存，Windows XP系統(tǒng),。原始圖像如圖1所示,，拍攝俯仰角順序如表1所示。
令拍攝次序t=1~36,，其中1~9表示相機光軸仰角為0°的圖像,，10~18表示仰角30°圖像……由此可見第t張圖像方位角θt仰角φt與t存在如式(1)所示的關系，其中［·］表示向下取整,，%表示求余,。
相機標定以標準網格角點作為匹配點進行多項式擬合［3］，將擬合關系式應用到整幅圖像并對坐標進行雙線性插值,，截去無效邊緣,，如圖2所示。
圖2相機標定
2圖像在球面上的配準
圖像拼接包括配準和融合兩個環(huán)節(jié),。本文圖像有確定的空間坐標,，能實現圖像在空間位置上的鑲嵌，從而找到配準點,。
以相機位置為原點,，拍第一幅圖像時相機光軸方向為0°方位角,，觀測點所在的大地水平面為xoy面，0°方位角為x軸的正向,，90°方位角為y軸正向，過原點xoy面指向天空的垂線為z軸,。0°方位角圖像在空間xoz平面上的剖面圖如圖3(a)所示,。其中E′表示相機光軸投影點E點在xoz面上的投影，在xoz平面上的坐標為(xE,zE),由圖3(a)可見,。
其中,，L表示物距，P為圖像中任意一點在xoz面上的投影,，如圖3(a)所示,，P在x軸上的坐標為xP=xE+Δx,P在z軸上的坐標為zP=zE－Δz，而
由此得出兩個關系式：
zp=Lsinφt－(j－a)cosφt(5)
xp=Lcosφt+(j－a)sinφt(6)
聯立式（5）,、式（6）解得：
zp=Lcscφt－xPcotφt
且
y∈［－b+1,－b+W］
x∈［Lcosφt+(1－a)sinφt,Lcosφt+(H－a)sinφt］(7)
0°方位角圖像在空間中xoy平面的剖面圖如圖3(b)所示,，可見
yp=－b+i(8)
數字圖像坐標是離散值，直接由（j,，i,，t）計算（x，y,，z）,，會導致有些空間坐標找不到與之對應的圖像坐標，投影后圖像出現空洞,。采用反算插值的方法可避免出現這種情況,，即根據投影面上每個點來反查其在原圖中的位置和像素值［4］。
當相機光軸θt=0°,，φt≠0°時,，聯立式（6）、式（8）反解j,、i得：
聯立式（7）,、式（9）即得到θt=0°，φt≠0°時球面各點（x,，y,，z）所對應的第t張圖的圖像坐標（j，i,，t）,。
同理，當θt=0°,、φt=0°時,，聯立公式（5）,、式（8），反解出
聯立式（6）,、式（10）即得到θt=0°,、φt=0°時球面各點對應的第t張圖的圖像坐標。
探測范圍方位角0°~360°,，仰角0°~180°,，半徑0~∞,該范圍的物體可投影到以相機為球心，給定值為半徑的半球面上,。假想景物已從空間投影到了天球面上,，相機拍攝將天球物體投影到CCD，因此應采用點光源在球心的球面投影方式［5］,。
設P(xP,yP,zP)是方位角0°,，仰角φt的圖像上任意一點，利用直角坐標與球坐標的轉換公式得到P點的球坐標為(θp,φp,rp),。
P點與投影點P′位于同一條直線上,，二者的θ和φ相同，而所有投影點的r都等于投影球面半徑R,，因此P點的投影坐標為P′(θp,φp,R),。
云臺順時針旋轉，那么球面投影坐標（θp－θt,，φp,，R）對應第t張圖像的圖像坐標（j，i）,。如式（1）所示,，改變θt就能通過方位角0°的圖像投影點求得其余8個方位對應的圖像投影坐標點，從而得到360°全景圖,。仰角90°時只取一張圖就能覆蓋拍攝范圍,，因此僅需用到前28張圖。為避免信息冗余和損失,，本文取物距與長邊弦長的中值為投影球面半徑,。
3垂直投影面的拼接
垂直投影光源為一束與z軸平行的平行光，投影面為xoy面,，投影公式如式（11）所示,，（x，y,，z）表示空間某點,，（x′，y′，z′）為其投影點,。
圖4全天空云圖垂直投影實驗結果MATLAB仿真垂直投影結果如圖4所示,，銜接無誤位置準確，但邊界亮度躍變,，這是相機自動曝光引起的,。數字圖像以矩陣形式存儲，反算下標會取整而產生誤差,，應盡量減少取整次數,，在最后一級投影才對下標取整和圖像融合。
加權平均融合是指在匹配點重合的n個子圖各自對應的3個顏色通道的亮度值進行加權平均,，新值賦給匹配點,，當權值相等時稱為均值融合,，均值融合結果如圖5所示,。
結果表明均值融合后圖像個體完整，但邊界明顯,，原因是沒有考慮權值差異,，距離匹配點近的子圖貢獻權值應更大。以距離比為權值,，將距離近的圖像在該點的值作為主元素進行融合,，稱為漸入漸出加權平均融合（后面簡稱加權融合）［6］。
子圖垂直投影之后的圖像是扇形的,，由于圖像是圍繞圓心旋轉的,，權值便不應用直線距離，而是該點到邊界的弧長,。對30°層進行加權融合實驗結果如圖6所示,，可見該方法解決了子圖邊緣突變，效果良好,。
對60°層進行一次加權融合,，僅顯示其重疊部分如圖7(a)所示，可見重疊部分還有重疊,，再次加權融合,，結果如圖7(b)所示，重疊消失,。最后將一級融合后的像元和無重疊像元放入,，形成60°層全景如圖7(c)所示。
先進行同仰角層內融合,，再進行層間融合,。層間融合用半徑的變化量來計算權值，圖8（a）是30°,、60°層間未融合圖像,，圖8（b）是30°,、60°層間加權融合，圖8（c）是融合后的全天空圖像,。
4等積柱面投影面的拼接
云量是按視面積來計算的,，等積投影［5］為云量的計算提供了借鑒。前人證明軸射柱面投影下,，任意形狀的球面與其在外切圓柱面上投影面積相等［7］,。設A是球面上任意一點，A的柱面坐標為(θA,ρA,zA),，A′是A的軸射柱面投影,，其柱面坐標為(θA′,ρA′,zA′)。投影公式為：
利用柱面坐標與直角坐標的轉換公式,，以直角坐標為媒介,，將球坐標轉換為柱面坐標。以柱面底圓的弧長為橫軸,，z軸為縱軸,，沿著0°方位角處的母線將柱面展開成投影平面。重合部位采用加權融合,，先同仰角層內融合,，再層間融合，60°仰角仍需二次融合,。
全天空圖加權融合實驗結果如圖9(a)所示,，完整呈現了低仰角圖像信息，云量無失真地得以保留,，但頂部（仰角90°附近）圖像拉伸產生了空洞,。采用漸入漸出加權平均插值算法利用空洞周邊像素對空洞插值。實驗結果如圖9(b)所示,，空洞消除且圖像呈連續(xù)自然變化,。
5結論
本文通過空間坐標轉換和地圖投影方法，用垂直投影和等積柱面投影實現了全天空云圖拼接,，采用漸入漸出加權平均融合算法解決了邊界融合,。該方法成本低廉，能夠保留子圖的信息量,，解決了低仰角圖像信息嚴重損失的弊端,，為識別云狀、計算云量提供了參考,。提高電機角速度可縮短子圖時差,，根據經緯度和時間計算位置，再插值可去除太陽，結合其他波段或雙目相機進行聯合觀測可進一步優(yōu)化,。
致謝：感謝中國氣象局大氣探測中心馬舒慶研究員,、陶法博士等對本課題的支持。
參考文獻
?。?］ 高太長,，劉磊，趙世軍,等.全天空測云技術現狀及進展［J］.應用氣象學報,，2010,，21(1):101108.
［2］ 孫曉鋼,，孫學金,，牛珍聰，等．全天空云圖獲取的一種方式及算法實現［J］．氣象科學,，2008,，28(3):338341.
［3］ 唐爍,繆源.基于Harris角點的圖像匹配算法［J］.微型機與應用,2013,32(2):4143.
?。?］ 李曉輝,周蔭清,王祖林,，等.基于曲面拼接的球面全景生成算法［J］.北京航空航天大學學報,2007,33(6):668681.
?。?］ 焦健,曾琪明.地圖學［M］.北京：北京大學出版社,，2005.
［6］ 鄭永安,宋建社,，曹繼平,，等.基于小波變換的區(qū)域特征加權自適應圖像融合算法［J］.計算機工程與應用，2004,，28(6):9495,，127.
［7］ 楊道富,彭明.球面和圓柱面相等對應關系的研究［J］.商丘師范學院學報,2005,21(5):6163.