세줄 요약
DSO/LDSO는 Photometric Energy Term을 Residual function으로 사용한다.
이는 두 프레임 $\textbf{I}_i, \textbf{I}_j$에서의 포인트 $\textbf{p}$에 대한 밝기 차이를 의미한다.
여기에 Huber norm과 Gradient에 대한 가중치를 주어 Robust 한 모델을 만들었다.
Introduction
Direct Approach를 사용하는 DSO/LDSO는 그에 걸맞게 Photometric Energy Term을 Residual Function으로 사용합니다. Model 생성 과정을 자세히 살펴보기 위해서는 LDSO보다는 DSO 논문을 추천합니다. 사실 LDSO Residual Function은 DSO와 같으며, 해당 내용에 대한 자세한 설명은 생략되어 있습니다. 따라서, 이번 포스트에서는 DSO 논문 기반으로 Photometric Energy Term을 모델링하는 방법을 살펴보도록 하겠습니다.
Notation
$\mathbf{c}$ | Intrinsic camera parameter $$\begin{pmatrix} f_x& 0 &c_x \\ 0 & f_y&c_y \\ 0 &0 & 1\\ \end{pmatrix}$$ |
${\Pi}_{\mathbf{c}}$ | Projection($\mathbb{R}^3 \mapsto \Omega $) $$\mathbf{c}d_{\mathbf{p}}$$ |
${\Pi}_{\mathbf{c}}^{-1}$ | Back-projection($\Omega\times\mathbb{R}\mapsto\mathbb{R}^3$) $$\mathbf{c}^{-1}/d_{\mathbf{p}}$$ |
$\mathbf{p}$ | Point $\mathbf{p}\in{\Omega}_i$ in $I_i$, observed in $I_j$ |
$d_{\mathbf{p}}$ | Inverse depth of $\mathbf{p}$ |
$\mathcal{N}_p$ | $E_{\mathbf{p}j}$를 구하기 위해 필요한 $\mathbf{p}$를 포함한 8개의 점 |
$\mathcal{F}$ | All Frames |
$\omega_{\mathbf{p}}$ | Gradient-dependent weighting |
$\left\|\cdot\right\|_\gamma$ | Huber norm |
$t_i, t_j$ | Exposure times of frame $I_i, I_j$ |
$a_i, a_j, b_i, b_j$ | Affine brightness transfer function parameters of frame $I_i, I_j$ |
이외에도 많은 표기법들이 사용되었지만, 보편적으로 사용되는 표기법이라 생략하였습니다. 대표적으로 rotation, translation, transformation을 나타내는 $\mathbf{R}, \mathbf{t}, \mathbf{T}$ 등이 있습니다.
Residual Function
Paper와 소개하는 순서를 조금 바꿔봤습니다. 우선 전체적인 Residual를 먼저 소개하고, 이를 해석하기 위해 필요한 요소들을 차근차근 살펴보죠. DSO에서 사용하는 Residual은 아래와 같습니다.
$$E_{\mathbf{p}j}:=\sum_{\mathbf{p}\in \mathcal{N}_{\mathbf{p}}}\omega_{\mathbf{p}} \left\|(I_j[\mathbf{p}']-b_j)-\frac{t_je^{a_j}}{t_ie^{a_i}}(I_i[\mathbf{p}]-b_i) \right\|_\gamma$$
단순하게 보면 $I_i$의 점 $\mathbf{p}$와 이를 Back-projection과 Projection을 거쳐 추론한 $I_j$의 점 $\mathbf{p}'$의 Intensity 차이를 통해 Photometric error를 계산하고 있습니다. 사실 이렇게만 이해하고 넘어가도 전체적인 흐름을 이해하는데 전혀 지장이 없습니다. 만약 어떤 Residual을 사용하고 있는지만 궁금하신 분은 여기까지만 읽으신 후 다시 코드나 논문으로 넘어가셔도 괜찮습니다.
하지만 이렇게만 알고 넘어가기에는 수식에 너무 많은 요소들이 대놓고 보란듯이 존재합니다. 이제부터는 그러한 요소들을 하나씩 살펴보도록 하겠습니다. 소개 순서는 모델을 설계할 때 꼭 필요하지는 않다고 생각되는 요소들을 먼저 소개하고, 꼭 필요하다고 생각되는 $\mathbf{p}'$는 마지막에 소개하도록 하겠습니다. 오히려 앞서 소개하는 요소들이 DSO만의 특장점이 될 수 있으니, DSO를 깊게 이해하고 싶으신 분들은 이 부분에 집중하시는 것도 좋을 것 같네요.
Residual Pattern
DSO에서는 Photometric error의 정확도를 올리기 위해 한 점이 아닌 Residual pattern을 이용해 Residual을 계산합니다. 즉, 점 $\mathbf{p}$에서의 error는 $\mathbf{p}$와 인접하는 점 7개에서의 energy term의 합이 됩니다. 이러한 과정을 SSD, Sum of Squared Distance라고 합니다. 이 패턴을 보면 너무나 당연하게도 모두 같은 질문을 생각할 수 있습니다.
왜 점 하나가 생략된거죠? Bottom-right의 점은 어디에 있나요?
처음 패턴을 보게 되면 인접 점들을 취했다고 추측할 수 있는데요. 한 곳이 비어 있어 마음을 불편하게 합니다. 대부분 이런 경우 어른들의 속사정으로 인해 생략되는데요. DSO에서도 비슷한 이유로 생략되고 말았습니다. 논문에서는 단순히 SSE-optimized processing을 하기 위해서 생략되었다고 언급되어 있습니다. 그리고 더 이상 이야기하지 않죠.
조금 더 이해가 가능하게 설명하자면, SSE2 레지스터를 사용한다면 네 번의 float 연산을 한 번에 수행할 수 있습니다. 즉, SSE2 레지스터를 활용해 속도를 가속하기 위해서는 4의 배수로 연산량이 끊겨야 합니다. 만약 패턴이 주인공 $\mathbf{p}$를 포함해 9개의 점이라면, 우리는 추가로 연산이 필요하게 됩니다. 그래서 결과에 크게 영향을 주지 않는 선에서 bottom-right의 점을 생략한 거죠. 여기까지가 어른의 속사정이었습니다. 사실 그렇게 크게 차이가 나는지 모르겠습니다.
아무튼 각설하고 다시 본론으로 돌아오면, 이렇게 여러 점들을 계산함으로서 얻을 수 있는 이점은 무엇일까요? 저자는 다음과 같이 Residual pattern을 사용해 얻을 수 있는 장점들을 소개하고 있습니다.
- Give a good trade-off between computations required for evaluation, robustness for motion blur, and providing sufficient information
- evaluating the SSD over such a small neighborhood of pixels is similar to adding first- and second-order irradiance derivative constancy terms (in addition to irradiance constancy) for the central pixel.
요약하면, 첫 번째 포인트는 계산량은 늘어나지만 증가하는 모델의 강건성과 정보량을 생각하면 꽤나 괜찮은 딜이었다는 주장입니다. 또한, 주인공 $\mathbf{p}$에 대한 결과 역시 좋아진다고 하고 있습니다. 이에 대한 보충 설명은 아래 블로그에서 자세히 다루고 있으니 시간 되실 때 보시는 걸 추천합니다.
Huber Normalization
DSO에서는 Huber norm을 이용해 모델의 Robustness를 높이고 있습니다. 여기서 잠깐 Huber normalization을 짚고 넘어가자면, Huber Norm은 L1 Norm과 L2 Norm의 장점들을 동시에 사용하고자 제안된 많은 기법 중 하나입니다. 우선 간단하게 L1 Norm과 L2 Norm을 살펴보도록 하겠습니다.
L1 Normalization
$$\left\|w \right\|_1$$
L2 Normalization
$$\frac{1}{2} \left\| w \right\|^2$$
L1 Norm은 residual에 따라 값이 linear 하게 증가하기 때문에 outlier에 대해 민감하게 반응하지 않습니다. 다른 정상적인 샘플들이 많다면, 자연스럽게 이상치는 잊어버리죠. 반면에 L2 Norm은 residual에 따라 값이 제곱으로 증가하기 때문에 outlier에 대해 매우 민감합니다. 그래서 이상치로 인해 발생한 Overfitting을 극복하기 위해서는 훨씬 더 많은 정상적인 샘플들이 필요합니다.
이번에는 미분 가능 측면에서 두 개의 Normalization 기법을 살펴보겠습니다. 위 그래프를 통해 알 수 있듯이 L1 Norm의 경우 residual이 0인 구간에서 미분이 불가능합니다. 이는 여러모로 생각해도 치명적일 수 있죠. 반면에 L2 Norm의 경우 모든 구간에서 미분이 가능합니다.
Huber Normalization
$$\left\{\begin{matrix}
\frac{1}{2}w^2 &,|w|\leq \alpha \\
\alpha(|w|-\frac{1}{2}\alpha)& ,otherwise \\
\end{matrix}\right.$$
Huber Norm은 L2 Norm처럼 모든 구간에서 미분 가능하며, L1 Norm처럼 이상치에 강건합니다. 즉, L1 Norm에서 미분 불가능한 구간은 L2 Norm으로 대체하고, L2 Norm에서 이상치에 예민해지는 구간은 L1 Norm으로 대체하는 형태를 띠고 있습니다. 두 개의 Norm이 전환되는 경계면의 경우 대부분 실험을 통해 최적값을 찾아 사용합니다. DSO에서는 9를 식의 $\alpha$로 사용 중입니다. 하나 중요한 점은 Huber Norm의 경우 경계면이 존재하기 때문에 Non-linear optimization을 통해 최적화를 진행합니다. 이를 기억하면 도움이 될 것입니다. 아래는 L1, L2, Huber Norm을 비교한 이미지입니다.
Gradient-dependent Weighting
DSO는 Huber Norm 이외에도 Gradient-dependent weighting factor $w_{\mathbf{p}}$를 가중치로 주고 있습니다. 기울기가 큰 지점에서 가중치가 적도록 설계되어 있으며, 식은 아래와 같습니다.
$$w_{\mathbf{p}}:=\frac{c^2}{c^2+\left\| \triangledown I_i(\mathbf{p})\right\|^2_2}$$
Tracking $\mathbf{p}'$
이제 남은 변수는 $\mathbf{p}'$입니다. 이는 $I_i$ 속 우리의 주인공 $\mathbf{p}$를 Back-projection해서 world coordinate에서의 $\mathbf{p}_w$를 구한 후, $I_j$ 속으로 Projection 시킨 좌표입니다. 매우 매우 어렵게 설명했지만, 요약하면 $I_j$에서의 $\mathbf{p}$를 찾는 것입니다. 식은 아래와 같습니다.
$$\mathbf{p}'={\Pi}_c(\mathbf{R}{\Pi}_c^{-1}(\mathbf{p}, d_{\mathbf{p}})+\mathbf{t})$$
이를 해석하기 위해서는 Pinhole camera model에 대한 배경 지식이 필요한데요. 자세한 내용은 앞으로 VSLAM 포스팅들을 작성하며 다룰 예정입니다. 이번 포스팅에서는 필요한 부분만 명제처럼 짚고 넘어가도록 하겠습니다. Pinhole model에서는 pixel coordinate 상의 점 $\mathbf{P}_{uv}$와 world coordinate 상의 점 $\mathbf{P}_w$는 아래의 식에 따라 서로 변환이 가능합니다.
$$Z\mathbf{P}_{uv}=Z\begin{bmatrix}
u \\
v\\1
\end{bmatrix}=\mathbf{K}(\mathbf{RP}_w+\mathbf{t})=\mathbf{KTP}_w$$
DSO의 식과 매우 유사한데요. 여기서 $\mathbf{K}$는 DSO에서의 $\mathbf{c}$와 매칭되고, $Z$를 반대편으로 넘겨 inverse depth로 표현해 주면 거의 다 왔습니다.
$$\mathbf{P}_{uv}=d_{\mathbf{P}}\mathbf{K}(\mathbf{RP}_w+\mathbf{t})$$
이제 $\mathbf{P}_w$를 구해야 합니다. 하지만 우리가 가진 것은 관측값으로 획득한 이미지 뿐이지요. 따라서 프레임에 관측된 $\mathbf{P}_{xy}$를 Back-projection 과정을 통해 $\mathbf{P}_w$를 구하도록 합시다. 원래 $(u, v)$를 사용해 표기하지만, 다른 이미지라는 구분을 위해 $(x, y)$로 표기하였습니다.
$$\mathbf{P}_{uv}=d_{\mathbf{P}}\mathbf{K}(\mathbf{RK}^{-1}\mathbf{P}_{xy}/d_{\mathbf{P}}+\mathbf{t})$$
이제, DSO 표기법에 맞게 식을 변형해 주면 ($\mathbf{P}_{xy}\to\mathbf{p}$, $\mathbf{P}_{uv}\to\mathbf{p}'$, $\mathbf{K}d_{\mathbf{p}}\to\Pi_{\mathbf{c}}$, $\mathbf{K}^{-1}/d_{\mathbf{p}}\to\Pi^{-1}_{\mathbf{c}}$) 최종적으로 식은 아래와 같습니다.
$$\mathbf{p}'=\frac{{\Pi}_{\mathbf{c}}}{d_{\mathbf{p}}}(\mathbf{R{\Pi}_{\mathbf{c}}^{-1}}\mathbf{p}+d_{\mathbf{p}}\mathbf{t})$$
$$=\Pi_{\mathbf{c}}(\mathbf{R}{\Pi}_{\mathbf{c}}^{-1}\mathbf{p}+\mathbf{t})$$
Conclusion
최종적으로 DSO의 Photometric Energy Term은 아래와 같은 sequence를 갖게 됩니다.
$$E_{photo}:=\sum _{i\in F}\sum_{\mathbf{p}\in P_i}\sum_{j\in \mathbf{obs(p)}}E_{\mathbf{p}j}$$
이를 우리에게 친숙한(?) Pseudo code로 표현하면, 아래와 같습니다.
// Reference: https://devbull.xyz/review-dso-original/
for i in F: # i over all frames F
for p in P[i]: # p over all points P[i]
for j in obs[p]: # j over all frame obs[p]
E_pj
마지막으로 논문에 소개된 이 모든 과정에 대한 그래프를 보겠습니다. 사실 저는 직관적으로 잘 이해가 되지는 않았지만, 자세히 보면 이해에 도움이 될 것입니다. 이번 포스팅은 DSO의 Residual function을 이해하고 식을 세우기 위해 꼭 필요한 요소들에 대한 소개를 했습니다. 더 많은 부분들을 소개하고 싶었지만 Affine Brightness Transfer Parameter라던가.. 시간과 능력이 부족해 다루지 못하였는데요. 기회가 된다면 좀 더 내용을 보강할 예정이니 가끔 생각나시면 방문해 주세요. 다음 포스팅은 오늘 살펴본 Photometric Energy Term이 코드의 어떤 부분에 어떤 방식으로 구현이 되었는지 살펴볼 예정입니다. 이 부분을 Ceres
를 이용해 최적화를 하고자 해 자세히 살펴보고 있는데요. 보면 볼수록 참 오묘한 코드인 것 같네요. 빠른 시일 내로 찾아뵙도록 하겠습니다.
Reference