Korean Journal of Remote Sensing (대한원격탐사학회지)

Korean Society of Remote Sensing (대한원격탐사학회)
권호 표지
DOI QR코드

DOI QR Code

Improvements for Atmospheric Motion Vectors Algorithm Using First Guess by Optical Flow Method

옵티컬 플로우 방법으로 계산된 초기 바람 추정치에 따른 대기운동벡터 알고리즘 개선 연구

  • Oh, Yurim (PhD Candidate, Department of Atmospheric Sciences, Division of Earth Environmental System, Pusan National University) ;
  • Park, Hyungmin (PhD Candidate, Department of Atmospheric Sciences, Division of Earth Environmental System, Pusan National University) ;
  • Kim, Jae Hwan (Professor, Department of Atmospheric Sciences, Division of Earth Environmental System, Pusan National University) ;
  • Kim, Somyoung (PhD Candidate, Department of Atmospheric Sciences, Division of Earth Environmental System, Pusan National University)
  • 오유림 (부산대학교 지구환경시스템학부 대기과학과 박사수료생) ;
  • 박형민 (부산대학교 지구환경시스템학부 대기과학과 박사수료생) ;
  • 김재환 (부산대학교 지구환경시스템학부 대기과학과 정교수) ;
  • 김소명 (부산대학교 지구환경시스템학부 대기과학과 박사수료생)
  • Received : 2020.08.19
  • Accepted : 2020.09.25
  • Published : 2020.10.31
Download PDF
⟨ Previous Next ⟩

Abstract

Wind data forecasted from the numerical weather prediction (NWP) model is generally used as the first-guess of the target tracking process to obtain the atmospheric motion vectors(AMVs) because it increases tracking accuracy and reduce computational time. However, there is a contradiction that the NWP model used as the first-guess is used again as the reference in the AMVs verification process. To overcome this problem, model-independent first guesses are required. In this study, we propose the AMVs derivation from Lucas and Kanade optical flow method and then using it as the first guess. To retrieve AMVs, Himawari-8/AHI geostationary satellite level-1B data were used at 00, 06, 12, and 18 UTC from August 19 to September 5, 2015. To evaluate the impact of applying the optical flow method on the AMV derivation, cross-validation has been conducted in three ways as follows. (1) Without the first-guess, (2) NWP (KMA/UM) forecasted wind as the first-guess, and (3) Optical flow method based wind as the first-guess. As the results of verification using ECMWF ERA-Interim reanalysis data, the highest precision (RMSVD: 5.296-5.804 ms-1) was obtained using optical flow based winds as the first-guess. In addition, the computation speed for AMVs derivation was the slowest without the first-guess test, but the other two had similar performance. Thus, applying the optical flow method in the target tracking process of AMVs algorithm, this study showed that the optical flow method is very effective as a first guess for model-independent AMVs derivation.

수치예보모델의 예측 바람장은 대기운동벡터 알고리즘의 표적 추적 과정에서 추적 정확도 향상이나 계산 시간 단축을 위해 초기 추정치로 사용된다. 대기운동벡터는 수치예보모델의 자료동화 시 활용가치가 높다고 알려졌으나, 초기 추정치로 사용된 수치예보모델 바람장이 대기운동벡터의 검증 과정에 참 값으로 사용된다는 모순이 있다. 이를 해결하기 위해서는 수치예보모델로부터 독립적인 초기 추정치가 필요하다. 본 연구에서는 Lucas and Kanade 옵티컬 플로우 방법을 적용하여 바람장을 도출한 후 이를 초기 추정치로 사용함으로써 표적 추적과정에서의 모델 의존성을 제거하고 계산 속도를 향상시키고자 하였다. 대기운동벡터 산출에는 2015년 8월 18일 ~ 9월 5일 00, 06, 12, 18시 동안의 정지궤도 위성 Himawari-8/AHI의 14번 채널 Level 1B 자료를 사용하였다. 옵티컬 플로우 방법이 대기운동벡터 산출에 미치는 영향을 평가하기 위하여 다음과 같은 세가지 방법으로 교차 검증을 수행 하였다. (1) 초기 추정치 없이, (2) KMA/UM 예보바람장을 초기 추정치로 사용하여, 그리고 (3) 옵티컬 플로우 방법으로 계산된 바람장을 초기 추정치로 사용하여 대기운동벡터를 산출하고 ECMWF ERA-Interim 재분석장과 비교 검증한 결과, 옵티컬 플로우 기반 바람장을 초기 추정치로 사용한 경우에 가장 높은 정밀도를 보였다(RMSVD: 5.296-5.804 ms-1). 계산 속도는 초기 추정치를 사용하지 않은 경우에 가장 느렸고, 나머지 테스트는 유사한 속도를 보였다. 그러므로 대기운동벡터 알고리즘의 표적 추적 과정에 옵티컬 플로우 방법을 적용하면, 모델 의존성 없는 고품질 바람벡터의 산출이 가능할 것으로 사료된다.

Keywords

1. 서론

대기운동벡터(Atmospheric Motion Vectors)는 정지궤도 또는 저궤도 위성으로 관측된 연속적인 영상으로부터 구름이나 수증기의 움직임을 추적하여 산출되는 바람정보이다(Forsythe et al., 2007). 대기운동벡터는 지상 관측이 제한적인 저위도 및 해양에서 시·공간적으로 균질한 바람정보를 제공하는 장점을 가지며 실시간 기상 감시 및 수치모델의 자료동화에 이용되어 기상 예측 정확도 향상에 기여한다(Velden et al., 2005; Wu et al., 2014; Lee et al., 2015; Kim et al., 2017; Le Marshall et al., 2016). 또한 정지궤도 위성 기반의 대기운동벡터는 짧은 주기의 관측이 가능하여 태풍 및 대류운의 추적 및 감시를 위한 자료로도 활용되고 있다(Mecikalski and Bedka, 2006; Langland et al., 2009; Kim et al., 2012; Apke et al., 2016; Oyama et al., 2018).

일반적으로 대기운동벡터 산출에는 특정 시간을 기준으로 이전과 이후 세 장의 연속된 위성영상이 이용되며, 산출과정은 표적선정, 고도할당, 표적추적, 품질관 리의 네 가지 단계로 다음과 같이 수행된다 (Holmlund, 1998; Shimoji, 2014; Oh et al., 2015; Oh et al., 2019). 기준 시간(t=t0)의 위성 관측 영상에서 추적하기 적당한 영역을 표적으로 선정한 후(표적선정 과정), 선정된 표적의 고도를 계산하고(고도할당 과정), 선정된 표적과 가장 유사한 영역을 이전 시간(t=t0-Δt) 영상과 이후 시간(t=t0+Δt) 영상을 이용하여 풍향 및 풍속을 계산한다(표적추적 과정). 마지막으로 산출된 대기운동벡터의 객관적인 품질 지표를 계산한다(품질관리 과정). 이 네 가지 과정 중 하나인 표적추적 과정에서, 계산 효율 증대를 위해 초기 추정치를 이용한다. 초기 추정치 사용으로 표적의 위치가 어느 정도 예측될 수 있기 때문에 추적 영역의 면적을 줄일 수 있고, 이로 인해 계산 시간 단축이 가능하다(Bedka and Mecikalski, 2005; Bresky et al., 2012). 또한 추적 영역이 축소됨에도 불구하고 그 위치가 유동적으로 결정되므로 기존의 고정된 추적 영역에서는 계산하기 어려운 빠른 풍속까지도 산출할 수 있다는 장점이 있다. 현재 대부분의 대기운동벡터 알고리즘에서는 초기 추정치로 모델 예측자료를 사용하고 있다. 3차원 바람장인 모델 예측자료에서 초기 추정치를 가져오기 위해서는 선정된 표적의 고도할당이 선행되어야 한다. 고도할당 과정은 대기운동벡터 산출과정 중 가장 큰 오차가 발생되는 부분이므로, 모델 예측자료에서 잘못 계산된 고도의 초기 추정치를 이용하면 문제가 발생할 수 있기 때문에 표적추적 과정의 오차를 불러올 수 있다(Borde and Doutriaux-Boucher, 2014; Salonen et al., 2015). 이러한 이유로 최근에는 표적추적 과정에서 초기 추정치를 사용하지 않고 대기운동벡터를 산출하는 시도가 늘어나고 있다(Borde and Garcia-Pereda, 2014; Garcia-Pereda and Borde, 2014). 그러나 초기 추정치를 사용하지 않으면, 추적 영역의 중심과 표적의 중심이 일치해야 하므로 추적 영역의 면적이 커야 하기 때문에 계산 시간이 증가한다. 또한 추적 영역의 면적이 커질수록 추적 영역 내에 표적과 유사한 패턴이 많아지게 되어 실제 표적의 이동을 추적하는데 어려움이 발생한다(Shimoji, 2014).

본 연구에서는 대기운동벡터 산출의 효율성과 정확성을 증대시킬 수 있는 초기 추정치를 사용하되, 예보 모델과 독립적인 바람장 산출을 위해 지금까지 시도된 적이 없었던 컴퓨터 비전에서 사용되고 있는 Lucas and Kanade 기반의 옵티컬 플로우(Optical Flow)로 계산된 바람장을 초기 추정치로 사용하고 이를 기존 방법과 비교하였다. 옵티컬 플로우는 컴퓨터 비전에서 널리 사용되는 기법으로, 연속된 이미지 영상에서 밝기 패턴이 동일한 곳을 찾아 표적의 이동을 찾는 방법이다(Horn and Schunck, 1981; Wu et al., 2016). 옵티컬 플로우는 위성에서 관측된 영상을 기반으로 계산되기 때문에 표적의 고도를 결정할 필요가 없고, 기존의 초기 추정치를 획득하는 과정에서 발생하는 벡터의 고도 오차를 무시할 수 있기 때문에 큰 장점이 있다. 그러나 옵티컬 플로우는 자동차나 사람 등 강체 (rigid body)의 움직임을 추적하는데 특화되어 있기 때문에 구름과 같이 발달과 소멸을 하면서 모양이 변하는 것을 추적하는 기존의 대기운동 벡터 측정 방법에 적용 시 한계가 있다. 따라서 본 연구에서는 타겟의 이동방향을 선정하는데 초기 추정치로의 사용 가능성을 처음으로 시도하였다.

옵티컬 플로우 기반 초기 추정치를 적용한 알고리즘을 현재 대기운동벡터 산출 알고리즘과 비교하여 평가하기 위해, 1) 초기 추정치 없이, 2) 기상예보장을 초기 추정치로 사용하여, 그리고 3) 옵티컬 플로우로 얻어진 바람장을 초기 추정치로 사용하여 대기운동벡터를 계산하고 그 결과를 ECMWF의 ERA-Interim 재분석장과 비교하였다. 대기운동벡터 산출에 사용된 위성 영상은 일본의 정지궤도 기상위성인 Himawari-8/AHI (Advanced Himawari Imager)의 주야간 사용 가능한 11.2 μm 적외채널(IR11.2) 영상이며, 표적추적은 교차상관 방법(Cross-Correlation method; CC)을, 고도할당에는 Equivalent Black-Body Temperature (EBBT) 방법을 사용하였다. 대기운동벡터 산출과 검증을 위해 사용된 자료는 2장에서 설명하였으며, 3장에는 옵티컬 플로우 알고리즘의 최적화를 위해 수행한 알고리즘 내 주요 인자에 대한 실험 결과와 대기운동벡터 산출 알고리즘에 대하여 기술하였다. 다양한 초기 추정치에 따라 산출된 대기운동벡터와 ERA-Interim 자료와의 검증 결과는 4장에서, 마지막으로 5장에서 요약과 결론을 기술하였다.

2. 연구자료

본 연구에서는 일본의 정지궤도 기상위성인 Himawari-8/AHI 자료를 이용하여 대기운동벡터를 산출하였다. 이 위성은 2014년 10월 7일 발사되어 동경 140.7도에서 관측을 수행하고 있으며, 10분 간격의 연속적인 전구 자료를 제공하고 있다(Bessho et al., 2016). Himawari-8 위성에 탑재된 AHI 센서는 0.5와 1 km의 공간 해상도를 갖는 4개의 가시채널과 2 km 공간 해상도를 갖는 3개의 단파적외채널, 3개의 수증기채널, 그리고 6개의 적외채널로 구성된다. 대기운동벡터는 전 채널에서 생산이 가능하지만, 본 연구에서는 시간 제약 없이 주·야간 대기운동벡터 산출이 가능하여 현업에서 활용도가 높은 적외채널(IR11.2)의 10분 간격 영상을 이용하여 대기운동벡터를 산출하였다. 산출에 이용된 자료의 기간은 2015년 8월 19일부터 9월 5일 0000, 0600, 1200, 1800 UTC이다. 표적추적을 위한 초기 추정치 결정과 고도할당을 위하여 수치모델 자료가 이용되었으며, 기상청에서 운영중인 통합모델(United Model; UM)의 6시간 간격 예보장을 이용하였다. 모델의 수평해상도는 약 25 km이며, 본 연구에서는 구름이 존재하는 위치의 모델자료를 얻기 위해 1013.25 hPa에서 85.18 hPa까지 연직 28층의 수평 바람과 온도, 상대습도, 기압 자료와 지표면 온도, 기압 자료를 사용하였다. 산출된 대기운동벡터의 검증을 위하여 European Centre for Medium-Range Weather Forecasts (ECMWF)의 ERA-Interim 자료를 이용하였다. 사용된 ERA-Interim 자료는 0.125°×0.125°의 공간해상도를 가지며, 1000 hPa에서 10 hPa까지 총 31개 연직층으로 이루어진 바람 성분 자료이다. 주로 육상의 관측만 가능한 라디오존데와는 달리, ERA-Interim 재분석 바람장 자료는 시간과 공간적으로 균질한 자료를 제공하므로, 육상과 해상 모두에서 정량적인 검증이 필요한 대기운동벡터 검증에 적합한 자료라고 할 수 있다. 검증에 사용된 자료는 대기운동벡터 산출과 동일 기간의 0000, 0600, 1200, 1800 UTC 자료를 사용하였다.

3. 연구방법

1) 이미지 피라미드를 활용한 옵티컬 플로우 바람장 산출

옵티컬 플로우 계산식은 영상의 밝기는 시간이 지나더라도 일정하다고 가정하며, 표적의 움직임이 1 화소 이내의 작은 움직임이라는 전제를 통해 편미분항의 고차항들을 생략하여 식 (1)과 같이 정의할 수 있다(Horn and Schunck, 1981).

\(\left(\frac{\partial I}{\partial x} \frac{\partial I}{\partial y}\right)\left(\frac{u}{v}\right)=-\frac{\partial I}{\partial t}\)       (1)

식 (1)의 I는 연속된 두 장의 영상 중에서 이후 시간 영상의 밝기를 의미하며, x, y는 표적의 위치를, t는 시간을, 그리고 u, v는 시간에 따른 표적의 이동을 나타낸다. 그러나 이 방정식은 방정식의 개수보다 미지수의 개수가 더 많기 때문에, 식을 풀기 위해 Lucas and Kanade(1981)는 탐색 영역이라는 개념을 도입했다. 탐색 영역은 정사각형인 영역 안에 포함된 각 화소가 모두 동일한 움직임을 가진다고 가정하는 공간으로, 추가된 제약조건을 통해 구한 최종 옵티컬 플로우 계산식은 식 (2)에 표현하였다.

\(\left(\begin{array}{cc} I_{x 1} & I_{y 1} \\ I_{x 2} & I_{y 2} \\ \vdots & \vdots \\ I_{x n} & I_{y n} \end{array}\right)\left(\begin{array}{l} u \\ v \end{array}\right)=-\left(\begin{array}{c} I_{t 1} \\ I_{Q} \\ \vdots \\ I_{t n} \end{array}\right)\)       (2)

AV=-d       (3)

식 (2)의 아래 첨자 x, y, t는 각 방향으로의 미분을, n은 탐색영역의 화소 개수를 의미하고, 식 (3)과 같이 간략하게 표현할 수 있다. 식 (3)의 행렬 A는 영상에서의 공간 기울기를 의미하며, 벡터 V는 표적의 움직임을, 벡터 d는 시간에 따른 변화를 의미한다. 탐색 영역을 통해 새롭게 정의된 선형방정식은 과결정된(over-determined) 시스템으로 최소자승법으로 최적의 해를 찾으면, 식 (4)와 같이 표현될 수 있으며, \(G=\left(\begin{array}{lll} \sum_{i=1}^{n} I_{x i} I_{x i} & \sum_{i=1}^{n} I_{x i} I_{y i} \\ \sum_{i=1}^{n} I_{x i} I_{y i} & \sum_{i=1}^{n} I_{y i} I_{y i} \end{array}\right)\), \(b=\left(\begin{array}{l} \sum_{i=1}^{n} I_{x i} I_{t i} \\ \sum_{i=1}^{n} I_{y i} I_{t i} \end{array}\right)\))이다.

V=(ATA)-1ATd=G-1b       (4)

기존의 옵티컬 플로우 알고리즘은 표적의 움직임이 1화소보다 적은 경우를 가정하지만, 본 연구에서는 1화소 이상의 강한 풍속까지도 계산하기 위해 이미지 피라미드 방법과 결합한 알고리즘을 사용하였다. 이미지 피라미드는 원본 영상의 크기를 단계적으로 축소시키면서 생산된 일련의 영상 집합을 의미한다. Fig. 1은 이미지 피라미드의 기본 개념을 나타낸 것으로 최하층 수준은 원본 영상을 나타내며, 상층 수준으로 올라갈수록 영상의 크기가 절반으로 줄어드는 것을 볼 수 있다. 그러므로 상층 수준에서의 1화소 이동은 그 아래 수준에서의 2화소 이동에 해당하며, 이미지 피라미드 수준의 개수가 많아질수록 더 큰 이동이 가능해짐을 알 수 있다. 이러한 성질을 이용하여 이미지 피라미드를 활용한 옵티컬 플로우 알고리즘은 기존 알고리즘의 한계를 극복하고 큰 움직임까지도 계산할 수 있다(Lee and Park, 2007). Fig. 2은 본 연구에서 사용한 이미지 피라미드 기반의 옵티컬 플로우 알고리즘 흐름도를 수식을 통해 나타낸 것이다. 연속된 두 장의 영상(H, I)이 입력되면 각 영상에서 이미지 피라미드를 생산하고, 최상위 수준에서부터 옵티컬 플로우 계산이 시작되어 최하위 수준인 원본 영상에 이를 때까지 계산이 진행된다. 피라미드 각 수준에서는 최적의 옵티컬 플로우를 구하기 위해 탐색 영역의 중심을 이동시키면서 반복적으로 해 Vk를 구하는 과정을 거친다. 이 전의 해와 현재 구해진 해의 차이가 특정 문턱값보다 작게 나타나는 조건을 만족하면, 해를 기반으로 하여 그 다음 피라미드 수준에서의 x, y 값을 갱신하고 옵티컬 플로우를 계산하는 과정이 반복된다. 다음 수준으로 넘어갈 때 각 점의 위치와 계산된 해는 이미지 피라미드가 하위 수준으로 갈수록 고해상도가 되므로 2배가 된다. 또한 이미지 피라미드의 수준은 식 (5)와 같이 제한된다.

OGCSBN_2020_v36n5_1_763_f0001.png 이미지

Fig. 1. Subsampling of the image pyramid using Himawari8 imagery.

OGCSBN_2020_v36n5_1_763_f0002.png 이미지

Fig. 2. Flowchart of the Optical Flow algorithm. L denotes the number of pyramid levels, and subscript n is the window size (2w + 1) × (2w + 1). The x, y of the image H means the position of the target. Vk is the optimal solution obtained by iterative computation and dL+1 is output displacement which updates the location of the target in the image I at L+1th level.

\(\begin{array}{} \text { max_pyramid_level = }\begin{array}{c} \text { Floor }\left(\log _{2}\left(\frac{\text { Min }(\text { Original domain size })}{64}\right)\right)= \operatorname{Floor}\left(\log _{2} \frac{5500}{64}\right)=6 \end{array} \end{array}\)       (5)

이는 이미지 피라미드를 구축하는 과정에서 수행되는 부 표본화 과정이 과하게 수행될 경우 이미지 가장자리 부분이 왜곡될 수 있기 때문이다. 따라서 본 연구에서는 이미지 피라미드 레벨의 개수를 여섯 개로 제한하였다.

기존의 Lucas-Kanade 방법을 이용하면, 1화소 이내의 작은 움직임만 계산할 수 있으며 이미지 피라미드를 활용한다고 하여도 최대 106 m/s의 바람밖에 계산하지 못한다는 단점이 있다. 따라서 본 연구에서는 이러한 한계를 보완하기 위해 기존의 Lucas-Kanade 방법으로 작은 움직임을 계산한 후, 움직인 위치로 기준점을 이동하여 다시 작은 움직임을 계산하고 이런 작은 움직임들을 합하는 과정을 반복하여 최종적으로 큰 움직임을 계산하는 개선된 Lucas-Kanade 방법을 활용하였다. 이러한 과정을 통해 피라미드 개수가 1개이어도 강한 풍속까지 계산될 수 있다(Fig. 2).

Fig. 2을 통해 알 수 있듯이 본 연구에서 사용한 옵티컬 플로우 알고리즘에서 정의되는 파라미터는 세 가지 (문턱값, 피라미드 개수, 탐색 영역 크기)이다. 첫 번째로 피라미드 각 수준에서 최적의 옵티컬 플로우를 계산하기 위해 해를 구하는 과정을 반복하여 수행하는데, 이 과정을 멈추기 위해 필요한 문턱값(threshold)을 설정해야 한다. 두 번째로 계산 가능한 최대 풍속 값에 영향을 미치는 이미지 피라미드의 개수(m)를 결정해야 하며, 마지막으로 옵티컬 플로우 방정식을 풀기 위해 추가된 제약 조건인 탐색 영역의 크기(n)를 정해야 한다.

본 연구에서는 옵티컬 플로우 계산을 위한 최적의 파라미터 선정을 위해 문턱값과 피라미드 수준의 개수, 그리고 탐색 영역을 각각 다르게 하여 대기운동벡터를 계산하고 그 결과를 비교하였다. 문턱값은 0.02부터 0.09까지 0.01씩 증가하도록, 피라미드 수준의 개수는 2개에서 6개까지 1개씩 증가하도록 설정하였다. 탐색 영역의 크기는, n=(2w+1)×(2w+1) 화소로 정의하고, w를 3에서 13 화소까지 2씩 증가하도록 설정하였다. 실험에 사용된 자료는 Himawari-8/AHI의 IR11.2 영상이며, 기간은 2015년 8월 20일의 0000, 0600, 1200, 1800 UTC이다. 계산된 옵티컬 플로우는 ECMWF의 ERA-Interim 재분석 바람장 자료와 비교하였으며, 그 결과를 Fig. 3에 나타내었다.

OGCSBN_2020_v36n5_1_763_f0003.png 이미지

Fig. 3. Experimental results for selecting optimal parameters defined by Optical Flow algorithm. The graph shows the results according to the number of pyramid levels, the threshold, and the window size. Solid line indicates the window size (w). The number represented by the color along solid line indicates the RMSVD (ms-1) at the highest accuracy and the window size at that time. The position of the number indicates the threshold.

Fig. 3에서 알 수 있듯이 피라미드 개수가 2개에서 5개까지는 유사한 경향을 보이는 것을 볼 수 있다. 피라미드 개수가 5개까지는 Root Mean Square Vector Difference(RMSVD)가 거의 14 ms-1로 유사하게 나타나지만, 6개일 때는 12.32 ms-1로 RMSVD가 현저하게 감소함을 볼 수 있다. 문턱값의 경우 증가할 수록 정밀도가 증가하는 추세를 보이며, 탐색 영역의 크기는 피라미드 개수가 6개일 때를 제외하고 작아질수록 정밀도가 증가함을 보인다. 이러한 결과에 따라, 본 연구에서는 옵티컬 플로우 알고리즘에서 사용되는 파라미터의 개수는 6개, 문턱값은 0.09, 탐색 영역의 크기는 9를 선정하였다.

2) 대기운동벡터 산출 알고리즘

표적 선정에는 총 세 장의 영상(t=t0-Δt, t0, t0+Δt) 중기준시간(t=t0)의 영상이 이용된다. 구름을 추적하는 적외채널 대기운동벡터 특성상 표적 내 구름화소가 일정 수 이상일 때, 추적이 가능하므로, 본 연구에서는 Level 2 구름탐지 자료를 이용하여(표적 내 구름화소 비율 10% 이상) 구름이 적은 표적은 제거하였다. 고도할당 방법에는 구름 표적의 고도할당에 널리 사용되는 EBBT 방법을 사용하었다. EBBT 방법은 표적의 대표 휘도온도와 수치 예측모델과 복사모델로 계산된 표적이 있는 위치의 온도프로파일을 매칭하여 최적의 고도를 얻어내는 방법이다. 본 연구에서 표적의 대표온도는 표적 내 구름 화소 중 휘도온도가 낮은 15%의 평균 휘도온도로 정하였다(NMSC, 2012). 표적추적은 표적을 이용하여 이전(t=t0-Δt)과 이후 시간(t=t0+Δt) 영상의 탐색영역에서 표적과 가장 유사한 영역을 찾는 과정이다. 기준시간에서 정의된 표적을 이미지 매칭 방법을 이용하여 이전과 이후 시간 영상에서 각각 유사한 영역을 찾아내며, 이미지 매칭 방법으로는 CC 방법을 적용하 였다. 연속된 세 장의 영상을 이용하여 벡터를 계산하 므로, 총 두 개의 벡터가 산출되며, 최종 벡터는 두 벡터의 평균으로 결정된다. 표적의 크기는 관측하고자 하는 바람의 시·공간 규모를 변화시킬 수 있으며, 산출 정확도에도 영향을 미치는 중요한 인자이다(Hillger and Vonder Harr, 1988; Garcia-Pereda and Borde, 2014; Oh et al., 2019). 그러므로 본 연구에서는 중규모 바람장을 산출하기에 적합한 대기운동벡터의 정확도를 비교하기 위해, 표적의 크기를 12×12, 24×24, 48×48 화소로 다양 하게 변화시키며 계산하였다. 대기운동벡터가 산출되는 격자 간격은 표적 크기의 절반이 되도록 정의하였다. 추적영역의 크기는 초기 추정치를 사용하지 않을 경우 최대 93.3 ms-1 까지 산출 가능하도록 정의하였다. 초기 추정치를 사용하는 경우에는 최대 30 ms-1 까지 산출 가능하도록 추적영역을 정의하고 초기 추정치가 20 ms-1 이상일 때 추적영역의 중심이 추정되는 지점으로 이동하도록 설정하였다.

선정된 표적에서 벡터가 계산되고, 고도가 결정되면 품질관리 과정을 통해 최종 대기운동벡터가 산출된다. 대기운동벡터의 품질을 나타내기 위해 사용되는 품질 지수인 Quality Indicator (QI)는 총 다섯 가지 항의 가중 평균으로 구성된다. 세 장의 위성 영상을 통해 구한 두개의 벡터 간의 풍향, 풍속, 벡터 일관성 그리고 최종 산출된 벡터의 공간 균질성, 수치 모델 바람 자료와의 일관성 검사를 나타내는 항의 가중 평균으로 정의된다 (Holmlund, 2000). QI 는 0에서 1 사이 값을 가지며, 1에가까울수록 산출된 대기운동벡터가 양질의 벡터임을 의미한다.

3) 실험방법

초기 추정치 적용 방법을 평가하기 위해, (1) 초기 추정치 없는 경우 (2) 기상예보장을 초기 추정치로 사용하는 경우, 그리고 (3) 옵티컬 플로우로 얻어진 바람장을 초기 추정치로 사용하여 대기운동벡터를 계산한 경우에 대해 ECMWF의 ERA-Interim 재분석장과 비교하였다. 초기 추정치 종류 및 유무에 따라 최적 표적 크기가 달라질 수 있으므로, 12×12, 24×24, 48×48 화소로 표적 크기를 변화시키며 계산하였다. 초기 추정치 사용을 제외한 모든 대기운동벡터의 산출 과정은 동일하며, 초기 추정치 사용 여부에 따라 탐색 영역의 중심 위치와 크기만 변경되었다. 초기 추정치를 사용할 경우에는 바람의 추정 위치를 탐색 영역을 중심으로 정의하므로 그 크기가 작은 반면, 초기 추정치를 사용하지 않는 경우에는 풍속이 강한 바람까지 계산하기 위해 탐색 영역의 크기가 보다 커야 한다. 본 연구에서는 초기 추정치를 사용하지 않는 경우에 대해서 최대 93.3 ms-1 까지 계산이 가능하도록 추적 영역의 크기를 결정하였다. Himawari-8/AHI 위성의 공간 해상도가 2 km 이므로, 초기 추정치를 사용하지 않는 경우에는 탐색 영역의 크기를 표적 크기가 12×12, 24×24, 48×48 화소로 변할 때 각각 68×68, 80×80, 104×104, 초기 추정치를 사용하는 경우에는 각각 30×30, 42×42, 66×66 화소로 정의하되, 추정치가 20 ms-1 이상일 때는 탐색 영역의 중심을 추정 위치로 이동 시켰다.

4. 결과 및 분석

3장에서 소개한 대기운동벡터 산출 알고리즘을 동일 하게 적용하고, 초기 추정치로 UM 자료와 옵티컬 플로우를 사용한 경우와 초기 추정치를 사용하지 않은 경우로 나누어 대기운동벡터를 산출하였다. 각 방법에 따라 하나의 장면에 대한 최종 대기운동벡터를 산출하는데 걸린 시간을 비교해보면, 초기 추정치로 UM 자료를 이용했을 때 13분 14.7초가 걸렸으며 동일한 조건에서 옵티컬 플로우를 사용한 경우에는 15분 53.2초로 거의 유사하게 나타났다. 반면, 초기 추정치를 사용하지 않은 경우 대기운동벡터 계산에 걸린 시간은 1시간 28분 38.9초로 초기 추정치를 사용한 경우보다 최대 569% 길어지는 것을 확인할 수 있다.

각 방법으로 산출된 대기운동벡터를 검증하기 위해, QI가 0.8 이상인 벡터를 ECMWF의 ERA-Interim 재분석 바람장과 비교하였다. 격자 자료인 ERA-Interim 재분석 바람장과 그렇지 않은 대기운동벡터의 1:1 검증을 위해서, 대기운동벡터를 기준으로 수평과 연직 방향으로 가장 가까운 재분석 바람장 8개의 지점을 선택하고, 수평은 이중 선형 보간법을, 연직은 기압에 대한 log-선형 보간법을 사용하였다. 검증에는 대기운동벡터의 평균 풍속을 나타내는 SPD, 대기운동벡터와 재분석 바람장 풍속의 차이를 나타내는 Bias, 벡터 차이를 나타내는 RMSVD, 그리고 Bias, RMSVD를 각각 재분석 바람장의 풍속으로 나누어 정규화한 Normalized BIAS (NBias), Normalized RMSVD (NRMSVD)가 사용되었다. Table. 1은 초기 추정치를 다르게 적용하여 산출된 QI가 0.80 이상인 대기운동벡터를 ERA-Interim 자료와 검증한 결과를 수치적으로 나타낸다. 대기운동벡터의 개수는 초기 추정치를 사용하지 않은 경우와 UM 자료를 사용한 경우 모두 유사하게 나타나며, 표적의 크기가 커질수록 벡터의 개수가 적어지는 경향을 보인다. 대기운동벡터 풍속은 대체로 표적의 크기가 작을 때 빠르며, 초기 추정치로 옵티컬 플로우를 이용하여 산출된 대기운동벡터의 풍속이 전체적으로 느리게 나타나는 것을 알 수있다. Bias와 NBias의 경우 표적이 작을 수록 높은 정확도를 보이며, 방법별 차이는 크지 않은 결과를 보였다.

Table 1. Validation results for ERA-Interim data and AMVs using IR11.2 images of Himawari-8 satellite with different target size between 12×12 and 48×48. No FG indicates the case not using the first guess to define the search area, FG (UM) indicates the case where the UM data is used as the first guess, and FG (OF) indicates the optical flow is used as the first guess

OGCSBN_2020_v36n5_1_763_t0001.png 이미지

RMSVD는 5-6 ms-1 의 분포를 보이고, NRMSVD는 0.3-0.4의 분포를 보였으며, 초기 추정치로 옵티컬 플로우를 사용하였을 때 가장 좋은 결과를 보였다. 즉, 대기운동 벡터의 정확도는 초기 추정치 사용 여부에 따른 큰 차이가 없었지만, 정밀도는 옵티컬 플로우를 초기 바람 추정치로 사용하였을 때 가장 좋은 결과를 보였다. Fig. 4는 초기 추정치를 다르게 적용하여 산출된 대기 운동벡터를 ERA-Interim 자료와 검증한 결과를 QI에따라 나타낸 것으로 QI가 각각 0.80, 0.85, 0.90, 0.95 이상일 때 표적 크기에 따른 대기운동벡터의 (a) 개수, (b) 풍속 Bias, (c) RMSVD를 각각 보여준다. Fig. 4(a)는 IR11.2 영상을 이용하여 산출된 대기운동벡터의 개수 분포이다. 초기 추정치 사용 유무에 따른 개수의 차이는 표적의 크기가 가장 작을 때 그 차이가 크게 나타나며, 표적의 크기가 24×24 이상부터는 초기 추정치를 사용하지 않고 산출된 경우와 UM 자료를 사용하여 산출된 대기운동 벡터의 개수가 유사하게 나타난다. 표적 크기가 12×12 일 때 전체적으로 초기 추정치로 UM 자료를 사용하여 산출된 대기운동벡터가 가장 많게 나타나지만, 표적 크기가 커질수록 초기 추정치를 사용하지 않고 산출된 대기운동벡터의 개수가 높게 나타나며 UM 자료를 사용 하여 산출된 대기운동벡터보다 최대 3.6% 많게 나타난다. 초기 추정치로 옵티컬 플로우를 사용하여 산출된 대기운동벡터는 다른 방법으로 산출된 대기운동벡터에 비해 적게 나타나는 특징을 보인다. UM 자료를 이용하여 산출된 대기운동벡터와 비교하면 QI가 0.80 이상이 면서 표적의 크기가 각각 12×12일 때 38%, 24×24일 때 25%, 48x48일 때 5% 적게 산출되었다.

OGCSBN_2020_v36n5_1_763_f0004.png 이미지

Fig. 4. The statistics of AMVs derived from No FG, FG(UM), and FG(OF) in IR11.2. (a) Number, (b) Windspeed Bias, (c) RMSVD.

Fig. 4(b)는 대기운동벡터의 풍속 Bias를 나타내고, 대체로 표적의 크기가 작을수록 풍속 bias가 작아짐을 알수 있다. 초기 추정치에 따라 산출된 대기운동벡터의 풍속 bias를 비교해보면, QI가 0.80 이상일 때는 옵티컬 플로우를 활용한 대기운동벡터의 정확도가 가장 낮게 나타나지만 QI가 증가할수록 정확도가 높아진다. 특히 QI가 0.95 이상이고 표적의 크기가 24×24일 때는 BIAS가 -0.001 ms-1 로 정확도가 가장 높게 나타나는 것을 알수 있다.

Fig. 4(c)는 대기운동벡터의 RMSVD를 나타낸다. QI 가 증가할수록 RMSVD가 감소하는 경향을 보이며 QI 가 0.80 이상일 때 각 산출 방법, 표적의 크기에 따른 RMSVD 차이가 가장 크게 나타난다. 표적의 크기가 24×24이면서 초기 추정치를 사용하지 않고 산출된 대기운동벡터의 RMSVD는 6.467 ms-1 로, 옵티컬 플로우를 이용한 경우는 5.296 ms-1 로 초기 추정치를 사용하지 않은 경우가 최대 22.1%까지 높게 나타나는 것으로 보인다. 가장 높은 정밀도는 모두 QI가 0.95 이상일 때이며, 그 중에서도 옵티컬 플로우를 활용하고 24×24인 표적 으로 산출된 대기운동벡터가 RMSVD 2.710 ms -1 로 가장 높은 정밀도를 보인다.

5. 요약 및 결론

본 연구에서는 초기 바람 추정치가 대기운동벡터의 산출 과정인 추적 과정에 미치는 효과를 검증하기 위하여 초기 추정치를 다르게 적용하여 그 결과를 비교 하였다. 적용된 방법은 총 세 가지로 표적 추적 과정에 초기 추정치를 사용하지 않은 경우, 그리고 초기 추정 치로 UM 자료 및 지금까지 시도된 적이 없었던 옵티컬 플로우를 활용한 경우이다. 표적 추적 이외의 나머지 산출 과정은 모두 동일하게 적용하여 대기운동벡터를 산출하였다. 표적 추적 과정의 변화에 따라 최적 표적 크기가 달라질 수 있으므로 표적의 크기를 12×12, 24×24, 48×48로 각각 변화시켜 대기운동벡터를 산출하였으며, 각 결과는 ERA-Interim 자료를 이용하여 검증하였다.

연구결과를 통해서 초기 추정치를 사용하지 않고 산출된 대기운동벡터의 정확도는 모든 방법에서 큰 차이가 없었지만, 정밀도는 초기 추정치로 옵티컬 플로우를 적용하였을 때 가장 높은 것을 확인할 수 있다. 즉 보편적으로 사용하고 있는 대기운동벡터 산출 방법인 UM 모델 자료를 초기 추정치로 사용하지 않고, 옵티컬 플로우 방법을 대신 사용하면, 현재보다 정확도가 향상될 수있으며, 모델 의존성 문제를 해결할 수 있음을 보였다. 또한 계산 속도는 UM 모델 자료를 초기 추정치로 사용 하는 경우와 비교하여 큰 차이 없었고, 초기바람장을 사용하지 않는 경우보다는 대폭 빨라짐을 확인하였다. 한편 옵티컬 플로우를 초기 추정치로 산출된 대기운동벡터 개수는 UM 자료를 이용한 결과와 비교하였을 때 적게 산출되는 경향을 보이며, 표적크기가 24×24일 때 25%의 비율로 비교적 크게 줄어들어 추후 자료동화 활용 측면에서 불리할 수 있다. 그러나 옵티컬 플로우를 이용하면 대기운동벡터 알고리즘 중 영상 추적과정에서 수치 모델의 개입을 배제할 수 있고, 이는 궁극적으로 대기운동벡터의 관측자료로서의 독립성을 높이는 것이다. 또한 수치모델 바람장을 초기 추정치로 사용하는 과정에서 발생하던 고도할당 과정의 오차가 추적과정까지 이어지는 문제가 옵티컬 플로우 기반 바람장을 사용함에 따라 제거되어 정확도 개선의 여지도 커진다. 따라서 옵티컬 플로우를 초기 추정치로 활용하면 대기 운동벡터 산출 알고리즘의 개선이 가능할 것으로 판단된다.

본 연구는 밤과 낮에 동시에 적용하기 위하여 적외선 영역에 한정하였으나 가시채널을 추가한다면 보다 향상된 대기운동벡터 도출이 가능할 것으로 기대된다. 또한 옵티컬 플로우의 한계성을 보완하여 구름과 같은 비강체에 바로 적용하여 대기운동벡터를 구하는 추가 연구를 지속할 계획이다.

사사

이 논문은 부산대학교 기본연구지원사업(2년)에 의하여 연구되었음.

References

  1. Apke, J. M., J. R. Mecikalski, and C. P. Jewett, 2016. Analysis of mesoscale atmospheric flows above mature deep convection using super rapid scan geostationary satellite data, Journal of Applied Meteorology and Climatology, 55(9): 1859-1887. https://doi.org/10.1175/JAMC-D-15-0253.1
  2. Bedka, K. M. and J. R. Mecikalski, 2005. Application of satellite-derived atmospheric motion vectors for estimating mesoscale flows, Journal of Applied Meteorology, 44(11): 1761-1772. https://doi.org/10.1175/JAM2264.1
  3. Bessho, K., K. Date, M. Hayashi, A. Ikeda, T. Imai, H. Inoue, Y. Kumagai, T. Miyakawa, H. Murata, and T. Ohno, 2016. An introduction to Himawari-8/9-Japan's new-generation geostationary meteorological satellites, Journal of the Meteorological Society of Japan, 94(2): 151-183. https://doi.org/10.2151/jmsj.2016-009
  4. Borde, R. and M. Doutriaux-Boucher, 2014. Extraction des vecteurs vents a partir d'images satellite, La Meteorologie, 87: 27-33.
  5. Borde, R. and J. Garcia-Pereda, 2014. Impact of wind guess on the tracking of atmospheric motion vectors, Journal of Atmospheric and Oceanic Technology, 31(2): 458-467. https://doi.org/10.1175/JTECH-D-13-00105.1
  6. Bresky, W. C., J. M. Daniels, A. A. Bailey, and S. T. Wanzong, 2012. New methods toward minimizing the slow speed bias associated with atmospheric motion vectors, Journal of Applied Meteorology and Climatology, 51(12): 2137-2151. https://doi.org/10.1175/JAMC-D-11-0234.1
  7. Forsythe, M., 2007. Atmospheric motion vectors: past, present and future, Presentation Paper of ECMWF Annual Seminar, Reading, UK, Sep. 3-7, pp. 1-79.
  8. Garcia-Pereda, J. and R. Borde, 2014. The impact of the tracer size and the temporal gap between images in the extraction of atmospheric motion vectors., Journal of Atmospheric and Oceanic Technology, 31(8): 1761-1770. https://doi.org/10.1175/JTECH-D-13-00235.1
  9. Hillger, D. W. and T. H. Vonder Haar, 1988. Estimating noise levels of remotely sensed measurements from satellites using spatial structure analysis, Journal of Atmospheric and Oceanic Technology, 5(2): 206-214. https://doi.org/10.1175/1520-0426(1988)005<0206:ENLORS>2.0.CO;2
  10. Holmlund, K., 1998. The utilization of statistical properties of satellite-derived atmospheric motion vectors to derive quality indicators, Weather and Forecasting, 13(4): 1093-1104. https://doi.org/10.1175/1520-0434(1998)013<1093:TUOSPO>2.0.CO;2
  11. Holmlund, K., 2000. The atmospheric motion vector retrieval scheme for meteosat second generation, Proc. of 5th International Winds Workshop, Lorne, Australia, Feb. 28-Mar. 3, pp. 201-208.
  12. Horn, B. K. and B. G. Schunck, 1981. Determining optical flow, Proc. of Techniques and Applications f Image Understanding Symposium, Washington, D.C., USA, Apr. 21-22, vol. 0281, pp. 185-203.
  13. Kim, M., H. M. Kim, J. Kim, S.-M. Kim, C. Velden, and B. Hoover, 2017. Effect of enhanced satellite-derived atmospheric motion vectors on numerical weather prediction in East Asia using an adjointbased observation impact method, Weather and Forecasting, 32(2): 579-594. https://doi.org/10.1175/WAF-D-16-0061.1
  14. Kim, S., J.-H. Park, M.-L. Ou, H. Cho, and E.-H. Sohn, 2012. Optimization of mesoscale atmospheric motion vector algorithm using Geostationary Meteorological Satellite Data, Atmosphere, 22(1): 1-12 (in Korean with English abstract). https://doi.org/10.14191/Atmos.2012.22.1.001
  15. Langland, R. H., C. Velden, P. M. Pauley, and H. Berger, 2009. Impact of satellite-derived rapid-scan wind observations on numerical model forecasts of Hurricane Katrina, Monthly Weather Review, 137(5): 1615-1622. https://doi.org/10.1175/2008MWR2627.1
  16. Le Marshall, J., Y. Xiao, D. Howard, C. Tingwell, J. Freeman, P. Gregory, T. Le, D. Margetic, T. Morrow, and J. Daniels, 2016. First results from the generation and assimilation of 10-minute atmospheric motion vectors in the Australian region, Journal of Southern Hemisphere Earth Systems Science, 66(1): 12-18. https://doi.org/10.22499/3.6601.003
  17. Lee, J.-K. and C.-J. Park, 2007. Algorithm for Arbitrary Point Tracking using Pyramidal Optical Flow, Journal of Korea Multimedia Society, 10(11): 1407-1416 (in Korean with English abstract).
  18. Lee, S., K. Salonen, and N. Bormann, 2015. Assessment of AMVs from COMS in the ECMWF system, ECMWF (European Centre for Medium-Range Weather Forecasts), Reading, Berkshire, UK.
  19. Lucas, B. D. and T. Kanade, 1981. An iterative image registration technique with an application to stereo vision, Proc. of 7th International Joint Conference on Artificial Intelligence, Vancouver, BC, Aug. 24-28, pp. 674-679.
  20. Mecikalski, J. R. and K. M. Bedka, 2006. Forecasting convective initiation by monitoring the evolution of moving cumulus in daytime GOES imagery, Monthly Weather Review, 134(1): 49-78. https://doi.org/10.1175/MWR3062.1
  21. National Meteorological Satellite Center (NMSC), 2012. Atmospheric Motion Vector (AMV) Algorithm Theoretical Basis Document, National Meteorological Satellite Center, Jincheon, South Korea, pp. 1-30.
  22. Oh, S. M., R. Borde, M. Carranza, and I.-C. Shin, 2019. Development and Intercomparison Study of an Atmospheric Motion Vector Retrieval Algorithm for GEO-KOMPSAT-2A, Remote Sensing, 11(17): 2054. https://doi.org/10.3390/rs11172054
  23. Oh, Y., J. H. Kim, H. Park, and K. Baek, 2015. Development and analysis of COMS AMV target tracking algorithm using gaussian cluster analysis, Atmosphere, 31(6): 531-548 (in Korean with English abstract).
  24. Oyama, R., M. Sawada, and K. Shimoji, 2018. Diagnosis of tropical cyclone intensity and structure using upper tropospheric Atmospheric Motion Vectors, Journal of the Meteorological Society of Japan, 96B: 3-26. https://doi.org/10.2151/jmsj.2017-024
  25. Salonen, K., J. Cotton, N. Bormann, and M. Forsythe, 2015. Characterizing AMV height-assignment error by comparing best-fit pressure statistics from the Met Office and ECMWF data assimilation systems, Journal of Applied Meteorology and Climatology, 54(1): 225-242. https://doi.org/10.1175/JAMC-D-14-0025.1
  26. Shimoji, K., 2014. Motion tracking and cloud height assignment methods for Himawari-8 AMV, Proc. of 12th International Winds Workshop, Copenhagen, Denmark, Jun. 16-20, pp. 80-90.
  27. Velden, C., J. Daniels, D. Stettner, D. Santek, J. Key, J. Dunion, K. Holmlund, G. Dengel, W. Bresky, and P. Menzel, 2005. Recent innovations in deriving tropospheric winds from meteorological satellites, Bulletin of the American Meteorological Society, 86(2): 205-224. https://doi.org/10.1175/BAMS-86-2-205
  28. Wu, Q., H.-Q. Wang, Y.-J. Lin, Y.-Z. Zhuang, and Y. Zhang, 2016. Deriving AMVs from geostationary satellite images using optical flow algorithm based on polynomial expansion, Journal of Atmospheric and Oceanic Technology, 33(8): 1727-1747. https://doi.org/10.1175/JTECH-D-16-0013.1
  29. Wu, T.-C., H. Liu, S. J. Majumdar, C. S. Velden, and J. L. Anderson, 2014. Influence of assimilating satellite-derived atmospheric motion vector observations on numerical analyses and forecasts of tropical cyclone track and intensity, Monthly Weather Review, 142(1): 49-71. https://doi.org/10.1175/MWR-D-13-00023.1