网友给发了一个天气雷达拼图数据文件，问MeteoInfo是否能打开并绘图。因为是二进制文件，没有数据文件格式说明文档是无法解析的。没找到正式的说明文档，不过在网上找到一个用Python解析数据文件的帖子（https://blog.csdn.net/weixin_45863084/article/details/140492948），试了一下应该是同一种数据文件。参考该贴写了一个MeteoInfoLab脚本来解析和绘制雷达拼图数据文件，MeteoInfoLab的脚本语言是Jython，和Python的语法一致，因此移植很方便。
雷达拼图二进制文件的格式比较简单，有一个256字节的文件头，后面就是等经纬度格点数据了。文件头里包含了数据文件产品描述信息、时间、格点分布等。用一个类 BinaryReader 来读取雷达拼图数据：
import structimport bz2class BinaryReader:"""Radar data reader"""def __init__(self, buffer, little_endian=True):self.buffer = bufferself.index = 0self.little_endian = little_endianself.read_header()def seek(self, offset, origin=0):if origin == 0: # Set absolute positionself.index = offsetelif origin == 1: # Set relative to current positionself.index += offsetelif origin == 2: # Set relative to end of fileself.index = len(self.buffer) + offsetdef read_bytes(self, length):data = self.buffer[self.index:self.index + length]self.index += lengthreturn datadef read_string(self, length):fmt = '{}s'.format(length)value, = struct.unpack_from(fmt, self.buffer, self.index)value.replace('\x00', '').replace(' ', '')self.index += struct.calcsize(fmt)return valuedef read_int(self):fmt = '<i' if self.little_endian else '>i'value, = struct.unpack_from(fmt, self.buffer, self.index)self.index += struct.calcsize(fmt)return valuedef read_int16(self):fmt = '<h' if self.little_endian else '>h'value, = struct.unpack_from(fmt, self.buffer, self.index)self.index += struct.calcsize(fmt)return valuedef read_header(self):"""Read data header"""self.label = self.read_string(4) # Data file IDself.version = self.read_string(4) # Data file versionself.file_bytes = self.read_int()self.mosaic_id = self.read_int16()self.coordinate = self.read_int16()self.varname = self.read_string(8)self.description = self.read_string(64)BlockPos = self.read_int()BlockLen = self.read_int()self.TimeZone = self.read_int()self.year = self.read_int16()self.month = self.read_int16()self.day = self.read_int16()self.hour = self.read_int16()self.minute = self.read_int16()self.second = self.read_int16()self.radartime = datetime.datetime(self.year, self.month,self.day, self.hour, self.minute, self.second)ObsSeconds = self.read_int()ObsDates = self.read_int16()GenDates = self.read_int16()GenSeconds = self.read_int()self.lat_min = self.read_int() / 1000self.lon_min = self.read_int() / 1000self.lat_max = self.read_int() / 1000self.lon_max = self.read_int() / 1000cx = self.read_int() / 1000cy = self.read_int() / 1000# Read data informationself.nX = self.read_int() # Number of data columnsself.nY = self.read_int() # Number of data rowsself.dx = self.read_int() / 10000self.dy = self.read_int() / 10000self.height = self.read_int16()self.compress = self.read_int16()self.num_of_radars = self.read_int()self.unZipBytes = self.read_int()self.scale = self.read_int16()unUsed = self.read_int16()RgnID = self.read_string(8)self.units = self.read_string(8)reserved = self.read_string(8)def read_data(self):"""Read data bytes by skipping 256 bytes of data header"""if self.compress == 0:return self.buffer[256:]else:return self.unbzip2_decompress(self.buffer[256:])def unbzip2_decompress(self, data):"""Decompress data with bz2"""try:return bz2.decompress(data)except Exception as e:logger.error("Error decompressing bz2 data: {}".format(e))return None
需要注意的是数据可能是bz2压缩过的，通过头文件中的 compress 信息判断，如果是数据有压缩需要进行解压。
BinaryReader类初始化时最重要的参数是数据文件中读出来的所有字节 buffer：
fn = 'D:/Temp/binary/radar/ACHN_CREF_20240718_234800.BIN'with open(fn, 'rb') as f:buffer = f.read()br = BinaryReader(buffer, little_endian=True)
然后用BinaryReader类对象br的read_data()函数读取数据（为Python字节数组），用numeric模块中的frombuffer函数将字节数组转为NDArray一维short数组（frombuffer函数在MeteoInfo 3.9.2版本才有），然后根据头文件中的nY, nX信息转为二维数据，还需要除以 scale并对Y轴反向得到雷达反射率数据，小于-20的值赋为nan不在图形中绘制。
br = BinaryReader(buffer, little_endian=True)data = br.read_data()data = np.frombuffer(data, dtype=dtype.short)data = data.reshape(br.nY, br.nX)data = data / br.scaledata = data[::-1]data[data<-20] = nan
头文件中包含了数据的经纬度范围，结合经纬度方向的格点数可以将数据的二维格点坐标lon, lat构建出来：
lon = np.linspace(br.lon_min, br.lon_max, br.nX)lat = np.linspace(br.lat_min, br.lat_max, br.nY)
剩下就是绘图了。
完整代码如下：
import structimport bz2class BinaryReader:"""Radar data reader"""def __init__(self, buffer, little_endian=True):self.buffer = bufferself.index = 0self.little_endian = little_endianself.read_header()def seek(self, offset, origin=0):if origin == 0: # Set absolute positionself.index = offsetelif origin == 1: # Set relative to current positionself.index += offsetelif origin == 2: # Set relative to end of fileself.index = len(self.buffer) + offsetdef read_bytes(self, length):data = self.buffer[self.index:self.index + length]self.index += lengthreturn datadef read_string(self, length):fmt = '{}s'.format(length)value, = struct.unpack_from(fmt, self.buffer, self.index)value.replace('\x00', '').replace(' ', '')self.index += struct.calcsize(fmt)return valuedef read_int(self):fmt = '<i' if self.little_endian else '>i'value, = struct.unpack_from(fmt, self.buffer, self.index)self.index += struct.calcsize(fmt)return valuedef read_int16(self):fmt = '<h' if self.little_endian else '>h'value, = struct.unpack_from(fmt, self.buffer, self.index)self.index += struct.calcsize(fmt)return valuedef read_header(self):"""Read data header"""self.label = self.read_string(4) # Data file IDself.version = self.read_string(4) # Data file versionself.file_bytes = self.read_int()self.mosaic_id = self.read_int16()self.coordinate = self.read_int16()self.varname = self.read_string(8)self.description = self.read_string(64)BlockPos = self.read_int()BlockLen = self.read_int()self.TimeZone = self.read_int()self.year = self.read_int16()self.month = self.read_int16()self.day = self.read_int16()self.hour = self.read_int16()self.minute = self.read_int16()self.second = self.read_int16()self.radartime = datetime.datetime(self.year, self.month,self.day, self.hour, self.minute, self.second)ObsSeconds = self.read_int()ObsDates = self.read_int16()GenDates = self.read_int16()GenSeconds = self.read_int()self.lat_min = self.read_int() / 1000self.lon_min = self.read_int() / 1000self.lat_max = self.read_int() / 1000self.lon_max = self.read_int() / 1000cx = self.read_int() / 1000cy = self.read_int() / 1000# Read data informationself.nX = self.read_int() # Number of data columnsself.nY = self.read_int() # Number of data rowsself.dx = self.read_int() / 10000self.dy = self.read_int() / 10000self.height = self.read_int16()self.compress = self.read_int16()self.num_of_radars = self.read_int()self.unZipBytes = self.read_int()self.scale = self.read_int16()unUsed = self.read_int16()RgnID = self.read_string(8)self.units = self.read_string(8)reserved = self.read_string(8)def read_data(self):"""Read data bytes by skipping 256 bytes of data header"""if self.compress == 0:return self.buffer[256:]else:return self.unbzip2_decompress(self.buffer[256:])def unbzip2_decompress(self, data):"""Decompress data with bz2"""try:return bz2.decompress(data)except Exception as e:logger.error("Error decompressing bz2 data: {}".format(e))return Nonefn = 'D:/Temp/binary/radar/ACHN_CREF_20240718_234800.BIN'with open(fn, 'rb') as f:buffer = f.read()br = BinaryReader(buffer, little_endian=True)data = br.read_data()data = np.frombuffer(data, dtype=dtype.short)data = data.reshape(br.nY, br.nX)data = data / br.scaledata = data[::-1]data[data<-20] = nanlon = np.linspace(br.lon_min, br.lon_max, br.nX)lat = np.linspace(br.lat_min, br.lat_max, br.nY)#Plotaxesm()geoshow('cn_province', edgecolor='gray')geoshow('country')levs = [10,15,20,25,30,35,40,45,50,55,60,65,70]imshow(lon, lat, data, levs, cmap='radar_2', vmin=0, vmax=75)colorbar(shrink=0.8)title('CREF - {}'.format(br.radartime))xlim(70, 140)ylim(15, 55)#Add south China Seaaxesm(position=[0.76,0.18,0.15,0.2], facecolor='w', axison=False)geoshow('cn_border')xlim(106, 123)ylim(2, 23)
