标签: spatial

将几何体从一个SRID转换/投影到另一个SRID

我有一个数据库表,目前在SRID 27700(英国国家电网)中保存几何数据.在检索数据时,我需要将其转换为SRID 4326(WGS84).有没有办法将PostGIS中找到的ST_Transform等函数应用到我的数据中以获得我需要的结果?

注意:解决方案需要能够使用T-SQL而不是存储过程等来实现.我必须能够构造一个语句并将其保存在表中作为字符串字段以便稍后检索.这是因为我的解决方案与数据库无关.

我目前在Oracle中这样做的方式如下:

select CLUSTER_ID, 
       NUM_POINTS, 
       FEATURE_PK, 
       A.CELL_CENTROID.SDO_POINT.X, 
       A.CELL_CENTROID.SDO_POINT.Y, 
       A.CLUSTER_CENTROID.SDO_POINT.X, 
       A.CLUSTER_CENTROID.SDO_POINT.Y, 
       TO_CHAR (A.CLUSTER_EXTENT.GET_WKT ()),  
       TO_CHAR (A.CELL_GEOM.GET_WKT ()), 
       A.CLUSTER_EXTENT.SDO_SRID 
from (SELECT CLUSTER_ID, 
             NUM_POINTS, 
             FEATURE_PK, 
             SDO_CS.transform (CLUSTER_CENTROID, 4326) cluster_centroid,
             CLUSTER_EXTENT, 
             SDO_CS.transform (CELL_CENTROID, 4326) cell_centroid, 
             CELL_GEOM FROM :0) a  
where sdo_filter( A.CELL_GEOM, 
                  SDO_CS.transform(mdsys.sdo_geometry(2003, :1, NULL, mdsys.sdo_elem_info_array(1,1003,3),mdsys.sdo_ordinate_array(:2, :3, :4, :5)),81989)) = 'TRUE'
Run Code Online (Sandbox Code Playgroud)

在使用PostGIS的PostgreSQL中,我这样做:

select CLUSTER_ID, 
       NUM_POINTS, 
       FEATURE_PK, ST_X(a.CELL_CENTROID), 
       ST_Y(a.CELL_CENTROID), 
       ST_X(ST_TRANSFORM(a.CLUSTER_CENTROID, 4326)),  
       ST_Y(ST_TRANSFORM(a.CLUSTER_CENTROID, 4326)), 
       ST_AsText(a.CLUSTER_EXTENT),  
       ST_AsText(a.CELL_GEOM), 
       ST_SRID(a.CLUSTER_EXTENT)  
FROM (SELECT CLUSTER_ID, 
      NUM_POINTS, 
      FEATURE_PK, 
      ST_TRANSFORM(ST_SetSRID(CLUSTER_CENTROID, 27700), 4326) cluster_centroid, 
      CLUSTER_EXTENT, 
      ST_TRANSFORM(ST_SetSRID(CELL_CENTROID, 27700), 4326) cell_centroid, 
      CELL_GEOM 
from …
Run Code Online (Sandbox Code Playgroud)

spatial sql-server-2008

6
推荐指数
1
解决办法
1万
查看次数

使用空间多边形更改栅格值

要更改下的栅格值,SpatialPoints您只需使用[.

r <- raster(system.file("external/test.grd", package="raster"))
rTp <- rasterToPoints(r, spatial = T)

set.seed(666)
rTpS <- rTp[sample(1:length(rTp), 500),]
plot(r)
plot(rTpS, add = TRUE, pch = ".")
Run Code Online (Sandbox Code Playgroud)

现在,空间点下方的栅格值更改如下:

r1 <- r
r1[rTpS] <- 99999
plot(r1)
Run Code Online (Sandbox Code Playgroud)

带有 SpatialPoints 的栅格(左)和具有更改值的栅格(右): 具有空间点的栅格(左)和具有更改值的栅格(右)

因此,我假设多边形的工作方式相同。因此,从此栅格数据集生成一些多边形:

rst <- as.integer(stretch(r, 0, 10))
rTpol <- rasterToPolygons(rst, dissolve = T)
rTpol <- rTpol[rTpol$layer > 3,]

plot(r)
plot(rTpol, add = T)

rst[rTpol[,]] <- 100
plot(rst)

# also tried 
# rst[rTpol] <- 100
Run Code Online (Sandbox Code Playgroud)

但是,当我尝试将栅格值更改为低于某个值时SpatialPolygons,不知何故更改了多边形的 ID(或其他值)。然而,它应该是100。

带有 SpatialPolygons 的栅格(左)和值发生奇怪变化的栅格(右): 带有 SpatialPolygons 的栅格(左)和值发生奇怪变化的栅格(右)

因此,我的问题是如何更改多边形内的栅格值。

plot r raster spatial r-raster

6
推荐指数
1
解决办法
2097
查看次数

将 POLYGON 聚合为 MULTIPOLYGON 并保留 data.frame

我有一个sfPOLYGON 几何类型的对象。我想使用分组属性 (group_attr) 将这些多边形聚合到 MULTIPOLYGON 中,并将新的 MULTIPOLYGON 对象与属性表连接起来。因此,作为结果,我将拥有一个sf两行三列的对象(group_attr、second_attr、geometry)。我已经尝试过使用st_cast- 它适用于sfc对象,但不适用于sf对象。是否可以使用sf包来做到这一点?

p1 <- rbind(c(0,0), c(1,0), c(3,2), c(2,4), c(1,4), c(0,0))
pol1 <-st_polygon(list(p1))
p2 <- rbind(c(3,0), c(4,0), c(4,1), c(3,1), c(3,0))
pol2 <-st_polygon(list(p2))
p3 <- rbind(c(4,0), c(4,1), c(5,1), c(5,0),c(4,0))
pol3 <-st_polygon(list(p3))
p4 <- rbind(c(3,3), c(4,2), c(4,3), c(3,3))
pol4 <-st_polygon(list(p4))

d <- data.frame(group_attr = c(1, 1, 2, 2), 
                second_attr = c('forest', 'forest', 'lake', 'lake'))
d$geometry <- st_sfc(pol1, pol2, pol3, pol4)
df<- st_as_sf(d)
plot(df)
df …
Run Code Online (Sandbox Code Playgroud)

aggregate r spatial r-sf

6
推荐指数
1
解决办法
2550
查看次数

基于最短地理距离匹配数据框

我有两个数据框,它们都包含纬度和经度坐标。第一个数据帧是对事件的观察,其中记录了位置和时间。第二个数据框是地理要素,其中记录了该要素的位置和信息。

my_df_1 <- structure(list(START_LAT = c(-33.15, -35.6, -34.08333, -34.13333, 
-34.31667, -47.38333, -47.53333, -34.08333, -47.38333, -47.15
), START_LONG = c(163, 165.18333, 162.88333, 162.58333, 162.76667, 
148.98333, 148.66667, 162.9, 148.98333, 148.71667)), row.names = c(1175L, 
528L, 1328L, 870L, 672L, 707L, 506L, 981L, 756L, 210L), class = "data.frame", .Names = c("START_LAT", 
"START_LONG"))

my_df_2 <- structure(list(latitude = c(-42.7984, -34.195, -49.81, -35.417, 
-28.1487, -44.657, -42.7898, -36.245, -39.1335, -31.8482), longitude = c(179.9874, 
179.526, -176.68, 178.765, -168.0314, 174.695, -179.9873, 177.7873, 
-170.0583, 173.2424), depth_top = c(935L, 2204L, 869L, 1973L, 
4750L, …
Run Code Online (Sandbox Code Playgroud)

r spatial latitude-longitude geographic-distance

6
推荐指数
1
解决办法
623
查看次数

R 中面板数据的空间回归

我有一个面板数据集,其中包含数百个区域、约 10 年以及这些区域的空间数据。我用包创建了一个权重矩阵spdep(通过标准方式,然后,nb2listw)。因此,我有一个矩阵,其中每个区域的权重(相对于其他区域) - 但每个区域仅表示一次。

我想运行spdep包 ( lagsarlm, errorsarlm) 中的一些空间回归,但出现错误:

Error in subset.listw(listw, subset, zero.policy = zero.policy) : Not yet able to subset general weights lists

Error in lagsarlm(y ~ x1 + x2: Input data and weights have different dimensions

我认为这是因为权重矩阵每个区域只有一行(然后只能计算一年)。您对如何解决这个问题有什么建议吗?我的想法围绕以下几点:

  • 扩展空间权重矩阵 OR
  • 告诉 spdep 这些区域将以相同的顺序重复(但如何重复?)

期待您的建议。

regression r panel spatial spdep

6
推荐指数
0
解决办法
750
查看次数

R - 如何围绕一个点绘制半径并使用该结果过滤其他点?

我希望在一个经纬度点周围绘制一个半径,然后使用该缓冲区来过滤适合其中的其他点。例如:

#stores datasets
stores = data.frame(store_id = 1:3,
                    lat = c("40.7505","40.7502","40.6045"),
                    long = c("-73.8456","-73.8453","-73.8012")
                    )

#my location
me  = data.frame(lat = "40.7504", long = "-73.8456")

#draw a 100 meter radius around me 


#use the above result to check which points in dataset stores are within that buffer
Run Code Online (Sandbox Code Playgroud)

不知道如何处理这个。我以前曾与over点和多边形相交,但不确定如何在孤立点上运行类似的场景。

r spatial r-sp r-sf

6
推荐指数
1
解决办法
2245
查看次数

EF Core Spatial Data:查询并获取以米为单位的距离

我可以查询并获取距客户最近的位置(有学位):

serviceQuery = serviceQuery.Where(s => s.Location.Distance(currentLocation) < <degree>)
    .OrderBy(s => s.Location.Distance(currentLocation));
Run Code Online (Sandbox Code Playgroud)

位置是NetTopologySuite.Geometries.PointC# 和 SQL Server 中的类型geography

如何在查询中提供米而不是度数,并获取米距离而不是度数?在客户端将结果度数转换为米是可以接受的。

sql-server spatial entity-framework-core asp.net-core

6
推荐指数
1
解决办法
5159
查看次数

是否有运行“空间向量自回归”的 R 包?

我正在寻找一个可以运行“空间向量自回归”的 R 包。\n tandfonline.com/doi/full/10.1080/17421770701346689

\n

根据Chen和Conley(2001)的说法,这是一个“向量自回归(VAR),其系数矩阵和冲击协方差矩阵是主体之间经济距离的函数。其他主体\xe2\x80\x99变量对条件均值的影响给定的主体\xe2\x80\x99s 变量是其与该主体的经济距离的函数。类似地,VAR 冲击的协方差是前一时期主体之间距离的函数,我们将这种属性称为各向同性。

\n

(Chen, X & Conley, TG (2001) 面板\n时间序列的新半参数空间模型,计量经济学杂志,105, 59\xe2\x80\x9383)

\n

然而,令人惊讶的是,我只能看到“空间自回归”,这仍然不是我的目的所需要的。我可以帮忙找到这个包吗?否则,我可以知道使用 R 编程运行此空间向量自回归模型的官方方法吗?

\n

r spatial package spatial-regression

6
推荐指数
1
解决办法
731
查看次数

空间 bam (gam) 中的时间自相关

我正在根据声学标签检测对河流中的鱼类深度进行建模(这意味着数据并不完全是完美间隔的连续时间序列)。我预测,深度会根据河流的空间位置而有所不同,因为不同的区域有不同的可用深度,一天中的不同时间,因为深度对光有反应,一年中的一天也有相同的原因,并且个体之间也有所不同。那么基本模型就是

深度 ~ s(lon, lat) + s(小时) + s(yday) + s(ID, bs="re")

有几百万次检测,所以模型很糟糕,所以

bam(深度 ~ s(经度、纬度) + s(小时) + s(yday) + s(ID, bs="re")

每个人的深度应该与之前的记录自相关(当然这取决于它最后一次注册的时间,但我不太知道如何解释时间上的离散间隔)。

我知道 rho 参数在 bam 中用作一种 corAR1 函数,我想这可以解释自相关。我还考虑将 lag(深度,by=ID) 作为预测变量,它的表现相当好,但我不确定这种方法的有效性。

我跟踪了几个面包屑发现 rho 可以从没有相关结构的模型中估计 rho<-acf(resid(m1),plot=FALSE)$acf 2 -

对于每个人,我添加了一个 ARSTART 变量来调用 AR.start = df$ARSTART 来考虑个体之间不同的时间序列 - 所以我的模型是

m2<-bam(depth~s(lon, lat)+s(yday)+s(hour, bs="cc")+s(fID, bs="re"), AR.start=df$ARSTART, discrete=T, rho=rho, data=df) 
Run Code Online (Sandbox Code Playgroud)

一切都进展顺利,根据 AIC,具有自相关结构的模型拟合得更好(更好),但效果的后验估计非常不准确(或缩放比例很差)。与没有结构的模型相比,根据 lon、lat 平滑器的空间效应变得极端(且均匀),其中空间平滑器似乎非常有效地捕获空间方差,表明它们被预测在较深的区域中更深且在较浅的区域较浅。

图 1. 深度的后验空间估计用于没有自相关结构的模型,但显然没有考虑深度和滞后深度将相关的现实(按 ID 嵌套)。注意 z 值在原始值范围内调整后验估计是有意义的(大多数探测深度为 0-4 m)

图 2. 当模型包含 rho 参数来解释时间自相关时,深度使用的后验空间估计。注意 z 的巨大尺度,它远远超出了我们测量的范围,并且根本无法解释

如果需要,我可以提供示例代码,但问题本质上是,与模型相比,自相关结构会如此显着地改变后验估计的值是否有意义,并且时间自相关结构是否吸收了所有方差是否与空间效应相关(在具有自相关结构的模型中似乎被否定)?

一些想法-我不知道什么是最好的:

  1. 盲目遵循 AIC,而没有真正理解为什么后验估计如此奇怪(巨大),或者为什么空间效应消失,尽管基于系统的生物学知识显然很重要
  2. 报告我们将自相关结构拟合到数据中,它拟合得很好,但没有改变关系的形状,因此我们呈现没有该结构的模型结果
  3. 没有自相关结构但使用 s(lagDepth) 变量作为固定效应的模型
  4. 模型的变化是深度而不是深度,这似乎消除了一些自相关。

非常感谢所有帮助-非常感谢

spatial smoothing gam mgcv autocorrelation

6
推荐指数
0
解决办法
723
查看次数

使用 R 进行 1 亿个点的空间连接

我有一个包含 1 亿个点的数据集。我需要在空间上将它们连接到它们相交的区域。我正在使用sfand st_join,但毫不奇怪,当我尝试执行此操作时,我的机器内存不足。

为了解决这个问题,我将数据分割成更易于管理的大小,并循环遍历它们以获得我想要的结果。然而,我想看看是否有更聪明的方法来实现我的最终目标(将区域 ID 附加到我的原始数据中)。

require("data.table")
require("dplyr")
require("sf")

DT <- data.table(X=runif(100000000, min=300000, max=500000),
                 Y=runif(100000000, min=200000, max=400000),
                 UID = 1:100000000,
                 Group = rep(1:50, each=100000000/50))
                             
Polygon <- st_read("https://ons-inspire.esriuk.com/arcgis/services/Administrative_Boundaries/Regions_December_2019_Boundaries_EN_BFE/MapServer/WFSServer?request=GetCapabilities&service=WFS") %>%
        select(-gml_id) %>%
        st_cast(., "GEOMETRYCOLLECTION") %>%
        st_collection_extract(., "POLYGON") %>% 
        group_by(across(1:(ncol(.)-1))) %>%
        summarise() %>%
        st_cast("MULTIPOLYGON") %>%
        ungroup()

DT_Split <- split(DT, by=c("Group"))

DT_Joined <- data.table()
    for (i in 1:length(DT_Split)){
    DT_Join <- st_join(st_as_sf(DT_Split[[i]], coords = c("X","Y"),crs=27700), Polygon, join = st_intersects)
    DT_Split_Join <- DT_Join %>% st_drop_geometry %>% as.data.table
    DT_Joined <- rbindlist(list(DT_Joined, DT_Split_Join))
    }
Run Code Online (Sandbox Code Playgroud)

r spatial data.table r-sf

6
推荐指数
0
解决办法
479
查看次数