- 积分
- 1009
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2011-8-31
- 最后登录
- 1970-1-1
|
发表于 2022-3-28 23:21:57
|
显示全部楼层
#-----------------------------------------------------# Author: Yaqiang Wang
# Date: 2014-11-16
# Purpose: Read FY2E hdf5 data
# Note: Sample
#-----------------------------------------------------
from org.meteoinfo.data.meteodata import MeteoDataInfo
from org.meteoinfo.data.meteodata import Dimension
from org.meteoinfo.data.meteodata import DimensionType
from org.meteoinfo.data.meteodata import DrawMeteoData
from org.meteoinfo.projection import ProjectionInfo
from org.meteoinfo.projection import KnownCoordinateSystems
from org.meteoinfo.legend import LegendManage
from org.meteoinfo.legend import LegendType
from org.meteoinfo.shape import ShapeTypes
import os.path
import jarray
#Set data directory
dataDir = 'D:/Temp/hdf/'
#Create MeteoDataInfo object
mdi = MeteoDataInfo()
projInfo = ProjectionInfo("+proj=geos +lon_0=104.5 +h=35785864")
#Read hdf data file
#fn = dataDir + 'FY2E_CTA_MLT_NOM_20141116_0200.hdf'
fn = dataDir + 'FY2E_FDI_ALL_NOM_20110713_0100.hdf'
if os.path.isfile(fn):
print fn
mdi.openNetCDFData(fn)
dataInfo = mdi.getDataInfo()
print mdi.getInfoText()
xmin = -5750000.0
ymin = -5750000.0
xnum = 2288
ynum = 2288
xmax = 5750000.0
ymax = 5750000.0
xdelt = (xmax - xmin) / xnum
ydelt = (ymax - ymin) / ynum
xlist = []
ylist = []
for i in range(0,xnum):
xlist.append(xmin + xdelt * i)
for i in range(0,ynum):
ylist.append(ymin + ydelt * i)
X = jarray.array(xlist, 'd')
Y = jarray.array(ylist, 'd')
xDim = Dimension(DimensionType.X)
xDim.setValues(X)
dataInfo.setXDimension(xDim)
yDim = Dimension(DimensionType.Y)
yDim.setValues(Y)
dataInfo.setYDimension(yDim)
#var = dataInfo.getVariable('FY2E_CLM_Hourly_Cloud_Amount')
var = dataInfo.getVariable('NOMChannelVIS')
dimList = [yDim, xDim]
var.setDimensions(dimList)
print 'Get grid data...'
gData = mdi.getGridData(var.getName())
gData.yReverse()
gData.projInfo = projInfo
#gData.missingValue = 0.0
print 'Reproject grid data to Lon/Lat...'
toProj = KnownCoordinateSystems.geographic.world.WGS1984
gData = gData.project(toProj)
print 'Ouput grid data to Sufer ASCII grid data file...'
outfn = 'D:/Temp/test/NOMChannelVIS.grd'
gData.saveAsSurferASCIIFile(outfn)
print 'Create raster layer from the grid data...'
aLS = LegendManage.createLegendSchemeFromGridData(gData, LegendType.GraduatedColor, ShapeTypes.Polygon)
aLayer = DrawMeteoData.createRasterLayer(gData, "Test_FY2E", aLS)
#aLayer = DrawMeteoData.createShadedLayer(gData, aLS, 'Test_HDF_Shaded', 'AOD', True)
aLayer.setProjInfo(toProj)
#aLayer.setProjInfo(projInfo)
print 'Add the layer...'
mf = miapp.getMapDocument().getActiveMapFrame()
mf.addLayer(aLayer)
mf.moveLayer(aLayer, 0)
print 'Finished!' |
|