4. 消息1001 解析

任务说明

本案例实现RTCM 3.2协议中消息1001（L1-only GPS RTK观测）的二进制报文解析。输入为符合RTCM 3.2标准的1001消息字节流，输出为结构化字段集合，包含消息头、卫星数量、以及每颗卫星的L1伪距、载波相位、多普勒和信噪比观测值。解析结果可直接用于后续RTK定位解算。

消息1001格式定义

消息1001承载GPS L1频段的单频观测数据，不包含载波相位平滑伪距。其报文结构如下：

代码实现

import struct

class RTCM1001Parser:
    """RTCM 3.2 消息1001解析器（L1-only GPS RTK观测）"""
    
    def __init__(self):
        self.msg_type = None
        self.station_id = None
        self.gps_epoch = None
        self.sync_flag = None
        self.sat_count = None
        self.smoothing = None
        self.satellites = []  # 每颗卫星的观测数据
    
    def _get_bits(self, data, start, length):
        """从字节数组提取位字段（大端序）"""
        if length == 0:
            return 0
        byte_start = start // 8
        bit_offset = start % 8
        # 读取足够字节
        needed = (bit_offset + length + 7) // 8
        chunk = data[byte_start:byte_start + needed]
        val = int.from_bytes(chunk, 'big')
        # 右移并掩码
        val >>= (8 * needed - bit_offset - length)
        return val & ((1 << length) - 1)
    
    def _sign_extend(self, value, bits):
        """有符号数扩展"""
        if value & (1 << (bits - 1)):
            return value - (1 << bits)
        return value
    
    def parse(self, payload: bytes) -> 'RTCM1001Parser':
        """
        解析1001消息负载（不含3字节前导和长度字段）
        payload: 消息体字节流（从消息号开始）
        """
        if len(payload) < 8:
            raise ValueError(f"负载过短: {len(payload)} 字节")
        
        # 消息头（前64位）
        self.msg_type = self._get_bits(payload, 0, 12)
        if self.msg_type != 1001:
            raise ValueError(f"消息类型错误: {self.msg_type}")
        
        self.station_id = self._get_bits(payload, 12, 12)
        self.gps_epoch = self._get_bits(payload, 24, 30)
        self.sync_flag = self._get_bits(payload, 54, 1)
        self.sat_count = self._get_bits(payload, 55, 5)
        self.smoothing = self._get_bits(payload, 60, 1)
        
        # 逐卫星解析（每颗卫星72位）
        self.satellites = []
        pos = 61  # 当前位偏移
        
        for _ in range(self.sat_count):
            if pos + 72 > len(payload) * 8:
                raise ValueError(f"数据不足: 需要{pos+72}位, 仅有{len(payload)*8}位")
            
            sat = {}
            sat['prn'] = self._get_bits(payload, pos, 6)
            pos += 6
            
            # L1伪距（24位无符号，单位0.02m，偏移299792.458m）
            raw_pseudo = self._get_bits(payload, pos, 24)
            sat['pseudo_range'] = raw_pseudo * 0.02 + 299792.458
            pos += 24
            
            # L1载波相位（20位有符号，单位0.0005/2π周）
            raw_phase = self._sign_extend(self._get_bits(payload, pos, 20), 20)
            sat['carrier_phase'] = raw_phase * 0.0005 / (2 * 3.141592653589793)
            pos += 20
            
            # 锁定时间（7位）
            sat['lock_time'] = self._get_bits(payload, pos, 7)
            pos += 7
            
            # 多普勒（8位有符号，单位0.0001Hz）
            raw_doppler = self._sign_extend(self._get_bits(payload, pos, 8), 8)
            sat['doppler'] = raw_doppler * 0.0001
            pos += 8
            
            # CNR（6位，单位0.25dB-Hz）
            sat['cnr'] = self._get_bits(payload, pos, 6) * 0.25
            pos += 6
            
            # 半周期标志（1位）
            sat['half_cycle'] = self._get_bits(payload, pos, 1)
            pos += 1
            
            self.satellites.append(sat)
        
        return self
    
    def summary(self) -> str:
        """生成解析摘要"""
        lines = [
            f"消息类型: {self.msg_type}",
            f"参考站ID: {self.station_id}",
            f"GPS历元: {self.gps_epoch}",
            f"同步标志: {self.sync_flag}",
            f"卫星数量: {self.sat_count}",
            f"平滑标志: {self.smoothing}",
            "--- 卫星观测 ---"
        ]
        for i, sat in enumerate(self.satellites):
            lines.append(
                f"卫星{i+1}: PRN={sat['prn']:2d} | "
                f"伪距={sat['pseudo_range']:.3f} m | "
                f"相位={sat['carrier_phase']:.6f} 周 | "
                f"多普勒={sat['doppler']:.4f} Hz | "
                f"CNR={sat['cnr']:.2f} dB-Hz | "
                f"锁定={sat['lock_time']} | "
                f"半周期={sat['half_cycle']}"
            )
        return "\n".join(lines)

# ===== 使用示例 =====
if __name__ == "__main__":
    # 构造一条包含2颗卫星的1001消息负载（十六进制）
    # 消息头: 类型1001(0x3E9), 站ID=1, 历元=100, 同步=0, 卫星数=2, 平滑=0
    # 卫星1: PRN=5, 伪距原始值=10000, 相位原始值=2000, 锁定=10, 多普勒=-5, CNR=40, 半周期=0
    # 卫星2: PRN=12, 伪距原始值=10050, 相位原始值=-1500, 锁定=3, 多普勒=8, CNR=35, 半周期=1
    
    # 手动构造位流（简化：使用struct打包）
    # 消息头64位 = 8字节
    header = 0
    header |= 1001 << 52   # 12位消息号
    header |= 1 << 40       # 12位站ID
    header |= 100 << 10     # 30位历元
    header |= 0 << 9        # 1位同步
    header |= 2 << 4        # 5位卫星数
    header |= 0 << 3        # 1位平滑（剩余3位填充0）
    
    payload = bytearray(header.to_bytes(8, 'big'))
    
    # 卫星1数据（72位 = 9字节）
    sat1 = 0
    sat1 |= 5 << 66         # PRN 6位
    sat1 |= 10000 << 42     # 伪距 24位
    sat1 |= (2000 & 0xFFFFF) << 22  # 相位 20位（正数）
    sat1 |= 10 << 15        # 锁定 7位
    sat1 |= ((-5) & 0xFF) << 7  # 多普勒 8位（有符号）
    sat1 |= 40 << 1         # CNR 6位
    sat1 |= 0                # 半周期 1位
    
    payload.extend(sat1.to_bytes(9, 'big'))
    
    # 卫星2数据
    sat2 = 0
    sat2 |= 12 << 66        # PRN
    sat2 |= 10050 << 42     # 伪距
    sat2 |= ((-1500) & 0xFFFFF) << 22  # 相位（负数补码）
    sat2 |= 3 << 15         # 锁定
    sat2 |= (8 & 0xFF) << 7 # 多普勒
    sat2 |= 35 << 1         # CNR
    sat2 |= 1                # 半周期
    
    payload.extend(sat2.to_bytes(9, 'big'))
    
    # 解析
    parser = RTCM1001Parser()
    parser.parse(bytes(payload))
    print(parser.summary())

验证方法

运行上述代码，观察输出。预期输出应包含：消息类型1001、参考站ID=1、GPS历元=100、卫星数量=2。卫星1的伪距应为 10000×0.02 + 299792.458 = 299992.458 m，相位为 2000×0.0005/(2π) ≈ 0.159155 周，多普勒为 -5×0.0001 = -0.0005 Hz，CNR为 40×0.25 = 10.00 dB-Hz。卫星2的伪距应为 10050×0.02 + 299792.458 = 299993.458 m，相位为 -1500×0.0005/(2π) ≈ -0.119366 周，多普勒为 0.0008 Hz，CNR为 8.75 dB-Hz，半周期标志为1。若输出与上述计算一致，则解析逻辑正确。