Comparison of Co-registration Algorithms for TOPS SAR Image

TOPS 모드 SAR 자료의 정합기법 비교분석

Abstract

For TOPS InSAR processing, high-precision image co-registration is required. We propose an image co-registration method suitable for the TOPS mode by comparing the performance of cross correlation method, the geometric co-registration and the enhanced spectral diversity (ESD) matching algorithm based on the spectral diversity (SD) on the Sentinel-1 TOPS mode image. Using 23 pairs of interferometric pairs generated from 25 Sentinel-1 TOPS images, we applied the cross correlation (CC), geometric correction with only orbit information (GC1), geometric correction combined with iterative cross-correlation (GC2, GC3, GC4), and ESD iteration (ESD_GC, ESD_1, ESD_2). The mean of co-registration errors in azimuth direction by cross correlation and geometric matching are 0.0041 pixels and 0.0016 pixels, respectively. Although the ESD method shows the most accurate result with the error of less than 0.0005 pixels, the error of geometric co-registration is reduced to 0.001 pixels by repetition through additional cross correlation matching between the reference and resampled slave image. The ESD method is not applicable when the coherence of the burst overlap areas is low. Therefore, the geometric co-registration method through iterative processing is a suitable alternative for time series analysis using multiple SAR data or generating interferogram with long time intervals.

TOPS InSAR 처리를 위해서는 고정밀도의 영상정합이 요구된다. Sentinel-1 TOPS 모드영상에 교차상관 영상정합, 기하정합, Spectral Diversity 정합에 기반 한 Enhanced Spectral Diversity 정합 알고리즘 성능 비교분석을 통해 TOPS모드에 적합한 영상정합 방법을 제시하고자 한다. 25개의 Sentinel-1 TOPS 영상으로부터 생성된 23개의 간섭쌍을 이용하여, 교차상관정합(CC), 기하보정(GC1), 기하보정 후 교차상관정합(GC2, GC3, GC4) 그리고 ESD를 이용한 정합(ESD_GC, ESD_1, ESD_2) 총 8가지 방법을 적용하였다. 교차상관정합과 기하정합에 따른 azimuth 방향 정합오차를 평균한 결과는 각각 0.0041 화소, 0.0016 화소이다. 비록 ESD 방법은 azimuth 정합오차가 0.0005 화소 이하로 가장 정확한 결과를 보이지만, 기하정합 결과는 추가적인 교차상관정합을 통한 반복 과정을 통해 0.001 화소 정도로 오차가 감소하였다. ESD 방법은 burst 중첩지역의 긴밀도가 낮은 경우 적용이 불가능하다. 따라서 반복 적용을 통한 기하정합 방법은 다수의 SAR 자료를 이용한 시계열 분석 또는 긴 시간간격을 갖는 간섭도 생성을 위해서 적합한 대안이 된다.

1.서론

최근에는 1 m급의 고해상도 SAR 위성영상을 이용한 연구뿐만 아니라 ALOS-2 ScanSAR 모드, Sentinel-1 TOPS (Terrain Observation by Progressive Scans) 모드와 같이 해상도보다는 넓은 지역을 짧은 방문주기로 촬영하여 자연재해, 기후변화 등 광역적 변위를 관측하는 연구가 활발히 이루어지고 있다(Himenmatsu and Furuya, 2016; Yague-Martinez et al., 2016). TOPS모드는 azimuth 방향으로 안테나를 스캔하여 짧은 시간 안에 광역 범위의 고해상도 SAR영상을 얻을 수 있기 때문에 한국항공 우주연구원에서 개발 중인 다목적실용위성6호에도 적용될 예정이다(Lim, 2017).

SAR 영상의 공간해상도가 점점 향상됨에 따라 영상 분석에 요구되는 영상정합의 정밀도도 점점 높아지고 있다. SAR 영상을 이용한 간섭기법에서 영상정합 오차는 간섭도의 긴밀도를 낮출뿐만 아니라 TOPS 영상에 서 burst 간의 위상차를 발생시킨다. Sentinel-1 TOPS 자료의 InSAR에서 burst 간의 위상차를 ±1.5° 이내로 유지하기 위해 요구되는 영상정합의 정밀도는 약 1 cm(~0.001 pixels) 정도이다. 하지만, 기존 진폭영상 상관계수 및 fringe 선명도를 이용한 방법(Grabriel and Goldstein, 1988; Lin et al., 1992)의 정밀도는 일반적으로 0.05 pixel 정도이다(Li and Goldstein, 1990). TOPS InSAR 처리를 위해서는 고정밀도의 영상정합, 부영상의 보간과 azimuth 방향의 spectral filtering이 요구된다(Yague-Marinez et al., 2016). 본 연구에서는 정확도 높은 간섭영상생성을 위해 Sentinel-1 TOPS 모드영상에 교차상관 영상정합, 기하정합, Spectral Diversity 정합에 기반한 Enhanced Spectral Diversity 정합 알고리즘을 적용해보고자 한다. 이를 통해 TOPS모드에 적합한 고정밀도 영상정합 방법을 제시하고자 한다.

2. TOPS 모드 자료를 이용한 영상정합

1) TOPS 모드 자료 특성

TOPS 모드는 전통적인 ScanSAR 모드와 달리 안테나 빔 방향을 ψ각도로 제어하면서 촬영한다(Fig. 1). ScanSAR모드와 같이 광역 관측을 수행하지만, 기존의 ScanSAR모드의 스캘러핑(scalloping) 현상과 azimuth에 따라 달라지는 signal-to-ambiguity 비를 안테나 빔 방향을 제어하는 것으로 극복하였다(Prats-Iraola et al., 2012). 하지만, TOPS 모드에서는 안테나 빔 방향의 변화로 인해 azimuth 방향으로 도플러 중심주파수가 변화한다. 또한, 도플러값 변화가 azimuth 방향 샘플링 주파수보다 높은 값을 가지고 있고 burst 모드 관측되기 때문에 각 burst 자료 간의 동기화가 요구된다(Rodriguez-Cassola et al., 2015).

Fig. 1. Sketch of TOPS acquisition geometry. T_B is the burst duration, rmid is the midrange slant-range distance in a given subswath, v_s is the satellite velocity, hsat is the satellite altitude, k_θ is the steering angle, and ψ_s and ψ_e are the staring and the end steering angles with the burst, respectively (Meta et al., 2010).

TOPS 영상은 Fig. 1과 같이 총 3개의 sub-swaths를 가지며 각 swaths내에 burst로 나뉜다. Swaths 내의 각각의 burst들은 SLC 영상 자료처리 시 각 burst들의 도플러 값을 이용하여 모자이크 처리가 가능하다. 따라서 TOPS InSAR 처리를 위해서는 burst 추출이 필요하며 각 burst 들의 도플러 값을 이용한 개별적 자료처리 후 모자이크 과정이 필요하다. Fig. 2(a)는 TOPS 자료를 이용한 InSAR 자료처리 시 영상 정합 단계에서 발생하는 오차로 인해 burst 간에 위상차가 발생하는 현상을 보여주고 있다. Azimuth 방향으로 정합오차는 0.05화소에 불과 하지만 burst 중첩지역에서 약 145°의 위상차가 발생하였다. ESD 알고리즘을 이용한 정밀 정합 결과(Fig. 2(b))에서는 burst 간 위상차가 거의 대부분 제거되었음을 확인할 수 있다.

Fig. 2. TerraSAR-X TOPS interferograms (a) with a coregistration error of 0.05 samples (b) after precise coregistration using ESD (De Zan et al., 2014).

2) 교차상관 영상정합

두 SAR 영상간의 정합과정과정은 크게 진폭상관계수를 이용하는 방법(Kwoh et al., 1994)과 fringe의 명확도를 이용하는 방법(Gabriel and Goldstein, 1988)이 사용된다. 진폭상관계수가 일반적으로 사용되는 방법이나, 경사가 완만한 지역에서 진폭영상에서의 밝기값의 변화가 거의 없는 경우 fringe 명확도를 이용하는 방법이 효과적이다. SAR 영상의 진폭상관계수를 이용한 영상정합은 incoherent cross-correlation방법이며,일반적으로coherent cross-correlation의 최대 정밀도는 Cramér-Rao boun에 의해 주어지는데, 오차의 표편차는 식 (1)과 같다(Bamler and Eineder, 2015).

\(\sigma_{\hat{\Delta x}}=\sqrt{\frac{3}{2 N}} \frac{\sqrt{1-\gamma^{2}}}{\pi \gamma}\)       (1)

여기서,N은 영상정합 추정을 위해 사용된 window의 크기이며,γ는 긴밀도이고 오차는 공간해상도 단위의 측정값이다.

3) 기하정합

위성의 관측 궤도에 따른 주영상과 부영상 SAR 자료간의 오프셋은 각 화소의 고도에 따라 차이가 발생한다. TerraSAR-X와 같은 고해상도 SAR 영상은 기존의 JERS-1, ERS-1/2, ENVISAT과 달리 두 SAR 영상간의 coregistration을 위해서는 지형에 의한 영향을 고려하여야 한다. ERS의 경우 단순한 affine 또는 고차다항식을 이용한 영상정합 모델링을 수행하여도 지형고도에 의한 coregistrtion 오차가 약 0.1 pixel 이내 이기 때문에 부영상을 재배열하는 과정에 크게 문제가 없으나, Fig. 3과 같이 TerraSAR-X coregistration의 경우 다항식 모델링에 의한 영상정합 모델과 기하영향을 고려한 정밀 정합 결과와의 차이가 1 화소에 이른다. Sansosti et al.(2006)는 ERS-1/2 SAR 자료중 수직 기선거리가 매우 긴(약 1 km) 간섭쌍의 경우 영상정합 과정에서 지형고도에 따른 왜곡 현상을 고려하지 않았을 때 긴밀도가 매우 떨어지는 현상을 보완하기 위한 기하정합 알고리즘을 제안하였다. 본 연구에서는 위성궤도 정보와 DEM에 기반하여 계산된 영상 오프셋을 반영하는 기하영상정합 과정을 적용하였다. 영상정보를 사용하지 않는 순수 기하영상정합에 요구되는 위성 궤도의 정밀도는 수 cm~수십 cm 정도이지만, 기하영상정합은 일반적으로 진폭영상의 상관계수를 이용한 정합과 병행하여 적용되며 이 경우 궤도의 정밀도가 다소 낮아도 충분히 적용이 가능하다. 또한 기하정합에 사용되는 DEM은 SRTM 정도의 정밀도로도 충분히 적용 가능하다(Sansosti et al., 2006).

Fig. 3. (a) Topography induced sample shift with the optimal polynomial modelling of the shifts (green line), (b) error caused by polynomial modelling (Eineder and Adam, 2007).

4) Spectral Diversity 정합

TOPS 자료에서 burst간에 중첩지역이 발생하지 않는 경우에 정밀정합을 수행하기 위해서는 SD (Spectral Diversity) 기반 알고리즘을 적용하여한다(Scheiber and Moreira, 2000). SD를 위해서 우선, 각 burst별 두 개의 full resolution SAR 영상 (Master와 Slave)으로 부터 각 각 forward 및 bakward 방향의 sub-look 영상을 Fig. 4와 같이 생성한다.

Fig. 4. Generation of sub-look images in forward and backward direction.

두 개의 sub-look으로 분해된 M_f, M_b, S_f, S_b 영상으로부터 forward 영상간의 간섭도 및 backward 영상간의 간 섭도를 생성하고, 이를 차분하여 식 (2)와 같이 ø_sd 간섭 영상을 생성한다.

\(\mathscr{\varnothing}_{s d}=\left(M_{f} \cdot S_{f}^{*}\right) \cdot\left(M_{b} \cdot S_{b}^{*}\right)^{*}\)       (2)

ø_sd는 식 (3)와 같이 두 sub-look 영상의 도플러 중심주파수 차이(Δf)와, azimuth 방향의 영상정합 오차(sec) Δt로 표현되므로 측정된 ø_sd를 이용하여 영상정합 오차를 구할 수 있다(Prats-Iraola et al., 2012).

\(\emptyset_{s d}=2 \pi \Delta f \Delta t\)       (3)

SD 알고리즘을 적용 시, Fig. 5와 같이 sub-look으로 분해된 영상으로 생성된 간섭도는 노이즈가 심하기 때문에 multilook을 통해 먼저 잡음을 감소시키고 차분해야 한다.

Fig. 5. Differential interferomgram (a) based on the 1-look forward/backward interferogram, and (b) based on the 10×2 look forward/backward interferogram.

5) Enhanced Spectral Diversity 정합

Enhanced Spectral Diversity(ESD) 알고리즘은 Scheiber and Moreira(2000)이 제안한 spectral diversity(SD)에 기반한 알고리즘으로,sub-look영상을 생성하지 않고 azimuth 방향의 미세한 영상정합 오차를 보정하기위해 Fig. 6과 같이 burst 간 중첩 지역의 위상을 탐색하여 영상정합 오차를 보정하는 방법이다. ESD 알고리즘은 계산이 간편하고 긴밀도가 낮은 지역에서도 안정적인 값을 추정 할 수 있다는 장점이 있다.

Fig. 6. Busted overlap of TOPS data.

주영상(master)과 부영상(slave)간에 생성된 i번째 burst 와 i+1번째 urst 간의 중첩지역간의 위상차는 식 (4)와 같다(Yague-Martinez et al., 2016).

\(\emptyset_{E S D}=\arg \left(m_{i} \cdot s_{i}^{*}\right) \cdot\left(m_{i+1} \cdot s_{i+1}^{*}\right)^{*}\)       (4)

화소 p에서 측정된 ø_ESD와 azimuth 방향 오프셋(Δy)간의 관계식은 식 (5)와 같다.

\(\emptyset_{E S D, p}=2 \pi \Delta f_{D C, p}^{o v l} \frac{\Delta y}{f_{a z}}\)       (5)

여기서, \(\Delta f_{D C, p}^{o v l}\)은 burst 중첩지역에서의 도플러중심주 파수 차이이며, f_az는 azimuth 방향 sampling 주파수이다. 따라서 azimuth 방향 오프셋(Δy)의 추정값은 식 (6)과 같이 정리될 수 있다.

\(\widehat{\Delta y}=\frac{f_{a z}}{2 \pi} \cdot \frac{\arg \left\{\left\langle e^{j \emptyset E S D, p\rangle}\right\}\right.}{\left\langle\Delta f_{D C, p}^{o v l}\right\rangle}\)       (6)

3. Sentinel-1 TOPS 자료를 이용한 실험 결과 및 분석

1) 비교분석에 사용된 자료

TOPS 모드에 적합한 고정밀도 영상정합 방법을 알아보기 위해 본 연구에서 분석한 자료는 2017년 3월 15일부터 2018년 8월 1일 동안 한반도 남동부 지역을 촬영한 Sentinel-1 TOPS 영상을 사용하였다. 두 번째 subswaths에서 4-6까지의 3개의 burst를 추출하여 분석에사용하였다(Fig. 7). 총 25개의 영상을 사용하였으며, 간섭쌍은 연속해서 촬영된 자료들 중 시간간격이 짧은 23개의 간섭쌍를 제작 분석하였다. 간섭쌍의 수직기선 거리는 최대 200 m 이내로 대부분 수십 m 이내의 비교적 짧은 기선거리를 가지고 있다(Fig. 8).

Fig. 7. Sentinel-1 TOPS preview image and processing area for analysis.

Fig. 8. Interferometric pair of Sentinel-1 TOPS data.

2) 자료처리 및 영상정합 성능 비교

Sentinel-1 영상으로부터 추출된 3개의 burst 영상을 이용하여 두 SAR 영상정합에 교차상관 영상정합, 기하 영상정합 및 ESD 영상정합 알고리즘을 각각 적용하여 간섭도 영상을 생성하였다. Fig. 9는 추출된 연구지역의 Sentinel-1 영상과 3개의 burst간의 중첩지역을 보여주고 있다. 또한 burst 중첩 지역에서의 위상차를 계산하기 위해 긴밀도를 기준으로 선택된 masking 지역의 예를 도시한 결과이다

Fig. 9. (a) Amplitude image, () burst overlap and (c) mask for phase shift estimation.

분석을 위해 생성된 23개의 간섭도 중 2017년 4월 8일 자료와 2017년 4월 20일 자료를 이용하여 생성된 간섭도는 Fig. 10과 같다. 전통적인 교차상관 영상정합 기법을 이용하여 생성된 간섭도(Fig. 10(a))에서 3개의 burst의 중첩 부분에서 간섭위상 오차가 불연속적임을 확인할 수 있다. 반복적인 기하정합을 수행한 결과를 생성된 간섭도(Fig. 10(b))에서는 burst 간의 경계에서 위상차가 거의 나타나지 않는 것으로 보인다. 교차상관정합, 기하정합, ESD 정합으로 생성된 차분간섭도 간의 비교를 위해, ESD 정합 간섭도를 기준으로 교차상관정합과 기하정합 간섭도를 차분한 결과는 Fig. 11과 같다. Fig. 10(a)에서 관측되었던 위상 불연속 현상이 ESD 간섭도와 차분한 결과(Fig. 11(a))에서 뚜렷하게 관측되는 반면, ESD 간섭도와 기하정합 간섭도를 차분한 Fig. 11(b)에서는 거의 잔여위상이 발생하지 않는 것으로 확인된다.

Fig. 10. Interferogram produced by (a) cross-correlation co-registration and (b) geometric co-registration.

Fig. 11. Difference phase image of (a) between ESD interferogram and cross-correlation interferogram, (b) between ESD interferogra and geometric interferogram.

기하정합을 위해 SRTM-1arc DEM을 이용하여 레이더좌표계 고도값을 생성하였다. 비록 Sentinel-1 위성자료의 궤도가 정밀하지만, 위성진행방향과 거리방향으로의 offset이 존재할 수 있다. 따라서 기하정합을 수행 하더라도 최소 1회 이상의 교차상관정합을 추가적으로 수행하는 것이 일반적이다. 본 연구에서는 순수하게 기하정합만을 실행한 결과와 1회에서 3회까지 반복과정을 통해 교차상관정합을 수행하여 영상정합의 정밀도를 갱신하였다. 또한 ESD 정합 알고리즘 역시 초기에는 기하정합을 수행하고, 이후 burst 중첩지역에서의 위 상차를 이용하여 azimuth 방향으로의 offset을 보정하는 과정을 수행하였으며, 두 번의 반복과정을 수행하였다. 따라서, 영상정합 방법에 따른 정밀도 분석을 위한 실험에서, 교차상관정합(CC), 기하보정(GC1), 기하보정 후 교차상관정합 횟수에 따라 GC2, GC3, GC4, 그리고 ESD에서 사용된 초기 기하보정(ESD_GC), 추가 ESD 보정 1회(ESD_1), 추가 ESD 보정 2회(ESD_2) 총 8가지 방법을 적용하였다. 각 정합방법의 정확도 분석을 위해 ESD_2의 정합 결과를 가장 우수한 것으로 가정할 수 있기 때문에 ESD_2를 기준으로 각 정합방법의 오차를 추정하였다. 교차상관정합 결과와 ESD_2 결과로 생성된 간섭도의 차분을 수행하고, burst 중첩지역에서의 위상 평균을 구한 후 식 (6)에 따라 azimuth 방향 화소 오차를 구하였다. 기하정합 정확도는 최종 기하정합 결과 (GC_4)의 간섭도와 ESD_2 최종 간섭도간의 burst중첩 지역 위상차를 이용하여 구한 후, 나머지 기하정합 실험 (GC_1, GC_2, GC_3)의 기하정합을 통해 보정된 azimuth 방향 offset을 이용하여 구하였다. 먼저 생성된 23개의 간섭도에서 교차상관정합(CC)와 기하정합(GC_4)로 생성된 간섭도를 ESD_2 간섭도와 비교하여 구한 azimuth 방향 평균 정합 오차는 Fig. 12와 같다. 예상하는 바와 같이 교차상관정합은 대부분 0.002 화소 이상의 오차를 보이고 있으며 일부 0.01 화소를 초과하는 간섭쌍도 존재 한다. 그에 반해 기하정합 결과는 대부분 0.002 이하의 정합오차를 보이고 있다. 23개 간섭도에서 구한 azimuth 방향 정합오차를 평균한 결과는 교차상관정합 0.0041 화소, 기하정합 0.0016 화소이다.

Fig.12. Azimuth co-registration error by CC and GC co-registration method.

ESD 정합방법을 이용해 생성된 최종 간섭도(ESD_2)를 기준으로 나머지 7개의 영상정합 방법에 의한 azimuth 방향 정합 오차를 도시한 결과는 Table 1, Fig. 13과 같다. 초기 기하정합 결과(GC_1)은 앞서 설명한 바와 같이 위성궤도 오차에 의한 offset이 재하기 때문에 다른 방법에 비해 높은 오차를 보이지만 추가적인 교차상관정 합을 통해 급속히 오차가 줄어들며, 4번 이상 반복을 통 해 0.001 화소 정도의 오차로 감소하는 것을 볼 수 있다. ESD 와 기하정합을 수행한 결과는 평균 0.002 화소 이하의 오차를 보이며 1회 burst 지역 위상값을 이용한 보정을 통해 0.0005 화소 이하의 정확도를 보임을 확인할 수 있다.

Table 1. Statistics of azimuth co-registration error by CC, GC, and ESD with respect to ESD_2 method

Fig. 13. Plot of mean and standard deviation of azimuth co-registration error by CC, GC, and ESD with respect to ESD_2 method.

비록 ESD 적용을 통해 매우 낮은 영상정합오차를 구현할 수 있지만, 실질적인 자료처리를 위해서는 burst 중첩지역에서 매우 높은 긴밀도가 요구된다. 따라서, 간섭쌍에 사용된 두 SAR 자료간의 시간 간격이 길 경우, ESD 방법은 적용이 불가능하다. 반면 기하정합과 교차 상관정합을 결합한방법은 SAR 영상의 진폭영상만을 이용하기 때문에 긴밀도에 영향을 거의 받지 않는다. 따라서, 다수의 SAR 자료를 이용한 시계열 분석 또는 시간간격을 갖는 간섭도 생성을 위해서는 기하정합 방법이 적합한 대안이 될 수 있다.

4.결론

TOPS 모드 SAR 자료의 InSAR 자료처리를 위한 영상정합 알고리즘 비교분석을 수행하였다. 24개의 Sentinel-1 TOPS 영상으로부터 생성된 23개의 간섭쌍을 이용하여, 교차상관정합(CC), 기하보정(GC1), 기하보정 후 교차상관정합(GC2, GC3, GC4) 그리고 ESD를 이용한 정합(ESD_GC, ESD_1, ESD_2) 총 8가지 방법을 적용하였다. 교차상관정합은 대부분 0.002 화소 이상의 오차를 보이고 있으며, 기하정합 결과는 대부분 0.002 이하의 정합오차를 보였다. Azimuth 방향 정합오차를 평균한 결과는 교차상관정합은 0.0041화소, 기하정합은 0.0016 화소이다. 초기 기하정합 결과(GC_1)는 위성궤도 오차로 인해 높은 오차를 보이지만 추가적인 교차상 관정합을 통해 급속히 오차가 줄어들며, 4번 이상 반복을 통해 0.001 화소 정도로 오차가 감소하였다. ESD 방법은 burst 중첩지역 위상값 보정을 통해 0.0005 화소 이하의 정확도를 보임을 확인할 수 있다. ESD 방법은 가장 우수한 정합결과를 제공하지만 burst 중첩지역의 긴밀도가 낮은 경우 적용이 불가능하다. 기하정합 방법의 경우 반복적인 기하정합 적용을 통해 약 0.001 화소 정도의 정확도 구현이 가능하다. 따라서 기하정합 방법의 정확도가 ESD 방법보다 낮지만, 다수의 SAR 자료를 이용한 시계열 분 또는 긴 시간간격을 갖는 간섭도 생성을 위해서는 기하정합 방법이 적합한 대안이 될 수 있다.