MeteoInfoLab脚本示例：MODIS AOD

MODIS的气溶胶光学厚度（AOD）产品应用很广，数据可以在Giovanni上下载：http://disc.sci.gsfc.nasa.gov/giovanni/overview/index.html。有HDF, netCDF和文本格式，这里示例打开netCDF格式数据，提取AOD数据并绘图。

脚本程序如下：

f = addfile('D:/Temp/nc/MOD08_D3.A2015121.051.2015122103938.pscs_000500931513.Optical_Depth_Land_And_Ocean_Mean.G3.nc')

v_aod = f['Optical_Depth_Land_And_Ocean_Mean']

aod = v_aod[[33,41],[113,122]]

axesm()

china = shaperead('D:/Temp/map/bou2_4p.shp')

geoshow(china, linecolor='gray')

world = shaperead('D:/Temp/map/country1.shp')

geoshow(world)

levels = arange(0, 1.0, 0.02)

layer = contourfm(aod, levels, interpolate=True)

title('Aerosol Optical Depth at 550 nm')

colorbar(layer)

show()

要和网站上的图形一致的话需要在等值线追踪前将格点数据插值成2倍的分辨率（layer = contourfm(aod, levels, interpolate=True)）。

HDF格式数据文件中的经纬度维有问题，需要在脚本中指定，示例如下：

f = addfile('D:/Temp/hdf/MOD08_D3.A2015121.051.2015122103938.pscs_000500931513.Optical_Depth_Land_And_Ocean_Mean.G3.hdf')

v_aod = f['Optical_Depth_Land_And_Ocean_Mean']

#Set x/y

x = linspace(-180.0,180.0, 360)

y = linspace(-90.0,90.0, 180)

#Set x/y dimensions

v_aod.setdim('Y', y)

v_aod.setdim('X', x)

#Get data

aod = v_aod[[33,41],[113,122]]

#Plot

axesm()

china = shaperead('D:/Temp/map/bou2_4p.shp')

geoshow(china, linecolor='gray')

world = shaperead('D:/Temp/map/country1.shp')

geoshow(world)

levels = arange(0, 1.0, 0.02)

layer = contourfm(aod, levels, interpolate=True)

title('Aerosol Optical Depth at 550 nm')

colorbar(layer)

show()

巴特西

MeteoInfoLab脚本示例：MODIS AOD

最新文章

热门文章