0

I'm trying to graph two shapefiles which should be overlapping, but they're nearly a hemisphere away from each other when graphed. The issue is not specific to graphing because I've explored it with things like st_nearest_feature and confirmed the issue is with my data somehow.

I noticed that the proj4string and epsb attributes of the two shapefiles, which RStudio displays after they're input by st_read, are different.

For my first file (which I'm reading from .shp):

epsg (SRID):    32645
proj4string:    +proj=utm +zone=45 +datum=WGS84 +units=m +no_defs

and for my second file (which I'm reading from .csv):

epsg (SRID):    NA
proj4string:    NA

Could this be the problem? If so, how do I fix it?

(I am also doing %>% st_set_crs(st_crs(my_first_shapefile)) because I saw an answer to a question like mine but it didn't seem to help. So the crs doesn't seem to be the problem.)

here are the messages when they're read in: (note the large bbox for first)

Reading layer X from data source `.shp' using driver `ESRI Shapefile'
Simple feature collection with 222 features and 12 fields
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 745400.7 ymin: 2402049 xmax: 998753.2 ymax: 2671392
epsg (SRID):    32645
proj4string:    +proj=utm +zone=45 +datum=WGS84 +units=m +no_defs
options:        X_POSSIBLE_NAMES=lat Y_POSSIBLE_NAMES=lon 
Reading layer Y from data source `.csv' using driver `CSV'
Simple feature collection with 5 features and 3 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 24.70942 ymin: 87.91372 xmax: 24.70942 ymax: 88.20073
epsg (SRID):    NA
proj4string:    NA

dputs

Here is the dput for the shapefile after it's read in:

structure(list(Transporta = 31159.5791334, Transpor_1 = 1, Transpor_2 = 2, 
    Transpor_3 = 0, TransportS = 1, Transpor_4 = 120.8226354, 
    Transpor_5 = 40.37383139, Transpor_6 = 0.334157844, Transpor_7 = 0.6603, 
    Transpor_8 = 0.8694, Transpor_9 = structure(NA_integer_, .Label = character(0), class = "factor"), 
    class = structure(1L, .Label = "0", class = "factor"), geometry = structure(list(
        structure(c(745439.842100001, 745546.5795, 745632.3839, 
        745788.981400001, 745917.343899999, 745985.294699999, 
        746066.299, 746154.675899999, 746368.8333, 746569.5261, 
        746772.2455, 746881.8819, 747004.401999999, 747073.7139, 
        747272.5434, 747470.2489, 747669.9781, 747812.230399999, 
        747975.2242, 748282.6014, 748497.1264, 748706.9185, 749043.816300001, 
        749292.6146, 749610.7259, 749925.135099999, 750001.7305, 
        750624.1876, 750825.9755, 750948.688199999, 751059.9489, 
        751163.9131, 751327.932799999, 751563.730499999, 751814.7569, 
        751971.032300002, 752126.386900001, 752195.9161, 752213.0906, 
        752386.998999999, 752474.422399999, 752562.3503, 752732.1897, 
        752984.445600001, 753125.1189, 753316.914400001, 753451.7143, 
        753804.5972, 753975.4802, 754070.401899998, 754153.675599999, 
        754287.556099999, 754405.755000002, 754531.939100001, 
        754559.107699999, 754588.9883, 754515.016200001, 754515.4296, 
        754545.5989, 754566.464399999, 754662.8412, 754914.5953, 
        755505.306800001, 756006.1144, 757647.819899998, 757818.6545, 
        757857.6266, 757878.117400001, 757888.392300001, 757979.6043, 
        758039.2786, 758104.023299999, 758273.8898, 758353.649, 
        758363.7606, 758425.7181, 758461.5542, 758412.883999999, 
        758361.6557, 758337.5731, 758369.8289, 758400.552799999, 
        758599.791, 758782.7463, 758948.5656, 759066.780899999, 
        759122.1718, 759140.5877, 759163.8683, 759310.9385, 759418.1858, 
        761587.9766, 761697.946600001, 761732.9896, 761735.1635, 
        761771.6095, 761834.645200001, 761880.7744, 761988.7883, 
        762117.364700001, 762363.861100001, 762461.7972, 762577.650000001, 
        2549587.0462, 2549520.9314, 2549439.3669, 2549255.2505, 
        2549117.1188, 2548990.8909, 2548851.1706, 2548730.2398, 
        2548549.1469, 2548369.9096, 2548184.485, 2548055.6101, 
        2547843.9775, 2547702.8617, 2547377.6804, 2547071.7279, 
        2546836.4809, 2546710.2235, 2546607.1267, 2546468.6995, 
        2546397.9558, 2546368.2405, 2546431.9041, 2546481.667, 
        2546520.1404, 2546550.2619, 2546534.1647, 2546540.4899, 
        2546584.3733, 2546647.1947, 2546721.98, 2546784.4899, 
        2546843.9544, 2546872.2121, 2546826.1245, 2546737.2423, 
        2546642.6158, 2546522.2403, 2546392.8828, 2545671.1662, 
        2545403.889, 2545106.76, 2544870.7301, 2544665.9495, 
        2544572.9174, 2544320.907, 2544112.5035, 2543574.9242, 
        2543279.2004, 2543016.8945, 2542750.9041, 2542376.3794, 
        2542042.1049, 2541679.6057, 2541467.5463, 2541190.3413, 
        2540905.4868, 2540780.0052, 2540631.2105, 2540525.5178, 
        2540373.1915, 2540070.0945, 2539482.8831, 2539046.8166, 
        2537839.0783, 2537724.5041, 2537640.0879, 2537539.1494, 
        2537377.2627, 2537091.1622, 2536893.6567, 2536619.2583, 
        2536123.8235, 2535849.6827, 2535679.4848, 2535292.3446, 
        2534815.1725, 2534426.1384, 2534186.3689, 2533976.9259, 
        2533708.7239, 2533530.0804, 2533205.0132, 2532931.801, 
        2532675.3421, 2532349.1915, 2532184.0268, 2531985.8156, 
        2531722.8626, 2531466.0844, 2531342.3258, 2531359.46, 
        2531296.537, 2531224.2141, 2531098.6506, 2530945.319, 
        2530557.4524, 2530278.6875, 2529895.6515, 2529626.4184, 
        2529270.0938, 2529090.8971, 2528794.2838), .Dim = c(103L, 
        2L), class = c("XY", "LINESTRING", "sfg"))), class = c("sfc_LINESTRING", 
    "sfc"), precision = 0, bbox = structure(c(xmin = 745439.842100001, 
    ymin = 2528794.2838, xmax = 762577.650000001, ymax = 2549587.0462
    ), class = "bbox"), crs = structure(list(epsg = 32645L, proj4string = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs"), class = "crs"), n_empty = 0L)), sf_column = "geometry", agr = structure(c(Transporta = NA_integer_, 
Transpor_1 = NA_integer_, Transpor_2 = NA_integer_, Transpor_3 = NA_integer_, 
TransportS = NA_integer_, Transpor_4 = NA_integer_, Transpor_5 = NA_integer_, 
Transpor_6 = NA_integer_, Transpor_7 = NA_integer_, Transpor_8 = NA_integer_, 
Transpor_9 = NA_integer_, class = NA_integer_), .Label = c("constant", 
"aggregate", "identity"), class = "factor"), row.names = 1L, class = c("sf", 
"data.frame"))

and here is the dput for the second shapefile (from the csv):

structure(list(width = structure(c(1L, 2L, 5L, 3L, 4L), .Label = c("2.4716396025460600", 
"2.6281478054445400", "3.0063089766321100", "3.3519104918569500", 
"4.263698279008270"), class = "factor"), lat = c(24.7094195311052, 
24.7094195311052, 24.7094195311052, 24.7094195311052, 24.7094195311052
), lon = c(87.9137151118854, 87.9139846064706, 88.1161055453975, 
88.1996488668206, 88.2007268451616), geometry = structure(list(
    structure(c(24.7094195311052, 87.9137151118854), class = c("XY", 
    "POINT", "sfg")), structure(c(24.7094195311052, 87.9139846064706
    ), class = c("XY", "POINT", "sfg")), structure(c(24.7094195311052, 
    88.1161055453975), class = c("XY", "POINT", "sfg")), structure(c(24.7094195311052, 
    88.1996488668206), class = c("XY", "POINT", "sfg")), structure(c(24.7094195311052, 
    88.2007268451616), class = c("XY", "POINT", "sfg"))), class = c("sfc_POINT", 
"sfc"), precision = 0, bbox = structure(c(xmin = 24.7094195311052, 
ymin = 87.9137151118854, xmax = 24.7094195311052, ymax = 88.2007268451616
), class = "bbox"), crs = structure(list(epsg = NA_integer_, 
    proj4string = NA_character_), class = "crs"), n_empty = 0L)), sf_column = "geometry", agr = structure(c(width = NA_integer_, 
lat = NA_integer_, lon = NA_integer_), .Label = c("constant", 
"aggregate", "identity"), class = "factor"), row.names = c(NA, 
5L), class = c("sf", "data.frame"))

setting crs or proj4string doesn't help

I've tried setting the crs of the second to the crs of the first and it didn't help. I also tried changing the proj4string of the second to the proj4string of the first and that didn't help either

# make a copy of second_sf
second_sf_transformed <- second_sf

# make EPSG of the copy of the second the same as the first
st_crs(second_sf_transformed) <- st_crs(first_sf)

# make the proj4string of the copy of the second the same as the first
second_sf_transformed <- st_transform(second_sf_transformed, '+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs')

Despite the first shapefile being the one with the large bbox, that is graphed correctly, and the one with correct X-Y coordinates isn't

JohnDoeVsJoeSchmoe
  • 671
  • 2
  • 8
  • 25
  • Hopefully, they're the same projection and the second one lost its geographic info somehow. What returns if you call `crs(secondshapefile)`? If NA, you could try defining the second's projection as the first's, with `crs(secondshapefile) <- crs(firstshapefile)` and then plot and see if it gets you closer. – dbo Jul 29 '19 at 02:26
  • Could you share the data? – agila Jul 29 '19 at 11:59
  • @agila I put dputs and input messages in the body of the question. Do those help? Please let me know if there's anything else I can do to make it easier to be helped. – JohnDoeVsJoeSchmoe Jul 29 '19 at 15:40
  • You can't simply change of CRS to the second shapefile until you set it's proper CRS. What's the CRS of the second shapefile? – agila Jul 29 '19 at 18:07
  • It's NA until I do `st_crs(second_sf) <- st_crs(first_sf)` but that doesn't seem to actually help at all. Once I do that it's the same as the first one: `EPSG: 32645 proj4string: "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs"` – JohnDoeVsJoeSchmoe Jul 29 '19 at 19:24
  • 1
    st_crs function is used to set CRS to a sf object with no CRS but, of course, your shapefiles do not have the same CRS. You should first identify CRS of the second shapefile, let's call it CRS_X, then you should set st_crs(second_sf) <- CRS_X and, after that, you can transform the second shapefile as the same CRS as the first one using st_transform. – agila Jul 29 '19 at 19:52
  • I don't know what's the CRS of the second shapefile but you should look for it where you got the shapefile. – agila Jul 29 '19 at 19:54

0 Answers0