st_contains 匹配:CRS 坐标系问题/不匹配

st_contains matching: CRS Coordinate System Problem / mismatch

我有两个数据框,df_points 和 df_polygon;我想根据给定多边形内的点进行匹配。但是,由于 CRS 不匹配,我 运行 遇到了一些错误。这就是我所做的:

首先,我需要这两个都是 sf class;我检查了 classes / CRS 并转换如下(运行 输出):

>class(df_polygon)
[1] "sf"         "data.frame"
>class(df_points)
[1] "data.frame"

>df_points <- st_as_sf(df_points, coords = c("X", "Y"))
>class(df_points)
[1] "sf"         "data.frame"

> st_crs(df_points)
Coordinate Reference System: NA
 st_crs(df_polygon)
Coordinate Reference System:
  User input: NAD83 / North Carolina (ftUS) 
  wkt:
PROJCRS["NAD83 / North Carolina (ftUS)",
    BASEGEOGCRS["NAD83",
        DATUM["North American Datum 1983",
            ELLIPSOID["GRS 1980",6378137,298.257222101,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4269]],
    CONVERSION["SPCS83 North Carolina zone (US Survey feet)",
        METHOD["Lambert Conic Conformal (2SP)",
            ID["EPSG",9802]],
        PARAMETER["Latitude of false origin",33.75,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8821]],
        PARAMETER["Longitude of false origin",-79,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8822]],
        PARAMETER["Latitude of 1st standard parallel",36.1666666666667,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8823]],
        PARAMETER["Latitude of 2nd standard parallel",34.3333333333333,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8824]],
        PARAMETER["Easting at false origin",2000000,
            LENGTHUNIT["US survey foot",0.304800609601219],
            ID["EPSG",8826]],
        PARAMETER["Northing at false origin",0,
            LENGTHUNIT["US survey foot",0.304800609601219],
            ID["EPSG",8827]]],
    CS[Cartesian,2],
        AXIS["easting (X)",east,
            ORDER[1],
            LENGTHUNIT["US survey foot",0.304800609601219]],
        AXIS["northing (Y)",north,
            ORDER[2],
            LENGTHUNIT["US survey foot",0.304800609601219]],
    USAGE[
        SCOPE["unknown"],
        AREA["USA - North Carolina"],
        BBOX[33.83,-84.33,36.59,-75.38]],
    ID["EPSG",2264]]

>df_points<- st_transform(df_points, crs = st_crs(df_polygon))
Error in st_transform.sfc(st_geometry(x), crs, ...) : 
  cannot transform sfc object with missing crs

>combine <- st_contains(df_polygon, df_points)
Error in st_geos_binop("contains", x, y, sparse = sparse, prepared = prepared,  : 
  st_crs(x) == st_crs(y) is not TRUE

我目前的理解是我已经加载了这两个文件;一个不是空间文件,所以我把它变成了一个 (st_as_sf),然后我检查了 CRS 并使它们匹配 (st_transform)。但是,无论出于何种原因,匹配都失败了,因为缺少 CRS。因此,当我 运行 st_contains 时,crs 不匹配并且也失败了。

对我做错了什么有什么建议吗?我在 答案中看到我可以做一个 st_transform,我将 crs 设置为一个数字,但我不确定哪些数字是有效的。

谢谢!

[编辑] 根据评论添加(一些)dput; df_points 和 df_polygon 都有接近 550 万的观察值,所以只贴到最后;不确定这是否足够有用:

dput(df_polygons)


5509796L, 5509797L, 5509798L, 5509799L, 5509800L, 5509801L, 5509802L, 
5509803L, 5509804L, 5509805L, 5509806L, 5509807L, 5509808L), class = c("sf", 
"data.frame"), sf_column = "Shape", agr = structure(c(PARNO = NA_integer_, 
ALTPARNO = NA_integer_, OWNNAME = NA_integer_, IMPROVVAL = NA_integer_, 
LANDVAL = NA_integer_, PARVAL = NA_integer_, MAILADD = NA_integer_, 
MUNIT = NA_integer_, MCITY = NA_integer_, MSTATE = NA_integer_, 
MZIP = NA_integer_, SITEADD = NA_integer_, SUNIT = NA_integer_, 
SCITY = NA_integer_, STNAME = NA_integer_, SZIP = NA_integer_
), .Label = c("constant", "aggregate", "identity"), class = "factor"))

dput(df_points)

5), class = c("XY", 
        "POINT", "sfg")), structure(c(-84.033419, 35.085265), class = c("XY", 
        "POINT", "sfg")), structure(c(-84.033419, 35.085265), class = c("XY", 
        "POINT", "sfg"))), n_empty = 0L, precision = 0, crs = structure(list(
        input = NA_character_, wkt = NA_character_), class = "crs"), bbox = structure(c(xmin = -106.399281, 
    ymin = 26.436531, xmax = -70.927672, ymax = 44.229347), class = "bbox"), class = c("sfc_POINT", 
    "sfc"))), row.names = c(NA, 5447703L), class = c("sf", "data.frame"
), sf_column = "geometry", agr = structure(c(PERSON_ID = NA_integer_, 
MOVE_ID = NA_integer_, address = NA_integer_, city = NA_integer_, 
st = NA_integer_, zip = NA_integer_, oaddress = NA_integer_, 
ocity = NA_integer_, ostate = NA_integer_, ozipcode = NA_integer_, 
MOVE_COUNT_TU = NA_integer_, PERSON_ADDRESS_ID = NA_integer_, 
FULL_ADD = NA_integer_, OFULL_ADD = NA_integer_, ADDRESS_ID_RAW = NA_integer_, 
FULL_ADDR_C = NA_integer_, ADDRESS_ID_CLEAN = NA_integer_, Loc_name = NA_integer_, 
Status = NA_integer_, Score = NA_integer_, Match_type = NA_integer_, 
Match_addr = NA_integer_, Addr_type = NA_integer_), class = "factor", .Label = c("constant", 
"aggregate", "identity")))

如果 sf 对象的当前 CRS 已知,则只能将其转换为新的 CRS。您的 df_points 没有关联的 CRS,因此转换失败。

例如,如果 df_points 在 WSG84 中(一个合理的猜测,但您应该确定其来源),那么您可以

df_points_with_crs <- st_set_crs(df_points, 4326)

在转换新的(引用的)数据集并将其与多边形组合之前。