Single-Channel Seismic Data Processing via Singular Spectrum Analysis

Woodon Jeong; Chanhee Lee; Seung-Goo Kang

doi:10.7582/GGE.2024.27.2.091

Preview

Research Article

Geophysics and Geophysical Exploration. 31 May 2024. 91-107
https://doi.org/10.7582/GGE.2024.27.2.091

Single-Channel Seismic Data Processing via Singular Spectrum Analysis

특이 스펙트럼 분석 기반 단일 채널 탄성파 자료처리 연구

Woodon Jeong¹^*

Chanhee Lee¹

Seung-Goo Kang²

정 우돈¹^*

이 찬희¹

강 승구²

¹Department of Geophysics, Kangwon National University, Chuncheon, South Korea

²Division of Glacier & Earth Science, Korea Polar Research Institute, Incheon, South Korea

¹강원대학교 지구물리학과

²한국해양과학기술원 부설 극지연구소 빙하지각연구본부

^{*Corresponding Author}

License (open-access, http://creativecommons.org/licenses/by-nc/3.0/):

This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (http://creativecommons.org/licenses/by-nc/3.0/) which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

ABSTRACT

Single-channel seismic exploration has proven effective in delineating subsurface geological structures using small-scale survey systems. The seismic data acquired through zero- or near-offset methods directly capture subsurface features along the vertical axis, facilitating the construction of corresponding seismic sections. However, substantial noise in single-channel seismic data hampers precise interpretation because of the low signal-to-noise ratio. This study introduces a novel approach that integrate noise reduction and signal enhancement via matrix rank optimization to address this issue. Unlike conventional rank-reduction methods, which retain selected singular values to mitigate random noise, our method optimizes the entire singular value spectrum, thus effectively tackling both random and erratic noises commonly found in environments with low signal-to-noise ratio. Additionally, to enhance the horizontal continuity of seismic events and mitigate signal loss during noise reduction, we introduced an adaptive weighting factor computed from the eigenimage of the seismic section. To access the robustness of the proposed method, we conducted numerical experiments using single-channel Sparker seismic data from the Chukchi Plateau in the Arctic Ocean. The results demonstrated that the seismic sections had significantly improved signal-to-noise ratios and minimal signal loss. These advancements hold promise for enhancing single-channel and high-resolution seismic surveys and aiding in the identification of marine development and submarine geological hazards in domestic coastal areas.

Keywords

Single channel seismic survey

Seismic data processing

Denoising

Singular spectrum analysis

단일 채널 탄성파 탐사는 소규모 자료획득 시스템으로 지하 지질구조를 파악하는 효과적인 방법이다. 영벌림거리 혹은 가까운 벌림거리를 사용하여 획득한 단일 채널 탄성파 자료는 연직 방향의 지하 지질구조를 직접 반영하므로 탄성파 단면도를 효과적으로 작성할 수 있다. 그러나 공통중간점 중합 과정을 적용할 수 없어 신호 대 잡음비가 매우 낮으므로 단면에 나타나는 반사 구조의 정밀한 해석에 있어 중합 단면 대비 불리함을 가진다. 본 연구에서는 단일 채널 탄성파 자료의 신호 대 잡음비를 향상시키기 위해 특이 스펙트럼 분석을 기반으로 한 잡음 제거 및 신호 향상 방법을 제안한다. 기존의 특이 스펙트럼 분석 방법은 행렬의 특정 특잇값을 임의로 추출하여 자료 내에 있는 무작위 잡음을 제거하는 방식으로 수행되었으나, 이는 낮은 신호 대 잡음비나 이상 잡음이 있는 자료에 적용할 수 없다. 따라서 본 연구에서는 행렬의 특잇값을 최적화하고 저계수 근사를 수행하여 무작위 및 이상 잡음을 동시에 효과적으로 제거한다. 또한, 잡음 제거로 인한 신호 손실을 보정하고 탄성파 이벤트의 수평적 연속성을 향상시키기 위해 행렬의 고유 영상에 기반한 가중치를 계산하여 탄성파 단면의 품질을 향상시킨다. 본 연구에서 제안하는 기술의 적용성 및 우수성을 확인하기 위해 북극해 척치해저고원에서 획득한 단일 채널 스파커 탄성파 자료에 대한 자료 처리 실험을 수행하였으며, 수치 예제를 통해 매우 높은 수준의 신호 대 잡음비와 최소의 신호 손실을 가진 탄성파 단면을 얻을 수 있었다. 본 연구에서 제안하는 단일 채널 탄성파 자료 처리 기술은 향후 국내 연안지역의 해양개발과 해저 지질재해를 규명하기 위한 단일 채널 및 초고해상도 탄성파 탐사에 매우 효과적으로 기여할 것으로 기대된다.

키워드

단일 채널 탄성파 탐사

탄성파 자료처리

잡음 제거

특이 스펙트럼 분석

MAIN

서 론
이 론
단일 채널 탄성파 자료취득
특이 스펙트럼 분석을 사용한 잡음 제거
견고한 특이 스펙트럼 분석을 사용한 잡음 제거
특이 스펙트럼 분석을 이용한 고유영상 기반 신호 향상
특이 스펙트럼 분석 기반 단일 채널 탄성파 자료처리
수치예제
결 론

서 론

단일 채널 해양 탄성파 탐사는 음파의 왕복주시로부터 해저 지질구조를 규명할 수 있는 탄성파 탐사 기법 중 하나로 천해저 지반조사 및 해양 엔지니어링을 위한 탐사에 사용되어 왔다. 특히 근해에서 실시하는 고해상도 천부 탄성파 탐사는 근해의 환경적 제약으로 인해 단일 채널 방식이 선호되며 천부에 대한 높은 수직 및 수평 분해능을 필요로 하므로 석유 및 가스 등 지하자원 탐사 목적의 방식과는 대비되는 특성을 갖는다. 단일 채널 해양 탄성파 탐사의 주요 장점은 빠르고 간편한 자료 획득이며, 이는 탐사 시간에 대한 제약을 크게 완화시킬 수 있다. 그러나 영 벌림거리나 가까운 벌림거리를 사용한 탐사 기법은 수평적인 변화가 크거나 급경사가 존재하는 지형에 대한 탐지가 어려울 수 있으며 탄성파 속도 분석을 적용할 수 없다는 한계가 있다. 이러한 한계를 극복하기 위해 단일 채널을 사용하여 효과적으로 고해상도 탄성파 단면을 획득하고 해석하기 위한 다양한 연구가 수행되었다. Kim and Kim (2004)은 단일 채널 탄성파 자료의 수평 해상도의 향상을 위한 수평 축 변환 연구를 수행하였으며, Ko et al. (2014)은 단일 채널 탄성파 자료의 속성 분석을 통해 지층 구조 및 내부 구성 물질 파악을 수행하였다. 또한, 속도 분석을 적용할 수 없다는 한계를 극복하기 위해 완전파형역산을 적용하여 지하 속도모델을 규명하는 연구가 수행되었다(Choi et al., 2014; Shin et al., 2017).

지하로부터 반사된 신호가 단일 채널에 기록되는 경우, 이는 지하 지질구조의 형태를 직접 반영하므로 자료에 존재하는 다량의 잡음을 제거하고 자료의 신호 품질을 향상시키는 것은 단일 채널 탄성파 자료를 사용하여 지하 지질구조를 효과적으로 규명하기 위한 매우 중요한 과정이다. 일반적으로, 단일 채널 탄성파 탐사 자료의 잡음 제거 및 신호 향상은 간소한 사용자 정의 작업흐름을 포함한다(Kim et al., 2010; Ko et al., 2014; Kang et al., 2021). 이 작업흐름은 자료 편집, 너울 필터, 정보정, 대역 통과 필터, 디컨볼루션, 이득 회수 등의 단계로 구성된다. 특히, Kluesner et al. (2019)는 단일채널 탄성파 자료의 취득 배열 및 소스 파형이 자료에 미치는 영향을 분석하여 잔차보정과 디컨볼루션 단계를 개선한 작업흐름을 제안하고, 이를 통해 스파커 음원을 사용한 단일채널 탄성파 자료의 품질이 개선될 수 있음을 보였다. 그러나 상대적으로 높은 수준의 잡음이 존재하는 경우, 기존의 자료처리 작업흐름을 통한 신호 향상은 매우 제한적일 수 있으며 단일 채널 탄성파 탐사의 장점을 효과적으로 활용할 수 없게 된다. 단일 채널 탄성파 탐사 자료는 다중 채널 탐사 자료에 비해 상대적으로 자료 밀도가 매우 낮기 때문에 현대의 진보된 잡음 제거 기술은 적용되지 않았으며, 이에 대한 연구가 이루어진 사례는 거의 없다. 따라서 본 연구에서는 단일 채널 탄성파 탐사 자료의 품질을 획기적으로 향상시킬 수 있는 자료처리 방법을 제안한다.

단일 채널 탄성파 탐사 자료의 품질을 저하시키는 잡음은 대표적으로 무작위 잡음과 이상 잡음으로 구분할 수 있다. 무작위 잡음은 자료 내에서 무작위로 분포하며 일관성이 없는 경향을 보이는 잡음으로, 이를 제거하여 자료의 신호 품질을 향상시키는 연구는 매우 다양한 접근 방식을 통해 수행되었다. 중합 전, 중합 후 탄성파 자료의 무작위 잡음 제거를 위한 연구는 예측 기반(prediction-based) 필터(Canales, 1984; Abma and Claerbout, 1995; Soubaras, 1995), 희소성 증진(sparsity-promoting) 필터(Donoho, 1995; Candès et al., 2006; Mallat, 2008; Fomel and Liu, 2010; Lv, 2020), 저계수 근사(low rank approximation) 필터 (Trickett, 2003; Candès and Plan, 2010; Zhou et al., 2010; Oropeza and Sacchi, 2011) 등의 방법이 있으며, 최근 다양한 기계학습 모델들이 제안되어 효과적으로 무작위 잡음을 제거하는 연구가 활발히 수행되고 있다(Si and Yuan, 2018; Kim et al., 2019; Zheng et al., 2020; Liu et al., 2018; Jun et al., 2021). 그러나 이러한 방법들은 주로 가우시안 잡음 모델을 전제로 하고 있어, 이상 잡음과 같은 비 가우시안 잡음이 존재할 경우에는 제대로 작동하지 않는 한계가 있다. 이상 잡음은 갑작스러운 불연속성을 보이는 샘플로, 주변의 자료와 비교하여 진폭의 절대치가 매우 크거나 명확히 구분되는 거동을 보이므로 위의 언급된 일관성 강조(coherence-enforcing) 기반 알고리즘의 성능을 저하시킨다. 이상 잡음을 제거하기 위한 가장 일반적인 방법은 중간값 필터이다. 특히, 견고한(robust) 형태의 중간값 필터는 자료의 방향성과 탄성파 단면의 구조를 고려하여 이상값을 효과적으로 제거한다(Huo et al., 2012; Huang et al., 2017; Chen et al., 2019; Huang et al., 2021). 이상 잡음 제거를 위한 다른 접근 방식으로 견고한 저계수 근사(low-rank approximation) 필터를 고려할 수 있다. Trickett et al. (2012), Chen and Sacchi (2015), Cheng et al. (2015) 등은 견고한 저계수 근사 필터를 제안하여 ℓ₁ 놈 목적함수 및 투영된 경사법(projected gradient method) 등에 기반한 반복 처리를 통해 이상 잡음을 효과적으로 제거하였다. 또한, 특정 제약조건을 만족시키는 반복적 희소성 증진 필터를 통한 이상잡음 제거 연구도 활발히 수행되고 있다(Chen et al., 2014; Gan et al., 2016; Zhao et al., 2018).

전통적인 저계수 근사 필터는 주성분 분석 혹은 고유영상(eigenimage) 필터(Ulrych et al., 1988; Trickett, 2003)의 형태로 무작위 잡음 제거 및 신호 향상 기술로써 매우 활발히 적용되었다. 이후 캐드조우(Cadzow) 필터링 및 특이 스펙트럼 분석(singular spectrum analysis; SSA)의 형태로 응용되어 탄성파 자료처리 분야에 적용한 연구가 이루어졌다(Cadzow, 1988; Trickett, 2008; Sacchi, 2009). 특이 스펙트럼 분석 방법은 신호와 잡음이 각각 특정 그리고 전체 행렬계수에 분포한다는 점에 기반하여, 특잇값 분해(singular value decomposition; SVD)를 통해 탄성파 자료의 낮은 행렬계수만을 추출하는 방식으로 적용된다. 소수의 특정 특잇값을 유지하고 나머지 특잇값을 0으로 대체하는 방식으로 무작위 잡음을 제거하고, ℓ₁ 놈 목적함수 및 투영된 경사법 등을 추가로 적용하여 이상 잡음을 제거한다. 그러나 실제 자료처리 작업흐름에서 사용될 경우 몇 개의 특잇값을 유지할지 결정하는 것은 사용자 정의에 의존하여 자료의 종류에 따라 최적의 변수를 찾기 위한 반복적인 실험이 필요하므로 자료처리 기술의 일관성 및 적용성이 저하될 수 있다.

본 연구에서는 단일 채널 탄성파 자료 내 무작위 잡음과 이상 잡음을 동시에 처리하는 견고한 특이 스펙트럼 분석(robust singular spectrum analysis) 방법을 제안한다. 트레이스 놈(trace norm) 혹은 누클리어 놈(nuclear norm)과 ℓ₁ 놈을 동시에 최소화하는 희소 및 저계수 공동 근사(joint sparse and low rank approximation; JSLR)를 적용하여 자료 내 무작위 잡음과 이상 잡음에 대한 동시 처리를 수행한다. 이는 가속화된 근접 경사법(accelerated proximal gradient; APG)이나 다중계수 대안방법(alternating direction method for multipliers; ADMM)과 특잇값에 대한 수축 및 임계 연산자(shrinkage and thresholding operators)에 기반하여 구성된 방법이다. 본 연구에서 제안하는 방법은 행렬의 특잇값을 선택하는 것이 아닌 최적화하는 방식으로 이루어져 사용자 정의 변수가 요구되지 않으며 자료에 대한 의존성 또한 낮다는 장점이 있다. 추가적으로, 잡음이 제거된 단면의 수평적 연속성을 향상시킬 수 있는 견고한 특이 스펙트럼 분석 기반 신호 향상 방법을 제안한다. 이상 잡음이 제거된 샘플에는 상대적으로 매우 작은 크기의 신호만이 남는 경우가 발생하며, 해당 샘플과 함께 동일한 탄성파 이벤트(event)를 구성하는 주변의 샘플들과 비교했을 때 불연속적인 성질을 보이게 된다. 따라서 본 연구에서는 잡음이 제거된 탄성파 단면에 대한 추가적인 특이 스펙트럼 분석을 수행하여 단면의 주성분을 추출하고, 이를 사용하여 생성한 보정계수를 통해 단면의 수평적 연속성 향상을 수행한다. 제안된 기술을 통해 단일 채널 탄성파 자료의 신호 대 잡음비 향상 뿐만 아니라 효과적인 천부 지질구조 해석이 가능할 것으로 기대되며, 현장의 데스크톱 컴퓨터 혹은 랩톱 컴퓨터에서 사용할 수 있는 가볍고 견고한 형태의 잡음 제거 및 신호 향상 방법을 제안하여 현장 품질관리의 역할까지 수행할 수 있는 단일 채널 탄성파 탐사 시스템의 새로운 요소를 제안하고자 한다.

본 논문은 다음과 같이 구성된다. 먼저 수치 실험에 사용할 단일 채널 탄성파 자료에 대해 간략하게 살펴본 뒤, 1) 견고한 특이 스펙트럼 분석을 사용한 잡음제거, 2) 특이 스펙트럼 분석 기반 신호 향상 알고리즘을 설명한다. 마지막으로, 현장 자료에 대한 적용을 통해 제안된 기술의 타당성과 효용성을 검증한다. 본 연구는 천부의 고해상도 영상을 효율적으로 획득하기 위해 단일 채널 탄성파 자료의 무작위 잡음 및 이상 잡음 제거, 수평적 연속성 향상에 초점을 맞추고 있으므로 자유면 다중반사파 및 지층 다중 반사파 등 특정 패턴이나 구조를 갖는 일관된 잡음에 대해 다루지 않는다.

이 론

단일 채널 탄성파 자료취득

본 연구는 북극해 척치해저고원 지역에서 쇄빙연구선 아라온호에 탑재된 스파커 음원(sparker source)과 단일 채널 수신기로 구성된 해양 탐사 시스템으로 획득한 자료를 사용하여 잡음 제거 및 신호 향상 기술의 적합성 및 타당성을 평가한다. 해당 자료는 5,000 J의 전기에너지 용량을 갖는 SIG Pulse L5 장비를 음원으로 사용하였으며 발파 간격은 3초로 설정하였다. 총 48개의 하이드로폰이 1 m 간격으로 배열되어 2.5초 동안 0.25 ms 간격으로 자료가 기록되며 각각의 신호는 기록 시스템에서 하나의 신호로 중합되어 자료로 저장된다. 자료 획득 시 돌발적으로 나타난 해빙에 의해 일부 구간에서 쇄빙에 의한 잡음과 회피 기동에 의한 잡음이 자료에 포함되었다. 탐사에 대한 세부 사항은 Kang et al. (2021)을 참고하여 확인할 수 있다.

Fig. 1은 취득한 현장 자료(Fig. 1a)와 정보정, 대역 통과 필터링, 디컨볼루션의 자료처리 과정을 거친 최종 탄성파 단면(Fig. 1b)을 나타낸다. Fig. 2a는 Fig. 1b의 특정 부분을 확대한 단면이며, 해당 단면에서는 해저산 구조와 주변지역 하부에 불연속적으로 나타나는 음향공백대 및 위상이 역전된 강한 반사파를 관찰할 수 있다. 특히, 이 반사파 신호는 해저면 모방 반사면(bottom simulating reflection; BSR)으로 추정되며 가스하이드레이트 안정영역(gas hydrate stability zone; GHSZ)과의 비교를 통해 가스하이드레이트의 부존 여부에 대한 평가를 수행할 수 있는 중요한 지표이다(Kang et al., 2021). 그러나 최종 단면과 해당 단면의 주파수-공간파수 스펙트럼(Fig. 2b)으로부터 기존의 자료처리 과정이 적용된 후에도 다량의 잡음이 존재하는 것을 알 수 있으며, 이로 인해 해저면 모방 반사면으로 추정되는 신호에 대한 정확한 해석과 이해가 어렵다. 따라서 본 연구에서 제안하는 방법을 통해 Figs. 1, 2에서 확인할 수 있는 잡음들을 제거하여 단일 채널 탄성파 탐사 자료 해석의 정확성을 향상하고자 한다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F1.jpg

Fig. 1.

Examples of single-channel seismic data showing (a) raw and (b) preprocessed sections, respectively. Note that the time axis is depicted on a relative scale.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F2.jpg

Fig. 2.

(a) An enlarged view of the preprocessed seismic section highlighted by the red-dashed box in Fig. 1b and (b) its corresponding frequency-wavenumber spectrum. Black solid arrows and dashed circles indicate ARAON mounds and bottom simulating reflections (BSRs), respectively.

특이 스펙트럼 분석을 사용한 잡음 제거

본 연구에서 잡음제거 기술은 특이 스펙트럼 분석(singular spectrum analysis; SSA) 혹은 캐드조우(Cadzow) 필터 이론에 기반한다(Cadzow, 1988; Trickett, 2008; Sacchi, 2009). 단일 기울기로 이루어진 2차원 시간-공간 영역 탄성파 자료 $D \in ℝ^{N_{t} \times N_{x}}$ 를 가정하면, 해당 자료의 주파수 영역 표현은 각주파수(angular frequency) 𝜔와 기울기 성분의 식으로 나타낼 수 있다.

(1)

d_{ω} (j) = A_{ω} (j) \exp (i ω p j Δ x), j = \{1, 2, 3, \dots, N_{x}\},

여기서 $i = \sqrt{- 1}$ 이며, $N_{t}$ 와 $N_{x}$ 는 각각 샘플과 트레이스의 개수를 의미한다. $p$ 는 특정 선형 기울기 성분을 뜻하며 $Δ x$ 는 트레이스 간의 공간적 간격을 의미한다. $d$ 와 $A$ 는 각각 주파수 영역 탄성파 자료와 이의 진폭을 의미하며, 아래 첨자 𝜔는 단일 주파수 성분을 나타낸다. 따라서 단일 선형 기울기 성분을 갖는 2차원 주파수-공간 영역의 자료는 단일 주파수 성분에 대하여

(2)

d_{ω} = [d_{ω} (1) d_{ω} (2) d_{ω} (3) \dots d_{ω} (N_{x} - 1) d_{ω} (N_{x})]^{T},

로 표현된다. 특이 스펙트럼 분석의 첫 번째 단계는 주어진 자료의 한켈화(Hankelization)이며, 이는 아래와 같이 한켈 행렬(Hankel matrix) $H_{ω} \in ℂ^{m \times n}$ 로 나타낼 수 있다.

(3)

H_{ω} = [\begin{array}{cccc} d_{ω} (1) & d_{ω} (2) & \dots d_{ω} (N_{x} - m + 1) \\ d_{ω} (2) & d_{ω} (3) & \dots d_{ω} (N_{x} - m + 2) \\ ⋮ & ⋮ & ⋱ & ⋮ \\ d_{ω} (m) & d_{ω} (m + 1) & \dots & d_{ω} (N_{x}) \end{array}],

여기서 $m = N_{x} / 2$ , $n = N_{x} - m + 1$ 로 정의된다. 특이 스펙트럼 분석은 식 (1)에 따라 각 인접한 자료 간의 재귀적 성질을 가정하고 있으므로 식 (3)의 행렬을 이루는 각 요소는

(4)

d_{ω} (j) = P d_{ω} (j - 1),

로 표현할 수 있다. 여기서 $P = \exp (i ω p ∆ x)$ 이다. 따라서 식 (4)를 사용하여 식 (3)은

(5)

H_{ω} = [\begin{array}{cccc} d_{ω} (1) & P d_{ω} (1) & \dots & P^{n - 1} d_{ω} (1) \\ d_{ω} (2) & P d_{ω} (2) & \dots & P^{n - 1} d_{ω} (2) \\ ⋮ & ⋮ & ⋱ & ⋮ \\ d_{ω} (m) & P d_{ω} (m) & \dots & P^{n - 1} d_{ω} (m) \end{array}],

와 같이 재구성될 수 있으며 이는 한켈 행렬 $H_{ω}$ 의 행렬계수(rank)가 1임을 명확하게 보여준다. 식 (1), (2), (3), (4)는 단일 기울기로 이루어진 자료를 가정하였으므로, 이를 일반화하여, r개의 평면파의 중첩으로 이루어진 신호를 고려한다면 한켈 행렬 $H_{ω}$ 는 최대 r의 행렬계수를 갖게 됨을 알 수 있다.

탄성파 자료에 존재하는 잡음은 행렬 $H_{ω}$ 의 모든 요소에 분포하게 되어 행렬의 계수를 증가시킨다. 따라서 특이 스펙트럼 분석은 행렬의 소수의 특잇값만 유지하고 나머지 특잇값을 0으로 대체하는 과정을 통해 유효한 신호 성분은 유지한 채 잡음 성분만 제거하는 저계수 근사 해를 제시한다. 식 (5)는 특잇값 분해를 통해 다음과 같이 쓸 수 있다.

(6)

H_{ω} = U \sum V^{T}

$U \in ℂ^{m \times m}$ 와 $V \in ℂ^{n \times n}$ 는 각각 고유벡터(eigenvector)의 집합으로 구성된 유니터리 행렬(unitary matrix)이며 $\sum \in ℝ^{m \times n}$ 은 특잇값으로 구성된 대각 행렬이다. 식 (6)의 특잇값 행렬의 요소 중 가장 큰 몇 개의 특잇값을 제외한 나머지를 0으로 대체하여 다시 행렬을 구성하면 한켈 행렬 $H_{ω}$ 의 저계수 근사를 통해 자료 내 잡음을 제거할 수 있다.

(7)

H_{ω}^{l} = U \sum_{l} V^{T},

$\sum_{l}$ 는 계수 $l$ 로 근사 된 특잇값 행렬을 의미한다. 그러나 저계수 근사 된 한켈 행렬 $H_{ω}^{l}$ 는 식 (2)의 2차원 주파수-공간 영역의 탄성파 자료와 다른 형태이므로 반 대각 평균화(anti-diagonal averaging)을 이용해 식 (2)와 동일한 형태의 벡터로 변환한다.

(8)

d_{w}^{l} (j) = \{\begin{array}{l} \begin{array}{cc} \begin{array}{cc} \frac{1}{j} \sum_{k = 1}^{j} H_{w}^{l} (j - k + 1, k), & j \leq m \end{array} \end{array} \\ \frac{1}{N_{x} - j + 1} \sum_{k = 1}^{N_{x} - j + 1} H_{w}^{l} (j - k + 1, j - m + k), j > m \end{array},

$H_{ω}^{l} (ㆍ, ㆍ)$ 는 계수 $l$ 로 근사 된 한켈 행렬 $H_{ω}^{l}$ 의 각 행과 열에 해당하는 성분을 의미하며, $d_{ω}^{l}$ 는 계수 $l$ 로 근사 되어 잡음이 제거된 탄성파 자료를 나타낸다. 최종 탄성파 단면, $D_{SSA} \in ℝ^{N t \times N_{x}}$ , 은 식 (8)로부터 재구성한 자료 벡터에 대해 역 푸리에 변환을 적용하여 구할 수 있다.

견고한 특이 스펙트럼 분석을 사용한 잡음 제거

주어진 자료가 소수의 선형 기울기를 가진 신호들의 중첩으로 구성되어 있다면, 식 (7)에서 소수의 특잇값을 선택하는 것은 잘린 특잇값 분해(truncated singular value decomposition; TSVD)에 해당하는 이차 목적함수의 해를 구하는 것과 같으므로 무작위 잡음에 대해 충분한 잡음 제거 효과를 얻을 수 있다. 하지만 해당 방법의 경우 특정 특잇값을 선택하기 위한 사용자 정의가 필요하며 이상값과 같은 비 가우시안 잡음이 동시에 존재할 경우 정상적으로 작동하지 않는 한계가 있다. 이와 같은 한계를 극복하기 위해, 본 연구에서는 희소 및 저계수 공동 근사(joint sparse and low rank approximation; JSLR)에 기반한 견고한 특이 스펙트럼 분석 방법을 제안한다. 희소 및 저계수 공동 근사는 트레이스 놈(trace norm) 혹은 누클리어 놈(nuclear norm)과 ℓ₁ 놈을 동시에 최소화하여 행렬계수와 이상값 구성요소를 최적화하는 잡음 제거 방법을 제안한다(Beck and Teboulle, 2009; Toh and Yun, 2010; Ji et al., 2010; Candès et al., 2011; Tao and Yuan, 2011; Sternfels et al., 2015).

단일 채널 탄성파 탐사로 획득한 자료를 유의미한 신호와 무작위 잡음, 이상 잡음의 중첩으로 가정하면 식 (6)의 한켈 행렬은 다음과 같이 표현할 수 있다.

(9)

H_{ω} = H_{s} + H_{n} + H_{e},

$H_{s} \in ℂ^{m \times n}$ 는 신호를 의미하며, $H_{n} \in ℂ^{m \times n}$ 와 $H_{e} \in ℂ^{m \times n}$ 는 각각 무작위 잡음과 이상 잡음 성분을 의미한다. 또한, 주변 샘플에 비해 상대적으로 높은 진폭을 갖거나 명확히 관련성이 없는 거동을 보이는 이상 잡음 $H_{e}$ 는 매우 적은 수의 0이 아닌 요소를 포함한 희소행렬로 가정한다. 희소 및 저계수 공동 근사를 위한 최적화 문제는 아래의 식을 최소화하는 신호성분 $H_{s} *$ 와 이상잡음 성분 $H_{e} *$ 를 동시에 구하여 해결할 수 있다.

(10)

H_{s}^{*}, H_{e}^{*} = \arg min_{H_{s}, H_{e}} {∥ H_{s} ∥}_{*} + λ {∥ H_{e} ∥}_{1} + \frac{1}{2 μ} {∥ H_{n} ∥}_{F}^{2},

$H_{n} = H - H_{s} - H_{e}$ 이며 $‖ \cdot ‖_{*}$ 는 특잇값의 총합을 측정하는 누클리어 놈, 그리고 $‖ \cdot ‖_{1}$ 과 $‖ \cdot ‖_{F}$ 는 각각 ℓ₁놈과 프로비니우스(Frobenius) 놈을 나타낸다. $λ = 1 / \sqrt{\max (m, n)}$ 와 $μ = (\sqrt{m} + \sqrt{n}) σ$ 는 정규화(regularization) 변수이다(Candès and Plan, 2010; Candès et al., 2011). 본 연구에서는 식 (10)의 목적함수를 해결하기 위하여 ℓ₁놈과 누클리어 놈의 동시 최적화 문제를 해결할 수 있는 수치해석 솔버인 가속화된 근접 경사법(accelerated proximal gradient; APG)을 사용하였다 (Ji et al., 2010). 가속화된 근접 경사법의 자세한 내용은 Algorithm 1과 같다. Algorithm 1의 $S_{τ} (x) = sgn (x) \max (| x | - τ)$ 는 부드러운 수축 연산자(soft shrinkage operator)이며 $D_{τ} (x) = U S_{τ} (\sum) V^{T}$ 는 특잇값 수축 연산자(singular value shrinkage operator)이다(Beck and Teboulle, 2009). 잡음이 제거된 최종 탄성파 단면, $D_{JSLR} \in ℝ^{N t \times N_{x}}$ 은 $H_{s} *$ 에 대한 반 대각 평균화 및 역 푸리에 변환을 적용하여 구할 수 있다.

Fig. 3a는 단일 채널 탄성파 자료(Fig. 1b)의 일부분을 나타낸다. 지하 지질구조를 나타내는 선형적인 신호들과 무작위 잡음, 불규칙하며 주변에 비해 큰 진폭을 보이는 이상 잡음 등을 관찰할 수 있다. Fig. 3b는 해당 자료에 한켈화 및 특잇값 분해를 적용하여 특잇값 행렬(식 (6)의 $\sum \in ℝ^{m \times n}$ )의 대각 성분을 모든 주파수 성분에 대해 도시한 결과이다. Fig. 3c와 3d는 해당 자료에 대한 특이 스펙트럼 분석(계수 2로 가정)을 적용한 결과와 이에 대응하는 특잇값 행렬의 대각 성분을 나타낸다. 기존의 특이 스펙트럼 분석은 무작위 잡음에 대한 효과적인 제거 능력을 보이지만, 신호를 지나치게 선형화 하거나 이상 잡음이 존재하는 부분에서 심한 왜곡이 발생하는 현상을 관찰할 수 있다. Fig. 3e와 3f는 동일한 자료에 대하여 견고한 특이 스펙트럼 분석을 적용한 결과이다. 무작위 잡음과 이상 잡음이 동시에 효과적으로 제거된 것을 관찰할 수 있으며, 행렬의 특잇값은 Algorithm 1에 의해 최적화된 형태로 표현되는 것을 관찰할 수 있다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F3.jpg

Fig. 3.

(Left panels) Single-channel seismic sections and (Right panels) their corresponding singular values with respect to the frequencies. The upper, middle, and lower rows hierarchically display seismic sections and singular values (a and b) for input, (c and d) with denoising applied, and (e and f) with signal enhancement applied, respectively.

Algorithm 1.

Accelerated proximal gradient (APG) method.

Input:	Hankel matrix $H_{ω}$
Output:	Hankel matrix $H_{JSLR}$
1.	Initialize: $H_{s} = H_{e} = 0; i = 0; t_{- 1} = t_{0} = 0; μ_{0} = σ_{0} = t o l =$ STANDARDDEVIATION( $H_{ω}$ );
2.	WHILE $μ_{i} \geq t o l$ DO
3.	$Y_{s}^{i} = H_{s}^{i} + \frac{t_{i - 1} - 1}{t_{i}} (H_{s}^{i} - H_{s}^{i - 1})$ ;
4.	$Y_{e}^{i} = H_{e}^{i} + \frac{t_{i - 1} - 1}{t_{i}} (H_{e}^{i} - H_{e}^{i - 1})$ ;
5.	$G_{s}^{i} = Y_{s}^{i} + \frac{1}{2} (Y_{s}^{i} - Y_{e}^{i} - H_{ω})$ ;
6.	$G_{e}^{i} = Y_{e}^{i} + \frac{1}{2} (Y_{s}^{i} - Y_{e}^{i} - H_{ω})$ ;
7.	$H_{s}^{i + 1} = D_{μ_{i} / 2} (G_{s}^{i})$ ; $H_{e}^{i + 1} = S_{λ μ_{i} / 2} (G_{e}^{i})$ ;
8.	$σ_{i + 1} =$ STANDARDDEVIATION $(H_{ω} - H_{s}^{i + 1} - H_{e}^{i + 1})$ ;
9.	$μ_{i + 1} = (\sqrt{m} + \sqrt{n}) σ_{i + 1}$ ;
10.	$t_{i + 1} = (1 + \sqrt{4 t_{i}^{2} + 1}) / 2$ ;
11.	$i = i + 1$ ;
12.	END WHILE
13.	$H_{JSLR} = H_{s}^{i}$ ;
14.	RETURN;

특이 스펙트럼 분석을 이용한 고유영상 기반 신호 향상

Fig. 3e의 잡음 제거 후의 탄성파 단면을 살펴보면, 주변의 명확한 탄성파 이벤트와 구조적 연속성을 보이지만 잡음 제거 과정을 거치며 신호 누락이 발생한 샘플들을 관찰할 수 있다. 이러한 샘플들은 동일한 탄성파 이벤트를 구성하는 주변 샘플에 비해 상대적으로 작은 진폭을 나타내므로 해당 탄성파 이벤트는 수평적 불연속 구조로 해석될 여지가 있다. 따라서, 단일 채널 탄성파 탐사를 사용하여 천부의 상세한 구조를 규명하는 것은 잡음 제거뿐만 아니라 신호 누락에 대한 후처리를 통해 보다 효과적으로 수행될 수 있다.

Ulrych et al. (1988)은 하나의 탄성파 단면을 여러 개의 고유영상(Eigenimage)의 합으로 나타낼 수 있음을 밝혔다.

(11)

X = \sum_{i = 1}^{R} σ_{i} u_{i} v_{i}^{T},

$X$ 는 잡음이 제거된 영상( $D_{JSLR}$ )을 의미하며 R은 영상의 계수이다. $u_{i}$ 와 $v_{i}$ 는 각각 ${XX}^{T}$ 와 $X^{T} X$ 의 $i$ 번째 고유벡터(eigenvector)를 나타내며, $σ_{i}$ 는 $X$ 의 $i$ 번째 특잇값이다. 식 (11)은 식 (6)과 동일한 형태임을 알 수 있으며, 특이 스펙트럼 분석의 한켈 행렬도 여러 개의 주파수 영역 고유영상의 합으로써 해석할 수 있음을 유추할 수 있다. 견고한 특이 스펙트럼 분석을 통해 비 가우시안 잡음이 제거된 탄성파 단면의 신호(Fig. 3e)는 특이 스펙트럼 분석에 적합한 형태로 변환되었다고 가정할 수 있으며, 따라서 Fig. 3e의 단면은 식 (6) 혹은 식 (11)과 같이 특정 주성분을 포함한 고유영상의 합으로써 근사할 수 있다. Fig. 3e에 대해 가장 큰 몇 개의 고유영상을 추출한 결과는 Fig. 4와 같다. Fig. 4a는 탄성파 단면의 계수를 1로 가정하고 가장 큰 특잇값만을 사용하여 생성한 고유영상으로 Fig. 3e에 포함된 단면의 주된 성분, 즉 대부분의 유의미한 성분들을 매우 잘 표현하고 있다. 동일한 단면(Fig. 3e)에 대해 두 번째, 세 번째, 그리고 네 번째 고유영상을 추출한 결과는 Fig. 4b, 4c 그리고 4d에 나타내었으며 매우 작은 구조들과 잡음만이 포함되어 있는 것을 알 수 있다. 따라서 견고한 특이 스펙트럼 분석을 통해 대부분의 잡음이 제거되었으며 수직적으로 큰 구조적 변화가 없는 단면의 경우, 1~2개의 고유영상으로 단면을 대표하는 주성분을 얻을 수 있다고 판단할 수 있다. 그러나 주어진 탄성파 단면의 가장 주된 성분만을 추출함으로써 탄성파 이벤트의 선형성만을 지나치게 강조하므로 구조의 과도한 선형화가 발생함을 관찰할 수 있다. 또한 주어진 단면에서 신호가 누락된 부분은 추출한 고유영상에도 누락되어 있음을 알 수 있다. 따라서 본 연구에서는 추출된 고유영상을 잡음 제거의 결과물로 사용하는 대신, 신호가 누락되거나 진폭이 매우 약한 부분의 신호 향상을 위한 모델(model)로 사용하여 가중치를 적용하는 방법을 제안한다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F4.jpg

Fig. 4.

Eigenimages generated using (a) the first, (b) the second, (c) the third, and (d) the fourth singular values, respectively.

신호향상을 위해 가중치를 적용하는 기준은 다음과 같다; 모델의 주된 구조를 대표하며 진폭의 크기가 큰 트레이스에는 1에 가까운 가중치를 적용하고, 신호가 누락되어 진폭이 매우 작은 부분에는 상대적으로 큰 값의 가중치를 적용한다. 즉, 각 트레이스의 에너지를 측정하기 위해 평균 제곱근 진폭(root mean square amplitude)을 구하고 정규화(normalization)한 뒤, 이의 역수를 가중치로 사용하면 진폭이 상대적으로 큰 뚜렷한 선형 구조에는 1에 가까운 보정을, 그리고 신호가 누락되어 진폭이 매우 작은 부분에는 상대적으로 큰 보정을 적용하여 신호를 향상시킬 수 있다. 이를 수식으로 나타내면 아래와 같다.

(12)

a = [a_{1}, a_{2}, \dots, a_{k}]^{T}

$a_{j} = {(\frac{1}{n} \sum_{i = 1}^{n} M (i, j)^{2})}^{\frac{1}{2}}$ 이며 $j$ 번째 트레이스의 평균 제곱근 진폭이다. M(i,j)는 잡음 제거 후 단면으로부터 생성된 고유영상 모델의 $j$ 번째 트레이스, $i$ 번째 샘플을 의미한다. $n$ 은 트레이스의 샘플 개수를 의미한다. 식 (12)를 사용하여 모델의 각 트레이스에 최종적으로 적용되는 가중치는 다음과 같다.

(13)

w_{j} = (a_{j} / ‖ a ‖_{\infty})^{- 1},

$w_{j}$ 는 $j$ 번째 트레이스에 적용되는 가중치를 의미하며, 가중치가 적용되어 보정된 모델은 $w_{j}$ 와 $M$ 의 곱으로 구할 수 있다. $‖ \cdot ‖_{\infty}$ 은 벡터를 구성하는 요소의 절댓값 중 최댓값을 의미한다. 한편, 입력된 원본 자료의 특성과 신호가 보정된 모델의 특성을 모두 반영하기 위해 최종 탄성파 단면은 두 자료의 평균으로 나타낼 수 있다.

(14)

\hat{D} (i, j) = \frac{(D_{JSLR} (i, j) + w_{j} M (i, j))}{2},

$\hat{D}$ 은 잡음이 제거된 자료와 가중치가 적용되어 신호가 보정된 새로운 모델의 중합으로 생성한 최종 단면이다. 식 (12), (13), (14)에 따르면, 모델 상에서 상대적으로 진폭이 큰 선형 구조에서는 가중치가 적용된 모델 $(w_{j} \cdot M (i, j))$ 과 입력 자료( $D_{JSLR}$ )가 거의 동일한 값을 가지므로( $w \approx 1$ ) 원본 영상에 선형성이 강조된 모델의 특성이 더해지는 형태로 최종 영상이 구성된다. 또한, 입력자료 상에서 신호의 누락으로 인해 향상이 필요한 부분에는 가중치가 적용된 새로운 모델의 진폭이 더해져 신호 향상의 효과를 얻을 수 있다. Fig. 5는 Fig. 3e에 나타난 잡음 제거 후의 탄성파 단면의 일부분과 이에 대한 첫 번째 고유영상을 추출한 모델, 그리고 계산된 가중치가 적용된 모델을 나타낸다. 각 트레이스에 적용된 가중치는 Fig. 6에 나타내었다. Fig. 7은 단일채널 탄성파 단면, 잡음 제거 후의 단면, 그리고 신호 향상 과정을 거친 단면이며 각 단계를 거치며 잡음 제거와 신호 향상이 매우 효과적으로 적용되었음을 관찰할 수 있다. 본 연구에서 제안하는 신호 향상 방법은 Algorithm 2에 나타내었다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F5.jpg

Fig. 5.

Local windows for (left panel) a portion of the denoised seismic section, its corresponding models (middle panel) illustrating the first eigenimage, and (right panel) with the weighting factor applied, respectively. The local window is extracted from the area delineated by the black-dashed box in Fig. 1b.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F6.jpg

Fig. 6.

The weighting factor generated from the middle panel depicted in Fig. 5.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F7.jpg

Fig. 7.

Local windows for (left panel) a portion of the input seismic section, its corresponding sections (middle panel) with denoising applied and (right panel) with signal enhancement applied, respectively. The local window is extracted from the area delineated by the black-dashed box in Fig. 1b.

Algorithm 2.

An eigenimage-based signal enhancement method.

Input:	Denoised seismic data $D_{JSLR}$ , An eigenimage model $M$
Output:	Signal-enhanced seismic data $\hat{D}$
1.	FOR $j \leftarrow 1, 2, \dots, n_{x}$ DO
2.	FOR $i \leftarrow 1, 2, \dots, n_{t}$ DO
3.	$a (j) = a (j) + M (i, j)^{2}$ ;
4.	END FOR
5.	$a (j) = \frac{1}{n_{t}} (a (j))^{\frac{1}{2}}$ ;
6.	END FOR
7.	$a_{\max} = \max \{a (1), a (2), \dots, a (n_{x})\}$ ;
8.	FOR $j \leftarrow 1, 2, \dots, n_{x}$ DO
9.	$w (j) = {\{a (j) / a_{\max}\}}^{- 1}$ ;
10.	FOR $i \leftarrow 1, 2, \dots, n_{t}$ DO
11.	$\hat{D} (i, j) = (D_{JSLR} (i, j) + w (j) M (i, j)) / 2$ ;
12.	END FOR
13.	END FOR
14.	RETURN;

특이 스펙트럼 분석 기반 단일 채널 탄성파 자료처리

본 연구에서 제안하는 단일 채널 탄성파 자료처리 작업흐름(workflow)은 잡음 제거와 신호 향상의 단계로 구성된다. 잡음 제거 단계는 희소 및 저계수 공동 근사를 통한 견고한 특이 스펙트럼 분석을 수행하여 자료 내의 무작위 및 이상 잡음에 대한 처리를 수행하고(Algorithm 1), 잡음이 제거된 단면에 특이 스펙트럼 분석을 적용하여 고유영상을 추출하고 가중치를 계산 및 적용하는 방식으로 수행된다(Algorithm 2). 추가적으로, 두 단계는 모두 국지 창에 기반하여 각 창에 대해 독립적으로 적용이 되므로 최초의 입력 자료를 국지 창으로 나누는 함수와( $W$ ) 각각의 창을 재구성하여 원본 자료와 같은 형태로 만드는 재구축 함수가( $M$ ) 필요하다. 각 함수는 각각의 국지 창의 경계 부분에서 발생할 수 있는 불연속적 특성 및 자료의 손실을 방지하기 위해 공간 및 시간 축의 방향으로 중첩되는 영역을 갖는 해닝(hanninng) 이동창을 기반으로 설계한다. 여기서 $\sum_{k = 1}^{n_{w}} M_{k} W_{k} ≅ I$ 이며, $n_{w}$ 는 주어진 자료의 크기에 따라 계산된 국지 창 전체의 개수를 의미한다. 본 연구에서 제안하는 행렬계수 최적화 기반 잡음 제거 및 신호 향상 작업흐름은 Algorithm 3에 나타내었다.

Algorithm 3.

The proposed workflow for the single-channel seismic data processing.

Input:	Single-channel seismic data $D$
Output:	Processed seismic data $D_{new}$
1.	FOR $k \leftarrow 1, 2, \dots, n_{w}$ DO
2.	${\tilde{D}}^{k} (x, ω) \leftarrow$ fft $[W_{k} D (x, t)]$ ;
3.	FOR $j \leftarrow 1, 2, \dots, n_{ω}$ DO
4.	$H_{ω}^{k} \leftarrow$ HANKELIZATION $[{\tilde{D}}^{k} (x, j)]$ ;
5.	$H_{ω}^{k} \leftarrow$ ALGORITHM1 $[H_{ω}^{k}]$ ;
6.	$M_{ω}^{k} \leftarrow U \sum_{l} V^{T}$ ; eq. (7) with an user-defined rank $l$
7.	${\tilde{D}}^{k} (x, j)$ ← ANTIDIAGONALAVERAGING $[H_{ω}^{k}]$ ; denoised data
8.	${\tilde{M}}^{k} (x, j)$ ← ANTIDIAGONALAVERAGING $[M_{ω}^{k}]$ ; eigenimage (rank $l$ )
9.	END FOR
10.	$D_{JSLR}^{k} (x, t) \leftarrow$ ifft $[{\tilde{D}}^{k} (x, ω)]$ ;
11.	$M^{k} (x, t) \leftarrow$ ifft $[{\tilde{M}}^{k} (x, ω)]$ ;
12.	${\hat{D}}^{k} (x, t)$ ALGORITHM2 $[D_{JSLR}^{k} (x, t), M^{k} (x, t)]$ ;
13.	END FOR
14.	$D_{new} = \sum_{k = 1}^{n_{w}} M_{k} {\hat{D}}^{k} (x, t)$ ;
15.	RETURN;

수치예제

본 연구에서 제안하는 작업흐름의 타당성을 검증하기 위해 북극 척치해저고원 지역에서 스파커 음원과 단일 채널 수신기로 획득한 자료를 사용하였다. 자료 취득 후 정보정, 대역 통과 필터링, 디컨볼루션의 전처리 과정이 적용된 탄성파 단면을 본 연구의 입력 자료로 사용하였다(Fig 1b). Fig. 8은 입력 자료, 잡음 제거 후의 단면(Algorithm 1), 신호 향상 후의 단면(Algorithm 2)을 나타낸다. 신호 향상 단계는 두 번째 고유영상까지 고려하여 모델을 생성한 뒤 알고리듬 2를 적용하였다. 보다 자세한 비교를 위해 단면을 확대한 결과와 이에 대응하는 주파수-공간파수 스펙트럼을 Fig. 9에 나타내었다. Fig. 9a와 9c의 비교를 통해 매우 효과적인 잡음 제거가 수행되었음을 알 수 있으며 이의 효과는 스펙트럼에서 매우 명확하게 나타난다. Fig. 9b에서 볼 수 있듯 이미 전처리가 적용되었음에도 불구하고 입력 자료의 스펙트럼은 신호의 희소성(sparsity)을 관찰하기 어려우며 잡음의 영향이 대부분의 주파수와 공간파수 영역에 분산되어 나타나고 있다. 반면, Fig. 9d의 스펙트럼은 잡음의 영향이 대부분 제거되어 스펙트럼의 희소성이 향상되었음을 관찰할 수 있다. 신호 향상 단계를 거친 단면(Fig. 9e)은 탄성파 이벤트의 수평적 연속성 및 신호 누락에 대한 보정에 있어 보다 향상된 결과를 보이며, 스펙트럼 또한 향상된 희소성과 일관성을 보이는 것을 관찰할 수 있다. 정량적 분석을 위해 신호 대 잡음비를 아래의 식을 사용하여 구하였다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F8.jpg

Fig. 8.

Single-channel seismic data for (a) the input (preprocessed as denoted in the main text), (b) with denoising applied, and (c) with denoising applied followed by signal enhancement, respectively.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F9.jpg

Fig. 9.

(Left panels) Enlarged versions of the single-channel seismic sections shown in Fig. 8 (highlighted by the red-dashed box) and (Right panels) their corresponding frequency-wavenumber spectra. The upper, middle, and lower rows hierarchically display seismic sections (a and b) for the input, (c and d) with denoising applied, and (e and f) with denoising applied followed by the signal enhancement application, respectively.

(15)

SNR = \frac{{\hat{D}}_{new}}{\hat{D} - {\hat{D}}_{new}} .

식 (15)를 통해 각 단면의 신호 대 잡음비를 구한 결과는 Fig. 10에 나타내었다. 대부분의 유의미한 신호가 분포해 있는 100~300 Hz의 영역에서 잡음 제거(Algorithm 1) 및 신호 향상(Algorithm 2) 단계를 거치며 수 데시벨(decibel; dB)의 신호 대 잡음비 향상을 관찰할 수 있으며, 특히 신호 향상 단계를 통해 보다 많은 정보를 복구하고 있다는 것을 알 수 있다. 또한, 신호가 거의 존재하지 않으며 상대적으로 잡음의 비율이 높은 100 Hz 이하의 대역에서는 5 데시벨 이상의 잡음이 제거되고 있는 것을 알 수 있다. 추가적으로, 국지 유사성(local similarity; Fomel, 2007)을 계산하여 신호의 누출(signal leakage)에 대한 정량 평가를 수행하였다. 입력자료(Fig. 8a)와 자료처리 최종 결과(Fig. 8c)의 잔차를 나타내는 Fig. 11a를 살펴보면, 대부분 잡음으로 구성되어 있으나 구조적인 경향을 띠는 부분이 관찰되므로 잡음 제거 및 신호 향상 단계에서 신호의 누출이 발생했을 가능성이 있다. Fig. 11b는 자료처리 최종 결과(Fig. 8c)와 잔차(Fig. 11a)와의 국지 유사성을 나타낸다. 국지 유사성의 정량적인 값들을 통해 유의미한 신호로 추정되는 최종 단면과 이에 대응하는 잡음 단면과의 관계성이 매우 낮음을 관찰할 수 있으며, 따라서 본 연구에서 제안하는 자료처리 방법은 신호와 잡음을 매우 잘 분리하였으며 신호의 누출이 거의 존재하지 않음을 알 수 있다. 또한, Fig. 11a의 구조적인 경향을 띠는 부분은 신호에 존재하는 산란(scattering) 형태의 비 선형적 이벤트들이 구조를 따라 제거된 결과로 판단된다. 결과적으로, 본 연구에서 제안하는 작업 흐름을 통해 단일 채널 탄성파 자료의 잡음 제거 및 신호 향상이 매우 성공적으로 수행되었음을 알 수 있다.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F10.jpg

Fig. 10.

A comparison of the signal-to-noise ratio (SNR). The red solid line indicates the SNR of the input data, whereas the blue and black solid lines represent the SNRs with denoising applied and with denoising applied followed by signal enhancement, respectively.

https://cdn.apub.kr/journalsite/sites/ggekseg/2024-027-02/N0540270201/images/ggekseg_27_02_01_F11.jpg

Fig. 11.

(a) The residual section between the input (Fig. 8a) and final output (Fig. 8c) sections and (b) a local similarity map between the residual (Fig. 11a) and final output (Fig. 8c) sections.

결 론

본 연구는 특이 스펙트럼 분석 기반의 단일채널 탄성파 자료처리 작업 흐름을 제안한다. 먼저, 선행 연구에서 수행한 자료처리 적용 후에도 존재하는 강한 잡음들을 제거하기 위해 견고한 특이 스펙트럼 분석을 적용하여 잡음을 제거하였다. 견고한 특이 스펙트럼 분석은 트레이스 놈(trace norm) 혹은 누클리어 놈(nuclear norm)과 ℓ₁ 놈을 동시에 최소화하여 자료 내 계수와 이상값 구성요소를 최적화하는 희소 및 저계수 공동 근사(joint sparse and low rank approximation; JSLR)에 기반하며, 가속화된 근접 경사법(accelerated proximal gradient method)를 사용하여 최적화 문제를 해결하였다. 이를 통해 무작위 잡음뿐만 아니라 비 가우시안 잡음에 해당하는 이상 잡음에 대한 제거를 수행할 수 있으므로 일반 탄성파 자료와 비교하여 신호 대 잡음비(signal-to-noise ratio; SNR)가 낮은 단일채널 탄성파 자료에 대해 매우 효과적인 잡음 제거를 수행할 수 있었다. 추가적으로, 특이 스펙트럼 분석을 사용하여 단일 채널 탄성파 단면의 고유영상을 추출하고, 이를 통해 지하 지질구조의 수평적 불연속성 및 잡음 제거 후 발생할 수 있는 신호 누락에 대한 보정을 수행하였다. 수치 예제를 통해 본 연구에서 제안하는 잡음 제거 및 신호 향상 방법의 견고함을 확인할 수 있었으며, 주파수-공간파수 스펙트럼, 신호 대 잡음비 및 국지 유사성 등의 정성적 및 정량적 분석을 통해 제안하는 작업흐름의 타당성을 확보할 수 있었다.

후속 연구로 잡음 제거 및 신호 향상이 적용된 단면의 해저면 모방 반사면(bottom simulating reflection; BSR) 규명 및 가스하이드레이트 안정영역(gas hydrate stability zone; GHSZ)과의 비교를 통해 가스하이드레이트의 부존 여부에 대한 평가가 수행된다면 북극 척치해저고원 남서부 해저사면의 가스하이드레이트 분포에 대한 보다 정확한 지표를 제시할 수 있을 것으로 예상된다. 또한, 본 논문에서 제안하는 기술을 원시자료에 적용하고 후속 공정의 결과 향상에 미치는 영향을 평가하여, 단일 채널 탄성파 자료처리의 성능을 제고할 수 있는 최적의 작업흐름을 도출할 수 있을 것으로 예상한다. 수치 예제에 사용된 자료뿐만 아니라, 본 연구에서 제안된 기술은 일반적인 단일 채널 탄성파 자료처리 시스템의 새로운 요소로서 사용될 수 있을 것으로 예상되어 향후 해양 단일 채널 탄성파 자료의 신호 대 잡음비 향상에 기여할 수 있을 것으로 기대된다.

Acknowledgements

이 논문은 2024년도 해양수산부 재원으로 해양수산과학기술진흥원의 지원을 받아 수행된 연구입니다(RS-2023-00259686). 또한 2024년도 정부(교육부)의 재원으로 한국연구재단의 지원을 받은 기초연구사업(No.2019R1A6A1A03033167) 및 2023년도 강원대학교 대학회계 학술연구 조성비의 지원을 받아 연구한 결과를 바탕으로 작성되었습니다. 저자 중 이찬희는 정부(산업통상자원부)의 재원으로 한국원자력환경공단의 지원을 받은 고준위방폐물관리 학계 전문인력양성 대학지원 사업(과제번호: 202308470002)의 지원을 받았습니다. 저자 중 강승구는 한국해양과학기술원 부설 극지연구소에서 해양수산부의 재원으로 수행하는 극지 해양환경 및 해저조사 연구 사업(R&D) ‘북극해 해저지질 조사 및 해저환경 변화 연구(과제번호: 20210632)’의 지원을 받았습니다.

References

Abma, R., and Claerbout, J., 1995, Lateral prediction for noise attenuation by t-x and f-x techniques, Geophysics, 60, 1887-1896.

10.1190/1.1443920

Beck, A., and Teboulle, M., 2009, A fast iterative shrinkage-thresholding algorithm for linear inverse problems, SIAM Journal of Imaging Sciences, 2, 183-202.

10.1137/080716542

Cadzow, J. A., 1988, Signal enhancement - a composite property mapping algorithm, IEEE Transactions on Acoustics, Speech and Signal Processing, 36, 49-62.

10.1109/29.1488

Canales, L., 1984, Random noise reduction, SEG Technical Program Expanded Abstracts 1984, 525-527.

10.1190/1.1894168

Candès, E. J., and Plan, Y., 2010, Matrix completion with noise, Proceedings of the IEEE, 98, 925-936.

10.1109/JPROC.2009.2035722

Candès, E. J., Romberg, J. K., and Tao, T., 2006, Stable signal recovery from incomplete and inaccurate measurements, Communications on Pure and Applied Mathematics, 59, 1207-1223.

10.1002/cpa.20124

Candès, E. J., Li, X., Ma, Y., and Wright, J., 2011, Robust principal component analysis?, Journal of Association for Computing Machinery, 58, 11:1-11:37.

10.1145/1970392.1970395

Chen, K., and Sacchi, M. D., 2015, Robust reduced-rank filtering for erratic seismic noise attenuation, Geophysics, 80, V1-V11.

10.1190/geo2014-0116.1

Chen, Y., Fomel, S., and Hu, J., 2014, Iterative deblending of simultaneous-source seismic data using seislet-domain shaping regularization, Geophysics, 79, V179-V189.

10.1190/geo2013-0449.1

Chen, Y., Zu, S.,Wang, T., and Chen, X., 2019, Deblending of simultaneous source data using a structure-oriented space-varying median filter, Geophysical Journal International, 216, 1214-1232.

10.1093/gji/ggy487

Cheng, J., Chen, K., and Sacchi, M. D., 2015, Application of robust principal component analysis (RPCA) to suppress erratic noise in seismic records, SEG Technical Program Expanded Abstracts2015, 4646-4651.

10.1190/segam2015-5869427.126642941

Choi, Y. J., Shin, S. R., Ha, J. H., Chung, W. K., and Kim, W. S., 2014, Velocity model building using waveform inversion from single channel engineering seismic survey, Geophysics and Geophysical Exploration, 17(4), pp. 231-241.

10.7582/GGE.2014.17.4.231

Donoho, D. L., 1995, De-noising by soft-thresholding, IEEE Transactions on Information Theory, 41, 613-627.

10.1109/18.382009

Fomel, S., 2007, Local seismic attributes, Geophysics, 72, A29-A33.

10.1190/1.2437573

Fomel, S., and Liu, Y., 2010, Seislet transform and seislet frame, Geophysics, 75, V25-V38.

10.1190/1.3380591

Gan, S.,Wang, S., Chen, Y., and Chen, X., 2016, Simultaneous-source separation using iterative seislet-frame thresholding, IEEE Geoscience and Remote Sensing Letters, 13, 197-201.

10.1109/LGRS.2015.2505319

Huang, G., Bai, M., Zhao, Q., Chen, W., and Chen, Y., 2021, Erratic noise suppression using iterative structure-oriented space-varying median filtering with sparsity constraint, Geophysical Prospecting, 69, 101-121.

10.1111/1365-2478.13032

Huang, W., Wang, R., Gong,X., and Chen, Y., 2017, Iterative deblending of simultaneous-source seismic data with structuring median constraint, IEEE Geoscience and Remote Sensing Letters, 15, 58-62.

10.1109/LGRS.2017.2772857

Huo, S., Luo, Y., and Kelamis, P. G., 2012, Simultaneous source separation via multidirectional vector-median filtering, Geophysics, 77, V123-V131.

10.1190/geo2011-0254.1

Ji, H., Huang, S., Shen, Z., and Xu, Y., 2010, Robust video restoration by joint sparse and low rank approximation, SIAM Journal on Imaging Sciences, 4, 1122-1142.

10.1137/100817206

Jun, H., Kim, C.-H., and Kim, H.-J., 2021, Machine-learning based noise attenuation of field seismic data using noise data acquisition, Journal of the Korean Society of Mineral and Energy Resources Engineers, 58(5), 408-417.

10.32390/ksmer.2021.58.5.408

Kang, S.-G., Jang, U., Kim, S., Choi, Y., Kim, Y.-G., Hong, J. K., and Jin, Y. K., 2021, Exploration of the gas hydrates on the southwestern continental slope of the Chuckchi plateau in the Arctic ocean, Journal of the Korean Society of Mineral and Energy Resources Engineers, 58(5), 418-432.

10.32390/ksmer.2021.58.5.418

Kim, H. D., and Kim, J. H., 2004, Horizontal distance correction of single channel marine seismic data, Geophysics and Geophysical Exploration, 7(4), 245-250.

Kim, W. S., Park, K., Kim, H., Cheong, S., Koo, N., Lee, H., and Park, E., 2010, Detailed processing and analysis on the single-channel seismic data for site survey of Caecheon-Wonsando subsea tunnel, Geophysics and Geophysical Exploration, 13(4), 336-348.

Kim, Y., Hardisty, R., and Marfurt, K., 2019, Seismic random noise attenuation in f-x domain using complex-valued residual convolutional neural network, SEG Technical ProgramExpanded Abstracts 2019, Society of Exploration Geophysicists, 2579-2583.

10.1190/segam2019-3216543.1

Kluesner, J., Brothers, D., Hart, P., Miller, N., and Hatcher, G., 2019, Practical approaches to maximizing the resolution of sparker seismic reflection data, Marine Geophysical Research, 40(3), 279-301.

10.1007/s11001-018-9367-2

Ko, H. K., Ha, J. H., Chung, W. K., and Shin, S. R., 2014, Seismic attribute analysis and interpretation with single channel seismic data of shallow marine, Journal of the Korean Society of Mineral and Energy Resources Engineers, 51(1), 55-67.

10.12972/ksmer.2014.51.1.55

Liu, D., Wang, W., Chen, W., Wang, X., Zhou, Y., and Shi, Z., 2018, Random noise suppression in seismic data: What can deep learning do?, SEG Technical Program Expanded Abstracts 2018, Society of Exploration Geophysicists, 2016-2020.

10.1190/segam2018-2998114.1

Lv, H., 2020, Single-channel and multi-channel seismic random noise suppression based on the regularized non-stationary decomposition, Journal of Applied Geophysics, 175, 103986.

10.1016/j.jappgeo.2020.103986

Mallat, S., 2008, A Wavelet Tour of Signal Processing: The Sparse Way, San Diego: Academic.

Oropeza, V., and Sacchi, M. D., 2011, Simultaneous seismic data denoising and reconstruction, Geophysics, 76, V25-V32.

10.1190/1.3552706

Sacchi, M. D., 2009, FX singular spectrum analysis, CSPG CSEG CWLS Convention, 392-395.

Shin, S. R., Chung, W. K., Ha, J. H., and Choi, Y. J., 2017, Seismic waveform inversion method and apparatus for underground velocity information derived from the single-channel seismic data, Korean Intellectual Property Office, KR101695080B1.

Si, X., and Yuan, Y., 2018, Random noise attenuation based on residual learning of deep convolutional neural network, SEG Technical Program Expanded Abstracts2018, Society of Exploration Geophysicists, 1986-1990.

10.1190/segam2018-2985176.1

Soubaras, R, 1995, Prestack random and impulsive noise attenuation by f-x projection filtering, SEG Technical Program Expanded Abstracts 1995, 525-527.

10.1190/1.1887383

Sternfels, R., Viguier, G., Gondoin, R., and Meur, D. L., 2015, Multidimensional simultaneous random plus erratic noise attenuation and interpolation for seismic data by joint low-rank and sparse inversion, Geophysics, 80, WD129-WD141.

10.1190/geo2015-0066.1

Tao, M., and Yuan, X., 2011, Recovering low-rank and sparse components of matrices from incomplete and noisy observations, SIAM Journal of Optimization, 21, 57-81.

10.1137/100781894

Toh, K. C. and Yun, S., 2010, An accelerated proximal gradient algorithm for nuclear norm regularized least squares problems, Pacific Journal of Optimization, 6, 615-640.

Trickett, S. R., 2003, F-xy eigenimage noise suppression, Geophysics, 68, 751-759.

10.1190/1.1567245

Trickett, S. R., 2008, F-xy Cadzow noise suppression: eigenimage processing of seismic sections, SEG Technical Program Expanded Abstracts 2008, 2586-2590.

10.1190/1.3063880

Trickett, S. R., Burroughs, L., and Milton, A., 2012, Robust rankreduction filtering for erratic noise, SEG Technical Program Expanded Abstracts 2012, 1-5.

10.1190/segam2012-0129.1

Ulrych, T. J., Sacchi, M. D., and Freire, S. L. M., 1988, Eigenimage processing of seismic sections, SEG Technical Program Expanded Abstracts 1988, 1261-1265.

10.1190/1.1892508

Zhao, Q., Du, Q., Gong, X., and Chen, Y., 2018, Signal-preserving erratic noise attenuation via iterative robust sparsity-promoting filter, IEEE Geoscience and Remote Sensing Letters, 56, 3547-3560.

10.1109/TGRS.2018.2802462

Zheng, Y., Yuan, Y., and Si, X., 2020, The improved DnCNN linear noise attenuation, SEG 2019 Workshop: Mathematical Geophysics: Traditional vs Learning, Beijing, China, Global Meeting Abstracts, 56-59.

10.1190/iwmg2019_14.1

Zhou, Z., Li, X., Wright, J., Candès, E. J., and Ma, Y., 2010, Stable principal component pursuit, 2010 IEEE International symposium on Information Theory, 1518-1522.

10.1109/ISIT.2010.5513535

Geophysics and Geophysical Exploration ISSN:1229-1064(Print) 2384-051X(Online) 지구물리와 물리탐사

Preview

Single-Channel Seismic Data Processing via Singular Spectrum Analysis

ABSTRACT

MAIN

Fig. 1.

Examples of single-channel seismic data showing (a) raw and (b) preprocessed sections, respectively. Note that the time axis is depicted on a relative scale.

Fig. 2.

(a) An enlarged view of the preprocessed seismic section highlighted by the red-dashed box in Fig. 1b and (b) its corresponding frequency-wavenumber spectrum. Black solid arrows and dashed circles indicate ARAON mounds and bottom simulating reflections (BSRs), respectively.

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

Fig. 3.

Algorithm 1.

Accelerated proximal gradient (APG) method.

(11)

Fig. 4.

Eigenimages generated using (a) the first, (b) the second, (c) the third, and (d) the fourth singular values, respectively.

(12)

(13)

(14)

Fig. 5.

Fig. 6.

The weighting factor generated from the middle panel depicted in Fig. 5.

Fig. 7.

Local windows for (left panel) a portion of the input seismic section, its corresponding sections (middle panel) with denoising applied and (right panel) with signal enhancement applied, respectively. The local window is extracted from the area delineated by the black-dashed box in Fig. 1b.

Algorithm 2.

An eigenimage-based signal enhancement method.

Algorithm 3.

The proposed workflow for the single-channel seismic data processing.

Fig. 8.

Single-channel seismic data for (a) the input (preprocessed as denoted in the main text), (b) with denoising applied, and (c) with denoising applied followed by signal enhancement, respectively.

Fig. 9.

(15)

Fig. 10.

A comparison of the signal-to-noise ratio (SNR). The red solid line indicates the SNR of the input data, whereas the blue and black solid lines represent the SNRs with denoising applied and with denoising applied followed by signal enhancement, respectively.

Fig. 11.

(a) The residual section between the input (Fig. 8a) and final output (Fig. 8c) sections and (b) a local similarity map between the residual (Fig. 11a) and final output (Fig. 8c) sections.

Acknowledgements

References