[ncl-talk] Conformal projection

[Carlos Peralta]

Hello,

I am following example 3 here:

https://www.ncl.ucar.edu/Applications/Scripts/lcnative_3.ncl

to plot a grb file where I am not sure of what the appropriate
parallels are.

I get the figure below, which seems correct, but not sure
why I see the axes displayed with "W" and "E" in both the vertical
and horizontal.

thanks for any help.


My code is below

;------------------------------------------------
; lcnative_3.ncl
;------------------------------------------------
;
; 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"

begin

;------------------------------------------------
; open file and read in data
;------------------------------------------------
  gfile= "./CW_Tskin.grib"
  f     = addfile (gfile, "r")
  printVarSummary(f)
  ;fili  = "ruc2.bgrb.20020418.i12.f00.grb"
  ;f     = addfile (fili, "r")

  ;vpt       = f->VPTMP_252_HYBL
  vpt       = f->TMP_GDS3_HTGL(g3_x_0|:,g3_y_1|:)
  ;lat2d     = f->gridlat_252
  ;lon2d     = f->gridlon_252
   lat2d=f->g3_lat_0
   lon2d=f->g3_lon_1


;------------------------------------------------
; create plot
;------------------------------------------------
  wks = gsn_open_wks("png","lcnative")      ; send graphics to PNG file

  res                   = True              ; plot mods desired
  res at tiMainString      = "Native Lambert Conformal"
  res at gsnAddCyclic = False
;------------------------------------------------
; lambert conformal projections are limited using
; the corners method rather than the latlon method
; seen for cylindrical equidistant projections
;------------------------------------------------
  res at mpLimitMode       = "Corners"         ; choose range of map
  res at mpLeftCornerLatF  = lat2d at corners(0)
  res at mpLeftCornerLonF  = lon2d at corners(0)
  res at mpRightCornerLatF = lat2d at corners(2)
  res at mpRightCornerLonF = lon2d at corners(2)
;------------------------------------------------
; there are many instances when you may not know the
; appropriate parallels. Here is a method to find out...
; 1) as a first guess, start with the parellels in the
;    center of the grid.
; 2) look at the lines plotted on the map (see below).
;    you know the projection is correct when the red
;    and blue lines form a right angle between themselves
;    and the boundary.
; 3) iterate as necessary to get the closest solution.
;------------------------------------------------
  res at mpProjection          = "LambertConformal"

  res at mpLambertParallel1F = (max(lat2d)+min(lat2d))*0.5
  res at mpLambertParallel2F = (max(lat2d)+min(lat2d))*0.5
  res at mpLambertMeridianF  = (max(lon2d)+min(lon2d))*0.5

;------------------------------------------------
  res at pmTickMarkDisplayMode = "Always"
  res at mpFillOn              =  False           ; turn off map fill
  res at mpOutlineDrawOrder    = "PostDraw"       ; draw continental outline
last
  res at mpOutlineBoundarySets = "GeophysicalAndUSStates" ; state boundaries

;------------------------------------------------
; usually, when data is placed onto a map, it is TRANSFORMED to the
specified
; projection. Since this model is already on a native lambert conformal
grid,
; we want to turn OFF the transformation.
;------------------------------------------------

  res at tfDoNDCOverlay = True
; res at tfDoNDCOverlay = "NDCViewport"          ; NCL V6.5.0 or later


;------------------------------------------------
; set additional resources and create plot
;------------------------------------------------
  res at cnFillOn              = True            ; color plot desired
  res at cnLinesOn             = False           ; turn off contour lines
  res at cnFillPalette         = "gui_default"   ; set color map
  res at gsnFrame              = False           ; don't advance frame yet

  ;plot = gsn_csm_contour_map(wks,vpt(9,:,:),res) ; draw plot
  plot = gsn_csm_contour_map(wks,vpt(:,:),res) ; draw plot

;------------------------------------------------
; add lines that help determine correct projection
;------------------------------------------------
  dims  = dimsizes(vpt)                       ; get dimension of grid
  ;nlat  = dims(1)
  ;nlon  = dims(2)
  nlat  = dims(0)
  nlon  = dims(1)

  gres                  = True                ; polyline mods desired
  gres at gsLineThicknessF = 3.0                 ; line thickness
  gres at gsLineColor      = "Red"               ; line color
  gsn_polyline (wks, plot, lon2d(nlat/2,:), lat2d(nlat/2,:), gres)

  gres at gsLineColor      = "Blue"
  gsn_polyline (wks, plot, lon2d(:,nlon/2), lat2d(:,nlon/2), gres)

  gres at gsLineColor      = "Black"
  gsn_polyline (wks, plot, lon2d(:,0)     , lat2d(:,0), gres)
  gsn_polyline (wks, plot, lon2d(:,nlon-1), lat2d(:,nlon-1), gres)
  gsn_polyline (wks, plot, lon2d(0,:)     , lat2d(0,:), gres)
  gsn_polyline (wks, plot, lon2d(nlat-1,:), lat2d(nlat-1,:), gres)

  frame(wks)

end




[image: image.png]
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20200915/1e799ad4/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 169869 bytes
Desc: not available
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20200915/1e799ad4/attachment.png>