[ncl-talk] Polar stereographic projection

David Brown dbrown at ucar.edu
Wed Oct 12 12:30:42 MDT 2016


The current algorithm for deciding whether to put latitude or
longitude along any side of a plot is based on which has a bigger
difference over the length of the side. On the left side of this plot
there is a difference of 50+ in longitude. I can't tell exactly what
the latitude difference is but it is presumably less than that amount.
Obviously it would be good if there was more user-level control of
these tick marks, but that is the way it is done now. Improving map
tickmark control is definitely on the to-do list, and will likely be a
priority after the next release.
 -dave


On Wed, Oct 12, 2016 at 9:49 AM, Rick Brownrigg <brownrig at ucar.edu> wrote:
> Hi,
>
> I don't know for certain what is going on there, but I suspect those labels
> are not inherently wrong given the apparent projection, but perhaps just not
> what you wanted. It would be interesting to know the values of
>
>   opts at mpLeftCornerLatF      = lat2d(0,0)
>   opts at mpLeftCornerLonF      = lon2d(0,0)
>   opts at mpRightCornerLatF     = lat2d(nlat-1,mlon-1)
>   opts at mpRightCornerLonF     = lon2d(nlat-1,mlon-1)
>
>   opts at mpCenterLonF          = f1 at CEN_LON
> ;  opts at mpCenterLatF          = f1 at CEN_LAT
>   opts at mpCenterLatF          = 90.         ; This is necessary to fix
> the wrong value on the WRF file.
>
>
> It may be that for the mpCenterLatF to the pole is the issue.  Also,
> temporarily setting this resource to True may help explain what projection
> is like:
>
>   opts at mpGridAndLimbOn       = False         ; turn on lat/lon lines
>
> Hope that helps some...
> Rick
>
>
> On Wed, Oct 12, 2016 at 1:20 AM, <am at bio.mie-u.ac.jp> wrote:
>>
>> Hello NCL users,
>>
>> I am having a problem with plotting figures on polar stereographic
>> projection.
>>
>> The labels on left, right and top sides of the figures seem somewhat
>> strange. They all show longitude. I expect the labels on left and right
>> show latitude as http://www.ncl.ucar.edu/Applications/wrfpo.shtml. Would
>> you please give me suggestions on my script attached below? If you need
>> more information, please let me know. Thank you very much.
>>
>> $ uname -a
>> Linux aofd165.bio.mie-u.ac.jp 2.6.18-398.el5 #1 SMP Tue Sep 16 20:50:52
>> EDT 2014 x86_64 x86_64 x86_64 GNU/Linux
>>
>> $ cat /etc/redhat-release
>> CentOS release 5.11 (Final)
>>
>> $ ncl -V
>> 6.3.0
>>
>>
>> $ ncl.run.sh
>> 1st Argument (script name) =h_sst.ncl
>> 2nd Argument (runname) =PL1101_00.02
>> 3rd Argument (domain) =d02
>> 4th Argument (nc_file) =wrfout_d02_2011-01-18_00:00:00.nc
>> 5th Argument (indir) =../../../../WRF.Result/
>> 6th Argument (nn) =0
>>
>> Input file: ../../../../WRF.Result/PL1101_00.02/WRF_PL1101_00.02/
>> wrfout_d02_2011-01-18_00:00:00.nc
>>
>> fig file = ../../../../WRF.Plots/PL1101_00.02/PL1101_00.02_h_sst.ncl/
>> PL1101_00.02.d02.h_sst.000
>> warning:PlotManagerSetValues: TickMark annotation cannot be added after
>> NhlCreate
>>
>>
>> ======================
>> ncl.run.sh
>> ======================
>> #!/bin/sh
>>
>> exe=./runncl.sh
>>
>> domain=d02
>>
>> #ncl="h_t2.ncl"
>> #ncl="h_t2.polarstereo.ncl"
>> #ncl="h_ice.ncl"
>> ncl="h_sst.ncl"
>>
>> indir="../../../../WRF.Result/"
>>
>> list_runname=" \
>> PL1101_00.02 \
>> "
>>
>> datetimelist="\
>> 2011-01-18_00:00:00 \
>> "
>>
>> nn=0
>> for runname in $list_runname; do
>>
>>   for datetime in $datetimelist; do
>>
>>     nc_file=wrfout_${domain}_${datetime}.nc
>>
>>     $exe $ncl "$runname" "$domain" "$nc_file" "$indir" "$nn"
>>
>>     nn=$(expr $nn + 1)
>>
>>   done
>>
>> done
>>
>>
>> exit 0
>>
>> ----------------------
>>  End of ncl.run.sh
>> ----------------------
>>
>>
>>
>> ======================
>> runncl.sh
>> ======================
>> #!/bin/bash
>> #
>> # Universal wrapper script for ncl.
>> # Pass arguments from the command line to environment variables
>> #
>> # version 0.1, Thierry Corti, C2SM ETH Zurich
>> #
>>
>> E_BADARGS=65
>>
>> if [ ! -n "$1" ]
>> then
>>   echo "Usage: `basename $0` script.ncl argument1 argument2 etc."
>>   exit $E_BADARGS
>> fi
>>
>> # save number of arguments to environment variable NCL_N_ARG
>> export NCL_N_ARGS=$#
>>
>> # save command line arguments to environment variable NCL_ARG_#
>> for ((index=1; index<=$#; index++))
>> do
>>   eval export NCL_ARG_$index=\$$index
>> done
>>
>> # run ncl
>> ncl -nQ $1
>>
>> ----------------------
>>  End of runncl.sh
>> ----------------------
>>
>>
>>
>> ======================
>> h_sst.ncl
>> ======================
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>> load "./add_lc_labels.ncl"
>>
>>
>> begin
>>   wallClock1 = systemfunc("date") ; retrieve wall clock time
>>
>>   a1 = getenv("NCL_ARG_1")
>>   print("1st Argument (script name) ="+a1)
>>   runname1 = getenv("NCL_ARG_2")
>>   print("2nd Argument (runname) ="+runname1)
>>   domain = getenv("NCL_ARG_3")
>>   print("3rd Argument (domain) ="+domain)
>>   nc_file = getenv("NCL_ARG_4")
>>   print("4th Argument (nc_file) ="+nc_file)
>>   indir = getenv("NCL_ARG_5")
>>   print("5th Argument (indir) ="+indir)
>>   nn = getenv("NCL_ARG_6")
>>   print("6th Argument (nn) ="+nn)
>>
>>   scriptname=systemfunc("basename " + a1 + " .ncl")
>>   scriptname2=a1
>>
>>   outdir= "../../../../WRF.Plots/" + runname1 + "/" \
>>   + runname1 + "_" + scriptname2 + "/"
>>
>>   command="mkdir -p " + outdir
>>   system(command)
>>   print("")
>>
>>   infile1= indir + runname1 + "/" + "WRF_" + runname1 + "/" + nc_file
>>   print("Input file: "+ infile1)
>>   print("")
>>
>> ;
>> ; The WRF ARW input file.
>> ; This needs to have a ".nc" appended, so just do it.
>>   f1 = addfile(infile1,"r")
>>
>> ; Get land mask data
>>   it=0
>>   oro = wrf_user_getvar(f1,"XLAND",it)
>>
>>   lat2d = wrf_user_getvar(f1,"XLAT",it)   ; latitude
>>   lon2d = wrf_user_getvar(f1,"XLONG",it)  ; longitude
>>   nlat = dimsizes(lat2d(:,0))
>>   nlon = dimsizes(lat2d(0,:))
>>
>>
>> ; What times and how many time steps are in the data set?
>>   times = wrf_user_getvar(f1,"times",-1)  ; get all times in the file
>>   ntimes = dimsizes(times)         ; number of times in the file
>>
>> ; Check dimension
>>   it=0
>>   lh  = wrf_user_getvar(f1, "LH",it)
>>   dsizes  = new( (/2/),integer)
>>   dsizes  = dimsizes(lh)
>> ;  print(dsizes)
>>   im=dsizes(0)
>>   jm=dsizes(1)
>>
>>   dx = f1 at DX
>>   dy = f1 at DY
>>
>>   lat1=50.
>>   lat2=83.
>>   lon1=-30.
>>   lon2=32.
>>
>> ; We generate plots, but what kind do we prefer?
>>   type = "png" ; ps x11  pdf  ngcm
>>
>>   res                       =  True
>>   res at gsnDraw      =  False                   ; do not draw the plot
>>   res at gsnFrame              =  False                   ;-- don't advance
>> frame
>>   res at gsnAddCyclic          =  False                   ;-- data are not
>> global, don't add lon cyclic point
>>
>> ;  Figure size
>>   res at gsnMaximize      =  False ;True           ; maximize plot in frame
>>   res at vpXF = 0.1
>>   res at vpYF = 1
>>   res at vpHeightF = 0.8
>>   res at vpWidthF  = 0.6
>>
>>   res at gsnDraw          = False                    ; turn off draw and
>> frame
>>   res at gsnFrame         = False                    ; b/c this is an
>> overlay plot
>>   res at pmTickMarkDisplayMode = "Always"                 ;-- draw nicer
>> tickmarks
>>   res at tiMainString    = runname1   ; add titles
>>   res at pmTickMarkDisplayMode = "Always"     ; turn on automatic tickmarks
>>
>>   res at tiMainString          = runname1
>>   res at gsnAddCyclic          = False
>>
>> ;   Plotting options
>>     opts = res
>>     opts at mpFillOn  = False
>>     opts at mpDataBaseVersion     = "MediumRes"              ;-- choose map
>> database
>>     opts at mpGeophysicalLineColor = "Gray" ;changes the outline line color.
>>     opts at mpGeophysicalLineThicknessF = 1 ;changes the thickness of
>> continental outlines.
>>     opts at gsnLeftString   = ""; sst1 at long_name
>>     opts at gsnRightString = times(it)
>>     opts at gsnAddCyclic = False
>>     opts at sfXArray = lon2d
>>     opts at sfYArray = lat2d
>>     opts at mpDataBaseVersion = "HighRes"
>>
>>   dimll  = dimsizes(lat2d)                ; get size of dimensions
>>   nlat   = dimll(0)
>>   mlon   = dimll(1)
>>
>>   opts at mpProjection          = "Stereographic"
>>   opts at mpLimitMode           = "Corners"
>>   opts at mpLeftCornerLatF      = lat2d(0,0)
>>   opts at mpLeftCornerLonF      = lon2d(0,0)
>>   opts at mpRightCornerLatF     = lat2d(nlat-1,mlon-1)
>>   opts at mpRightCornerLonF     = lon2d(nlat-1,mlon-1)
>>
>>   opts at mpCenterLonF          = f1 at CEN_LON
>> ;  opts at mpCenterLatF          = f1 at CEN_LAT
>>   opts at mpCenterLatF          = 90.         ; This is necessary to fix
>> the wrong value on the WRF file.
>>
>>   opts at cnFillOn = True
>>   opts at mpGridAndLimbOn       = False         ; turn on lat/lon lines
>>   opts at cnLevelSelectionMode = "ManualLevels"     ; set manual contour
>> levels
>>   opts at cnMinLevelValF       =  0.                ; set min contour level
>>   opts at cnMaxLevelValF       =  14.                ; set max contour
>> level
>>   opts at cnLevelSpacingF      =  1.               ; set contour spacing
>>
>>   opts2=res
>>   opts2 at cnFillOn = False
>>   opts2 at sfXArray = lon2d
>>   opts2 at sfYArray = lat2d
>>
>>
>>   nt=0
>>   nstart=0 ; Time loop starts a day after the model initial time
>>   nstep=1
>>   nend=1; ntimes-1
>> ;  do it = nstart, nend, nstep              ; TIME LOOP
>>
>>   it=0
>>
>>     lh = wrf_user_getvar(f1,"LH",it) ;
>>     t2 = wrf_user_getvar(f1,"SST",it) ;
>>     t2=t2-273.15
>>
>>     slp=wrf_user_getvar(f1,"slp",it)
>>
>>     figfile= outdir + runname1 + "." + domain + "." + scriptname + \
>>     "." + sprinti("%0.3i",stringtointeger(nn))
>>
>>     print("fig file = " + figfile )
>>     wks = gsn_open_wks(type, figfile)
>>     gsn_define_colormap(wks,"WhBlGrYeRe")
>>
>>     plot = gsn_csm_contour_map(wks,t2,opts)
>>
>>     plot2 = gsn_csm_contour(wks,slp,opts2)
>>
>>     overlay(plot,plot2)
>>
>>
>>    ;MAKE PLOTS
>>     draw(plot)
>>     frame(wks)
>>
>> ;    nt=nt+1
>> ;  end do
>>
>>   wallClock2 = systemfunc("date") ; retrieve wall clock time
>>   print("Started  at " + wallClock1)
>>   print("Finished at " + wallClock2)
>>
>> end
>>
>> ----------------------
>>  End of h_sst.ncl
>> ----------------------
>>
>> Atsuyoshi
>>
>> ---
>> Atsuyoshi Manda, PhD
>> Mie University, Japan
>> am at bio.mie-u.ac.jp
>>
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>


More information about the ncl-talk mailing list