Re: (no subject)

From: Mary Haley <haley_at_nyahnyahspammersnyahnyah>
Date: Tue Aug 16 2011 - 11:54:37 MDT

Hi,

I haven't tried this, but what about making as many copies that you need of "map" after you attach the shapefile polylines, and then just use these copies to overlay the additional data on?

--Mary

On Aug 16, 2011, at 3:45 AM, ugo merlini wrote:

> Hi
>
> I take 50 seconds
>
> Is possibile to create a shapefile map only one time an after use it as background for ploting over an input data for example 2m temperature at time t, t+1 (other loop)? I file for each time.
>
> this way should be more fast
>
> now I'm trying with overlay function so before I'm creating two maps one with the shapefile only the other one the contour only here the code but I'm getting thes warning messages
>
> warning:mpLimitMode is not a valid resource in JAL_contour at this time
> warning:mpMinLonF is not a valid resource in JAL_contour at this time
> warning:mpMaxLonF is not a valid resource in JAL_contour at this time
> warning:mpMinLatF is not a valid resource in JAL_contour at this time
> warning:mpMaxLatF is not a valid resource in JAL_contour at this time
>
> code
>
> ;----------------------------------------------------------------------
> ; This example shows how to read geographic data
> ; from Natural Earth shapefiles
> ; and plot them as polylines and polygons.
> ;----------------------------------------------------------------------
> ; This particular example plots data for Switzerland.
> ;----------------------------------------------------------------------
> ; Download the shapefiles from http://www.naturalearthdata.com/
> ; Unzip to a directory
> ;----------------------------------------------------------------------
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>
> begin
>
> ;--- Open workstation and define colormap.
> wks_type = "png" ; use "newpdf" instead of "pdf" for smaller files
> wks = gsn_open_wks(wks_type,"JAL")
> gsn_define_colormap(wks,"BlAqGrYeOrRe") ; choose colormap
>
> res = True ; plot mods desired
> res@gsnMaximize = True
> ;res@cnFillOn = False ; turn on color fill
> res@gsnTickMarksOn = False
> res@mpLimitMode="LatLon"
> res@mpMinLatF = 35.0
> res@mpMaxLatF = 48.0
> res@mpMinLonF = 6.5
> res@mpMaxLonF = 19.0
> res@gsnDraw = False ; do not draw the plot
> res@gsnFrame = False ; do not flip the page
>
> map = gsn_csm_map(wks,res)
>
>
>
>
>
>
> ;---Shapefiles
> fnames = "/mnt/internetserver/map/shapefile/europa/aaa_full/" + (/"europa"/) + ".shp" ; Files to be open
> linecolors = (/"Black"/) ; color definition for each shapefile
> fillcolors = (/"Transparent"/) ; Fill color definition for each shapefile
> thicks = (/2/) ; Thickness for each shapefile
> lnres = True ; resources for polylines
> plres = True
>
> ;---Loop through files that we want to read geographic information from.
> prims = True
> lines = True
> do n=0,dimsizes(fnames)-1 ; Loop through files that we want to read geographic information from.
> f = addfile(fnames(n),"r") ; Open the shapefile number n.
> segments = f->segments ; Read data off the shapefile
> geometry = f->geometry
> segsDims = dimsizes(segments)
> geomDims = dimsizes(geometry)
> geom_segIndex = f@geom_segIndex ; Read global attributes
> geom_numSegs = f@geom_numSegs
> segs_xyzIndex = f@segs_xyzIndex
> segs_numPnts = f@segs_numPnts
> geometry_type = f@geometry_type
> numFeatures = geomDims(0)
> lon = f->x
> lat = f->y
> if (geometry_type.eq."polygon") then ; Put if statement outside the loop
> plres@gsFillColor = fillcolors(n)
> plres@gsEdgesOn = True ; Draw border around polygons
> plres@gsEdgeColor = linecolors(n)
> plres@gsEdgeThicknessF = thicks(n)
> do i=0, numFeatures-1 ; Section to draw polygons on map.
> 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
> dumstr = unique_string("lines") ; This call adds the polygon.
> map@$dumstr$ = gsn_add_polygon(wks, map , lon(startPT:endPT), lat(startPT:endPT), plres)
> end do
> end do
> else
> lnres@gsLineThicknessF = thicks(n)
> lnres@gsLineColor = linecolors(n)
> do i=0, numFeatures-1 ; Section to draw polylines on map.
> 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
> dumstr = unique_string("primitive") ; This call adds the line segment.
> map@$dumstr$ = gsn_add_polyline(wks, map, lon(startPT:endPT), lat(startPT:endPT), lnres)
> end do
> end do
> end if
> delete(lat) ; Clean up before we read in same variables again.
> delete(lon)
> delete(segments)
> delete(geometry)
> delete(segsDims)
> delete(geomDims)
> end do
>
> ;***************************
> ; read in data
> ;***************************
>
> url = "http://nomad1.ncep.noaa.gov:9090/dods/gfs_master/gfs20110815/"
> filename = url + "gfs_master_00z"
>
> exists = isfilepresent(filename)
> if(.not.exists) then
> print("OPeNDAP isfilepresent test unsuccessful.")
> print("Either file doesn't exist, or NCL does not have OPeNDAP capabilities on this system")
> else
> print("OPeNDAP isfilepresent test successful.")
> a = addfile(filename,"r")
> vars = getfilevarnames(a)
> end if
>
> TEMP2M = a->tmp2m(0,:,:) ; temperature
> TEMP2M = TEMP2M - 273.15
> TEMP2M@units = "(C)"
> TEMP2M@long_name = "Temperatura a 2 metri dal suolo"
>
> res1= True ; plot mods desired
> res1@cnFillOn = True ; turn on color
> res1@cnLinesOn = True ; turn contour lines
> res1@gsnSpreadColors = True ; use entire color map
> res1@lbLabelAutoStride = True ; let NCL determine label spacing
>
> ;res1@gsnDraw = False ; do not draw the plot
> ;res1@gsnFrame = False ; do not flip the page
> ;res1@cnLineDrawOrder = "PostDraw"
> ;res1@gsnSpreadColors = True ; use full range of color map
> ;res1@gsnMaximize = True
>
> plot = gsn_csm_contour(wks,TEMP2M,res)
> overlay(map,plot)
>
> ;---We're done adding primitives, draw everything and advance frame.
>
> draw(map)
> frame(wks)
> end
>
>
> ugo
>
>
>
>
>
>
>
>
>
>
>
>
>
> > Date: Mon, 15 Aug 2011 13:15:22 -0600
> > From: dave.allured@noaa.gov
> > Subject: Re: (no subject)
> > To: ugomerlini@hotmail.com
> > CC: ncl-talk@ucar.edu
> >
> > Drawing shape files in NCL is slow because of nested loops. On my
> > computer I get these speeds to make only one plot:
> >
> > 3.3 million points X11 output 24 seconds
> > 3.3 million points PDF output 36 seconds
> >
> > What speed are you getting?
> >
> > --Dave
> >
> > On 8/15/2011 11:44 AM, ugo merlini wrote:
> > > i made my first script which use a shapefile as source but it is
> > > very slow to produce the output
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Tue Aug 16 11:54:46 2011

This archive was generated by hypermail 2.1.8 : Mon Aug 22 2011 - 08:13:37 MDT