Monitoring QZSS CLAS-based VRS-RTK Positioning Performance

Abstract

The Centimeter Level Augmentation Service (CLAS) is the Precise Point Positioning (PPP) - Real Time Kinematic (RTK) correction service utilizing the Quasi-Zenith Satellite System (QZSS) L6 (1278.65 MHz) signal to broadcast the Global Navigation Satellite System (GNSS) error corrections. Compact State-Space Representation (CSSR) corrections for mitigating GNSS measurement error sources such as satellite orbit, clock, code and phase biases, tropospheric error, ionospheric error are estimated from the ground segment of QZSS CLAS using the code and carrier-phase measurements collected in the Japan's GNSS Earth Observation Network (GEONET). Since the CLAS service begun on November 1, 2018, users with dedicated receivers can perform cm-level precise positioning using CSSR corrections. In this paper, CLAS-based VRS-RTK performance evaluation was performed using Global Positioning System (GPS) observables collected from the refence station, TSK2, located in Japan. As a result of performing GPS-only RTK positioning using the open-source software CLASLIB and RTKLIB, it took about 15 minutes to resolve the carrier-phase ambiguities, and the RTK fix rate was only about 41%. Also, the Root Mean Squares (RMS) values of position errors (fixed only) are about 4cm horizontally and 7 cm vertically.

1. INTRODUCTION

Real-Time Kinematic (RTK) 또는 Network RTK는 대표적인 Observation Space Representation (OSR) 방식의 반송파 보강항법시스템으로, 사용자는 기준국-사용자 간 측정치 차분을 통해 사용자 측정치에 포함된 Global Navigation Satellite System (GNSS) 오차를 제거하고, 반송파 미지정수 결정함으로써 정밀 위치 정보를 산출할 수 있다 (Park & Kee 2010, Park et al. 2014). 이러한 RTK 측위 기법은 수 초 이내에 cm급 고정밀 위치 정보를 획득할 수 있다는 장점이 있지만, 서비스 범위가 70 km - 100 km 수준으로 매우 제한적이고, 국가 단위의 광역 서비스를 위해서는 수 만에서 수십 만 bps의 bandwidth가 필요하므로 OSR 방식은 위성을 통한 보정정보 방송 서비스에 적합하지 않다. 반면, State Space Representation (SSR) 방식은 GNSS 오차를 분리하여 추정하고, 각 오차 성분에 대한 보정정보를 오차 특성에 맞춰 스케줄링 할 수 있으므로, OSR 대비 상대적으로 적은 bandwidth로 위성 기반 nation-wide 보정정보 방송 서비스가 가능하다 (Lim 2022). 현재 이용 가능한 SSR 보정정보 기반의 정밀 측위 서비스로는 Trimble RTX (Trimble 2022), u-blox PointPefect (U-blox 2022), Galileo High Accuracy Service (HAS) (Susi et al. 2021), Quasi-Zenith Satellite System (QZSS) Centimeter Level Augmentation System (CLAS) (Miya et al. 2015) 및 Multi-GNSS Advanced Demonstration tool for Orbit and Clock Analysis (MADOCA) (Harima et al. 2017) 등이 있다.

CLAS는 일본 내 위치한 다수의 GNSS 상시 기준국(Continuously Operating Reference Station; CORS)에서 수집된 GNSS 관측 데이터를 활용하여 GNSS 코드 및 반송파 측정치 오차 보상을 위한 보정정보가 포함된 Compact SSR (CSSR) 메시지를 생성하고 (Hirokawa et al. 2019), 이를 QZSS 위성 L6 신호를 통해 2,000 bps의 전송속도로 방송하는 nationwide Precise Point Positioning (PPP)/PPP-RTK 서비스로, 2018년 11월 1일부터 공식적으로 서비스를 시행하였다 (Hirokawa & Hernández 2020). CSSR은 SSR 보정정보의 효율적인 방송을 위해 개발된 공개 포맷으로 SSR 관련 국제 표준 포맷인 Radio Technical Commission for Maritime Services (RTCM) SSR (RTCM 2013)과 호환이 가능하다는 특징을 가지며, 현재 QZSS CLAS와 Galileo HAS가 CSSR 포맷을 활용하여 정밀 측위 서비스를 제공하고 있다.

CLAS 보정정보가 QZSS L6 신호를 통해 실시간으로 방송됨에 따라 이를 이용하여 다양한 분야에 적용하여 정밀 측위를 수행하고, 그 성능을 평가하기 위한 많은 연구들이 진행되고 있다. 기존 연구에서는 QZSS 위성을 통해 제공되는 CLAS 및 MADOCA 보정정보를 활용한 로봇 콤바인의 정밀 측위 시스템을 구성하고, 정적 및 동적 환경에 대해서 각 보정정보 서비스의 정밀 측위 성능을 비교 분석하였고 (Udompant et al. 2021), CLAS 기반 PPP-RTK 측위를 위해 자체 개발한 GNSS 수신기와 u-blox ZED-F9P와 같은 상용 수신기를 활용하여 동적 환경 하에서의 정밀 측위 성능을 평가하였다 (Hirokawa et al. 2019). 뿐만 아니라 일본 내 위치한 다수의 GNSS 상시 관측소에서 장시간 수집된 GNSS 관측 데이터를 활용하여 정적 환경에서의 CLAS 서비스의 측위 성능을 분석함으로써 QZSS가 제시하는 목표 성능을 만족하는 지 평가하였다 (Motooka et al. 2019). 이러한 기존 연구들은 QZSS CLAS 서비스의 성능을 측위 정확도 측면에서만 분석하였다.

따라서 본 논문에서는 일본 내 위치한 GNSS 기준국 관측 데이터에 CLAS 보정정보를 적용하여 수행된 RTK 기반 정밀 측위 결과를 통해 QZSS CLAS 서비스의 위치 정확도를 평가하였을 뿐만 아니라 보정정보와 함께 제공된 무결성 정보를 활용하여 사용자 위치 오차에 대한 보호수준을 산출하고, 이를 통해 QZSS CLAS 기반 정밀 측위 성능 감시를 수행하였다. 일반적으로 네트워크 보정정보 기반의 정밀 측위 서비스의 경우는 지상 기준국 네트워크 중심부에 위치한 사용자의 측위 성능이 가장 좋을 것으로 예상되므로 본 논문에서도 CLAS 네트워크에 중심부(일본 내)에 위치한 GNSS 기준국을 사용자로 가정하고, 해당 기준국에서 수집된 Global Positioning System (GPS) 관측 데이터에 CSSR 포맷의 CLAS 보정정보 적용 후, 반송파 측정치 잔차 분석을 통해 QZSS CLAS 보정정보의 반송파 미지정수해 도출 가능성을 분석하였다. 또한, CLAS 보정정보를 기반으로 사용자 위치에서의 가상 기준국 (Virtual Reference Station; VRS) 관측 데이터를 생성한 후, 이를 사용자 관측 데이터와 함께 RTK 측위에 특화된 오픈 소스 소프트웨어인 RTKLIB (RTKLIB 2022)에 입력하여 출력된 RTK 측위 결과를 바탕으로 CLAS 보정정보의 GNSS 측정치 오차 제거 성능을 평가하였다. 또한, 각 CLAS 보정정보에 대한 신뢰도를 의미하는 quality indicator를 활용하여 RTK 측위에 사용된 위성의 코드 및 반송파 측정치 잔차에 대한 공분산 행렬을 구성하였고, 이를 기반으로 사용자 위치 오차에 대한 보호수준을 산출하였다. 최종적으로 산출된 보호수준이 위치 오차를 적절하게 bound 하는지 여부를 확인함으로써 CLAS 기반 VRS-RTK 측위 성능 감시를 수행하였다.

2. CENTIMETER LEVEL AUGMENTATION SERVICE (CLAS)

2.1 Overview of CLAS Service

CLAS는 QZSS 위성의 L6 신호를 통해 GNSS 코드 및 반송파 측정치 오차 보상을 위한 보정정보를 일본 전역에 방송함으로써 서비스 영역 내 사용자의 PPP-RTK 측위를 가능케 하는 cm급 정밀 측위 서비스이다. CLAS는 GNSS Earth Observation Network (GEONET)이라 불리는 약 1,200개소의 일본 내 상시 기준국에서 GPS, Galileo, QZSS 등 다중 위성군 관측 데이터를 수집하고, 이를 중앙 처리국으로 전송한다. 중앙 처리국에서는 각 기준국에서 수집한 GNSS 측정치를 활용하여 위성 궤도 및 시계오차, 위성 신호 코드 및 반송파 바이어스, 대류층 및 전리층 오차에 대한 보정정보를 생성하고, 이를 QZSS 위성을 통해 CSSR 포맷으로 사용자에게 방송한다. Fig. 1은 QZSS CLAS 시스템 구조를 도식화한 것이다.

Fig. 1. System schematics of QZSS CLAS (Motooka et al. 2019).

CLAS 서비스 영역의 크기는 대략 800,000 km2이며 (Hirokawa et al. 2019), Fig. 2와 같이 12개의 네트워크와 212개의 격자점(grid)으로 구성되어 있고, 인접한 격자점 간의 간격은 약 60 km 정도로 설정되어 있다 (Cabinet Office 2021). 현재 CLAS는 GPS, Galileo, QZSS를 포함한 다중 위성군에 대한 보정정보 및 무결성 정보를 방송하고 있으며, 향후 GLONASS CDMA 신호 (G1a, G2a)를 지원할 예정이다. Table 1은 Performance Standard (PS)-QZSS 문서 (Cabinet Office 2022)에 정의된 CLAS 서비스의 PPP-RTK 측위 성능, 보정 신호 등을 정리한 것이다.

Fig. 2. Local network and grid distribution of CLAS.

Table 1. Specifications of CLAS (Cabinet Office 2022).

Parameter		Specification
Augmented GNSS signals		GPS L1 C/A, L1C, L2P, L2C, L5 QZSS L1 C/A, L1C, L2C, L5 Galileo E1b, E5a
Accuracy	Static	Horizontal ≤ 6 cm (95%), Vertical ≤12 cm (95%)
	Kinematic	Horizontal ≤12 cm (95%), Vertical ≤24 cm (95%)
Availability		≥0.99 (at least 3 QZS satellites)
Continuity		≥1-2×10-4/hour (each QZS satellite)
Integrity (Time to alert)		9.2 sec
Time to first fix		≤60 sec (95%)

 

2.2 Compact SSR Message Format

Fig. 3은 QZSS L6 메시지 구조를 나타낸 것이며, 각 메시지는 49 bits header, 1695 bits data, 256 bits Reed-Solomon, 총 2000 bits로 구성되어 있다. 또한, data 파트에는 CSSR의 sub-type 메시지가 포함되어 있으며, 각 sub-type 메시지는 표준 RTCM 3 메시지와 유사한 구조를 갖는다.

Fig. 3. QZSS L6 message structure (Cabinet Office 2021).

CSSR 메시지는 기존 RTCM 3 포맷이 호환이 가능한 메시지로 proprietary message type 4073를 통해 방송되고 있으며, Table 2와 같이 GNSS 측정치 오차 보정정보 또는 무결성 정보 관련 파라미터가 포함된 12개의 sub-type 메시지가 존재한다. 또한, CSSR 보정정보는 크게 모든 네트워크 및 격자점에서 동일하게 적용할 수 있는 common correction과 각 네트워크 또는 격자점마다 제공되는 local correction으로 구분할 수 있다. 먼저, common correction은 각 위성에 대한 궤도 및 시계 오차, 신호 코드 및 반송파 바이어스, User Range Accuracy (URA) 정보로 구성되어 있으며, 이는 RTCM SSR 1단계 및 2단계 메시지와 호환이 가능한 형태로 제공된다 (Hirokawa et al. 2019). CSSR은 RTCM SSR와 마찬가지로 다중 위성군을 지원하며, 효율적인 메시지 전송을 위해 RTCM Multiple Signal Messages (MSM)와 유사하게 satellite mask를 활용하고 있다. Local correction은 서비스 영역을 구성하는 각 네트워크 및 격자점에 대한 대류층 오차 및 전리층 오차로 구성되어 있다. 전리층 오차 보정정보는 네트워크 단위로 생성되는 area correction과 가상의 격자점마다 생성되는 residual correction으로 구성되어 있다. 특히, area correction은 위도 및 경도의 함수로 표현되는 다항식 (polynomial)의 형태로 제공된다. 대류층 오차 보정정보는 모든 위성에 대하여 동일하게 적용할 수 있는 천정 방향 지연 값을 격자점별로 제공한다.

Table 2. Compact SSR Messages (Cabinet Office 2021).

Group name	Sub-group name	Message type ID, sub-type ID
Common correction	Compact SSR mask Compact SSR orbit correction Compact SSR clock correction Compact SSR satellite code bias Compact SSR satellite phase bias Compact SSR satellite code and phase bias Compact SSR URA Compact SSR combined correction	4073, 1 4073, 2 4073, 3 4073, 4 4073, 5 4073, 6 4073, 7 4073, 11
Local correction	Compact SSR STEC correction Compact SSR gridded correction Compact SSR atmospheric correction	4073, 8 4073, 9 4073, 12
Service	Compact SSR service information	4073, 10

3. RANGE ERROR CORRECTION

CSSR은 GNSS 측정치에 포함된 오차 요소의 특징에 따라 보정 정보의 종류가 구별되고 각 오차 성분에 대하여 모델링 된 파라미터를 제공한다. 오차 성분별로 모델링하기 때문에 오차 요소의 변화율을 고려하여 전송주기 설정하며, 다른 오차 성분에 비해 시변화율이 큰 위성 시계 오차에 대한 보정정보의 전송주기는 약 5초인 반면, 위성 궤도 오차 및 바이어스, 대류층 및 전리층 오차 등 비교적 천천히 변하는 오차 성분에 대한 보정 정보의 전송 주기는 약 30초이다.

3.1 Satellite Orbit and Clock Corrections

현재 CLAS에서 제공하는 위성 궤도 오차에 대한 보정정보는 sub-type IDs 2, 11에 포함되어 있다. Sub-type 2 메시지는 위성 궤도 오차 성분, 방송 궤도력 Issue of Data, Ephemeris (IODE) 등을 포함하고 있으며, sub-type 11 메시지는 위성 궤도 오차 보정정보 뿐만 아니라 시계 오차 보정정보를 함께 제공한다. CSSR은 RTCM SSR과의 호환성을 중점을 두고 있기 때문에 위성 궤도 오차 보정정보를 Fig. 4와 같이 radial, along-track, cross-track 방향 성분으로 제공하고 있다. 따라서 radial, along-track, cross-track 방향의 궤도 오차 성분을 지구중심 지구고정 좌표계 (Earth-Centered Earth-Fixed; ECEF)로 변환하는 과정이 필요하다 (Cabinet Office 2021). 또한, CLAS는 사용자가 항법 메시지에 포함된 파라미터를 이용하여 계산된 위성 시계에 대한 오차를 보상하기 위한 보정정보를 제공하며, 이는 CSSR 메시지 sub-type IDs 3, 11에 포함되어 있다. Sub-type 3 메시지에는 GPS, GLONASS, Galileo, QZSS 위성 시계 오차 보정정보가 포함되어 있으며, 사용자는 이를 이용하여 항법 메시지 기반으로 산출된 위성 시계에 대한 오차를 보상할 수 있다 (Cabinet Office 2021).

Fig. 4. Radial, along-track and cross-track components of satellite orbit.

3.2 Satellite Code and Phase Biases

CLAS를 포함한 대부분의 SSR 기반 보강항법시스템은 위성 궤도 및 시계오차 보정정보 외에도 코드 및 반송파 측정치에 포함된 바이어스에 대한 보정정보 또한 함께 제공한다. CSSR 메시지 sub-type 4 및 5 메시지는 각각 위성 코드 바이어스와 반송파 바이어스를 포함하고 있으며, sub-type 6은 위성 코드 및 반송파 바이어스를 함께 제공하는데, 사용자는 각 위성 바이어스에 대한 flag를 통해 위성 코드 및 반송파 바이어스 보정정보의 제공 여부를 확인할 수 있다. 일반적으로 sub-type 6 메시지의 갱신 주기와 sub type 2 메시지에 포함된 위성 궤도 오차 보정정보의 갱신 주기는 서로 다르기 때문에 임의의 epoch t와 t-1에서 SSR 보정정보의 불연속이 발생할 수 있다. 이러한 보정정보 불연속은 사용자 반송파 측정치의 급격한 변화를 유발하고, 이로 인해 편향된 사용자 위치 결과가 산출될 수 있다.

3.3 Ionospheric Correction

CLAS는 전리층 오차 보정정보를 다양한 주파수 측정치에 적용 가능하도록 Slant Total Electron Content (STEC) 형태로 제공하고 있다. 또한, 앞서 언급한 바와 같이 CSSR 전리층 오차 보정정보는 각 네트워크별로 생성되는 polynomial 성분과 사전 정의된 격자점을 기준으로 생성된 residual 성분으로 구분된다. 전리층 오차 보정정보의 다항식 성분은 CSSR 메시지 sub-type 8, 잔차 성분은 sub-type 9에 포함되어 있으며, sub-type 12 메시지의 경우는 다항식 성분과 잔차 성분을 함께 제공한다. STEC 보정정보를 계산하기 위해 필요한 다항식 계수는 sub-type 8 또는 12 메시지에 포함된 “STEC Correction Type” 파라미터에 의해 결정된다. CLAS 사용자가 코드 및 반송파 측정치 전리층 오차를 보상하기 위해서는 사용자를 둘러싼 최소 3개 이상의 격자점이 필요하며, 각 격자점에서의 전리층 보정정보들의 가중치 합으로 사용자 측정치의 전리층 오차에 대한 보정정보를 산출한다.

3.4 Tropospheric Correction

CSSR 대류층 오차 보정정보는 격자점을 기준으로 제공되므로 보정정보 모델링 파라미터는 전리층 오차 보정정보와 함께 sub-type 9, 12 메시지에 포함되어 있다. CLAS가 제공하는 건조 및 습윤 대류층 오차 보정정보는 상수 (constant) 값과 변동 (variation) 값으로 구분되어 제공되며, 건조 지연에 대한 상수 값으로는 2.3 m, 습윤 지연의 경우는 sub-type 9 메시지에 정의된 0.252 m 또는 sub-type 12에 포함된 “Troposphere Residual Offset” 값을 사용한다. 전리층 오차 보정정보와 마찬가지로 사용자를 둘러싼 3개 또는 4개의 격자점 대류층 오차 보정정보를 보간하여 사용자 위치에서의 천정 방향 건조 및 습윤 지연 값을 계산할 수 있다. 그 다음, 각 위성에 대한 대류층 오차 사상함수 (mapping function)를 산출된 천정 방향 대류층 오차에 적용함으로써 최종적으로 경사 대류층 지연 오차 보정정보를 산출할 수 있다.

3.5 User Dependent Error Correction

IS-QZSS 문서에 따르면 CLAS는 사용자 위치에 영향을 주는 고체지구 조석(solid earth tide), 해양 조석 (ocean tide), 극 조석 (pole tide)과 반송파 측정치의 위상 말림 현상 (phase windup effect), 안테나 위상 중심 오프셋 (phase center offset) 및 안테나 위상 중심 변동 (phase center variation) 등 사용자 위치 및 안테나 모델에 따라 달라질 수 있는 오차 성분에 대한 보정정보는 CSSR 메시지를 통해 제공하고 있지 않으므로, 사용자는 자신의 위치 정보를 반영하여 이러한 오차를 보상해야 한다 (Cabinet Office 2021). 본 논문에서 코드 및 반송파 측정치 오차 보정정보 생성 시 적용한 오차 모델은 다음과 같다.

    \(\cdot\) Phase windup
    \(\cdot\) Shapiro delay (Relativistic path range effect)
    \(\cdot\) Solid earth tide
    \(\cdot\) Ocean tide
    \(\cdot\) Pole tide

4. INTEGRITY INFORMATION

4.1 Compact SSR Integrity

CLAS는 GNSS 측정치 오차 보정정보 이외에도 alert flag와 quality indicator로 구성된 무결성 정보를 함께 제공하고 있다. 1-bit alert flag는 L6 메시지 header에 포함되어 있으며, 이 값이 1로 설정된 경우, 사용자는 L6 신호를 통해 방송되는 CSSR 보정정보를 사용하지 말아야 한다. CSSR 메시지에는 오차 보정정보 뿐만 아니라 보정정보의 불확실성 (uncertainty)을 의미하는 quality indicator를 함께 제공하고 있다. CSSR 메시지 sub-type ID 7은 각 위성에 대한 URA 값이 포함되어 있으며, 대류층 오차 보정정보에 대한 quality indicator는 sub-type IDs 9 및 12, 전리층 오차 보정정보에 대한 quality indicator는 sub-type IDs 8 및 12를 통해 제공된다. 각 보정정보에 대한 quality indicator는 식 (1)과 같이 3-bit CLASS와 3-bit VALUE의 조합으로 표현된다 (Cabinet Office 2021).

\( \rm Quality~Indicator~[mm]≤3^{CLASS}\left(1+{VALUE \over 4}\right)-1~[mm]\)                                             (1)

CSSR URA 값과 대류층 오차 보정정보에 대한 quality indicator는 식 (1)을 이용하여 변환 가능하며, 전리층 오차 보정정보의 경우는 IS-QZSS 문서에 수록된 CSSR STEC quality indicator와 STEC 불확실성 간의 관계표를 활용하여 quality indicator를 산출한다. 사용자가 CSSR 보정정보를 적용한 후 측정치에 남아 있는 잔여 오차는 식 (2)와 같이 모델링 된다 (Cabinet Office 2021).

\(\sigma_{CLAS,s}=\sqrt{(\sigma_{s,user} )^2+(\sigma_{s,sis}/10)^2+\left({40.3 \over f^2}  \sigma_{s,iono}×100\right)^2+(\sigma_{s,tropo}/10/\sin⁡ E_s )^2} \)                   (2)

여기서 \(E_s\)는 위성 앙각, \(\sigma_{s,user}\)는 수신기 잡음, 다중경로오차 등 사용자 수신기 및 환경에 의한 오차, \(\sigma_{s,sis}\)는 위성 궤도 및 시계오차, 신호 바이어스의 불확실성을 의미하는 Signal In Space (SIS) 오차를 의미한다. 또한, \(\sigma_{s,iono}\)는 STEC 보정정보 quality indicator로부터 계산된 전리층 오차 보정정보에 대한 불확실성, \(\sigma_{s,tropo}\)는 대류층 오차 보정정보에 대한 불확실성을 의미한다.

CLAS 사용자는 식 (2)를 통해 계산한 오차 값을 이용하여 CSSR 보정정보 적용 후의 측정치에 대한 innovation 및 residual 검사 또는 산출된 위치 결과에 대한 다양한 무결성 검사를 수행할 수 있다.

4.2 Protection Level

최근 위성 또는 인터넷을 기반으로 PPP 및 PPP-RTK 측위를 위한 보정정보를 방송하는 다양한 정밀 측위 서비스가 등장함에 따라 이러한 서비스를 드론 및 차량과 같은 이동체의 정밀 항법 시스템 구축에 활용하려는 연구가 다각도로 이루어지고 있다. 특히, 무인비행체 또는 자율주행차량의 경우는 GNSS 기반 항법 성능이 안전과 직결되기 때문에 정확한 위치 정보뿐만 아니라 항법해의 신뢰도를 의미하는 파라미터 중 하나인 보호수준 (Protection Level; PL)을 함께 제공해주는 것이 매우 중요하다. 이에 따라 Trimble, u-blox 등과 같은 GNSS 시장의 선두 기업들은 자사 RTK-GNSS 수신기 및 정밀 측위 서비스 (Trimble RTX, u-blox Point Perfect) 활용 시, cm급 정밀 항법해와 그에 대한 보호수준을 함께 제공하고 있다. 위성을 통해 GNSS 오차 보정정보를 무료로 제공하는 PPP-RTK 측위 서비스 중 하나인 QZSS CLAS는 현재 사용자가 직접 산출된 항법해의 무결성 검사를 수행할 수 있도록 식 (2)와 같이 CSSR 보정정보의 표준편차 값을 제공하고 있지만, 이를 활용한 사용자 보호 수준 산출 기법은 IS-QZSS 문서 상에 명시되어 있지 않다. 이에 본 논문에서는 SBAS 사용자 보호 수준 계산 식을 기반으로 CLAS 사용자 보호 수준을 계산함으로써 CSSR 보정정보 적용에 따른 측위 오차가 적절하게 bound 되고 있는 지 여부를 확인하였다.

Satellite-Based Augmentation System (SBAS)는 GPS 코드 측정치 오차 보상을 위한 보정정보를 위성을 통해 방송하는 광역보강시스템이다 (Lim et al. 2018, 2020). SBAS는 QZSS CLAS와 동일하게 각 오차 성분별 SSR 보정정보를 제공한다. 또한, SBAS는 보정정보와 함께 해당 보정정보의 confidence bound를 의미하는 표준편차를 제공하며, SBAS 보정정보 적용 후 코드 측정치 잔차의 분산 값 (\(\sigma_{SBAS,i}^2\))은 식 (3)과 같이 계산할 수 있다 (RTCA 2006).

\(\sigma_{SBAS,i}^2=\sigma_{flt,i}^2+\sigma_{UIRE,i}^2+\sigma_{tropo,i}^2+\sigma_{air,i}^2\)                                                               (3)

여기서 \(\sigma_{flt,i}^2\)는 SBAS 위성 궤도 및 시계오차 보정정보에 대한 분산, \(\sigma_{UIRE,i}^2\)는 SBAS 전리층 오차 보정정보에 대한 분산, \(\sigma_{tropo,i}^2\)는 SBAS 대류층 오차 보정정보에 분산, \(\sigma_{air,i}^2\)는 수신기 잡음 및 다중경로 오차에 대한 분산을 의미한다. SBAS 사용자는 식 (3)을 통해 계산된 가시 위성의 코드 측정치 분산 값을 이용하여 보호수준 (Protection Level; PL)을 계산함으로써 SBAS 보정정보 적용 후 산출된 위치 정보에 대한 가용성을 판단할 수 있다. SBAS 시스템에서 수평 보호수준 (Horizontal Protection Level; HPL) 및 수직 보호수준 (Vertical Protection Level; VPL)은 식 (4)와 (5)를 이용하여 계산할 수 있다 (RTCA 2006).

\(HPL=K_h \sigma_{major}\)                                                                                    (4)

\(VPL=K_v d_v \)                                                                                            (5)

with

\(K_v=5.33, \left\{ \begin{array}. K_h=6.18& {\rm for~en~route~through~LNAV} \\ K_h=6.00& {\rm for~LNAV, VNAV, LP, LPV~approach} \end{array}\right.\)

where

\(\mathbf{S}=(\mathbf{G}^T \mathbf{WG})^{-1} \mathbf{G}^T \mathbf{W}= \left[ \begin{array}{cccc}s_{e,1}&s_{e,2}&⋯&s_{e,k}\\s_{n,1}&s_{n,2}&⋯&s_{n,k}\\s_{u,1}&s_{u,2}&⋯&s_{u,k}\\s_{t,1}&s_{t,2}&⋯&s_{t,k} \end{array}\right],\)

\(\mathbf{W}^{-1}=\left[ \begin{array}{cccc}\sigma_{SBAS,1}^2&0&⋯&0\\0&\sigma_{SBAS,2}^2&⋯&0\\⋮&⋮&⋱&0\\0&0&⋯&\sigma_{SBAS,k}^2 )\end{array}\right]\)

\(\mathbf{G}=\left[ \begin{array}{ccc}-\cos ⁡E_1 \sin⁡ A_1 &-\cos ⁡E_1 \cos ⁡A_1 &-\sin⁡ E_1 &1\\-\cos ⁡E_2 \sin⁡ A_2 &-\cos ⁡E_2 \cos ⁡A_2 &-\sin⁡ E_2 &1\\⋮&⋮&⋯&⋮\\-\cos ⁡E_k \sin⁡ A_k &-\cos ⁡E_k \cos ⁡A_k &-\sin⁡ E_k &1\end{array}\right]\)

\(d_e^2=\sum\limits_{i=1}^k s_{e,i}^2  \sigma_{SBAS,i}^2,~d_n^2=\sum\limits_{i=1}^k s_{n,i}^2  \sigma_{SBAS,i}^2,~d_en=\sum\limits_{i=1}^k s_{e,i} s_{n,i} \sigma_{SBAS,i}^2~,d_v^2=\sum\limits_{i=1}^k s_{u,i}^2  \sigma_{SBAS,i}^2\)

\(\sigma_{major}=\sqrt{ {d_e^2+d_n^2 \over 2}+\sqrt{\left({d_e^2-d_n^2 \over 2}\right)^2+d_en^2}} \)

여기서 \(\mathbf{S}\)는 측정치 영역에서 위치 영역으로의 투영 행렬 (projection matrix), \(\mathbf{G}\)는 기하 행렬 (geometry matrix), \(\mathbf{W}\)는 가중치 행렬 (weighting matrix)를 의미한다.

5. TEST RESULTS

5.1 Test Configurations

QZSS CLAS 기반 VRS-RTK 측위 성능 평가 실험을 위해 Fig. 5와 같이 일본 Tsukuba에 위치한 International GNSS Service (IGS) 기준국 (TSK2)을 정적 사용자로 선정하였다. 실험에 사용할 사용자 GNSS 측정치는 UTC 기준 2022년 6월 15일 07시부터 약 3시간동안 TSK2 기준국에서 수집된 관측 데이터를 MIRAI (Go! GNSS 2022) NTRIP 서버를 통해 RTCM 포맷으로 취득하였고, 이를 Receiver INdependent Exchange (RINEX) 포맷으로 변환하여 후처리 과정에 사용하였다. 또한, 사용자 관측 데이터 날짜에 방송된 QZSS L6 메시지는 QZSS CLAS 웹사이트 (QZSS 2022)를 취득하였으며, 이를 자체 개발한 MATLAB 기반 프로그램으로 디코딩하여 CSSR 사용자 보정정보를 생성하였다.

Fig. 5. QZSS CLAS user selection (IGS 2022).

5.2 Range Domain Analysis

측정치 영역 성능 분석에서는 사용자 반송파 측정치에 CSSR 보정정보를 적용한 후 위성 간 단일 차분 측정치의 잔차를 분석하였으며, 측정치 잔차 계산 과정은 다음과 같다. 먼저, 각 위성에 대한 CSSR 보정정보를 측정치 영역에서 재조합하여 Carrier-Phase Correction (CPC)를 생성한 다음, 식 (6)과 같이 CPC를 사용자 반송파 측정치에 적용하여 오차를 보상한 후, 위성 간 측정치 단일 차분(\(\triangledown\))을 수행한다.

\(\triangledown\phi_f+CPC_f-\triangledown d=\triangledown( \delta T-\delta I_f + ε_{\phi_f} )+\triangledown N_f \lambda_f\)                                                 (6)

단일 차분된 측정치에는 반송파 미지정수가 남아 있으므로, 이를 제거해야만 잔차를 계산할 수 있다. 여기서, CSSR 보정정보를 적용 후의 대류층 지연 오차의 residual(\(\delta T\)) 및 전리층 지연 오차의 residual(\(\delta I\))과 반송파 측정치 잡음(\(ε_{\phi_f}\))이 Gaussian 분포를 따른다고 가정하면, 반송파 측정치의 cycle slip이 발생하지 않은 구간에 대하여 단일 차분 반송파 측정치의 평균을 취한 후, 반올림하여 후처리 미지정수(\(\triangledown\hat{N}\))를 결정할 수 있다.

\(∇\hat{N}_f={\rm round}\left ( {∇\phi_f-∇d \over \lambda_f} \right )\)                                                                               (7)

식 (7)을 통해 결정된 후처리 미지정수를 식 (8)와 같이 단일 차분 측정치에 적용함으로써 최종적으로 CSSR 보정정보 적용에 따른 사용자 반송파 측정치 잔차를 산출하였다.

\(\triangledown_{res_{Φ_f}}=(\triangledown Φ_f+CPC_f-\triangledown d)-\triangledown \hat{N}_f \lambda_f \)                                                                (8)

Fig. 6은 GPS L1/L2 반송파 측정치 단일 차분 잔차의 시계열을 도시한 것이며, 각 그래프에서 파란색 실선은 반송파 파장 길이의 절반을 의미한다. 단일 차분 잔차 계산 시, 전체 구간 중 앙각이 가장 높게 유지되었던 GPS PRN 2 위성을 기준 위성으로 선정하였다. 일반적으로 측정치 영역 성능 검증에서는 반송파 측정치 잔차 값이 파장의 1/2을 초과하는 경우는 미지정수가 잘못 결정될 확률이 매우 높으므로 계산된 보정정보 적용 후의 잔차 값이 임계 값(파장 1/2)을 초과하는 지 여부를 분석함으로써 해당 보정정보의 미지정수해 도출 가능성을 예측해볼 수 있다 (Song 2016, Kim 2019). Fig. 6의 결과를 보면 실험 구간 전반에 걸쳐 대부분의 잔차 값이 파장의 절반 이하로 유지됨을 확인할 수 있는데, GPS PRN 5, 7, 13, 15 등 일부 위성의 경우에는 잔차가 파장 1/2을 초과하는 구간이 다수 발생하였음을 확인하였다. Fig. 7은 위성 앙각에 따른 반송파 측정치 단일 차분 잔차를 도시한 것인데, 두 그래프의 결과를 보면 L1/L2 측정치 모두 위성 앙각이 약 30° 수준임에도 불구하고, 잔차 값이 임계 값을 초과하는 위성이 다수 존재함을 알 수 있다. 뿐만 아니라 50° 이상의 고앙각 위성임에도 불구하고 잔차의 최대 값이 약 0.12 m (L1), 0.14 m (L2) 수준인 경우가 발생하는 것을 확인하였다. 이상의 결과를 종합해볼 때, 사용자 측정치에 CSSR 보정정보 적용 후에도 제거되지 않고 남아 있는 오차가 존재함을 확인할 수 있고, 이는 반송파 미지정수 결정 성능에 상당한 영향을 미칠 것으로 예상된다.

Fig. 6. Time series of GPS single-differenced phase residuals, (a) L1, (b) L2.

 

Fig. 7. Elevation angle vs. GPS single-differenced phase residuals. (a) L1, (b) L2.

5.3 Position Domain Analysis

앞서 언급한 바와 같이 본 논문에서는 RTKLIB 프로그램을 이용하여 CSSR 보정정보의 VRS-RTK 측위를 수행하고, 그 결과를 바탕으로 위치 오차, 반송파 미지정수 결정율 등 RTK 측위 성능 지표를 분석함으로써 일본 내 위치한 사용자에 대한 QZSS CLAS의 성능을 평가하였다. 이를 위해서는 정밀 상대 측위를 위한 기준국 RINEX 관측 데이터 파일이 필요한데, 본 논문에서는 RTKLIB 기반의 오픈 소스 소프트웨어인 CLAS library (CLASLIB) (Motooka et al. 2019)를 이용하여 사용자(TSK2) 위치에서의 GNSS 가상 관측 데이터를 생성하였다. CLASLIB의 SSR2OBS 프로그램은 SSR 보정정보를 OSR 보정정보로 변환하는 모듈로써 QZSS CLAS L6 메시지를 디코딩하여 CSSR 보정정보를 추출하고, 이를 활용하여 사용자에 의해 입력된 3차원 좌표에 대한 가상 측정치를 생성한 후, RINEX 포맷으로 출력한다. Fig. 8은 CLASLIB 및 RTKLIB 프로그램을 활용한 VRS-RTK 측위 과정을 도식화한 것이다.

Fig. 8. VRS-RTK positioning process using CLASLIB and RTKLIB software.

앞서 언급한 바와 같이 본 연구에서는 SBAS 사용자 보호 수준 산출 기법을 CLAS 사용자에 적용하였는데, RTK 측위 과정에서는 위성 간 단일 차분 기법을 활용하므로 식 (5)의 geometry matrix(\(\mathbf{G}\))와 weighting matrix(\(\mathbf{W}\))를 식 (9)와 (10) 같이 변경하여 CLAS 사용자 보호 수준을 산출하였다.

\(\mathbf{G}_{CLAS}=\left[ \begin{array}{cccc} \left(\mathbf{G}_{(1,1)}-\mathbf{G}_{(2,1)} \right)&\left(\mathbf{G}_{(1,2)}-\mathbf{G}_{(2,2)} \right)&\left(\mathbf{G}_{(1,3)}-\mathbf{G}_{(2,3)} \right)&0\\\left(\mathbf{G}_{(1,1)}-\mathbf{G}_{(3,1)} \right)&\left(\mathbf{G}_{(1,2)}-\mathbf{G}_{(3,2)} \right)&\left(\mathbf{G}_{(1,3)}-\mathbf{G}_{(3,3)} \right)&0\\⋮&⋮&⋯&⋮\\\left(\mathbf{G}_{(1,1)}-\mathbf{G}_{(k,1)} \right)&\left(\mathbf{G}_{(1,2)}-\mathbf{G}_{(k,2)} \right)&\left(\mathbf{G}_{(1,3)}-\mathbf{G}_{(k,3)} \right)&0 \end{array}\right]\)                               (9)

\(\mathbf{W}_{CLAS}^{-1}=\left[ \begin{array}{cccc}\sigma_{CLAS,1}^2+\sigma_{CLAS,2}^2&\sigma_{CLAS,1}^2&⋯&\sigma_{CLAS,1}^2\\ \sigma_{CLAS,1}^2&\sigma_{CLAS,1}^2+\sigma_{CLAS,3}^2&⋯&\sigma_{CLAS,1}^2\\ ⋮&⋮&⋱&⋮\\ \sigma_{CLAS,1}^2&\sigma_{CLAS,1}^2&⋯&\sigma_{CLAS,1}^2+\sigma_{CLAS,k}^2 \end{array}\right]\)                (10)

Fig. 9는 GPS L1/L2/L5 측정치를 이용한 VRS-RTK 측위 결과 (Case 1)와 사용자 (TSK2)로부터의 기선 거리 약 16 km인 지점에 위치한 ISHI 기준국의 관측 데이터를 활용하여 단일 기준국 RTK 측위를 수행한 결과 (Case 2)의 사용자 측위 오차를 도시한 것이며, 측위 오차 계산 시, IGS에서 제공하는 Solution INdependent Exchange (SINEX) 데이터 (igs22P22143.snx)에 포함된 TSK2 기준국의 3차원 좌표를 참값으로 활용하였다. Table 3은 반송파 미지정수가 결정된 RTK 고정해 (fixed solution)에 대한 수평 및 수직 위치 오차의 통계치, Time-To-First-ambiguity-Fix (TTFF) 및 RTK fix rate를 정리한 것이다. 여기서, TTFF는 반송파 미지정수가 결정되어 RTK 고정해가 산출되기까지 소요되는 시간을 의미하는데, CLAS 기반 VRS-RTK와 단일 기준국 RTK 각각 약 900초, 1초가 소요되었다. 측위 정확도 비교 결과, 95% 오차는 Case 1과 Case 2 각각 수평 방향 5.1 cm, 6.9 cm, 수직 방향 7.1 cm, 11.1 cm로 SSR 보정정보의 오차 제거 성능이 OSR 대비 다소 떨어지기는 하지만 현격한 차이를 보이지는 않는 것을 확인하였다. 다만, 전체 구간 중 RTK 고정해가 출력된 시간의 비율을 의미하는 RTK fix rate, 즉 반송파 미지정수 결정 성능은 CLAS의 경우, 약 41% 수준으로 전통적인 OSR 보정정보 기반의 단일 기준국 RTK 결과 대비 현저히 떨어지는 것으로 나타났다.

Fig. 9. Single base RTK vs. CLAS-based VRS-RTK, (a) horizontal, (b) vertical.

Table 3. RTK positioning error statistics (fixed only), TTFF and fix rate.

		Error statistics		TTFF [s]	Fix rate [%]
		Horizontal	Vertical
CLAS-based VRS-RTK (Case 1)	RMS (m) Mean (m) 95% (m)	0.0402 0.0341 0.0686	0.0670 -0.0335 0.1113	897	40.6
Single base RTK (Case 2)	RMS (m) Mean (m) 95% (m)	0.0268 0.0231 0.0510	0.0349 0.0081 0.0710	1	100

 

Fig. 10은 CLAS 보정정보 기반 VRS-RTK 측위 결과의 수평 및 수직 오차와 보호수준을 함께 도시한 것이며, 붉은색 그래프와 파란색 그래프는 각각 RTK 유동해 (float solution)와 고정해의 수평 위치 오차 (Horizontal Position Error; HPE) 및 수직 위치 오차 (Vertical Position Error; VPE)이고, 초록색 그래프는 RTK 고정해에 대한 수평 및 수직 보호수준을 의미한다. 실험 전체 구간에서 RTK 고정해의 위치 오차가 보호수준에 의해 적절히 bound되는 것을 확인할 수 있지만, 붉은 색 원으로 표시된 구간 (GPS TOW 288375s – 288636s)에서는 미지정수가 결정되었음에도 불구하고 수평 및 수직 보호 수준이 각각 약 5 m, 8 m 수준으로 다른 epoch 대비 최대 약 8배가량 큰 것을 확인하였다. 해당 문제에 대한 원인을 파악하기 위해 Fig. 11과 같이 측위 전 구간에서의 Position Dilution of Precision (PDOP)과 RTK 측위에 사용된 위성 수 시계열을 분석하였다. 먼저, 위 문제가 발생한 구간에서의 PDOP 값이 약 14로 정상적인 구간 대비 최대 7배 정도 큰 수준으로 나타났다. 해당 구간에서의 PDOP 값이 급격하게 증가한 것은 RTK 고정해 산출에 사용된 위성 수가 급감하여 위성의 기하학적 배치가 매우 좋지 않았던 것이 주요 원인으로 파악된다. 위성의 기하학적 배치는 보호 수준에도 상당한 영향을 끼치므로 지나치게 높은 PDOP으로 인해 미지정수가 정상적으로 결정된 구간 대비 비정상적으로 큰 보호 수준이 산출된 것으로 예상할 수 있다.

Fig. 10. CLAS-based VRS-RTK position errors and protection levels, (a) horizontal, (b) vertical.

 

Fig. 11. PDOP and number of satellites.

이와 같이 위성 수가 급격히 감소한 정확한 원인 (CLAS 보정정보 제공 여부, 반송파 측정치 이상 여부 등)을 RTKLIB 소프트웨어의 debug trace 기능을 통해 분석한 결과, 비정상 구간 시점 (GPS TOW 288375s)에서 RTKLIB 소프트웨어에 내장된 단일 차분 geometry-free 반송파 측정치 변화율 기반 cycle slip 검출 알고리즘에 의해 GPS PRN 7, 9, 30 위성의 반송파 측정치에 cycle slip flag가 발생한 것으로 확인되었다. 다만, 이러한 cycle slip flag 발생 원인이 CLAS 보정정보 적용 여부 때문인지, 아니면 반송파 측정치 자체에 이상 현상이 발생한 것인지 확인을 위해 CLAS 보정정보를 적용하기 전후의 반송파 측정치를 활용하여 테스트를 진행해보았다. 그 결과, CLAS 보정정보를 적용하지 않은 경우에는 해당 epoch에서 cycle slip flag가 발생하지 않았다. 즉, CLAS 보정정보의 적용으로 인해 단일 차분 geometry-free 반송파 측정치 변화율이 임계 값을 초과하여 GPS PRN 7, 9, 30 위성에 대한 cycle slip flag가 발생하였고, 결과적으로 해당 위성들이 RTK 고정해 산출 과정에서 제외됨으로써 위성 수가 급격히 감소하고, 비정상적으로 큰 PDOP 값이 산출된 주요 원인으로 작용하였음을 확인하였다.

6. CONCLUSIONS

본 논문에서는 일본 내 위치한 사용자에 대하여 현재 QZSS 위성 L6 신호를 통해 방송 중인 CLAS 서비스의 측정치 영역 및 위치 영역 성능 분석을 수행하였다. 이 과정에서 CSSR 보정정보 적용에 의한 사용자 측정치 오차 제거 성능을 보다 명확하게 확인하기 위해 RTK 정밀 측위 시, 다중 위성군이 아닌 GPS 단일 위성군 측정치만을 활용하였다. 또한, SBAS 사용자 보호 수준 산출 기법을 활용하여 CSSR 보정정보 적용 후 산출된 정밀 위치 정보에 대한 보호 수준을 계산하고, 계산된 보호 수준이 측위 오차를 적절하게 bound 하는지 여부를 확인하였다.

먼저, CSSR 보정정보의 오차 보정 성능을 확인하기 위해 사용자 반송파 측정치에 CSSR 기반 CPC 보정정보 적용 후, 위성 간 단일 차분 측정치의 잔차를 분석하였다. 이 과정에서 반송파 미지정수는 반송파 측정치의 cycle slip이 발생하지 않은 구간에 대하여 시간 평균을 취한 후, 그 값을 반올림하여 결정하였다. GPS L1/L2 측정치 영역 분석 결과, 대부분의 위성 단일 차분 잔차는 반송파 파장의 절반 이하였으나, 일부 위성은 앙각이 약 30° 임에도 불구하고, 잔차 값이 파장의 1/2을 초과하는 경우가 다수 발생하였다. 뿐만 아니라 50° 이상의 고앙각 위성의 경우에도 보정정보 적용 후 측정치 잔차가 약 10 cm 이상 발생함을 확인하였다.

다음으로 위치 영역 성능 분석에서는 CSSR 보정정보를 활용하여 사용자 위치에서의 가상 기준국 관측 데이터를 생성하고, RTKLIB 오픈 소스 프로그램 기반 VRS-RTK 측위를 수행하였다. GPS L1/L2/L5 측정치를 활용한 정밀 측위 결과, CLAS 기반 VRS-RTK 위치해의 95% 오차는 수평 및 수직 각각 약 7 cm, 11 cm 수준으로 단일 기준국 RTK에 비해 다소 떨어지는 것으로 나타났다. 또한, 실험 전체 구간 중 RTK 고정해, 즉, 반송파 미지정수가 결정된 측위 해가 산출되기까지 약 15분이 소요되었으며, 미지정수 결정율은 약 41% 수준에 그쳤으며, 이는 전통적인 OSR 보정정보 기반의 단일 기준국 RTK 측위 성능 (TTFF 1초, Fix Rate 100%)에 비해 현저히 떨어지는 결과이다.

CLAS 무결성 정보를 이용하여 산출된 RTK 고정해에 대한 보호수준을 분석한 결과, 전반적으로 고정해의 위치오차가 보호수준에 의해 bound 되었지만, 일부 epoch에서는 미지정수가 결정되었음에도 불구하고 수평 및 수직 보호 수준이 각각 약 5 m, 8 m 수준으로 다른 epoch 대비 최대 약 8배가량 큰 것을 확인하였다. 해당 문제의 원인을 다양한 측면에서 분석한 결과, CLAS 보정정보의 적용으로 인해 단일 차분 geometry-free 반송파 측정치 변화율이 임계 값을 초과하여 일부 위성에 대하여 cycle slip flag가 발생하였고, 결과적으로 해당 위성들이 RTK 고정해 산출 과정에서 제외됨으로써 위성 수가 급격히 감소하고, 비정상적으로 큰 보호수준이 산출된 주요 원인임을 확인하였다.

위의 결과를 종합해볼 때, QZSS CLAS 기반 GPS-only RTK 측위 시, 미지정수가 결정된 항법 구간에 대해서는 cm급 위치해가 산출되었지만, 반송파 미지정수 결정 성능 측면에서 볼 때, SSR 보정정보의 측정치 오차 제거 성능이 전통적인 OSR 보정정보에 수준에 미치지 못하는 것을 확인하였다. 결과적으로 SSR 보정정보 기반의 cm급 정밀 측위 시, 보다 안정적인 성능 (측위 정확도, 반송파 미지정수 결정, 사용자 무결성 등)을 확보하기 위해서는 다중 위성군 측정치 활용이 필요할 것으로 보인다. 또한, 각 위성별 SSR 보정정보의 신뢰도를 의미하는 표준편차 뿐만 아니라 위성의 기하학적 배치에 직접적인 영향을 받는 보호 수준은 GPS 단독 측위의 경우, 최대 약 2 m 수준 (GPS TOW 288375s – 288636s 구간 제외)으로 정밀 측위에 적용하기에는 적합하지 않을 것으로 판단된다.