【发布时间】：2018-06-26 15:05:06
【问题描述】：
我正在尝试使用 akima 包 interp 函数将一些网格值插入到路段上的一个点。由于那些时候数据丢失，数据点并不总是在那里。例如，请参见下图。
显示的 lat/lng 点是具有来自我的数据集的 dBz 值的点（如果没有值，我们的数据库不会记录条目）。我的目标是使用此网格将 dBz 值插入到特定的纬度/经度点（即纬度 = 39.95045 和经度 = -86.3574）
我的数据如下...
mrms
tstamp: POSIXct: format: "2018-06-09 16:56:00" ...
dBz: int 27 26 50 ...
latitude: num 39.9 39.9 39.9 ...
longitude: num -86.4 -86.4 -86.4
如果我对特定时间的数据进行子集化（因为每个时间戳有多个条目。我的数据如下所示。
tstamp dBz latitude longitude
1 6/9/2018 16:58 50 39.92501 -86.37501
2 6/9/2018 16:58 52 39.93498 -86.36501
3 6/9/2018 16:58 29 39.93498 -86.355
4 6/9/2018 16:58 38 39.91499 -86.37501
5 6/9/2018 16:58 44 39.92501 -86.36501
6 6/9/2018 16:58 23 39.92501 -86.355
7 6/9/2018 16:58 52 39.93498 -86.37501
8 6/9/2018 16:58 50 39.94498 -86.36501
9 6/9/2018 16:58 23 39.90501 -86.36501
10 6/9/2018 16:58 51 39.94498 -86.37501
11 6/9/2018 16:58 18 39.90501 -86.37501
12 6/9/2018 16:58 37 39.91499 -86.36501
在对数据进行子集化后，我运行以下命令...
interpolated <- interp(x = mrms_sub$longitude, y = mrms_sub$latitude, z = mrms_sub$dBz, xo = road$longitude, yo = road$latitude, linear = TRUE)
产生 x = num -86.4、y = num 40 和 z = num NA
我不知道为什么 interp 不断产生 NA 值。当我通过 csv 文件读取所有数据时，我以前使用过这种方法。我已经检查了我的 csv 方法的数据类型与将数据从 SQL Server 拉入数据帧，它们是相同的。非常感谢任何想法或建议！
更新
require(tidyverse)
require(RODBC)
require(akima)
##### Road Details #####
cnx_wx <- odbcConnect("dbName", uid = "", pwd = "")
cnx_rd <- odbcConnect("dbName", uid = "", pwd = "")
queryString <-
"SELECT
XDSegID,
RoadNumber,
RoadName,
Bearing,
geog.EnvelopeCenter().STAsText() as roadCenter
FROM [itsdb1].[inrix_xd].[dbo].[__xd] as inrix_xd
WHERE XDSegID = 1363484955 AND version = '2018-04-25'"
road <- sqlQuery(cnx_rd, queryString)
#Massage road center lat long points
road$roadCenter <- gsub("POINT|\\(|\\)", "", road$roadCenter)
road = separate(road, roadCenter, into = c("junk", "longitude", "latitude"), sep = " ")
road$longitude <- as.numeric(road$longitude)
road$latitude <- as.numeric(road$latitude)
##### Speed Data #####
queryString <-
"SELECT
CAST(tstamp as smalldatetime) as tstamp,
xdid,
speed,
score,
cvalue
FROM itsdb1.inrix_xd.dbo.xdspeeds
WHERE tstamp >= '2018-06-09 15:00:00' AND tstamp <= '2018-06-09 18:00:00'
AND xdid = 1363484955
ORDER BY tstamp"
speed <- sqlQuery(cnx_rd, queryString)
#attr(speed$tstamp, "tzone") <- "GMT"
odbcClose(cnx_rd)
##### MRMS Data #####
queryString <-
"SELECT
CAST(tstamp as smalldatetime) as tstamp,
dBz,
latitude,
longitude
FROM weather_mrms.dbo.v_seamless_hsr
WHERE tstamp >= '2018-06-09 15:00:00' AND tstamp <= '2018-06-09 18:00:00'
AND latitude BETWEEN 39.90 AND 39.95 AND longitude BETWEEN -86.38 AND -86.29
ORDER BY tstamp"
mrms <- sqlQuery(cnx_wx, queryString)
odbcClose(cnx_wx)
##### Combined Data #####
numOfwxLevels <- nlevels(as.factor(mrms$tstamp))
i <- 1
obsTime <- mrms$tstamp[2]
output <- data.frame(0,0)
colnames(output) <- c("tstamp", "Interpolated dBz")
while (i <= numOfwxLevels) {
mrms_sub <- filter(mrms, mrms$tstamp == obsTime)
interpolated <- akima::interp(x = mrms_sub$longitude, y = mrms_sub$latitude, z = mrms_sub$dBz,
xo = road$longitude, yo = road$latitude, linear = TRUE)
val <- c(obsTime, interpolated$z)
newTable <- rbind(output, val)
obsTime <- obsTime + 120
i <- i + 1
}
【问题讨论】：
-
您的预期输出是什么？
road数据是什么？ -
我期待某种 dBz 值在 40-50 左右。我很抱歉一开始没有说清楚。还有什么我遗漏的吗？编辑->道路数据只是一个数据框，其中包含要插入的经纬度点，它还有一些其他信息，但这些对插值都不重要。
-
您用来创建输出网格的点是什么？我只是想复制它。
-
我想我明白你在问什么，我用于输出网格的点是道路 $ 经度，道路 $ 纬度。我只对插值感兴趣的一点是我的“网格”。也许我误解了如何使用 interp 功能，如果是这样，我很抱歉。我仍在努力真正掌握 R。
-
你感兴趣的经纬度值是多少？
标签： r interpolation