Hello everyone,
I am working on a project where I map conflict data (csv format, geometry type = points) on polygons (kml files). Both layers are data.frames, share the same CRS, and the polygons are single feature.
The mapping (with st_join) works well before I apply st_buffer and st_union to the objects (see conflictMap_1225). However, after these two functions, st_join somehow does not connect the two layers (see conflictMap1225_buf). Please find the reprex attached.
Does someone know why st_join behaves like this and knows how to fix it? Or alternatively, could point me to another solution?
Many thanks, Dominique
library(dplyr)
#>
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#>
#> filter, lag
#> The following objects are masked from 'package:base':
#>
#> intersect, setdiff, setequal, union
library(sf)
#> Linking to GEOS 3.7.2, GDAL 2.4.2, PROJ 5.2.0
library(lwgeom)
#> Linking to liblwgeom 2.5.0dev r16016, GEOS 3.6.1, PROJ 4.9.3
#> Warning in fun(libname, pkgname): GEOS versions differ: lwgeom has 3.6.1 sf
#> has 3.7.2
#> Warning in fun(libname, pkgname): PROJ versions differ: lwgeom has 4.9.3 sf
#> has 5.2.0
library(tmap)
library(raster)
#> Loading required package: sp
#>
#> Attaching package: 'raster'
#> The following object is masked from 'package:dplyr':
#>
#> select
library(pryr)
#>
#> Attaching package: 'pryr'
#> The following object is masked from 'package:raster':
#>
#> subs
library(ggplot2)
library(rnaturalearth)
u = "https://arkham.enterprises/r-script/KML_1225.kml"
download.file(url = u, destfile = "/tmp/KML_1225.kml")
# load KML files (project polygons)
project1225_raw <- st_read("/tmp/KML_1225.kml")
#> Reading layer `KACP_AREAS' from data source `/private/tmp/KML_1225.kml' using driver `KML'
#> Simple feature collection with 28 features and 2 fields
#> geometry type: MULTIPOLYGON
#> dimension: XY
#> bbox: xmin: 34.24501 ymin: -0.400429 xmax: 34.58954 ymax: 0.81301
#> epsg (SRID): 4326
#> proj4string: +proj=longlat +datum=WGS84 +no_defs
# Read CSV file (conflict data, points)
conflictData <- read.csv("https://arkham.enterprises/r-script/GED_20200103_copy2.csv")
# Convert the CSV data frame to an sf object
conflictData_sf4326 <- st_as_sf(conflictData, coords = c("longitude", "latitude"), crs = 4326)
# Project CRS of conflict data and KML
conflictData_sf <- st_transform(conflictData_sf4326, 3395)
project1225 <- st_transform(project1225_raw, 3395)
# Buffer KML polygons
project1225_buf <- st_buffer(project1225, 3000)
# Dissolve buffered polygon into a single polygon
project1225_union_buf <- st_union(project1225_buf, by_feature = FALSE)
project1225_union_buf_sf <- st_as_sf(project1225_union_buf)
# Map conflict data with project data
conflictMap1225_buf <- st_join(conflictData_sf, project1225_union_buf_sf)
# That's how the mapping should look like
conflictMap_1225 <- st_join(conflictData_sf, project1225)
Created on 2020-01-07 by the reprex package (v0.3.0)