[ncl-talk] Polymarker on lon/lat grid does not appear
Orestis Speyer
ospeyer at noa.gr
Fri Oct 28 00:37:17 MDT 2016
Dear all,
I am trying to plot a dot (a station) in a map but the dot does not appear.
I looked into the archive search but could not find a solution. I attach the
relevant code below. I am not sure if the problem is a
map/transform/projection issue, an issue of Draw priority or simply NCL does
not know where to plot the dot (as the main plot coordinates are created
"dynamically"). The frustrating/interesting part is that with overlay
(ovelay countour lines to a basic countour map plot in another script) the
dot appears.
Here is the code (sorry for the length). Thank you for your time.
;================================================;
; horizontal contour plot with station
;================================================;
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/csm/contributed.ncl"
; ================================================;
begin
;=================================================;
; choose wks
;=================================================;
nbeglon = (/-0.20/)
nbeglat = (/-0.2/) ; chooses begin of lon/lat(rotated coordinates of the
data)
nendlon = (/0.55/)
nendlat = (/0.58/)
starthour = 0
endhour = 2 ;for three days
ndata = 72
aa = starthour+1;
bb = endhour;
name = "1month"
PRM = "totDRE_surface" ;SOBS_RAD";"T_2M";
UNIT = " [delta W/m^2]"
SCL1 = -2.6;2.5;.4;0.2 ;-0.5;-1
SCL2 = 2.6;2.5;.4;0.2 ;0.5;1
SCL3 = 0.2;0.05;0.25 ;0.1
IYR = 2013;
MON = "12" ; MONTH DEC OR JAN (01)
filename = ("winter."+IYR+"."+PRM+name+"")
;=================================================;
; choose level
;=================================================;
lvl = 39
;=================================================;
; open file and read in data, I cut this out as It does not matter
;=================================================;
;================================================;
; read in data I cut this out as It does not matter
;================================================;
;vardumbs is the variable I'm going to countour (The basic countour plot to
which I want to add a dot)
dum1 = dim_avg_n(vardumbs,0)
plotvar = (dum1) ; in ppb
aveplotvar = dim_avg_n(plotvar,0) ;
aveplotvar2 = dim_avg_n(aveplotvar,0)
printVarSummary (plotvar)
printMinMax (plotvar, True)
delete(dum1)
lat = f1[0]->lat({nbeglat:nendlat},{nbeglon:nendlon}) ;to retrieve lat
lon = f1[0]->lon({nbeglat:nendlat},{nbeglon:nendlon}) ;to retrieve lon
nlat = dimsizes(lat) ;this gives nlat=330
nlon = dimsizes(lon) ;this gives nlon=384
;================================================;
; Create wks according to definition
;================================================;
wks = gsn_open_wks ("png",filename)
;================================================;
; Define Colormap
;================================================;
cmap = (/"white","black","white","lemonchiffon1",\
"lightblue1", "lightblue2","deepskyblue1","deepskyblue3","orchid3",\
"darkorchid","blueviolet","darkorchid4","deeppink2"/)
gsn_define_colormap(wks,"temp_19lev") ;meteoro
;===================================;
; Manage ressources and title
;===================================;
title = ""+PRM+" "+UNIT+""
tdate = "20.12.2013 - 21.1.2014 "
;================================================;
; Set some resources that will apply to the base
; contour/map plot
;================================================;
res = True
res at tiMainString = " "
res at cnLevelSelectionMode = "ManualLevels" ;
set manual contour levels
res at lbLabelAngleF = 45
; angle labels
res at cnMinLevelValF = SCL1
res at cnMaxLevelValF = SCL2
res at cnLevelSpacingF = SCL3
res at cnLabelBarEndStyle = "ExcludeOuterBoxes"
res at cnFillPalette = "temp_19lev" ;meteor
res at gsnLeftString = " "
res at gsnRightString = " "
;===============================================;
; standard stuff
;===============================================;
res at gsnDraw = False ; Do not draw plot now
(helps adding stuff to the plot or producing labels etc)
res at gsnFrame = False ; Do not advance frame
res at cnFillOn = True ; color Fill
res at cnFillMode = "AreaFill" ; raster mode
res at cnRasterSmoothingOn = True
res at cnRasterMinCellSizeF = 0.0005
res at cnLinesOn = True ; Turn off contour lines
res at cnLineLabelsOn = False ; Turn off label of contour
lines
res at cnMaxLevelCount = 100
res at cnInfoLabelOn = False ; do not plot info label
res at gsnAddCyclic = False ; data is not cyclic
res at mpDataBaseVersion = "MediumRes" ; map resolution
res at mpProjection = "CylindricalEquidistant"
res at mpOutlineBoundarySets = "National"
res at mpFillOn = True
res at mpPerimOn = False
res at mpMaskAreaSpecifiers = (/"land"/)
res at mpFillDrawOrder = "PreDraw"
res at mpGridAndLimbOn = True ; en- or disable lon lat
lines
res at mpGeophysicalLineThicknessF = 2.5
res at pmTickMarkDisplayMode = "Always"
res at tmXTMinorOn = False
res at tmYRMinorOn = False
;=================================================;
; map projection
;=================================================;
res at tfDoNDCOverlay = True
res at mpDataBaseVersion = "HighRes" ; map resolution
res at mpProjection = "CylindricalEquidistant" ; map projection
res at mpOutlineBoundarySets = "National"
res at mpCenterLonF = lon({rlon|0},{rlat|0})
res at mpCenterLatF = lat({rlon|0},{rlat|0})
res at mpLimitMode = "Corners"
res at mpLeftCornerLatF = lat(0,0) ; latitude left lower
corner from model output
res at mpLeftCornerLonF = lon(0,0) ; longitude left lower
corner from model output
res at mpRightCornerLatF = lat(nlat(0)-1,nlat(1)-1) ; latitude right upper
corner from model output
res at mpRightCornerLonF = lon(nlon(0)-1,nlat(1)-1) ; longitutde right
upper corner from model output
res at sfXCStartV = nbeglon ; minimal length
(rotated): embedd data into map
res at sfXCEndV = nendlon ; maximal length
(rotated): embedd data into map
res at sfYCStartV = nbeglat ; minimal length
(rotated): embedd data into map
res at sfYCEndV = nendlat ; maximal length
(rotated): embedd data into map
;=================================================;
; Scaling
;=================================================;
res at vpWidthF = 0.6 ; change the aspect ratio
res at vpHeightF = 0.6
res at vpXF = .1 ; location of where plot starts
;res at vpYF = 1
bot_plot = gsn_csm_contour_map(wks,plotvar(:,:),res)
;-----------add polymarker-----------------;
light_gray = NhlNewColor(wks,0.85,0.85,0.85) ; Add light gray
mkres =True
thissio_lat = 0.2732
thissio_lon =0.1720
mkres at gsMarkerIndex = 17 ; Filled circle
mkres at tfPolyDrawOrder = "PreDraw"
mkres at gsMarkerSizeF = 0.03
Station = gsn_add_polymarker(wks,bot_plot,thissio_lon,thissio_lat,mkres)
draw (bot_plot)
;*****************************************************
; Manually create and attach titles
;*****************************************************
;
; Attach some titles at the top.
;
res_text = True
res_text at txFontHeightF = 0.03 ; change font size
txid_top = gsn_create_text(wks, title, res_text)
amres = True
amres at amJust = "BottomCenter"
amres at amParallelPosF = 0.0 ; This is the center of the plot.
amres at amOrthogonalPosF = -0.72 ; This is above the top edge of the plot.
annoid_top = gsn_add_annotation(bot_plot, txid_top, amres)
res_text at txFontHeightF = 0.02 ; change font size
txid_mid = gsn_create_text(wks, tdate,res_text)
amres at amOrthogonalPosF = -0.62 ; This is just below the previous title.
annoid_mid = gsn_add_annotation(bot_plot, txid_mid, amres)
pres = True
maximize_output(wks,pres)
frame (wks)
delete(plotvar)
end
Once again thank you,
Orestis Speyer,
Research Fellow, National Observatory of Athens
Institute for Environmental Research and Sustainable Development
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20161028/dd6dfe94/attachment.html
More information about the ncl-talk
mailing list