从 R 中的 Radolan 数据进行 GPS Lat/Long 查找 Table,将立体投影转换为 R 中的经纬度网格
Make a GPS Lat/Long Lookup Table from Radolan Data in R, convert stereographic projection into lat long grid in R
我正在使用此示例在 R 中加载 Radolan 数据
Example,
RADOLAN Specification
binary_filepath <- paste0("raa01-rw_10000-1708091950-dwd---bin")
header_end <- regexpr("[=12=]3", readLines(binary_filepath, 1))
rb_stream <- file(binary_filepath, "rb")
skip_temp <- readBin(rb_stream, "raw", n = header_end, endian = "little")
rb_data <- readBin(rb_stream, "raw", size = 1, n = 900*900*2, endian = "big")
rb_data <- rawToBits(rb_data)
rbi <- lapply(seq(0, 900*900-1, 1), function(x){
sum(as.integer(rb_data[(x*16+1):(x*16+12)]) * 2**seq(0, 11, 1))
})
rbi <- unlist(rbi)
rbi[1:16]
rb <- raster(t(matrix(rbi, ncol = 900, nrow = 900)))
rb <- flip(rb, "y")
radolan_proj <-
CRS("+proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +to_meter=1000 +no_defs")
extent(rb) <- extent(-523.4622, 376.5378, -4658.645, -3758.645)
projection(rb) <- radolan_proj
RADOLAN 数据是立体网格。
此网格中的每个点相距 1 公里。
class : RasterLayer
dimensions : 900, 900, 810000 (nrow, ncol, ncell)
resolution : 1, 1 (x, y)
extent : -523.4622, 376.5378, -4658.645, -3758.645 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +to_meter=1000 +no_defs
data source : in memory
names : layer
values : 0, 2500 (min, max)
我想对 GPS 数据进行 Lat/Long 查找 table。
我如何将此 rb 对象转换为 Lat/Long 网格?
谢谢!
这个rb
对象是一个光栅。如果你想要非投影地理坐标(wgs84),你只需要重新投影它:
library(rgdal)
library(raster)
library(maps) # to add country limits
# raster with Lat/Long non projected coordinates
rb_wgs84 <- projectRaster(rb, crs = "+proj=longlat +datum=WGS84 +no_defs")
# Show the raster in a plot
plot(rb_wgs84)
maps::map("world", add = TRUE)
# Get coordinates and values in a dataframe
xyz <- data.frame(coordinates(rb_wgs84), z = values(rb_wgs84))
这是您要找的吗?
我正在使用此示例在 R 中加载 Radolan 数据
Example, RADOLAN Specification
binary_filepath <- paste0("raa01-rw_10000-1708091950-dwd---bin")
header_end <- regexpr("[=12=]3", readLines(binary_filepath, 1))
rb_stream <- file(binary_filepath, "rb")
skip_temp <- readBin(rb_stream, "raw", n = header_end, endian = "little")
rb_data <- readBin(rb_stream, "raw", size = 1, n = 900*900*2, endian = "big")
rb_data <- rawToBits(rb_data)
rbi <- lapply(seq(0, 900*900-1, 1), function(x){
sum(as.integer(rb_data[(x*16+1):(x*16+12)]) * 2**seq(0, 11, 1))
})
rbi <- unlist(rbi)
rbi[1:16]
rb <- raster(t(matrix(rbi, ncol = 900, nrow = 900)))
rb <- flip(rb, "y")
radolan_proj <-
CRS("+proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +to_meter=1000 +no_defs")
extent(rb) <- extent(-523.4622, 376.5378, -4658.645, -3758.645)
projection(rb) <- radolan_proj
RADOLAN 数据是立体网格。 此网格中的每个点相距 1 公里。
class : RasterLayer
dimensions : 900, 900, 810000 (nrow, ncol, ncell)
resolution : 1, 1 (x, y)
extent : -523.4622, 376.5378, -4658.645, -3758.645 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +to_meter=1000 +no_defs
data source : in memory
names : layer
values : 0, 2500 (min, max)
我想对 GPS 数据进行 Lat/Long 查找 table。 我如何将此 rb 对象转换为 Lat/Long 网格?
谢谢!
这个rb
对象是一个光栅。如果你想要非投影地理坐标(wgs84),你只需要重新投影它:
library(rgdal)
library(raster)
library(maps) # to add country limits
# raster with Lat/Long non projected coordinates
rb_wgs84 <- projectRaster(rb, crs = "+proj=longlat +datum=WGS84 +no_defs")
# Show the raster in a plot
plot(rb_wgs84)
maps::map("world", add = TRUE)
# Get coordinates and values in a dataframe
xyz <- data.frame(coordinates(rb_wgs84), z = values(rb_wgs84))
这是您要找的吗?