鱼眼参数的数值计算优化方法
? ? ? ? 在下面兩篇文章中,我們大概介紹了開源圖像處理軟件 Hugin 中魚眼圖像的矯正方法,即如何將目標全景圖上的坐標映射到魚眼圖像上,從而獲取相應的像素信息。
對于以上前面兩篇文章,我們忽略了一個問題,即以上矯正所用到的參數大部分都是不準確或者未知的,這樣會導致不同的魚眼圖像經過矯正后不能很好地拼接在一起,最直接的后果是在拼接縫處出現斷層,從而使得全景圖像的畫面不連續,影響觀感。因此我們需要對一些參數進行優化。這里主要還是參考 Hugin 中的方法,一種基于數值計算迭代優化的方法。
圖1 魚眼圖像矯正流程? ? ? ? 首先,我們對魚眼矯正的流程做個總結,如圖 1 所示。那么,我們可優化的一些參數如表 1 所示。一般來說,為了覆蓋整張全景圖像即達到 360°×180°360 \degree \times 180 \degree360°×180° 的視場,需要至少兩張魚眼圖像,例如典型的用兩個視場大于 180°180 \degree180° 的魚眼鏡頭分別在偏航為 0°0 \degree0° 和 180°180 \degree180° 的進行拍攝,這樣就可以保證覆蓋整張全景圖,同時保留一定的重疊區域以方便參數的優化以及保證圖像融合的效果。當然使用更多的魚眼鏡頭可以使得畫面更加精細,畢竟魚眼鏡頭越遠離光軸的景象越被壓縮,而且重疊區面積增大也可以提供更多的控制點對,提升參數優化和圖像融合的質量,但也會帶來成本增加,機械設備誤差或者故障發生的風險增大的問題。另外,增加一個相機意味著增加一組參數,那么參數優化的復雜度也會相應提升。
表1 可優化參數列表| 鏡頭參數 | 鏡頭有效視場 | FOV | Field of View,圓形 |
| 鏡頭徑向畸變參數 | a, b, c | 徑向畸變矯正參數 | |
| 傳感器傾斜參數 | tilt, slant, rotate | 較少使用 | |
| 傳感器中心偏移 | offset_x, offset_y | – | |
| 剪切效應參數 | shear_x, shear_y | 較少使用 | |
| 方位參數 | 鏡頭姿態 | yaw, pitch, roll | 偏航,俯仰,旋轉 |
| 光心偏移 | X, Y, Z | 詳見文章 2 | |
| 平面偏航、俯仰 | plane_y, plane_p | 詳見文章 2 |
? ? ? ? 在系統參數優化問題中,參數的數量越多,所需要的約束條件就越多,同時優化問題收斂的條件就越苛刻。為了減少需要優化的參數的數量,通常假設這些相機的固有鏡頭參數是相同的,即包括視場、徑向畸變參數、傳感器傾斜和中心偏移等這些外部一般不能改變的參數,因此只需要優化一組鏡頭參數。另外,通常也會選定一個相機作為全景圖像的基準,稱為主相機,這樣就可以減少一組需要優化的方位參數,其他的相機則可以主相機為基準進行優化。注意,這里是參數優化,而不是純粹的參數從無到有的求解,所以每個相機都應該首先提供一個合理的參數初始值,否則會大大地提升參數優化的難度,導致無法收斂。
? ? ? ? 為了進行參數優化,我們首先得確定一個優化的目標。我們知道,物理空間的一個點以不同的方位參數拍攝會落在不同的魚眼圖像的不同位置,把這些對應著物理空間同一個點的不同位置稱為控制點,一般來說以兩個作為一對。因為目標全景圖像只有一張,那么該物理空間點可以通過球面透視表示為唯一的球面坐標。因此,利用這些魚眼圖像各自的參數將控制點對根據圖 1 從下至上的過程逆推回全景圖像對應的球面時,理論上這些坐標應該是重合的。所以,我們的目標就是在不同的魚眼圖像上找出一系列的控制點對,然后優化出最佳的參數,使得這些控制點對逆映射到球面上時,彼此之間的距離最小。
? ? ? ? 假設圖 1 從下至上的逆映射過程可定義為一個非線性函數,即
xp=g(xf;θf),xf∈R2,xp∈R3,∥xp∥=1.(1){{\mathbf{x}}_p} = {\mathbf{g}}\left( {{{\mathbf{x}}_f};{{\mathbf{\theta }}_f}} \right),{\text{ }}{{\mathbf{x}}_f} \in {\mathbb{R}^2},{\text{ }}{{\mathbf{x}}_p} \in {\mathbb{R}^3},{\text{ }}\left\| {{{\mathbf{x}}_p}} \right\| = 1.\tag{1}xp?=g(xf?;θf?),?xf?∈R2,?xp?∈R3,?∥xp?∥=1.(1)
其中 xf{\mathbf{x}}_fxf? 代表魚眼圖像上的坐標,xp{\mathbf{x}}_pxp? 代表全景圖像對應的球面坐標,θf{\mathbf{\theta}}_fθf? 代表表 1 中的參數,不同魚眼圖像對應著不同的參數。我們假設這個函數是連續可導的。雖然在開頭兩篇文章的分析中,我們主要討論了圖 1 從上至下的推導過程,但我們也能很容易推導出其逆向映射過程,除了鏡頭徑向畸變校正這一步。因為我們采用的畸變校正模型中,實際落點距離 rrr 是關于入射角 θ\thetaθ 的多項式函數,那么在逆映射的過程中,我們需要從實際落點距離 rrr 逆推回入射角 θ\thetaθ,即使用校正模型的反函數。然而,該反函數我們很難推導,因此可以采用數值近似的方法,一種常用的方法是牛頓迭代法,這里不細說。
? ? ? ? 如果我們有一對控制點 (xf1,xf2)({\mathbf{x}}_{f1}, {\mathbf{x}}_{f2})(xf1?,xf2?),可以得到它們的球面坐標對 (xp1,xp2)({\mathbf{x}}_{p1}, {\mathbf{x}}_{p2})(xp1?,xp2?),那么定義這一對控制點的代價函數為
l(x;θ)=ψ(xp1,xp2),x=[xf1T,xf2T]T,θ=[θf1T,θf2T]T.(2)l\left( {{\bf{x}};{\bf{\theta }}} \right) = \psi \left( {{{\bf{x}}_{p1}},{\rm{ }}{{\bf{x}}_{p2}}} \right),{\rm{ }} {\bf{x}} = {\left[ {{\bf{x}}_{f1}^T,{\rm{ }}{\bf{x}}_{f2}^T} \right]^T},{\bf{\theta }} = {\left[ {{\bf{\theta }}_{f1}^T, {\rm{ }}{\bf{\theta }}_{f2}^T} \right]^T}.\tag{2}l(x;θ)=ψ(xp1?,xp2?),x=[xf1T?,xf2T?]T,θ=[θf1T?,θf2T?]T.(2)
圖2 球面距離示意圖其中 ψ{\psi}ψ 是兩點的距離函數,可根據實際情況選擇。通常來說,因為這兩個點是在球面上的,因此一般使用球面距離,而不使用歐幾里德距離。如圖 2 所示,假如球面上有兩點 A(φ1,θ1)A({\varphi}_1, {\theta}_1)A(φ1?,θ1?), B(φ2,θ2)B({\varphi}_2, {\theta}_2)B(φ2?,θ2?),其中 φ\varphiφ, θ\thetaθ 分別是緯度和經度。那么我們可定義其球面距離為 A, B 之間的大圓劣弧長度,即
ψ(a,b)=R?∠AOB=R?arccos?(aTb∥a∥∥b∥)=R?arcsin?(∥a×b∥∥a∥∥b∥).(3)\psi \left( {{\bf{a}},{\rm{ }}{\bf{b}}} \right) = R \cdot \angle AOB = R \cdot \arccos \left( {\frac{{{{\bf{a}}^T}{\bf{b}}}}{{\left\| {\bf{a}} \right\|\left\| {\bf{b}} \right\|}}} \right) = R \cdot \arcsin \left( {\frac{{\left\| {{\bf{a}} \times {\bf{b}}} \right\|}}{{\left\| {\bf{a}} \right\|\left\| {\bf{b}} \right\|}}} \right).\tag{3}ψ(a,b)=R?∠AOB=R?arccos(∥a∥∥b∥aTb?)=R?arcsin(∥a∥∥b∥∥a×b∥?).(3)
那么,根據式(3)我們就得到了一條約束方程。另外,我們也可以分別計算 A, B 兩點的南北和東西方向距離,其中東西方向距離可定義為平均緯度上的小圓劣弧長度,南北方向距離則可定義為緯度之間沿著經線的弧線長度,即
ψφ(a,b)=R?(∣φ1?φ2∣),(4.1){\psi _\varphi }\left( {{\bf{a}},{\rm{ }}{\bf{b}}} \right) = R \cdot \left( {\left| {{\varphi _1} - {\varphi _2}} \right|} \right), \tag{4.1}ψφ?(a,b)=R?(∣φ1??φ2?∣),(4.1)
ψθ(a,b)=R?cos?(φ1+φ2)?(∣θ1?θ2∣).(4.2){\psi _\theta }\left( {{\bf{a}},{\rm{ }}{\bf{b}}} \right) = R \cdot \cos \left( {{\varphi _1} + {\varphi _2}} \right) \cdot \left( {\left| {{\theta _1} - {\theta _2}} \right|} \right).\tag{4.2}ψθ?(a,b)=R?cos(φ1?+φ2?)?(∣θ1??θ2?∣).(4.2)
這樣,我們實際上就得到了兩個約束方程。一般來說,在優化過程中,我們不會同時使用式(3)和式(4)。
? ? ? ? 由此,我們可以定義參數優化所用到的代價函數為
F(Θ)=12∑i=1m(fi(Θ))2=12∥f(Θ)∥2,,(5.1)F\left( {\mathbf{\Theta }} \right) = \frac{1}{2}\sum\limits_{i = 1}^m {{{\left( {{f_i}\left( {\mathbf{\Theta }} \right)} \right)}^2}} = \frac{1}{2}{\left\| {{\mathbf{f}}\left( {\mathbf{\Theta }} \right)} \right\|^2},{\tag{5.1}}, F(Θ)=21?i=1∑m?(fi?(Θ))2=21?∥f(Θ)∥2,,(5.1)
fi(Θ)=li(xi;θi),Θ=θ1∪θ2∪?∪θm∈Rn.(5.2){f_i}\left( {\mathbf{\Theta }} \right) = {l_i}\left( {{{\mathbf{x}}_i};{{\mathbf{\theta }}_i}} \right),{\text{ }}{\mathbf{\Theta }} = {{\mathbf{\theta }}_1} \cup {{\mathbf{\theta }}_2} \cup \cdots \cup {{\mathbf{\theta }}_m} \in {\mathbb{R}^n}. \tag{5.2}fi?(Θ)=li?(xi?;θi?),?Θ=θ1?∪θ2?∪?∪θm?∈Rn.(5.2)
其中 mmm 代表總共有 mmm 個約束方程,注意其不一定等于控制點的對數,因為在式(4)中一個控制點對可以定義兩個約束方程。nnn 代表總共有 nnn 個需要優化的參數,主要每個控制點對用到的魚眼圖像可能不同,所以其參數也就不一定相同,但也有一些公共的參數,例如鏡頭參數,因此這里使用了它們的并集。盡管在式(3)和(4)中我們使用的是球面距離,但在求所有控制點對的整體距離時還是使用了平方和,這相當于把球面距離當作一種樣本誤差,我們的目標是要最小化整體的誤差平方和。對于 F(Θ)F\left( {\mathbf{\Theta }} \right)F(Θ) 來說我們其實并不知道其真正的函數形式,但為了方便一般假設其是一個高度非線性的連續函數,這時可將問題形式化為
Θ?=arg?min?Θ{F(Θ)}.(6){{\bf{\Theta }}^ * } = \arg {\min _{\bf{\Theta }}}\left\{ {F\left( {\bf{\Theta }} \right)} \right\}.\tag{6}Θ?=argΘmin?{F(Θ)}.(6)
? ? ? ? 由于 F(Θ)F\left( {\mathbf{\Theta }} \right)F(Θ) 為平方和形式,因此也將式(6)稱為最小二乘問題,相關的優化算法也有很多,比較常用的是 LM 算法,其原理可參考我的另一篇博客。注意,在 LM 算法優化過程中需要用到 F(Θ)F\left( {\mathbf{\Theta }} \right)F(Θ) 的 Jacobian 矩陣,也就是各 fi(Θ){f_i}\left( {\mathbf{\Theta }} \right)fi?(Θ) 關于各個參數的一階偏導組成的矩陣,由于我們并不知道各 fi(Θ){f_i}\left( {\mathbf{\Theta }} \right)fi?(Θ) 的具體形式,因此可通過一些數值近似方法,例如比較簡單的有限差分法,即給某一個參數增加非常小的步長 δ\deltaδ,其他的參數保持不變,求解出新的函數值,那么函數值的差值除以 δ\deltaδ 即為關于該參數的近似導數,即
f(x+δ)=f(x)+f′(x)δ+o(δ),(7.1)f(x + \delta ) = f(x) + f'(x)\delta + o(\delta) , \tag{7.1} f(x+δ)=f(x)+f′(x)δ+o(δ),(7.1)
f′(x)≈f(x+δ)?f(x)δ,δis??sufficient??samll.(7.2)f'(x) \approx \frac{{f(x + \delta ) - f(x)}}{\delta },{\rm{ }}\delta \text{{\rm{ is {\rm{ }}sufficient {\rm{ }}samll}}{\rm{.}}} \tag{7.2}f′(x)≈δf(x+δ)?f(x)?,δ?is??sufficient??samll.(7.2)
不過這種方法計算量比較大,而且 δ\deltaδ 應該足夠小,否則不一定能滿足有限差分的收斂性條件。
? ? ? ? 在 LM 算法中,我們需要保證 F(Θ)F\left( {\mathbf{\Theta }} \right)F(Θ) 的 Jacobian 矩陣是列滿秩的,一般來說也就要求約束方程要比需要優化的參數多,所以我們應該提供足夠多的質量好的控制點對。PanoTools 全景工具庫中使用了兩種優化策略,也就是式(3)和(4)的區別。式(3)中一個控制點對只能提供一個約束方程,由于約束性相對不那么強,因此收斂速度會慢一些,但是穩定性相對較好。而式(4)中一個控制點對提供了兩個約束方程,約束性要強一些,因此在初始參數設置較好時會有較快的收斂速度,但如果初始參數較差時穩定性會差一些。因此 PanoTools 中也采用了混合的優化策略,即前期使用式(3),只要每次迭代 F(Θ)F\left( {\mathbf{\Theta }} \right)F(Θ) 的下降比例都達到一定要求,否則切換到式(4)繼續優化,這樣就可以提高收斂速度和穩定性。
總結
以上是生活随笔為你收集整理的鱼眼参数的数值计算优化方法的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: C# 色系表配色 颜色表 美工必备
- 下一篇: indexeddb_深入IndexedD