IFFT 기반의 바다 시뮬레이션

IFFT을 통해 구현한 바다

컴퓨터 그래픽에서 바다를 구현하는 건 마냥 쉬운 일은 아닙니다. 시각적으로 그럴싸한 바다를 만드는 것 자체도 어렵지만, 사용자의 성능 또한 고려해야 한다면 더더욱 어려워지죠.

사인파 (Sine Wave)

가장 단순한 형태로 "바다"를 그린다면 어떤 모양일까요? 아마 사인파 (Sine Wave)에서 출발하게 될 겁니다.

제가 이전에 구현했었던 Wave 아트워크 역시 이 사인파를 기반으로 만든 메쉬들을 이래저래 겹쳐 움직이도록 한 것이 전부입니다. 물론 이건 추상적인 형태로 바다를 그려낸 것이라 실제 바다에 가깝게 묘사했다고 보기엔 어려움이 있죠.

사인파로 바다를 구현한다고 했을 때에 눈에 띄는 문제점은, 사인파로는 단순히 상하 운동만 처리할 수 있다는 것입니다. 반면 실제 해수면의 물 입자는 단순히 위 아래로 흔들리는 것이 아니라 원형 궤도를 그리며 움직인다는 차이가 있죠.

거스트너 웨이브 (Gerstner Wave)

기존 사인파가 정점의 높이($y$)만 변화시키는 것과 달리, 거스트너 웨이브는 정점을 파도의 진행방향 $x$와 $z$로도 함께 밀어냅니다. 그 결과 파도의 마루(Peak)는 뾰족해지고 골(Trough)은 넓고 평평해지는 바다 특유의 파형이 만들어집니다.

$i$번째 파도가 정점을 얼마나 이동시키는지를 나타내는 변위 벡터 $\Delta \mathbf{P}_i(\mathbf{p}_0, t)$는 아래와 같이 간략하게 나타낼 수 있습니다.

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

\Delta \mathbf{P}_i(\mathbf{p}_0, t) = \begin{pmatrix} d_{i,x} \, a_i \cos(f_i) \\ a_i \sin(f_i) \\ d_{i,z} \, a_i \cos(f_i) \end{pmatrix}

제 경우는 이전에 Godot으로 아래와 같은 레트로 스타일의 Ocean 아트워크를 작업했던 적이 있는데, 이 때 구현한 파도는 3개의 거스트너 웨이브를 합산한 형태입니다.

영상에서 보이는 바다는 합산하는 파도의 개수 자체가 그렇게 많지 않기 때문에, 다소 인위적이고 단순한 형태입니다. 거스트너 웨이브의 개수가 늘어나면 늘어날수록 점점 더 디테일한 형태의 바다를 구현할 수 있겠지만, 그만큼 계산 과정이 무거워집니다.

거스트너 웨이브 기반의 바다는 컴퓨트 셰이더나 렌더타깃 텍스처를 사용하지 않고 단일 셰이더 하나만으로도 즉시 동작하는 바다를 구현할 수 있기 때문에, 그 간편함에 강점이 있습니다.

IFFT(Inversed Fast Fourier Transform) 기반의 바다

IFFT 기반의 바다는 이번 포스트를 통해 다루고자 했던 핵심적인 주제입니다. IFFT 기반 바다의 핵심 원리는 해양학에서 실제 바다의 파랑 에너지를 계측하여 정립된 통계적 파랑 스펙트럼 (Statistical Wave Spectrum)을 주파수 도메인에 생성한 후, 이를 역 푸리에 변환(IFFT)하여 공간 도메인 상의 높이장(Height field)와 변위장(Displacement field)을 한번에 얻어내는 것입니다.

파수 벡터 $\mathbf{k}$

버퍼의 해상도가

N \times N

이고, 실제 바다 타일의 물리적 크기가 $L \times L$ 일 때, 그리드 인덱스

n, m \in [-N/2, N/2)

에 대응하는 파수 벡터(Wave Vector) $\mathbf{k}$ 는 다음과 같습니다.

\mathbf{k} = (k_x, k_z) = \left( \frac{2\pi n}{L}, \frac{2\pi m}{L} \right), \quad k = \Vert{}\mathbf{k}\Vert{} = \sqrt{k_x^2 + k_z^2}

\mathbf{k} = (k_x, k_z) = \left( \frac{2\pi n}{L}, \frac{2\pi m}{L} \right), \quad k = \Vert{}\mathbf{k}\Vert{} = \sqrt{k_x^2 + k_z^2}

\mathbf{k} = (k_x, k_z) = \left( \frac{2\pi n}{L}, \frac{2\pi m}{L} \right), \quad k = \Vert{}\mathbf{k}\Vert{} = \sqrt{k_x^2 + k_z^2}

중앙 쪽으로 갈수록 파수(Wavenumber) $k$가 낮으며, 파수가 낮다는 것은 파장(Wavelength)가 길다는 것을 의미합니다. 다시 말해 가운데로 갈수록 파장이 긴 평평한 물결(저주파)가 되며, 테두리 쪽으로 갈수록 파장이 짧은 잔물결(고주파)이 됩니다.

초기 스펙트럼 생성 - $\tilde{h}_0(\mathbf{k})$

각 주파수 픽셀이 가져야 할 초기 에너지와 위상을 결정하는 과정입니다. 이는 시뮬레이션 시작 시, 또는 바람이 바뀔 때마다 최초 1회 실행됩니다.

Phillips 스펙트럼

해양학의 Phillips 스펙트럼 $P_h(\mathbf{k})$ 는 풍속 $V$ , 중력 $g$ , 풍향 단위벡터 $\hat{\mathbf{w}}$ 가 주어졌을 때 각 파수 벡터 $\mathbf{k}$ 가 가지는 에너지 크기를 정의합니다.

P_h(\mathbf{k}) = A \frac{\exp\left(-1 / (k L_{max})^2\right)}{k^4} \vert{}\hat{\mathbf{k}} \cdot \hat{\mathbf{w}}\vert{}^2, \quad \left(L_{max} = \frac{V^2}{g}\right)

단, 여기서 $\hat{\mathbf{k}} \cdot \hat{\mathbf{w}} \le 0$ 인 경우 $P_h(\mathbf{k}) = 0$으로 둡니다. 이를 통해 역풍에 동일한 에너지가 할당되는 것을 막습니다.

$\hat{\mathbf{w}} = (1, 1), V = 120$ 일 때의 $P_h(\mathbf{k})$

하지만 이 스펙트럼을 그대로 사용할 경우, 지나치게 규칙적인 형태의 바다를 얻게 됩니다. 이를 해결하기 위해 평균이

0

이고 표준편차가

1

인 가우시안 난수 쌍 $(X_r, X_i)$ 를 곱합니다.

Box-Muller Transform

가우시안 난수쌍을 이용해 복소수를 생성하는 가장 대표적인 방법은 박스-뮬러 변환(Box-Muller Transform)입니다. 구간 $(0, 1)$에서 균일하게 추출된 두 개의 독립 난수를 $U_1$과 $U_2$라고 할 때, 가우시안 난수쌍 $(X_r, X_i)$는 다음과 같이 계산됩니다.

$$ X_r=\sqrt{-2\ln U_{1}}\cos (2\pi U_{2}) \\ X_i=\sqrt{-2\ln U_{1}}\sin (2\pi U_{2}) $$

이를 앞서 구한 스펙트럼에 곱해 복소수의 초기 진폭 $\tilde{h}_0(\mathbf{k})$를 만듭니다.

\tilde{h}_0(\mathbf{k}) = \frac{1}{\sqrt{2}} (\xi_r + i \xi_i) \sqrt{P_h(\mathbf{k})}

\tilde{h}_0(\mathbf{k}) = \frac{1}{\sqrt{2}} (\xi_r + i \xi_i) \sqrt{P_h(\mathbf{k})}

\tilde{h}_0(\mathbf{k}) = \frac{1}{\sqrt{2}} (\xi_r + i \xi_i) \sqrt{P_h(\mathbf{k})}

\tilde{h}_0(\mathbf{k}) = \frac{1}{\sqrt{2}} (\xi_r + i \xi_i) \sqrt{P_h(\mathbf{k})}

시간 진행에 따른 스펙트럼 갱신

스펙트럼 계산 이후 실시간으로 움직이는 파도를 얻기 위해선 해당 스펙트럼을 시간 $t$에 따라 갱신해주어야 합니다.

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

\tilde{h}(\mathbf{k}, t) = \tilde{h}_0(\mathbf{k}) e^{i\omega(k)t} + \tilde{h}_0^*(-\mathbf{k}) e^{-i\omega(k)t}

IFFT

이제 파도가 주파수 도메인에서 잘 표현되었지만, 이것을 실제 화면으로 렌더링하기 위해서는 공간 도메인으로 돌아가야 합니다. 이 계산에서 필요한 것이 바로 IFFT 입니다.

해당 작업은 스펙트럼의 해상도에 따라 연산량이 크게 달라지기 때문에, 컴퓨트 셰이더를 통해 GPU에서 연산하는 것이 안전합니다.

FFT 알고리즘을 이 포스트에서 다루기엔 다소 부담스럽기 때문에, 여기서는 우리가 앞서 계산했던 실시간 스펙트럼에 IFFT를 적용하여 아래와 같이 Displacement Map을 얻을 수 있다는 점만 짚고 가도록 하겠습니다.

IFFT를 거쳐 얻은 Displacement Map

결과

이로부터 얻은 Displacement Map을 바다 메쉬의 버텍스 상에 적용하면, 실제 바다처럼 출렁거리는 표면을 확인할 수 있습니다.

저는 이를 제 아트워크인 "피난처"의 구현에 활용했습니다. 사실 제 작업에는 노멀값의 활용이 필요없었기 때문에 Displacement Map으로도 충분했으나, 좀 더 현실적인 형태의 렌더링을 요구하는 상황이라면 노멀에 대한 계산도 동일하게 IFFT를 거쳐 이루어질 수 있습니다.

이 블로그 검색

Shubidumdu Blog