Analysis of the Differences Between CLAS L6 Message from the QZSS Public Archives and CLAS Correction Messages Received by Actual GNSS Receivers

QZSS 공개 아카이브의 CLAS L6 메시지와 실제 GNSS 수신기에서 수신하는 CLAS 보정정보 차이 분석

Abstract

This study compares and analyzes CLAS correction messages received by an actual GNSS receiver with CLAS L6 data from the QZSS public archive. The Quasi-Zenith Satellite System (QZSS) broadcasts Centimeter Level Augmentation Service (CLAS) correction messages via the L6 band, providing correction information for satellite orbit, clock, code and phase biases, as well as ionospheric and tropospheric errors for GPS, Galileo, and the QZSS constellations. In this study, data collected from the Septentrio AsteRx4 GNSS receiver and the QZSS public archive from October to December 2023 were used for a comparative analysis, with a focus on the variations in the Message Generation Facility (MGF) ID, which indicates the institution generating the CLAS message, and their impact on correction values and positioning accuracy. The results showed that changes in the MGF ID led to variations in correction information, causing differences of about 1-2 cm in the PPP-RTK positioning results.

1. 서론

Global Navigation Satellite System (GNSS)은 전 세계적으로 위치 정보를 제공하는 핵심 기술로, 다양한 분야에서 널리 활용되고 있다. 그러나 GNSS는 위성 궤도 및 시계 오차, 대류권 및 전리층 오차, 다중 경로 오차 등의 요소로 인해 위치 결정 정확도가 저하되는 한계가 존재한다 (Misra & Enge 2011). 이를 해결하기 위해 다양한 오차 모델링 및 보정 기술이 연구되었으며, 크게 Observation Space Representation (OSR) 방식과 State Space Representation (SSR) 방식으로 분류할 수 있다. 그 중에서도 SSR은 GNSS 상시관측소 네트워크에서 관측된 데이터를 취합하여 이를 기반으로 오차 요소별 보정정보를 생성 및 제공하는 방식이다 (Wübbena et al. 2005, Teunissen & Khodabandeh 2015, Hirokawa et al. 2021). SSR 방식은 GNSS 상시관측소 네트워크 내부에서 특정 기준국과의 거리와 관계없이 일관된 정확도를 제공하며, 양방향 통신이 필요하지 않은 장점이 있다. 이러한 이유로 SSR은 GNSS 측위의 신기술로 평가받고 있다.

다양한 SSR 전송 규격 중 Compact SSR은 일본의 지역위성항법시스템인 Quasi-Zenith Satellite System (QZSS)에서 제공하는 센티미터 수준의 보강서비스인 Centimeter Level Augmentation Service (CLAS)를 위한 전송 규격이다 (Hirokawa et al. 2019, Cabinet office 2022). CLAS는 2018년 11월 18일부터 정식 운용을 시작했으며, GPS, Galileo, QZSS 위성에 대한 위성 궤도, 시계, 코드 및 위상 바이어스 오차는 물론, 일본 지역에 대한 전리층 및 대류권 오차 보정정보를 제공한다 (Hirokawa et al. 2023). 이러한 보정 메시지는 QZSS 위성의 L6 밴드를 통해 사용자에게 전송되며, L6 밴드 수신이 가능한 안테나와 수신기를 통해 CLAS 메시지를 직접 수신할 수 있다. 한편, QZSS는 공개 아카이브(public archives)를 통해 과거의 CLAS 메시지를 제공하고 있고 (QZSS 2024a), 대부분의 선행 CLAS 관련 연구는 이러한 공개 아카이브에서 배포하는 후처리된 CLAS 메시지를 기반으로 수행되었다 (Hao et al. 2020, Kim et al. 2022, Lim et al. 2022, Lee et al. 2023).

현재 GNSS 수신기 업계는 QZSS CLAS를 지원하는 기능을 제공하는 추세이며, 대표적으로 u-blox와 Septentrio에서 CLAS를 지원하는 수신기 모듈을 개발 및 보급하고 있다 (QZSS 2024b). u-blox는 CLAS 메시지에 접근할 수 있는 QZSS L6 전용 수신기인 NED-D9C 모듈을 개발하였고, 자사의 ZED-F9P/F9R 등의 모듈과 결합하여 사용함으로써 CLAS 기반 정밀측위 기능을 지원한다. Septentrio는 AsteRx4 등의 수신기와 mosaic-go CLAS 등의 모듈에 QZSS CLAS 메시지를 수신 및 저장할 수 있는 기능을 탑재하였고, Septentrio Binary Format (SBF)라는 Septentrio 사의 독자적인 저장 형식으로 CLAS 메시지를 저장할 수 있다 (Septentrio 2023, 2024). 이를 기반으로 CLAS를 적용한 GNSS 정밀측위 알고리즘을 개발 및 검증한 연구가 수행된 바 있다 (Lee & Park 2020, Lee 2021).

CLAS 메시지를 활용한 Precise Point Positioning – Real-Time Kinematic (PPP-RTK) 측위 알고리즘을 개발하는 과정에서, QZSS 아카이브에서 제공되는 후처리된 CLAS 메시지와 실시간으로 GNSS 수신기를 통해 직접 수신한 CLAS 메시지를 적용한 측위 결과에 차이가 있음을 확인하였다. 이러한 차이의 원인을 규명하기 위해 국내외의 선행 연구를 조사하였으나, Lee & Park (2020)과 Lee (2021)를 제외하고 실시간 환경에서 GNSS 수신기를 통해 CLAS 메시지를 직접 수신 및 사용하는 연구는 확인되지 않았다. 또한, CLAS 메시지 간의 차이를 비교하고, 그 차이가 GNSS 측위 성능에 미치는 영향을 분석한 연구는 수행되지 않았다. 따라서 본 연구에서는 실시간으로 수신된 CLAS 메시지와 QZSS 공개 아카이브의 CLAS 메시지 간의 차이를 직접 비교하여, 그 차이의 특성을 분석하고 발생 원인을 규명하였다. 그리고 이러한 차이가 실제 PPP-RTK 측위 성능에 어떠한 영향을 미치는지 파악하였으며, 실시간 환경에서 CLAS 메시지를 사용함에 따라 발생하는 문제점과 사용자 단에서 유의해야 할 사항을 제시하고자 한다.

본 논문의 구성은 다음과 같다. 2장에서는 CLAS 보정 메시지의 이론적 배경을 설명한다. 이후 3장에서 QZSS 공개 아카이브의 CLAS 메시지와 실시간으로 수신한 CLAS 메시지 간의 차이를 분석하고, 그 원인을 규명한다. 마지막으로 결론에서는 연구 결과를 종합하며 실시간 환경에서 CLAS 메시지를 사용할 시 유의할 사항을 언급하고, 향후 연구 방향을 논의한다.

2. CLAS 보정메시지 개요

CLAS L6 메시지는 헤더(header), 데이터 파트(data part), 그리고 reed-Solomon 코드(code)로 구성된다. Fig. 1은 CLAS L6 메시지의 구성을 나타낸 것이다.

Fig. 1. CLAS L6 message structure (Cabinet office 2022).

2.1 Header of CLAS L6 Message

CLAS 메시지의 헤더는 Fig. 2와 같이 Preamble 32비트, PRN 8비트, L6 message type ID 8비트 그리고 Alert Flag 1비트로 총 49비트로 구성되어 있다. Preamble은 메시지의 시작을 알리며, PRN은 L6 메시지를 송출한 QZSS 위성의 Pseudo Random Noise (PRN) 번호를 나타낸다. L6 message type ID는 Table 1에 기술한 것과 같이 7-5비트에서 Vendor ID, 4-3비트에서 Message Generation Facility ID (MGF ID)를 표현하고, 마지막 비트에는 서브프레임(subframe)의 시작 여부를 나타내는 Subframe indicator가 위치해 있다.

Fig. 2. Header of CLAS L6 message (Cabinet Office 2022).

Table 1. L6 message type ID (Cabinet Office 2022).

Bit field	Date name	Note
7-5	Vendor ID	“101” for CLAS
4-3	Message Generation Facility ID	“00”, “01”: Hitachi-Ota “10”, “11”: Kobe
2-1	Reserved	N/A
0	Subframe indicator	“1”: first data part of a subframe “0”: others

 

그 중에서 특히 MGF ID는 CLAS 메시지를 생성한 기관을 식별할 수 있는 중요한 정보로, CLAS 메시지의 출처를 나타낸다. Table 1에 따르면, MGF ID가 "00" 또는 "01"일 경우 히타치 오타(Hitachi-Ota) 주관제소에서 생성된 메시지이며, "10" 또는 "11"은 고베(Kobe) 주관제소에서 생성된 메시지임을 의미한다. 이처럼 MGF ID는 CLAS 메시지 생성 기관을 구분하고, CLAS 메시지를 생성한 주관제소에 따른 보정정보의 특성 및 경향의 차이를 분석하는 데 중요한 역할을 한다. 실제로 QZSS의 지상국(ground system)은 Fig. 3과 같이 2개의 주관제소, 7개소의 Telemetry, Tracking and Command (TTC) 스테이션으로 구성되어 있으며 (Numuta 2024), 이 중 히타치 오타 주관제소는 QZSS의 메인 관제소로서 시스템 운영과 위성의 실시간 모니터링 및 제어를 담당하며, 고베 주관제소는 백업 관제소로 히타치 오타 주관제소의 역할을 보조하거나 필요 시 대신한다.

Fig. 3. QZSS Ground system (Numuta 2024).

2.2 Data Part of CLAS L6 Message

CLAS의 데이터 파트에는 Compact SSR 형식의 GPS, Galileo, QZSS 위성의 궤도(orbit), 시계(clock), 코드 및 위상 바이어스(code & phase bias), 그리고 일본 전역의 대기권(atmospheric) 오차에 대한 보정정보가 포함된다. CLAS 보정정보의 종류와 Compact SSR 메시지 제공 간격은 Table 2와 같다. 이러한 정보는 QZSS 위성의 L6 밴드를 통해 방송되며, L6 밴드를 수신할 수 있는 GNSS 수신기는 CLAS 메시지를 수신할 수 있다.

Table 2. Messages of compact SSR and update interval (Cabinet Office 2022).

Message name	Sub type ID	Nominal update interval (s)
Compact SSR mask Compact SSR GNSS orbit correction Compact SSR GNSS clock correction Compact SSR GNSS satellite code bias Compact SSR GNSS satellite phase bias Compact SSR GNSS satellite code and phase bias Compact SSR GNSS URA Compact SSR STEC correction Compact SSR Gridded correction Compact SSR Service Information Compact SSR GNSS combined correction Compact SSR atmospheric correction Null message	1 2 3 4 5 6 7 8 9 10 11 12 N/A	30 30 5 30 30 30 30 30 30 N/A 5 or 30 30 N/A

 

2.2.1 Compact SSR GNSS orbit correction (sub type 2)

Compact SSR의 GNSS Orbit correction message에는 위성 별 Radial/Along-Track/Cross-Track 방향의 위성 위치 보정량과, 위성의 위치를 계산하는데 사용해야 할 방송궤도력의 Issue of Data Ephemeris (IODE)가 포함되어 있다. 위성 궤도 보정량은 매 30초마다 갱신된다. 사용자는 궤도 보정량과 함께 제공된 IODE에 해당하는 궤도력을 사용하여 위성의 위치를 계산하고, Radial/Along-Track/Cross-Track 방향의 위성 위치 보정량을 Earth-Centered Earth-Fixed (ECEF) 좌표계 상의 보정량으로 환산하여 위성 궤도 오차를 보정한다. 그 방법은 Eqs. (1-5)와 같다.

\(e_{along}=\frac{\dot{r} (t) }{\left|\dot{r} (t)\right|} \)                                                                                         (1)

\(e_{cross}=\frac{r (t) \times\dot{r} (t)}{\left|r (t) \times\dot{r} (t)\right|}\)                                                                                    (2)

\(e_{radial}=e_{along}(t)\times e_{cross} (t)\)                                                                           (3)

\(\delta X=\left[e_{radial}\ e_{along}\ e_{cross}\right] \left[\begin{matrix}\delta O_{radial} \\ \delta O_{along} \\ \delta O_{cross} \end{matrix}\right]\)                                                                (4)

\(X_{orbit}=r(t)-\delta X\)                                                                                         (5)

여기서 \(r(t)\)와 \(\dot{r} (t)\)는 각각 위성의 위치 및 속도벡터를 나타내고 \(\delta O_{radial}, \delta O_{along}, \delta O_{cross}\)는 각각 radial, along-track, cross-track 방향의 미터 단위 보정량을 의미한다. 또한 \(\delta X\)는 ECEF 좌표계 상에서의 위성 궤도 보정량, \(X_{orbit}\)은 궤도 보정이 적용된 위성 위치 벡터이다.

2.2.2 Compact SSR GNSS clock correction (sub type 3)

5초 간격으로 갱신되는 Compact SSR GNSS clock correction message에서는 방송궤도력으로 모델링한 위성 시계 오차(offset)에 대한 보정량을 제공한다. 보정 방법은 Eq. (6)과 같다.

\(T_{sat}=T_{broadcast}+\frac{\delta C}{C}\)                                                                                 (6)

여기서 \(T_{sat}\)은 위성 시계 보정이 적용된 위성 시계 오차, \(T_{broadcast}\)는 방송궤도력으로 모델링한 위성 시계 오차이며 \(\delta C\)는 미터 단위의 위성 시계 오차 보정량, \(C\)는 진공에서의 빛의 속도(299,792,458 m/s)이다.

2.2.3 Compact SSR GNSS code bias & phase bias correction (sub type 4 & 6)

CLAS는 Sub Type 4에서 코드 바이어스를, Sub Type 6에서 위상 바이어스를 제공하고 있다. 특히, Sub Type 6의 위상 바이어스는 코드 바이어스와 달리 CLAS 격자망의 각 네트워크에 대한 위상 바이어스 값이 제공된다. 따라서 CLAS 사용자는 본인이 위치한 지점의 네트워크 번호를 정확히 인지하고, 해당 네트워크 번호에 대한 위상 바이어스 보정정보를 적용해야 함을 유의해야 한다. Fig. 4는 CLAS 격자 네트워크 및 격자점 분포를 나타내며, 코드 및 위상 바이어스의 보정방법은 Eqs. (7-8)에 제시하였다.

\(P{R_j^s}^\prime=PR_j^s-{BIAS}_{PR_j^s}\)                                                                           (7)

\({{CP}_j^s}^\prime={CP}_j^s-{BIAS}_{{CP}_j^s}\)                                                                           (8)

여기서 \(PR_j^s, {CP}_j^s\)는 사용자가 받은 \(s\) 위성의 \(j\) 채널 코드의사거리와 반송파위상 관측값을 의미하고, \({BIAS}_{PR_j^s}\)는 \(PR_j^s\)에 대한 코드 바이어스를, \({BIAS}_{{CP}_j^s}\)는 \({CP}_j^s\)에 대한 위상 바이어스를 의미한다. \(P{R_j^s}^\prime\)과 \({{CP}_j^s}^\prime\)는 코드 및 위상 바이어스 보정이 적용된 GNSS 관측값을 의미한다.

Fig. 4. Definition of CLAS network and grid points.

2.2.4 Compact SSR atmospheric correction (sub type 12)

CLAS 서브타입 12(compact SSR atmospheric correction)는 대류권과 전리층 오차에 대한 보정정보를 제공하기 위한 메시지로, 일본 전역 214개 격자점에 대한 대류권 및 전리층 다항 계수(polynomial coefficients)와 잔차(residuals)로 구성되어 있다.

Compact SSR의 Sub Type 12에서는 대류권 건조지연을 모델링할 수 있는 대류권 다항계수(troposphere polynomial coefficients)와 습윤지연을 모델링할 수 있는 대류권 잔차(troposphere residual) 및 대류권 잔차 오프셋(troposphere residual offset)을 제공하며, 사용자 위치에서의 천정방향 대류권 지연 계산 방법은 Eqs. (9-11)과 같다.

\(\mathrm{\Delta}T_{hs}=T_{00}+T_{01}\left(\phi-\phi_0\right)+T_{10}\left(\lambda-\lambda_0\right)+T_{11}\left(\phi-\phi_0\right)\left(\lambda-\lambda_0\right)\)                                       (9)

\(T_{hs}={\mathrm{\Delta}T_{hs}+\delta}_{hs0}\)                                                                                     (10)

\(T_{wet}={\mathrm{\Delta}T_{wet}+\delta}_{wet0}\)                                                                                   (11)

여기서, \(T_{00}  \sim  T_{11}\)은 대류권 다항계수를, \(\phi\)와 \(\phi_0\)은 각각 사용자 위치의 위도, 격자점의 위도이며 \(\lambda\)와 \(\lambda_0\)는 각각 사용자 위치의 경도, 격자점의 경도이다. \(T_{hs}\)와 \(T_{wet}\)는 각각 계산된 대류권 건조지연 및 습윤지연 보정량을 의미하고, \(\delta_{hs0}\)은 건조지연 변동 상수(hydrostatic variation constant value)로 2.3 m이며 \(\mathrm{\Delta}T_{wet}\)는 대류권 잔차, \(\delta_{wet0}\)은 대류권 잔차 오프셋이다. 즉, 대류권 건조지연 보정량은 특정 네트워크에 대해 제공되는 다항계수를 사용해 사용자의 위치에 대한 보정량을 계산하는 방식으로 계산된다. 반면, 대류권 습윤지연 보정량은 사용자가 특정 격자점에 대해 계산한 값을 적용하는 방식을 따른다.

이후, 대류권 오차 모델링의 정확도를 향상시키기 위해 인접한 여러 격자점으로부터 사용자 위치에서의 대류권 오차를 모델링한다. 이때 주로 주변 4개의 격자점을 활용하며 (Cabinet Office 2022), 인접 4개 격자점으로부터 모델링한 대류권 오차를 보간하기 위해 각 격자점과 사용자의 위치 간의 거리를 사용하여 각 격자점에 대한 가중치를 Eq. (12)와 같이 계산한다. 그리고 이를 바탕으로 Eqs. (13-14)를 통해 대류권 건조지연 및 습윤지연량을 선형보간하고, Eq. (15)와 같이 보간된 값에 사상함수를 적용하여 시선방향 대류권 총 지연량을 계산한다.

\(W_k=\frac{\frac{1}{R_{user,\ k}}}{\frac{1}{R_{user,\ 1}}\ +\ \frac{1}{R_{user,\ 2}}\ +\ \frac{1}{R_{user,\ 3}}\ +\ \frac{1}{R_{user,\ 4}}}(k=1,\ 2,\ 3,\ 4)\)                                                  (12)

\(T_{user,hs}=\sum_{k=1}^{4}{W_kT_{k,hs}}\)                                                                              (13)

\(T_{user,wet}=\sum_{k=1}^{4}{W_kT_{k,wet}}\)                                                                              (14)

\(T_{user}^{PRN}=T_{user,hs}\ast M_{hs}^{PRN}+ T_{user,wet}\ast M_{wet}^{PRN}\)                                                             (15)

여기서 \(W_k\)는 특정 격자점 \(k\)에 대한 가중치이고, 격자점 \(k\)와 사용자 위치 간의 거리 \(R_{user, k}\)의 역수와 인접한 4개 격자점(격자점 \(k\) 포함)과 사용자 위치 간의 거리 \(R_{user, 1}\)부터 \(R_{user, 4}\)까지의 역수의 총합의 비율으로 계산된다. \(T_{user, hs}\)와 \(T_{user, wet}\)는 인접한 4개 격자점으로부터 모델링한 천정방향의 대류권 건조지연 및 습윤지연량을 보간한 값이며, \(T_{user}^{PRN}\)은 사상함수 \(M_{hs}^{PRN}, M_{wet}^{PRN}\)를 적용하여 계산된 시선방향 대류권 지연량을 의미한다.

전리층 오차 보정량은 인접한 각 격자점에서의 위성 별 Slant Total Electron Content (STEC)을 계산하고, 그 STEC 값들을 보간하여 사용자 위치에서의 STEC 값을 산출하는 방식으로 계산된다. 특정 격자점에서의 STEC 계산 방법은 Eq. (16)과 같다.

    \(STEC_k^{PRN}=C_{00}^{PRN}+C_{01}^{PRN}\left(\phi-\phi_0\right)+C_{10}^{PRN}\left(\lambda-\lambda_0\right)\)                                                                                          
\({+\ C}_{11}^{PRN}\left(\phi-\phi_0\right)\left(\lambda-\lambda_0\right) + C_{02}^{PRN}\left(\phi-\phi_0\right)^2+C_{20}^{PRN}\left(\lambda-\lambda_0\right)^2+\delta STEC_k^{PRN}\)       (16)

여기서 \(C_{00}^{PRN}\sim C_{20}^{PRN}\)은 STEC 다항 계수(polynomial coefficients)이고, \(\delta STEC_k^{PRN}\)는 STEC 잔차(residual)를 의미한다. 또한 \(\phi, \lambda와 \phi_0, \lambda_0\)는 각각 특정 격자점의 위도와 경도, 해당 격자점이 속한 네트워크에서의 1번 격자점의 위도와 경도이다. 예를 들어 네트워크 4의 3번 격자점에서의 STEC을 계산할 때 \(\phi, \lambda\)는 네트워크 4의 3번 격자점의 위도/경도, \(\phi_0, \lambda_0\)는 네트워크 4의 1번 격자점의 위도/경도가 된다.

이후 각 격자점에서 구한 STEC을 미터 단위의 전리층 오차 보정량으로 환산한 뒤 Eqs. (12, 17)에서 구한 각 격자점에 대한 가중치를 활용하여 이를 보간함으로써 사용자 위치에서의 전리층 오차 보정량을 Eq. (18)과 같이 계산한다.

\(I_k^{PRN}=\frac{40.31\ast{10}^{16}}{f^2}\cdot STEC_k^{PRN}\)                                                                      (17)

\(I_{user}^{PRN}=\sum_{k=1}^{4}{W_kI_k^{PRN}}\)                                                                            (18)

여기서 \(f\)는 위성 신호의 주파수(Hz), \(I_k^{PRN}\)는 미터 단위의 각 격자점에서의 전리층 오차 보정량이며 \(I_{user}^{PRN}\)는 인접 4개 격자점에서의 전리층 오차 보정량을 보간한 값이다. 이렇게 계산된 전리층 오차 보정량 \(I_{user}^{PRN}\)은 위상 관측값에서는 앞섬 효과를, 코드 관측값에서는 지연 효과를 유발한다. 따라서, 코드의사거리 관측방정식에는 양(+)의 부호로, 반송파위상 관측방정식에는 음(-)의 부호로 각각 적용해야 한다.

3. QZSS 공개 아카이브의 CLAS 메시지와 QZSS 위성으로부터 직접 수신한 CLAS 메시지 비교

이번 장에서는 QZSS 공개 아카이브에서 제공되는 CLAS 메시지(이하 편의상 archive)와, 실시간 환경에서 GNSS 수신기로 직접 수신한 CLAS 메시지(이하 편의상 receiver) 간의 차이를 탐색 및 분석하였다. 이를 위해 QZSS 웹사이트의 QZSS 공개 아카이브에서 배포하는 후처리된 CLAS 메시지를 다운로드하였으며 (QZSS 2024a), 직접 GNSS 수신기로 CLAS 메시지를 수신하기 위하여 (주)피피솔에서 운영 및 관리하고 있는 울산 GNSS 상시관측소인 PPBY에 설치된 Septentrio AsteRx4 수신기를 사용하였다. Septentrio의 AsteRx4는 CLAS 메시지 처리와 저장이 가능한 수신기이며, CLAS 메시지를 직접 취득 및 SBF 형식으로 저장이 가능하다 (Septentrio 2023). (주)피피솔에서 운영하는 PPBY의 주변 관측 환경은 Fig. 5, 상세 정보는 Table 3과 같다 (PP-Solution Inc. 2024).

Fig. 5. Surrounding environment of PPBY (PP-Solution Inc. 2024).

Table 3. Details of PPBY (PP-Solution Inc. 2024).

Item	Note
Reference station	PPBY (ByeongYeong)
IGb14 Coordinates @ Sep 1, 2023	Latitude: 35° 34' 9.34608" N  Longitude: 129° 20' 41.35202" E Ellipsoid height: 59.324 m
Antenna	Tallysman TWIVC6050 (SCIS)
Receiver	Septentrio AsteRx4

 

3.1 Receiver와 Archive의 CLAS 송출 위성 및 MGF ID 변화 양상 확인

2023년 10월 28일(DOY 300)부터 2023년 12월 31일(DOY 365)까지 archive와 receiver 간의 비교 및 분석을 수행하였으며, Fig. 6은 이 기간 동안 archive와 receiver에 명시된 CLAS 송출 위성 정보와 MGF ID를 나타낸 것이다.

Fig. 6. PRN & MGF ID from CLAS of Archive and Receiver during DOY 300-365, 2023.

이번 연구에서 사용한 Septentrio 사의 GNSS 수신기에서는 QZSS L6 밴드를 처리할 때 관측된 모든 QZSS 위성의 L6 밴드 관측값을 처리하지 않고, 임의의 한 위성의 L6 밴드 관측값만 처리한다. 따라서 QZSS 위성으로부터 CLAS 메시지를 수신할 때, QZSS 위성의 위치에 따라 CLAS 메시지를 받아오는 QZSS 위성이 달라지고 CLAS 메시지 송출 위성 정보가 변경된다. 그 결과 Fig. 6과 같이 Archive의 CLAS 송출 위성은 J01으로 고정되어 있는 반면, Receiver의 CLAS 송출 위성 번호는 일정한 패턴을 따라 변화한다.

MGF ID의 경우, archive와 receiver 모두에서 변동 양상이 나타났다. archive의 MGF ID 변동 빈도는 receiver에 비해 현저히 낮았지만, 공개 아카이브에서 제공되는 메시지의 생성 기관 또한 간헐적으로 변경될 수 있다는 사실을 확인하였다. 일반적으로 한 시간 또는 수 시간 동안 MGF ID는 "00"과 "01" 사이에서 변동하지만, DOY 309에서는 MGF ID가 "10"과 "11"로 나타나기도 하며 이는 고베 주관제소에서 생성된 보정정보가 공개 아카이브에서 배포될 수 있음을 의미한다.

3.2 CLAS 송출 위성 및 MGF ID 변화에 따른 CLAS 보정량 차이 분석

Receiver의 CLAS 송출 위성의 변화 또는 MGF ID 변화에 따른 CLAS 보정량의 변화 여부와 양상을 확인하기 위해 다음 두 가지 사례의 archive와 receiver 데이터를 비교 및 분석하였다. 첫 번째 사례(이하 Case 1)는 MGF ID가 변동하는 2023년 5월 27일(DOY 086) UTC 기준 06시에서 07시 사이의 데이터를 사용한 것이며, 두 번째 사례(이하 Case 2)는 MGF ID는 동일하지만 CLAS 송출 위성만 변동하는 같은 일자 09시에서 10시 사이의 데이터를 사용한 것이다. Fig. 7은 Case 1과 Case 2에서의 archive와 receiver의 CLAS 송출 위성과 MGF ID 변화를 나타낸다. Fig. 7에 따르면 Case 1에서 receiver의 CLAS 송출 위성은 변하지 않았으며, MGF ID만 01에서 00으로 변경되었다. 이는 수신기에서 QZSS L6 밴드 관측값을 받아오는 위성은 변함이 없으나, 메시지를 생성한 기관만 변경되었음을 의미한다. 반면, Case 2에서는 MGF ID는 변경되지 않았으나 receiver의 CLAS 송출 위성은 J02에서 J03으로 달라졌다.

Fig. 7. PRN & MGF ID of CLAS messages from archive and receiver at Case 1 (left) & Case 2 (right).

Case 1과 Case 2에서 archive와 receiver 간의 정량적인 보정정보 차이 여부를 확인하기 위해 궤도, 시계, 신호 바이어스, 대류권 및 전리층 오차 보정량의 총합인 Pseudo Range Correction (PRC) 및 Carrier Phase Correction (CPC)의 시계열 변화를 분석하였다. Figs. 8과 9는 각각 Case 1과 Case 2에서 archive와 receiver의 보정량 총합과 차이를 나타낸다. GPS L1/L2, Galileo E1/E5a 그리고 QZSS L1/L2 주파수에 해당하는 PRC 및 CPC를 산출하였고, 그 계산 방법은 Eqs. (18, 19)와 같다 (Cabinet office 2022). 이때 궤도 보정량 \(\delta X^{PRN}\cdot los_{user}^{PRN} (t)\), 대류권 오차 보정량 \(T_{user}^{PRN} (t)\) 그리고 전리층 오차 보정량 \(I_{user}^{PRN} (t)\)은 일본에 위치한 International GNSS Service (IGS) 상시관측소 중 하나인 STK2를 기준으로 계산되었다.

\({CPC}_{user}^{PRN} (t)=-\delta C^{PRN} (t)-\delta X^{PRN}\cdot los_{user}^{PRN} (t)+BIAS_{CP^{PRN}}-I_{user}^{PRN} (t)+T_{user}^{PRN} (t)\)           (19)

\({PRC}_{user}^{PRN} (t)=-\delta C^{PRN} (t)-\delta X^{PRN}\cdot los_{user}^{PRN} (t)+BIAS_{{PR}^{PRN}}+I_{user}^{PRN} (t)+T_{user}^{PRN} (t)\)           (20)

Fig. 8. Total of CLAS correction and its difference between archive and receiver in Case 1.

Fig. 9. Total of CLAS correction and its difference between archive and receiver in Case 2.

Fig. 8은 receiver의 MGF ID가 바뀌었을 때, GPS, Galileo, QZSS 시스템 별 임의의 한 위성의 보정량 및 보정량 차이의 시계열 변화를 보여준다. Receiver의 MGF ID가 01에서 00으로 변경되는 06시 12분부터 보정량 간의 차이가 사라지는 것을 확인할 수 있다. 이를 통해 MGF ID가 동일한 두 CLAS 메시지는 보정량 차이가 존재하지 않는다는 것을 알 수 있다. 또한, 임의로 선정한 G14, E15, J02 외의 CLAS 보정정보가 제공된 다른 모든 위성에서도 보정량 차이가 발생할 때 그 크기만 미세하게 달라질 뿐, receiver의 MGF ID가 01에서 archive와 동일한 00으로 변경될 때 보정량 간 차이가 사라지는 양상은 동일하게 나타났다. Case 1에서 PRC는 ±10 cm 이내의 차이가 발생했고, CPC에서는 최대 40~50 cm의 차이가 발생했으며, MGF ID가 다른 모든 경우에서도 유사한 크기의 차이가 일관되게 관찰되었다.

반면 receiver의 CLAS 송출 위성이 J02에서 J03으로 변경되는 Case 2에서는 receiver와 archive 간의 보정량 간 차이는 확인되지 않았다. 이를 통해 CLAS 송출 위성의 변동이 CLAS 메시지의 내용에 직접적인 영향을 미치지 않는다는 사실을 알 수 있었다. 단, 송출 위성이 변경되는 시점(09시 03분)에 receiver에서 25초 동안 일시적인 데이터 누락이 발생했으며, 이는 새로운 QZSS 위성으로부터 L6 밴드 관측값을 수신하는 과정에서 발생한 현상으로 파악된다. Case 1과 마찬가지로 임의로 선정한 위성 외의 모든 위성에서도 동일한 양상이 나타났으며, MGF ID가 동일한 다른 모든 경우에서 일관된 결과가 나타났다.

3.3 MGF ID 차이에 따른 PPP-RTK 측위 결과 비교

본 연구에서는 자체 개발한 PPP-RTK 측위 알고리즘을 적용하여 MGF ID에 따른 측위 결과를 비교 및 분석하였다. 사용된 PPP-RTK 측위 알고리즘의 자료처리 전략은 Table 4에 요약하였다. 이 알고리즘은 GPS, Galileo, QZSS의 결합되지 않은 이중 주파수 코드의사거리와 반송파위상 관측값을 처리하고, 이에 따른 inter system bias 및 수신기 단의 differential code bias를 함께 추정한다. 오차 보정은 CLAS에서 제공되는 위성 궤도, 시계, 신호 바이어스와 전리층 및 대류권 보정정보를 활용하여 수행된다. 특히, PPP-RTK 측위 알고리즘의 안정도 향상을 위해 전리층과 대류권 오차 보정정보를 가상의 추가 관측치로 취급하여, 추정 변수에 포함시키고 모델에 제약 조건을 설정한다 (An et al. 2023). 또한, 단독측위에서 고려되어야 할 각종 조석에 의한 지구 지각 변위 (Petit & Luzum 2010), Phase Wind Up (Wu et al. 1993), 상대성 효과 그리고 수신기 안테나 PCO/PCV 등의 정밀오차를 모델링하고 보정한다. 모호정수 결정은 Least-squares AMBiguity Decorrelation Adjustment (LAMBDA) 기법을 적용하여 수행되었고 (Teunissen 1995), 모호정수 결정 성공 여부의 판정(accept test)은 모호정수 결정 성공률(fix success rate) 확인, 비율 검정(ratio test)을 거쳐 이루어졌다 (Teunissen & Verhagen 2008). PPP-RTK 측위 알고리즘의 추정기(estimator)로는 확장 칼만 필터(extended Kalman filter)를 사용하였다.

Table 4. Processing strategies of PPP-RTK positioning algorithm.

Item

Descriptions

System Signal selection (RINEX 3) Inter system bias Satellite orbit & clock Satellite signal bias Receiver code bias Tropospheric delay Ionospheric delay Earth displacements Phase wind up Receiver PCO/PCV Relativity effect Ambiguity Estimator

GPS + Galileo + QZSS GPS, 1C/2W; Galileo, 1X/5X; QZSS, 1C/2X Estimated Broadcast ephemeris + CLAS correction CLAS Correction Differential code bias estimation for each GNSS CLAS correction: estimated and treated as pseudo-observations (An et al. 2023) CLAS correction: estimated and treated as pseudo-observations (An et al. 2023) Solid Earth Tides, Ocean Tides Loading and Pole Tides (Petit & Luzum 2010) Wu et al. (1993) igs20.atx Considered: Satellite signal & clock Fixed (Teunissen 1995, Teunissen & Verhagen 2008) Extended Kalman filter

 

MGF ID가 달라질 때 측위 정확도에 유의미한 영향을 미치는지 확인하기 위해, 2023년 11월 8일(DOY 312) 00:00–02:00와 17:00–19:00 시간대의 archive 및 receiver의 CLAS 메시지를 활용하고 STK2의 GNSS 관측값을 처리하여 모호정수 결정 적용 후의 위치결정 결과와 모호정수 결정률의 차이를 중점적으로 분석하였다. Fig. 10은 DOY 312의 archive와 receiver의 MGF ID 변화를 나타낸 그림으로, 2023년 DOY 312에서 receiver의 MGF ID가 00시 40분과 17시 52분에 변경되는 반면, archive의 MGF ID는 항상 00으로 고정되어 변화하지 않는 것을 알 수 있다.

Fig. 10. MGF ID of CLAS messages from archive and receiver at DOY 312, 2023.

한편, 측위 결과 산출 과정에서 매 정각마다 필터를 초기화했으며, 갑작스러운 CLAS 보정량 변화로 인한 측위 알고리즘의 비정상적인 동작을 방지하기 위해 receiver에서 MGF ID가 변동될 때도 필터를 초기화하였다. Fig. 11과 Table 5는 archive 및 receiver의 CLAS 메시지를 적용한 2023년 DOY 312 00시부터 02시까지의 PPP-RTK 측위 결과를 시각적 및 수치적으로 나타낸 것이다.

Fig. 11. Positioning results and their differences during 00:00-02:00, DOY 312, 2023.

Table 5. Positioning results of archive and receiver during 00:00-02:00 at DOY 312, 2023.

Time	CLAS type	MGF ID	RMSE (cm)			Fix rate (%)
			ΔN	ΔE	ΔV
00:00 – 00:40	Archive Receiver	00 01	1.8 2.9	1.8 1.7	5.4 5.9	78.8 80.6
00:41 – 02:00	Archive Receiver	00 00	2.4 2.2	1.9 1.9	9.0 9.1	98.6 98.6

 

Fig. 11과 Table 5에 따르면 receiver의 MGF ID가 01에서 00으로 변동되기 전, 00시 40분에는 archive와 receiver의 측위 결과 간에 미세한 차이가 존재하였다. 이 차이는 주로 약 1~2 cm 크기로 무시할 수 있는 수준이었으나, 간헐적으로 약 10 cm 수준의 차이가 발생하는 경우도 있었다. MGF ID가 변경되는 시점에는 receiver의 PPP-RTK 측위 알고리즘이 초기화되면서 측위 결과가 일시적으로 발산하였고, 그에 따라 차이가 증가하였다. 그러나 이후 수렴 과정에서 측위 오차의 차이는 감소하였고, archive와 receiver의 PPP-RTK 측위 알고리즘이 모두 초기화된 01시 정각 이후 측위 결과는 완전히 일치하였다.

동일한 일자의 17시부터 19시까지의 archive 및 receiver CLAS 메시지를 적용한 PPP-RTK 측위 결과는 Fig. 12와 Table 6과 같다.

Fig. 12. Positioning results and their differences during 17:00-19:00, DOY 312, 2023.

Table 6. Positioning results of archive and receiver during 17:00-19:00 at DOY 312, 2023.

Time	CLAS type	MGF ID	RMSE (cm)			Fix rate (%)
			ΔN	ΔE	ΔV
17:00 – 17:52	Archive Receiver	00 00	2.0 2.0	3.0 3.0	9.3 9.3	99.4 99.4
17:53 – 19:00	Archive Receiver	00 01	2.1 3.2	5.5 5.6	6.0 7.2	97.9 97.0

 

Fig. 12와 Table 6을 통해 17시 52분까지는 receiver의 MGF ID가 00으로 유지되며 archive와 receiver의 측위 결과는 완전히 동일함을 확인할 수 있다. 그러나 receiver의 MGF ID가 01로 변경되는 시점(17시 52분)에서, receiver의 PPP-RTK 측위 알고리즘이 초기화되면서 일시적으로 archive와 receiver의 측위 결과 간에 큰 차이가 발생하였다. 이후 수렴 과정에서 차이는 점차 줄어들었으나, 수렴 이후에도 그 차이가 완전히 사라지지는 않았음을 확인할 수 있다. 이는 17시 52분 이후에는 receiver의 MGF ID가 01으로 archive의 MGF ID와 다르기 때문으로 파악된다.

위의 결과를 요약하면, 00시부터 02시는 receiver와 archive의 MGF ID가 달랐다가 동일해지는 시간대로, receiver의 MGF ID가 변경되는 00시 40분 이전에는 1-2 cm 수준의 측위 정확도 차이와 2.2%p의 모호정수 결정률 차이가 발생하였다. 또한, 01시 정각 이후 측위 결과는 완전히 일치하였다. 반면, 17시부터 19시는 receiver의 MGF ID가 archive와 동일했다가 달라지는 시간대로, receiver의 MGF ID가 변경되는 시점인 17시 52분 이전에는 측위 결과가 완전히 동일하였다. MGF ID가 변경된 후에는 1-2 cm 수준의 측위 정확도 차이 그리고 0.9%p의 모호정수 결정률 차이가 확인되었다. 즉, MGF ID가 동일할 때는 보정정보가 같아 측위 결과가 일치하는 반면, MGF ID가 다를 경우 측위 결과가 일치하지 않는다는 점이 명확하게 드러났다. 이 결과를 통해 receiver와 archive의 MGF ID가 다를 경우 보정정보 간의 차이로 인해 측위 결과에도 차이가 발생함을 규명하였다.

4. 결론

본 연구는 실시간 환경에서 CLAS 메시지를 수신할 때, GNSS 수신기로 직접 수신한 CLAS 메시지와 QZSS 공개 아카이브에서 배포한 CLAS 메시지 간의 차이를 비교 및 분석하였다. 2023년 10월 28일부터 12월 31일까지 데이터를 수집하여 분석한 결과, CLAS 메시지를 송출하는 위성의 변화는 CLAS 보정정보와 직접적인 관련이 없음을 확인하였다. 단, QZSS L6 밴드 관측값을 수신해오는 위성이 변경될 시 L6 밴드 관측값을 다시 추적하는 과정에서 일시적인 메시지 수신 누락(5-30초)이 발생할 수 있어, 실시간 환경에서 CLAS 메시지를 사용하는 경우 이를 감안해야 할 필요가 있을 것으로 보인다. 또한, CLAS 메시지의 생성기관을 의미하는 MGF ID가 간헐적으로 변경되고, 특히 실시간 환경에서 GNSS 수신기로 직접 수신한 CLAS 메시지와 QZSS 공개 아카이브의 후처리된 CLAS 메시지 간에 MGF ID의 차이가 존재할 수 있음을 입증하였다. MGF ID에 따라 보정량의 차이가 발생하지만, 그 크기에 비해 PPP-RTK 측위 결과에 미치는 직접적인 영향은 1-2 cm 수준에 불과했다. 다만, 이러한 결과는 MGF ID가 변경되는 시점에 측위 알고리즘 초기화를 적용한 경우에 해당하며, 그렇지 않을 경우 특히 실시간 환경에서는 보정량의 급격한 변화로 인해 측위 알고리즘에 오류나 불안정이 발생할 수 있음을 유의해야 한다.

본 연구는 실시간 환경에서의 안정적이고 원활한 CLAS 메시지 활용을 위한 유의미한 결과를 도출하였으나, 보다 심층적인 분석을 위해 향후 연구에서는 1년 이상의 장기간 데이터를 활용한 분석과 CLAS 메시지 생성 기관에 따른 측위 성능 차이에 대한 추가 연구가 수행될 수 있을 것이다. 또한, 실시간 환경에서 MGF ID가 변경되더라도 필터를 초기화하지 않고 연속적인 정밀 위치 결정을 가능하게 하는 연구를 추가적으로 수행할 수 있을 것이다. 이러한 후속 연구를 통해 CLAS의 신뢰성과 실시간 환경에서의 활용성을 더욱 향상시킬 수 있을 것으로 기대된다.