• • Gravity measurements and forward modeling analysis are performed to estimate the vertical and horizontal structure of the Yangsan Fault core line

• • A prominent gravity difference of more than 1.5 mGal across fault core, indicating density difference of the crust

• • The fault core is elongated vertically to a depth of about 2,000m and extended horizontally 3,000m to the NNE direction

2016년 9월 12일 한반도에서 발생한 계기지진 중 최대 규모의 경주지진(그림 1A의 2016M5.8_GY)과 2017년 11월 15일, 두 번째 규모의 포항지진(그림 1A의 2017M5.4_PO)이 일어나자 대규모의 지진을 일으킨 단층이 어디인가에 대한 관심이 일어났다. 일반적으로 경주 및 포항지진과 같은 큰 규모의 지진들은 양산단층과 같은 대형 단층의 움직임에 의한 것으로 알려져 있다. 그래서, 큰 규모의 지진이 일어나기 위해서는 관련된 활성 단층이 있어야하며, 이런 활성단층의 위치를 알아내는 것은 현재의 지진활동뿐만 아니라 미래의 지진활동을 예측하는데 매우 중요하다는 것이 틀림없다.

Figure 1. (A) The study area (black solid rectangle) locates in the northern part of the Yangsan fault line (white dotted line). 2016M5.8_GY: Gyeongju earthquake epicenter, 2017M5.4_PO: Pohang earthquake epicenter. (B) Geologic map indicates that the Yangsan fault line in the study area is the boundary between the Cretaceous biotite Granite and the Cretaceous grey Sandstone. (C) presents measured gravity stations (dots) and electrical resistivity survey line (DD-08) in the Naengsuri. 2D electrical resistivity structure presented in (D) will be used as one of the important constraints to model density structure.

경주 및 포항지진 이후 양산단층의 활성화 여부를 판단하고자 다양한 지구과학적인 연구가 진행되고 있다(예, Kim et al., 2020a; Kim et al., 2020b). 단층은 지각이 갈라져 깨진 부분이라 암석이 주위에 비해 상대적으로 약하기 때문에 지표면에서는 보통 선상의 저 지형을 형성하는 경우가 많다. 따라서 일반적으로 다음과 같은 방식으로 단층의 활성화 여부에 대한 연구를 실시한다. 먼저, 위성사진이나 LiDAR (Light Detection and Ranging) 등을 이용해서 이러한 지형적 특징과 제 4기 지질시대에 발생한 지질학적인 변이 정도를 확인한다. 다음으로는 야외 조사를 실시하여 이 단층의 활동이력과 시기를 파악하는 작업을 진행한다. 마지막으로 굴착(트렌치) 조사를 통해 단층의 최근 활동이력, 활동시기, 재발주기 그리고 변위량 등 고지진에 대한 정보를 얻어내어 활성 단층 여부를 판단할 수 있지만(kim et al., 2020b), 한반도에서 발생한 거의 모든 지진의 진원(Hypocenter)은 최소 5 km 이상의 심부에 위치하고 있다(예, 경주 및 포항 지진의 진원지는 각각 약 12 km, 5 km, Choi et al., 2020). 그렇기 때문에, 위에 언급한 방법 이외의 적절한 연구 방법이 필요하다.

이 연구는 중력장을 측정하고 그 데이터를 분석하여 지질학적으로 확인된 천부 단층의 위치를 확인하고 중력 이상의 원인이 되는 지각 밀도 분포를 모델링하여 단층의 3차원적인 규모를 결정하는 해석 과정을 제시하고자 한다. 연구를 위해 선정한 지역은 경상북도 포항시 북구 신광면 냉수 1리 지역(그림 1B)인데 최근 실시된 LiDAR 판독, 트렌치 조사 및 전기 비저항 탐사(그림 1D)를 종합한 지질학적인 연구에서 단층이 활성화 되었다고 보고된 곳으로(kim et al., 2020a), 본 연구에서 제시할 중력장 데이터 분석을 통한 단층의 수직 및 수평적 규모를 평가하는데 매우 적합한 지역이라고 판단된다.

연구 지역은 경상북도 포항시 북구 신광면 냉수리 일대이다(그림 1B). 이 지역은 2017년부터 시작된 국가활성단층지도 제작(Kim et al., 2020b)의 하나로 양산단층(그림 1A와 1B의 파선)을 따라서 수행된 여러 지표 조사 지역 중의 한 지역으로, 위성사진과 LiDAR 분석을 통해서 양산단층이라고 추정되는 선형 구조(그림 1C의 파선)가 뚜렷하게 나타나는 곳이다(Kim et al., 2020a; Kim et al., 2020b). 특히 냉수1리(그림 1B의 사각형) 지역에서 실시된 트렌치 조사에서는 선형 구조가 나타난 부분에서 아주 뚜렷하게 변형된 제4기 사질층의 미구조와 맞물린 부분적인 파쇄구조가 발견되었는데, 이와 같은 파쇄층은 지진성 제4기 단층운동이 발생하여 지표 근처의 기존 단층핵과 미고결 퇴적층 내에 빠른 미끌림 현상이 일어났다는 것을 보여준다(Kim et al., 2020a). 이 파쇄층은 서쪽의 불국사 화강암류에 속하는 흑운모 화강암(그림 1B와 그림 1C의 Bulkuksa Granite)과 동쪽의 경상분지 하양층군 대구층(그림 1B와 1C에서 Cretaceous Sediments)의 경계에 해당하는데 이는 제4기에 발생했던 활성 단층을 의미한다(Kim et al., 2020a). 트렌치 조사와 더불어 실시된 전기비저항 측정 결과(그림 1D)도 이 파쇄층에 아주 뚜렷한 전기 비저항 이상대가 나타나고 있으며, 이 비저항대가 지표면으로부터 50 m 깊이까지 연장되어 있다는 것을 보여주고 있다. 결론적으로 트렌치 조사 및 전기비저항 탐사는 이 선형구조가 활성단층이라는 것을 강력히 시사하고 있다(Kim et al., 2020a). 한편 이 지역에서 채취한 불국사 화강암의 밀도는 약 2.60 g/cm³였고 경상계 퇴적암은 약 2.67 g/cm³ 이었다(Park et al., 2009; Kim et al., 2019). 이 암석 밀도는 아래의 중력장 해석에서 매우 중요한 역할을 한다.

연구는 단층 지역에서 측정된 중력장 데이터의 역산과 순산을 수행한다. 첫째로 오일러 디콘볼루션(Euler deconvolution) 방법을 이용하여, 천부에서 확인된 단층의 수직적 분포를 확인한다. 둘째, 곡률분석(Curvature analysis) 역산 식을 이용하여, 지표면에서 확인된 단층의 수평적 분포를 가시화한다. 셋째, 3차원 순산 모델링 프로그램(예, IGMAS+, Schmidt et al., 2004)을 이용하여 단층의 수평 및 수직적인 밀도 구조를 모델링 하여 중력장의 변화를 수반하는 지각 내부의 3차원적인 밀도 구조를 판단하고 지질 및 지구조학적인 의미를 찾고자 한다.

3.1. 오일러 디콘볼루션(Euler deconvolution) 방법

중력장 변화의 원인이 되는 밀도의 경계 및 깊이(source depth)를 3차원적으로 계산하기 위해서 주로 많이 사용되는 계산 방법으로 오일러디콘볼루션 방식이 있다. 이 식의 수학적 표현은 다음과 같다(Thomson, 1982; Reidet al., 1990; Pašteka et al., 2009).

여기서 x, y, z는 측정 점의 좌표이며, x₀, y₀, z₀는 얻고자 하는 밀도 변화 경계 지점의 좌표이고, g는 측정된 중력 값, Δg은 중력값의 1차 미분값이다. N은 위 식에서 가장 중요한 변수로서 구조지수 (structural index, SI) 라고 하며, 위 식을 계산하기 위해 임의로 주어진 변수다. 이 구조지수(SI) 개념을 이해해야만 지질학적인 밀도 구조의 깊이와 위치를 계산 할 수 있다. 연구에서 사용한 오일러디콘볼루션 알고리듬(REDGER, Pašteka et al., 2009)에서는 지각 내의 밀도 변화의 평균 깊이(center of source depth)를 계산하기 위해서는 SI 값을 2로 주어야 하며(그림 2A의 사각형 표시), 지표면과 가장 가까운 깊이(Minimum source depth)를 계산하기 위해서는 0 에 가까운 값(그림 2A의 원 표시)을, 그리고 단층면과 같은 곳(그림 2A의 십자 표시)을 찾아내기 위해서는 1의 값을 주도록 설계되었다(Pašteka et al., 2009). 식(1)에서 2차원 중력 측정 데이터를 위해서는 y 축에 대한 모든 값을 0 으로 처리해야 하므로 식(2)와 같이 정리된다.

Figure 2. Simplified schema of Euler deconvolution method (A), Curvature analysis (B), and 3D modeling algorithm to understand theoretical background of gravity inverse interpretations and forward modelling (C). SI: Structure Index, RHO: Density

3.2. 곡률 분석(Curvature analysis) 역산

중력장의 곡률 분석은 지표면 가까운 곳에 존재하는 밀도의 경계면을 알아내기 위해서 사용되는 방법이다. 이 방법은 초기 3차원 음파탐사 자료(Roberts, 2001)의 분석을 위해서 개발되었으나, 중력장 자료를 이용한 지표면 밀도 분포에도 다양한 적용이 가능해졌다(예, Choi et al., 2020). 그림 2B에 방법의 이론적인 배경을 제시하였다. 밀도가 같은 2 개의 지표면 하부 구조 (RHO1, RHO3) 사이에 밀도가 낮은 한 개의 지각 구조(RHO 2)가 끼어있다고 가정하면, 그로 인해서 나타나는 이론적인 중력 이상값은 가운데가 음의 값을 나타내는 것은 당연하다. 중력장의 각 측점의 곡면(Curvature, 식(3)의 k)은 그 측점에 연결된 원의 반경(그림 2B에서 R)의 일차 미분식 즉 곡면의 법선(Tangent, 그림 2B에서 T)로도 나타낼 수가 있다. 다른 식으로 표현한다면, 한 곡면의 점은 곡선의 각(그림 2B에서 dω(ds))과 길이에 좌우되므로, 곡면의 기울기는 k = dω/ds 이다. 이로부터 다음과 같은 식이 유도된다.

여기서 dx, dy, dz는 수평 및 수직 방향의 1차 미분값들이다. 위의 식을 좀 더 간단하게 2차 미분식으로도 표현할 수 있다(Roberts, 2001). 만약 이차 미분식에 의해서 계산된 값(k)이 0에 가까이 접근을 하면, 이는 이 지역의 하부에 중력장에 영향을 미치는 밀도의 변화가 없다는 것을 뜻한다(그림 2B에서 노란색과 초록색으로 표시된 지점). 그러나 밀도가 음의 값에서 양의 값으로 변하는 경계 지역의 2차 미분 값(k)은 당연히 양의 값을 갖게 되고(그림 2B에서 빨간색으로 표시된 지점), 반대의 경우에는 음의 값을 갖게 된다(그림 2B에서 파란색으로 표시된 지점). 그러므로 전체 중력 측정 지역을 Curvature analysis를 이용하여 계산했을 때, 그림 2B의 아래 그림처럼, 빨간색과 파란색이 교차되어 나타나는 곳은 밀도의 변화가 지각 내에 있다는 것을 나타내며, 이는 곧 지질의 경계 내지는 단층일 가능성을 의미한다.

3.3. 순산 모델링(IGMAS+)

위에서 언급한 방식들은 측정된 중력 값들을 해석하여 역으로 원인이 되는 지각 물질의 지구물리적 성질(예, 밀도, 대자율)을 규명하는 것이다. 반대로 중력 및 자기장을 이용한 순산 모델링 해석은 다음과 같은 방식으로 진행된다. 첫째로, 원인이 되는 지각 물성(예, 밀도, 대자율)의 분포를 먼저 가정(그림 2C의 모델링 참조)한 다음에 순산 방식의 알고리듬(예, Talwani et al., 1959; Götze & Lahmeyer, 1988)을 이용하여 지각 물성의 중력 및 자기장 효과를 계산하는 것이다. 둘째로, 이렇게 계산된 값들을 실제로 측정된 중력 및 자기장 값들과 비교한다. 마지막으로 계산된 효과 값들이 측정된 실제 값들과 만족할 정도로 근접했을 때 제시된 지각 모델을 지질학적으로 해석한다.

중력장을 이용한 2차원 순산 모델링 알고리듬은 1959년 Talwani에 의해서 처음으로 그 이론이 제시되었으며, 3차원 모델링의 이론은 1988년 Götze에 의해서 정립되었다. 2차원 방식은 2차원 평면구조 모델의 좌에서 우로 주어진 밀도와 지각구조의 이론적인 중력 및 자기장 값들을 적분 방식으로 계산하는 것이다(Talwani et al., 1959). 이 경우 모델 라인에 존재하는 밀도 구조에 대해서는 매우 효과적으로 그 중력 및 자기장의 효과 값들을 계산할 수 있으나, 모델 라인에 포함되지 않은 곳에 위치하는 중력 및 자기장의 효과는 고려하지 못하는 단점을 갖고 있다. 이와 같은 단점을 보완하기 위해서 제시된 방법이 3차원 계산 방식이다. 그림 2C에서 보이는 것처럼, 몇 개의 2차원적 모델(그림 2C의 plane)들이 삼각형(Triangulation)으로 연결되어 있다고 가정할 경우, 각각의 2차원 모델은 Green 2차 적분 식으로 계산되어 질 것이다. 이렇게 각각 계산된 2차원 모델은 가우스 3차 적분 식으로 연결되어 3차원적인 중력 및 자기장 효과 값들을 계산해 낼 수 있다(Götze & Lahmeyer, 1988). 모델링 과정에서 가장 중요한 것은 초기에 주어져야 할 지각 모델의 정확성이다. 그렇기 때문에, 모델을 결정할 수 있는 가능한 많은 연구 자료가 종합적으로 비교 분석되어야 한다. 지난 30년 동안 베를린 자유대학과 킬 대학 중력측정 연구소에서 개발한 3차원 중력장 및 자기장 순산 모델링 프로그램(IGMAS+)은 3D 지반정보(geoinformation system)와 객체 및 응답형 기능(interactive function)이 탑재되어 다양한 형태의 연구 결과들을 3차원적으로 비교 분석하여 초기 지각구조 모델을 발전시킬 수 있다. 그런 기능으로 인해서 IGMAS+에 의해서 제시된 최종 모델은 3차원 밀도 및 자기장 분포에 대한 신뢰성을 높일 수 있으며, 지질학적인 해석에 지구물리학적인 의미를 부여할 수 있을 것이다.

4.1. 중력 측정

활성 단층에 대한 확실한 검증 및 단층의 깊이 확인을 위해서 트렌치 조사와 전기비저항 탐사를 수해한 측선(그림 1C의 DD_08 측선)으로부터 약 100 m 북쪽에 위치하고 있는 도로를 따라서 중력 측정하였다(그림 3A 참조). 사용된 중력계는 Scintrex사의 CG5이며 이 장비의 정밀도는 0.001 mGal 이나, 현지 측정시 평균 오차 값은 약 0.01 mGal 정도였다. 중력 보정을 위해서는 중력 측정 지점의 위치, 특히 고도를 정확히 측정해야 한다. 고도 측위를 위해 사용한 장비는 Trimble 사의 4000SSE이며 평균 오차는 약 5 mm이며, 이 오차로 인해서 발생할 수 있는 중력 오차 값은 약 0.001 mGal 으로 중력 해석에 영향을 줄 수 있는 정도보다 낮은 값이다. 중력 측정 시의 간격은 최소 15 m에서 최대 40 m로 현지 상황에 따라서 조정하였으며 평균 측점 간격은 약 20 m였다. 측정 결과로서 길이 약 580 m의 구간(그림 3A 참조)에 걸쳐서 총 30점의 중력 및 고도 측정값을 확보하였다. 지각 내의 밀도 분포를 계산해 내기 위해서는 획득된 중력 측정값들을 부게이상 값으로 전환하여야 한다. 이를 위해서 사용한 부게판(Bouguer plate)에 전 지구 평균 밀도인 2.67 g/cm³을 적용하였다.

Figure 3. (A) Gravity measurement line parallel to the electrical survey line (ER_DD-08). (B) The complete Bouguer anomaly along the WE profile with a mean value of approximately 28.5 mGal and a range from 27.7 to 29.5 mGal. The lower Bouguer gravity anomalies are generally observed in the western part of the Yangsan fault line, while the area to the east is characterized by higher anomalies. A prominent anomaly change from 27.5 to 29.0 mGal is well coincided with the low resistivity area (LRZ in A).

그림 3B는 냉수1리 지역을 서쪽에서 동쪽으로 측정한 부게이상 분포이다. 평균 부게이상은 약 28.5 mGal이며 최소 약 27.7 mGal에서 최대 약 29.5 mGal의 변이를 보이고 있다. 상대적으로 낮은 중력 값들은 양산단층 서쪽에 분포하는데, 이는 지질도에서 나타난 것처럼 주변에 비해서 밀도가 낮은 불국사 화강암류에 속하는 흑운모 화강암(그림 1B, 1C에서 Bulkuksa granite, Kim et al., 2020b)의 분포와 밀접한 관계가 있는 것으로 추측한다. 한편, 평균보다 높은 부게이상은 주로 양산단층의 동쪽에서 나타나는데 이는 이 지역에 불국사 화강암보다 밀도가 높은 물질이 존재한다는 것을 추측하게 한다. 지질도(그림 1B)와 비교하면 양산단층 동쪽의 높은 부게이상은 백악기 경상분지 하양층군 대구층(그림 1B와 1C의 Cretaceous sediments)을 포함한 동쪽 지층이 흑운모 화강암이 주를 이루고 있는 서쪽의 지층보다 상대적으로 높은 밀도를 갖는다는 것을 의미한다. 또한, 트렌치 조사에서 파쇄대로 나타난 지역과 전기비저항 탐사에서 저비저항 이상대(그림 3의 Low resistivity zone)가 나타난 150 m 이내의 구간에서 중력장은 약 28.0 mGal에서 약 29.5 mGal로 급격하게 상승한다. 이는 이 구간을 경계로 서쪽과 동쪽의 밀도 분포가 뚜렷하게 변화한다는 것을 나타내며, 이런 경계 구간(그림 3B에서 서쪽으로부터 250–400 m 사이)의 지하에는 서쪽의 불국사 화강암층과 동쪽의 대구층이 혼합된 경계 지층(파쇄대)이 존재한다는 것을 의미한다. 연구지역에서 실시된 약 600 m의 중력장 데이터 분석 결과는 지금까지 제시된 트렌치 조사와 전기비저항 탐사가 보여준 파쇄대 부분(그림 3A에서 LRZ로 표시된 부분)이 아주 뚜렷한 밀도의 경계면에 위치하고 있다는 것을 보여준다.

다음 장에서는 전기비저항 탐사에서 파쇄대로 표현된 양산단층의 깊이를 계산하기 위해서 2차원 Euler deconvolution 방식(식 2)을 사용하여 밀도의 경계면(양산단층)의 깊이를 계산하고자 한다(이론적 배경은 그림 2A 참조).

4.2. 단층의 수직적 분포

위에서 언급한 것처럼, 냉수리 지역의 서쪽에서 동쪽으로 측정한 부게이상의 평균값은 약 28.5 mGal, 그 변이는 최소 약 27.7 mGal에서 최대 약 29.5 mGal 이며(그림 4A), 양산단층의 파쇄대(그림 4의 Fault core) 지역은 중력 변이가 가장 심하게 나타나는 구간이다. 이 파쇄대의 최대 깊이를 계산하기 위해서 2차원 Euler deconvolution 방법(식 2)을 사용하였다. x축의 미분 값(그림 4B에서 xderivative로 표시)은 파쇄대 지역에서 뚜렷한 변곡선을 나타내고 있으며, 수직 방향에 대한 미분 값(그림 4B에서 z-derivative로 표시)의 변이 역시 파쇄대 구간에서 다른 지역에 비해 뚜렷하게 상승 곡선을 나타내고 있다. 중력의 변이를 야기한 원인 밀도물질의 평균 깊이(source depth) 계산을 위해서, 위의 두 미분 값들과 SI=2의 값(Pašteka et al., 2009)을 수식 E2에 대입하여 평균 깊이를 계산하였다. 그 결과(그림 4C)는 양산단층에 해당하는 파쇄대의 깊이가 전기비저항 탐사에서 확인된 50 m보다 더 깊은 100 m 이상이었다. 그러나 이 깊이가 양산단층에 존재하는 파쇄대의 평균 깊이를 의미하는 것은 아니다. 그 이유는 중력 측정 거리에 비례하여 깊이 분포에 대한 한계가 있기 때문인데, 실험 결과에 의하면(Pašteka et al., 2009) 대략 중력장 변이 거리의 절반에 해당하는 깊이까지 해석이 가능하다고 보고되었다. 그렇기 때문에 더 깊은 곳까지 파쇄대가 연결되어 있을 가능성을 확인하기 위해서는 연구 지역의 범위를 확장하여 광역의 중력장 데이터의 해석이 필요하다.

Figure 4. (A) Bouguer anomalies across the Yangsan fault line are shown (refer to Fig. 1C and Fig. 3A for the location of the profile). (B) Application of the Euler deconvolution technique by using the complete Bouguer anomaly field (A) to calculated x and z-derivative. (C) Distribution of source points. The Yangsan fault core is characterized by the prominent disturbances of xand z-derivatives, source depths are about 120 meters, definitely deeper than the mean depth (50 meters).

5.1. 부게이상도

더 넓은 지역의 단층의 수평 및 수직적 도달 범위를 계산하기 위해서 냉수리 측량 지역을 중심으로 동서 방향으로 약 5 km, 남북 방향으로 약 9 km 지역에 대한 부게 이상 지도를 제작하였다(그림 5A). 그림에서 포항분지의 서북쪽 지역은 2011년부터 지열발전소 및 CO₂ 저장시설 설치를 목적으로 한국 지질자원연구소 주관으로 중력장이 측정되었으며(Lim et al., 2019), 그 중에서 신광면 일대에는 약 80여점의 중력장 데이터가 획득되었다(그림 5A와 5B의 Gravity station). 평균 부게이상은 약 30 mGal이며 최소 약 20 mGal에서 최대 약 40 mGal의 범위를 보인다. 평균보다 높은 부게이상은 주로 신광면을 중심으로 북동쪽에서 나타나는데 이는 이 지역의 지각의 수평 평균 밀도가 남쪽보다 높다는 것을 의미한다.

Figure 5. A shows that the Bouguer anomaly map over the study area (Sinkwangmyen, Pohang) has a mean value of approximately 30 mGal with a range from 20 to 40 mGal. The higher Bouguer gravity anomalies in the NE area from Sinkwangmyen, while the SW is characterized by lower anomalies. (B) The mapped Dip-curvature shows a prominent density contrast along the Yangsan fault line. Several density discontinuities along the Yangsan fault line (marked with S1-S5) are anticipated to be the fault segmentations.

한편 평균값보다 낮은 부게이상 값들은 신광면으로부터 남쪽에 주로 분포하는데, 이는 이 지역에 넓게 분포하는 주변의 경상계 퇴적암(2.67 g/cm³, Park et al., 2009)에 비해서 밀도가 낮은 불국사 화강암(2.60 g/cm³, Park et al., 2009)의 분포와 밀접한 관계가 있을 것이다. 양산단층(그림 5의 파선) 주변의 부게이상도 신광면을 중심으로 남쪽과 북쪽이 아주 뚜렷한 차이가 있다. 신광면 남쪽 지역의 양산단층은 부게이상이 평균보다 낮은 지역을 관통하지만, 주변에 비해서 상대적으로 높은 부게이상 값을 보여주고 있기 때문에 중력이상도 상에서 양산단층의 위치가 비교적 잘 나타나고 있다. 반면 북쪽 양산단층 지역은 부게이상 지도에서 뚜렷한 위치를 찾기 어려울 정도로 주변과 차이가 없다. 그 원인에 대한 것은 앞으로 지질학 및 지구물리학적인 연구가 이루어져야 할 것이다.

이 연구에서 아쉬운 점은 신광면 남쪽 지역(그림 5A에서 사각형으로 표시된 지역) 중에서 양산 단층의 서쪽과 남쪽에는 비교적 고르게 중력 데이터가 분포하고 있지만, 양산단층과 포항 분지 사이의 약 1 km 지역은 획득된 중력 데이터가 없다는 점이다. 이런 이유로 이 지역에서의 중력장 해석의 신뢰성은 다른 지역에 비해서 떨어지는 것은 당연하다.

5.2. 양산단층의 수평적 분포

앞의 3.2에서 설명한 것처럼, 곡률분석(Curvature analysis, 관계식 3)를 이용하여 부게이상(예, 그림 5A)의 곡면 기울기를 계산한 결과를 색깔로 나타내었을 때, 빨간색(maximum)과 파란색(minimum)이 교차되어 나타나는 곳(그림 2B 참조)은 밀도의 수평적인 경계면 즉 단층일 가능성을 의미한다. 그림 5B는 신광면 지역의 부게이상의 곡면 기울기를 계산한 지도이다. 냉수리 일대(그림 5B에서 사각형 영역)에서 확인된 파쇄대는 서쪽의 파란색과 동쪽의 빨간색으로 구분되어 양산단층이 뚜렷한 밀도의 경계에 위치하고 있음을 명확하게 보여주고 있다. 이 파쇄대에 해당되는 색깔의 경계면은 남쪽으로는 약 1,000 m 정도, 북쪽으로는 약 3,000 m 정도 연장되어 있는 것으로 보여 진다. 이는 냉수리 지역에서 확인된 파쇄대(양산단층)의 수평적인 분포를 보여준다. 북쪽으로 연결된 색깔의 경계는 현재 지도상에 표시된 양산 단층보다 서쪽으로 약간 치우쳐 나타나고 있다. 이것은 양산 단층이 지표면에서 보이는 것보다 지각 내에서는 서쪽으로 기울어져 존재할 가능성, 혹은 해석에 사용된 중력 데이터가 충분하지 않아서 나타난 현상일 가능성도 있다. 또한 지도상에 표시된 양산단층은 전 구간이 연속된 것이 아니라 구간별로 분절(segmentation)된 것처럼 보인다(그림 5B의 S1-S5). 정확한 단층의 수평적 연장을 확인하기 위해서는 냉수리에서 실시한 것처럼 좁은 간격의 중력 측정이 양산단층을 전체 지역을 따라서 진행되어야 할 것이다.

5.3. 양산단층의 수직적 분포

양산단층의 수직적 연장을 계산하기 위해 3차원 Euler deconvolution 방식(식 1)을 사용하였다. 중력의 변화를 유발한 원인 밀도물질의 최대 깊이(source depth) 계산을 위해서, 그림 5A에 제시한 부게이상, SI=2의 값(Pašteka et al., 2009)을 식 1에 대입한 후, 3차원 Euler deconvolution software(REDGER, Pašteka et al., 2009)를 이용하여 깊이를 계산한 결과(그림 6A), 연구 지역을 남남서에서 북북동으로 관통하는 양산단층의 평균 깊이는 약 1,500 m이라는 것과 최대 깊이는 약 3,000 m라는 것을 확인하였다(그림 6A에서 색깔로 표시된 네모들).

Figure 6. (A) Results of the application of the Euler deconvolution technique to the Bouguer anomalies in Fig. 5A get information on the 3D distribution of source points. (B) The Bouguer anomalies along the WE profile varies from about 27.0 mGal to about 29.2 mGal. (D) The distribution of the mean source depth showing the fault core extends to a depth of 2 km. Surface details in (C).

그림 6B는 양산단층을 서쪽에서 동쪽으로 가로지르는 라인(그림 6A의 W-E)을 따라서 나타난 중력 값이며, 그림 6C와 6D는 중력의 변이를 야기한 원인 물질의 깊이를 Euler deconvolution 방식으로 계산한 것이다. 그림 6D에서 평균 깊이는 양산단층의 서쪽 지역은 약 1,200 m이며, 양산단층이 있는 곳으로 갈수록 평균 깊이가 약 2,500 m까지 뚜렷하게 하강하다가 더 동쪽으로 가면서 다시 1,200 m로 올라가는 것을 보여준다(그림 6D에서 회색 점). 이는 냉수리에서 계산한 파쇄대(약 150 m, 그림 4C와 그림 6C 참조)가 지표면으로부터 지각 내 2,500 m까지 연장되어 있을 가능성이 높다는 것을 시사한다. 다만, 냉수리 측선에서는 파쇄대의 넓이를 정확하게 확정할 수 있지만, 그림 6D에서 나타난 파쇄대의 넓이는 약 1,000 m 정도로 상당히 불확실한 것을 알 수 있다. 그 원인은 획득한 중력 데이터의 측점 간격에 있는 것이기 때문에, 냉수리 측선처럼 측점 간격을 매우 좁게 그리고 전 지역에 걸치며 비슷한 측점 간격으로 중력 데이터를 획득한다면, 양산 단층대를 따라서 존재하는 파쇄대의 동서 방향의 규모도 확정할 수 있을 것이다.

그림 7B는 양산단층을 따라서 남남서에서 북북동 쪽으로 이어진 라인(그림 7A의 SSW-NNE)을 따라서 나타난 중력 값이며, 그림 7C는 중력의 변이를 야기한 원인 물질의 깊이를 Euler deconvolution 방식으로 계산한 것이다. 이 라인의 부게이상은 최소 28.5 mGal에서 최대 29.2 mGal의 범위(그림 7B)를 보이며, 중력 변화를 수반한 원인 물질의 평균 깊이는 냉수리 남부 지역에서 약 800 m로 가장 낮게 나타나며, 동-서 라인(그림 7B에서 WE-profile로 표시된 지역)이 지나가는 지역에서는 깊이가 약 2,000 m로 가장 깊게 나타난다. 양산단층을 따라서 나타나는 이와 같은 평균 깊이의 변이는 여러 가지 원인이 있을 수 있지만, 양산 단층에 가해진 각 지질시대별 남북 방향의 수평 응력에 의한 것이 아닐까 추측해 본다.

Figure 7. (A) Index map is the same as Fig. 6(A). (B) the Bouguer anomalies along the SSW-NNE profile ranging from about 28.5 mGal to about 29.2 mGal. (C) The mean source depth (black dotted line) shows that the depth of the Yangsan fault core is ranging from 2000 to 800 meters, which can be explained by different depths of fault segments.

6.1. 냉수리 지역

냉수리 지역 트렌치 조사와 전기비저항 탐사 지역(그림 3A)에서 실시한 중력 측정 데이터를 이용하여 2차원 순산 지각구조 모델링을 실시하였다. 그림 8A의 실선은 냉수리 지역을 서에서 동쪽으로 횡단한 부게이상 분포이다(측정 위치는 그림 3A 참조). 중력 데이터를 이용한 모델링에서는 지층이나 암석의 구조(polygon) 그리고 각 구조들의 밀도(density)가 주어져야만 측정 점들에서의 중력효과를 계산할 수 있다(이론적인 배경은 그림 2C 및 3.3절 참조). 지하 구조 모델을 위한 참고 데이터(constraints)로 지하 50 m까지의 전기 비저항 탐사 결과(그림 8B의 사각형 영역) 그리고 역산을 통해 계산해 낸 심도 분포(source depth, 그림 8B의 사각형)를 활용하였다. 또한, 실험실에서 측정한 불국사 화강암의 밀도(2.60 g/cm³, Park et al., 2009)와 하양층군 퇴적암의 밀도(2.67 g/cm³, Park et al., 2009)는 기반암 밀도 데이터로 활용하였다. 다만, 지표면 가까운 곳에는 공극율이 대략 10%에 이른다고 판단하여, 실제 모델링에 입력한 밀도는 화강암, 2.40 g/cm³, 퇴적암, 2.50 g/cm³ 로 하였다. 지표면 가까운 지층들의 밀도는 포항분지에서 획득한 시추자료(Choi et al., 2020)들을 참고하여 20 m까지는 평균 밀도 약 1.50 g/cm³, 20–50 m 사이에서는 평균 약 1.80 g/ccm³으로 하였다(그림 8C). 또한 트렌치 조사에서 나타난 파쇄대(그림 8C에서 Fault core로 표시된 회색 부분)는 서쪽과 동쪽의 암반들이 서로 혼합된 수직 지층(Kim et al., 2020a)으로 판단하여, 양쪽 밀도들의 중간 값들을 깊이에 따라 차등하여 모델링하였다.

Figure 8. The final gravity model along the E-W profile in Naengsuri. The subsurface structure (C) is modeled with iterative modifications by using all available constraints from the interpretation of the gravity field (A), results of electrical resistivity survey (B), and the source depths analysis (B).

그림 8A의 점선은 위에서 언급한 지층 및 암석의 분포와 밀도 값들에 의해서 계산한 부게이상이며, 이 계산값들은 전 부분에 걸쳐서 측정된 부게이상과 오차 범위 내에서 일치하는 것으로 나타났다. 이는 그림 8C의 지하구조와 밀도 분포가 신빙성이 있다는 것을 의미한다. 즉 양산단층의 서쪽 지역에 나타나는 평균값보다 낮은 중력 이상은 주로 밀도가 상대적으로 낮은 기반암인 불국사 화강암에 의한 것이며, 동쪽 지역에 나타나는 높은 중력이상은 반대로 밀도가 상대적으로 높은 기반암인 백악기 퇴적암에 의한 것으로 판단된다. 또한 양산단층을 지나며 급격하게 상승하는 부게이상 값은 두 기반암의 밀도차에 의한 중력 효과인 것으로 설명할 수 있다.

6.2. 신광면 지역

냉수리 모델링에서 파쇄대의 깊이를 확인하기 위해서 양산 단층을 서쪽에서 동쪽으로 가로지르는 약 4km의 중력 모델링을 실시하였다(그림 9A에서 W-E로 표시). 그림 9B에서 보이는 실선(Measured로 표시)은 신광면 지역을 서쪽에서 동쪽으로 측정한 중력의 부게이상이다. 냉수리 모델링(그림 8C)의 경우와 같이 지하 구조 모델(polygon)들을 위한 참고데이터(constraints)로 역산을 통해서 계산해 낸 심도 분포(source depth, 그림 9C)를 활용하였다. 또한, 실험실에서 측정한 불국사 화강암의 밀도(2.60 g/cm³, Park et al., 2009)와 하양층군 퇴적암의 밀도(2.67 g/cm³, Park et al., 2009)는 기반암의 밀도 데이터로 활용하였다. 다만, 지표면으로부터 약 1,000 m까지는 대략 5%의 공극을 포함(Kim et al., 2020a)하고 있다고 판단하여, 실제 모델링에 주어진 밀도는 화강암 2.50 g/cm³, 퇴적암 2.60 g/cm³로 정하였다. 지표면 가까운 지층들의 밀도는 포항분지에서 획득한 시추자료(Choi et al., 2020)들을 참고하여 150 m까지는 평균 밀도 약 1.90 g/cm³으로 하였다(그림 9C 참조). 또한 트렌치 조사에서 나타난 파쇄대(그림 9C에서 Fault core로 표시된 회색 부분)는 서쪽과 동쪽의 암반들이 서로 혼합된 수직 지층으로 취급하여, 양쪽 밀도들의 중간 값을 깊이에 따라 차등 적용하였다. 그림 9B의 파란색 파선은 위에서 언급한 지층 및 암석의 분포와 밀도 값들로 계산한 부게이상인데, 이는 측정된 부게이상과 전 부분에 걸쳐 오차 범위 내에서 일치하는 것으로 나타났다. 이는 그림 9C에서 제시한 지하구조 및 밀도 분포 모델이 어느 정도의 신빙성이 있다는 것을 의미한다. 즉 양산단층의 서쪽 지역에 나타나는 평균값보다 낮은 중력 이상은 지표면 가까이로부터 지하 2,000 m 층까지 분포하고 있으며, 밀도가 상대적으로 낮은 기반암인 불국사 화강암(그림 9C에서 GR)에 의한 것이라고 판단된다. 또한 동쪽 지역에 나타나는 높은 중력 이상은 화강암과 같은 두께를 가지면서도 밀도가 상대적으로 높은 기반암인 백악기 퇴적암(그림 9C에서 CS)에 의한 것으로 판단된다. 한편, 양산단층을 지나 동쪽으로 급격하게 상승하는 부게이상 값은 두 기반암의 밀도 차이에 의한 중력 효과인 것으로 설명되며, 두 기반암 사이에 위치한 파쇄대의 깊이는 최소 2,000 m에 이르는 것으로 모델링 하였다.

그림 9D는 양산단층을 따라서 나타난 부게이상 변화이다. 남남서에서는 매우 낮은 중력 이상(약 28.5 mGal)이 나타나며 북북동쪽으로 갈수록 중력 값이 증가하다가 최대값인 29.2 mGal에 다다른 뒤에는 다시 급격하게 하강하는 모습을 보이고 있다. 모델링 결과(그림 9E), 이와 같은 중력 변화의 원인은 양산단층을 따라서 분포하는 암상들의 깊이 차이에 의한 것으로 판단된다(예로 냉수리에서는 약 800 m, 서-동 라인에서는 약 2,000 m). 이러한 깊이 차이의 원인은 양산 단층에 가해진 남북 방향의 응력에 의해서 형성된 분절 현상(그림 5B)에 의한 것으로 추측되나, 더 정확한 결론을 도출해 내기 위해서는 신광면 지역에서도 냉수리에서 실시한 정밀 중력 측정이 선행되어야 할 것이다.

Figure 9. (A) Geological map with locations of E-W and SSW-NNE profiles. Final gravity models (C and E) along the E-W and SSW-NNE profiles (locations refer to Fig. A) are shown, which are constrained by the interpretation of the gravity field (B and D) and calculated source depths (Yellow rectangles in Fig. C and E).

이 연구에서는 지표면에서 확인된 활성단층의 3차원적인 규모를 파악하기 위해서 실시한 중력장 측정과 이의 해석(예, Euler deconvolution, Curvature analysis) 그리고 지각구조 모델링에 대한 이론 및 실제 적용사례를 제시하였다.

연구 지역은 경상북도 포항시 신광면 냉수 1리 지역이다. 이곳의 양산 단층을 따라서 실시된 LiDAR 탐사, 트렌치 조사 및 전기비저항 탐사 결과를 통하여 이 지역에서 남남서에서 북북동으로 이어진 양산 단층의 파쇄대를 확인하였다. 연구 지역에서 실시된 중력장 측정 및 해석 결과는 오일러디콘볼루션에 의한 불연속면의 깊이가 지표면으로부터 수직적으로 지각 내 약 2,000 m까지 이르며, 수평적으로는 남남서에서 북북서 방향으로 최소 3,000 m 정도 연장되어 있는 것으로 분석되었는데. 저자들은 단층 파쇄대의 공간적 분포 규모도 이와 관련이 있을 것으로 추측한다. 이러한 결과는 지금까지 주로 지표면에서 실시된 2차원적인 활성단층 확인 절차와 방식을 3차원 규모로 확대하는 중력장 해석 및 순산 모델링 방식이 매우 효과적이라는 것을 의미한다. 아울러 이 연구는 활성단층의 정확한 규모를 파악하기 위해서는 냉수리에서 실시한 것과 같은 정밀 중력 측정과 해석 연구가 필요하다는 것을 시사하고 있다. 다만 중력장 해석의 신뢰성을 높이기 위해서는 다양한 지질 및 지구물리학적인 연구 결과를 공유하는 것이 매우 중요하다.

이 연구는 한국연구재단 이공분야기초연구사업(2020R1F1A1054863)의 지원으로 수행되었습니다.