爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 5013|回复: 5

[作图] ncl数据提取作图

[复制链接]

新浪微博达人勋

发表于 2019-2-22 16:13:00 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
代码出图如下,只想提取中国区域数据。 c3dif0.png
;These files are loaded by default in NCL V6.2.0 and newer
; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"  
; load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"  
load "$NCARG_ROOT/lib/ncarg/nclscripts/shapefile/shapefile_utils.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/cnmap/cnmap.ncl"

;************************************************
f = addfile("/test/data/20181204data/c3PftFrac_Lmon_MIROC-ESM_esmControl_r1i1p1_185001-210012.nc","r")
;sname = "/shapefile/gadm36_CHN_0.shp"
var = f->c3PftFrac
    var@missing_Value = 1e+20
    var@_FillValue = 1e+20
    varyear = month_to_annual(var, 0)
vardif = var(250,:,:)-var(0,:,:)

varyear = lonFlip(varyear)
vardif = lonFlip(vardif) ; 重新排列经纬度坐标reorder
;printVarSummary(c3dif)
t=max(vardif)
;print (c3avg)
printVarSummary(varyear)
printVarSummary(vardif)
copy_VarAtts(var, vardif)   
;printMinMax(c3dif , opt)
print (t)
;************************************************
; create default plot
;************************************************
wks =gsn_open_wks("png","c3dif0")
lon1 = 72
lon2 = 136
lat1 = 17
lat2 = 55

outline_areas = (/"China","Taiwan","Disputed area between India and China","India:Arunachal Pradesh"/)

res = True

res@vpXF = 0.15
res@vpYF = 0.8
res@vpHeightF = 0.6
res@vpWidthF = 0.8

res@mpMinLatF = lat1
res@mpMaxLatF = lat2
res@mpMinLonF = lon1
res@mpMaxLonF = lon2
res@cnLinesOn = False

  res@gsnDraw             = False         ; don't draw yet
  res@gsnFrame            = False         ; don't advance frame yet
  res@gsnMaximize = True

res@mpOutlineOn = False ;-- outline land (default: False)
res@mpFillOn              = True   
res@mpOutlineSpecifiers = outline_areas ;-- which country to be outlined
res@mpPerimOn             = False               ; don't draw box around map
res@mpGridAndLimbOn       = False               ; Don't draw lat/lon lines.
res@mpLandFillColor         = "white"
res@mpInlandWaterFillColor  = "white"
res@mpOceanFillColor        = "white"
res@mpOutlineBoundarySets   = "NoBoundaries"



res@tiMainString = "Total C3 PFT Cover Fraction change of 1850-2100"
res@tiMainFontHeightF = 0.025
res@gsnMajorLatSpacing = 10
res@gsnMinorLatSpacing = 1
res@tmYLLabelStride = 1
res@tmYLLabelFontHeightF = 0.015
res@tmYLMajorLengthF = 0.02
res@tmYLMinorLengthF = 0.01
res@tmYLMajorLineColor = "black"
res@tmYLMinorLineColor = "grey20"
res@tmYLLabelFontColor = "black"

res@gsnMajorLonSpacing = 10
res@gsnMinorLonSpacing = 1
res@tmXBLabelStride = 1
res@tmXBLabelFontHeightF = 0.015
res@tmXBMajorLengthF = 0.015
res@tmXBMinorLengthF = 0.01
res@tmXBMajorLineColor = "black"
res@tmXBMinorLineColor = "grey20"
res@tmXBLabelFontColor = "black"

res@gsnRightString = "%"

String = "MIROC-ESM"
ndcres = True
ndcres@txFontColor = "black"
ndcres@txFontHeightF = 0.015
gsn_text_ndc(wks, String, 0.8, 0.12, ndcres)
res@cnFillOn = True
res@cnFillPalette = "BlRe"
res@lbBoxLinesOn = False

res@lbBoxMinorExtentF = 0.2 ;-- decrease height of labelbar boxes and vp
res@lbBoxLinesOn = False ;-- don't draw lines around labelbar boxes
res@lbLabelFontHeightF = 0.015 ;-- label font height
res@lbLabelFont = "helvetica-bold";-- label font
res@lbLabelOffsetF = 0.07 ;-- move the labelbar labels downwards
res@pmLabelBarWidthF = 0.8 ;-- labelbar width; default is shorter
res@pmLabelBarHeightF = 0.1 ;-- labelbar height; default is taller
res@pmLabelBarOrthogonalPosF = 0.07 ;-- y-position (positive: downward; def: 0.02)
res@pmLabelBarParallelPosF = 0.5 ;-- x-position (CenterCenter); default: 0.5
res@cnLevelSelectionMode = "ManualLevels"
res@cnLevelSpacingF = 5
plot = gsn_csm_contour_map(wks,vardif,res)   

    ;.....................
    ;add China map
     cnres           = True
cnres@china     = True       ;draw china map or not
cnres@river     = True       ;draw changjiang&huanghe or not
cnres@province  = True       ;draw province boundary or not
cnres@nanhai    = True       ;draw nanhai or not
cnres@diqu      = False       ; draw diqujie or not
chinamap = add_china_map(wks,plot,cnres)
draw(plot)
  frame(wks)


密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-2-22 16:40:58 | 显示全部楼层
这个得利用中国边界的shp文件去判断把点筛出来,我目前只能提供这一种思路。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-2-22 16:47:26 | 显示全部楼层
mpAreaMaskingOn   mpMaskAreaSpecifiers   用mp中mask相关的属性试一下
详细参考http://www.ncl.ucar.edu/Applications/mask.shtml
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-2-22 19:07:43 | 显示全部楼层
使用中国区域shapefile文件做一下遮盖就可以了
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-2-24 16:37:54 | 显示全部楼层
mask中国区域的以外部分
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2019-2-26 10:03:56 | 显示全部楼层
我不喜欢吃生姜 发表于 2019-2-22 16:47
mpAreaMaskingOn   mpMaskAreaSpecifiers   用mp中mask相关的属性试一下
详细参考http://www.ncl.ucar.edu ...

谢谢,好的,很挺有帮助的
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表