[R-sig-Geo] Mapping my own polygon?

[Michael Sumner]

I took a guess that your coordinates are not in EPSG:6487 but in global
Mercator (EPSG:3857) which seems to give a reasonable region from online
image servers.

https://gist.github.com/mdsumner/e6997c2f4a54c743e078aca8401537a0

If that looks ok?     Then, simply replace 6487 with 3857 in your code, but
also please take care to track down where your neighbourhood coordinates
for the x,y range came from.  3857 is the infamous global google mercator,
and that might be what you're using or simply close to it, you should make
sure you know. :)

I'm using in-dev code in the example, it's not because I want you to use
that or advise you to - it's just to keep a record of what I did. I
experimented with the code and scale to guess at what might be the problem.

If the pic in my gist is not on the right track I'm happy to follow up,
best of luck!

HTH, Mike



On Sun, Jun 11, 2023 at 5:27 AM Kevin Zembower via R-sig-Geo <
r-sig-geo using r-project.org> wrote:

> In my continuing work on reporting on US Census data for my
> neighborhood, I'd like to draw a map of the boundaries of it. I was
> successful in creating and printing an OSM basemap, with the US Census
> blocks that make up my neighborhood on it.
>
> Now, I'd like to create my own polygon, of the boundaries of my
> neighborhood, because the census blocks don't line up exactly with the
> neighborhood boundaries. I need help creating a polygon that I can
> submit to read_osm() that will correctly return an OSM map of my area.
>
> Here's what I've tried so far:
> ## Reproducible simple example:
> library(tidyverse)
> library(sf)
>
> rw <- data.frame( ## Simplified neighborhood rectangle
>      Longitude = c(-8528150, -8528500, -8528500, -8528150),
>      Latitude  = c( 4771475,  4771475,  4771880,  4771880)
> )
>
> rw ## Returns (as expected):
>
>  > rw
>    Longitude Latitude
> 1  -8528150  4771475
> 2  -8528500  4771475
> 3  -8528500  4771880
> 4  -8528150  4771880
>  >
>
> rw %>%
>      st_as_sf(coords = c("Longitude", "Latitude"), dim = "XY") %>%
>      st_set_crs(6487) %>% ## CRS 6487 is NAD83 (2011) Maryland in meters
>      st_transform(crs = 4269) ## CRS 4269 is NAD83
>
> ## Returns:
>
> Simple feature collection with 4 features and 0 fields
> Geometry type: POINT
> Dimension:     XY
> Bounding box:  xmin: 171.777 ymin: 24.65904 xmax: 171.7818 ymax: 24.66314
> Geodetic CRS:  NAD83
>                     geometry
> 1 POINT (171.7818 24.66192)
> 2 POINT (171.7806 24.65904)
> 3  POINT (171.777 24.66026)
> 4 POINT (171.7781 24.66314)
>  >
>
> I expected the POINTS to look like:
> 1 POINT (-76.61246 39.35010)
>
> Can anyone suggest what I'm doing wrong? Thanks so much in advance. I've
> worked on it all day today, without making much progress.
>
> -Kevin
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>


-- 
Michael Sumner
Software and Database Engineer
Australian Antarctic Division
Hobart, Australia
e-mail: mdsumner using gmail.com

	[[alternative HTML version deleted]]