反軌跡線通過網格機率分布圖

Table of contents

背景

terr.ncl與pm10.ncl差異說明

對照表

項目pm10.nclterr.ncl說明
輸入檔cmaq標準輸出檔含有經緯度的m3nc檔前者需另由GRIDCRO2D檔案讀取經緯度;後者模板為tmplateD1_27km.nc
時間序列 
外加軌跡線後者有軌跡方向之迴圈
取log10等值界線也隨之差異
中國省界底圖bou2_4p.shpbou2_4p.shp 
台灣縣市底圖小範圍需要
國界後者範圍未涉及國界

輸入檔案差異

  • prob.nc檔案內含經緯度,不需另外給定。
kuang@125-229-149-182 /Users/Data/cwb/e-service/btraj_WRFnests
$ diff ./kmean_FG123/terr.ncl ~/NCL_scripts/contour_with_basemap/pm10.ncl
155,167c155,166
<   
< ; a = addfile("/nas1/backup/data/cwb/e-service/btraj_WRFnests/tmplateD1_27km.nc","r")
<   str=(/"prob",path,".nc"/)
<   a = addfile(str_concat(str),"r")
<   lt = a->LAT(:,:)
<   ln = a->LON(:,:)
<   u = a->O(0,0,:,:)  
<   u@lat2d = lt
<   u@lon2d = ln
< 
<   wks   = gsn_open_wks ("png", path )  ; send graphics to PNG file
---
>   cmaqfile = addfile("PM10.nc","r")
>   a = addfile("GRIDCRO2D_1804_run5.nc","r")
>   lt = a->LAT(0,0,:,:)
>   ln = a->LON(0,0,:,:)
> 
>   do t = 0,2 ;15 
>   pm10 = log10(cmaqfile->PM10(t,0,:,:))  
>   pm10@lat2d = lt
>   pm10@lon2d = ln
>   st=tostring_with_format(t,"%3.3d")
>   fname=str_concat((/"pm10",st/))
>   wks   = gsn_open_wks ("png", fname )  ; send graphics to PNG file

標題、界限

181c180
<   res@tiMainString   = str_concat((/"Trajectory Cluster from ",path/))
---
>   res@tiMainString   = "2018/3/31-4/8 PM10 Episode, Log10 ug/m3"
185,198c184,188
<   pth=(/"SH", "BJ", "BH", "SC", "SW", "LOCAL", "MC"/)
<   idx=get1Dindex(pth,path)
<   bds=(/(/20,40,100,130/),(/15,max(u@lat2d),105,135/),(/20,max(u@lat2d),110,135/),(/10,35,100,130/),\
< 	(/10,30,110,135/),(/24.8,25.5,120.8,122/),(/20,40,100,130/)/)
< ;  res@mpMinLatF            = min(u@lat2d)        ; zoom in on map
< ;  res@mpMaxLatF            = max(u@lat2d)
< ;  res@mpMinLonF            = min(u@lon2d)
< ;  res@mpMaxLonF            = max(u@lon2d)
<   res@mpMinLatF            = bds(idx,0);min(u@lat2d)        ; zoom in on map
<   res@mpMaxLatF            = bds(idx,1);max(u@lat2d)
<   res@mpMinLonF            = bds(idx,2);min(u@lon2d)
<   res@mpMaxLonF            = bds(idx,3);max(u@lon2d)
< 
<   res@cnLinesOn          = False        ; turn off contour lines
---
>   res@mpMinLatF            = min(pm10@lat2d)        ; zoom in on map
>   res@mpMaxLatF            = max(pm10@lat2d)
>   res@mpMinLonF            = 60;min(pm10@lon2d)
>   res@mpMaxLonF            = max(pm10@lon2d)

等值界線的劃分方式與底圖

>   res@cnLevels           = (/-1,-0.5,0,0.5,1,1.5,2,2.5,3/)
200,218c190,195
<   if (path .eq. "LOCAL") then
<     res@mpOutlineOn           = False
<     filename="/nas1/backup/data/cwb/e-service/btraj_WRFnests/ncl_scripts/shp/TWN_COUNTY.shp"
<     res@cnLevels    = (/0.001, .002,.003,.004, .005, 0.01,.1/);(/0.001, .005, 0.01,.02,.03,.04, .05, .06,.07,.08,.09,0.1/)   ; set levels
<     res@pmTickMarkDisplayMode = "Always"            ; turn on built-in tickmarks
<   else
<     res@mpDataBaseVersion       = "MediumRes"
<     res@mpDataSetName           = "Earth..4"
<     res@mpAreaMaskingOn         = True
<     res@mpOutlineBoundarySets = "National"
<     res@mpMaskAreaSpecifiers    = (/"China","Taiwan","Disputed area between India and China","India:Arunachal Pradesh"/)
<     filename="/nas1/backup/data/cwb/e-service/btraj_WRFnests/ncl_scripts/shp/bou2_4p.shp"
<     res@cnLevels    = (/0.001, 0.01, 0.1,.12,.14,.16, .18, 0.2, 0.5,1,1.5,2/)   ; set levels
<   end if
< ;read the polylines segments
<   lonlat=asciiread(str_concat((/path,".csv"/)),(/10,3/),"float")
<   gres                  = True                ; polyline mods desired
<   gres@gsLineThicknessF = 4.0                 ; line thickness
<   gres@gsLineColor      = "Red"               ; line color 
---
>   res@cnLinesOn          = False        ; turn off contour lines
>   res@mpDataBaseVersion       = "MediumRes"
>   res@mpDataSetName           = "Earth..4"
>   res@mpAreaMaskingOn         = True
>   res@mpOutlineBoundarySets = "National"
>   res@mpMaskAreaSpecifiers    = (/"China","Taiwan","Disputed area between India and China","India:Arunachal Pradesh"/)

輸出圖形

220,221c197,198
<   plot = gsn_csm_contour_map(wks, u, res)   
<   gsn_polyline (wks, plot, lonlat(:,0), lonlat(:,1), gres)
---
>   filename="bou2_4p.shp"
>   plot = gsn_csm_contour_map(wks, pm10, res)   
233c210

富貴角與金門反軌跡NCL程式差異

  • 經調整後2者有相同的邊界與等值界線劃分

邊界調整

kuang@125-229-149-182 /Users/Data/cwb/e-service/btraj_WRFnests
$ diff ./kmean_FG123/terr.ncl ./kmean_JM123/terr.ncl
188,192c188,192
< 	(/10,30,110,135/),(/24.8,25.5,120.8,122/),(/20,40,100,130/)/)
< ;  res@mpMinLatF            = min(u@lat2d)        ; zoom in on map
< ;  res@mpMaxLatF            = max(u@lat2d)
< ;  res@mpMinLonF            = min(u@lon2d)
< ;  res@mpMaxLonF            = max(u@lon2d)
---
> 	(/10,30,110,135/),(/24.3,25,117.8,118.7/),(/20,40,100,130/)/)
>   res@mpMinLatF            = min(u@lat2d)        ; zoom in on map
>   res@mpMaxLatF            = max(u@lat2d)
>   res@mpMinLonF            = min(u@lon2d)
>   res@mpMaxLonF            = max(u@lon2d)

機率值界線調整

202,204c204
<     filename="/nas1/backup/data/cwb/e-service/btraj_WRFnests/ncl_scripts/shp/TWN_COUNTY.shp"
<     res@cnLevels    = (/0.001, .002,.003,.004, .005, 0.01,.1/);(/0.001, .005, 0.01,.02,.03,.04, .05, .06,.07,.08,.09,0.1/)   ; set levels
---
>     res@cnLevels    = (/0.001, .002,.003,.004, .005, 0.01,.1/)   ; set levels

Table of contents