Re: shapefile mask

Hi Ciara,

I haven't been following the rest of the thread, but a quick glance over
your code suggests you're falling prey to the need to delete(wks) at the
end of each loop. Insert that exact line 'delete(wks)' after the calls to
draw and frame and before 'end do' and you should see the images pop out.
Let the list know if that doesn't work...

Kyle

----------------------------------------
Kyle S. Griffin
Department of Atmospheric and Oceanic Sciences
University of Wisconsin - Madison
Room 1421
1225 W Dayton St, Madison, WI 53706
Email: ksgriffin2@wisc.edu

On Fri, Mar 14, 2014 at 7:10 AM, Ciara O'Hara <cohara1985@hotmail.com>wrote:

> Hi everyone,
>
> I have tried to put into practise what Rick said below, and I am not
> getting any errors, but I'm also getting any graphs. The script just keeps
> running as if it is doing something, but doesn't stop or output anything.
> If I remove the lines corresponding to the mask it works again. Can anyone
> see any problem with the script below? This is the text that appears on the
> command line when I run the script:
>
> (0) ==================================================
> (0) Shapefile: /home/cohara/IRL_adm/IRL_adm0.shp
> (0) Areas of interest: the whole shapefile
> (0) min_lat_chk: 51.4196
> (0) max_lat_chk: 55.4504
> (0) min_lon_chk: -10.6631
> (0) max_lon_chk: -5.99361
> (0) min_lat_data: 50
> (0) max_lat_data: 56.99
> (0) min_lon_data: -11
> (0) max_lon_data: -5
> (0) 47034 data values originally
> (0) Will keep data values inside given shapefile areas
>
> This is the script:
>
> --------------------------------------------------------------------------------------------------------------------------
> 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 "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
> load "./scripts/rainfall/for_map/shapefile_mask_data.ncl"
>
> begin
> ;8
> n=0
> i=1981
> print("Hi")
> y=systemfunc("ls
> /home/cohara/RainfallData/for_map/timsum_*.nc")
> ;12
> x=addfiles(y,"r")
> do while(i.le.2013)
> data=x[n]->var61(0,0,:,:)
>
> shpfile="/home/cohara/IRL_adm/IRL_adm0.shp"
> ;16
> opt=True
> opt@debug=True
> ;opt@shape_names=(/"Ireland"/)
>
> data_mask=shapefile_mask_data(data,shpfile,opt)
> ;20
> plres=True
> plres@gsEdgesOn=True
> plres@gsLineColor="black"
> ;24
> plres@gsEdgeColor="black"
> plres@PanelPlot=True
> resources=True
> resources@gsnMaximize=True
> ;28
> resources@gsnDraw=False
> resources@gsnFrame=False
> resources@tiMainString="Total Annual Rainfall " +
> sprinti("%0.4i",i)
> resources@cnMonoLineColor=False ;Turn of the drawing of
> contour lines in 1 colour ;32
> resources@tiMainFontColor="black"
> resources@tiMainFont="times-roman"
> resources@cnFillOn=True ;Turn on contour line fill
> resources@cnMonoFillPattern=True ;Turn on single fill
> pattern (default) ;36
> resources@cnMonoFillColor=False ;(default)
> resources@cnMonoLineColor=False ;(default)
> resources@cnLinesOn=False ;(default)
> resources@tmYLOn=False ;Removes tick marks from
> left of y axis ;40
> resources@tmYROn=False
> resources@tmXTOn=False
> resources@tmXBOn=False
> resources@cnLineLabelsOn=False ;Removes contour line
> labels ;44
> resources@cnInfoLabelOn=False ;Removes the
> informational label
> resources@pmLabelBarDisplayMode="Always";Turns on the label bar
> resources@lbLabelBarOn=True
> resources@lbBoxLinesOn=False ;Removes lines from
> between the label bar colours ;48
> resources@lbTitleOn=True ;Turns on label bar title
> resources@lbTitleString="Rainfall (mm)"
> resources@lbPerimOn=False ;Removes perimeter from
> label bar
> resources@lbTitleFontColor="black"
> ;52
> resources@lbTitleFont="times-roman"
> resources@lbOrientation="Horizontal"
> resources@lbLabelPosition="Bottom"
> resources@lbTitlePosition="Top"
> ;56
> ;resources@cnExplicitLegendLabelsOn=True
> ;resources@lgLabelStrings
> =(/"4","5","6","7","8","9","10","11","12","13","14"/)
> ;resources@cnLineLabelStrings=resources@lbLabelStrings
> resources@lbTitleDirection="Across"
> ;60
> resources@lbTitleFontHeightF=0.02
> resources@lbLabelFont="times-roman"
> resources@lbLabelFontColor="black"
> ;resources@cnLevelSelectionMode="ManualLevels" ;Set the label
> bar levels manually ;64
> ;resources@cnLevelSpacingF=200
> resources@cnLevelSelectionMode="ExplicitLevels"
> resources@cnLevels
> =(/400,500,600,700,800,900,1000,1100,1200,1300,1400,1500,1600,1700,1800,1900,2000,2100,2200,2300,2400,2500,2600,2700,2800,2900,3000,3100,3200,3300,3400,3500,3600,3700,3800,3900,4000/)
> resources@mpProjection="Mercator" ;Map
> Projection ;68
> resources@gsnAddCyclic=False
> resources@mpLimitMode="LatLon" ;The viewable portion of
> the map determined by latitude and longitude
> resources@mpMinLatF=51
> resources@mpMaxLatF=55.7
> ;72
> resources@mpMinLonF=-11
> resources@mpMaxLonF=-5.2
> resources@mpGridAndLimbOn=True
> resources@mpDataSetName="Earth..4" ;The data set used (most
> up to date) ;76
> resources@mpDataBaseVersion="MediumRes"
> resources@mpDataResolution="Finest" ;Map render resolution
> resources@mpOutlineOn=False
> resources@mpFillOn=False
> ;80
> ;resources@vpWidthF=1
> ;resources@vpHeightF=0.7
> wtype="png"
> ;wtype@wkWidth=2500
> ;84
> ;wtype@wkHeight=2500
>
> ;sfiles=(/"/home/cohara/IRL_adm/IRL_adm0","/home/cohara/GBR_adm/GBR_adm0"/)+".shp"
> xwks=gsn_open_wks(wtype,"Total_" + sprinti("%0.4i",i))
>
> gsn_define_colormap(xwks,"precip4_11lev")
> ;88
> plot=gsn_csm_contour_map(xwks,data,resources)
>
> poly0=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/IRL_adm/IRL_adm1.shp",plres)
>
> poly1=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/GBR_adm/GBR_adm2.shp",plres)
>
> ;poly2=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/IRL_adm/IRL_adm0.shp",plres)
> ;92
> contour_mask=wrf_contour(x,xwks,data_mask,resources)
> plot_mask=wrf_map_overlays(x,xwks,contour_mask,plres,resources)
> draw(plot)
>
> frame(xwks)
> ;96
> n=n+1
> i=i+1
> end do
> end
> ;100
>
> Thanks for any help,
>
> Ciara
>
>
>
>
>
>
> > From: brownrig@ucar.edu
> > Subject: Re: fatal:mask: dimension sizes of parameter 0 and
> parameter 1 do not match
> > To: cohara1985@hotmail.com; ncl-talk@ucar.edu
> > Date: Wed, 5 Mar 2014 08:37:41 -0700
> >
> > Hi Ciara,
> >
> > On line 17:
> >
> > data_mask=mask(data,shpfile,opt)
> >
> > I think you want to use function shapefile_mask_data(). This is an
> > external function that can be obtained by looking at example #14 at
> > the link below; that example is also one of performing masking in the
> > fashion you are trying to do:
> >
> > http://www.ncl.ucar.edu/Applications/shapefiles.shtml
> >
> > I hope that helps...
> >
> >
> > Rick
> >
> > On Wed, 5 Mar 2014 12:40:17 +0000
> > Ciara O'Hara <cohara1985@hotmail.com> wrote:
> > > Hi,
> > >
> > > I am running the script below, and getting the following error
> > >message:
> > >
> > >
> ------------------------------------------------------------------------------------------------------------------------
> > > fatal:mask: dimension sizes of parameter 0 and parameter 1 do not
> > >match
> > > fatal:["Execute.c":8128]:Execute: Error occurred at or near line 17
> > >in file /home/cohara/scripts/rainfall/for_map/contour_map_script.ncl
> > >
> ------------------------------------------------------------------------------------------------------------------------
> > >
> > > I am trying to use a shapefile to mask out the sea, so that my plot
> > >has nice clean edges around the coast. I understand that it is a
> > >latitude and longitude dimension size problem, but I'm not sure how
> > >to fix it. I've run print_Var_Summary on one of my data variables
> > >(they are all the same except for the date attached), here is the
> > >output:
> > >
> > > ------------------------------------------------
> > > ncl 7> printVarSummary(data)
> > >
> > > Variable: data
> > > Type: float
> > > Total Size: 188136 bytes
> > > 47034 values
> > > Number of Dimensions: 2
> > > Dimensions and sizes: [lat | 234] x [lon | 201]
> > > Coordinates:
> > > lat: [ 50..56.99]
> > > lon: [ -11.. -5]
> > > Number Of Attributes: 4
> > > height : 0
> > > time : 19810930
> > > table : 1
> > > _FillValue : -9e+33
> > > ------------------------------------------------
> > >
> > > Any suggestions would be appreciated,
> > >
> > > Best Regards,
> > > Ciara
> > >
> > > ---------------------------------
> > > script:
> > > ---------------------------------
> > >
> > > 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"
> > >
> > > begin
> > > n=0
> > > i=1981
> > > print("Hi")
> > > y=systemfunc("ls
> > >/home/cohara/RainfallData/for_map/timsum_*.nc")
> > > x=addfiles(y,"r")
> > > do while(i.le.2013)
> > > data=x[n]->var61(0,0,:,:)
> > > shpfile="/home/cohara/IRL_adm/IRL_adm0.shp"
> > > opt=True
> > > opt@return_mask=True
> > > data_mask=mask(data,shpfile,opt)
> > > data_land_mask=where(data_mask.eq.1,data,data@_FillValue)
> > > copy_VarMeta(data,data_land_mask)
> > > plres=True
> > > plres@gsEdgesOn=True
> > > plres@gsLineColor="black"
> > > plres@gsEdgeColor="black"
> > > resources=True
> > > resources@gsnMaximize=True
> > > resources@gsnDraw=False
> > > resources@gsnFrame=False
> > > resources@tiMainString="Total Annual Rainfall " +
> > >sprinti("%0.4i",i)
> > > resources@cnMonoLineColor=False ;Turn of the drawing
> > >of contour lines in 1 colour
> > > resources@tiMainFontColor="black"
> > > resources@tiMainFont="times-roman"
> > > resources@cnFillOn=True ;Turn on contour line
> > >fill
> > > resources@cnMonoFillPattern=True ;Turn on single fill
> > >pattern (default)
> > > resources@cnMonoFillColor=False ;(default)
> > > resources@cnMonoLineColor=False ;(default)
> > > resources@cnLinesOn=False ;(default)
> > > resources@tmYLOn=False ;Removes tick marks
> > >from left of y axis
> > > resources@tmYROn=False
> > > resources@tmXTOn=False
> > > resources@tmXBOn=False
> > > resources@cnLineLabelsOn=False ;Removes contour line
> > >labels
> > > resources@cnInfoLabelOn=False ;Removes the
> > >informational label
> > > resources@pmLabelBarDisplayMode="Always";Turns on the label
> > >bar
> > > resources@lbLabelBarOn=True
> > > resources@lbBoxLinesOn=False ;Removes lines from
> > >between the label bar colours
> > > resources@lbTitleOn=True ;Turns on label bar
> > >title
> > > resources@lbTitleString="Rainfall (mm)"
> > > resources@lbPerimOn=False ;Removes perimeter
> > >from label bar
> > > resources@lbTitleFontColor="black"
> > > resources@lbTitleFont="times-roman"
> > > resources@lbOrientation="Horizontal"
> > > resources@lbLabelPosition="Bottom"
> > > resources@lbTitlePosition="Top"
> > > ;resources@cnExplicitLegendLabelsOn=True
> > > ;resources@lgLabelStrings
> =(/"4","5","6","7","8","9","10","11","12","13","14"/)
> > > ;resources@cnLineLabelStrings=resources@lbLabelStrings
> > > resources@lbTitleDirection="Across"
> > > resources@lbTitleFontHeightF=0.02
> > > resources@lbLabelFont="times-roman"
> > > resources@lbLabelFontColor="black"
> > > ;resources@cnLevelSelectionMode="ManualLevels" ;Set the
> > >label bar levels manually
> > > ;resources@cnLevelSpacingF=0.2
> > > resources@cnLevelSelectionMode="ExplicitLevels"
> > > resources@cnLevels
> =(/0,200,400,600,800,1000,1200,1400,1600,1800,2000,2200,2400,2600,2800,3000,3200,3400,3600,3800,4000/)
> > > resources@mpProjection="Mercator" ;Map Projection
> > > resources@gsnAddCyclic=False
> > > resources@mpLimitMode="LatLon" ;The viewable portion
> > >of the map determined by latitude and longitude
> > > resources@mpMinLatF=51
> > > resources@mpMaxLatF=55.7
> > > resources@mpMinLonF=-11
> > > resources@mpMaxLonF=-5.2
> > > resources@mpGridAndLimbOn=True
> > > resources@mpDataSetName="Earth..4" ;The data set used
> > >(most up to date)
> > > resources@mpDataBaseVersion="MediumRes"
> > > resources@mpDataResolution="Finest" ;Map render
> > >resolution
> > > resources@mpOutlineOn=False
> > > resources@mpFillOn=False
> > > ;resources@vpWidthF=1
> > > ;resources@vpHeightF=0.7
> > > wtype="png"
> > > ;wtype@wkWidth=2500
> > > ;wtype@wkHeight=2500
> > >
> ;sfiles=(/"/home/cohara/IRL_adm/IRL_adm0","/home/cohara/GBR_adm/GBR_adm0"/)+".shp"
> > > xwks=gsn_open_wks(wtype,"Total_" + sprinti("%0.4i",i))
> > > gsn_define_colormap(xwks,"precip4_11lev")
> > > plot=gsn_csm_contour_map(xwks,data_land_mask,resources)
> > >
> poly0=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/IRL_adm/IRL_adm1.shp",plres)
> > >
> poly1=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/GBR_adm/GBR_adm2.shp",plres)
> > >
> poly2=gsn_add_shapefile_polylines(xwks,plot,"/home/cohara/IRL_adm/IRL_adm0.shp",plres)
> > > draw(plot)
> > > frame(xwks)
> > > n=n+1
> > > i=i+1
> > > end do
> > > end
> > >
> > >
>
> _______________________________________________
> 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 Fri Mar 14 07:05:16 2014