RTK reference station coordinates determination using shell script

category: gnss
tags: rtk

I think that the reference stations of the Geographical Survey Institute of Japan (GSI) are often used to determine the coordinates of the RTK reference station in Japan. However, even though the reference point to be used is fixed, it is a hassle to search for and download the reference point every time. However, GSI also provides a way to download these data via ftp (file transfer protocol).

In addition, if we collect observation data oneself and process it interactively with RTKPOST application of RTKLIB, we will deepen our understanding of RTK. However, the work of calculating each time the antenna is changed, extracting the coordinates fixed with, say Microsoft Excel, and then calculating the average is a hassle. However, RTKLIB also ships with an application that can be run non-interactively called rnx2rtkp.

Therefore, I created a shell script that performs a series of processes from downloading the observation data of GSI to determining the coordinates of the RTK reference station.

I assume:

  1. You have obtained an ftp account (ID, password) of GSI.
  2. You know the reference station number, latitude, longitude, and ellipsoidal height of the reference station.
  3. curl, convbin, rnx2rtkp, awk are in a directory in your path.

Since the satellite constellation has a 12-hour cycle or a 24-hour cycle, it is necessary to observe for 24 hours or more to determine the coordinates of the RTK reference station. On the other hand, the RINEX format is often used for archiving the reference station, but the RINEX format cannot combine multiple files as they are. In order to combine multiple observation files, u-blox format files, which are raw data, are used here. The raw data is in stream format, so you can combine the files as is.

Please rewrite as appropriate according to your environment for YEAR, MONTH, DAY, the reference station number POS, coordinates LAT, LON, and HEIGHT , GSI account and password GSIID, GSIPW, and the URI of calculating the coordinate of the reference station RTKURI.


YEAR=2021      # date
POS=1157        # GSI reference station number
LAT=34.39704732 # position
GSIID= GSI ID # change me
GSIPW= GSI password # change me

T1=`date +%s -d "$YEAR$MONTH$DAY"`
T2=`date +%s -d "${YEAR}0101"`
DAYNO=`echo $T1 $T2 | awk '{printf "%03d", ($1-$2)/60/60/24+1}'`
echo Getting GSI observation data on $YEAR/$MONTH/$DAY at $POS
curl -Ss -u $GSIID:$GSIPW -O $GSIURI/$YEAR/$DAYNO/$POS${DAYNO}0.${YEAR:2:2}o.gz

echo -n "Getting RTK observation data: "   # an example of obtaining GNSS raw data
for hourcode in {a..x}; do
    echo -n "$YEAR$MONTH$DAY$hourcode "
    curl -Ss -O $RTKURI/$YEAR$MONTH/$YEAR$MONTH$DAY$hourcode.ubx
cat $YEAR$MONTH${DAY}?.ubx > $YEAR$MONTH$DAY.ubx && rm $YEAR$MONTH${DAY}?.ubx

echo Converting RTK raw to RINEX
convbin $YEAR$MONTH$DAY.ubx -v 3.02 -ti 1 && rm $YEAR$MONTH$DAY.ubx

echo Executing RTK
rnx2rtkp -p 2 -h -t -u \
    $YEAR$MONTH$DAY.obs $YEAR$MONTH$DAY.nav $POS${DAYNO}0.${YEAR:2:2}o.gz \
    -l $LAT $LON $HEIGHT \
    > $YEAR$MONTH$DAY.pos

echo Averaged position:
awk '$6==1{n+=1; lat+=$3; lon+=$4; height+=$5} END {printf "%.7f %.7f %.3f\n", lat/n, lon/n, height/n}' $YEAR$MONTH$DAY.pos


On the day when the observation data of the RTK reference station was obtained, the 7th digit after the decimal point of the longitude did not change even if it was processed on another day, so I think that the coordinates can be determined on the order of centimeters.