Announcement

Collapse
No announcement yet.
X
  • Filter
  • Time
  • Show
Clear All
new posts

  • geoplot with geoframe rclip - suppressing plotting of clip boundaries?

    Dear Statalist,

    I am using the excellent tools provided in Ben Jann's geoplot and related commands. I have used the geoframe rclip tool to focus on my geographic area of interest. Everything works beautifully except that the graph plots the boundaries of the rectangle selected and clipped by geoframe clip. I would like to draw the graph without these boundaries.

    Is this possible? Or possible but not advisable for some reason?

    Many thanks,
    Raymond

    geoplot and geoframe installed from the geoplot github site
    https://github.com/benjann/geoplot

    geoframe.ado
    *! version 1.3.5 12oct2024 Ben Jann
    geoplot.ado
    *! version 1.3.5 17sep2024 Ben Jann

    StataNow 19.5 MP on Windows 11

    Code:
    geoframe query bbox if dist_sadhanw_KM<10, padding(5)
    mat box_sadhanw_10KM = r(limits)
    
    foreach toClip in state district subdistrict roads {
      frame `toClip': geoframe rclip box_sadhanw_10KM ///
                        , into(`toClip'_sadhanw_10KM)
    }
    
    loc graph_name villages_LT10KM_sadhanw_clipped
    geoplot  ///
      (line district_sadhanw_10KM if inlist(pc11_d_id,`levels_distID') ///
        , lwidth(medium) lcolor(black)) ///
      (line subdistrict_sadhanw_10KM ///
        , lpattern(dot) lwidth(medthin) lcolor(gs8)) ///
      (point village ///
         if inlist(pc11_id,"`sadhanw_pc11_id'") ///
         , msymbol(diamond)) ///
      (point village if dist_sadhanw_KM<10 ///
        & !inlist(pc11_id,"`sadhanw_pc11_id'") ///
        & town==0 ///
        , msize(medsmall) mcolor(navy)) ///
      (point village if dist_sadhanw_KM<10 ///
        & !inlist(pc11_id,"`sadhanw_pc11_id'") ///
        & town==1,     ///
        msize(small) msymbol(triangle)) /// ///
      (label village place_name ///
        if (dist_sadhanw_KM<10 & town==1) | inlist(pc11_id,"`sadhanw_pc11_id'") ///
        , position(12)) ///
      (line roads_sadhanw_10KM if inlist(code,5115), lcolor(gs12)) ///
      (line roads_sadhanw_10KM if inlist(code,5112,5113,5114), lcolor(gs4)) ///
      , ///
      glegend(layout(1 "District boundary" 2 "Subdistrict boundary" 8 "Major road" 7 "Local road") ///
        position(ne) box) ///
      glegend(layout(5 "Town" 4 "Village") position(se) box) ///
        saving("output/fig/gph/`graph_name'.gph", replace) noisily

    graph twoway command returned by "noisily" option:

    Code:
    graph twoway
        (line Y X in 3/330, cmissing(n) lpattern(solid) lcolor(black)
            lwidth(medium))
        (line Y X in 331/658, cmissing(n) lcolor(gs8) lpattern(dot)
            lwidth(medthin))
        (scatter Y X in 659/659, msymbol(diamond))
        (scatter Y X in 660/713, mcolor(navy) msize(medsmall))
        (scatter Y X in 714/714, mlabel(MLAB) msize(0) msymbol(i)
            mlabposition(12))
        (line Y X in 715/1124, cmissing(n) lpattern(solid) lwidth(.15)
            lcolor(gs12))
        (line Y X in 1125/2177, cmissing(n) lpattern(solid) lwidth(.15)
            lcolor(gs4))
        (area Y X in 2178/2183, cmissing(n) nodropbase lalign(center)
            lpattern(solid) lwidth(.15) finten(100) fcolor(none)
            lcolor(gray))
        (line Y X in 2184/2185, cmissing(n) lcolor(black) lpattern(solid)
            lwidth(medium) pstyle(p1))
        (line Y X in 2186/2187, cmissing(n) lcolor(gs8) lpattern(dot)
            lwidth(medthin) pstyle(p2))
        (line Y X in 2188/2189, cmissing(n) lcolor(gs4) lpattern(solid)
            lwidth(.15) pstyle(p7))
        (line Y X in 2190/2191, cmissing(n) lcolor(gs12) lpattern(solid)
            lwidth(.15) pstyle(p6))
        (scatter Y X in 2192/2197, mlabel(MLAB) mlabvposition(MLABPOS)
            mlabcolor(black) mlabsize(vsmall) msize(0) msymbol(i) mlabgap(0))
        (area Y X in 2198/2203, cmissing(n) nodropbase lalign(center)
            lpattern(solid) lwidth(.15) finten(100) fcolor(none)
            lcolor(gray))
        (scatter Y X in 2204/2204, mcolor(navy) msize(medsmall) pstyle(p4))
        (scatter Y X in 2205/2207, mlabel(MLAB) mlabvposition(MLABPOS)
            mlabcolor(black) mlabsize(vsmall) msize(0) msymbol(i) mlabgap(0))
        ,
        legend(off)
        aspectratio(1, units)
        xscale(range(74.52512929562504 74.77266475337508) off)
        yscale(range(30.70949541310006 30.91568831890003) off)
        ylabels(none, labsize(vsmall))
        xlabels(none, labsize(vsmall))
        xtitle("", size(small))
        ytitle("", size(small))
        plotregion(style(none) istyle(none) margin(zero))
        graphregion(style(none) istyle(none) margin(small))
        bgcolor(white)
        saving("output/fig/gph/villages_LT10KM_sadhanw_clipped.gph", replace)
    Outputs

    Click image for larger version

Name:	villages_LT10KM_sadhanw_clipped.png
Views:	1
Size:	153.0 KB
ID:	1780402
    Click image for larger version

Name:	villages_LT10KM_kanech_clipped.png
Views:	1
Size:	278.4 KB
ID:	1780403
    Last edited by Raymond Guiteras; 28 Jul 2025, 14:32.

  • #2
    Please provide a reproducible example.

    Comment


    • #3
      Adapted from Asjad Naqvi 's "Maps in Stata III: geoplot"
      https://medium.com/the-stata-guide/m...t-a764cf42688a

      Code:
      local filelist nuts3 nuts3_shp 
      
      foreach x of local filelist {
       
       copy "https://github.com/asjadnaqvi/The-Stata-Guide/raw/master/GIS/`x'.dta" "`x'.dta", replace
       
      }
      
      geoframe create nuts3 , replace 
      frame change nuts3  
      
      geoframe query bbox if NUTS_ID=="DE300", pad(60) n(50) circle
      mat zoom1 = r(limits)  
      
      frame nuts3: geoframe rclip zoom1, into(nuts3_zoom1) replace
      
      geoplot (line nuts3_zoom1)

      Comment


      • #4
        Here is the output from the noisily option

        Code:
        graph twoway
            (line Y X in 3/338, cmissing(n) lpattern(solid) lwidth(.15) lcolor(gray))
            ,
            legend(off)
            aspectratio(1, units)
            xscale(range(4516441.960746481 4589811.945297712) off)
            yscale(range(3231111.311695597 3304336.517342095) off)
            ylabels(none, labsize(vsmall))
            xlabels(none, labsize(vsmall))
            xtitle("", size(small))
            ytitle("", size(small))
            plotregion(style(none) istyle(none) margin(zero))
            graphregion(style(none) istyle(none) margin(small))
            bgcolor(white)

        Comment


        • #5
          Code:
          geoplot line nuts3_zoom1  /// 
              ||  scatteri 4516442 3231111 4589812 3231111 4589812 3304337 /// 
              4516442 3304337 4516442 3231111 , recast(line) lc(white) lw(thick)

          Click image for larger version

Name:	stata.png
Views:	1
Size:	47.8 KB
ID:	1780419

          Comment


          • #6
            In addition to #5, you can use the -if- qualifier to restrict the range of the x-dimension. Since the y-dimension spans the entire graph region, you can simply manipulate the margins to remove the top and bottom lines. Changes are highlighted.

            Code:
            clear
            cap frame change default
            cap frame drop n*
            local filelist nuts3 nuts3_shp 
            
            foreach x of local filelist {
             
             copy "https://github.com/asjadnaqvi/The-Stata-Guide/raw/master/GIS/`x'.dta" "`x'.dta", replace
             
            }
            
            geoframe create nuts3 , replace 
            frame change nuts3  
            
            geoframe query bbox if NUTS_ID=="DE300", pad(60) n(50) circle
            mat zoom1 = r(limits)
            local xmin = zoom1[1, 1]  
            local xmax = zoom1[1, 2]  
            
            frame nuts3: geoframe rclip zoom1, into(nuts3_zoom1)  replace
            geoplot (line nuts3_zoom1 if inrange(_CX, `xmin', `xmax')) , graphregion(margin( b-4 t-4))



            Click image for larger version

Name:	Graph.png
Views:	1
Size:	42.1 KB
ID:	1780422

            Comment


            • #7
              Dear Raymond Guiteras,
              Even though I have moved all my workflow to geoplot, I do have this standalone package for clipping shapefiles from the pre-geoplot days: https://github.com/asjadnaqvi/stata-clipgeo


              You can use this to create clipped shapefiles and just geoframe them.

              Hope this helps!




              Comment


              • #8
                Thank you all for your helpful suggestions.

                Scott Merryman, if I understand what you're doing correctly, you are covering the border with an opaque white line. Is that correct?

                Andrew Musau, your proposal works perfectly but I can't understand why it works for the x-dimension., I understood _CX to be the centroid of the polygon, not the border, so I would have expected your qualifier (if inrange(_CX, `xmin', `xmax')) to mean "if the x-coordinate of the polygon's centroid is > xmin and < xmax", not "if the x-coordinate of this particular point is > xmin and < xmax". Do you have any insight on why it appears to mean the latter rather than the former?

                Asjad Naqvi, I may be misunderstanding but it appears to me from the examples on your GitHub site as though clippolygon adds the boundaries of the box to polygons, while clippolylines does not add boundaries to line shapes. Is this correct?

                Comment


                • #9
                  You are correct that _CX refers to the x-coordinate of the polygon’s centroid. So yes, it is not the coordinate of each vertex along the polygon's border. Thus, when I write:

                  geoplot (line nuts3_zoom1 if inrange(_CX, `xmin', `xmax'))
                  I am instructing geoplot to include only those polygons whose centroids fall within the specified x-range. I’m not clipping individual border points; rather, I’m filtering out entire shapes based on where their centroids lie. If you look at the right side of the figure in #6, you’ll notice that a border line is missing — that’s because its centroid does not fall within the specified range (compare this with #5).

                  If you do not want the selection (via the -if- qualifier), which applies at the level of each polygon rather than each coordinate point (and I would recommend this), you can follow the approach in #5, which introduces an opaque white line to cover the outer box, as you correctly deduced.

                  Comment


                  • #10
                    No need for tricks - simply add option -line- when clipping.

                    Code:
                    frame nuts3: geoframe rclip zoom1, into(nuts3_zoom1_lines) replace line
                    geoplot (line nuts3_zoom1_lines), tight

                    Click image for larger version

Name:	g1.png
Views:	1
Size:	99.9 KB
ID:	1780501

                    It may be convenient to have both variants of clipped shapes available:

                    Code:
                    geoplot (area nuts3_zoom1 i._ID) ///
                            (line nuts3_zoom1_lines, lc(black) lw(.5)), ///
                        tight legen(off)
                    Click image for larger version

Name:	g2.png
Views:	1
Size:	121.6 KB
ID:	1780502

                    ben

                    Comment


                    • #11
                      Thanks Ben Jann , I had misunderstood the description of the line option in the geoframe help file.

                      Comment

                      Working...
                      X