ploygon shapefile subset

Hi all,
I've been trying to fill ploygons in a shapefile different colors depending
on a range of values in the "gr" column of the attribute table of the
shapefile. Admittedly I am very new to NCL. I've tried using examples from
the website and the new one Mary just posted yesterday for subsetting
polylines. I keep getting the bellow error.
Dimension sizes on right hand side of assignment do not match dimension
sizes of left hand side
fatal:["Execute.c":8128]:Execute: Error occurred at or near line 120 in
file project2.ncl

If I comment out the section causing the error I get it to draw but the
entire shapefile is filled in with blue (the first color).
Here is my code...

 in_grib1 = addfile("batch-mars-web240-20130319140652-TnMhMt.grib","r")
  in_grib2 = addfile("batch-mars-web240-20130319145852-7dOvP2.grib","r")
  in_grib3 = addfile("batch-mars-web240-20130319174030-8SXgmI.grib","r")

function create_map(wks,BBS)
local a, res
begin

  res = True

  res@gsnDraw = False ; don't draw yet
  res@gsnFrame = False ; don't advance frame yet
  res@gsnMaximize = True ; maximize plot in frame

  res@mpProjection = "LambertConformal" ; choose projection
  res@mpLambertParallel1F = 33 ; first parallel
  res@mpLambertParallel2F = 45 ; second parallel
  res@mpLambertMeridianF = -98 ; meridian

  res@mpLimitMode = "Corners" ; corner method of zoom
  res@mpLeftCornerLatF = 30 ; left corner
  res@mpLeftCornerLonF = -125 ; left corner
  res@mpRightCornerLatF = 65 ; right corner
  res@mpRightCornerLonF = -35 ; right corner

  res@pmTickMarkDisplayMode = "Always" ; turn on tickmarks

  res@tiMainString = "BBS data"

;---Create map.
  map = gsn_csm_map(wks,res)

return(map)
end

;*************************************************
; main code
;*************************************************
begin

  wks = gsn_open_wks("pdf","project2")

;---Create the map
  map = create_map(wks,"2009")

;---Attach two sets of polylines to the same map.
  lnres = True
  lnres@gsLineThicknessF = 0.5
  lnres@gsLineOpacityF = 0.5
  lnres@gsLineColor = "gray55"
  poly1 = gsn_add_shapefile_polylines(wks,map,"s_06se12_Project.shp",lnres)

  lnres@gsLineColor = "gray55"
  lnres@gsLineThicknessF = 0.5
  lnres@gsLineOpacityF = 0.5
  poly0 = gsn_add_shapefile_polylines(wks,map,"PROVINCE.shp",lnres)

  pres = True
  pres@gsFillColor = "Transparent"
  pres@gsEdgesOn = True
  pres@gsLineThicknessF = 0.06
  pres@gsLineOpacityF = 0.5

  f = addfile("dom_2009_krig.shp","r")
  gr = f->GRIDCODE
  lon = f->lon
  lat = f->lat
  ; Read data off shapefile
;
  segments = f->segments
  geometry = f->geometry
  segsDims = dimsizes(segments)
  geomDims = dimsizes(geometry)
;
; Read global attributes
;
  geom_segIndex = f@geom_segIndex
  geom_numSegs = f@geom_numSegs
  segs_xyzIndex = f@segs_xyzIndex
  segs_numPnts = f@segs_numPnts
  numFeatures = geomDims(0)

  colors = (/"blue","green","yellow","red"/)
  npoly = sum(geometry(:,geom_numSegs))
  poly = new(npoly,graphic) ; Array to hold polygons
  segNum = 0

    do i=0, numFeatures-1

  if (gr(i).ge.290 .and. gr(i).lt.400) then
pres@gsFillColor = colors(0)
end if

if (gr(i).ge.482 .and. gr(i).lt.552) then
pres@gsFillColor = colors(1)
  end if

if (gr(i).ge.600 .and. gr(i).lt.710) then
pres@gsFillColor = colors(2)
end if

if (gr(i).ge.720 .and. gr(i).lt.857) then
pres@gsFillColor = colors(3)
  end if

startSegment = geometry(i, geom_segIndex)
     numSegments = geometry(i, geom_numSegs)
     do seg=startSegment, startSegment+numSegments-1
        startPT = segments(seg, segs_xyzIndex)
        endPT = startPT + segments(seg, segs_numPnts) - 1
        poly(segNum) = gsn_add_shapefile_polygons(wks, map, "dom_2009_krig.
shp", pres)
        segNum = segNum + 1
     end do
  end do

  poly2 = gsn_add_shapefile_polygons(wks,map,"dom_2009_krig.shp",pres)

;---Drawing the map will also draw the attached polygons.
  draw(map)
  frame(wks)
end

-- 
Jennifer McCabe
PhD Student, Ecology and Environmental Science
204 Roger Clapp Greenhouse, UMaine
Orono, ME 04469
p:307-413-5512
jennifer.mccabe@maine.edu <jennifer.mccabe@maine.edu>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Sun Apr 7 09:25:15 2013