Improvement of 2-pass DInSAR-based DEM Generation Method from TanDEM-X bistatic SAR Images

TanDEM-X bistatic SAR 영상의 2-pass 위성영상레이더 차분간섭기법 기반 수치표고모델 생성 방법 개선

Abstract

The 2-pass DInSAR (Differential Interferometric SAR) processing steps for DEM generation consist of the co-registration of SAR image pair, interferogram generation, phase unwrapping, calculation of DEM errors, and geocoding, etc. It requires complicated steps, and the accuracy of data processing at each step affects the performance of the finally generated DEM. In this study, we developed an improved method for enhancing the performance of the DEM generation method based on the 2-pass DInSAR technique of TanDEM-X bistatic SAR images was developed. The developed DEM generation method is a method that can significantly reduce both the DEM error in the unwrapped phase image and that may occur during geocoding step. The performance analysis of the developed algorithm was performed by comparing the vertical accuracy (Root Mean Square Error, RMSE) between the existing method and the newly proposed method using the ground control point (GCP) generated from GPS survey. The vertical accuracy of the DInSAR-based DEM generated without correction for the unwrapped phase error and geocoding error is 39.617 m. However, the vertical accuracy of the DEM generated through the proposed method is 2.346 m. It was confirmed that the DEM accuracy was improved through the proposed correction method. Through the proposed 2-pass DInSAR-based DEM generation method, the SRTM DEM error observed by DInSAR was compensated for the SRTM 30 m DEM (vertical accuracy 5.567 m) used as a reference. Through this, it was possible to finally create a DEM with improved spatial resolution of about 5 times and vertical accuracy of about 2.4 times. In addition, the spatial resolution of the DEM generated through the proposed method was matched with the SRTM 30 m DEM and the TanDEM-X 90m DEM, and the vertical accuracy was compared. As a result, it was confirmed that the vertical accuracy was improved by about 1.7 and 1.6 times, respectively, and more accurate DEM generation was possible with the proposed method. If the method derived in this study is used to continuously update the DEM for regions with frequent morphological changes, it will be possible to update the DEM effectively in a short time at low cost.

2-pass DInSAR (Differential Interferometric SAR) 기법을 이용한 DEM 생성 방법은 SAR 영상쌍 간의 정합, 간섭도 생성, 위상 unwrapping, DEM오차 계산, 좌표계 변환 등 복잡한 단계가 필요하다. 이러한 각 단계별 자료처리 정확도는 최종적으로 생성되는 DEM의 성능에 영향을 미치게 된다. 본 연구에서는 TanDEM-X bistatic SAR 영상의 2-pass DInSAR 기법 기반 DEM 생성 방법에 대한 성능 향상을 위한 개선된 방법을 개발하였다. 개발된 DEM 생성 방법은 unwrapped 위상 내의 DEM 오차와 좌표계 변환 시 발생할 수 있는 DEM오차를 현저히 줄일 수 있는 방법이다. 개발된 알고리즘의 성능 분석은 GPS 측량으로부터 생성한 지상기준점(Ground Control Point, GCP)을 이용하여 기존의 방법과 새로 제안된 알고리즘 적용 결과의 수직정확도(Root Mean Square Error, RMSE)를 비교하여 수행하였다. 실제 unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하지 않고 생성한 DInSAR 기반 DEM의 고도 오차는 수직정확도는 39.617 m로 관측되었고, 제안한 방법을 통하여 생성한 DEM의 수직정확도는 2.346 m로 향상됨을 확인할 수 있었다. 제안하는 2-pass DInSAR 기법을 통해 reference로 사용한 SRTM 30 m DEM(수직정확도 5.567 m)에 DInSAR로 관측한 SRTM DEM 오차를 보상하여 최종적으로 공간해상도는 약 5배, 수직 정확도는 약 2.4배 향상된 DEM을 생성할 수 있었다. 또한, 제안한 방법을 통하여 생성한 DEM의 공간해상도를 SRTM 30 m DEM과 TanDEM-X 90 m DEM과 일치시키고 수직정확도를 비교한 결과 각각 약 1.7배 및 1.6배 향상되어 제안하는 2-pass DInSAR 기반 DEM 생성 방법으로 보다 정확한 DEM 생성이 가능함을 확인할 수 있었다. 빈번한 형태학적 변화를 갖는 지역에 대한 DEM의 지속적인 업데이트를 위하여 본 연구에서 도출한 방법을 이용한다면 저비용으로 빠른 시간 내에 효과적으로 DEM을 갱신할 수 있을 것이다.

1. 서론

위성영상레이더(Synthetic Aperture Radar, SAR)는 위성체로부터 마이크로파 신호를 보내고 지표물질에 반사되어 되돌아온 신호를 통해 영상화하는 원격탐사 시스템이다. 위성영상레이더가 취득하는 신호는 복소수 형태(Single Look Complex, SLC)로 기록되어지며, 지표의 거칠기 및 유전율에 의하여 결정되는 반사도 정보의 진폭 값과 지표까지의 거리정보를 나타내는 위상 값을 포함한다. 이러한 위성영상레이더를 통해 관측되는 대상에 따른 반사 신호 특성과 그 차이를 이용하여 지표면에 대한 수치표고모델(Digital Elevation Model, DEM) 을 생성할 수 있다. SAR 영상을 이용하여 DEM을 생성하는 방법은 아래와 같이 네 가지로 나눌 수 있다.

1) 위성영상레이더 측량기법(Radargrammetry), 인접한 위성궤도를 따라 서로 다른 각도로 촬영된 동일 지역의 Stereo SAR 입체영상에서 발생하는 시차(parallax)를 이용하여 DEM 생성하는 방법;

2) 위성영상레이더 간섭기법(Interometric SAR, InSAR), 동일한 지역에 대하여 다른 시기에 촬영된 두 개의 SAR 영상의 위상정보를 이용하여 DEM 생성하는 방법;

3) 위성영상레이더 Radarclinometry, 주로 행성탐사에 활용하는 방법으로 한 장의 SAR 영상의 밝기 값을 이용하여 지표면의 경사각을 계산하고(Shape-from-shading) 이를 이용하여 DEM 생성하는 방법;

4)위성영상레이더편파간섭기법(PolInSAR,Polarimetric InSAR), full-polarization을 가지는 InSAR 영상쌍에 대하여 지표의 고도에 대한 각 편광별 산란 행렬 및 산란 벡터의 공분산 행렬 등을 이용하여 적절한 물리적 역산을 통해 DEM 생성하는 방법;

이 중, InSAR 기법은 지난 1990년대 이후 급속하게 발전한 마이크로파 원격탐사의 핵심 기술 중 하나로 동일한 지역에 대하여 다른 시기에 촬영된 두 개의 SAR 영상의 위상정보를 이용하여 지표면에 대한 DEM 생성 및 지표면의 변위를 관측할 수 있다(Rogers and Ingals, 1969; Graham, 1974; Rodriguez and Martin, 1992; Massonnet and Feigl, 1998). InSAR 기법을 이용한 DEM 생성에 관한 연구는 1980년대 후반부터 본격적으로 시작되어(Gabriel and Goldstein, 1988), 이후 DEM 생성을 위한 다양한 방법이 제안되었다. InSAR 기반 DEM 생성 방법은 1) SAR 영상쌍 정합(Co-registration), 2) 간섭도(Interferogram) 생성, 3) 위상 필터링(Phase Filtering), 4) 위상 불구속화(Phase Unwrapping), 5) 고도변환(Phase to Height Conversion), 6) 좌표계 변환(Geo-coding) 등 매우 복잡한 단계를 필요로 한다(Geudtner and Schwabisch, 1996; Kim, 2012). 이러한 알고리즘의 각 단계별 자료 처리 정확도는 최종적으로 생성되는 DEM의 성능에 영향을 미치게 된다. 따라서 각 단계별 자료 처리 정확도를 향상하기 위하여 정밀 영상정합, 위상 필터링, 간섭도 내 이온층 및 대기 오차 저감 등 다양한 연구가 수행되고 있다(Zebker et al., 1994; Rufino et al., 1998; Liu et al., 2020; Saied et al., 2020; Toutin and Gray, 2000).

특히 DEM 생성 성능 향상을 위하여 2-pass 차분간섭 기법(Differential InSAR, DInSAR)을 이용하는 방법이 제안되었다(Seymour, 1999; Yoon et al., 2001; Kim et al., 2005). 2-pass DInSAR(Massonnet et al., 1993)를 이용한 DEM 생성은 InSAR 조건을 만족하는 두 장의 SAR영상과 외부 DEM을 사용하여 차분간섭도를 생성하고 이차분간섭도에 남아 있는 잔여 위상에 의한 DEM오차를 다시 외부 DEM에 보상하여 정밀 DEM을 생성하는 방법이다. 그러나 기존의 연구에서는 위상 불구속화 및 좌표계 변환에서 발생할 수 있는 오차에 대한 정량적 분석이 수행되지 않았으며, 또한 이러한 오차에 대한 저감을 고려한 2-pass DInSAR 기반 DEM 생성 방법이 정립되지 않았다. 이러한 오차가 야기하는 고도 오차는 수십 m 이상으로 일반적인 허용 한계보다 커지게 된다. 따라서 본 연구에서는 위상 불구속화 및 좌표계 변환 시 발생할 수 있는 DEM오차를 고려한 개선된 DEM 생성 방법을 개발하였다. 개발된 알고리즘의 성능 분석은 GPS 측량으로부터 생성한 지상기준점(Ground Control Point, GCP)을 이용하여 기존의 방법과 새로 제안된 알고리즘 적용 결과의 수직정확도(Root Mean Square Error, RMSE)를 비교하여 수행하였다.

2. 연구지역 및 자료

Fig. 1은 본 연구의 연구지역인 에트나(Etna) 산을 나타낸다. 에트나 산의 위치는 북위 37.755064°, 동경 14.995246°으로 이탈리아 시칠리아 섬에 존재하는 해발 3,350미터의 유럽에서 가장 높은 활화산이다. 지난 15년 동안 에트나 화산은 빈번한 분화활동을 보였다. 2005년 이후 주목할 만한 분화는 2006-2007년과 2008-2009년 사이(1991-1993년 분화 이후 Etna의 시간적으로 가장 긴 화산분출 활동), 2011년, 2014년, 2015년, 2017년 및 2018년 12월에 기록되었다(Palaseanu et al., 2019). 이러한 빈번한 분화활동은 에트나 화산 정상 지역의 급격한 형태학적 변화를 일으켰다. 특히 분출 현상의 빈도는 지난 몇 년 동안 증가하였기 때문에 에트나 화산 지형의 지속적인 업데이트는 중요하다할 수 있다. 이러한 중요성에 Bisson et al.(2015)의 연구에서는 airborne LiDAR 원격탐사 자료, 그리고 Palaseanu et al.(2019)의 연구에서는 Pleiades 광학위성의 사진측량기법(photogrammetry)을 이용하여 에트나 화산에 대한 정밀 지형정보를 생성하였다.

Fig. 1. Study Area location with the acquired SAR Image footprint (red outline) and location and classification of the GPSs from the geodetic monitoring network of Mt Etna (Bisson et al., 2015).

본 연구에서는 에트나 화산 지역을 대상으로 2-pass DInSAR 기반 DEM을 정밀하게 생성하는 방법을 제안하고자 하며, 이를 위하여 연구지역에 대한 TanDEM-X(TerraSAR adds-on Digital Elevation Measurement) SAR 영상쌍을 획득하였다(Fig. 1의 붉은색 박스). Table 1은 획득한 TanDEM-X 영상쌍에 대한 정보를 나타낸다. TanDEM-X는 두 대의 X 밴드 위성으로 구성되어, 일정 거리를 유지하며 동시에 SAR 영상을 취득 가능한 singlepass interferometry 미션이다(Krieger et al., 2007). 영상은 하나의 위성이 신호를 지표로 보내고, 두 개의 위성이 동시에 신호를 받는 bistatic 모드로 거의 동시에 관측되기 때문에 interferogram 생성 시, temporal decorrelation을 완벽히 제거 가능하므로, 높은 긴밀도(coherence)가 측정되어, 최종 결과물인 DEM의 수직정밀도를 획기적으로 향상 시킬 수 있었다(Lee, 2020). 획득한 TanDEM-X 영상쌍의 한 fringe(위상의 변화)를 변화시키는데 필요한 고도의 변화인 고도 민감도는 약 65.45 m이다. 따라서 해당 TanDEM-X 영상쌍은 높은 긴밀도와 고도 민감도를 가진 자료로 정밀 DEM 제작이 가능하다.

Table 1. Summary of TanDEX-X SAR pair

Fig. 2는 획득한 외부 DEM을 나타낸다. 2-pass DInSAR 기반 DEM 생성을 위하여 외부 DEM이 필요하며, 이를 위하여 본 연구에서는 오픈소스로 자료 수급이 용이하고 글로벌 DEM으로 DInSAR 처리에 가장 널리 활용되는 SRTM(Shuttle Raar Topography Mission) DEM을 이용하였다(Fig. 2(a)). SRTM은 2000년 초반에 전지구의 약 80%의 DEM을 구축하여 데이터를 오픈소스로 제공하고 있다. SRTM DEM은 2000년에 구축한 자료로서 다양한 오차를 포함하고 있으며, 대표적으로 산등성이나 협곡 지역의 Void(or Hole) 오차를 포함한다. 이에 최신에 구축된 외부의 DEM을 이용하여 SRTM DEM의 Void를 보간(interpolation)하여 DEM을 업데이트 하고 있다. 본 연구에서는 사용한 SRTM DEM은 Version 3.0의 해상도 30 m 자료이다. 이 자료는 미국 지질 조사 웹사이트(https://earthexplorer.usgs.gov/)에서 획득할 수 있다. 추가적으로 제안한 방법으로 생성한 DInSAR 기반 DEM의 비교 검증을 위하여 TanDEM-X 90 m DEM을 획득하였다(Fig. 2(b)). TanDEM-X 90 m DEM은 TerraSARX 및 TanDEM-X 영상으로 생성한 공간해상도 12 m의 WorldDEM™(수직정확도 2 m)을 공간해상도 90 m로 다운샘플링하여 독일 항공우주센터 웹 사이트(https://geoservice.dlr.de/web/dataguide/tdm90/)에서 무료로 배포하고 있다. 공식적인 TanDEM-X 90 m DEM의 수직 정확도는 발표되지 않았으나 최근 연구에 따르면 TanDEM-X DEM의 범람지(floodplain)에서의 수직정확도는 3.10 m, 수목지역에서는 5.68 m, 경사가 가파른(경사도 10-15°) 산림지역에서 수직정확도는 6.07 m. 경사가 매우 가파른(경사도 15° 이상) 산림지역에서 수직정확도는 6.42 m로 관측되었다(Hawker et al., 2019). SRTM DEM 및 TanDEM-X DEM의 정량적 및 정성적 평가에 대해서는 Grohmann(2018)의 연구를 참고하도록 한다.

Fig. 2. External DEM used in the study.

TanDEM-X SAR영상을 이용한 2-pass DInSAR 기반 DEM을 검증하기 위하여 에트나 산에 운용 중인 56개의 측지 관측용 GPS 네트워크(geodetic monitoring GPS network)에서 2005년에 관측한 GCP(Ground Control Point)를 이용하였다(Bonforte and Puglisi, 2006). 본 연구에서는 1) TanDEM-X 영상 범위 안에 존재하며 균질한 GCP의 공간 분포, 2) GCP를 획득한 2005년과 TanDEM-X 영상 촬영한 2012년 사이의 화산활동으로 인한 형태학적 변화 또는 토지 사용 변화의 영향을 받지 않는 지점 상에 존재하는지를 고려하여 총 41개(Fig. 1의 별표)의 GCP를 선별하여 검증에 활용하였다(Palaseanu et al., 2019).

3. 2-pass DInSAR 기반 DEM 생성 방법 개선

2-pass DInSAR 기반 DEM 제작을 위하여 필요한 기초 이론은 Yoon et al. (2001) 및 Kim et al. (2005)의 연구를 참고하도록 한다. 본 연구에서는 제안하는 알고리즘의 자료처리 과정을 중심으로 설명을 하고자 한다. Fig. 3은 제안하는 알고리즘의 자료처리 과정을 나타낸다. 주요 과정은 크게 5가지로 요약된다. 1) 외부 DEM 준비, 2) InSAR 처리, 3) DInSAR 처리, 4) 위상 불구속화, 5) DInSAR 기반 DEM 생성이다.

Fig. 3. Detailed processing work flow for the proposed DInSAR based DEM generation method.

1) 외부 DEM 준비

DInSAR 처리 및 DEM 생성을 위하여 외부 DEM을 준비하였다. 본 연구에서는 SRTM DEM을 이용하였다. 획득한 SRTM DEM은 표고(Orthometric height) 자료로 WGS84 타원체고(Ellipsoid height)로 변환하기 위하여 지오이드고(Geoid height)를 더해주어야 한다. 미국 UNAVCO 웹 사이트(https://www.unavco.org/)에서 제공하는 geoid calculator를 이용하여 지오이드고를 획득하였으며, 이 오프셋을 원 SRTM DEM에 보상하여 WGS84 타원체고인 SRTM DEM을 계산하였다. 이렇게 준비된 SRTM WGS84 타원체고 DEM에 DInSAR 방법을 활용하여 정밀고도에 대한 업데이트를 수행한다.

2) InSAR 처리

두 SAR 영상에 대해 InSAR 처리를 통하여 interferogram을 생성하기 위해서는 두 SAR 영상쌍을 정합(coregistration)해야 한다. 서로 다른 궤도를 지나면서 관측한 SAR 영상쌍에서는 지표 위의 동일한 산란체에 대하여 서로 다른 좌표의 화소에 위치하게 되며, 이러한 다른 두 궤도에 의해서 생성된 영상을 완벽하게 동일한 위치에 놓이게 하는 것이 co-registration이다. 본 연구에서는 진폭 상관계수를 이용하여 Coarse-to-Fine 방법으로 co-registration을 수행하였다(Kwoh et al., 1994). 정합 결과를 바탕으로 부영상을 주영상에 맞추어 재배열(resampling)하고 두 SAR 영상 간의 위상차를 계산함으로써 interferogram을 획득할 수 있다. 이때, 긴밀도를 높이기 위하여 range 및 azimuth방향으로 각각 4룩(look)의 멀티룩 처리를 수행하였다. 이에 interferogmram의 픽셀 크기는 range 및 azimuth방향으로 약 5.46 m×8.82 m이다.

이러한 interferogram에는 지구 곡률에 의한 위상, 지형고도에 의한 위상, 지표변위에 의한 위상, 대기 및 이온층에 의한 위상, 노이즈 위상 등을 포함한다(Zebker and Goldstein, 1986). 이러한 위상 성분들에 대한 정확한 분리를 통한 지형고도에 의한 위상 성분을 추출하는 것이 정밀 DEM 생성을 위한 핵심이다.

3) 2-pass DInSAR 처리

전술한 바와 같이 생성된 interferogram에는 지형과 지표변위에 대한 위상을 모두 포함하고 있다. 2-pass DInSAR 방법에서는 준비된 DEM을 이용하여 지형정보 및 관측각도를 이용한 고도에 의한 간섭영상을 모사할 수 있으며, 이를 이용해 지형에 대한 위상을 제거할 수 있다. 간섭영상 모사를 위하여 준비한 지리좌표계의 외부 DEM을 SAR영상 좌표계로의 역 좌표계변환(Geocoding)이 필요하다. 이를 위해서 지리 좌표계의 표고자료와 SAR 영상 기하를 이용하여 각 픽셀에 해당하는 좌표를 생성하고 이를 매칭시키기 위한 수단인 geocoding lookup table을 이용한다. geocoding lookup table의 초기값 설정은 궤도 정보를 이용하여 계산하며 이때, 해상도 30 m의 SRTM DEM을 interferogmram의 픽셀 크기와 맞추기 위하여 SRTM DEM 기준 위도 및 경도 방향으로 각각 5 및 4의 비율로 오버샘플링(oversampling)하여 유사한 ground dimension(5 m×8 m)을 갖게 하였다. 또한, DEM과 SAR 영상의 해상도 차이로 인한 오차에 대해서는 registration offset을 계산/보정하여 geodocing lookup table을 개선(refine)한다. 생성한 geocoding lookup table을 활용하여 외부 DEM을 SAR영상 좌표계로 geocoding한 후, 위성 영상의 궤도 파라미터를 고려하여 고도 및 지구곡률에 의한 위상을 시뮬레이션한다. 모사한 위상을 원 interferogram에서 차분함으로써 차분간섭도(Differential interferogram)를 생성한다. 이렇게 생성된 Diff. interferogram에 남아 있는 위상 성분 Φ_diff.은

Φ_diff. = Φ_DEMerr. + Φ_def. + Φ_iono. + Φ_noise       (1)

이다. 여기서 Φ_DEMerr.는 Diff. interferogram 생성 시 고도에 대한 잔여위상, Φ_def.는 지표변위에 의한 위상, Φ_atm.는 대기에 의한 위상, Φ_iono.는 이온층 효과로 인한 위상, Φ_noise.는 노이즈 위상 성분이다. 일반적인 repeat pass SAR 영상쌍의 DInSAR 기반 DEM 생성 시에는 지표변위, 대기 및 이온층에 의한 위상성분의 고려가 필요하다. 그러나 본 연구에서는 bistatic 모드로 거의 동시에 획득된 TanDEM-X 영상쌍을 연구 자료로 활용하기 때문에 temporal decorrelation을 완벽히 제거할 수 있어 지표변위, 대기 및 이온층 효과로 인한 위상은 무시할 수 있다. 이로써 TanDEM-X Diff. interferogram에 남아있는 위상은 차분에 사용된 외부 DEM의 오차에 의한 성분만이 남아있게 된다.

4) 위상 불구속화

Diff. interferogram의 잔여 위상은 (-π, π] 범위로 구속화되어 있기 때문에 이에 대한 절대 위상값으로 변환하는 위상 불구속화(phase unwrapping)를 수행이 필요하다. MCF(Minimum Cost Flow)기법을 이용하여 phase unwrapping을 수행하였으며, 이때 긴밀도(coherence) 영상을 가중치로 사용하였다. phase unwrapping을 통해 생성한 unwrapped 영상은 DEM오차에 의한 위상성분을 나타낸다.

본 연구에서는 보다 정확한 unwrapped 영상(=DEM 오차로 인한 잔여 위상)을 획득하기 위하여 다음의 후처리과정을 제안한다. 우선 baseline 효과로 인하여 영상 전체에 오차로 포함되어 있는 trend를 제거하기 위하여 coherence가 0.3 이하인 지역을 마스킹한 후 2차원 다항식 모델(a0+a1_*y+a2_*x)을 계산하고 이를 unwrapped 영상에서 제거하였다. 추가적으로 unwrapped 영상에서 DEM오차의 특이 값(outlier)에 대한 보정을 수행하였다. DEM오차 dh(m)는

\(d h=d \Phi \times\left(\frac{d \Phi}{d h}\right)^{-1} \approx d \Phi \times \frac{\lambda \rho \sin \theta}{4 \pi B_{\perp}}\)       (2)

이다. 여기서 dΦ(rad)는 unwrapped DInSAR 영상의 DEM오차로 인한 잔여 위상, h(m)는 SRTM DEM의 고도, dΦ/dh(rad/m)는 고도에 따른 위상의 변화량이며, dΦ/dh는 B_⊥의 수직기선(Perpendicular baseline) 거리 및 SARP(SAR Processing) 및 궤도 상태 벡터(state vector)의 파라미터(λ의 SAR 마이크로파의 파장 길이, ρ는 위성체에서 지표까지의 관측거리, θ는 관측 입사각)를 기반으로 각 패스에 대한 SAR 위성체의 위치를 결정하여 지형과 관련된 위상을 시뮬레이션한다. 이렇게 생성한 DEM 오차 영상에서 (-10 m, 10 m] 이외의 값에 대하여 마스킹한 후, 2차원 다항식 모델에 fitting하여 trend를 제거하였다. 보정 후에 unwrapped 위상 영상에서 unwrapping 오차를 재차 확인하며, 문제가 있을 시 임계값을 조정해가며 제안하는 후처리 과정을 반복한다. 최종 보정된 unwrapped 위상 영상에는 보정을 통하여 마스킹된 부분이 존재하며 이를 내삽(interpolation)을 통해 빈값을 채워주었다.

5) DInSAR 기반 DEM 생성

보정된 unwrapped 영상은 DEM오차로 인한 잔여 위상을 나타내며, 전술한 바와 같이 unwrapped 영상과 고도에 따른 위상의 변화량을 이용하여 DEM오차 영상을 생성할 수 있다. 이렇게 생성된 DEM오차를 외부 DEM인 저해상도의 SRTM DEM에 더해주어 DInSAR 기반 고해상 DEM을 생성한다. 이러한 DInSAR 기반 DEM은 SAR 좌표계에서의 고도값을 나타내며, 이 값을 지리 좌표계로의 변환하여 최종 결과를 생성해야 한다.

본 연구에서는 보다 정확한 최종 DEM 생성을 위한 좌표계 변환 과정을 제안한다. 기존의 연구에서는 2-pass DInSAR 처리 시, 저해상도의 외부 DEM을 이용하여 생성한 geocoding lookup table을 이용하여 좌표계 변환을 수행하였다. SRTM DEM에 비하여 보정되는 DEM오차가 작은 경우에는 기존의 외부 DEM에서 생성한 geocoding lookup table을 이용한 좌표변환 방법은 문제가 되지 않는다. 그러나 본 연구의 연구지역인 화산과 같은 지역은 지표변화가 큰 지역으로 SRTM DEM에 비하여 보정되는 DEM오차가 크다(>50 m). 이러한 경우 기존의 SRTM DEM 기반 lookup table은 부정확하여 좌표계 변환 오차가 발생하여 정밀한 DEM 생성에 적합하지 못하다. 따라서 이러한 좌표변환 문제를 극복하기 위해서 앞서 생성한 DInSAR 기반 DEM을 활용하여 geocoding lookup table을 생성하고, 이 DEM과 SAR 영상의 offset을 계산/보정하여 geocoding lookup table을 개선하였다. 이렇게 생성한 고정밀 geocoding lookup table을 활용하여 지리좌표계 상의 DInSAR 기반 DEM을 생성한다. 이렇게 지리좌표계로 변환된 DInSAR 기반 DEM의 오차를 재차 확인하며, 문제가 있을 경우 생성한 DInSAR 기반 DEM을 활용하여 ‘3) 2-pass DInSAR 처리’ 부터 재수행한다. 새로운 단계로써 반복 수행을 통하여 재생성한 고도 오차는 직전 단계에서 생성하였던 DInSAR 기반 DEM에 더해주어 DInSAR 기반 DEM을 업데이트한다. 또한, 생성한 DInSAR 기반 DEM의 지리 좌표로의 좌표계 변환은 DEM오차에 따라 다음과 같이 수행한다. 1) DEM오차가 작다면(<50 m), 직전 단계에서 DInSAR 기반 DEM을 활용하여 생성하였던 geocoding lookup table을 활용, 2) 재수행 단계임에도 DEM오차가 크다면(>50 m), 해당 단계에서 생성한 DInSAR 기반 DEM을 활용하여 geocoding lookup table을 다시 생성하고 이를 이용하여 좌표계 변환을 수행한다.

지리좌표계로 변환된 DInSAR 기반 DEM에 문제가 없을 시, 영상 내 빈값에 대하여 interpolation하고 5×5 평균값필터(average filter)를 적용하여 노이즈를 저감하는 후처리 과정을 거쳐 최종 DInSAR 기반 DEM을 생성하였다.

4. DEM 생성 방법의 성능 검증

Fig. 4는 본 연구에서 생성한 DEM으로, (a)는 2-pass DInSAR 처리에 사용한 SRTM 30m DEM, (b)는 제안하는 unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하지 않은 2-pass DInSAR 기반 DEM, (c)는 제안하는 unwrapped 위상 오차 보정을 수행하고 좌표계변환 오차에 대한 보정을 수행하지 않은 2-pass DInSAR 기반 DEM, (d)는 제안하는 unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하여 생성한 2-pass DInSAR 기반 DEM을 나타낸다. 저해상의 SRTM DEM(공간해상도 30 m×30 m)에 DInSAR 기반 DEM오차(공간해상도 약 5 m×8 m)를 보상하여 고해상의 DEM을 획득하였다.

Fig. 4. Results of 2-pass DInSAR based-DEMs and SRTM DEM.

Fig. 5는 Fig. 4의 2-pass DInSAR 기반 DEM의 보다 명확한 육안 비교분석을 위하여 식(2)에서 생성한 보정된 DEM오차(SRTM DEM과의 차이)를 나타내는 것으로, (a)는 unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하지 않고 생성한 DEM오차, (b)는 제안하는 unwrapped 위상 오차 보정을 수행하고 좌표계변환 오차에 대한 보정을 수행하지 않고 생성한 DEM오차, (c)는 제안하는 unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하여 생성한 DEM오차, (d)는 Fig. 4(d)의 제안하는 방법으로 생성한 DEM을 기반으로 제안하는 2-pass DInSAR 기반 DEM 생성 방법 반복 수행하여 생성한 DEM오차를 나타낸다. Fig. 5(a)는 unwrapped 위상오차에 대하여 저감하지 않아 영상 전체적으로 대략 -50 m의 DEM오차가 분포한다. 또한, 산악지역에 따른 레이더신호 왜곡(layover, shadow, foreshortening) 지역 및 분화구 및 용암류(lava flow) 지역과 같이 형태학적 변화가 큰 지역에서 unwrapping 오차로 인하여 최대 -170 m의 DEM오차를 보인다. 이러한 unwrapping에 의한 오차는 위상 unwrapping 시, reference 처리를 어떻게 하느냐에 따라 그 값이 달라질 수 있다. Fig. 5(b)는 제안하는 unwrapped 위상 오차를 보정을 적용하여 영상 전체적인 DEM오차가 보정되었으며, 시간에 따른 고도 변화가 거의 없는 산의 저고도 평지에서 DEM오차가 0에 가깝게 보정되었다. 제안하는 unwrapped 위상 오차를 보정을 통하여 unwrapped 영상은 range 및 azimuth 방향으로 1,000 픽셀 당 각각 0.09 m 및 -2.2 m의 DEM오차 위상 성분의 offset이 보정되어 보다 정확한 DEM오차 획득이 가능함을 알 수 있다. 그러나 좌표계변환 과정에서 저해상도 SRTM DEM으로부터 계산한 geocoding lookup table을 사용하기 때문에 DEM과 SAR 해상도 차이에서 발생하는 오차 및 레이더신호 왜곡 지역에 대한 처리(interpolation) 오차 등이 포함된다. Fig. 5(c)는 unwrapped 위상 오차 및 좌표변환 오차를 보정하여 Fig. 5(a) 및 Fig. 5(b)에서 보인 오차들이 저감된 것을 확인할 수 있다. Fig. 5(d)는 Fig. 4(d)의 제안한 방법에 의해 생성한 DEM을 입력자료로 재활용하여 생성한 DEM오차를 나타내며, 제안한 DEM 생성방법의 반복 수행을 통하여 최대 약 6m의 DEM오차가 추가적으로 보정이 되었다.

Fig. 5. Results of 2-pass DInSAR based DEM errors (height corrections).

Table 2는 총 41점의 고도에 대한 GPS 측량성과(Bonforte and Puglisi, 2006)와 GPS 위치에 해당하는 생성한 DEM의 높이 값과 GPS 고도를 기준으로 하여 계산한 각 DEM의 수직정확도(RMSE)를 나타낸다. 또한, Fig. 6은 GPS에서 관측한 고도에 대한 생성한 각 DEM의 오차를 나타낸다. 2-pass DInSAR DEM 생성에서 reference로 사용한 SRTM 30 m DEM의 고도 오차(Fig. 6의 보라색 선)는 최대 16.665 m, 최소 -8.486 m이고 평균과 표준편차는 각각 2.010 m, 5.256 m 이며, SRTM DEM 의 수직정확도는 5.567 m이다. unwrapped 위상 및 좌표계변환 오차에 대한 보정을 수행하지 않고 생성한 DInSAR 기반 DEM의 고도 오차(Fig. 6의 하늘색 선)는 최대 -30.120 m, 최소 -8.486 m이고 평균과 표준편차는 각각 -41.372 m, 5.169 m 이며, 수직정확도는 39.617 m이다. 이러한 편향된 고도오차는 unwrapped 위상 및 좌표계변환 오차이다. 참고로 관측된 수직정확도(39.617 m)는 위상 unwrapping 시 reference 처리를 어떻게 하느냐에 따라 달라질 수 있는 값이다. 보다 정확한 DEM 생성을 위하여 이에 대한 보정이 필요하다. 이에 제안하는 unwrapped 위상 오차 보정을 수행하여 생성한 DInSAR 기반 DEM의 고도 오차(Fig. 6의 주황색 선)는 최대 5.264 m, 최소 -8.486 m이고 평균과 표준편차는 각각 1.384 m, 2.101 m 이며, 수직정확도는 2.495 m이다. 제안하는 unwrapped 위상 오차 보정을 통하여 DInSAR 기반 DEM의 수직정확도는 39.617 m에서 2.495 m로 향상되었다. 생성한 unwrapped 영상에 50 m 이상의 DEM오차가 포함되어 있는 것을 고려하여 추가적으로 좌표계변환 오차 보정을 통하여 생성한 최종 제안하는 DInSAR 기반 DEM의 고도 오차(Fig. 6의 빨간색 선)는 최대 4.185 m, 최소 -8.486 m이고 평균과 표준편차는 각각 0.659 m, 2.279 m 이며, 수직정확도는 2.346 m이다. 제안하는 unwrapped 위상 오차 및 좌표계변환 오차 보정을 통하여 DInSAR 기반 DEM의 수직정확도는 39.617 m에서 2.346 m로 향상되었으며, unwrapped 위상 오차 보정 후 추가적으로 좌표계변환 오차의 보정을 통하여 약 15 cm의 수직 정확도가 향상되었다. 제안한 방법으로 생성한 DInSAR 기반 DEM을 활용하여 DInSAR 처리 과정부터 제안한 DEM 생성 방법을 재수행하였고 이에 생성한 DEM의 고도오차(Fig. 6의 노란색 선)는 최대 4.440 m, 최소 -8.486 m이고 평균과 표준편차는 각각 0.563 m, 2.360 m 이며, 수직정확도는 2.399 m이다. 반복수행을 통하여 오차에 대한 과보정이 이루어져 수직정확도가 약 12 cm 저하된 것을 확인할 수 있다. Fig. 5(d)에서 확인할 수 있듯이 DEM오차가 미소량(최대 약 6 m) 포함되어 있을 시, 오차 저감을 위한 과정의 반복수행은 필요하지 않다.

Table 2. DEM vertical errors estimated from GPS surveys

Fig. 6. Result of vertical errors assessment (Generated DEM-GPS elevation).

제안하는 방법을 통하여 reference로 사용한 SRTM 30 m DEM(수직정확도 5.567 m)에 DInSAR로 관측한 SRTM DEM 오차를(Fig. 5(c)) 보상하여 최종적으로 수직정확도가 2.346 m인 DEM을 획득하였다. 다시 말해 제안하는 방법을 통하여 SRTM DEM에 비하여 공간해상도는 약 5배, 수직 정확도는 약 2.4배 향상된 DEM을 생성할 수 있다. SRTM DEM과 공간해상도를 일치하여 비교하기 위하여 DInSAR 기반 DEM을 다운샘플링하여 공간해상도가 30 m인 DEM을 생성하였으며, 이 DEM의 고도 오차(Fig. 6의 초록색 선)는 최대 9.807 m, 최소-8.486 m이고 평균과 표준편차는 각각 -0.098 m, 3.310 m이며, 수직정확도는 3.271 m이다. SRTM DEM과 공간해상도를 일치시키고 비교하여도 수직방향의 고도 정확도는 약 1.7배 향상된 것을 확인할 수 있다. 추가적으로 TanDEM-X 90 m DEM과의 수직 정확도 비교를 수행하였다. TanDEM-X 90 m DEM의 고도 오차(Fig. 6의 검은색 선)는 최대 21.079 m, 최소 -8.486 m이고 평균과 표준편차는 각각 -3.441 m, 9.183 m 이며, 수직정확도는 9.701 m이다. TanDEM-X 90 m DEM과 공간해상도를 일치하여 비교하기 위하여 DInSAR 기반 DEM을 다운샘플링하여 공간해상도가 90 m인 DEM을 생성하였으며, 이 DEM의 고도 오차(Fig. 6의 갈색 선)는 최대 12.423 m, 최소 -8.486 m이고 평균과 표준편차는 각각 -1.426 m, 6.007 m 이며, 수직정확도는 6.102 m이다. TanDEM-X 90 m DEM과 공간해상도를 일치시키고 수직방향의 고도 정확도를 비교하였을 때, 제안한 방법을 통해 생성한 DInSAR 기반 DEM의 수직정확도는 약 1.6배 높은 것을 확인할 수 있다. 추가적으로 SRTM 30 m DEM의 수직정확도가 TDX 90 m DEM의 수직정확도 보다 좋게 측정된 것은 공간해상도의 차이, GPS 관측 시기 및 화산 활동 등의 영향인 것으로 사료되어진다.

5. 결론

TanDEM-X bistatic SAR 영상의 2-pass DInSAR 기법을 이용한 DEM 생성 방법에 대한 성능 향상을 위한 개선된 DEM 생성 방법을 개발하였다. 개발된 DEM 생성 방법은 unwrapped 위상 내의 DEM 오차와 좌표계 변환 시 발생할 수 있는 DEM오차를 현저히 줄일 수 있는 방법이다. 실제 unwrapped 위상 및 좌표계 변환 오차에 대한 보정을 수행하지 않고 생성한 DInSAR 기반 DEM의 고도 오차는 수직정확도는 39.617 m이며, 제안한 방법으로 오차에 대한 보정을 통하여 생성한 DInSAR 기반 DEM의 수직정확도는 2.346 m로 제안한 방법을 통하여 생성한 DEM의 정확도가 향상됨을 확인할 수 있었다. 제안하는 2-pass DInSAR 기법을 통해 reference로 사용한 SRTM 30 m DEM(수직정확도 5.567 m)에 DInSAR로 관측한 SRTM DEM 오차를 보상하여 최종적으로 공간해상도는 약 5배, 수직 정확도는 약 2.4배 향상된 DEM을 생성할 수 있었다. 제안하는 방법으로 생성한 DInSAR 기반 DEM의 공간해상도를 SRTM 30 m DEM과 공간해상도를 일치시키고 성능을 비교 검증한 결과 수직방향의 고도 정확도가 약 1.7배 향상된 것을 확인할 수 있었다. 추가적으로 공간해상도를 90m로 다운샘플링하여 TanDEM-X 90 m DEM과의 수직 정확도를 비교하였으며, 제안한 방법을 통해 생성한 DInSAR 기반 DEM의 수직정확도가 TanDEM-X 90m DEM보다 약 1.6배 높은 것을 확인할 수 있었다. 따라서 제안하는 2-pass DInSAR 기반 DEM 생성 방법으로 보다 정확한 DEM 생성이 가능함을 확인하였다.

빈번한 형태학적 변화를 갖는 지역에 대한 DEM의 지속적인 업데이트를 위하여 본 연구에서 도출한 방법을 이용한다면 저비용으로 빠른 시간 내에 효과적으로 DEM을 갱신할 수 있을 것이다. 그러나 본 연구에서 제안한 방법을 일반적인 repeat pass SAR 영상쌍에 적용하기 위해서는 추가 연구가 필요하다. 추가적 연구과제로 일반적인 repeat pass SAR 영상쌍의 interferogram에서 지표변위, 대기, 이온층에 의한 위상성분에 대한 정확한 분리를 통한 지형고도에 의한 위상 성분을 추출 및 이를 이용한 정밀 DEM 생성을 하는 방법론을 정립하는 것이다.