4

我正在尝试绘制一个地图小型多重网格，该网格显示自 1900 年以来与佛罗里达州相交的飓风/热带风暴。我使用一些空间查询对该项目的所有大西洋风暴的数据库进行了子集化。

我现在正在佛罗里达州、一些毗邻的州、佛罗里达州的一些主要城市，当然还有奥基乔比湖的多边形顶部绘制我有限数量的飓风轨迹的线shapefile 。这是简单的代码：

library(maptools)
library(gdata)
library(RColorBrewer)

setwd("~/hurricanes")

# read shapefiles
florida <- readShapePoly("florida.shp")
south <- readShapePoly("south.shp")
hurricanes <- readShapeLines("hurricanes-florida.shp")

cities <- readShapePoints("cities.shp")
lakes <- readShapePoly("lakes.shp")

# miami, orlando and tallahassee (in FL)
cities <- subset(cities, ST == "FL")

# don't need ALL the 'canes
hurricanes1900 <- subset(hurricanes, Season >= 1900)

mycolors <- brewer.pal(5, "YlOrRd")

pdf(file = "hurricanemaps.pdf", ,width=8,height=20)
    par(mfrow=c(15,5), mar=c(1,1,1,1))
    for(i in 1:nrow(hurricanes1900)) 
        { 
        plot(south, col="#e6e6e6", border = "#999999")
        plot(florida, col="#999999", border = "#999999", add = TRUE)
        plot(lakes, col="#ffffff", border = "#999999", add = TRUE)
        plot(cities, pch=10, cex=.1,col="#000000", bg="#e38d2c", lwd=1, add = TRUE)
        plot(hurricanes1900[i,], col = mycolors[cut(hurricanes$MAX_Wind_W, breaks = 5)],
             lwd=3, add = TRUE); title(hurricanes1900$Title[i])
     }
dev.off()

我遇到的三个大问题：

1）循环给了我每场风暴的地图。我希望代码在网格中为每年（即使在没有风暴的年份）和当年的所有风暴生成佛罗里达/南部地图，最好带有标签。

2）我想在所有风暴中设置风速的颜色，而不仅仅是循环中每个特定行中的那些。这样，即使是一年中唯一的风暴，强风暴（如 1992 年的安德鲁）也会显得更暗。也许我可以通过重新编码一个可以相应设置样式的分类（H1、H2 等）变量来解决这个问题。

3) 假设我能找出第一条，我很难在每条风暴路径上渲染标签。maptools 文档不是很好。

无论如何，这是目前的输出（标题是 shapefile 中两个字段的串联）：

在此处输入图像描述

真正的问题是第 1 号。提前感谢您的帮助。

4

1 回答 1

0

鉴于没有可重复的数据，我为此演示收集了一些数据。请从下次开始为 SO 用户提供最小的可重现数据。这将帮助您获得更多帮助。

您想要实现的目标可以使用 ggplot2 来完成。如果你想为每一年绘制一张地图，你可以使用facet_wrap(). 如果要根据风添加颜色，可以在aes()绘制路径时执行此操作。如果你想添加飓风的名字，你可以使用 ggrepel 包，它可以让你轻松地提供注释。注意，如果要绘制平滑的路径，还需要数据处理。

library(stringi)
library(tibble)
library(raster)
library(ggplot2)
library(ggthemes)
library(ggrepel)
library(RColorBrewer)
library(data.table)

# Get some data. Credit to hmbrmstr for a few lines in the following code.

mylist <- c("http://weather.unisys.com/hurricane/atlantic/2007H/BARRY/track.dat",
            "http://weather.unisys.com/hurricane/atlantic/2007H/TEN/track.dat",
            "http://weather.unisys.com/hurricane/atlantic/2006H/ERNESTO/track.dat",
            "http://weather.unisys.com/hurricane/atlantic/2006H/ALBERTO/track.dat")

temp <- rbindlist(
            lapply(mylist, function(x){
                foo <- readLines(x)
                foo <- read.table(textConnection(gsub("TROPICAL ", "TROPICAL_",
                                  foo[3:length(foo)])), 
                                  header=TRUE, stringsAsFactors=FALSE)

                year <- stri_extract_first(str = x, regex = "[0-9]+")
                name <- stri_extract_first(str = x, regex = "[A-Z]{2,}")

                cbind(foo, year, name)

            }
        ))


### Add a fake row for 2017
temp <- temp %>%
        add_row(ADV = NA, LAT = NA, LON = NA, TIME = NA, WIND = NA,
                PR = NA, STAT = NA, year = 2017, name = NA)

### Prepare a map
usa <- getData('GADM', country = "usa", level = 1)
mymap <- subset(usa, NAME_1 %in% c("Florida", "Arkansas", "Louisiana",
                                   "Alabama", "Georgia", "Tennessee",
                                   "Mississippi",
                                   "North Carolina", "South Carolina"))
mymap <- fortify(mymap)


# Create a data.table for labeling hurricanes later.
label <- temp[, .SD[1], by = name][complete.cases(name)]

g <- ggplot() +
     geom_map(data = mymap, map = mymap,
              aes(x = long, y = lat, group = group, map_id = id),
                  color = "black", size = 0.2, fill = "white") +
    geom_path(data = temp, aes(x = LON, y = LAT, group = name, color = WIND), size = 1) +
    scale_color_gradientn(colours = rev(brewer.pal(5, "Spectral")), name = "Wind (mph)") +
    facet_wrap(~ year) +
    coord_map() +
    theme_map() +
    geom_text_repel(data = label,
                    aes(x = LON, y = LAT, label = name),
                    size = 2, 
                    force = 1,
                    max.iter = 2e3,
                    nudge_x = 1,
                    nudge_y = -1) +
    theme(legend.position = "right")

在此处输入图像描述

于 2016-09-02T02:40:18.180 回答