- 积分
- 55945
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2011-6-21
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
本帖最后由 MeteoInfo 于 2016-5-16 10:36 编辑
FY2C-E的gpf格式是二进制数据,由数据头、定标表、通道数据三部分组成,这里给出一个FY2E gpf文件(经纬度投影)读取、绘图的示例程序。先利用Python读取二进制数据文件的功能读取文件头中相关信息,然后利用MeteoInfoLab的binread函数读取通道数据至二维数组,再绘图。
- import struct
- fn = 'D:/Temp/binary/FY2E_2011_06_25_00_01_E_PJ3.gpf'
- #Read data header parameters
- f = open(fn, 'rb')
- fileid, = struct.unpack('2s', f.read(2))
- version, = struct.unpack('<h', f.read(2))
- satid, = struct.unpack('<h', f.read(2))
- year, = struct.unpack('<h', f.read(2))
- month, = struct.unpack('<h', f.read(2))
- day, = struct.unpack('<h', f.read(2))
- hour, = struct.unpack('<h', f.read(2))
- minute, = struct.unpack('<h', f.read(2))
- chnums, = struct.unpack('<h', f.read(2))
- pjtype, = struct.unpack('<h', f.read(2))
- width, = struct.unpack('<h', f.read(2))
- height, = struct.unpack('<h', f.read(2))
- clonres, = struct.unpack('<f', f.read(4))
- clatres, = struct.unpack('<f', f.read(4))
- stdlat1, = struct.unpack('<f', f.read(4))
- stdlat2, = struct.unpack('<f', f.read(4))
- earthr, = struct.unpack('<f', f.read(4))
- minlat, = struct.unpack('<f', f.read(4))
- maxlat, = struct.unpack('<f', f.read(4))
- minlon, = struct.unpack('<f', f.read(4))
- maxlon, = struct.unpack('<f', f.read(4))
- ltlat, = struct.unpack('<f', f.read(4))
- ltlon, = struct.unpack('<f', f.read(4))
- rtlat, = struct.unpack('<f', f.read(4))
- rtlon, = struct.unpack('<f', f.read(4))
- lblat, = struct.unpack('<f', f.read(4))
- lblon, = struct.unpack('<f', f.read(4))
- rblat, = struct.unpack('<f', f.read(4))
- rblon, = struct.unpack('<f', f.read(4))
- stdlon, = struct.unpack('<f', f.read(4))
- centerlon, = struct.unpack('<f', f.read(4))
- centerlat, = struct.unpack('<f', f.read(4))
- chindex = []
- for i in range(chnums):
- chindex.append(struct.unpack('b', f.read(1))[0])
- f.read(128 - chnums)
- plonres, = struct.unpack('<f', f.read(4))
- platres, = struct.unpack('<f', f.read(4))
- f.read(1808)
- #Read calibration table
- f.read(32768)
- f.close()
- #Read one channel data
- cn = 1 #Infrared channel 1
- skipn = 2048 + 32768
- for i in range(1, cn):
- if chindex[i-1] < 5: #Infrared channel
- byten = 2
- else:
- byten = 1 #Visible light channel
- skipn += width * height * byten
- if cn < 5:
- data = binread(fn, [height, width], 'short', skip=skipn)
- else:
- data = binread(fn, [height, width], 'byte', skip=skipn)
- data = data[::-1,:]
- lon = linspace(ltlon, rtlon, width)
- lat = linspace(lblat, ltlat, height)
- #Plot
- axesm()
- lworld = shaperead('D:/Temp/Map/country1.shp')
- lchina = shaperead('D:/Temp/Map/bou2_4p.shp')
- geoshow(lchina, edgecolor='b')
- geoshow(lworld, edgecolor='b')
- cmap = 'MPL_gist_gray'
- layer = imshowm(lon, lat, data, 20, cmap=cmap)
- colorbar(layer)
- t = datetime.datetime(year, month, day, hour, minute)
- title('FY2E (' + t.strftime('%Y-%m-%d') + ')')
|
|