The site data draw contour needs to first interpolate the site data to the lattice point data, meteoinfo the inverse distance weights (IDW) and Cressman two methods, wherein IDW method can have the option of interpolation radius. Here The example reads a micaps first class of data (ground full-feature observations), obtains 6-hour accumulated precipitation data (precipitation6h), and then uses the Griddata function of the site data to interpolate the site data into the grid data, Use the CONTOURFM function to create contour fill layers (contour intervals and colors can be customized).
Script (latitude and longitude projection):
#Set Data FoldersBasedir ='D:/myprogram/distribution/java/meteoinfo/meteoinfo'DataDir= Os.path.join (Basedir,'Sample/micaps') Mapdir= Os.path.join (Basedir,'Map')#Read shape FilesBou2_layer = Shaperead (Os.path.join (Mapdir,'bou2_4p.shp')) Bou1_layer= Shaperead (Os.path.join (Mapdir,'bou1_4l.shp')) China_layer= Shaperead (Os.path.join (Mapdir,'china.shp')) City_layer= Shaperead (Os.path.join (Mapdir,'res1_4m.shp'))#Read Station Dataf = addfile_micaps (Os.path.join (DataDir,'10101414.000')) PR= F.stationdata ('precipitation6h')#Griddata function-interpolatex = Arange (75, 135, 0.5) y= Arange (18, 55, 0.5) PRG= Pr.griddata ((x, y), method='IDW', radius=3)#Plotaxesm () geoshow (Bou2_layer, Edgecolor='Lightgray') geoshow (Bou1_layer, Facecolor= (0,0,255)) Geoshow (City_layer, Facecolor='R', Size=4, labelfield='NAME', Fontname=u'italics', Fontsize=16, yoffset=15) geoshow (china_layer, Visible=False) Levs= [0.1, 1, 2, 5, 10, 20, 25, 50, 100]cols= [(255,255,255), (170,240,255), (120,230,240), (200,220,50), (240,220,20), (255,120,10), (255,90,10), (240,40,0), (180,10,0), (120,10, 0)] Layer= Contourfm (PRG, Levs, colors=cols) Masklayer (china_layer, [layer]) Colorbar (layer) Xlim (72, 136) Ylim (16, 55) Text (,,, U'The National Precipitation Fact map', Fontname=u'blackbody', fontsize=16) Text (Up, Up, u'(2010-10-14 08:00 to 2010-10-14)', Fontname=u'blackbody', fontsize=14)#ADD South China SeaSc_layer =Bou1_layer.clone () axesm (position=[0.14,0.18,0.15,0.2], axison=False) geoshow (Sc_layer, Facecolor= (0,0,255)) Xlim (106, 123) Ylim (2, 23)
Script (Lambert projection):
#Set Data FoldersBasedir ='D:/myprogram/distribution/java/meteoinfo/meteoinfo'DataDir= Os.path.join (Basedir,'Sample/micaps') Mapdir= Os.path.join (Basedir,'Map')#Read shape FilesBou2_layer = Shaperead (Os.path.join (Mapdir,'bou2_4p.shp')) Bou1_layer= Shaperead (Os.path.join (Mapdir,'bou1_4l.shp')) China_layer= Shaperead (Os.path.join (Mapdir,'china.shp')) City_layer= Shaperead (Os.path.join (Mapdir,'res1_4m.shp'))#Read Station Dataf = addfile_micaps (Os.path.join (DataDir,'10101414.000')) PR= F.stationdata ('precipitation6h')#Griddata function-interpolatex = Arange (75, 135, 0.5) y= Arange (18, 55, 0.5) PRG= Pr.griddata ((x, y), method='IDW', radius=3)#Plotproj = Projinfo (proj='LCC', lon_0=105, lat_1=25, lat_2=47) Axesm (Projinfo=proj, position=[0, 0, 0.9, 1], Axison=false, Gridlabel=false, frameon=False) geoshow (Bou2_layer, Edgecolor='Lightgray') geoshow (Bou1_layer, Facecolor= (0,0,255)) Geoshow (City_layer, Facecolor='R', Size=4, labelfield='NAME', Fontname=u'italics', Fontsize=16, yoffset=15) geoshow (china_layer, Visible=False) Levs= [0.1, 1, 2, 5, 10, 20, 25, 50, 100]cols= [(255,255,255), (170,240,255), (120,230,240), (200,220,50), (240,220,20), (255,120,10), (255,90,10), (240,40,0), (180,10,0), (120,10, 0)] Layer= Contourfm (PRG, Levs, colors=cols) Masklayer (china_layer, [layer]) Colorbar (layer, shrink=0.5, aspect=15) axism ([78, 130, 14, 53]) text (------'The National Precipitation Fact map', Fontname=u'blackbody', fontsize=18) Text (A., Wuyi, U'(2010-10-14 08:00 to 2010-10-14)', Fontname=u'blackbody', fontsize=16)#ADD South China SeaSc_layer =Bou1_layer.clone () axesm (position=[0.1,0.05,0.15,0.2], axison=False) geoshow (Sc_layer, Facecolor= (0,0,255)) Xlim (106, 123) Ylim (2, 23)
Operation Result:
Meteoinfolab Script Example: Site data draw Contour Line