[ncl-talk] Fwd: How to add shapefiles to nesting domain

Alan Brammer abrammer at albany.edu
Thu Jan 28 08:29:53 MST 2016


See reply below from Tao Lu.  I don’t use the WRF functions enough to know whats going on in this cake. 

Unless obvious, always reply to ncl-talk not individual responders.
 
Hopefully someone familiar with wrf_wps_dom() can comment on adding extra polylines to the output. 



> 
> Thank you very much for replying to me.
> And I am sorry for my bad description.
> 
> In the code I attached, only 
> ; Test to add a shapefile
> ; I add the following contents
>   dir         = "/home/taolu/WRF/CASES/kinugawa/shp/"
>   shp_kinu = dir + "kinugawa.shp"
>   shp_tone = dir + "tonegawa.shp"
> 
>   preres                  = True
>   preres at gsLineColor      = "gray25"
>   preres at gsLineThicknessF = 2.0  
> 
>   kinu_id = gsn_add_shapefile_polylines(wks,mp,shp_kinu,preres)
>   tone_id = gsn_add_shapefile_polylines(wks,mp,shp_tone,preres)
>   
>   draw(mp);
> is added by me.
> 
> Run the original code, I could get the domains as showed below:
> 
> 
> But after I added my code, I could get this one:
> 
> The one I want is the combination of the two.
> 
> I also tried type =“png” and mpres at gsnDraw = True, but also failed.
> Do you have any clue about this?
> 
> On Thu, Jan 28, 2016 at 4:22 AM, Alan Brammer <abrammer at albany.edu <mailto:abrammer at albany.edu>> wrote:
>> So, I edit the plotgrids_new.ncl as showed below, but I could not get what I wanted.
> What do you mean by this? 
> Is a plot produced? 
> You’re outputting to x11 so it will only be temporary.  At the top of the script you might want to change type =“png” or “pdf” then you can send the output to the script and detail what you actually wanted.
> 
> mpres at gsnDraw = False ;  may prevent the wrf_wps_dom() function from drawing allowing your to overlay your shape file polylines.   I’m not familiar with wrf_wps_dom() though, and the documentation I found online was particularly lacking. 
> 
> 
> An example of why the plot isn’t what you want, will enable someone to probably either answer or forward this onto wrfhelp at ucar.edu <mailto:wrfhelp at ucar.edu>
> 
> 
> Alan
> 
> 
> 
> ##############################
> Alan Brammer,
> Post-Doc Researcher
> 
> Department of Atmospheric and Environmental Sciences,
> University at Albany, State University of New York, Albany, NY, 12222
> abrammer at albany.edu <mailto:abrammer at albany.edu>
> ##############################
> 
>> On 27 Jan 2016, at 03:21, Tao Lu <hakufu.asano at gmail.com <mailto:hakufu.asano at gmail.com>> wrote:
>> 
>> I try to show nesting domain and some shape files at the same time.
>> So, I edit the plotgrids_new.ncl as showed below, but I could not get what I wanted.
>> 
>> I found the reason is the command: draw(mp), But I do not know how to handle this problem.
>> Please help me.
>> 
>> 
>> ;   Script display location of model domains
>> ;   Only works for ARW domains
>> ;   Only works for NCL versions 6.2 or later
>> ;   Reads namelist file directly
>> 
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>> 
>> ;----------------------------------------------------------------------------
>> begin
>> ;
>> 
>> ; Check the version of NCL
>>   version = systemfunc("ncl -V")
>>   if(version.lt.6.2) then
>>     print("You need NCL V6.2 or later to run this script. Try running plotgrids.ncl. Stopping now...")
>>     return
>>   end if
>> 
>> ; We generate plots, but what kind do we prefer?
>>   type = "x11"
>> ; type = "pdf"
>> ; type = "ps"
>> ; type = "ncgm"
>>   wks = gsn_open_wks(type,"wps_show_dom")
>> 
>> ; read the following namelist file
>>   filename = "namelist.wps"
>> 
>> ; Set the colors to be used
>>   colors = (/"white","black","White","ForestGreen","DeepSkyBlue","Red","Blue"/)
>>   gsn_define_colormap(wks, colors)  
>> 
>> 
>> ; Set some map information ; line and text information
>>   mpres = True
>>   mpres at mpFillOn = True
>>   mpres at mpFillColors  = (/"background","DeepSkyBlue","ForestGreen","DeepSkyBlue", "transparent"/)
>>   mpres at mpDataBaseVersion           = "Ncarg4_1"
>>   mpres at mpGeophysicalLineColor      = "Black"
>>   mpres at mpGridLineColor             = "Black"
>>   mpres at mpLimbLineColor             = "Black"
>>   mpres at mpNationalLineColor         = "Black"
>>   mpres at mpPerimLineColor            = "Black"
>>   mpres at mpUSStateLineColor          = "Black"
>> ;  mpres at mpOutlineBoundarySets       = "AllBoundaries"
>>   ;mpres at mpGridSpacingF              = 45
>>   mpres at tiMainString                = " WPS Domain Configuration  "
>> 
>>   lnres = True 
>>   lnres at gsLineThicknessF = 2.5
>>   lnres at domLineColors    = (/ "white", "Red" , "Red" , "Blue" /)
>> 
>>   txres = True
>>   txres at txFont = "helvetica-bold"
>>   ;txres at txJust = "BottomLeft"
>>   txres at txJust = "TopLeft"
>>   txres at txPerimOn = False
>>   txres at txFontHeightF = 0.015
>> 
>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>> ; Do not change anything between the ";;;;;" lines
>> 
>>   maxdom = 21
>>   nvar = 19
>>   parent_idn = new (maxdom,integer)
>>   parent_grid_ration = new (maxdom,integer)
>>   i_parent_startn = new (maxdom,integer)
>>   j_parent_startn = new (maxdom,integer)
>>   e_wen = new (maxdom,integer)
>>   e_snn = new (maxdom,integer)
>>   plotvar = new((/maxdom,nvar/),float)
>>   plotvar at _FillValue = -999.0
>> 
>>   plotvar = wrf_wps_read_nml(filename)
>> 
>>   mpres at max_dom = floattointeger(plotvar(0,0))
>>   mpres at dx = plotvar(0,1)
>>   mpres at dy = plotvar(0,2)
>>   if (.not.ismissing(plotvar(0,3))) then
>>     mpres at ref_lat = plotvar(0,3)
>>   else
>>     mpres at ref_lat = 0.0
>>   end if
>>   if (.not.ismissing(plotvar(0,4))) then
>>     mpres at ref_lon = plotvar(0,4)
>>   else
>>     mpres at ref_lon = 0.0
>>   end if
>>   if (.not.ismissing(plotvar(0,5))) then
>>     mpres at ref_x = plotvar(0,5)
>>   end if
>>   if (.not.ismissing(plotvar(0,6))) then
>>     mpres at ref_y = plotvar(0,6)
>>   end if
>>   mpres at truelat1 = plotvar(0,7)
>>   mpres at truelat2 = plotvar(0,8)
>>   mpres at stand_lon = plotvar(0,9)
>>   mproj_int = plotvar(0,10)
>>   mpres at pole_lat = plotvar(0,11)
>>   mpres at pole_lon = plotvar(0,12)
>> 
>>   do i = 0,maxdom-1
>>     parent_idn(i) = floattointeger(plotvar(i,13))
>>     parent_grid_ration(i) = floattointeger(plotvar(i,14))
>>     i_parent_startn(i) = floattointeger(plotvar(i,15))
>>     j_parent_startn(i) = floattointeger(plotvar(i,16))
>>     e_wen(i) = floattointeger(plotvar(i,17))
>>     e_snn(i) = floattointeger(plotvar(i,18))
>>   end do
>> 
>> 
>>   if(mpres at max_dom .gt. 1) then
>>     do i = 1,mpres at max_dom-1
>> 
>>       ;Making sure edge is nested grid is at least 5 grid points from mother domain.
>>       if(i_parent_startn(i) .lt. 5) then
>>         print("Warning: Western edge of grid must be at least 5 grid points from mother domain!")
>>       end if
>>       if(j_parent_startn(i) .lt. 5) then
>>         print("Warning: Southern edge of grid must be at least 5 grid points from mother domain!")
>>       end if
>>       pointwe = (e_wen(i)-1.)/parent_grid_ration(i)
>>       pointsn = (e_snn(i)-1.)/parent_grid_ration(i)
>>       gridwe = e_wen(parent_idn(i)-1)-(pointwe+i_parent_startn(i))
>>       gridsn = e_snn(parent_idn(i)-1)-(pointsn+j_parent_startn(i))
>>       if(gridwe .lt. 5) then
>>         print("Warning: Eastern edge of grid must be at least 5 grid points from mother domain!")
>>       end if
>>       if(gridsn .lt. 5) then
>>         print("Warning: Northern edge of grid must be at least 5 grid points from mother domain!")
>>       end if
>> 
>>       ;Making sure nested grid is fully contained in mother domain.
>>       gridsizewe = (((e_wen(parent_idn(i)-1)-4)-i_parent_startn(i))*parent_grid_ration(i))-(parent_grid_ration(i)-1)
>>       gridsizesn = (((e_snn(parent_idn(i)-1)-4)-j_parent_startn(i))*parent_grid_ration(i))-(parent_grid_ration(i)-1)
>>       if(gridwe .lt. 5) then
>>         print("Warning: Inner nest (domain = " + (i+1) + ") is not fully contained in mother nest (domain = " + parent_idn(i) + ")!")
>>         print("For the current setup of mother domain = " + parent_idn(i) + ", you can only have a nest of size " + gridsizewe + "X" + gridsizesn + ". Stopping Program!")
>>         exit
>>       end if
>>       if(gridsn .lt. 5) then
>>         print("Warning: Inner nest (domain = " + (i+1) + ") is not fully contained in mother nest (domain = " + parent_idn(i) + ")!")
>>         print("For the current setup of mother domain = " + parent_idn(i) + ", you can only have a nest of size " + gridsizewe + "X" + gridsizesn + ". Stopping Program!")
>>         exit
>>       end if
>> 
>>       ;Making sure the nest ends at a mother grid domain point.
>>       pointwetrunc = decimalPlaces(pointwe,0,False)
>>       pointsntrunc = decimalPlaces(pointsn,0,False)
>>       if((pointwe-pointwetrunc) .ne. 0.) then
>>         nest_we_up = (ceil(pointwe)*parent_grid_ration(i))+1
>>         nest_we_dn = (floor(pointwe)*parent_grid_ration(i))+1
>>         print("Nest does not end on mother grid domain point. Try " + nest_we_dn + " or " + nest_we_up + ".")
>>       end if
>>       if((pointsn-pointsntrunc) .ne. 0.) then
>>         nest_sn_up = (ceil(pointsn)*parent_grid_ration(i))+1
>>         nest_sn_dn = (floor(pointsn)*parent_grid_ration(i))+1
>>         print("Nest does not end on mother grid domain point. Try " + nest_sn_dn + " or " + nest_sn_up + ".")
>>       end if
>> 
>>     end do
>>   end if
>> 
>>   mpres at parent_id = parent_idn(0:mpres at max_dom-1)
>>   mpres at parent_grid_ratio = parent_grid_ration(0:mpres at max_dom-1)
>>   mpres at i_parent_start = i_parent_startn(0:mpres at max_dom-1)
>>   mpres at j_parent_start = j_parent_startn(0:mpres at max_dom-1)
>>   mpres at e_we = e_wen(0:mpres at max_dom-1)
>>   mpres at e_sn = e_snn(0:mpres at max_dom-1)
>> 
>>   if(mproj_int .eq. 1) then
>>     mpres at map_proj = "lambert"
>>     mpres at pole_lat = 0.0
>>     mpres at pole_lon = 0.0
>>   else if(mproj_int .eq. 2) then
>>     mpres at map_proj = "mercator"
>>     mpres at pole_lat = 0.0
>>     mpres at pole_lon = 0.0
>>   else if(mproj_int .eq. 3) then
>>     mpres at map_proj = "polar"
>>     mpres at pole_lat = 0.0
>>     mpres at pole_lon = 0.0
>>   else if(mproj_int .eq. 4) then
>>     mpres at map_proj = "lat-lon"
>>   end if
>>   end if
>>   end if
>>   end if
>> 
>> ; Deal with global wrf domains that don't have dx or dy
>> 
>>   if (mpres at dx.lt <mailto:mpres at dx.lt>.1e-10 .and. mpres at dx.lt <mailto:mpres at dx.lt>.1e-10) then
>>     mpres at dx = 360./(mpres at e_we(0) - 1)
>>     mpres at dy = 180./(mpres at e_sn(0) - 1)
>>     mpres at ref_lat = 0.0
>>     mpres at ref_lon = 180.0
>>   end if
>> 
>>   mp = wrf_wps_dom (wks,mpres,lnres,txres)
>>   
>> 
>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>> 
>> ; Now you can add some information to the plot. 
>> ; Below is an example of adding a white dot over the DC location.
>>   ;pmres = True
>>   ;pmres at gsMarkerColor = "White"
>>   ;pmres at gsMarkerIndex = 16
>>   ;pmres at gsMarkerSizeF = 0.01
>>   ;gsn_polymarker(wks,mp,-77.26,38.56,pmres)
>>   
>> ; Test to add a shapefile
>> ; I add the following contents
>>   dir         = "/home/taolu/WRF/CASES/kinugawa/shp/"
>>   shp_kinu = dir + "kinugawa.shp"
>>   shp_tone = dir + "tonegawa.shp"
>> 
>>   preres                  = True
>>   preres at gsLineColor      = "gray25"
>>   preres at gsLineThicknessF = 2.0  
>> 
>>   kinu_id = gsn_add_shapefile_polylines(wks,mp,shp_kinu,preres)
>>   tone_id = gsn_add_shapefile_polylines(wks,mp,shp_tone,preres)
>>   
>>   draw(mp);
>> 
>> 
>>   frame(wks)           ; lets frame the plot - do not delete
>> 
>> end
>> 
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu <mailto:ncl-talk at ucar.edu>
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
> 
> 
> 
> 
> -- 
> ******************************************************
> 盧 涛 (TAO LU)
> 〒112-8551 東京都文京区春日1-13-27
> 中央大学理工学研究科都市環境学専攻
> 河川・水文研究室(山田正教授) 修士課程1年
> 
> TEL: 03-3817-3406;   Phone: 070-2188-7509
> E-mail1: hakufu.asano at gmail.com <mailto:mail%3Amet.yamos at gmail.com>
> E-mail2: lutao at civil.chuo-u.ac.jp <mailto:mail%3Ayamoto at civil.chuo-u.ac.jp>
> *******************************************************

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160128/9c5fce36/attachment.html 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 57648 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160128/9c5fce36/attachment.png 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 50169 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160128/9c5fce36/attachment-0001.png 


More information about the ncl-talk mailing list