[R-sig-Geo] reading into R IPCC data files

Robert Hijmans r.hijmans at gmail.com
Tue Jan 18 00:04:04 CET 2011


I agree with Andy, although dealing with ncdf files has its own problems
((although you can open most with raster: brick(filename))

So, for the record, reading these ASCII files is not that hard if you have
R. See below for an example that works on the one file I looked at. You
could easily extend this to automate downloading and processing all the
files with functions download.file() and uzip() 

Robert


library(raster)

filename <- HADCM3_A1F_DSWF_1980.mea

x <- readLines(filename)
x <- apply(matrix(x),1, function(x) gsub("  ", " ", x))
x <- apply(matrix(x),1, function(x) trim (gsub("  ", " ", x)))

start = which(substr(x,1,5)=="7008 ")+1
end = c(which(substr(x,1,4)=="IPCC")[-1] - 1, length(x))

nc = as.integer(substr(x[2], 8, 11))
nr = as.integer(substr(x[2], 13, 16))

# "Grid is 96 * 73  Month is Jan"         

if (length(start) != 12) {
	stop('something wrong')
} else {
	d = matrix(nrow=nc *nr, ncol=12)
	month = NULL
	for (m in 1:12) {
		d[,m] <- as.numeric(unlist(strsplit(x[start[m]:end[m]], " ")))
		mm <- x[start[m]-4]
		mm <- trim(substr(mm, nchar(mm)-3, nchar(mm)))
		month <- c(month, mm)
	}
}

b <- brick(nrow=nr, ncol=nc)
b <- setValues(b, d)
b <- rotate(b)
ln <- paste(month, ";", x[1], ";", x[3], ";", x[4])
layerNames(b) <- ln

plot(b, 1)
plot(b)

b <- writeRaster(b, 'filename.tif')







-- 
View this message in context: http://r-sig-geo.2731867.n2.nabble.com/reading-into-R-IPCC-data-files-tp5925687p5933842.html
Sent from the R-sig-geo mailing list archive at Nabble.com.



More information about the R-sig-Geo mailing list