Skip to Content
| GAGE


Trimble NetR9 Receiver - GSOF Messages

General Serial Output Format (GSOF) messages are a Trimble proprietary format and can be used to send information such as position and status to a third-party device.

More detailed information regarding the structure of this binary message format can be found at Trimble’s GSOF Specification Page.

To enable a GSOF stream on a Trimble NetR9 receiver:

  1. Go to ’I/O Configuration’ on the web interface.
  2. Either select an unused serial port or ’add TCP/IP or UDP Port.
  3. Select the ’GSOF’ message type from the pull down menu.
  4. If using ’add TCP/IP or UDP Port’, specify the port number that you’d like the stream to broadcast from.
  5. If using a serial port, configure the baud, parity and flow to match your third-pary device.
  6. Select which messages you would like the receiver to output.
  7. Click ’OK’ to enable the stream.

See below for an example python scripts ( getGsofStream.py  and utilities.py) connect to a TCP socket and parses GSOF message.

Sample output:

Reference not specified, using first epoch: X =-1282458.7777 Y=-4718377.8575 Z=4084176.9870
WN - GPS WEEK
SOW - GPS time
FLAG1 FLAG2
SATS - # SVS
X[m] Y[m] Z[m]
eX[m] eY[m] eZ[m] - Position VCV Position RMS
E[m] N[m] U[m]
# WN  SOW     FLAG1 FLAG2 SATS  X[m]          Y[m]            Z[m]       eX[m]   eY[m]     eZ[m]   E[m]   N[m]   U[m]
2035 319837.000 191  35 16  -1282458.7759  -4718377.8519   4084176.9822  0.0144  0.0232  0.0174  0.0002  0.0001 -0.0075
2035 319838.000 191  35 16  -1282458.7719  -4718377.8502   4084176.9807  0.0146  0.0234  0.0176  0.0037  0.0007 -0.0106
2035 319839.000 191  35 16  -1282458.7801  -4718377.8574   4084176.9830  0.0148  0.0237  0.0178 -0.0024 -0.0034 -0.0021

 

 

 

############### getGsofStream.py  ###################

from struct import unpack
from utilities import *
import socket
import math
import argparse

__author__ = "Henry T. Berglund"


class Gsof(object):
    """ A class to parse GSOF messages from a TCP stream. """

    def __init__(self, sock=None):
        if sock is None:
            self.sock = socket.socket(
                socket.AF_INET, socket.SOCK_STREAM)
        else:
            self.sock = sock
        self.msg_dict = {}
        self.msg_bytes = None
        self.checksum = None
        self.rec_dict = {}

    def connect(self, host, port):
        self.sock.connect((host, port))

    def get_message_header(self):
        data = self.sock.recv(7)
        msg_field_names = (’STX’, ’STATUS’, ’TYPE’, ’LENGTH’,
                           ’T_NUM’, ’PAGE_INDEX’, ’MAX_PAGE_INDEX’)
        self.msg_dict = dict(zip(msg_field_names, unpack(’>7B’, data)))
        self.msg_bytes = self.sock.recv(self.msg_dict[’LENGTH’] - 3)
        (checksum, etx) = unpack(’>2B’, self.sock.recv(2))

        def checksum256(st):
            """Calculate checksum"""
            return reduce(lambda x, y: x+y, map(ord, st)) % 256
        if checksum-checksum256(self.msg_bytes+data[1:]) == 0:
            self.checksum = True
        else:
            self.checksum = False

    def get_records(self):
        while len(self.msg_bytes) > 0:
            # READ THE FIRST TWO BYTES FROM RECORD HEADER
            record_type, record_length = unpack(’>2B’, self.msg_bytes[0:2])
            self.msg_bytes = self.msg_bytes[2:]
            self.select_record(record_type, record_length)

    def select_record(self, record_type, record_length):
        if record_type == 1:
            rec_field_names = (’GPS_TIME’, ’GPS_WEEK’, ’SVN_NUM’,
                               ’FLAG_1’, ’FLAG_2’, ’INIT_NUM’)
            rec_values = unpack(’>LH4B’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 2:
            rec_field_names = (’LATITUDE’, ’LONGITUDE’, ’HEIGHT’)
            rec_values = unpack(’>3d’, self.msg_bytes[0:record_length])
            rec_values = list(rec_values)
            rec_values = (math.degrees(rec_values[0]), math.degrees(rec_values[1]), rec_values[2])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 3:
            rec_field_names = (’X_POS’, ’Y_POS’, ’Z_POS’)
            rec_values = unpack(’>3d’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 4:
            rec_field_names = (’LOCAL_DATUM_ID’, ’LOCAL_DATUM_LAT’,
                               ’LOCAL_DATUM_LON’, ’LOCAL_DATUM_HEIGHT’, ’OPRT’)
            rec_values = unpack(’>8s3dB’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 5:
            rec_field_names = (’LOCAL_DATUM_ID’, ’LOCAL_ZONE_ID’,
                               ’LOCAL_ZONE_NORTH’, ’LOCAL_ZONE_EAST’, ’LOCAL_DATUM_HEIGHT’)
            rec_values = unpack(’>2s3d’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 6:
            rec_field_names = (’DELTA_X’, ’DELTA_Y’, ’DELTA_Z’)
            rec_values = unpack(’>3d’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 7:
            rec_field_names = (’DELTA_EAST’, ’DELTA_NORTH’, ’DELTA_UP’)
            rec_values = unpack(’>3d’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 8:
            rec_field_names = (’VEL_FLAG’, ’VELOCITY’, ’HEADING’, ’VERT_VELOCITY’)
            rec_values = unpack(’>B3f’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 9:
            rec_field_names = (’PDOP’, ’HDOP’, ’VDOP’, ’TDOP’)
            rec_values = unpack(’>4f’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 10:
            rec_field_names = (’CLOCK_FLAG’, ’CLOCK_OFFSET’, ’FREQ_OFFSET’)
            rec_values = unpack(’>B2d’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 11:
            rec_field_names = (’POSITION_RMS_VCV’, ’VCV_XX’, ’VCV_XY’, ’VCV_XZ’,
                               ’VCV_YY’, ’VCV_YZ’, ’VCV_ZZ’, ’UNIT_VAR_VCV’, ’NUM_EPOCHS_VCV’)
            rec_values = unpack(’>8fh’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 12:
            rec_field_names = (’POSITION_RMS_SIG’, ’SIG_EAST’, ’SIG_NORT’, ’COVAR_EN’, ’SIG_UP’,
                               ’SEMI_MAJOR’, ’SEMI_MINOR’, ’ORIENTATION’, ’UNIT_VAR_SIG’,
                               ’NUM_EPOCHS_SIG’)
            rec_values = unpack(’>9fh’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 15:
            rec_field_names = ’SERIAL_NUM’
            rec_values = unpack(’>l’, self.msg_bytes[0:record_length])
            self.rec_dict.update({rec_field_names: rec_values[0]})
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 16:
            rec_field_names = (’GPS_MS_OF_WEEK’, ’CT_GPS_WEEK’, ’UTC_OFFSET’, ’CT_FLAGS’)
            rec_values = unpack(’>l2hB’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 26:
            rec_field_names = (’UTC_MS_OF_WEEK’, ’UTC_GPS_WEEK’, ’UTC_SVS_NUM’, ’UTC_FLAG_1’, ’UTC_FLAG_2’,
                               ’UTC_INIT_NUM’)
            rec_values = unpack(’>lh4B’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 34:
            NUM_OF_SVS = unpack(’>B’, self.msg_bytes[0])
            self.msg_bytes = self.msg_bytes[1:]
            rec_field_names = (’PRN’, ’SV_SYSTEM’, ’SV_FLAG1’, ’SV_FLAG2’, ’ELEVATION’, ’AZIMUTH’,
                               ’SNR_L1’, ’SNR_L2’, ’SNR_L5’)
            for field in xrange(len(rec_field_names)):
                self.rec_dict[rec_field_names[field]] = []
            for sat in xrange(NUM_OF_SVS[0]):
                rec_values = unpack(’>5Bh3B’, self.msg_bytes[0:10])
                self.msg_bytes = self.msg_bytes[10:]
                for num in xrange(len(rec_field_names)):
                    self.rec_dict[rec_field_names[num]].append(rec_values[num])
        elif record_type == 37:
            rec_field_names = (’BATT_CAPACITY’, ’REMAINING_MEM’)
            rec_values = unpack(’>hd’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        elif record_type == 41:
            rec_field_names = (’GPS_TIME_REF’, ’GPS_WEEK_REF’, ’LAT_REF’, ’LONG_REF’,
                               ’HEIGHT_REF’, "QI")
            rec_values = unpack(’>LH3dB’, self.msg_bytes[0:record_length])
            self.rec_dict.update(dict(zip(rec_field_names, rec_values)))
            self.msg_bytes = self.msg_bytes[record_length:]
        else:
            """Unknown record type? Skip it for now!"""
            #print record_type
            self.msg_bytes = self.msg_bytes[record_length:]


def main():
    parser = argparse.ArgumentParser(prog=’getGsofStream’, usage=’%(prog)s -i IP_ADDRESS -p PORT -r X Y Z ’)
    parser.add_argument(’-i’, ’--IP’, action="store", default="69.44.87.212", type=str)
    parser.add_argument(’-p’, ’--PORT’, action="store", default=9999, type=int)
    parser.add_argument(’-r’, ’--ref’, nargs=3, action=’append’)
    results = parser.parse_args()
    # OPEN GSOF STREAM
    s = Gsof()
    s.connect(results.IP, results.PORT)
    try:
        # Use user provided reference coordinates if given
        ref_x, ref_y, ref_z = results.ref[0]
        print "# Reference specified, using X =%.3f Y=%.3f Z=%.3f" % (ref_x, ref_y, ref_z)
    except:
        # Get reference coordinates from first epoch if not provided by user
        s.get_message_header()
        s.get_records()
        ref_x, ref_y, ref_z = s.rec_dict[’X_POS’], s.rec_dict[’Y_POS’], s.rec_dict[’Z_POS’]
        print "# Reference not specified, using first epoch: X =%.4f Y=%.4f Z=%.4f" % (ref_x, ref_y, ref_z)
    print ""
    print "WN - GPS WEEK"
    print "SOW - GPS time"
    print "FLAG1 FLAG2"
    print "SATS - # SVS"
    print "X[m] Y[m] Z[m]"
    print "eX[m] eY[m] eZ[m] - Position VCV Position RMS"
    print "E[m] N[m] U[m]"
    print ""
    print "# WN  SOW     FLAG1 FLAG2 SATS  X[m]          Y[m]            Z[m]       eX[m]   eY[m]     eZ[m]   E[m]   N[m]   U[m]"
    while 1:
        # READ GSOF STREAM
        s.get_message_header()
        s.get_records()
        # PRINT GSOF STREAM
        x = s.rec_dict[’X_POS’]
        y = s.rec_dict[’Y_POS’]
        z = s.rec_dict[’Z_POS’]
        neu = xyz2neu(ref_x, ref_y, ref_z, x, y, z)
        output = "%04d %.3f %3d %3d %2d %14.4f %14.4f %14.4f %7.4f %7.4f %7.4f %7.4f %7.4f %7.4f" % (
            s.rec_dict[’GPS_WEEK’],
            s.rec_dict[’GPS_TIME’]/1000.0,
            s.rec_dict[’FLAG_1’],
            s.rec_dict[’FLAG_2’],
            s.rec_dict[’SVN_NUM’],
            s.rec_dict[’X_POS’],
            s.rec_dict[’Y_POS’],
            s.rec_dict[’Z_POS’],
            math.sqrt(s.rec_dict[’VCV_XX’]),
            math.sqrt(s.rec_dict[’VCV_YY’]),
            math.sqrt(s.rec_dict[’VCV_ZZ’]),
            neu[1],
            neu[0],
            neu[2],
        )
        print(output)


if __name__ == ’__main__’:
    try:
        main()
    except KeyboardInterrupt:
        print "\n#: Stopping, you hit ctl+C. "
        #traceback.print_exc()

 

###### utilities.py #################

__author__ = ’berglund’

"""
Utilities for various geodetic conversions
"""
from math import sin, cos, sqrt, atan2, degrees, radians

def llh2xyz(lat, lon, h, a=6378137.0, finv=298.257223563, input_degrees=True):
    """
    Converts ellipsoidal coordinates to cartesian.

    :param lat:
    :param lon:
    :param h:

    :return: A tuple containing cartesian coordinates in ECEF [X, Y, Z] meters.
    """
    # convert degrees to radians if necessary
    if input_degrees:
        lat = radians(lat)
        lon = radians(lon)
    # Default (a=6378137.0, finv=298.257223563) is for WGS84
    f = 1/finv
    b = a*(1 - f)
    e2 = 1 - (1 - f)**2
    # Compute the Cartesian coordinates
    v = a/sqrt(1-e2*sin(lat)*sin(lat))
    x = (v+h)*cos(lat)*cos(lon)
    y = (v+h)*cos(lat)*sin(lon)
    z = (v*(1-e2)+h)*sin(lat)
    return x, y, z

def xyz2llh(x, y, z, a=6378137.0, finv=298.257223563, output_degrees=True):
    """
    Converts cartesian coordinates to ellipsoidal. Uses an iterative algorithm.

    :param x: x coordinate in meters
    :param y: y coordinate in meters
    :param z: z coordinate in meters
    :param a: major semi-axis (Default: WGS84)
    :param finv: inverse of flattening (Default: WGS84)
    :param output_degrees: A boolean to select output in degrees (default is True)
    :return: A tuple containing geodetic coordinates in WGS84 (lat, lon, h).

    WGS84: a=6378137.0, finv=298.257223563
    """
    # Default (a=6378137.0, finv=298.257223563) is for WGS84
    f = 1/finv
    b = a*(1 - f)
    e2 = 1 - (1 - f)**2
    # Latitude and height convergence criteria
    elat = 1.0e-12
    eht = 1.e-5
    # Initial values for iteration
    p = sqrt(x*x + y*y)
    lat = atan2(z, p*(1-e2))
    h = 0
    dh = 1
    dlat = 1
    # Iterate until lat & h converge to elat & eht
    while dlat > elat or dh > eht:
        lat0 = lat
        h0 = h
        v = a/sqrt(1-e2*sin(lat)*sin(lat))
        h = p/cos(lat)-v
        lat = atan2(z, p*(1-e2*v/(v+h)))
        dlat = abs(lat-lat0)
        dh = abs(h-h0)
    lon = atan2(y, x)
    if output_degrees:
        # output degrees
        geodetic_coordinates = (degrees(lat), degrees(lon), h)
    else:
        # output radians
        geodetic_coordinates = (lat, lon, h)
    return geodetic_coordinates


def xyz2neu(X, Y, Z, x, y, z):
    """
    Rotates the <(x-X),(y-Y),(z-Z)> vector into a topocentric coordinate system
    tangential to the reference point XYZ.

    Assumes coordinates are given in meters.
    :param X: Reference Coordinate X (meters)
    :param Y: Reference Coordinate Y (meters)
    :param Z: Reference Coordinate Z (meters)
    :param x: Observed location x (meters)
    :param y: Observed location y (meters)
    :param z: Observed location z (meters)
    :return: A tuple containing topocentric coordinates in n, e, u (meters)

    WGS84: a=6378137.0, finv=298.257223563
    """
    lat, lon, h = xyz2llh(X, Y, Z, a=6378137.0, finv=298.257223563, output_degrees=False)
    e = -sin(lon) * (x - X) + cos(lon) * (y - Y)
    n = -sin(lat) * cos(lon) * (x - X) - sin(lat) * sin(lon) * (y - Y) + cos(lat) * (z - Z)
    u = cos(lat) * cos(lon) * (x - X) + cos(lat) * sin(lon) * (y - Y) + sin(lat) * (z - Z)
    return n, e, u


def xyz_to_globk_neu(X, Y, Z, x, y, z):
    """
    Finds the arc-length between two LLH points in dN and dE.
    The parameter dU is defined as the difference between the observed and reference ellipsoidal height.

    Assumes coordinates are given in meters.
    :param X: Reference Coordinate X (meters)
    :param Y: Reference Coordinate Y (meters)
    :param Z: Reference Coordinate Z (meters)
    :param x: Observed location x (meters)
    :param y: Observed location y (meters)
    :param z: Observed location z (meters)
    :return: A tuple containing n, e, u (meters)
    """

    def sub(a, b):
        """
        Subtracts elements in two lists
        """
        c = [i - j for i, j in zip(a, b)]
        return c

    def mul(a, b):
        """
        Multiplies elements in two lists
        """
        c = [i * j for i, j in zip(a, b)]
        return c

    P1 = xyz2llh(X, Y, Z, a=6378137.0, finv=298.257223563, output_degrees=False)
    P2 = xyz2llh(x, y, z, a=6378137.0, finv=298.257223563, output_degrees=False)
    # Subtract the reference coordinate from the observed coordinate
    P_diff = sub(P2, P1)
    # Convert the dLLH to dNEU
    # Assumes a sphere
    rot = [6378137.0, 6378137.0*cos(P1[0]), 1]
    NEU = mul(P_diff, rot)

    return tuple(NEU)


def covrot(X, Y, Z, cxx, cyy, czz, cxy, cxz, cyz):
    """
    Rotates the sigmas and the correlations from xyz to enu
    :param cxx:
    :param cyy:
    :param czz:
    :param cxy:
    :param cxz:
    :param cyz:
    :return:
    """
    def matmult(a, b):
        zip_b = zip(*b)
        return [[sum(ele_a*ele_b for ele_a, ele_b in zip(row_a, col_b)) for col_b in zip_b] for row_a in a]

    lat, lon, h = xyz2llh(X, Y, Z, a=6378137.0, finv=298.257223563, output_degrees=False)

    rxx = -sin(lat)*cos(lon)
    rxy = -sin(lat)*sin(lon)
    rxz = cos(lat)
    ryx = -sin(lon)
    ryy = cos(lon)
    ryz = 0
    rzx = cos(lat)*cos(lon)
    rzy = cos(lat)*sin(lon)
    rzz = sin(lat)

    rotation_matrix = [[rxx, rxy, rxz],
                       [ryx, ryy, ryz],
                       [rzx, rzy, rzz]]
    rotation_transpose = [[rxx, ryx, rzx],
                          [rxy, ryy, rzy],
                          [rxz, ryz, rzz]]
    covariance_matrix = [[cxx**2, cxy*cxx*cyy, cxz*cxx*czz],
                         [cxy*cxx*cyy, cyy**2, cyz*cyy*czz],
                         [cxz*cxx*czz, cyz*cyy*czz, czz**2]]

    tmp = matmult(rotation_matrix, covariance_matrix)

    result = matmult(tmp, rotation_transpose)

    cxx = result[0][0]**0.5
    cyy = result[1][1]**0.5
    czz = result[2][2]**0.5
    cxy = result[0][1]/(cxx*cyy)
    cxz = result[0][2]/(cxx*czz)
    cyz = result[1][2]/(cyy*czz)

    return cxx, cyy, czz, cxy, cxz, cyz


Code disclaimer information

This document contains programming examples. UNAVCO grants you a nonexclusive copyright license to use all programming code examples from which you can generate similar function tailored to your own specific needs. All sample code is provided by UNAVCO for illustrative purposes only. These examples have not been thoroughly tested under all conditions. All programs contained herein are provided to you "AS IS" without any warranties of any kind.

 

 

Last modified: 2019-12-27  16:36:35  America/Denver