绘制R中地理位置的LAT/LON坐标

Question

我有一张海冰地图作为"geotiff".最终的目标是在特定的纬度/经度坐标处提取海冰浓度.

地理信息可在以下网址找到:http://www.iup.uni-bremen.de:8084 /amsr2data /asi_daygrid_swath/n6250/2015/jan/as-AMSR2-n6250-20150101-v5.tif

我想要做的是使用raster()加载geotiff,然后通过我的位置覆盖它,然后使用extract()函数从特定位置的栅格文件中获取值.

但是我的纬度/经度积累在地图的中心.我哪里出错了？非常感谢任何帮助或输入!

library(raster)
library(sp)

r1 = raster("test.tif")

##check plot
plot(r1)

## check projection
projection(r1)

mydf <- structure(list(longitude = rep(22,7), latitude = seq(60,90,5)),.Names = c("longitude","latitude"), class = "data.frame", row.names = c(NA, -7L)) 

### Get long and lat from data.frame.
xy <- mydf[,c(1,2)]
spdf <- SpatialPointsDataFrame(coords = xy, data = mydf,
                           proj4string = CRS("+proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"))

 points(spdf, col="red", lwd=2)

 ##extract RGB values as reference for sea-ice concentration
 seaice_conc = extract(r1, spdf)

Answer 1

Geo-sp的解决方案可行,但不是最佳的(缓慢且不精确).您应该始终(重新)投影矢量(在这种情况下为点)数据,而不是栅格数据.投影栅格数据会更改值,而矢量数据则不是这种情况.投影栅格数据的计算量也要大得多.

因此,你应该这样做:

library(raster)

r <- raster("asi-AMSR2-n6250-20150101-v5.tif")
crs(r)
# CRS arguments:
# +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 

df <- data.frame(longitude = rep(22,7), latitude = seq(60,90,5), ID=1:7)

spdf <- SpatialPointsDataFrame(coords = df[,1:2], data = df,
         proj4string = CRS("+proj=longlat +datum=WGS84"))


library(rgdal)
p <- spTransform(spdf, crs(r))       

extract(r, p)

重要的是要注意你犯了一个很常见的错误:

spdf <- SpatialPointsDataFrame(coords = xy, data = mydf, proj4string =
          CRS("+proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m"))

您创建SpatialPointsDataFrame并分配错误的坐标参考系统(crs).您必须分配与您的数据实际匹配的crs,即"+proj=longlat +datum=WGS84".之后,您可以将数据转换为您希望拥有的crs(使用spTransform).