station plot

I read in a dataset with the following format:

lat1, lon1, treecover1
lat2, lon2, treecover2
...

I tried following station_1.ncl's example for plotting
station data as contour/shaded. But the map always
turns out empty. Any idea why my data doesn't appear?

Mike

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/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin

a = asciiread("fpc.out",(/10320,12/),"float")
lat = a(:,0)
lon = a(:,1)
fpc = a(:,2:11)
tree = dim_sum(a(:,2:8))
grass = dim_sum(a(:,9:11))

; 12.73N to 58.16N, 146.94W to 42.46W

   tlat1 = 60.0
   tlat2 = 30.0
   clon = -98.5
   clat = 36.3

   latlon2d=new((/2,71,151/),float)
   lat1=fspan(0.,70.,71)
   lon1=fspan(-170.,-20.,151)
   do ilat=0,70
     do ilon=0,150
       latlon2d(0,ilat,ilon)=lat1(ilat)
       latlon2d(1,ilat,ilon)=lon1(ilon)
     end do
   end do
   lat2d=latlon2d(0,:,:)
   lon2d=latlon2d(1,:,:)

   wks = gsn_open_wks("ps","plot_fpc") ; open a ps file
   gsn_define_colormap(wks,"WhViBlGrYeOrRe")
   gsn_reverse_colormap(wks)

   res = True
   res_at_gsnFrame = False ; So we can draw markers
   res_at_gsnMaximize = True
   res_at_gsnSpreadColors = True
   res_at_gsnSpreadColorEnd = 90
   res_at_pmTickMarkDisplayMode = "Always"
   res_at_trGridType = "TriangularMesh" ; The
default if you
                                                         ; have 1D data
   res_at_cnLineLabelPlacementMode = "Constant"
   res_at_cnLineDashSegLenF = 0.3
   res_at_cnLevelSelectionMode = "ManualLevels"
   res_at_cnMinLevelValF = 0.
   res_at_cnMaxLevelValF = 1.
   res_at_cnLevelSpacingF = 0.05
   res_at_cnFillOn = True
   res_at_cnLinesOn = True
   res_at_cnLineLabelsOn = True
   res_at_cnLevelFlags = new(139,"string")
   res_at_cnLevelFlags(:) = "NoLine"
   res_at_cnLevelFlags(0::20) = "LineAndLabel"
; res_at_lbLabelStride = 20
   res_at_lbBoxLinesOn = False
   res_at_tiMainString = "Forest Cover"
   res_at_sfXArray = lon
   res_at_sfYArray = lat
   res_at_mpProjection = "LambertConformal"
   res_at_mpLambertParallel1F = tlat1
   res_at_mpLambertParallel2F = tlat2
   res_at_mpLambertMeridianF = clon
   res_at_mpLimitMode = "Corners"
   res_at_mpLeftCornerLatF = lat2d(10,20)
   res_at_mpLeftCornerLonF = lon2d(10,20)
   res_at_mpRightCornerLatF = lat2d(60,130)
   res_at_mpRightCornerLonF = lon2d(60,130)
   res_at_mpFillOn = False
   res_at_mpOutlineDrawOrder = "PostDraw"
   res_at_mpFillDrawOrder = "PreDraw"
   res_at_mpOutlineBoundarySets = "GeophysicalAndUSStates"
   res_at_mpUSStateLineColor = "Gray10"
   res_at_mpUSStateLineDashPattern = 2
   map = gsn_csm_contour_map(wks,tree,res)

   frame(wks) ; Now advance the frame.

end

_______________________________________________
ncl-talk mailing list
ncl-talk_at_ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Aug 24 2006 - 13:23:05 MDT