Original Article

Journal of Korean Society of Disaster and Security. 30 June 2019. 73-82
https://doi.org/10.21729/ksds.2019.12.2.73

ABSTRACT


MAIN

  • 1. 서 론

  • 2. 본 론

  •   2.1 연구 방법론

  •   2.2 모형 적용

  • 3. 모의 결과 및 분석

  • 4. 결 론

1. 서 론

최근 이상기후로 인하여 국지성 집중호우의 발생빈도가 증가함에 따라 하천수위가 급격하게 상승하고 이로 인한 하천변 사회기반시설의 침수피해가 증가하고 있다(Korea Meteorological Administration, 2017; Seoul Metropolitan Government, 2016). 따라서 사회기반시설의 침수 가능성 여부에 대한 신속한 예·경보가 필요한 실정이다. 국내의 경우, 「하천법시행규칙 제23조 제2항」에 의하면 홍수 예·경보에 하천수위를 이용하고 있으며 주로 Navier-Stokes 방정식 등 다양한 미분방정식을 해석하는 수치모형을 통하여 하천수위를 예측하였다. Lee et al.(2005)은 1차원 모형DWOPER을 적용하여 하상변화와 지류 및 본류의 유량과 조위의 변화에 따른 한강의 수위변화를 분석하였으며, Lee and Lee(2010)은 1차원 부정류모형인 FLDWAV를 이용하여 조석 및 팔당댐 방류에 따른 수위변화를 분석하였다. 또한 Song et al.(2014)은 유량에 따른 배수(Backwater)흐름을 해석하기 위하여 천수방정식을 해석한 수치모형을 개발하고 한강본류구간(영동대교~행주대교)에서 유량을 변화시켜가며 수위를 검토한 바 있다. 그러나 CFD에 의한 수치해석은 정확하고, 정교한 수위예측 결과를 제공하는 반면, 수치해석시간이 다소 오래 소요됨에 따라 홍수 예·경보 시 필요한 골든타임 확보에 한계점이 있다. 따라서 기존의 수치모형의 한계점을 극복하면서, 예측의 정확도를 확보하는 새로운 방법의 적용이 필요하다.

이에 최근에는 방대한 양의 자료 수집이 가능해지고, 빅데이터 처리 기술이 발달함에 따라 인공신경망(Artificial Neural Network)등을 적용한 자료기반의 수위예측 모델이 많은 연구에서 사용되었다. Thirumalaiah and Deo(1998)은 신경망 모형의 3가지 알고리즘을 이용하여 홍수기 시 홍수위를 예측하고 분석을 통하여 종속상관알고리즘이 우수한 예측결과를 도출함을 확인하였으며, Chen et al.(2012), Kim and Tachikawa(2017)은 일반하천에서 ANN 기반의 예측 모형을 구축하여 하천의 수위예측을 수행하고 적용성을 검토하였다. 이외에도 Coulibaly and Anctil(1999)은 효율적 저수지운영을 위한 단기 유입량 예측을 위해 순환신경망(Recurrent Neural Network, RNN)을 적용하였으며 RNN이 효과적인 실시간 예측결과를 도출하는 것을 확인하였고, Kim et al.(2018)은 도시 홍수 예측 기법을 제시하기 위하여 동적 인공신경망을 통하여 최소의 입력 자료만을 활용하여 맨홀의 홍수유출량을 예측하였다. 그러나, 수문 및 수자원분야에서 다루는 자료가 대부분 시계열 형태의 자료이고 예측 시 입력자료와 결과자료의 시간지연(Time delay)이 존재함에도 불구하고, 다양한 시간적 매개변수를 고려한 인공신경망 모형을 사용한 연구는 많지 않은 실정이다. 시간적 매개변수는 예측 결과의 정확도에 영향을 미칠수 있기때문에, 이에 시간적 매개변수를 고려할 수 있는 인공신경망을 이용한 하천의 수위예측 연구가 필요하다.

따라서 본 연구에서는 시간 지연(Time delay)의 매개변수를 고려한 비선형 자기회귀 신경망 모형(Nonlinear Auto- Regressive eXogenous model, NARX)을 사용하였다. 또한 다양한 수문 입력 자료를 사용하기 보다는 강수량 및 하천의 수위만을 사용하여 예측 모형을 구축하였다. 수위예측 시 NARX 모형의 적합성을 판단하기 위하여 인공신경망(Artificial Neural Network, ANN), 순환신경망(Recurrent Neural Network, RNN) NARX model을 이용하여 한강대교의 3시간 후 수위를 예측한 후 각각의 모형의 정확도를 비교하였다.

2. 본 론

2.1 연구 방법론

2.1.1 NARX model

본 연구에서는 하천의 수위예측을 위한 인공신경망으로서 외생변수를 활용한 비선형 자기회귀 신경망 모형(Nonlinear Auto-Regressive eXogenous model, NARX)을 사용하였다. 기존의 인공신경망(Artificial Neural Network) 모형은 Fig. 1과 같이, 입력층(Input Layer)에 자료를 입력하게 되면 입력 값과 가중치(Weights)를 곱한 값을 전부 합산한 후 활성 함수(Activation Function)를 통하여 최종적으로 결과 값이 출력되는 전방 전달 신경망(Feedforward Neural Network)구조로 구성되어 있다.

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F1.jpg
Fig. 1.

Structure of ANN model

반면에, NARX 모형은 일반적인 인공신경망 구조와 달리 학습 시 도출되는 결과 출력 값이 다시 입력층의 입력 값으로 들어가 학습에 이용되는 피드백 구조를 가지고 있는 특징이 있어 자기상관성을 가지는 시계열(Time series)형태의 자료를 예측하는데 적합하다(Shen et al., 2013). NARX 모형의 구조는 Fig. 2와 같으며, Eq. (1)과 같이 나타낼 수 있다.

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F2.jpg
Fig. 2.

Structure of NARX model

$$o_t=F(o_{t-1},o_{t-2},\;...\;,o_{t-d},x_t,x_{t-1},x_{t-2},x_{t-3},\;...\;)+\varepsilon_t$$ (1)

여기서, ot는 시점 t에서의 시스템의 출력 값이고 xt는 시점 t에서의 시스템의 입력 값이고 εt는 오차의 항이다. d는 피드백 시킬 과거 출력 값의 개수인데, 이는 시계열의 자기상관성의 형태에 따라서 민감도 분석(Sensitivity Analysis)을 통해 최적의 값을 선택한다.

2.1.2 정확도 분석 기법

본 연구에서는 구축된 입력 자료를 ANN, RNN 및 NARX 모형을 통하여 학습시킨 후 예측 출력 결과 자료에 대한 통계분석을 통하여 모형의 정확도를 비교하고자 한다. 모형의 정확도 판단 지표로는 평균제곱근오차(Root Mean Square Error, RMSE)와 평균절대오차(Mean Absolute Error, MAE)를 사용하였으며, 각각의 오차 산정 식은 Eqs. (2), (3)과 같다. RMSE와 MAE 모두 실제 관측 값과 예측 결과의 오차를 정량화하여 나타낸 지표로 값이 ‘0’에 가까울수록 실측값과 유사하게 재현된다는 것을 판단할 수 있다. 본 연구에서는 하천수위 예측 시 ‘m’ 단위를 사용하였다.

$$RMSE=\sqrt{\frac{{\displaystyle\sum_{i=1}^n}(x_{p,i}-x_{o,i})^2}n}$$ (2)
$$MAE=\frac{{\displaystyle\sum_{i=1}^n}\left|x_{p,i}-x_{o,i}\right|}n$$ (3)

여기서, n은 데이터 총 개수, xp,i, xo,i는 모형의 i번째 예측값 및 실제 관측값을 의미한다.

정량적인 오차분석 외에도 결정계수(R2)에 대한 분석도 실시하였다. 결정계수는 상관계수의 제곱으로 나타내며 0과 1사이의 값을 갖는다. 결정계수 값이 1에 가까울수록 관측 값과 모형의 예측 결과 값이 선형관계를 갖는다고 판단할 수 있으며, 이를 통하여 예측 값과 실측값의 경향성을 확인할 수 있다.일반적으로 수문모형에서는 결정계수가 0.5보다 크면 유의미한 상관관계가 있다고 간주한다(Moriasi et al., 2007).

2.2 모형 적용

2.2.1 입력 자료간의 상관성 분석

하천변 사회기반시설 침수에 가장 큰 영향을 미치는 하천수위는 조위, 강수량, 하천 유량 등 다양한 인자들에 의해 영향을 받을 수 있으며, 본 연구에서는 다양한 영향인자 중 한강대교 수위와 직접적으로 관련이 있는 인자는 인공신경망 모델 학습 시 노이즈를 발생 시킬수 있다고 판단하였다. 따라서 노이즈를 발생시킬 수 있는 하천 유량 자료는 제외하고 조위 자료와 강수량 자료에 대하여 상관성 분석을 수행하였다.

위 그래프는 2011년 7월의 조위와 한강대교의 수위를 도시한 것이다. 해당 시기는 중부권에 집중호우가 있었던 시기로, Fig. 3을 통하여 수위와 조위가 진동하고 있는 양상을 확인하였으며 정량적인 주파수 분석을 위하여 이산푸리에 변환(Discrete Fourier Transform, DFT)을 수행하였다.

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F3.jpg
Fig. 3.

The tidal level at Incheon and water surface elevation at Hangang RIver Bridge in July 2011

이산푸리에변환 결과 Fig. 4을 통하여 조위 자료는 8.6 cycle/month의 주파수에서 첨두(Peak)값이 발견되었다. 수위자료와 비교하였을 때 조위자료는 비교적 뚜렷한 첨두값을 보여주었는데, 이는 조위는 조석현상에 지배적으로 영향을 받기 때문에 일정한 주기가 존재한다고 판단된다. 그러나 조위자료와 달리 수위자료에서는 매우 낮은 주파수(2cycle/month)이하에서 다양한 첨두값이 관찰되었고, 강수량 자료 역시 대부분의 주파수에서 다양한 첨두값을 확인할 수 있었다. 이는 해당시기에 댐의 방류, 집중호우, 지류나 상류로부터의 유입량 증가 등 비주기적인 현상들이 수위변화에 더욱 큰 영향을 주는 것임을 알 수 있다. 추가적으로 두 인자를 평균0, 분산1로 표준화 후 Pearson’s 상관계수(correlation coefficient)를 구하여 정량적으로 상관관계를 검토하였다. 그 결과 수위와 조위인자의 Pearson 상관계수는 0.023로 두 인자는 상관관계가 없다는 것을 확인할 수 있었다. 반면, 강수량인자의 경우 수위와 Pearson 상관계수가 0.25로 약한 양의 상관관계를 가지는 것을 확인할 수 있어 최종적으로 수위예측을 위하여 강수량 자료를 활용하기로 하였다. 그러나 홍수피해의 조기예측을 위해서는 일강수량이 아닌 시강수량을 사용해야하고 이는 매끄러운 변화가 아닌 거친 변화를 보인다. 따라서 이를 그대로 사용해서 미분방정식 기반 모델에 의거한 예측에 심각한 오류를 만들 수 있다. 보통 제어시스템이나 통계분야에선 이를 필터링하여서 사용하는데, 본 연구에서는 3시간동안의 누적강우량을 사용하는 방식으로 필터링하였다.

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F4.jpg
Fig. 4.

Discrete fourier transform of water surface elevation, tidal level and precipitation in July 2011

2.2.2 입력 자료 구축

본 연구에서 NARX model을 통하여 수위예측 시 학습에 사용된 자료는 과거 10년간(2009~2018년) 한강대교의 수위자료(관측지점: 경도 126° 57' 24'', 위도 37° 30' 52'')와 서울시의 강수량자료(관측지점: 경도 126° 57' 57.6'', 위도 37° 34' 15.6'')로, 수위자료는 한강홍수통제소에서, 강수량자료는 기상청의 기상정보공개포털에서 확보하였다. 각각의 자료에는 결측치(Missing value)가 존재하였는데, 수위자료의 결측치는 주로 동절기에 있었다. 이 시기는 수위의 변동이 조위변화밖에 없어 선형보간(Linear Interpolation)을 수행하였다. 강수량 자료의 결측치는 비강우시기에 존재하여 이는 ‘0’으로 간주하였다. 각각의 자료에 대하여 출처 및 특징에 대하여 정리하면 Table 1과 같다.

Table 1. Input data characteristics

Data Observation Duration
(Year)
Unit Source Interpolation
Water elevation Hangan bridge 2009~2018 1 hour Han River Flood Control Office ∙Winter season
∙Linear Interpolation
Precipitation Seoul Korea Meteorological Administration ∙Non-Rainfall period
∙Regarding value (zero)

2.2.3 Machine Learning Tool 및 민감도 분석(Sensitivity Analysis)

NARX model을 이용하여 수위예측을 수행하기 위하여 Matlab(Mathworks, version R2019a)에서 제공하는 Deep Learning Toolbox를 활용하였다. Deep Learning Toolbox는 GUI(Graphic User Interface)를 제공하여 입력자료와 예측변량을 설정한 후 모형 변수(Variable)로 은닉층(Hidden Layer) 개수와 시간지연(Time Delay)값을 설정하여 학습 및 예측을 수행한다. 본 연구에서의 예측 변량은 3시간 후의 수위이며, 민감도 분석을 통하여 은닉층 개수와 시간지연의 최적 값을 도출하였다(Fig. 5, Fig. 6).

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F5.jpg
Fig. 5.

Hidden layer sensitivity analysis result

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F6.jpg
Fig. 6.

Time delay sensitivity analysis result

은닉층의 경우, 7개 이상의 값에서 평균제곱근 오차 값은 수렴하는 경향을 보인 반면, 예측 소요시간(CPU Time)은 약 3배 증가되는 것을 확인할 수 있다. 따라서 은닉층의 최적값은 7개로 판단하였다. 시간지연(Time Delay)의 경우, 2시간 이후부터는 평균제곱근 오차 값은 수렴하는 경향을 보인 반면, 예측 소요시간은 약 2배 증가되는 것을 확인할 수 있다. 따라서 시간지연의 최적값은 2시간으로 판단하여 본 모의에서 최적값을 적용하여 한강대교 수위 예측을 수행하였다. Table 2는 본 모의 수행 시에 사용하였던 변수 및 인자 들을 정리하였다.

Table 2. Model setup

Variable Value
Hidden layer 7
Time delay 2 hour
Iterations 1,000
Learning rate 0.02
Prediction of elevation 3 hour

3. 모의 결과 및 분석

전술한 바와 같이 본 연구에서는 NARX 모형의 수위예측 적합성 여부를 판단하기 위해서 기존에 주로 사용되는 인공신경망 모형(Artificial Neural Network model, ANN)과 순환신경망 모형(Recurrent Neural Network model, RNN)의 결과와 NARX 모형의 결과를 비교하였다. 입력자료를 분할하여 학습에 70%, 검정과 평가에 각각 15%씩 사용하여 한강대교 수위를 예측하고 그 성능을 평가하였다. Table 3은 각 모형의 평균제곱근오차(RMSE)와 평균절대오차(MAE), 첨두수위 오차(Peak Error) 및 결정계수(R2) 분석결과를 나타내고 있다. 그 결과 외부입력을 재귀시키는 순환신경망 모형이 인공신경망 모형보다 높은 정확도를 보였고, 출력결과도 함께 재귀시키는 NARX 모형이 가장 높은 예측 정확도를 보여주는 것을 확인하였다.

Table 3. Comparison of accuracy about ANN, RNN, and NARX model

Neural network model RMSE (m) MAE (m) Peak error (m) R2
ANN 0.20 0.12 1.56 0.81
RNN 0.11 0.06 0.55 0.90
NARX 0.09 0.05 0.10 0.94

Fig. 7은 각각의 모형을 이용한 2018년 한강대교 수위 예측결과를 나타내고 있으며, 수위비교 결과에서도 나타났듯이 NARX 모형이 2018년 예측 기간에 대해 저수위 뿐 만아니라 고수위까지 정확하게 예측하는 것을 확인할 수 있다(첨두수위 오차 0.10 m). 반면에 인공신경망 모형과 순환신경망 모형은 관심수위(EL. 3.9 m)이상에서 실측수위와 예측수위의 오차가 NARX 모형보다 크게 발생하는 것을 확인할 수 있다. 이는 인공신경망 모형의 경우, 시간지연(Time delay)등 과거의 입력 자료를 고려하지 못하고 단일 입력자료를 통한 가중치(weight)학습을 통하여 출력 예측 값을 도출하기 때문에 시계열 자료의 변동 추세를 반영할 수 없다. 또한 활성함수(Activation Function)를 sigmoid 함수(출력 값 0 또는 1)를 사용하는 것이 고수위 예측 시 수위가 과소 산정되는 원인이라 판단된다. 이와 달리 순환신경망 모형 및 NARX 모형은 과거의 입력 자료를 고려함으로써 시계열 자료의 변동 추세를 학습 할 수 있으며, 또한 활성함수(Activation Function)를 쌍곡선탄젠트(hyperbolic tangent) 함수(출력 값 –1 또는 1) 또는 Rectified Linear Unit(ReLU) 함수(출력값 0 또는 ∞)를 사용하기 때문에 고수위 예측 시에도 정확한 예측을 수행할 수 있다고 판단된다. 다만 NARX 모형의 경우, RNN 모형에서 출력결과를 추가적으로 재귀시키기 때문에 RNN 모형보다 예측 정확도가 뛰어나다고 판단된다. 따라서 한강대교 수위 예측 시 NARX 모형을 사용하는 것이 적합하다고 판단된다.

http://static.apub.kr/journalsite/sites/ksds/2019-012-02/N0240120208/images/ksds_12_02_08_F7.jpg
Fig. 7.

The water surface elevation distribution (left) and scatter plot (right) according to neural network model

4. 결 론

본 연구에서는 자료기반의 수위예측 알고리즘을 개발하기 위하여 인공신경망(ANN), 순환신경망(RNN) 및 NARX 모형을 활용하였다. 예측변량은 한강대교의 3시간 후의 수위이며, 과거 10년간의 서울시 강수량, 인천 조위, 한강대교 수위자료를 수집하여 입력자료로 활용하였다. 수위예측 알고리즘 적용에 앞서 강수량, 수위 자료는 결측치가 존재하나, 결측치가 발생한 기간이 대부분 비강우시기 또는 동절기여서 각각의 자료 별 유의성을 위해 보간작업을 수행하였으며, 인자간의 상관관계를 분석하였다. 상관관계를 분석한 결과, 조위자료는 수위자료와 상관관계가 없는 것으로 확인된 반면, 강수량은 수위와 양의 상관관계를 가지는 것을 확인하여 수위 예측 시 강수량 자료만 이용하였다. NARX 모형에 사용되는 매개변수는 시간지연(Time delay)과 은닉층(Hidden Layer)의 개수가 있어 각각의 변수에 대하여 민감도 분석을 수행하였다. 최적의 값을 도출하기 위하여 매개변수 별 정확도 및 모의 소요시간을 비교하였고, 시간지연 2시간, 은닉층 7개일 때 가장 정확한 수위 예측이 가능하다고 확인하였다. 민감도 분석을 통하여 산정된 최적 값을 통하여 인공신경망(ANN), 순환신경망(RNN) 및 NARX 모형별 수위예측을 수행하였다. 입력자료의 경우, 70%를 학습에 검정과 평가에 각각 15%를 사용하여 각각의 모형의 평균제곱근오차와 평균절대오차를 비교하였다. 그 결과 현재의 값만을 사용하는 인공신경망 모형의 경우 3시간 후 수위의 예측정확도는 평균제곱근오차는 0.20 m, 평균절대오차는 0.12 m이며, 첨두수위의 오차(Peak error)는 1.56 m로 실측수위와 다소 상이한 예측결과를 나타내었다. 이에 반면 과거의 입력자료를 고려하는 순환신경망 모형의 경우 예측 정확도는 평균제곱근오차는 0.20 m에서 0.11 m로 평균절대오차는 0.12 m에서 0.06 m로 오차가 감소하였으며, 첨두수위의 오차도 0.55 m로 인공신경망 모형보다 정확도 높은 것을 확인하였다. 그리고 과거의 출력을 함께 재귀시키는 NARX 모형이 평균제곱근오차 0.09 m, 절대평균오차 0.05 m 첨두수위 오차 0.10 m로 실측수위와 가장 유사한 예측결과를 나타내었다. 따라서 과거의 입력자료를 고려함으로써 시계열 자료의 변동 추세를 학습하며, 쌍곡선탄젠트(Hyperbolic tangent) 함수 및 Rectified Linear Unit(ReLU) 함수를 활성함수로 사용하는 NARX 모형이 본 연구대상인 한강대교의 수위를 예측할 때 가장 적합한 인공신경망 모델이라고 판단하였다. 다만, 본 연구에서 사용한 순환신경망 모형 및 NARX 모형은 입력자료의 길이(Sequence Length)가 길 경우 기울기 소실(Vanishing gradient)문제가 발생하여 정확한 예측을 수행할 수 없다는 한계점이 있다. 이에 향후에는 NARX 모형의 기울기소실(Vanishing gradient) 문제를 해결한 Long-Short Term Method(LSTM)을 이용하고자 한다. 또한 다양한 활성함수(Activation function) 및 비용함수(Cost Function)를 고려하여 최적의 모형을 구축한다면 보다 정확도 높은 예측결과를 도출 할 수 있을 것으로 기대된다.

Acknowledgements

This work is supported by the Korea Agency for Infrastructure Technology Advancement(KAIA) grant funded by the Ministry of Land, Infrastructure and Transport (Grant 19AWMP-C140010-02).

References

1
Chen, W. B., Liu, W. C., and Hsu, M. H. (2012). Comparison of ANN Approach with 2D and 3D Hydrodynamic Models for Simulating Estuary Water Stage. Advances in Engineering Software. 45(1): 69-79.
10.1016/j.advengsoft.2011.09.018
2
Coulibaly, P. and Anctil, F. (1999). Real-time Short-term Natural Water Inflows Forecasting using Recurrent Neural Networks. Neural Networks. 1999. IJCNN'99. International Joint Conference on, IEEE: 3802-3805.
3
Kim, S. and Tachikawa, Y. (2017). Real-time River-stage Prediction with Artificial Neural Network based on Only Upstream Observation Data. Annual Journal of Hydraulic Engineering. JSCE. 62: 1375-1380.
10.2208/jscejhe.74.I_1375
4
Lee, E.R., Kim W., and Kim, S.H. (2005). Effect of Flood Stage by Hydraulic Factors in Han River. Journal of Korea Water Resources Association. 38(2): 121-131.
10.3741/JKWRA.2005.38.2.121
5
Lee, J. K. and Lee, J. H. (2010). A Study on Water Level Rising Travel Time due to Discharge of Paldang Dam and Tide of Yellow Sea in Downstream Part of Paldang Dam. Journal of the Korean Society of Hazard Mitigation. 10(2): 111-122.
6
Moriasi, D. N., Arnold, J. G., Van Liew, M. W., Bingner, R. L., Harmel, R. D. and Veith, T. L. (2007). Model Evaluation Guidelines for Systematic Quantification of Accuracy in Watershed Simulations. ASABE. 50(3): 885-900.
10.13031/2013.23153
7
Shen, H. Y. and Chang, L. C. (2013). Online Multistep-ahead Inundation Depth Forecasts by Recurrent NARX Networks. Hydrol Earth Syst. Sci. 17: 935-945.
10.5194/hess-17-935-2013
8
Song, C. G., Kim, H. J., and Rhee, D. S. (2014). Analysis of Flow Reversal by Tidal Elevation and Discharge Conditions in a Tidal River. Journal of the Korean Society of Safety. 29(6): 104-110.
10.14346/JKOSOS.2014.29.6.104
9
Thirumalaiah, K. and Deo, M.C. (1998). Real‐Time Flood Forecasting Using Neural Networks. Computer‐Aided Civil and Infrastructure Engineering. 13(2): 101-111.
10.1111/0885-9507.00090

Korean References Translated from the English

1
기상청 (2017). 이상기후 보고서.
2
김현일, 금호준, 한건연 (2018). 도시침수 해석을 위한 동적 인공신경망의 적용 및 비교. 대한토목학회논문집. 38(5): 671-683.
3
한강홍수통제소 (2016). 한강하천예보연감.
페이지 상단으로 이동하기