""" 大探标准格式解析代码 """ import gc import numpy as np import struct __all__ = ['get_radar_data'] def get_radar_data_new(f): """ :param filename: :return: dbz, vel, w, zdr, cc, dp, kdp, el_n, azm, el """ f.seek(0, 2) block_length_all = f.tell() f.seek(0) GENERIC_HEADER = { # noqa # 共32bit,保留字节16bit "Magic_Number": struct.unpack('I', f.read(4)), # 魔术字 - --固定标志,用来指示雷达数据文件 # "Major_Version": struct.unpack('h', f.read(2)), # 主版本号 # "Minor_Version": struct.unpack('h', f.read(2)), # 次版本号 "Major_Version": struct.unpack('H', f.read(2)), # 主版本号 "Minor_Version": struct.unpack('H', f.read(2)), # 次版本号 "Generic_Type": struct.unpack('I', f.read(4)), # 文件类型 - --1–基数据文件;2–气象产品文件;3 - --体扫配置文件;4 - --保留;5 - --体扫调度配置文件;6 - --雷达状态数据文件 "Product_Type": struct.unpack('I', f.read(4)), # 产品类型 - --仅在文件类型为2时生效 # 保留字节16bit "Reserved": f.read(16) } # f.seek(32) SITE_CONFIG = { # 共128bit,保留字节46bit "Site_Code": struct.unpack("8s", f.read(8))[0], # 站号 - --站号具有唯一性,用来区别不同的雷达站,如Z9010 "Site_Name": struct.unpack("32s", f.read(32))[0], # 站点名称 - --站点名称,如BeiJing "Latitude": struct.unpack('f', f.read(4))[0], # 纬度 - --雷达站天线所在位置纬度 "Longitude": struct.unpack('f', f.read(4))[0], # 经度 - --雷达站天线所在位置经度 "Antenna_Height": struct.unpack('I', f.read(4))[0], # 天线高度 - --天线馈源水平时海拔高度 "Ground_Height": struct.unpack('I', f.read(4)), # 地面高度 - --雷达塔楼地面海拔高度 "Frequency": struct.unpack('f', f.read(4)), # 工作频率 # TODO gyj modify # "Beam_Width_Hori": struct.unpack('f', f.read(4)), # 水平波束宽度 # "Beam_Width_Vert": struct.unpack('f', f.read(4)), # 垂直波束宽度 "Antenna_Type": struct.unpack('I', f.read(4)), # 天线类型 "TR_Number": struct.unpack('I', f.read(4)), # 收发单元数量 "RDA_Version": struct.unpack('I', f.read(4)), # RDA版本号 - --雷达数据采集软件版本号 # TODO gyj modify # 't':struct.unpack('h', f.read(2))[0], "Radar_Type": {1: "SA", 2: "SB", 3: "SC", 4: "SAD", 7:'SPAR',8:'SPARD', 33: "CA", 34: "CB", 35: "CC", 36: "CCJ", 37: "CD", 38: "CAD", 43:'CPAR',44:'CPARD', 65: "XA", 66: "XAD",69:'XPAR',70:'XPARD'}.get(struct.unpack('h', f.read(2))[0]), # 雷达类型 - --1–SA;2–SB;3–SC;4–SAD;33–CA;34–CB;35–CC;36–CCJ;37–CD;38–CAD;65–XA;66–XAD # 7–SPAR 8–SPARD 43–CPAR 44–CPARD 69–XPAR 70–XPARD # "Antenna_Gain": struct.unpack('h', f.read(2)), # 天线增益 - --编码为实际100倍 # "Transmitting_feeder_loss": struct.unpack('h', f.read(2)), # 发射馈线损耗 - --编码为实际损耗100倍 # "Receiving_feeder_loss": struct.unpack('h', f.read(2)), # 接收馈线损耗 - --编码为实际损耗100倍 # "Other_Loss": struct.unpack('h', f.read(2)), # 其他损耗 - --编码为实际损耗100倍 # 保留字节46bit "Reserved": f.read(54) } # f.seek(160) TASK_CONFIG = { # 共256bit,保留字节40bit "Task_Name": struct.unpack("32s", f.read(32)), # 任务名称 - --任务名称,如VCP21 "Task_Description": struct.unpack("128s", f.read(128)), # 任务描述 "Polarization_Type": struct.unpack('I', f.read(4))[0], # 极化方式 - --1 – 水平极化;2 – 垂直极化;3 – 水平 / 垂直同时;4 – 水平 / 垂直交替 "Scan_Type": struct.unpack('I', f.read(4) ), # 扫描任务类型 - --0 – 体扫;1–单层PPI;2 – 单层RHI;3 – 单层扇扫;4 – 扇体扫;5 – 多层RHI;6 – 手工扫描 # "Pulse_Width": struct.unpack('I', f.read(4) # ), # 脉冲宽度 - --发射脉冲宽度 # "Scan_Start_Time": struct.unpack('I', f.read(4) # ), # 扫描开始时间 - --扫描开始时间为UTC标准时间计数, 1970年1月1日0时为起始计数基准点 #TODO gyj modify 'Scan_Beam_Number':struct.unpack('I', f.read(4)), #扫描配置中的发射波束数量 "Cut_Number": struct.unpack('I', f.read(4) ), # 扫描层数 - --根据扫描任务类型确定的扫描层数,与基数据保存的层数一致 "Ray_Order": struct.unpack('I', f.read(4)), #径向数据排序,0 – 按径向采集时间排序,1 – 按先方位后俯仰方式排序 "Scan_Start_Time": struct.unpack('q', f.read(8) ), # 扫描开始时间 s - --扫描开始时间为UTC标准时间计数, 1970年1月1日0时为起始计数基准点 # "Horizontal_Noise": struct.unpack('f', f.read(4) # ), # 水平通道噪声 - --水平通道的噪声电平 # "Vertical_Noise": struct.unpack('f', f.read(4) # ), # 垂直通道噪声 - --垂直通道的噪声电平 # "Horizontal_Calibration": struct.unpack('f', f.read(4) # ), # 水平通道标定值 - --水平通道的反射率标定常数 # "Vertical_Calibration": struct.unpack('f', f.read(4) # ), # 水平通道噪声温度 - --垂直通道的反射率标定常数 # "Horizontal_Noise_Temperature": struct.unpack('f', f.read(4) # ), # 水平通道噪声温度 # "Vertical_Noise_Temperature": struct.unpack('f', f.read(4) # ), # 垂直通道噪声温度 # "ZDR_Calibration": struct.unpack('f', f.read(4) # ), # ZDR标定偏差 # "PHIDP_Calibration": struct.unpack('f', f.read(4) # ), # 差分相移标定偏差 # "LDR_Calibration": struct.unpack('f', f.read(4) # ), # 系统LDR标定偏差 "Reserved": f.read(68) # 保留字节68bit } # CUT_CONFIG = { # 每个仰角共256bit,保留72bit # } # print(np.frombuffer(f.read(128),'S128')[0].decode('Windows-1252').encode("utf-8")) cut_number = TASK_CONFIG.get("Cut_Number")[0] # 扫描层数 Scan_Beam_Number = TASK_CONFIG.get("Scan_Beam_Number")[0] #扫描配置中的发射波束数量 # print(cut_number) SCAN_CONFIG = {} for i in range(1, Scan_Beam_Number + 1): # f.seek(416) # 扫描波束配置块(一个) SCAN_CONFIG[i] = { 'BEAM_INDEX': struct.unpack('I', f.read(4)), # 按扫描波束依次编号,一个体扫可由多个扫描波束组成,一个扫描波束可输出多个扫描仰角数据 'BEAM_TYPE': struct.unpack('I', f.read(4)), # 扫描波束类型 1 – 宽波束 2 – 窄波束 'SubPulse_Number': struct.unpack('I', f.read(4)), # 子脉冲数量,主工作脉冲和补盲子脉冲数量总和,以下按作用距离由远及近排列1~4 'Tx_Beam_Direction': struct.unpack('f', f.read(4)), # 发射波束中心指向,发射波束中心的俯仰角指向-2.00~90.00 'Tx_Beam_Width-H': struct.unpack('f', f.read(4)), # 发射波束水平宽度 'Tx_Beam_Width-V': struct.unpack('f', f.read(4)), # 发射波束垂直宽度 'Tx_Beam_Gain': struct.unpack('f', f.read(4)), # 发射波束中心增益 'Reserved': f.read(100) # 保留字节100bit } # cut_number = SCAN_CONFIG['0'].get("SubPulse_Number")[0] # cut_number = 4 # print('cut_number: ',cut_number) # 子脉冲参数块(4个) for j in range(1, 5): SCAN_CONFIG[i]['%d' % j] = { # 共128bit,保留70bit 'SubPulse_Strategy': struct.unpack('I', f.read(4)), # 子脉冲策略 0=频分 1=时分 'SubPulse_Modulation': struct.unpack('I', f.read(4)), # 子脉冲调制方式 0=单载频 1=线性正调频ULFM 2=线性负调频DLFM 3=非线性正调频UNLFM 4=非线性负调频DNLFM 'SubPulse_Frequency': struct.unpack('f', f.read(4)), # 子脉冲频率,点频脉冲指工作频率,调频脉冲指中心频率, 1.00~999,000.00 'SubPulse_BandWidth': struct.unpack('f', f.read(4)), # 子脉冲带宽,对调频模式有效,表示双边带宽, MHz兆赫 1.00~99.00 'SubPulse_Width': struct.unpack('I', f.read(4)), # 子脉冲宽度,单位为纳秒ns,信号脉冲宽度(时长) 'Horizontal_Noise': struct.unpack('f', f.read(4)), # 水平通道噪声,dBm分贝毫瓦,水平通道的噪声电平-100.00~0.00 'Vertical_Noise': struct.unpack('f', f.read(4)), # 垂直通道噪声,dBm分贝毫瓦,垂直通道的噪声电平-100.00~0.00 'Horizontal_Calibration': struct.unpack('f', f.read(4)), # 水平通道标定偏差,dB分贝,水平通道的标定偏差0.00~200.00 'Vertical_Calibration': struct.unpack('f', f.read(4)), # 垂直通道标定偏差,dB分贝,垂直通道的标定偏差0.00~200.00 'Horizontal_Noise_Temperature': struct.unpack('f', f.read(4)), # 水平通道噪声温度,单位为K开氏温标,水平通道的噪声温度0.00~800.00 'Vertical_Noise_Temperature': struct.unpack('f', f.read(4)), # 垂直通道噪声温度,单位为K开氏温标,垂直通道的噪声温度0.00~800.00 'ZDR_Calibration': struct.unpack('f', f.read(4)), # ZDR标定偏差,dB分贝,水平通道的标定偏差 -10.00~10.00 'PHIDP_Calibration': struct.unpack('f', f.read(4)), # PHIDP标定偏差,度,水平通道的标定偏差 -180.00~180.00 'LDR_Calibration': struct.unpack('f', f.read(4)), # LDR标定偏差,dB分贝,水平通道的标定偏差 -60~0 'Pulse_Points': struct.unpack('h', f.read(2)), # 脉冲距离库数, 0~10000,该子脉冲探测的距离库数 'Reserved': f.read(70), # 保留字节70bit } # # # ========== 接收仰角配置块 ======== BEAM_CONFIG = {} for i in range(1, cut_number + 1): BEAM_CONFIG['%d' %i] = { 'Cut_Index': struct.unpack('h', f.read(2)), #扫描层索引, 1~256, 基数据仰角编号 'Tx_Beam_Index': struct.unpack('h', f.read(2)), # 发射波束索引, 1~256, 对应的发射波束索引 'Rx_Beam_Elevation': struct.unpack('f', f.read(4)), #接收波束指向,度, -2.00~90.00, PPI模式的俯仰角 'Tx_Beam_Gain': struct.unpack('f', f.read(4)), #发射波束增益,,dB分贝, 在本接收方向上,发射波束的增益 'Rx_Beam_Width_H': struct.unpack('f', f.read(4)), #接收波束水平宽度, 度 'Rx_Beam_Width_V': struct.unpack('f', f.read(4)), #接收波束垂直宽度, 度 'Rx_Beam_Gain': struct.unpack('f', f.read(4)), #接收波束垂直宽度, dB分贝, "Process_Mode": {1: "PPP", 2: "FFT"}.get(struct.unpack('I', f.read(4))[0]), # 处理模式---1---PPP;2---FFT "Wave_Form": {0: "CS连续监测", 1: "CD连续多普勒", 2: "CDX多普勒扩展", 3: "Rx Test", 4: "BATCH批模式", 5: "Dual PRF双PRF", 6: "Staggered PRT 参差PRT"}.get(struct.unpack('I', f.read(4))[0]), # 波形类别---0 – CS连续监测;1 – CD连续多普勒;2 – CDX多普勒扩展;3 – Rx Test;4 – BATC H批模式;5 – Dual PRF双PRF ;6---Staggered PRT 参差PRT "N1PRF1": struct.unpack('f', f.read(4) ), # 第一组脉冲重复频率1---单位赫兹Hz---对于Batch、双PRF和参差PRT模式,表示高PRF值。对于其它单PRF模式,表示唯一的PRF值。 "N1PRF2": struct.unpack('f', f.read(4) ), # 第一组脉冲重复频率2---单位赫兹Hz---对Batch、双PRF和参差PRT模式,表示低PRF值。对其它单PRF模式,无效。 "N2PRF1": struct.unpack('f', f.read(4) ), # 第二组脉冲重复频率1---单位赫兹Hz---对于Batch、双PRF和参差PRT模式,表示高PRF值。对于其它单PRF模式,表示唯一的PRF值。 "N2PRF2": struct.unpack('f', f.read(4) ), # 第二组脉冲重复频率2---单位赫兹Hz---对Batch、双PRF和参差PRT模式,表示低PRF值。对其它单PRF模式,无效。 "Dealiasing_Mode": struct.unpack('I', f.read(4)), # 速度退模糊方法---1 – 单PRF;2 –双PRF3:2模式;3 –双PRF4:3模式;4 –双PRF 5:4模式 "Azimuth": struct.unpack('f', f.read(4)), # 方位角---RHI模式的方位角, 度,0.00~360.00, PPI模式的俯仰角同参数接收波 # "Elevation": struct.unpack('f', f.read(4)), # 俯仰角---PPI模式的俯仰角 "Start_Angle": struct.unpack('f', f.read(4)), # 起始角度----度,-10.00~360.00---PPI扇扫的起始方位角,或RHI模式的高限仰角 "End_Angle": struct.unpack('f', f.read(4)), # 结束角度---度,-10.00~360.00---PPI扇扫的结束方位角,或RHI模式的低限仰角 "Angular_Resolution": struct.unpack('f', f.read(4)), # 角度分辨率--度,---径向数据的角度分辨率,仅用于PPI扫描模式 "Scan_Speed": struct.unpack('f', f.read(4)), # 扫描速度---度/秒---0.00~100.00---PPI扫描的方位转速,或RHI扫描的俯仰转速 "Log_Resolution": struct.unpack('f', f.read(4)), # 强度分辨率---米---1~5,000---强度数据的距离分辨率,支持浮点分辨率 "Doppler_Resolution": struct.unpack('f', f.read(4)), # 多普勒分辨率---米---1~5,000---多普勒数据的距离分辨率 "Maximum_Range1": struct.unpack('I', f.read(4)), # 最大距离1---米---1~500,000---对应脉冲重复频率1的最大可探测距离 "Maximum_Range2": struct.unpack('I', f.read(4)), # 最大距离2---米---1~500,000---对应脉冲重复频率2的最大可探测距离 "Start_Range": struct.unpack('I', f.read(4)), # 起始距离---米---1~500,000---数据探测起始距离 "Sample1": struct.unpack('I', f.read(4)), # 采样个数1---2~512---对应于脉冲重复频率1的采样个数 "Sample2": struct.unpack('I', f.read(4)), # 采样个数2---2~512---对应于脉冲重复频率2的采样个数 "Phase_Mode": struct.unpack('I', f.read(4)), # 相位编码模式---1 – 固定相位;2 – 随机相位;3 – SZ编码 "Atmospheric_Loss": struct.unpack('f', f.read(4)), # 大气衰减---分贝/千米--0.000000~10.000000---双程大气衰减值,精度为小数点后保留6位 "Nyquist_Speed": struct.unpack('f', f.read(4)), # 最大不模糊速度---米/秒---0~100---理论最大不模糊速度 "Moments_Mask": struct.unpack("q", f.read(8)), # 数据类型掩码---0~0xFFFFFFFFFFFFFFFF---以掩码的形式表示当前允许获取的数据类型,其中:0–不允许获取数据;1 –允许获取数据。(json1) "Moments_Size_Mask": struct.unpack("q", f.read(8)),# 数据大小掩码---0~0xFFFFFFFFFFFFFFFF---以掩码形式表示每种数据类型字节数,其中:0–1个字节;1 – 2个字节(json1) "Misc_Filter_Mask": struct.unpack('I', f.read(4)), # 滤波设置掩码--0~0xFFFFFFFF---0–未应用;1–应用(json2) "SQI_Threshold": struct.unpack('f', f.read(4)), # SQI门限---0.00~1.00---- "SIG_Threshold": struct.unpack('f', f.read(4)), # SIG门限--dB分贝---0.00~20.00 "CSR_Threshold": struct.unpack('f', f.read(4)), # CSR门限--dB分贝--0.00~100.00 "LOG_Threshold": struct.unpack('f', f.read(4)), # LOG门限--dB分贝--0.00~20.00 "CPA_Threshold": struct.unpack('f', f.read(4)), # CPA门限--0.00~100.00-- "PMI_Threshold": struct.unpack('f', f.read(4)), # PMI门限--0.00~1.00 "DPLOG_Threshold": struct.unpack('f', f.read(4)), # DPLOG门限--0.00~100.00 "Thresholds_r": struct.unpack('4s', f.read(4)), #阈值门限保留---水平通道的反射率标定常数 =======保留字段====== # 这个时保留字段 "Thresholds_r":data.slice(548,552),#阈值门限保留---水平通道的反射率标定常数 "dBT_Mask": struct.unpack('I', f.read(4)), # dBT质控掩码---0~0xFFFFFFFF---dBT数据使用的质控门限掩码,其中:0–未应用,1–应用(json3) "dBZ_Mask": struct.unpack('I', f.read(4)), # dBZ质控掩码---dBZ数据使用的质控门限掩码, 其中:0–未应用,1–应用(json3) "Velocity_Mask": struct.unpack('I', f.read(4)), # 速度质控掩码----0~0xFFFFFFFF---速度数据使用的质控门限掩码, 其中:0–未应用,1–应用(json3) "Spectrum_Width_Mask": struct.unpack('I', f.read(4)), # 谱宽质控掩码----0~0xFFFFFFFF--谱宽数据使用的质控门限掩码,其中:0–未应用,1–应用(json3) "DP_Mask": struct.unpack('I', f.read(4)), # 偏振量质控掩码----0~0xFFFFFFFF--偏振量数据使用的质控门限掩码,其中:0–未应用,1–应用(json3) "Mask_Reserved": f.read(12), # 保留字节12bit 'Reserved': f.read(4), # 保留字节4bit # 保留12bit # "Scan_Sync": struct.unpack('I', f.read(4)), # 扫描同步标志---用于多部雷达同步扫描标识 "Direction": struct.unpack('I', f.read(4)), # 天线运行方向---仅对PPI模式有效1 – 顺时针;2 – 逆时针 "Ground_Clutter_Classifier_Type": struct.unpack('h', f.read(2)), # 地物杂波图类型---1 – 所有数据不滤波;2 – 全程滤波;3 – 使用实时动态滤波图;4 – 使用静态滤波图 "Ground_Clutter_Filter_Type": struct.unpack('h', f.read(2)), # 地物滤波类型---0 –不滤波;1 – 频域自适应滤波;2 - 固定宽带频域滤波器;3 - 可变宽带频域滤波器;4 - 可变最小方差频域滤波器;5 – IIR时域滤波 "Ground_Clutter_Filter_Notch_Width": struct.unpack('h', f.read(2)), # 地物滤波宽度--0.1 米/秒--0.1~10.0 "Ground_Clutter_Filter_Window": struct.unpack('h', f.read(2)), # 滤波窗口类型---滤波算法FFT窗口类型;0 – 矩形窗;1 – 汉明窗;2 – Blackman窗;3 – 自适应窗口;4 – 无 "Reserved": f.read(44), # 保留字节44bit } #TODO # radar_range = np.arange(BEAM_CONFIG['1']['Sample1'][0]) * BEAM_CONFIG['1']['Log_Resolution'][0] data_block = f.tell() # print(data_block) # return data_shape = {} # for J in range(radial_number): RadialData = { # 径向头,共128字节 保留13 字节 "Radial_State": struct.unpack('I', f.read(4)), # 径向数据状态 - --0–仰角开始;1–中间数据;2–仰角结束;3–体扫开始;4–体扫结束;5–RHI开始;6–RHI结束 "Spot_Blank": struct.unpack('I', f.read(4)), # 消隐标志 - --0–正常;1–消隐 "Sequence_Number": struct.unpack('I', f.read(4)), # 序号 - 1~65536--每个体扫径向从1计数 "Radial_Number": struct.unpack('I', f.read(4)), # 径向数 - 1~1000--每个扫描从1计数 "Elevation_Number": struct.unpack('I', f.read(4)), # 仰角编号 - --仰角编号,每个体扫从1计数 "Azimuth": struct.unpack('f', f.read(4)), # 方位角 ---度 0.00~360.00 --扫描的方位角度 "Elevation": struct.unpack('f', f.read(4)), # 仰角 ---度 -2.00~90.00--扫描的俯仰角度 "Seconds": struct.unpack('q', f.read(8)), # 秒 - --径向数据采集的时间,UTC计数的秒数, 从1970年1月1日0时开始计数 "Microseconds": struct.unpack('I', f.read(4)), # 微秒 - --径向数据采集的时间除去UTC秒数后,留下的微秒数 "Length_of_data": struct.unpack('I', f.read(4)), # 数据长度 --1~100000 --仅本径向数据块所占用的长度如有压缩,长度为压缩后数据长度 "Moment_Number": struct.unpack('I', f.read(4)), # 数据类别数量 ---1~64 --径向数据类别(如Z,V,W等各占一种)的数量 "Scan Beam Index": struct.unpack('h', f.read(2)), # 波束编号 ---1-32 --波束编号 "Horizontal_Estimated_Noise": struct.unpack('h', f.read(2)), # 径向的水平估计噪声 --dB分贝 --编码为实际噪声的 - 100倍 "Vertical_Estimated_Noise": struct.unpack('h', f.read(2)), # 径向的垂直估计噪 - --编码为实际噪声的 - 100 倍 "PRF_FLAG": struct.unpack('I', f.read(4)), #重频标志, 11 – N1 PRF #1 ; 12 – N1 PRF #2; 21 – N2 PRF #1; 22 – N2 PRF #2; 排列方式:11128888:表示使用N1 PRF #1和PRF #2 #######################uint # 88882122:表示使用N2 PRF #1和PRF #2 11122122:表示使用N1和N2的PRF #1和PRF #2 (对应表2-7中09-13的描述)(8表示无效填充数据,弥补用0表示会出现的问题) # "Zip_Type": struct.unpack('h', f.read(2)), # 压缩类型 - --0 - 不压缩;1 - LZO压缩 "Reversed": f.read(70) # 保留70bit } radial_number = RadialData.get("Radial_Number")[0] # 径向数 J type_number = RadialData.get("Moment_Number")[0] # 数据类别数量 K typeHeader = {} # for J in range(radial_number): # for K in range(type_number): for i in range(type_number): # typeHeader[J][K] = { typeHeader = { # 共32字节 "Data_Type": struct.unpack('I', f.read(4)), # 数据类型 ---1~64 --具体径向数据类型见表2 - 6 "Scale": struct.unpack('I', f.read(4)), # 比例 --0~32768 --数据编码的比例 "Offset": struct.unpack('I', f.read(4)), # 偏移 --0~32768 --数据编码的偏移 "Bin_Length": struct.unpack('h', f.read(2)), # 库字节长度 - --保存一个距离库值用的字节数 "Flags": struct.unpack('h', f.read(2)), # 标志 - --数据标志位,暂不使用 "Length": struct.unpack('I', f.read(4)), # 长度 --1~32768 --距离库数据的长度,不包括当前的径向数据头大小 "Reserved": f.read(12) } block_length = typeHeader.get("Length")[0] bin_length = typeHeader.get("Bin_Length")[0] offset = typeHeader.get("Offset")[0] scale = typeHeader.get("Scale")[0] data_body = f.read(block_length) raw = np.frombuffer(data_body, 'u' + str(bin_length)).astype(np.float64) #TODO if i ==0: radar_range = np.arange(len(raw))*30 raw[raw <= 5.] = np.nan # 对于保存的编码值来说,5以下的值表示特殊意义,不应该被解码。 value = (raw - offset) / scale data_type = typeHeader.get("Data_Type")[0] data_shape[data_type] = value.shape[0] block_index = f.tell() roll_num = int((block_length_all - data_block) / (block_index - data_block)) azm_num = int(roll_num / cut_number) gc.collect() # print("=========",data_shape) # return get_radial_data(f, cut_number, azm_num, data_shape, data_block, SITE_CONFIG, radar_range, TASK_CONFIG['Polarization_Type']) return get_radial_data(f, cut_number, azm_num, data_shape, data_block, SITE_CONFIG, radar_range, TASK_CONFIG['Polarization_Type'], BEAM_CONFIG) def get_radial_data(f, cut_number, azm_num, data_shape, data_block, site_config, radar_range, polar_type, BEAM_CONFIG): """ :param f: :param cut_number: cut number :param azm_num: number of azm :param data_shape: length of block :return: dbz, vel, w, zdr, cc, dp, kdp, el_n, azm, el dbz:滤波后反射率(Reflectivity) dBZ vel:径向速度(Doppler Velocity) V w:谱宽(Spectrum Width) W zdr_:差分反射率(Differential Reflectivity ZDR cc_:协相关系数(Cross Correlation Coefficient) CC dp_:差分相移(Differential Phase) DP kdp_差分相移率(Specific Differential Phase) KDp azm:方位角 维度1 el:仰角 维度 1 """ _data_dict = {} for k, v in data_shape.items(): if k == 2: dbz = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['reflectData'] = dbz elif k == 3: vel = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['velData'] = vel elif k == 4: w = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['widData'] = w elif k == 7: zdr = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['zdrData'] = zdr elif k == 9: cc = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['ccData'] = cc elif k == 10: dp = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['dpData'] = dp elif k == 11: kdp = np.empty((cut_number, azm_num, v), dtype=np.float64) _data_dict['kdpData'] = kdp el_n = np.empty((cut_number, azm_num), np.int32) azm = np.empty((cut_number, azm_num), np.float32) el = np.empty((cut_number, azm_num), np.float32) f.seek(data_block) for j in range(cut_number): for k in range(azm_num): RadialData = { # 共64字节 保留13 字节 "Radial_State": struct.unpack('I', f.read(4)), # 径向数据状态 - --0–仰角开始;1–中间数据;2–仰角结束;3–体扫开始;4–体扫结束;5–RHI开始;6–RHI结束 "Spot_Blank": struct.unpack('I', f.read(4) ), # 消隐标志 - --0–正常;1–消隐 "Sequence_Number": struct.unpack('I', f.read(4) ), # 序号 - --每个体扫径向从1计数 "Radial_Number": struct.unpack('I', f.read(4) ), # 径向数 - --每个扫描从1计数 "Elevation_Number": struct.unpack('I', f.read(4) ), # 仰角编号 - --仰角编号,每个体扫从1计数 "Azimuth": struct.unpack('f', f.read(4) ), # 方位角 - --扫描的方位角度 "Elevation": struct.unpack('f', f.read(4) ), # 仰角 - --扫描的俯仰角度 "Seconds": struct.unpack('q', f.read(8) ), # 秒 - --径向数据采集的时间,UTC计数的秒数, 从1970年1月1日0时开始计数 "Microseconds": struct.unpack('I', f.read(4) ), # 微秒 - --径向数据采集的时间除去UTC秒数后,留下的微秒数 "Length_of_data": struct.unpack('I', f.read(4) ), # 数据长度 - --仅本径向数据块所占用的长度如有压缩,长度为压缩后数据长度 "Moment_Number": struct.unpack('I', f.read(4) ), # 数据类别数量 - --径向数据类别(如Z,V,W等各占一种)的数量 "Scan Beam Index": struct.unpack('h', f.read(2) ), # 波束编号 ---1-32 --波束编号 "Horizontal_Estimated_Noise": struct.unpack('h', f.read(2) ), # 径向的水平估计噪声 - --编码为实际噪声的 - 100倍 "Vertical_Estimated_Noise": struct.unpack('h', f.read(2) ), # 径向的垂直估计噪 - --编码为实际噪声的 - 100 倍 "PRF_FLAG": struct.unpack('I', f.read(4)), #重频标志, 11 – N1 PRF #1 ; 12 – N1 PRF #2; 21 – N2 PRF #1; 22 – N2 PRF #2; 排列方式:11128888:表示使用N1 PRF #1和PRF #2 #######################uint # 88882122:表示使用N2 PRF #1和PRF #2 11122122:表示使用N1和N2的PRF #1和PRF #2 (对应表2-7中09-13的描述)(8表示无效填充数据,弥补用0表示会出现的问题) # "Zip_Type": struct.unpack('h', f.read(2) # ), # 压缩类型 - --0 - 不压缩;1 - LZO压缩 "Reversed": f.read(70) # 保留13bit } # print(f.tell()) # f.seek(8672) el_n_ = RadialData.get("Elevation_Number")[0] azm_ = RadialData.get("Azimuth")[0] # el_ = RadialData.get("Elevation")[0] el_ = BEAM_CONFIG['%d' %(j+1)]["Rx_Beam_Elevation"][0] type_number = RadialData.get("Moment_Number")[0] for i in range(type_number): typeHeader = { # 共32字节 "Data_Type": struct.unpack('I', f.read(4)), # 数据类型 - --具体径向数据类型见表2 - 6 "Scale": struct.unpack('I', f.read(4) ), # 比例 - --数据编码的比例 "Offset": struct.unpack('I', f.read(4) ), # 偏移 - --数据编码的偏移 "Bin_Length": struct.unpack('h', f.read(2) ), # 库字节长度 - --保存一个距离库值用的字节数 "Flags": struct.unpack('h', f.read(2) ), # 标志 - --数据标志位,暂不使用 "Length": struct.unpack('I', f.read(4) ), # 长度 - --距离库数据的长度,不包括当前的径向数据头大小 "Reserved": f.read(12) } # print(typeHeader) data_type = typeHeader.get("Data_Type")[0] scale = typeHeader.get("Scale")[0] offset = typeHeader.get("Offset")[0] bin_length = typeHeader.get("Bin_Length")[0] block_length = typeHeader.get("Length")[0] data_body = f.read(block_length) dtype = "u%d" % bin_length # print(body) raw = np.frombuffer(data_body, dtype).astype(np.float64) # raw = np.asarray(struct.unpack(str(block_length) + "B", data_body)).astype(np.float64) raw[raw <= 5.] = np.nan value = (raw - offset) / scale # print(np.nanmin(value)) # print(np.nanmax(value)) if data_type == 2: # dbz[el_n_ - 1, k, :] = value dbz[el_n_ - 1, k, :len(value)] = value elif data_type == 3: # vel[el_n_ - 1, k, :] = value vel[el_n_ - 1, k, :len(value)] = value elif data_type == 4: # w[el_n_ - 1, k, :] = value w[el_n_ - 1, k, :len(value)] = value elif data_type == 7: # zdr[el_n_ - 1, k, :] = value zdr[el_n_ - 1, k, :len(value)] = value elif data_type == 9: # cc[el_n_ - 1, k, :] = value cc[el_n_ - 1, k, :len(value)] = value elif data_type == 10: # dp[el_n_ - 1, k, :] = value dp[el_n_ - 1, k, :len(value)] = value elif data_type == 11: # kdp[el_n_ - 1, k, :] = value kdp[el_n_ - 1, k, :len(value)] = value el_n[el_n_ - 1, k] = el_n_ azm[el_n_ - 1, k] = azm_ el[el_n_ - 1, k] = el_ # print(f"the No.{el_n_} elevation is {el_}") gc.collect() result = { 'radarNumber': site_config['Site_Code'], 'longitude': site_config['Longitude'], 'latitude': site_config['Latitude'], 'altitude': site_config['Antenna_Height'], 'radarRange': radar_range, 'azimuthData': azm[1, :], 'elevationData': el[:, 1], 'polar_type': polar_type } result.update(_data_dict) f.close() return result