[Raster-commits] r342 - in pkg/raster: . R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Mar 11 02:56:31 CET 2009
Author: rhijmans
Date: 2009-03-11 02:56:31 +0100 (Wed, 11 Mar 2009)
New Revision: 342
Added:
pkg/raster/R/read.R
pkg/raster/R/readBlock.R
pkg/raster/R/readCells.R
pkg/raster/R/readRaster.R
Removed:
pkg/raster/R/read.generic.R
pkg/raster/R/read.raster.R
Modified:
pkg/raster/NAMESPACE
pkg/raster/R/bilinearValue.R
pkg/raster/R/calc.R
pkg/raster/R/calcStack.R
pkg/raster/R/modal.R
pkg/raster/R/stackRead.R
pkg/raster/R/values.R
pkg/raster/man/read.Rd
pkg/raster/man/writeRaster.Rd
Log:
Modified: pkg/raster/NAMESPACE
===================================================================
--- pkg/raster/NAMESPACE 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/NAMESPACE 2009-03-11 01:56:31 UTC (rev 342)
@@ -1,8 +1,8 @@
importFrom("methods", Ops, Math)
importFrom("graphics", hist, plot, lines)
importFrom("stats", median, aggregate)
-importFrom("utils", stack)
+importFrom("utils", stack, unstack)
importFrom("sp", overlay, bbox, Spatial, SpatialPixels, SpatialPixelsDataFrame, SpatialGrid, SpatialGridDataFrame)
exportClasses(BoundingBox, BasicRaster, Raster, RasterLayer, RasterStack)
-exportMethods(calc, overlay, bbox, aggregate, stack, show, summary, plot, hist, ncol, nrow, ncell, dim, lines, median)
+exportMethods(calc, overlay, bbox, aggregate, stack, unstack, show, summary, plot, hist, ncol, nrow, ncell, dim, lines, median)
exportPattern("^[^\\.]")
\ No newline at end of file
Modified: pkg/raster/R/bilinearValue.R
===================================================================
--- pkg/raster/R/bilinearValue.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/bilinearValue.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -35,19 +35,18 @@
.bilinear <- function(x,y, x1,x2,y1,y2, q11,q12,q21,q22) {
div <- (x2-x1)*(y2-y1)
- if (all(div > 0)) {
- return( (q11/div)*(x2-x)*(y2-y) + (q21/div)*(x-x1)*(y2-y) + (q12/div)*(x2-x)*(y-y1) + (q22/div)*(x-x1)*(y-y1) )
- } else {
- print('oops, it happend')
- bil <- vector(length=length(div))
- bil[div>0] <- (q11/div)*(x2-x)*(y2-y) + (q21/div)*(x-x1)*(y2-y) + (q12/div)*(x2-x)*(y-y1) + (q22/div)*(x-x1)*(y-y1)
- bil[(x1==x2 && y1==y2)] <- q11
- div <- y2-y1
- bil[(x1==x2 && y1!=y2)] <- (q11/div)*(y2-y) + (q12/div)*(y-y1)
- div <- x2-x1
- bil[(x1!=x2 && y1==y2)] <- (q11/div)*(x2-x) + (q21/div)*(x-x1)
- }
+
+ bil <- ( (q11/div)*(x2-x)*(y2-y) + (q21/div)*(x-x1)*(y2-y) + (q12/div)*(x2-x)*(y-y1) + (q22/div)*(x-x1)*(y-y1) )
return(bil)
+
+# bil <- vector(length=length(div))
+# bil[div>0] <- (q11/div)*(x2-x)*(y2-y) + (q21/div)*(x-x1)*(y2-y) + (q12/div)*(x2-x)*(y-y1) + (q22/div)*(x-x1)*(y-y1)
+# bil[(x1==x2 && y1==y2)] <- NA
+# div <- y2-y1
+# bil[(x1==x2 && y1!=y2)] <- (q11/div)*(y2-y) + (q12/div)*(y-y1)
+# div <- x2-x1
+# bil[(x1!=x2 && y1==y2)] <- (q11/div)*(x2-x) + (q21/div)*(x-x1)
+
}
Modified: pkg/raster/R/calc.R
===================================================================
--- pkg/raster/R/calc.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/calc.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -10,10 +10,12 @@
}
setMethod('calc', signature(x='RasterLayer', fun='function'),
+
function(x, fun, filename="", overwrite=FALSE, filetype='raster', datatype='FLT4S', track=-1) {
if (length(fun(5)) > 1) {
stop("function 'fun' returns more than one value")
}
+
filename <- trim(filename)
outraster <- setRaster(x, filename)
outraster <- setDatatype(outraster, datatype)
Modified: pkg/raster/R/calcStack.R
===================================================================
--- pkg/raster/R/calcStack.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/calcStack.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -8,10 +8,13 @@
#mCalc <- function(...) { stop('mCalc has been replaced by generic function "calc"')}
setMethod('calc', signature(x='RasterStack', fun='function'),
+
function(x, fun, filename="", overwrite=FALSE, filetype='raster', datatype='FLT4S', track=-1) {
+
if (length(fun(seq(1:5))) > 1) {
stop("function 'fun' returns more than one value")
}
+
filename <- trim(filename)
outraster <- setRaster(x, filename)
outraster <- setDatatype(outraster, datatype)
Modified: pkg/raster/R/modal.R
===================================================================
--- pkg/raster/R/modal.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/modal.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -9,7 +9,7 @@
setGeneric("modal", function(x, ..., ties='random', na.rm=FALSE)
standardGeneric("modal"))
-
+
setMethod('modal', signature(x='ANY'),
function(x, ..., ties='random', na.rm=FALSE) {
#partly based on http://wiki.r-project.org/rwiki/doku.php?id=tips:stats-basic:modalvalue
@@ -27,19 +27,24 @@
return(z)
} else {
freq <- table(z)
- if (is.logical(z)){logic <- TRUE} else {logic <- FALSE}
- if (logic) {
- w <- as.logical(names(freq[max(freq)==freq]))
+ if (is.numeric(z)){
+ w <- as.numeric(names(freq[max(freq)==freq]))
+ } else if (is.logical(z)) {
+ w <- as.logical(freq[max(freq)==freq])
} else {
- w <- as.numeric(names(freq[max(freq)==freq]))
+ w <- names(freq[max(freq)==freq])
}
if (length(w) > 1) {
if (ties == 'lowest') {
w <- min(w)
- if (logic) { w <- as.logical(w) }
+ if (is.logical(z)) {
+ w <- as.logical(w)
+ }
} else if (ties == 'highest') {
w <- max(w)
- if (logic) { w <- as.logical(w) }
+ if (is.logical(z)) {
+ w <- as.logical(w)
+ }
} else if (ties == 'NA') {
w <- NA
} else { # random
@@ -52,4 +57,3 @@
}
)
-
Added: pkg/raster/R/read.R
===================================================================
--- pkg/raster/R/read.R (rev 0)
+++ pkg/raster/R/read.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -0,0 +1,82 @@
+# R function for the raster package
+# Author: Robert J. Hijmans
+# International Rice Research Institute. Philippines
+# contact: r.hijmans at gmail.com
+# Date : November 2008
+# Version 0.8
+# Licence GPL v3
+
+
+### readAll ###
+
+if (!isGeneric("readAll")) {
+ setGeneric("readAll", function(object)
+ standardGeneric("readAll"))
+}
+setMethod('readAll', signature(object='RasterLayer'),
+ function(object){ return(.rasterRead(object, -1))}
+)
+setMethod('readAll', signature(object='RasterStack'),
+ function(object){ return(.stackRead(object, -1))}
+)
+
+
+### readRow ###
+
+if (!isGeneric("readRow")) {
+ setGeneric("readRow", function(object, rownr)
+ standardGeneric("readRow"))
+}
+setMethod('readRow', signature(object='RasterLayer'),
+ function(object, rownr){ return(.rasterRead(object, rownr))}
+)
+setMethod('readRow', signature(object='RasterStack'),
+ function(object, rownr){ return(.stackRead(object, rownr))}
+)
+
+
+### readRows ###
+
+if (!isGeneric("readRows")) {
+ setGeneric("readRows", function(object, startrow, nrows=3)
+ standardGeneric("readRows"))
+}
+
+setMethod('readRows', signature(object='RasterLayer'),
+ function(object, startrow, nrows=3) {
+ #read multiple rows
+ return(.rasterReadBlock(object, startrow, nrows))
+ }
+)
+
+
+### readBlock ###
+
+if (!isGeneric("readBlock")) {
+ setGeneric("readBlock", function(object, startrow, nrows=3, startcol=1, ncolumns=(ncol(object)-startcol+1))
+ standardGeneric("readBlock"))
+}
+
+setMethod('readBlock', signature(object='RasterLayer'),
+ function(object, startrow, nrows=3, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
+ return(.rasterReadBlock(object, startrow, nrows, ncolumns))}
+)
+
+
+### readPartOfRow ###
+
+if (!isGeneric("readPartOfRow")) {
+ setGeneric("readPartOfRow", function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1))
+ standardGeneric("readPartOfRow"))
+}
+
+setMethod('readPartOfRow', signature(object='RasterLayer'),
+ function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
+ return(.rasterRead(object, rownr, startcol, ncolumns))}
+)
+
+setMethod('readPartOfRow', signature(object='RasterStack'),
+ function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
+ return( .stackRead(object, rownr, startcol, ncolumns) ) }
+)
+
Deleted: pkg/raster/R/read.generic.R
===================================================================
--- pkg/raster/R/read.generic.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/read.generic.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -1,82 +0,0 @@
-# R function for the raster package
-# Author: Robert J. Hijmans
-# International Rice Research Institute. Philippines
-# contact: r.hijmans at gmail.com
-# Date : November 2008
-# Version 0.8
-# Licence GPL v3
-
-
-### readAll ###
-
-if (!isGeneric("readAll")) {
- setGeneric("readAll", function(object)
- standardGeneric("readAll"))
-}
-setMethod('readAll', signature(object='RasterLayer'),
- function(object){ return(.rasterRead(object, -1))}
-)
-setMethod('readAll', signature(object='RasterStack'),
- function(object){ return(.stackRead(object, -1))}
-)
-
-
-### readRow ###
-
-if (!isGeneric("readRow")) {
- setGeneric("readRow", function(object, rownr)
- standardGeneric("readRow"))
-}
-setMethod('readRow', signature(object='RasterLayer'),
- function(object, rownr){ return(.rasterRead(object, rownr))}
-)
-setMethod('readRow', signature(object='RasterStack'),
- function(object, rownr){ return(.stackRead(object, rownr))}
-)
-
-
-### readRows ###
-
-if (!isGeneric("readRows")) {
- setGeneric("readRows", function(object, startrow, nrows=3)
- standardGeneric("readRows"))
-}
-
-setMethod('readRows', signature(object='RasterLayer'),
- function(object, startrow, nrows=3) {
- #read multiple rows
- return(.rasterReadBlock(object, startrow, nrows))
- }
-)
-
-
-### readBlock ###
-
-if (!isGeneric("readBlock")) {
- setGeneric("readBlock", function(object, startrow, nrows=3, startcol=1, ncolumns=(ncol(object)-startcol+1))
- standardGeneric("readBlock"))
-}
-
-setMethod('readBlock', signature(object='RasterLayer'),
- function(object, startrow, nrows=3, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
- return(.rasterReadBlock(object, startrow, nrows, ncolumns))}
-)
-
-
-### readPartOfRow ###
-
-if (!isGeneric("readPartOfRow")) {
- setGeneric("readPartOfRow", function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1))
- standardGeneric("readPartOfRow"))
-}
-
-setMethod('readPartOfRow', signature(object='RasterLayer'),
- function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
- return(.rasterRead(object, rownr, startcol, ncolumns))}
-)
-
-setMethod('readPartOfRow', signature(object='RasterStack'),
- function(object, rownr, startcol=1, ncolumns=(ncol(object)-startcol+1)) {
- return( .stackRead(object, rownr, startcol, ncolumns) ) }
-)
-
Deleted: pkg/raster/R/read.raster.R
===================================================================
--- pkg/raster/R/read.raster.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/read.raster.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -1,220 +0,0 @@
-# R code for reading raster (grid) data
-# Author: Robert J. Hijmans, r.hijmans at gmail.com
-# International Rice Research Institute
-# Date : June 2008
-# Version 0.8
-# Licence GPL v3
-
-
-#read a block of data (a rectangular area of any dimension)
-.rasterReadBlock <- function(raster, startrow, nrows=3, startcol=1, ncolumns=(ncol(raster)-startcol+1)) {
- if (startrow < 1 ) { stop("startrow too small") }
- if (startrow > nrow(raster) ) { stop("startrow too high") }
- if (nrows < 1) { stop("nrows should be > 1") }
- if (startcol < 1) { stop("startcol < 1") }
- if (startcol > ncol(raster)) { stop("startcol > ncol(raster)") }
- if (ncolumns < 1) { stop("ncolumns should be > 1") }
- if ((startcol + ncolumns - 1) > ncol(raster) ) {
- warning("ncolumns too high, truncated")
- ncolumns <- ncol(raster)-startcol }
-
- endrow <- startrow+nrows-1
- if (endrow > nrow(raster)) {
- warning("Rows beyond raster not read")
- endrow <- nrow(raster)
- nrows <- endrow - startrow + 1
- }
- raster <- .rasterRead(raster, startrow, startcol, ncolumns)
- blockvalues <- values(raster)
- if (nrows > 1) {
- for (r in (startrow+1):endrow) {
- raster <- .rasterRead(raster, r, startcol, ncolumns)
- blockvalues <- c(blockvalues, values(raster))
- }
- }
- raster at data@values <- blockvalues
- raster at data@content <- 'block'
- firstcell <- cellFromRowCol(raster, startrow, startcol)
- lastcell <- cellFromRowCol(raster, endrow, (startcol+ncolumns-1))
- raster at data@indices <- c(firstcell, lastcell)
-
- return(raster)
-}
-
-
-#read part of a single row
-.rasterRead <- function(raster, rownr, startcol=1, ncolumns=(ncol(raster)-startcol+1)) {
- rownr <- round(rownr)
- if (rownr == 0) { stop("rownr == 0. It should be between 1 and nrow(raster), or < 0 for all rows") }
- if (rownr > nrow(raster)) { stop("rownr too high") }
- if (startcol < 1) { stop("startcol < 1") }
- if (startcol > ncol(raster)) { stop("startcol > ncol(raster)") }
- if (ncolumns < 1) { stop("ncols should be > 1") }
-
- endcol <- startcol + ncolumns - 1
- if (endcol > ncol(raster)) {
- endcol <- ncol(raster)
- ncolumns <- ncol(raster) - startcol + 1
- }
-
- if (dataSource(raster)=='ram') {
- result <- valuesRow(raster, rownr)[startcol:endcol]
- } else if (.driver(raster) == 'raster') {
-
- rastergri <- .setFileExtensionValues(filename(raster))
- if (!file.exists( filename(raster))) {
- stop(paste(filename(raster)," does not exist"))
- }
- con <- file(rastergri, "rb")
-
- dtype <- .shortDataType(raster at file@datanotation)
- if (dtype == "INT" | dtype == "LOG" ) {
- dtype <- "integer"
- } else {
- dtype <- "numeric"
- }
- dsize <- dataSize(raster at file@datanotation)
- dsign <- dataSigned(raster at file@datanotation)
-
- if (rownr > 0) {
- seek(con, ((rownr-1) * ncol(raster) + (startcol-1)) * dsize)
- result <- readBin(con, what=dtype, n=ncolumns, dsize, dsign, endian=raster at file@byteorder) }
- else {
- result <- readBin(con, what=dtype, n=ncell(raster), dsize, dsign, endian=raster at file@byteorder)
- }
- close(con)
-# result[is.nan(result)] <- NA
- if (dtype == 'numeric') {
- result[result <= (0.999 * .nodatavalue(raster)) ] <- NA
- } else {
- result[result == raster at file@nodatavalue ] <- NA
- }
- if (dtype == 'logical') {
- result <- as.logical(result)
- }
- }
- else { #use GDAL
- if (is.na(raster at file@band)) { result <- NA }
- else {
- if (rownr <= 0) {
- offs <- c(0, 0)
- reg <- c(nrow(raster), ncol(raster)) # reg <- dim(raster at file@gdalhandle[[1]])
- }
- else {
- offs= c((rownr-1), (startcol-1))
- reg <- c(1, ncolumns)
- }
- }
- result <- getRasterData(raster at file@gdalhandle[[1]], offset=offs, region.dim=reg, band = raster at file@band)
-
-# if (!is.vector(result)) { result <- as.vector(result) }
-
- # if setNAvalue() has been used.....
- if (raster at file@nodatavalue < 0) {
- result[result <= raster at file@nodatavalue ] <- NA
- } else {
- result[result == raster at file@nodatavalue ] <- NA
- }
-
- }
-
- raster at data@values <- as.vector(result)
- if (rownr < 0) {
- raster at data@indices <- c(1, ncell(raster))
- raster at data@content <- "all"
- raster <- setMinMax(raster)
- } else if (startcol==1 & ncolumns==(ncol(raster)-startcol+1)) {
- raster at data@indices <- c(cellFromRowCol(raster, rownr, startcol), cellFromRowCol(raster, rownr, endcol))
- raster at data@content <- "row"
- } else {
- raster at data@indices <- c(cellFromRowCol(raster, rownr, startcol), cellFromRowCol(raster, rownr, endcol))
- raster at data@content <- "block"
- }
-
- return(raster)
-}
-
-
-
-#read data on the raster for cell numbers
-.rasterReadCells <- function(raster, cells) {
- uniquecells <- na.omit(unique(cells[order(cells)]))
- uniquecells <- uniquecells[(uniquecells > 0) & (uniquecells <= ncell(raster))]
- res <- cbind(cells, NA)
- if (length(uniquecells) > 0) {
- if (dataContent(raster) == 'all') {
- vals <- cbind(uniquecells, values(raster)[uniquecells])
- } else if (dataSource(raster) == 'disk') {
- if (.driver(raster) == 'gdal') {
- vals <- .readCellsGDAL(raster, uniquecells)
- } else {
- vals <- .readCellsRaster(raster, uniquecells)
- }
- } else {
- vals <- cbind(uniquecells, NA)
- }
- if (length(vals) == 2) {
- res[res[,1]==vals[1],2] <- vals[2]
- } else {
- for (i in 1:length(vals[,1])) {
- res[res[,1]==vals[i,1],2] <- vals[i,2]
- }
- }
- }
- return(res[,2])
-}
-
-
-.readCellsGDAL <- function(raster, cells) {
- colrow <- matrix(ncol=5, nrow=length(cells))
-# valuename <- raster at file@shortname
-# if (valuename == "") {valuename <- "value" }
-# colnames(colrow) <- c("id", "colnr", "rownr", "cell", valuename)
- for (i in 1:length(cells)) {
- colrow[i,1] <- colFromCell(raster, cells[i])
- colrow[i,2] <- rowFromCell(raster, cells[i])
- colrow[i,3] <- cells[i]
- colrow[i,4] <- NA
- }
- rows <- na.omit(unique(colrow[order(colrow[,2]), 2]))
- for (i in 1:length(rows)) {
- raster <- .rasterRead(raster, rows[i])
- thisrow <- subset(colrow, colrow[,2] == rows[i])
- for (j in 1:length(thisrow[,1])) {
- colrow[colrow[,3]==thisrow[j,3],4] <- raster at data@values[thisrow[j,1]]
- }
- }
- return(colrow[,3:4])
-}
-
-
-
-.readCellsRaster <- function(raster, cells) {
-# cells <- cbind(cells, NA)
-# valuename <- raster at file@shortname
-# if (valuename == "") {valuename <- "value" }
-# colnames(cells) <- c("id", "cell", valuename)
-# uniquecells <- na.omit(unique(cells[order(cells[,2]),2]))
-
- rastergri <- .setFileExtensionValues(filename(raster))
- if (!file.exists(filename(raster))) { stop(paste(filename(raster)," does not exist")) }
- con <- file(rastergri, "rb")
-
- res <- vector(length=length(cells))
- res[] <- NA
- dsize <- dataSize(raster at file@datanotation)
- dtype <- .shortDataType(raster at file@datanotation)
- for (i in 1:length(cells)) {
- seek(con, (cells[i]-1) * dsize)
- if (dtype == "FLT") {
- dtype <- "numeric"
- } else {
- dtype <- "integer"
- }
- res[i] <- readBin(con, what=dtype, n=1, size=dsize, endian=raster at file@byteorder)
- }
- close(con)
- res[res <= max(-3e+38, .nodatavalue(raster))] <- NA
- return(cbind(cells,res))
-}
-
Added: pkg/raster/R/readBlock.R
===================================================================
--- pkg/raster/R/readBlock.R (rev 0)
+++ pkg/raster/R/readBlock.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -0,0 +1,35 @@
+# Author: Robert J. Hijmans, r.hijmans at gmail.com
+# Date : June 2008
+# Version 0.8
+# Licence GPL v3
+
+
+#read a block of data (a rectangular area of any dimension)
+.rasterReadBlock <- function(raster, startrow, nrows=3, startcol=1, ncolumns=(ncol(raster)-startcol+1)) {
+ if (startrow < 1 ) { stop("startrow too small") }
+ if (nrows < 1) { stop("nrows should be >= 1") }
+
+ endrow <- startrow+nrows-1
+ if (endrow > nrow(raster)) {
+ if (options('verbose')[[1]]) {
+ warning("Rows beyond end of raster not read")
+ }
+ endrow <- nrow(raster)
+ nrows <- endrow - startrow + 1
+ }
+
+ blockvalues <- vector()
+ for (r in (startrow:endrow)) {
+ raster <- .rasterRead(raster, r, startcol, ncolumns)
+ blockvalues <- c(blockvalues, values(raster))
+ }
+
+ raster at data@values <- blockvalues
+ raster at data@content <- 'block'
+ firstcell <- cellFromRowCol(raster, startrow, startcol)
+ lastcell <- cellFromRowCol(raster, endrow, (startcol+ncolumns-1))
+ raster at data@indices <- c(firstcell, lastcell)
+
+ return(raster)
+}
+
Added: pkg/raster/R/readCells.R
===================================================================
--- pkg/raster/R/readCells.R (rev 0)
+++ pkg/raster/R/readCells.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -0,0 +1,89 @@
+# Author: Robert J. Hijmans, r.hijmans at gmail.com
+# Date : June 2008
+# Version 0.8
+# Licence GPL v3
+
+
+#read data on the raster for cell numbers
+.rasterReadCells <- function(raster, cells) {
+ uniquecells <- na.omit(unique(cells[order(cells)]))
+ uniquecells <- uniquecells[(uniquecells > 0) & (uniquecells <= ncell(raster))]
+ res <- cbind(cells, NA)
+ if (length(uniquecells) > 0) {
+ if (dataContent(raster) == 'all') {
+ vals <- cbind(uniquecells, values(raster)[uniquecells])
+ } else if (dataSource(raster) == 'disk') {
+ if (.driver(raster) == 'gdal') {
+ vals <- .readCellsGDAL(raster, uniquecells)
+ } else {
+ vals <- .readCellsRaster(raster, uniquecells)
+ }
+ } else {
+ vals <- cbind(uniquecells, NA)
+ }
+ if (length(vals) == 2) {
+ res[res[,1]==vals[1],2] <- vals[2]
+ } else {
+ for (i in 1:length(vals[,1])) {
+ res[res[,1]==vals[i,1],2] <- vals[i,2]
+ }
+ }
+ }
+ return(res[,2])
+}
+
+
+.readCellsGDAL <- function(raster, cells) {
+ colrow <- matrix(ncol=5, nrow=length(cells))
+# valuename <- raster at file@shortname
+# if (valuename == "") {valuename <- "value" }
+# colnames(colrow) <- c("id", "colnr", "rownr", "cell", valuename)
+ for (i in 1:length(cells)) {
+ colrow[i,1] <- colFromCell(raster, cells[i])
+ colrow[i,2] <- rowFromCell(raster, cells[i])
+ colrow[i,3] <- cells[i]
+ colrow[i,4] <- NA
+ }
+ rows <- na.omit(unique(colrow[order(colrow[,2]), 2]))
+ for (i in 1:length(rows)) {
+ raster <- .rasterRead(raster, rows[i])
+ thisrow <- subset(colrow, colrow[,2] == rows[i])
+ for (j in 1:length(thisrow[,1])) {
+ colrow[colrow[,3]==thisrow[j,3],4] <- raster at data@values[thisrow[j,1]]
+ }
+ }
+ return(colrow[,3:4])
+}
+
+
+
+.readCellsRaster <- function(raster, cells) {
+# cells <- cbind(cells, NA)
+# valuename <- raster at file@shortname
+# if (valuename == "") {valuename <- "value" }
+# colnames(cells) <- c("id", "cell", valuename)
+# uniquecells <- na.omit(unique(cells[order(cells[,2]),2]))
+
+ rastergri <- .setFileExtensionValues(filename(raster))
+ if (!file.exists(filename(raster))) { stop(paste(filename(raster)," does not exist")) }
+
+ con <- file(rastergri, "rb")
+
+ res <- vector(length=length(cells))
+ res[] <- NA
+ dsize <- dataSize(raster at file@datanotation)
+ dtype <- .shortDataType(raster at file@datanotation)
+ for (i in 1:length(cells)) {
+ seek(con, (cells[i]-1) * dsize)
+ if (dtype == "FLT") {
+ dtype <- "numeric"
+ } else {
+ dtype <- "integer"
+ }
+ res[i] <- readBin(con, what=dtype, n=1, size=dsize, endian=raster at file@byteorder)
+ }
+ close(con)
+ res[res <= max(-3e+38, .nodatavalue(raster))] <- NA
+ return(cbind(cells,res))
+}
+
Added: pkg/raster/R/readRaster.R
===================================================================
--- pkg/raster/R/readRaster.R (rev 0)
+++ pkg/raster/R/readRaster.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -0,0 +1,101 @@
+# R code for reading raster (grid) data
+# Author: Robert J. Hijmans, r.hijmans at gmail.com
+# International Rice Research Institute
+# Date : June 2008
+# Version 0.8
+# Licence GPL v3
+
+
+.rasterRead <- function(raster, rownr, startcol=1, ncolumns=(ncol(raster)-startcol+1)) {
+ rownr <- round(rownr)
+ ncolums <- round(ncolumns)
+ if (rownr == 0) { stop("rownr == 0. It should be between 1 and nrow(raster), or < 0 for all rows") }
+ if (rownr > nrow(raster)) { stop("rownr too high") }
+ if (startcol < 1) { stop("startcol < 1") }
+ if (startcol > ncol(raster)) { stop("startcol > ncol(raster)") }
+ if (ncolumns < 1) { stop("ncols should be > 1") }
+
+ endcol <- startcol + ncolumns - 1
+ if (endcol > ncol(raster)) {
+ warning("ncolumns too high, truncated")
+ endcol <- ncol(raster)
+ ncolumns <- ncol(raster) - startcol + 1
+ }
+
+ if (dataSource(raster)=='ram') {
+ result <- valuesRow(raster, rownr)[startcol:endcol]
+ } else if (.driver(raster) == 'raster') {
+# if dataContent(raster=='all')
+
+ rastergri <- .setFileExtensionValues(filename(raster))
+ if (!file.exists( filename(raster))) {
+ stop(paste(filename(raster)," does not exist"))
+ }
+ con <- file(rastergri, "rb")
+
+ dtype <- .shortDataType(raster at file@datanotation)
+ if (dtype == "INT" | dtype == "LOG" ) {
+ dtype <- "integer"
+ } else {
+ dtype <- "numeric"
+ }
+ dsize <- dataSize(raster at file@datanotation)
+ dsign <- dataSigned(raster at file@datanotation)
+
+ if (rownr > 0) {
+ seek(con, ((rownr-1) * ncol(raster) + (startcol-1)) * dsize)
+ result <- readBin(con, what=dtype, n=ncolumns, dsize, dsign, endian=raster at file@byteorder) }
+ else {
+ result <- readBin(con, what=dtype, n=ncell(raster), dsize, dsign, endian=raster at file@byteorder)
+ }
+ close(con)
+# result[is.nan(result)] <- NA
+ if (dtype == 'numeric') {
+ result[result <= (0.999 * .nodatavalue(raster)) ] <- NA
+ } else {
+ result[result == raster at file@nodatavalue ] <- NA
+ }
+ if (dtype == 'logical') {
+ result <- as.logical(result)
+ }
+ }
+ else { #use GDAL
+ if (is.na(raster at file@band)) { result <- NA }
+ else {
+ if (rownr <= 0) {
+ offs <- c(0, 0)
+ reg <- c(nrow(raster), ncol(raster)) # reg <- dim(raster at file@gdalhandle[[1]])
+ }
+ else {
+ offs= c((rownr-1), (startcol-1))
+ reg <- c(1, ncolumns)
+ }
+ }
+ result <- getRasterData(raster at file@gdalhandle[[1]], offset=offs, region.dim=reg, band = raster at file@band)
+
+# if (!is.vector(result)) { result <- as.vector(result) }
+
+ # if setNAvalue() has been used.....
+ if (raster at file@nodatavalue < 0) {
+ result[result <= raster at file@nodatavalue ] <- NA
+ } else {
+ result[result == raster at file@nodatavalue ] <- NA
+ }
+
+ }
+
+ raster at data@values <- as.vector(result)
+ if (rownr < 0) {
+ raster at data@indices <- c(1, ncell(raster))
+ raster at data@content <- "all"
+ raster <- setMinMax(raster)
+ } else if (startcol==1 & ncolumns==(ncol(raster)-startcol+1)) {
+ raster at data@indices <- c(cellFromRowCol(raster, rownr, startcol), cellFromRowCol(raster, rownr, endcol))
+ raster at data@content <- "row"
+ } else {
+ raster at data@indices <- c(cellFromRowCol(raster, rownr, startcol), cellFromRowCol(raster, rownr, endcol))
+ raster at data@content <- "block"
+ }
+
+ return(raster)
+}
Modified: pkg/raster/R/stackRead.R
===================================================================
--- pkg/raster/R/stackRead.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/stackRead.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -13,14 +13,16 @@
}
return(rstack)
}
- rstack at data@values <- matrix(nrow=length(values(raster)), ncol=nlayers(rstack))
for (i in seq(nlayers(rstack))) {
- raster <- .rasterRead(rstack at layers[[i]], rownumber, startcol, ncolumns)
- rstack at data@values[,i] <- values(raster)
+ r <- .rasterRead(rstack at layers[[i]], rownumber, startcol, ncolumns)
+ if (i==1) {
+ rstack at data@values <- matrix(nrow=length(values(r)), ncol=nlayers(rstack))
+ rstack at data@content <- dataContent(r)
+ rstack at data@indices <- dataIndices(r)
+ }
+ rstack at data@values[,i] <- values(r)
}
- rstack at data@content <- dataContent(raster)
- rstack at data@indices <- dataIndices(raster)
return(rstack)
}
Modified: pkg/raster/R/values.R
===================================================================
--- pkg/raster/R/values.R 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/R/values.R 2009-03-11 01:56:31 UTC (rev 342)
@@ -8,10 +8,14 @@
values <- function(object, format='', names=FALSE) {
format <- trim(format)
if (class(object) == 'RasterLayer') {
- if (format=='') {format <- 'vector'}
+ if (format=='') {
+ format <- 'vector'
+ }
return(.rasterValues(object, format, names=names))
} else {
- if (format=='') {format <- 'matrix'}
+ if (format=='') {
+ format <- 'matrix'
+ }
return(.stackValues(object, format, names=names))
}
}
Modified: pkg/raster/man/read.Rd
===================================================================
--- pkg/raster/man/read.Rd 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/man/read.Rd 2009-03-11 01:56:31 UTC (rev 342)
@@ -1,4 +1,4 @@
-\name{Read Raster* data}
+\name{read}
\alias{readAll,RasterLayer-method}
\alias{readAll,RasterStack-method}
@@ -15,14 +15,13 @@
\alias{readPartOfRow}
\alias{readBlock}
-
\title{Read values from disk}
\description{
Read values from a raster file associated with a RasterLayer or RasterStack object into memory.
-You can either read all the data (memory permitting), or read data by row, part of row, block, by cellnumber,
-or for (x,y) coordinates. Data can be read from 'native raster' binary files, as well as for other formats
-that are supported by the rgdal package.
+You can either read all the data (memory permitting), or read data for a single or multiple rows, part of row,
+or a block (rectangular area). Data can be read from 'native raster' binary files, as well as for other formats
+that are supported by the rgdal package \code{\link[rgdal]{gdalDrivers()}}.
}
@@ -34,6 +33,7 @@
readBlock(object, startrow, nrows=3, startcol=1, ncolumns=(ncol(object)-startcol+1))
}
+
\arguments{
\item{object}{a Raster* object}
\item{rownr}{the row number of the row to read (between 1 and nrows(raster))}
@@ -46,7 +46,7 @@
\note{After read* the values are accessible with the "values(object)" function}
-\seealso{ \code{\link[raster]{writeRaster}}, \code{\link[raster]{xyValues}} }
+\seealso{ \code{\link[raster]{xyValues}}, \code{\link[raster]{cellValues}}, \code{\link[raster]{writeRaster}}}
\author{Robert J. Hijmans}
@@ -66,6 +66,10 @@
# read a block
block <- readBlock(rs, startrow=50, nrows=5, startcol=30, ncolumns=10)
+
+# read a row from a stack
+ st <- stack(rs, rs, rs)
+ st <- readRow(st, 50)
}
\keyword{classes}
Modified: pkg/raster/man/writeRaster.Rd
===================================================================
--- pkg/raster/man/writeRaster.Rd 2009-03-10 06:45:05 UTC (rev 341)
+++ pkg/raster/man/writeRaster.Rd 2009-03-11 01:56:31 UTC (rev 342)
@@ -45,7 +45,7 @@
\author{Robert J. Hijmans}
-\seealso{ \code{\link[rgdal]{readGDAL}} }
+\seealso{ writeFormats \code{\link[raster]{writeFormats}} }
\examples{
rs <- rasterFromFile(system.file("external/test.ag", package="sp"))
More information about the Raster-commits
mailing list