Converting a shapefile from easting/northing to lat/long

Question

I am working with an open data shapefile that uses easting/northing coordinates. My geocoded data set uses lat/long, so I am trying to reproject/transform the shapefile to use lat/long.

```{r}
nwk_shp <- st_read('http://data.ci.newark.nj.us/dataset/61038922-a789-4abe-9928-79f9fe96c307/resource/90ead04f-dffd-4b3a-93c7-39cd4b190d4f/download/superneighborhoodswgs84.geojson') %>% 
  as('Spatial')
```

Reading layer `OGRGeoJSON' from data source `http://data.ci.newark.nj.us/dataset/61038922-a789-4abe-9928-79f9fe96c307/resource/90ead04f-dffd-4b3a-93c7-39cd4b190d4f/download/superneighborhoodswgs84.geojson' using driver `GeoJSON'
Simple feature collection with 21 features and 8 fields
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 561009.6 ymin: 671183.6 xmax: 598979.6 ymax: 712226.2
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs

The shapefile uses northing/easting coordinates:

```{r}
nwk_shp@data[1:3, ] %>% select(-AREA, -PERIMETER, -CATEGORY)
```
1   0     18         LOWER CLINTON HILL   686419.4  575672.5
2   1      8 MT PLEASANT/LOWER BROADWAY   701184.5  583709.9
3   2     12            LOWER ROSEVILLE   701167.2  577849.7

Here is the projection string on my shapefile:

```{r}
proj4string(nwk_shp)
```
[1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

Based on this tutorial, and several others, I have tried several variations on spTransform - none of which seem to transform my data into expected lat/longitude coordinates.

For instance:

```{r}

nwk_shp_transformed <- spTransform(nwk_shp, CRS("+proj=longlat +datum=WGS84"))

nwk_shp_transformed@data[1:3, ] %>% select(-AREA, -PERIMETER, -CATEGORY)
```
 FID MSLINK                       NAME NORTHING_Y EASTING_X
1   0     18         LOWER CLINTON HILL   686419.4  575672.5
2   1      8 MT PLEASANT/LOWER BROADWAY   701184.5  583709.9
3   2     12            LOWER ROSEVILLE   701167.2  577849.7

Depending on what common lat/long CRS I use, I get a 'non-finite transformation detected':

```{r}

nwk_shp_again <- spTransform(nwk_shp, CRS("+init=epsg:4238"))

nwk_shp_again@data[1:3, ] %>% select(-AREA, -PERIMETER, -CATEGORY)
```
35 projected point(s) not finitenon finite transformation detected:
[1] 574389.9 685357.4      Inf      Inf
Error in .spTransform_Polygon(input[[i]], to_args = to_args, from_args = from_args,  : 
  failure in Polygons 1 Polygon 1 points 1

All of this leads me to suspect that there's something about either

• the proj string of the original shapefile?
• something omitted from the CRS I am using in my spTransform command?
  at work here.

I've read the 'Overview of Coordinate Reference Systems in R' but am still stuck.

Desired output would be a shapefile / SpatialPolygonsDataFrame with coordinates in lat / longitude - eg, Latitudes around 40 degrees, and Longitudes around -74 degrees.

Answer 1

It looks like you're correct with your first suspicion, the original proj4string is incorrectly defined as "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0". You can tell because the coordinates don't look like longitude/latitude decimal degrees.

If you knew the correct coordinate system of the original data, you could correct the proj4string like this

new_shp@proj4string <- CRS("+proj=the correct proj parameters")

and then call spTransform like you did

nwk_shp_transformed <- spTransform(nwk_shp, CRS("+proj=longlat +datum=WGS84"))

I tried this, guessing that the correct coordinate system of the original was web mercator, but that didn't seem to be correct either.

Looking at the original data source - http://data.ci.newark.nj.us/dataset/newark-neighborhoods/resource/90ead04f-dffd-4b3a-93c7-39cd4b190d4f - the problem is pretty clearly on their end.

But this might be a working version of the same? http://data.ci.newark.nj.us/dataset/newark-neighborhoods/resource/fcbbb04d-4684-439c-886f-929cb4e7c200

Answer 2

The coordinate values in the nwk_shp appear to be using a New Jersey State Plane zone, units = US survey feet, based on NAD 1983 or a later re-adjustment (like NAD 1983 (2011), which one is impossible to determine by inspection). The EPSG well-known ID for the NAD 1983 version is 3424.