作者AndrewShi (没有你的我)
看板R_Language
标题[问题] 选取资料特定区域作图
时间Thu Nov 8 16:08:11 2018
[问题类型]:
程式谘询(我想用R 做某件事情，但是我不知道要怎麽用R 写出来)
[软体熟悉度]:
入门(写过其他程式，只是对语法不熟悉)
[问题叙述]:
这笔资料为2010年全球的降雨资料，资料可从这里取得:
http://0rz.tw/JI056，下面的程
式码是已经取特定的经纬度范围(东亚)做月平均的降雨计算所画出来的图，我想请教的是
要怎麽做才能画出经纬度范围只有东亚这块(也就是把旁边地图空白的部分都去除)??
我试过改image的经纬度范围，但是会出现变数(precip1)的长度不等於经度乘上纬度的错
误讯息，而把变数(precip1)的长度设成跟经纬度一样或是在计算月平均时就放进跟选取
特定经纬度的范围一样大的矩阵里也不行，还是一开始在读降雨这个变数时就只读进特定
的区域也不行，下面加上#号的程式码是我上述试的方法，如果这些方法是有可行的，只
是程式码写的不对，还烦请各位大大给予指点，谢谢。
[程式范例]:
rm(list=ls())
library(ncdf4)
library(data.table)
TRMM_data <- "C:\\Users\\TOM\\Desktop\\R(资料库)\\TRMM资料\\trmm_2010.nc"
obs <- nc_open(TRMM_data)
print(obs)
lon <- ncvar_get(obs,"lon")
nlon <- dim(lon)
lat <- ncvar_get(obs,"lat",verbose = F)
nlat <- dim(lat)
#lon <- seq(60,150,by=0.25)
#lat <- seq(-15,60,by=0.25)
time <- ncvar_get(obs,"time")
tunits <- ncatt_get(obs, "time", "units")
precip <- ncvar_get(obs,"r")
#precip <- ncvar_get(obs,"r",start=c(240,140,1),count=c(361,261,365))
nc_close(obs)
time <- as.Date(time / 24, format = "%Y-%m-%d", origin = "2010-01-01")
#precip1 <- matrix(NA_real_,361,261)
precip1 <- matrix(NA_real_,1440,400)
for(i in c(240:600)){
for(j in c(140:400)){
precip1[i,j] <- mean(precip[i,j,time=(1:31)])
}
}
library(sp)
library(RColorBrewer)
library(maps)
#lon <- seq(60,150,by=0.25)
#lat <- seq(-15,60,by=0.25)
#nlon <- 361L
#nlat <- 261L
image(lon,lat,precip1,col=rev(brewer.pal(10,"RdBu")))
#image(lon=c(60,150),lat=c(-15,60),precip1,col=rev(brewer.pal(10,"RdBu")))
#precip2 <- matrix(precip1,ncol=length(lat),nrow=length(lon),byrow=F)
#image(lon,lat,precip2,col=rev(brewer.pal(10,"RdBu")))
map(database='world2',add= TRUE)
#map(database='world2',xlim=c(60,150),ylim=c(-15,60),add= TRUE)
[环境叙述]:
https://imgur.com/a/RW1YjV0
[关键字]:
画图 特定区域 nc档
--
※ 发信站: 批踢踢实业坊(ptt.cc), 来自: 220.136.193.182
※ 文章网址: https://webptt.com/cn.aspx?n=bbs/R_Language/M.1541664493.A.4EE.html
1F：→ andrew43: 看一下image()中xlim和ylim这二个参数怎麽用 11/08 16:43
2F：→ andrew43: 简单来说，你不用预先为了特定区域而改变计算任何code 11/08 16:44
3F：→ andrew43: 只要限制image()的xlim和ylim就好了。或许再加上asp=1。 11/08 16:46
4F：→ AndrewShi: andrew大~好像画不出来，画出来是下图这样子，且有等了 11/08 18:35
5F：→ AndrewShi: 超过10分钟，图还是这样(空白)。 11/08 18:35
7F：→ andrew43: 我在image()中加了xlim和ylim後，直到map()都一切正常 11/08 18:39
8F：→ AndrewShi: 阿...我知道了，lon和lat不能删掉，xlim和ylim加在後 11/08 18:41
9F：→ AndrewShi: 面就好，谢谢andrew大^^ 11/08 18:41