[ncl-talk] Creating Canada Mask...too slow ):

[Michael Notaro]

The script below reads in a Canadian GADM shapefile and
creates a 0.1-degree gridded mask where canada=1 for grid cells within Canada
and canada=0 for grid cells outside of Canada.  It works but probably will take
an estimated 15 hours to run.  Can anyone recommend a way to speed it up?
Thanks, Michael


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

a=addfile("create_gridded_snwd.nc","r")
lat=a->lat
lon=a->lon

b=addfile("gadm36_CAN_shp/gadm36_CAN_0.shp","r") ; num_segments=24519
shape_lon=b->x
shape_lat=b->y
geometry=b->geometry
segments=b->segments
segsDims = dimsizes(segments)
geomDims = dimsizes(geometry)
geom_segIndex = b at geom_segIndex
geom_numSegs  = b at geom_numSegs
segs_xyzIndex = b at segs_xyzIndex
segs_numPnts  = b at segs_numPnts
numFeatures = geomDims(0)

canada=new((/201,301/),integer)
canada=0

segNum = 0
do i=0, numFeatures-1
  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
    lons=shape_lon(startPT:endPT)
    lats=shape_lat(startPT:endPT)
    do ilat=60,200
      do ilon=0,300
        inout=gc_inout(lat(ilat),lon(ilon),lats,lons)
        if (inout) then
          canada(ilat,ilon)=1
        end if
      end do
    end do
    segNum = segNum + 1
    delete(lons)
    delete(lats)
  end do
end do

canada!0="lat"
canada!1="lon"
canada&lat=lat
canada&lon=lon

system("rm create_canada_mask.nc")
out=addfile("create_canada_mask.nc","c")
out->canada=canada
out->lat=lat
out->lon=lon

end





Michael Notaro
Associate Director
Nelson Institute Center for Climatic Research
University of Wisconsin-Madison
Phone: (608) 261-1503
Email: mnotaro at wisc.edu
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20190605/9a710023/attachment.html>