MeteoInfo脚本示例：GrADS to netCDF

这里给出一个将GrADS数据文件转为netCDF数据文件的脚本示例程序，其它格式数据转netCDF可以参考：

#-----------------------------------------------------
# Author: Yaqiang Wang
# Date: 2015-3-12
# Purpose: Convert CUACE grads data to netCDF (CUACE/Dust)
# Note: Sample
#-----------------------------------------------------
print ‘Loading classes...‘
from org.meteoinfo.data import GridData
from org.meteoinfo.data import DataMath
from org.meteoinfo.data.meteodata import MeteoDataInfo
from org.meteoinfo.geoprocess.analysis import ResampleMethods
from org.meteoinfo.data.meteodata.netcdf import NetCDFDataInfo
from org.meteoinfo.projection import ProjectionInfo
from org.meteoinfo.projection import ProjectionNames
from org.meteoinfo.projection import KnownCoordinateSystems
from ucar.nc2 import NetcdfFileWriter
from ucar.nc2 import Attribute
from ucar.ma2 import DataType
from ucar.ma2 import Array
import os.path
import jarray
import datetime
from java.util import Date
from java.text import SimpleDateFormat

#Set date
year = 2014
month = 4
day = 23
hour = 0
sdate = datetime.datetime(year, month, day, hour)
print sdate
#Set directory
dataDir = ‘U:/data/cuace_dust/dust_example/2014_case‘
outDir = dataDir
infn = os.path.join(dataDir, ‘dust_post_‘+ sdate.strftime(‘%Y%m%d%H‘) + ‘.ctl‘)
outfn = os.path.join(dataDir, ‘WMO_SDS-WAS_Asian_Center_Model_Forecasting_CUACE-Dust_CMA_‘     + sdate.strftime(‘%Y-%m-%d‘) + ‘.nc‘)
#Set output X/Y coordinates and projection
toProjInfo = KnownCoordinateSystems.geographic.world.WGS1984
sx = 70.0
xnum = 161
sy = 20
ynum = 71
delta = 0.5
xlist = []
ylist = []
for i in range(0, xnum):
    xlist.append(sx + delta * i)
for i in range(0, ynum):
    ylist.append(sy + delta * i)
X = jarray.array(xlist, ‘d‘)
Y = jarray.array(ylist, ‘d‘)

#Read GrADS data file
print ‘Open GrADS data file...‘
mdi = MeteoDataInfo()
mdi.openGrADSData(infn)
dataInfo = mdi.getDataInfo()
dataInfo.setBigEndian(True)
fromProjInfo = mdi.getProjectionInfo()
tnum = dataInfo.getTimeNum()
mvalue = dataInfo.getMissingValue()

#Set output nc data file
print ‘Create output NC file: ‘ + outfn
ncfile = NetcdfFileWriter.createNew(NetcdfFileWriter.Version.netcdf3, outfn)
#Add dimensions
print ‘Add dimensions...‘
xDim = ncfile.addDimension(None, ‘lon‘, xnum)
yDim = ncfile.addDimension(None, ‘lat‘, ynum)
tDim = ncfile.addDimension(None, ‘time‘, tnum)

#Add global attributes
print ‘Add global attributes...‘
ncfile.addGroupAttribute(None, Attribute(‘Conventions‘, ‘CF-1.6‘))
ncfile.addGroupAttribute(None, Attribute(‘Title‘, ‘Sand and dust storm forecasting‘))
ncfile.addGroupAttribute(None, Attribute(‘Model‘, ‘CUACE/Dust‘))
ncfile.addGroupAttribute(None, Attribute(‘Center‘, ‘WMO SDS-WAS Asian Center‘))
ncfile.addGroupAttribute(None, Attribute(‘Agency‘, ‘China Meteorological Administration‘))

#Add variables
xvar = ncfile.addVariable(None, ‘lon‘, DataType.FLOAT, [xDim])
xvar.addAttribute(Attribute(‘units‘, ‘degrees_east‘))
xvar.addAttribute(Attribute(‘long_name‘, ‘Longitude‘))
xvar.addAttribute(Attribute(‘standard_name‘, ‘longitude‘))
xvar.addAttribute(Attribute(‘axis‘, ‘X‘))
yvar = ncfile.addVariable(None, ‘lat‘, DataType.FLOAT, [yDim])
yvar.addAttribute(Attribute(‘units‘, ‘degrees_north‘))
yvar.addAttribute(Attribute(‘long_name‘, ‘Latitude‘))
yvar.addAttribute(Attribute(‘standard_name‘, ‘latitude‘))
yvar.addAttribute(Attribute(‘axis‘, ‘Y‘))
tvar = ncfile.addVariable(None, ‘time‘, DataType.INT, [tDim])
tvar.addAttribute(Attribute(‘units‘, ‘hours since 1900-01-01 00:00:0.0‘))
tvar.addAttribute(Attribute(‘long_name‘, ‘Time‘))
tvar.addAttribute(Attribute(‘standart_name‘, ‘time‘))
tvar.addAttribute(Attribute(‘axis‘, ‘T‘))
#Data variables
vnames = [‘load‘,‘con‘,‘dry‘,‘wet‘,‘aod‘]
varlist = []
var = ncfile.addVariable(None, ‘LOAD_DUST‘, DataType.FLOAT, [tDim, yDim, xDim])
var.addAttribute(Attribute(‘long_name‘, ‘Dust load‘))
var.addAttribute(Attribute(‘units‘, ‘kg/m2‘))
var.addAttribute(Attribute(‘missing_value‘, -9999.0))
varlist.append(var)
var = ncfile.addVariable(None, ‘SCONC_DUST‘, DataType.FLOAT, [tDim, yDim, xDim])
var.addAttribute(Attribute(‘long_name‘, ‘Surface dust concentration‘))
var.addAttribute(Attribute(‘units‘, ‘ug/m3‘))
var.addAttribute(Attribute(‘missing_value‘, -9999.0))
varlist.append(var)
var = ncfile.addVariable(None, ‘DDEPO_DUST‘, DataType.FLOAT, [tDim, yDim, xDim])
var.addAttribute(Attribute(‘long_name‘, ‘3-hour accumulated dry deposition‘))
var.addAttribute(Attribute(‘units‘, ‘kg/m2‘))
var.addAttribute(Attribute(‘missing_value‘, -9999.0))
varlist.append(var)
var = ncfile.addVariable(None, ‘WDEPO_DUST‘, DataType.FLOAT, [tDim, yDim, xDim])
var.addAttribute(Attribute(‘long_name‘, ‘3-hour accumulated wet deposition‘))
var.addAttribute(Attribute(‘units‘, ‘kg/m2‘))
var.addAttribute(Attribute(‘missing_value‘, -9999.0))
varlist.append(var)
var = ncfile.addVariable(None, ‘AOD550_DUST‘, DataType.FLOAT, [tDim, yDim, xDim])
var.addAttribute(Attribute(‘long_name‘, ‘Dust optical depth at 550nm‘))
var.addAttribute(Attribute(‘units‘, ‘-‘))
var.addAttribute(Attribute(‘missing_value‘, -9999.0))
varlist.append(var)

#Write nc file
ncfile.create()
#Write x,y,z,t variables
print ‘Write x variable...‘
shape = jarray.array([xnum], ‘i‘)
data = Array.factory(DataType.FLOAT, shape)
for i in range(0,xnum):
    data.set(i, X[i])
ncfile.write(xvar, data)

print ‘Write y variable...‘
shape = jarray.array([ynum], ‘i‘)
data = Array.factory(DataType.FLOAT, shape)
for i in range(0,ynum):
    data.set(i, Y[i])
ncfile.write(yvar, data)

print ‘Write time variable...‘
format = SimpleDateFormat(‘yyyy-MM-dd‘)
sdate = format.parse(‘1900-01-01‘)
tvalues = dataInfo.getTimeValues(sdate, ‘hours‘)
shape = jarray.array([tnum], ‘i‘)
data = Array.factory(DataType.INT, shape)
for i in range(0,tnum):
    data.set(i, tvalues[i])
ncfile.write(tvar, data)

#Write data variables
print ‘Write data variable...‘
for vname, var in zip(vnames, varlist):
    for t in range(0, tnum):
        print ‘Time: ‘ + str(t + 1)
        mdi.setTimeIndex(t)
        gData = mdi.getGridData(vname)
        ngData = gData.project(fromProjInfo, toProjInfo, X, Y, ResampleMethods.Bilinear)
        origin = jarray.array([t, 0, 0], ‘i‘)
        ncfile.write(var, origin, NetCDFDataInfo.gridToArray3D(ngData))

#Close nc file
ncfile.flush()
ncfile.close()

print ‘Finished‘

　　

上面转换的netCDF文件绘制模式结果和地面天气现象观测叠加动画图的示例脚本：

# coding=utf-8
#-----------------------------------------------------
# Author: Yaqiang Wang
# Date: 2015-3-13
# Purpose: Read CUACE/Dust netCDF data and MICAPS observation data to plot figures
# Note: Sample
#-----------------------------------------------------
print ‘Loading classes...‘
from org.meteoinfo.layout import MapLayout
from org.meteoinfo.data import GridData
from org.meteoinfo.data.meteodata import MeteoDataInfo, DrawMeteoData
from org.meteoinfo.legend import LegendScheme
from org.meteoinfo.shape import ShapeTypes
from org.meteoinfo.global.image import AnimatedGifEncoder
import os.path
import jarray
import datetime
import sys
from java.util import Date, Calendar, Locale
from java.text import SimpleDateFormat
from java.awt import Color

#Set date
year = 2013
month = 3
day = 1
hour = 0
sdate = datetime.datetime(year, month, day, hour)

#sdate = datetime.date.today()
#if len(sys.argv) >= 2:
#    sdate = sdate - datetime.timedelta(days=int(sys.argv[1]))
#    sdate = sdate + datetime.timedelta(days=1)
print sdate
dformat = SimpleDateFormat(‘HH dd MMM yyy‘, Locale.ENGLISH)
dformat1 = SimpleDateFormat(‘yyMMddHH‘)
cal = Calendar.getInstance()

#Set model
#model = ‘CUACE-DUST_CMA‘
model = ‘ADAM2_KMA‘

#Set directory
dataDir = ‘D:/Working/2015/International/SDS_Asian_Region_Center/Model_Verification‘
obsDir = ‘U:/data/micaps/2014/plot‘
obsDir = ‘E:/MetData/micaps‘
runDir = dataDir
outDir = os.path.join(dataDir, ‘figure‘)
if not os.path.exists(outDir):
    os.mkdir(outDir)
#Set input/output file names
infn = os.path.join(dataDir, ‘WMO_SDS-WAS_Asian_Center_Model_Forecasting_‘ + model + ‘_‘     + sdate.strftime(‘%Y-%m-%d‘) + ‘.nc‘)
projfn = os.path.join(runDir, ‘sds_asian.mip‘)

#Plot data
print ‘Plot data...‘
mapLayout = MapLayout()
mapLayout.loadProjectFile(projfn)
mf = mapLayout.getActiveMapFrame()
title = mapLayout.getTexts().get(2)
legend = mapLayout.getLegends()[0]

#---- Set weather list - sand and dust storm
weathers = [6, 7, 8, 9, 30, 31, 32, 33, 34, 35]
#---- Set weather list - sand and dust storm and haze
#weathers = [5, 6, 7, 8, 9, 30, 31, 32, 33, 34, 35]

#---- Create MeteoDataInfo object
mdi = MeteoDataInfo()
omdi = MeteoDataInfo()

#---- Plot loop
mdi.openNetCDFData(infn)
lsfn = os.path.join(runDir,‘dust_conc.lgs‘)
print ‘Read data file: ‘ + infn
aLS = LegendScheme(ShapeTypes.Polygon)
aLS.importFromXMLFile(lsfn)
tnum = mdi.getDataInfo().getTimeNum()
#tnum = 3
s = ‘SCONC_DUST‘
giffn = os.path.join(outDir, ‘V_‘ + s + ‘_‘ + model + ‘_‘ + sdate.strftime(‘%Y%m%d‘) + ‘--loop-.gif‘)
print giffn
encoder = AnimatedGifEncoder()
encoder.setRepeat(0)
encoder.setDelay(1000)
encoder.start(giffn)
sTime = mdi.getDataInfo().getTimes().get(0)
for t in range(1, tnum):
    mdi.setTimeIndex(t)
    aTime = mdi.getDataInfo().getTimes().get(t)
    cal.setTime(aTime)
    cal.add(Calendar.HOUR, 8)
    bjTime = cal.getTime()
    #---- Open observation weather data
    obsfn = os.path.join(obsDir, dformat1.format(bjTime) + ‘.000‘)
    print obsfn
    if not os.path.exists(obsfn):
        continue
    omdi.openMICAPSData(obsfn)
    wData = omdi.getStationData(‘WeatherNow‘)
    weatherLayer = DrawMeteoData.createWeatherSymbolLayer(wData, weathers, ‘Weather‘)
    #for lb in weatherLayer.getLegendScheme().getLegendBreaks():
    #    lb.setColor(Color.red)
    weatherLayer.setAvoidCollision(False)
    mf.removeMeteoLayers()
    mf.addLayer(weatherLayer)
    #---- Get grid data and create a shaded layer
    gData = mdi.getGridData(s)
    aLayer = DrawMeteoData.createShadedLayer(gData, aLS, ‘Forecasting_‘ + s, ‘Data‘, True)
    aLayer.setProjInfo(mdi.getProjectionInfo())
    mf.addLayer(aLayer)
    mf.moveLayer(aLayer, 0)
    #---- Set title
    title.setLabelText(‘Run: ‘ + dformat.format(sTime) + ‘    Valid: ‘ + dformat.format(aTime)         + ‘(H+‘ + str(t * 3) + ‘)‘)
    #---- Set legend
    legend.setLegendLayer(aLayer)
    mapLayout.paintGraphics()
    encoder.addFrame(mapLayout.getViewImage())
    figurefn = os.path.join(outDir, ‘V_‘ + model + ‘_‘ + s + ‘_‘ + dformat1.format(aTime) + ‘.png‘)
    print ‘Output figure: ‘ + figurefn
    mapLayout.exportToPicture(figurefn)

encoder.finish()
print ‘Finished‘

　　

时间： 2024-07-30 10:12:01

MeteoInfo脚本示例：GrADS to netCDF的相关文章

MeteoInfo脚本示例：读取FY3A AOD HDF文件

FY3A卫星有AOD产品数据,HDF格式,这里示例用MeteoInfo脚本程序读取和显示该类数据. 脚本程序如下: #----------------------------------------------------- # Author: Yaqiang Wang # Date: 2015-3-18 # Purpose: Read FY3A AOD hdf5 data # Note: Sample #---------------------------------------------

MeteoInfo脚本示例：Trajectory

示例读取HYSPLIT模式输出的气团轨迹数据文件,生成轨迹图层,并显示轨迹各节点的气压图. 脚本程序: f = addfile_hytraj('D:/MyProgram/Distribution/java/MeteoInfo/MeteoInfo/sample/HYSPLIT/tdump') tlayer = f.trajlayer() stlayer = f.trajsplayer() axesm(position=[0.1,0.4,0.8,0.6]) mlayer = shaperead('D

MeteoInfoLab脚本示例：创建netCDF文件（合并文件）

在MeteoInfoLab中增加了创建netCDF文件并写入数据的功能,这里利用合并多个netCDF文件为一个新的netCDF文件为例. 1.创建一个可写入的netCDF文件对象(下面用ncfile表示),用addfile函数,第一个参数是文件名,第二次参数'c'表示创建新的netCDF文件.ncfile = addfile(outfn, 'c') 2.添加维(Dimensions),用ncfile的adddim函数,两个参数分别是维名称和维长度.stn = 26564stdim = ncfil

linux添加开机自启动脚本示例详解

来源: linux添加开机自启动脚本示例详解 linux下(以RedHat为范本)添加开机自启动脚本有两种方法,先来简单的; 一.在/etc/rc.local中添加如果不想将脚本粘来粘去,或创建链接什么的,则:step1. 先修改好脚本,使其所有模块都能在任意目录启动时正常执行;step2. 再在/etc/rc.local的末尾添加一行以绝对路径启动脚本的行;如:$ vim /etc/rc.local#!/bin/sh## This script will be executed *after*

SqlServer 审核（脚本示例）

此文章主要是脚本示例,更多说明看官方文档:审核(数据库引擎) -- 必须在 master 数据库中创建审核 USE master; GO -- 创建服务器审核对象 -- https://msdn.microsoft.com/zh-cn/library/cc280448(v=sql.100).aspx CREATE SERVER AUDIT [Audit_ToFile] TO FILE ( --目标类型:FILE(文件)/APPLICATION_LOG(应用程序日志)/SECURITY(安全日志)

MeteoInfoLab脚本示例：FY-3C全球火点HDF数据

FY-3C全球火点HDF数据包含一个FIRES二维变量,第一维是火点数,第二维是一些属性,其中第3.4列分别是火点的纬度和经度.下面的脚本示例读出所有火点经纬度并绘图. 脚本程序: #Add data file fn = 'D:/Temp/hdf/FY3C_VIRRX_GBAL_L2_GFR_MLT_GLL_20150811_POAD_1000M_MS.HDF' f = addfile(fn) #Get data variable v = f['FIRES'] #Get data array d

shell脚本程序中的部分常用环境变量和参数变量的说明以及简单shell脚本示例

环境变量 $HOME 当前用户的家目录 $PATH 以冒号分隔的用来搜索命令的目录列表 $PS1 命令提示符,通常是$字符,但在bash中,可以使用一些更复杂的值.例如,字符串[\[email protected]\h\w]$就是一个流行的默认值,它给出用户名/机器名和当前的目录名,当然也包括一个$提示符. $PS2 二级提示符,用来表示后续的输入,通常是 > 字符. $IFS 输入域分隔符.当shell读取输入时,它给出用来分隔单词的一组字符,他们通常是空格,制表符和换行符. $0 shell

Linux shell脚本示例（四）

Example No. 2作为shell编写人员,经常面对数据格式不一致的问题,将数据标准话后输出是一个需要解决问题,本示例以MySQL的时间为例,脚本输入月.日.年三个参数,将其标准化后输出,月份以英文标准输出,年份如果是4个数字,直接输出,如果是0~69,则年份为2000~2069,如果是70~99,则年份为1970~1999:脚本示例: #!/bin/bash #shell name:shell_format.sh #author by woon #数据标准化输出示例 Month_to_n

MeteoInfoLab脚本示例：站点数据绘制等值线

站点数据绘制等值线需要首先将站点数据插值为格点数据,MeteoInfo中提供了反距离权法(IDW)和cressman两个方法,其中IDW方法可以有插值半径的选项.这里示例读取一个MICAPS第一类数据(地面全要素观测),获取6小时累积降水数据(Precipitation6h),然后用站点数据的griddata函数将站点数据插值为格点数据,再利用contourfm函数创建等值线填色图层(等值线间隔和颜色可以自定义). 脚本程序(经纬度投影): #Set data folders basedir =