首页 > 解决方案 > 如何在 R 中重新缩放和变换阿拉斯加的美国地图上投影阿拉斯加的空间点

问题描述

我正在尝试用 ALASKA 绘制美国地图,然后在其上投影一组点。我使用发布的功能fixupfix1以下链接重新定位阿拉斯加和夏威夷。

使用 ggplot2 在美国专题地图上重新定位阿拉斯加和夏威夷

然后投影具有纬度和经度值的我的点,dput如下所示

my_points= structure(c(44.334567, 44.209571, 44.049845, 44.752622, 44.791511, 
44.391792, 44.540124, 46.075669, 46.007611, 45.836781, 46.595384, 
45.703999, 45.47594, 44.715117, 44.385954, 42.804004, 43.349842, 
43.252618, 42.891499, 56.212978, 55.391604, 55.395194, 60.476929, 
61.709743, 62.76729, 62.346439, 59.93483, 64.472359, 64.88569, 
-122.047007, -122.256733, -123.42621, -122.128965, -122.578974, 
-122.497582, -122.435915, -121.998698, -122.347319, -122.466208, 
-122.459556, -123.755405, -123.725121, -123.887335, -123.831778, 
-123.610907, -122.728942, -123.026172, -124.070652, -131.638358, 
-131.195583, -132.408627, -151.081668, -149.231938, -149.693379, 
-150.019201, -158.190006, -146.926265, -147.249648), .Dim = c(29L, 
2L))

我的代码是:

require(maptools)
require(rgdal)

fixup <- function(usa,alaskaFix,hawaiiFix){
  alaska=usa[usa$STATE_NAME=="Alaska",]
  alaska = fix1(alaska,alaskaFix)
  proj4string(alaska) <- proj4string(usa)
  hawaii = usa[usa$STATE_NAME=="Hawaii",]
  hawaii = fix1(hawaii,hawaiiFix)
  proj4string(hawaii) <- proj4string(usa)
  usa = usa[! usa$STATE_NAME %in% c("Alaska","Hawaii"),]
  usa = rbind(usa,alaska,hawaii)
  return(usa)
}

fix1 <- function(object,params){
  r=params[1];scale=params[2];shift=params[3:4]
  object = elide(object,rotate=r)
  size = max(apply(bbox(object),1,diff))/scale
  object = elide(object,scale=size)
  object = elide(object,shift=shift)
  object
}

## download the USA shapefile to the computer
    setwd(tempdir())
download.file("https://dl.dropbox.com/s/wl0z5rpygtowqbf/states_21basic.zip?dl=1", 
              "usmapdata.zip", 
              method = "curl")

 #This is a mirror of http://www.arcgis.com/home/item.html?
 #id=f7f805eb65eb4ab787a0a3e1116ca7e5
unzip("usmapdata.zip")

然后读入我的shapefile. 使用rgdal

us = readOGR(dsn = "states_21basic",layer="states")

现在转换为等面积,并运行 fixup 函数:

usAEA = spTransform(us,CRS("+init=epsg:2163"))
usfix = fixup(usAEA,c(-35,1.5,-2800000,-2600000),c(-35,1,6800000,-1600000))
plot(usfix)

这些参数分别是阿拉斯加和夏威夷的旋转、缩放、x 和 y 位移,是通过反复试验获得的。仔细调整它们。甚至将夏威夷的比例参数更改为 0.99999 也将其从地球上发射出去,因为涉及的人数众多。

将其转回经纬度:

usfixLL = spTransform(usfix,CRS("+init=epsg:4326"))
plot(usfixLL)

现在,处理我的观点,将它们添加到生成的地图中

require(sp)
my_NEW_S<-as.data.frame(my_points)
colnames(my_NEW_S)<-c("y","x")
coordinates(my_NEW_S)=~x+y #column names of the lat long cols
proj4string(my_NEW_S)<- CRS("+init=epsg:2163")

# add my points to the produced map
points(my_NEW_S,col="red",pch=16,cex=0.5)

美国本土的所有点都投影在地图上,除了那些与阿拉斯加和 HI 有关的点。如果有人帮助我如何使我的阿拉斯加点投影与重新安置的阿拉斯加或任何其他解决方案相同,将不胜感激。谢谢你。

标签: rggmaprgdalmaptools

解决方案


推荐阅读