benjann / geoplot

Stata module to draw maps
MIT License
28 stars 3 forks source link

Is it possible to zoom in on an aspect of a line? #17

Open hhiihh opened 8 months ago

hhiihh commented 8 months ago

Hi Ben,

Thanks for providing this package. I was wondering if it was possible to zoom in to focus on a line at a certain section?

For example, let's say that we have accidents mapped against roads, but want to zoom in on a certain section of road where most of the accidents cluster.

Using zoom it seems to want to grab the entire road as a layer. I am guessing that this is due to how the line is coded via the shapefile as a single object and not one that can be segmented via _CX and _CY coordinates.

Any thoughts on how this might be accomplished?

benjann commented 8 months ago

Yes, zoom() operates at the level of layers, i.e. the zoom will include everything that is part of the selected layers. If you only want to zoom in on a section of a road and the corresponding accidents, you would need to clip the road and accidents and then generate layers from the clipped data. Here is an example based on the data that I used in [https://ideas.repec.org/p/boc/lsug23/23.html](the talk at the UK Stata Conference); frame Borough contains the shapefile for boroughs in greater London, frame Accidents contains geolocated data on accidents. There are no roads, so I apply the clipping to boroughs. The same command can be used to clip roads.

local limits 531000 536000 183000 187000
frame Borough:   geoframe rclip `limits', into(Bclip) replace
frame Accidents: geoframe rclip `limits', into(Aclip) replace split
geoplot (area Borough) (point Accidents, msymbol(p) mcolor(gray)) ///
    (area Bclip, fc(white)) (point Aclip, msymbol(p) mcolor(cranberry)) ///
    , zoom(3 4: 4 130 170, box lc(blue))

image

ben

hhiihh commented 8 months ago

Hi Ben,

Thanks for the reply. Yes, I was working off an older version of geoplot and have now updated to the latest build on Github (v1.1.6 30oct2023).

I can't seem to find the rclip options in the help file, is that something that might be available soon?

I have tried to adopt your example for my data and it generally works well. There are just some things which I wished to bring to your attention as potential feature requests:

  1. Can clipped points be used with weights? It seems that when I tried something like geoplot (point Aclip [w=number]) to get different proportioned points by say the number of incidents it did not quite work. I imagine that the layer has to be created somehow with the weighted points first before it can be clipped and not after the fact?
  2. I would like to add marker labels to the zoomed points. However, marker labels do not seem to be correct for the zoomed image if using something like geoplot (point Bclip , ml(name)). Some of the marker labels are correct, but others are wildly off. I think this might have to do with how the local `limits' is used to clip. I am guessing that you use some coordinate of lat long min max to define the clip, but I might be hitting some issues and not getting the right points during the clip.

For the sake of troubleshooting I will set up your worked example for the UK Stata Conference to narrow down any weird issues that might be present in my datasets.

benjann commented 8 months ago

geoframe rclip is available since the 05oct2023 update. Type help geoframe and then under "Manipulation" in the table of subcommands look for rclip.

ad 1: Weights should work without problem, but note that weights will be ineffective with msymbol(p); this is default behavior of Stata graphs; marker symbol p is a small dot that will not be scaled, for alternative symbols see help symbolstyle. Furthermore, you might need to play around with msize() or wmax() to get a pleasing result. Here is an example:

frame Aclip: generate W = runiform()
geoplot (area Borough) (point Accidents, msymbol(p) mcolor(gray)) ///
    (area Bclip, fc(white)) ///
    (point Aclip [w=W], msym(o) msize(.3) mcol(cranberry%30) mlcol(%0)) ///
    , zoom(3 4: 4 130 170, box lc(blue))
image

ad 2: If you apply clipping to an attribute frame that is linked to a shape frame, then the polygons in the shape frame will be clipped, but the centroids that might exist in the attribute frame (and which are used to place the labels) will not be updated automatically. You will need to do that yourself. In the above example, frame Borough is the attribute frame of boroughs, containing information such as ID, Name, and centroids, and the frame is linked to Borough_shp which contains the polygons of the boroughs. So geoframe rclip applied to Borough wil clip the shapes in Borough_shp and remove units from Borough for which no shape is left in Borough_shp after the clipping, but it will not update the centroids in Borough. However, you can update the centroids by applying command geoframe generate centroid. Here is how my example graph looks like without fixing the centroids:

geoplot (area Borough) (point Accidents, msymbol(p) mcolor(gray)) ///
    (area Bclip, fc(white)) ///
    (point Aclip, msym(p) mcol(cranberry)) ///
    (point Bclip, ml(NAME)) ///
    , zoom(3/5: 4 130 170, box lc(blue))
image

Now I update the centroids and plot again:

frame Bclip: geoframe generate centroids, replace
geoplot (area Borough) (point Accidents, msymbol(p) mcolor(gray)) ///
    (area Bclip, fc(white)) ///
    (point Aclip, msym(p) mcol(cranberry)) ///
    (point Bclip, ml(NAME)) ///
    , zoom(3/5: 4 130 170, box lc(blue))
image

ben

hhiihh commented 8 months ago

Hi Ben,

Thanks for the reply.

Apologies, for some reason I was thinking that the rclip documentation was in the geoplot help, not the geoframe one.

This is an incredibly well featured package that you have created and has drastically simplified my charting workflow in Stata.

I will work through your reply and hopefully can recreate your charts above.

benjann commented 8 months ago

Unfortunately I did not yet have the time to produce much applied documentation for geoplot. In case you are interested, some brief tutorials ("Challenges") can found at doi.org/10.48350/188248 (clipping is not covered in these tutorials, though).

hhiihh commented 8 months ago

Hi Ben,

I have played around a bit and found out why I was having issues with the clipping.

My points are places that I have tried to geocode and some of them were not successfully found so it returned a missing value. I have left those missing values in to remind myself that some locations still need to be geocoded by hand or otherwise.

With the point geoframe containing missing _X _Y coordinates, when trying to use rclip it sometimes changes the coordinates in the new frame from the original. This causes things to generally line up in columns.

This might be useful to know.

benjann commented 8 months ago

I see. Yes, missing coordinates have special meaning in the geoplot system. They separate shape items. So if there are missings, at least in some contexts, this makes geoframe think that the coordinates are lines or polygons. It may be best to remove the missing points after loading the data. Alternatively, give each point has a unique ID, e.g. by typing

generate _ID = _n

This should also resolve the issue.

Would it be possible for you to share your code and data with me? This will maybe help me improve geoplot. ben

hhiihh commented 8 months ago

Hi Ben,

Apologies for the delay. I just decided to drop the missing observations for now, but it is useful to know that spatial data treats missing values differently.

What would be the best way to share information with you?

benjann commented 8 months ago

If total size is less than, say, 5mb you can send the files directly to me (ben.jann@gmail.com). Otherwise send me a download link. ben