我想提取那些与多边形重叠的网格。我使用了“crop”和“mask”,但它无法正常工作并且没有正确检测到所有网格。
#US shapefile
Us <- readOGR(dsn ="C:/shapefile",layer="boundaries") # read shapefile
#ozone ncdf file
data <- brick(path="C:/",pattern = "nc", varname = "o3")
data_us <- crop(data, Us)
data_us2 <- mask(data_us, Us)
#urban shapefile
urban <- readOGR("C:/shapefile_urban") # read shapefile
urban_T <- st_as_sf(spTransform(urban, CRS("+proj=longlat +datum=WGS84 +no_defs
+ellps=WGS84 +towgs84=0,0,0"))) # Transform projection
#mask with urban
data_us3 <- crop(data_us2, urban_T)
data_us4 <- mask(data_us3, urban_T)
df_i <- as.data.frame(data_us4, xy = TRUE)
names(df_i) <- c("Longitude", "Latitude", "value")
df_i <- na.omit(df_i)
#Plot
ggplot(df_i, aes(x = Longitude, y = Latitude, fill = value)) +
geom_tile() +
#add US map layer with state boundaries
geom_map(data = map_data("state"), map = map_data("state"),
aes(long, lat, map_id = region),
color = "black", fill = NA,size = 0.5) +
expand_limits(x = map_data("state")$long, y = map_data("state")$lat)