[ncl-talk] extra profile in plot with shapefile, how can i get rid of them?
Mary Haley
haley at ucar.edu
Mon Oct 26 09:48:38 MDT 2015
I think these are map outlines that you see. Try setting:
res at mpOutlineOn = False
--Mary
On Mon, Oct 26, 2015 at 9:21 AM, dyjbean at 163.com <dyjbean at 163.com> wrote:
> hi,
> i made a plot with shapefile , but there exist some extra lines in
> figure, as follows:
>
>
> there some extra lines on figure above. how can i get rid of them?
>
>
> the following is my 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 "./shapefile_mask_data.ncl"
>
> begin
> ;forgname="../../globecover_hubei_liuhuan_usgs_ronghe_snow-ice1.txt"
> forgname="beijing_globecover_liuhuan30_usgs_changed.txt"
> ;fconname="../../globecover_huabei_liuhuan_usgs_new_allocate1.txt"
> fconname="beijing_globecover_liuhuan30_usgs_changed_lir-mir-hir.txt"
>
> onlat=5403
> omlon=5403
>
> nnlat=1801
> nmlon=1801
>
> odat=new((/onlat,omlon/),"integer","No_FillValue")
> ndat=new((/nnlat,nmlon/),"integer","No_FillValue")
>
> odat=readAsciiTable(forgname,omlon,"integer",5)
> ndat=asciiread(fconname,(/nnlat,nmlon/),"integer")
> delete(odat at _FillValue)
> delete(ndat at _FillValue)
>
> ;; confused lc
> olat=new((/onlat/),"float","No_FillValue")
> olon=new((/omlon/),"float","No_FillValue")
>
> do i=0,onlat-1
> olat(i)=45.00413067-i*0.0027778
> end do
>
>
> do j=0,omlon-1
> olon(j)=109.99606533+j*0.0027778
> end do
>
> olat!0="lat"
> olat at units="degrees_north"
> olat at long_name="latitude"
>
> olon!0="lon"
> olon at units="degrees_east"
> olon at long_name="longitude"
>
> odat!0="lat"
> odat!1="lon"
> odat&lat=olat
> odat&lon=olon
> odat at units=""
> odat at long_name="land cover confused with globaland30 and globercover"
>
> ;; new classified lc
> nlat=new((/nnlat/),"float","No_FillValue")
> nlon=new((/nmlon/),"float","No_FillValue")
>
> do i=0,nnlat-1
> nlat(i)=45.00413067-i*0.0083334
> end do
>
>
> do j=0,nmlon-1
> nlon(j)=109.99606533+j*0.0083334
> end do
>
> nlat!0="lat"
> nlat at units="degrees_north"
> nlat at long_name="latitude"
>
> nlon!0="lon"
> nlon at units="degrees_east"
> nlon at long_name="longitude"
>
> ndat!0="lat"
> ndat!1="lon"
> ndat&lat=nlat
> ndat&lon=nlon
> ndat at units=""
> ndat at long_name="land cover new categories after confused"
>
> ;
> ; If you already have the mask NetCDF file, set this to False.
> ; Creating the mask can be slow!
> ;
> CREATE_MASK = True
>
> ;---Whether to draw lat/lon points and shapefile outlines
> ADD_LATLON_POINTS = False
> ADD_SHAPEFILE_OUTLINES = True
>
> ;---Name of shapefile containing USA outlines
> shp_fname = "jingjinji.shp"
>
> ;---Name of file to write mask to or to read mask from.
> mask_oname = "jingjinji_mask_300m.nc"
> mask_nname = "jingjinji_mask_1km.nc"
>
>
> ;---Rough area we are interested in. Everything outside this will be masked.
> minlon = 113.
> maxlon = 120.
> minlat = 36.0
> maxlat = 43.0
>
> ;; odat
> if(CREATE_MASK) then
> print("Creating the mask file...")
>
> ;---Create a new mask using a shapefile of USA
> odims = dimsizes(odat)
> opto = True
> opto at return_mask = True
> ; opto at debug = True
> opto at minlon = minlon ; Makes the shapefile masking
> opto at maxlon = maxlon ; go faster.
> opto at minlat = minlat
> opto at maxlat = maxlat
> odat_mask = shapefile_mask_data(odat,shp_fname,opto)
>
> ;---Write new mask to file
> system("rm -f " + mask_oname)
> oout = addfile(mask_oname,"c")
> oout->odat_mask = odat_mask
> else
> print("Reading mask off file.")
>
> ;---Read the new mask from the NetCDF file
> omask = addfile(mask_oname,"r")
> odat_mask = omask->odat_mask
> end if
>
> ;; ndat
> if(CREATE_MASK) then
> print("Creating the mask file...")
>
> ;---Create a new mask using a shapefile of USA
> ndims = dimsizes(ndat)
> optn = True
> optn at return_mask = True
> ; optn at debug = True
> optn at minlon = minlon ; Makes the shapefile masking
> optn at maxlon = maxlon ; go faster.
> optn at minlat = minlat
> optn at maxlat = maxlat
> ndat_mask = shapefile_mask_data(ndat,shp_fname,optn)
>
> ;---Write new mask to file
> system("rm -f " + mask_nname)
> nout = addfile(mask_nname,"c")
> nout->ndat_mask = ndat_mask
> else
> print("Reading mask off file.")
>
> ;---Read the new mask from the NetCDF file
> nmask = addfile(mask_nname,"r")
> ndat_mask = nmask->ndat_mask
> end if
>
> ;---Create masked data array
> omask = where(odat_mask.eq.1,odat,odat at _FillValue)
> copy_VarMeta(odat,omask)
>
> ;---Create masked data array
> nmask = where(ndat_mask.eq.1,ndat,ndat at _FillValue)
> copy_VarMeta(ndat,nmask)
>
>
> print(odat(onlat-1,omlon-1))
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ;; create plot
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> plot=new(2,"graphic")
> info=(/ \
> "1 urban_and_built-in" ,\
> "2 dryland_cropland" ,\
> "3 irrigated_cropland" ,\
> "4 mixed_dry_cropland" ,\
> "5 cropland_grassland" ,\
> "6 cropland_woodland" ,\
> "7 grassland" ,\
> "8 shrubland" ,\
> "9 mixed_shrub_grassland",\
> "10 savanna" ,\
> "11 deciduous_broadleaf" ,\
> "12 deciduous_needleleaf" ,\
> "13 evergreen_broadleaf" ,\
> "14 evergreen_needleleaf" ,\
> "15 mixed_forest" ,\
> "16 water" ,\
> "17 herbaceous_wetland" ,\
> "18 wooded_wetland" ,\
> "19 barren_or_sparse" ,\
> "31 low_intensity_res" ,\
> "32 mid_intensity_res" ,\
> "33 high_intensity_res"/)
>
> ninfo=dimsizes(info)
>
> ninfo = dimsizes(info) ; # of classifications
>
>
> ;************************************************
> ; create plot
> ;************************************************
> pltType = "png" ; ps, pdf, png, x11, eps
>
>
>
> colorscheme =(/"red","SpringGreen2","SpringGreen3","springgreen",\
> "limegreen","darkolivegreen3","green","lightsalmon","green3",\
>
> "LightSeaGreen","palegreen","olivedrab3","chartreuse2","chartreuse3", \
> "darkgreen","blue","hotpink2","hotpink","light cyan", \
> "khaki1","violetred1","purple3"/)
>
> ncolors = dimsizes(colorscheme)
>
> if (ninfo.ne.ncolors) then ; make sure # of colors match categories (classes)
> print("size mismatch: ninfo="+ninfo+" ncolors="+ncolors)
> exit
> end if
>
> wks = gsn_open_wks(pltType,"plot_lc_jingjinji_new")
> plot=new(2,graphic)
> res = True ; plot mods desired
> res at gsnDraw = False
> res at gsnFrame = False
> ;res at gsnMaximize = True ; ps, pdf
> res at gsnAddCyclic = False
>
> ;res at mpMinLatF = min(odat&lat)
> ;res at mpMaxLatF = max(odat&lat)
> ;res at mpMinLonF = min(odat&lon)
> ;res at mpMaxLonF = max(odat&lon)
>
> res at cnFillOn = True ; color Fill
> res at cnFillMode = "RasterFill"
> res at cnLinesOn = False ; Turn off contour lines
> res at cnLineLabelsOn = False
> res at cnLevelSelectionMode = "ExplicitLevels" ; set explict contour levels
> res at cnLevels
> = tobyte((/2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,31,32,33/)) ; tobyte(ispan(1,ninfo-1,1)+1)
> res at cnFillPalette
> = colorscheme ; distinct colors for categories
> res at gsnSpreadColors = False ; use each color sequentially
> res at cnInfoLabelOn = False
>
> ;res at lbLabelPosition = "Center" ; label position
> ;res at lbLabelAlignment = "BoxCenters" ; label orientation
> ;res at lbLabelStrings
> = (/"1","2","3","4","5","6","7","8","9","10","11","12","13","14","15","16",\
> ;"17","18","19","31","32","33"/)
> ;res at lbLabelStrings = ispan(0,ninfo,1)
> ;res at lbLabelStride = 1
> ;res at lbLabelAutoStride
> = False ; in V6.1.0 and up, this is defaulted
>
> ; to True, and overrides lbLabelStride
>
> ; default. Set this to ensure that
>
> ; your label scheme is preserved.
> ;res at pmLabelBarHeightF = 0.075
> ;res at pmLabelBarWidthF = 0.60 ; default is 0.6
> ;res at pmLabelBarOrthogonalPosF = 0.08 ; move up smidge
>
> res at mpFillOn = False
>
> ;res at mpGridAndLimbOn = True
> ;res at mpGridLatSpacingF = 5.0
> ;res at mpGridLonSpacingF = 5.0
> ;res at mpGridLineDashPattern = "2"
>
> res at gsnAddCyclic = False ; regional data
> ;res at trGridType = "TriangularMesh"
> res at tmXBLabelFontHeightF = 0.01 ; Make lon & lat text smaller
> res at tmYLLabelFontHeightF = res at tmXBLabelFontHeightF
>
> res at gsnLeftString = ""
> res at gsnRightString = ""
>
> res1 = True
> res1 = res
> res1 at mpMinLatF = minlat
> res1 at mpMaxLatF = maxlat
> res1 at mpMinLonF = minlon
> res1 at mpMaxLonF = maxlon
>
> res1 at tmXBMode = "Explicit"
> res1 at tmXBValues = (/113.0,114.0,115.0,116.0,117.0,118.0,119.0,120.0/)
> res1 at tmXBLabels
> = (/"113.0~S~o~N~E","114.0~S~o~N~E","115.0~S~o~N~E","116.0~S~o~N~E","117.0~S~o~N~E",\
> "118.0~S~o~N~E","119.0~S~o~N~E","120.~S~o~N~E"/)
>
>
> res1 at tmYLMode = "Explicit"
> res1 at tmYLValues = (/36.0,37.0,38.0,39.0,40.0,41.0,42.0,43.0/)
> res1 at tmYLLabels
> = (/"36.0~S~o~N~N","37.0~S~o~N~N","38.0~S~o~N~N","39.0~S~o~N~N","40.0~S~o~N~N",\
> "41.0~S~o~N~N","42.0~S~o~N~N","43.0~S~o~N~N"/)
>
> res1 at tiMainString = "confused LC"
> res1 at lbLabelBarOn = False
> plot(0) = gsn_csm_contour_map_ce(wks, omask, res1) ; create plot
>
> res2 = True
> res2=res
> res2 at mpMinLatF = minlat
> res2 at mpMaxLatF = maxlat
> res2 at mpMinLonF = minlon
> res2 at mpMaxLonF = maxlon
>
> res2 at tmXBMode = "Explicit"
> res2 at tmXBValues = (/113.0,114.0,115.0,116.0,117.0,118.0,119.0,120.0/)
> res2 at tmXBLabels
> = (/"113.0~S~o~N~E","114.0~S~o~N~E","115.0~S~o~N~E","116.0~S~o~N~E","117.0~S~o~N~E",\
> "118.0~S~o~N~E","119.0~S~o~N~E","120.~S~o~N~E"/)
>
> res2 at tmYLMode = "Explicit"
> res2 at tmYLValues = (/36.0,37.0,38.0,39.0,40.0,41.0,42.0,43.0/)
> res2 at tmYLLabels
> = (/"36.0~S~o~N~N","37.0~S~o~N~N","38.0~S~o~N~N","39.0~S~o~N~N","40.0~S~o~N~N",\
> "41.0~S~o~N~N","42.0~S~o~N~N","43.0~S~o~N~N"/)
>
> res2 at lbLabelBarOn = False
> res2 at tiMainString = "reclassified LC"
> plot(1) = gsn_csm_contour_map_ce(wks,nmask,res2)
>
> pres = True
> pres at gsnFrame = False
> pres at gsnPanelLabelBar = True
> pres at lbLabelFontHeightF = 0.01
> ;pres at gsnPanelBottom = 0.015
> pres at lbLabelPosition = "Center" ; label position
> pres at lbLabelAlignment = "BoxCenters" ; label orientation
> pres at lbLabelStrings
> = (/"1","2","3","4","5","6","7","8","9","10","11","12","13","14","15","16",\
> "17","18","19","31","32","33"/)
> ;pres at lbLabelStrings = ispan(0,ninfo,1)
> pres at lbLabelStride = 1
> pres at lbLabelAutoStride = False
>
> gsn_panel(wks,plot,(/1,2/),pres)
>
> rtxt = True
> rtxt at txJust = "TopLeft"
> rtxt at txFontHeightF = 0.0125
>
> ; Add text: rows x columns of text (arbitrary)
> ; Usually must play with xx, yy and txFontHeightF
>
> nrow = 5 ; # rows
> ncol = 6 ; # columns
>
> n = -1
> xx = 0.015 ; arbitrary
> do nc=0,ncol-1
> yy = 0.20
> do nr=0,nrow-1
> n = n+1
> if (n.le.(ninfo-1)) then
> gsn_text_ndc (wks,info(n),xx,yy,rtxt)
> yy = yy - 2*rtxt at txFontHeightF
> end if
> end do
> xx = xx + 0.20
> end do
>
> ;draw(plot)
> frame(wks)
>
> end
> ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
>
> what kind of properties should i add in the script?
> or something is wrong with the shapefile?
>
>
> thanks for your advices
>
>
>
> ------------------------------
> dyjbean at 163.com
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20151026/7e4d1cee/attachment-0001.html
-------------- next part --------------
A non-text attachment was scrubbed...
Name: bg.jpg
Type: image/jpeg
Size: 20365 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20151026/7e4d1cee/attachment-0001.jpg
-------------- next part --------------
A non-text attachment was scrubbed...
Name: InsertPic_.png
Type: image/png
Size: 247761 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20151026/7e4d1cee/attachment-0001.png
More information about the ncl-talk
mailing list