High-quality Stitching Method of 3D Multiple Dental CT Images

3차원 다중 치과 CT 영상의 고화질 스티칭 기법

Abstract

In this paper, we propose a high-quality stitching method of 3D multiple dental CT images. First, a weighted function is generated using the difference of two distance functions that calculate a distance from the nearest edge of an overlapped region to each position. And a blending ratio propagation function for two gradient vectors is parameterized by the difference and magnitude of gradient vectors that is also applied by the weighted function. When the blending ratio is propagated, an improved region growing scheme is proposed to decide the next position and calculate the blending intensity. The proposed method produces a high-quality stitching image. Our method removes the seam artifact caused by the mean intensity difference between images and vignetting effect. And it removes double edges caused by local misalignment. Experimental results showed that the proposed method produced high-quality stitching images for ten patients. Our stitching method could be usefully applied into the stitching of 3D or 2D multiple images.

1. 서 론

영상 스티칭 기법은 중첩되는 영역을 가진 여러 장의 영상들을 정합한 후, 하나의 영상으로 합성하는 기법이다[1]. 영상 스티칭 기법은 다양한 분야에 활용되고 있다. 특히 의학 영상 분야에서는 치과 CT(computed tomography) 장비는 기계적으로 화각의 한계가 있어서 치아 전체를 한 번에 촬영할 수가 없기 때문에 여러 번 영역을 나누어서 촬영하게 된다. 이 경우에 촬영된 다수의 3차원 영상들을 정합한 후 하나의 3차원 영상으로 합성해야 하는데, 영상 간의 평균 밝기값 차이나 영상 모서리나 주변부가 어두워지는 비네팅 효과(vignetting effect) 때문에 영상의 중첩 부분에서 심 결함(seam artifact)이 생기게 된다. 또한 지역적 정합 오류가 있는 경우에 이중 에지(double edge)가 생기게 된다.

기존에 제안된 스티칭 기법은 다음과 같다. 일반적으로 우선 영상들을 정합하기 위해 밝기값 기반 정합(intensity-based registration)이나 특징 기반 정합(feature-based registration) 방법을 이용한다[2]. 다음으로 중첩되는 영역에서 어떤 픽셀을 최종 합성 영상에 사용할지 선택하는 방법과 가시적인 심(visible seam), 흐려짐(blurring), 고스팅(ghosting) 등을 최소화하기 위해 픽셀을 최적으로 블렌딩(blending)하는 방법을 이용하여 영상을 합성하게 된다[1]. Uyttendaele 등은 페더링(feathering) 기법을 제안하였다[3]. 페더링 기법은 설정한 윈도우(window)의 경계선으로부터의 거리에 따라서 두 영상 비율에 가중치를 준다. 윈도우의 크기를 잘 설정하면 좋은 스티칭 결과를 얻을 수 있지만 결과 영상의 품질이 윈도우 크기에 영향을 많이 받고, 노출 차이(exposure difference)로 인하여 가시적인 심이 나타날 수 있다. 또한 정합 오류로 인해 고스팅이나 이중 에지가 발생하기도 쉽다. Adelson 등은 피라미드 블렌딩(pyramid blending) 기법을 제안하였다[4]. 이 기법은 각기 다른 주파수 대역(frequency bands)들을 그에 알맞은 알파 마스크(alpha masks)를 이용하여 블렌딩한다. 이 기법도 정합 오류가 있는 경우 이중 에지가 발생하게 된다. Davis는 중첩되는 영역에서 어떤 픽셀을 최종 합성 영상에 사용할지를 선택하기 위하여 최적 심 기법을 제안하였다[5]. 이 기법은 영상 간의 차이를 최소화하는 심 커브를 찾아서, 이를 기준으로 각 이미지를 대응되는 면에 복사한다. 이 최적의 심을 찾기 위해 그래프 컷 최적화(graph cut optimization)를 이용하는 방법이 연구되기도 하였다[6]. 이 방법들 역시 노출 차이로 인해 가시적인 심이 나타나기도 하며, 사용자의 입력을 요구하기도 한다. Zomet 등은 기울기 영역에서 영상 스티칭의 품질을 평가할 수 있는 비용 함수(cost function)를 제안하였고, 이 함수를 최적화시킴으로써 노출 차이와 정합 오류를 극복할 수 있다는 것을 보여주었다[7]. 이 알고리즘은 비교적 좋은 스티칭 결과를 보여주었지만 수행 속도가 느리고 함수 값이 항상 최적으로 수렴하지 않는다는 단점이 있는데, 쿼드 트리를 사용함으로써 수행 속도를 개선하려는 연구도 진행되었다[8]. Jia 등은 최적의 심을 계산한 후, 심 위에 존재하는 두 영상의 1차원 특징점들을 정합시키고 변형 벡터(deformation vector)를 구해서 이를 전파(propagation)시키는 방법을 제안하였다. 이 기법은 정합 오류를 정정하는데 좋은 성능을 보여주었다[9]. 김 등은 전립선 병리조직 영상들을 기하 보정 및 강체 정합 기법으로 자동 스티칭하는 방법을 제안하였다[10]. 소벨 연산자로 에지를 검출하고, 해리스 코너 검출자로 코너 점을 검출하여 에지와 코너 특징을 이용하여 정합을 수행하였다. 안 등은 영상에서 특징점을 추출하고, 이를 정합하는 과정에서 의미 있는 특징점을 분류하고, 이를 사용하는 고속 영상 스티칭 기법을 제안하였다[11]. 하지만, 김 등의 기법[10]과 안 등의 기법[11]은 특징 추출 오차에 정합 정확성이 크게 영향을 받는 문제점이 있다.

Fig. 1.Proposed method.

본 논문에서는 3차원 치과 CT 영상의 스티칭에서 발생할 수 있는 기존 기법들의 문제점들을 해결하며 빠르게 수행되는 고화질 영상 스티칭 알고리즘을 제안한다. 먼저 두 영상을 정합하기 위하여 밝기값 기반 방식을 사용한다. 유사도 척도(similarity measure)로는 정규상호상관(normalized cross correlation, NCC) 계수를 사용한다. 다음으로 물체 영역과 외곽선을 배경으로부터 추출한 후, 중첩되는 영역에서 각 영상의 외곽선을 기준으로 거리 맵(distance map)을 생성한다. 각 위치에서 두 거리 맵의 차이 값을 이용하여 가중치 함수를 구성한 후, 두 영상의 기울기 벡터(gradient vector)의 차이와 크기에 앞서 계산한 가중치를 적용함으로써 기울기 벡터를 블렌딩할 비율을 계산한다. 그리고 이 비율과 현재 위치에서의 밝기 값을 개선된 시드 기반 영역 확장(seeded region growing, SRG) 방식을 사용하여 3차원으로 전파시키며 다음 위치의 밝기 값을 계산한다. 이로써 두 개의 3차원 영상을 스티칭할 때 영상 간의 평균밝기 값 차이나 비네팅 효과로 발생하는 심 결함이나, 지역적 정합 오류로 인한 이중 에지를 제거할 수 있게 된다.

본 논문은 다음과 같은 순서로 구성되어 있다. 2장에서는 본 논문에서 제안한 3차원 다중 치과 CT 영상의 고화질 스티칭 기법에 대하여 설명한다. 3장에서 는 실험 결과를 기술하고, 4장에서는 결론을 맺는다.

 

2. 치과 CT 영상의 스티칭 기법

본 논문에서 제안한 치과 CT 영상의 스티칭 기법은 크게 네 가지 단계로 구성된다. 먼저 NCC를 이용하여 두 CT 영상을 정합하고, 다음으로 객체 영역을 분할하고, 객체 경계를 추출한다. 다음으로 객체 경계에 대하여 거리 맵을 생성하고, 마지막으로 심 결함과 이중 에지를 제거하기 위하여 블렌딩 비율 전파함수를 이용하여 두 영상을 스티칭한다.

2.1 NCC 기반 3차원 영상 정합

치과 CT 영상에서 해부학적 구조물을 특징으로 추출하기 위해서는 반자동 분할 기법의 도움이 있더라도 많은 시간이 소요되기 때문에 본 논문에서는 특징 추출이 필요가 없는 밝기값 기반 방식을 사용하였다[2]. 이 방식은 두 영상 간의 밝기값에 기반을 둔 유사성 척도를 최대화시키는 변환 인자를 찾는다. 유사성 척도로는 단일 모달리티의 영상 정합에 적합한 식 (1)의 정규상호상관(normalized cross correlation, NCC) 계수를 사용하였다.

이 식에서 n은 복셀의 개수이고, i1(x, y, z)는 첫 번째 CT 영상의 복셀 밝기값, i2(x, y, z)는 두 번째 CT 영상의 복셀 밝기값을 나타내고, 는 첫 번째 CT 영상의 복셀 밝기값의 평균, 는 두 번째 CT 영상의 복셀 밝기값의 평균을 나타낸다. 또한, σi1과 σi2는 가각 첫 번째와 두 번째 CT 영상의 복셀의 밝기값의 표준편차를 나타낸다.

치과 CT 영상 내의 객체는 강체 구조물로 가정할 수 있고, 다만 두 CT 영상의 촬영 인자의 차이에 의한 객체의 확대/축소는 존재할 수 있다. 따라서 본 논문에서는 변환 인자는 x, y, z축 방향의 이동 변환, 회전 변환, 확대/축소 변환의 인자를 사용하였다. 최적화 기법으로 브렌트 방법(Brent's method)을 사용하였고, 다중 해상도 기법을 적용함으로써 정합 속도를 향상시켰다. 정합 후 단순하게 두 영상을 스티칭하면 영상 간의 밝기 값 차이나 비네팅 효과 때문에 영상의 중첩 부분에서 심 결함이, 지역적 정합 오류로 인해 이중 에지가 나타나게 된다(Fig. 2).

Fig. 2.Cases of problems of previous methods. (a) Average intensity differences (b) Registration error (c) Vignetting effect

2.2 객체 경계 추출 및 거리 맵 생성

치과 CT 영상에서 공기로 구성된 배경을 분할하여 제거하면, 남는 영역이 객체가 된다. 따라서 본 논문에서는 CT 영상에서 각 슬라이스의 최외곽 모서리 네 점을 씨앗점으로 설정하고, 3차원 SRG 기법을 사용하여 배경을 분할하여 제거하여 객체를 분할한다[12]. 객체의 경계는 canny 경계 검출자를 이용하여 추출한다. 두 CT 영상에서 객체 영역과 경계를 추출한 후 객체 영역 내에 각 복셀에서 가장 가까운 경계까지의 거리를 나타내는 거리맵을 각각 생성한다. 이렇게 계산된 거리 값은 두 영상의 스티칭 과정에서 두 영상 내의 각 복셀 들의 기울기 벡터를 블렌딩하기 위한 비율을 결정할 때 가중치로 활용하게 된다.

본 논문에서 거리 맵 생성은 챔퍼 거리 변환[12]을 이용하여 근사한다. 챔퍼 거리 맵은 전체 영상을 인접 거리 정보를 갖고 있는 전방향 매스크(forward mask)와 역방향 매스크(backward mask)로 각각 영상을 순회하면서 생성한다. 제안 기법의 거리맵은 3차원에서 생성되지만, 2차원의 예를 들면, 다음과 같이 생성된다. 전방향 스캔에서는 위에서 아래 방향과 왼쪽에서 오른쪽 방향으로 영상 내의 모든 점 p에 대하여 f1(p)를 다음과 같이 계산한다.

이 식에서 B(p)는 점 p에 바로 인접한 점들 중 스캔 방향에 반대 방향에 위치한 (x, y+1), (x-1, y), (x-1, y+1), (x+1, y+1)로 구성된 점들의 집합을 나타낸다. 역방향 스캔에서는 아래 방향에서 위 방향과 오른쪽에서 왼쪽 방향으로 모든 점 p에 대하여 f2(p)를 다음과 같이 계산한다.

이 식에서 A(p)는 점 p에 바로 인접한 점들 중 스캔 방향에 반대 방향에 위치한 (x, y-1), (x+1, y), (x+1, y-1), (x-1, y-1)로 구성된 점들의 집합을 나타낸다.

2.3 스티칭 기법

정합된 두 영상을 단순 합성하게 되면, Fig. 2와 같이 심 결함이나 이중 에지, 흐려지는 현상 등이 발생하게 된다. 이런 결함을 제거하기 위하여 본 논문에서는 두 영상 각각의 기울기 벡터와 거리 맵을 이용하여 블렌딩 비율을 계산하고 이를 전파하는 함수를 제안하였다. 정합 후 두 영상의 중첩되는 영역에서 두 영상의 기울기 벡터의 차이와 크기에 거리 맵을 이용하여 가중치를 부여한 값으로 블렌딩할 비율을 계산한다. 이 비율은 개선된 영역 확장 방식으로 주변으로 전파되어 가면서 주변의 블렌딩 비율에 영향을 주게 된다. 이 비율에 따라 두 영상의 기울기 벡터를 적절하게 블렌딩한 후 이를 이용하여 영상의 밝기 값을 계산하게 된다. 아래 식 (4)는 두 영상의 기울기 벡터를 섞는 비율 P(q)를 구하는 식이다.

중첩 영역에서 P(q)는 1로 초기화되어 있으며, 이것은 두 영상 I1(q)와 I2(q)의 기울기 벡터를 1:0 비율로 블렌딩한다는 의미이다. q는 현재 위치이며 q′은 확장이 진행될 방향의 위치이다. P(q′) 값은 개선된 영역 확장 기법의 확장 진행 방향에 따라 같이 계산되어 나가기 때문에 최대 4번 씩 주변의 블렌딩 비율값으로부터 영향을 받게 되며, 이미 방문되어 값이 저장되어 있는 경우 식 (4)와 같이 평균을 내어 저장한다. 이로써 블렌딩 비율이 변하는 위치를 어느 정도 서로 모을 수 있기 때문에 원 영상의 특징을 유지하는데 도움이 된다.

또한, 두 영상의 기울기 벡터의 차이 값 ｜∇IA(q)-∇IB(q)｜ 이 0에 가까울수록 즉, 차이가 적을수록 두 영상이 부드럽게 연결될 수 있는 적절한 위치라고 가정할 수 있다. 따라서 I2(q)의 기울기 벡터를 블렌딩하는 비율을 높여가게 된다.

두 영상의 기울기 벡터 크기 또한 블렌딩 비율에 영향을 주게 되며 식 (5)의 β값이 이에 해당된다. β값은 식 (6)에 의해서 결정된다.

식 (4), 식 (5)에서 두 영상의 기울기 벡터의 차이와 크기에 관련된 두 항에는 가중치 값 ωdistance(q)가 곱해진다. 앞 단계에서 생성한 두 영상의 거리 맵을 이용하여 각 위치에서의 가중치 값을 식 (7)에서와 같이 다르게 할당한다. D1(q)와 D2(q)는 두 영상 I1과 I2의 중첩 영역에 대한 거리 맵의 해당 위치에 거리 값을 나타낸다.

제안 기법은 두 영상의 거리 맵 값의 차이를 이용하여 중첩 영역의 외곽선에 가까운 영역에서는 가중치를 매우 작게 주어서, 이 영역에서 블렌딩 비율이 변화하는 것을 억제시킨다. 중앙으로 갈수록 가중치는 빠른 속도로 커지며 기울기 벡터의 차이와 크기에 따라 블렌딩 비율을 크게 변화시키게 된다. 이러한 비율 조절은 비네팅 효과로 인한 심 결함을 제거하는데 매우 효과적이다.

식 (4)에서 계산된 P(q′) 값을 확장이 진행되는 방향의 위치에 저장한 후, 식 (8)과 같이 두 영상의 기울기 벡터를 P(q′) : 1-P(q′) 비율로 블렌딩하여 최종 영상에서 현재 위치의 밝기 값 I(q)에 더해줌으로써 새로운 위치의 밝기 값 I(q′)을 계산한다. I(q′)를 계산하는 식은 다음과 같다.

개선된 영역 확장 방식으로 영상 I1의 물체 영역 좌측 상단의 한 점으로부터 확장해 나가며 I(q′)를 계산하기 때문에 한 점당 최대 4번 씩 주변에 위치한 점들로부터 영향을 받게 된다. 이 과정으로 두 기울기 벡터의 사용이 전환되는 위치가 각각 다르기 때문에 발생할 수 있는 노이즈를 최소화할 수 있다.

 

3. 치과 CT 영상 스티칭 기법의 설계 및 구현

3.1 시스템 구성

본 논문에서 개발한 치과 CT 영상 스티칭 기법을 위하여 Fig. 3과 같은 세부 모듈을 갖는 시스템을 설계 및 구현하였다. 제안 시스템은 Dicom 영상 입력 모듈, 영상 정합 모듈, 객체 분할 모듈, 영상 스티칭 모듈, 스티칭 결과 출력 모듈로 구성되고, 각 모듈은 이전 단계의 모듈의 출력물을 입력으로 받아 순차적으로 실행된다.

Fig. 3.Flowchart of system organization.

3.2 시스템 설계 및 구현

제안 시스템을 구성하는 각 모듈은 다음과 같은 기능 명세를 갖고, 실제 해당하는 기능을 수행하게 된다. 먼저 Dicom 영상 입력 모듈은 치과 CT 영상 파일을 입력으로 받아 dicom 헤더를 분석하여 두 CT 영상의 실제 영상 데이터를 입력으로 받는다. 다음으로 영상 정합 모듈은 두 CT 영상 데이터를 정합하는 최적의 정합 인자를 탐색한다. 객체 분할 모듈은 영상에서 배경과 객체를 구분하고, 객체의 최외곽 경계를 추출한다. 다음으로 영상 스티칭 모듈은 정합된 영상들에서 스티칭된 새로운 영상을 생성하는 역할을 한다. 마지막으로 스티칭 결과 출력 모듈에서 스티칭된 영상을 입력 영상과 함께 시각화한다.

 

4. 실험 결과

제안 기법에 대한 실험은 Intel Core i5 760 2.8GHz CPU와 4GB 메모리를 갖는 시스템에서 수행되었다. 개발 환경은 Microsoft Visual Studio 2010을 이용하여 C++ 언어로 개발하였다. 제안 기법은 두 영상(기준(reference) 영상, 부유(floating) 영상)이 각각 248 × 248 × 248의 크기를 갖는 10명의 치과 CT 영상 데이터에 적용이 되었다. 이 중 첫 번째 데이터에 대하여 α 값을 5부터 50까지 5의 간격으로 변화시키고, T 값을 50부터 500까지 50의 간격으로 변화시키면서 각각의 조합에 대하여 시각적으로 가장 정확한 결과를 보여주었던 인자 값을 찾았다. 그 결과 식 (5)의 α 값은 10으로 식 (6)의 T 값은 150으로 설정하였다.

Fig. 4(a)의 기준 치과 CT 영상에 Fig. 4(b)의 부유 치과 CT 영상을 본 논문에서의 제안 기법으로 정합 후 스티칭하면, Fig. 4(c)와 같이 스티칭 결과 영상을 얻을 수 있다. Fig. 4(c)에서와 같이 제안 기법으로 밝기 값 차이와 심 결함을 제거한 좋은 품질의 합성 영상을 얻을 수 있음을 확인할 수 있다.

Fig. 4.Result image of proposed method. (a) Reference image (b) Floating image (c) Stitching result image

다음으로 제안 기법의 스티칭 결과와 기존 기법의 스티칭 결과를 비교하였다. 기존 기법으로는 합성이 될 두 영상 중 하나의 영상에서 밝기값을 선택하는 단순 합성 기법과 블렌딩 비율을 1:1로 하는 기법 두 가지에 대하여 비교를 수행하였다. 먼저 Fig. 5에 단순 합성 기법과 제안 기법과의 스티칭 결과를 비교하였다. Fig. 5(a)와 5(b)에서 볼 수 있듯이 두 영상의 평균 밝기 값 차이가 크기 때문에 기존 기법으로는 Fig. 5(c)와 같이 심 결함이 매우 심하게 나타났다. 본 논문의 제안 기법을 사용하면 Fig. 5(d)와 같이 평균 밝기 값 차이를 복원하고 심 결함을 완벽하게 제거한 스티칭 결과 영상을 얻을 수 있었다.

Fig. 5.Comparison image. (a) Reference image (b) Floating image (c) Previous method (d) Proposed method

다음으로 Fig. 6에 블렌딩 비율을 1:1로 합성한 기존 기법과 제안 기법과의 스티칭 결과를 비교하였다. Fig. 5의 영상 데이터와는 달리 Fig. 6에 제시된 영상 데이터는 두 영상을 정합 후 지역적인 정합 오류가 발생하였다. Fig. 6(a)와 Fig. 6(b)에 나타난 영상을 정합하게 되면 Fig. 6(c)와 같이 전역적으로는 두 영상이 정합된 것처럼 보이지만, 두 영상의 중첩 부분에서 1:1의 비율로 블렌딩하는 방식으로 합성하면 고스팅 결함이 나타났다. 각 영상에서 하나의 선이었던 어둡고 길게 늘어진 노이즈가 Fig. 6(c)에서는 흐려진 두 개의 선으로 나타나는 것을 볼 수 있었다. 또한, 중첩 영역의 경계선 부분에서 미세한 밝기 값 차이로 인해 심 결함이 나타나는 것도 볼 수 있었다. Fig. 6(d)는 본 논문의 제안 기법으로 스티칭한 결과 영상이다. 제안 기법은 중첩 영역의 경계선으로부터 거리와 주변의 블렌딩 비율로 계산된 가중치를 두 영상의 기울기 벡터의 차이와 크기에 적용하여 각 위치에서 새로운 블렌딩 비율을 계산해 낸다. 그리고 이 비율에 따라 두 기울기 벡터를 블렌딩하기 때문에 Fig. 6(d)에서 볼 수 있듯이 중첩 영역의 경계선에서 생기는 심 결함을 제거하였고, 비록 노이즈이긴 하지만 길게 늘어진 노이즈가 연결되어 하나로 표현될 수 있는 최적의 심을 블렌딩 비율 계산을 통해 찾음으로써 지역적 정합 오류를 해결하였다.

Fig. 6.Comparison image. (a) Reference image (b) Floating image (c) Previous method (d) Proposed method

10개의 영상 데이터에 대하여 중첩되기 전의 두 CT 영상에서 의사가 각각 5개의 landmark point를 정의하고, 이들 사이의 거리를 의사가 반자동으로 정합한 결과에 대하여 제안 기법으로 스티칭을 수행한 후 평균 거리 오차를 측정하였다. 10개의 영상 데이터에 대하여 평균 거리 오차는 0.37mm 였고, 픽셀로 환산할 경우 평균 1픽셀 정도의 정합 오차로 매우 정확한 결과를 보여주었다.

마지막으로 10개의 영상 데이터에 대하여 제안 기법의 수행 시간의 평균을 측정하였다. 영상 정합 과정의 경우 평균 3.8초의 시간이 소요되었고, 영상 합성의 경우 평균 12.9초의 시간이 소요되었다.

 

5. 결 론

본 논문에서는 3차원 다중 치과 CT 영상의 고화질 스티칭 기법을 제안하였다. 제안 기법은 먼저 중첩 영역의 경계선을 기준으로 두 영상의 거리 맵을 생성한 뒤, 각 위치에서 두 거리 맵의 차이 값을 이용하여 가중치 함수를 구성하였다. 다음으로 두 영상의 기울기 벡터의 차이와 크기에 앞서 계산한 가중치 함수를 적용함으로써 기울기 벡터를 블렌딩할 비율을 계산하고 이를 개선된 영역 확장 방식으로 전파시키며 블렌딩할 영상의 밝기 값을 계산하였다. 기울기 벡터를 이용하기 때문에 영상 간의 평균 밝기 값 차이를 보정할 수 있었고, 거리 맵의 차이 값을 이용한 가중치 함수가 중첩 영역의 경계선 부분에서 블렌딩 비율이 변하는 것을 억제시키기 때문에 심 결함 없이 영상을 스티칭할 수 있었다. 하지만, 치과 CT 영상에는 잡음이 많은 편이고, 최종 결과 영상에도 많은 잡음이 존재하게 된다. 향후 연구로 치과 CT 중첩 영상에 잡음을 최소화하는 방안을 연구할 계획이다.