以高分辨率重新投影栅格(WGS84 到 BNG)
reproject raster (WGS84 to BNG) with large resolution
我正在做一个足够简单的操作,将栅格从 WGS84 重新投影到英国国家网格,但我想知道一些结果 post 重新投影。结果栅格的总和是完全不同的;这是由于 'projectRaster' 中的分辨率和双线性插值的工作方式造成的吗?
我从一个 csv 文件开始,它具有覆盖 -180 到 180 度范围的纬度和经度(仅限整数值)的全局数据,以及一些 z 值。这是子集,制作成投影 WGS84 的栅格并转换为 BNG(下面的子集):
x <- rep(c(-10:3), times = 10)
n <- 14
y <- rep(48:57, each=n)
z <- rnorm(n=140, mean=20, sd=5)
ind <- which(z %in% sample(z, 45))
z[ind]<-NA
df <- data.frame("x"=x,"y"=y,"value"=z)
bng <- '+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs'
wgs84 <- '+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs'
coordinates(df) = ~x + y
ras = rasterFromXYZ(df, crs=wgs84)
cellStats(ras,sum)
plot(ras)
ras_bng <- projectRaster(ras,crs=bng)
plot(ras_bng)
cellStats(ras_bng,sum)
总和从 1919 到 2625(在我的例子中),相当大的一部分。
是否纯粹是由于重投影在 'edges' 周围创建了这么多额外的单元格?如果我要重新投影到分辨率小得多(5 公里)的栅格,这会大大减少不同的总和吗?
谢谢,S
> ras
class : RasterLayer
dimensions : 10, 14, 140 (nrow, ncol, ncell)
resolution : 1, 1 (x, y)
extent : -10.5, 3.5, 47.5, 57.5 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0
data source : in memory
names : layer
values : 3.243918, 32.21532 (min, max)
> ras_bng
class : RasterLayer
dimensions : 12, 19, 228 (nrow, ncol, ncell)
resolution : 67900, 111000 (x, y)
extent : -375677.8, 914422.2, -343553.2, 988446.8 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894
data source : in memory
names : layer
values : 8.174968, 28.31331 (min, max)
这是因为 'ras' 中的 NA
个值。
取出z[ind]<-NA
后差别很小。这是因为 projectRaster
隐含地使用了 "na.rm=TRUE";也许需要争论来改变它。这不像在其他情况下那么简单,例如,只有一个从中插入的单元格可能是 NA
,在这种情况下,可能应该计算它。
在实践中很难找到这么多的 NA 值,通常它们仅限于(土地)的边缘。
顺便说一句WGS84
是一个数据,不是投影。您使用的投影将是 'longitude/latitude'(也称为其他名称),除了这也不是投影,因为投影的整个点是从这样的 angular 坐标到平面坐标。因此,您正在将一个坐标参考系 ('long/lat, WGS84') 转换为另一个 (BNG)。
我正在做一个足够简单的操作,将栅格从 WGS84 重新投影到英国国家网格,但我想知道一些结果 post 重新投影。结果栅格的总和是完全不同的;这是由于 'projectRaster' 中的分辨率和双线性插值的工作方式造成的吗?
我从一个 csv 文件开始,它具有覆盖 -180 到 180 度范围的纬度和经度(仅限整数值)的全局数据,以及一些 z 值。这是子集,制作成投影 WGS84 的栅格并转换为 BNG(下面的子集):
x <- rep(c(-10:3), times = 10)
n <- 14
y <- rep(48:57, each=n)
z <- rnorm(n=140, mean=20, sd=5)
ind <- which(z %in% sample(z, 45))
z[ind]<-NA
df <- data.frame("x"=x,"y"=y,"value"=z)
bng <- '+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs'
wgs84 <- '+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs'
coordinates(df) = ~x + y
ras = rasterFromXYZ(df, crs=wgs84)
cellStats(ras,sum)
plot(ras)
ras_bng <- projectRaster(ras,crs=bng)
plot(ras_bng)
cellStats(ras_bng,sum)
总和从 1919 到 2625(在我的例子中),相当大的一部分。
是否纯粹是由于重投影在 'edges' 周围创建了这么多额外的单元格?如果我要重新投影到分辨率小得多(5 公里)的栅格,这会大大减少不同的总和吗?
谢谢,S
> ras
class : RasterLayer
dimensions : 10, 14, 140 (nrow, ncol, ncell)
resolution : 1, 1 (x, y)
extent : -10.5, 3.5, 47.5, 57.5 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0
data source : in memory
names : layer
values : 3.243918, 32.21532 (min, max)
> ras_bng
class : RasterLayer
dimensions : 12, 19, 228 (nrow, ncol, ncell)
resolution : 67900, 111000 (x, y)
extent : -375677.8, 914422.2, -343553.2, 988446.8 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894
data source : in memory
names : layer
values : 8.174968, 28.31331 (min, max)
这是因为 'ras' 中的 NA
个值。
取出z[ind]<-NA
后差别很小。这是因为 projectRaster
隐含地使用了 "na.rm=TRUE";也许需要争论来改变它。这不像在其他情况下那么简单,例如,只有一个从中插入的单元格可能是 NA
,在这种情况下,可能应该计算它。
在实践中很难找到这么多的 NA 值,通常它们仅限于(土地)的边缘。
顺便说一句WGS84
是一个数据,不是投影。您使用的投影将是 'longitude/latitude'(也称为其他名称),除了这也不是投影,因为投影的整个点是从这样的 angular 坐标到平面坐标。因此,您正在将一个坐标参考系 ('long/lat, WGS84') 转换为另一个 (BNG)。