[Raster-commits] r160 - pkg/raster/R

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Jan 19 18:31:08 CET 2009


Author: rhijmans
Date: 2009-01-19 18:31:06 +0100 (Mon, 19 Jan 2009)
New Revision: 160

Modified:
   pkg/raster/R/exportGDAL.R
   pkg/raster/R/group.generic.functions.R
Log:


Modified: pkg/raster/R/exportGDAL.R
===================================================================
--- pkg/raster/R/exportGDAL.R	2009-01-19 13:29:22 UTC (rev 159)
+++ pkg/raster/R/exportGDAL.R	2009-01-19 17:31:06 UTC (rev 160)
@@ -18,13 +18,11 @@
 }
 
 
-exportGDAL <- function(raster, filename, gdalfiletype = "GTiff", overwrite=FALSE, ForceIntOutput = FALSE ) {
-	mvFlag = NA
-	options = NULL
+.getGDALtransient <- function(raster, filename, gdalfiletype, mvFlag, options, overwrite, ForceIntOutput)  {
 	.isSupportedGDALdriver(gdalfiletype)
 	
-# this is a RasterLayer hence:
-    nbands = 1
+# this is a RasterLayer hence nbands = 1:
+    nbands = nlayers(raster)
 # but we keep this for later (stack, brick)
 
 	if (file.exists(filename)) {
@@ -56,7 +54,18 @@
     p4s <- projection(raster)
     .Call("RGDAL_SetProject", transient, p4s, PACKAGE = "rgdal")
 	
+	return(transient)
+}
+
+
+
+exportGDAL <- function(raster, filename, gdalfiletype = "GTiff", overwrite=FALSE, ForceIntOutput = FALSE ) {
+	mvFlag = NA
+	options = NULL
+	nbands = nlayers(raster)
 	
+	transient <- .getGDALtransient(raster, filename, gdalfiletype, mvFlag, options, overwrite, ForceIntOutput)
+	
     for (band in 1:nbands) {
 
 		if (dataContent(raster)=='all') {
@@ -65,15 +74,15 @@
 #			x <- putRasterData(transient, t(values(raster, format='matrix')), band, c(0, 0)) 
 			for (r in 1:nrow(raster)) {
 				x <- putRasterData(transient, valuesRow(raster, r), band, c((r-1), 0)) 
-			}	
+			}
 		} else {
 			if (dataSource(raster)=='ram') {
 				stop("No data on disk, and not all values in memory. Cannot write the file")
 			}
 			for (r in 1:nrow(raster)) {
 				x <- putRasterData(transient, values(readRow(raster, r)), band, c((r-1), 0)) 
-			}	
-		}	
+			}
+		}
 
 #        if (!is.na(mvFlag)) {
 #            transient_b <- getRasterBand(dataset = transient, band = band)
@@ -94,3 +103,36 @@
 	return(outras)
 }
 
+
+writeGDALrow <- function(raster, filename, gdalfiletype = "GTiff", rownr, overwrite=FALSE, ForceIntOutput = FALSE ) {
+	if (rownr == 1) {
+		mvFlag = NA
+		options = NULL
+		transient <- .getGDALtransient(raster, filename, gdalfiletype, mvFlag, options, overwrite, ForceIntOutput)
+		attr(raster, "transient") <- transient
+	}	
+	
+    for (band in 1:nlayers(raster)) {
+		x <- putRasterData(raster at transient, values(raster, rownr), band, c((rownr-1), 0)) 
+	}	
+
+#        if (!is.na(mvFlag)) {
+#            transient_b <- getRasterBand(dataset = transient, band = band)
+#            .Call("RGDAL_SetNoDataValue", transient_b, as.double(mvFlag), PACKAGE = "rgdal")
+#       }
+    
+	if (rownr == nrow(raster)) {
+		saveDataset(raster at transient, filename)
+		GDAL.close(raster at transient) 
+#  do NOT do this, it removes the driver for future use!!! ????	
+#	GDAL.close(driverobj) 
+		outras <- rasterFromFile(filename)
+		outras at data@min <- raster at data@min
+		outras at data@max <- raster at data@max
+		if (!is.na((outras at data@min))) { outras at data@haveminmax <- TRUE }
+		return(outras)
+	} else {
+		return(raster)
+	}
+}
+

Modified: pkg/raster/R/group.generic.functions.R
===================================================================
--- pkg/raster/R/group.generic.functions.R	2009-01-19 13:29:22 UTC (rev 159)
+++ pkg/raster/R/group.generic.functions.R	2009-01-19 17:31:06 UTC (rev 160)
@@ -12,6 +12,7 @@
 .CanProcessInMemory <- function(raster, n=2, datasize=16) {
 #	memalloc <- n * ncell(raster) * 8
 #	memavailable <- memory.limit()-memory.size()
+# For now something simplistic :
 	maxalloc <- 10^8
 	if ( (ncell(raster) * n * datasize) > maxalloc ) { 
 		return( FALSE )
@@ -100,9 +101,19 @@
 		if (!isTRUE(is.atomic(e2) & length(e2)==1)) {
 			stop('second argument should be a single number')
 		}
-		rs <- setRaster(e1, values=callGeneric(.getRasterValues(e1), rep(e2, ncell(e1)) ) )
-		rs <- setDatatype(rs, datatype='integer', datasize=2)
-		return(rs)
+		if (.CanProcessInMemory(e1, 2)) {
+			raster <- setRaster(e1, values=callGeneric(.getRasterValues(e1), rep(e2, ncell(e1)) ) )
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		} else {
+			raster <- setRaster(e1, filename=tempfile())
+			rowrep <- rep(e2, ncol(e1))
+			for (r in 1:nrow(e1)) {
+				raster <- setValues(raster, callGeneric( .getRowValues(e1, r), rowrep ), r)
+				raster <- writeRaster(raster)
+			}
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		}
+		return(raster)
 	}
 )	
 
@@ -111,9 +122,19 @@
 		if (!isTRUE(is.atomic(e2) & length(e2)==1)) {
 			stop('first argument should be a single number')
 		}
-		rs <- setRaster(e2, values=callGeneric(.getRasterValues(e2), rep(e1, ncell(e2)) ) ) 
-		rs <- setDatatype(rs, datatype='integer', datasize=2)
-		return(rs)
+		if (.CanProcessInMemory(e2, 2)) {
+			raster <- setRaster(e2, values=callGeneric(.getRasterValues(e2), rep(e1, ncell(e2)) ) )
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		} else {
+			raster <- setRaster(e2, filename=tempfile())
+			rowrep <- rep(e1, ncol(e2))
+			for (r in 1:nrow(e2)) {
+				raster <- setValues(raster, callGeneric( .getRowValues(e2, r), rowrep ), r)
+				raster <- writeRaster(raster)
+			}
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		}
+		return(raster)
 	}
 )	
 
@@ -123,9 +144,18 @@
 		if (!cond) {
 			stop("Cannot compare RasterLayers that have different BasicRaster attributes. See compare()")
 		}	
-		rs <- setRaster(e1, values=callGeneric(.getRasterValues(e1), .getRasterValues(e2) ) ) 
-		rs <- setDatatype(rs, datatype='integer', datasize=2)
-		return(rs)
+		if (.CanProcessInMemory(e1, 2)) {
+			raster <- setRaster(e1, values=callGeneric(.getRasterValues(e1), .getRasterValues(e2) ) ) 
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		} else {
+			raster <- setRaster(e1, filename=tempfile())
+			for (r in 1:nrow(e1)) {
+				raster <- setValues(raster, callGeneric( .getRowValues(e1, r), .getRowValues(e2, r) ), r)
+				raster <- writeRaster(raster)
+			}
+			raster <- setDatatype(raster, datatype='integer', datasize=2)
+		}
+		return(raster)
 	}
 )	
 
@@ -196,19 +226,39 @@
 	}
 )
 
+
 setMethod("Arith", signature(e1='RasterLayer', e2='numeric'),
     function(e1, e2){ 
-		return(setRaster(e1, values=callGeneric(.getRasterValues(e1), e2)))
+		if (.CanProcessInMemory(e1, 2)) {
+			return(setRaster(e1, values=callGeneric(.getRasterValues(e1), e2)))
+		} else {
+			raster <- setRaster(e1, filename=tempfile())
+			for (r in 1:nrow(e1)) {
+				raster <- setValues(raster, callGeneric( .getRowValues(e1, r), e2) , r)
+				raster <- writeRaster(raster)
+			}
+			return(raster)
+		}		
 	}
 )
 
 setMethod("Arith", signature(e1='numeric', e2='RasterLayer'),
     function(e1, e2){ 
-		return(setRaster(e2, values=callGeneric(.getRasterValues(e2), e1)))
+		if (.CanProcessInMemory(e2, 2)) {
+			return(setRaster(e2, values=callGeneric(.getRasterValues(e2), e1)))
+		} else {
+			raster <- setRaster(e2, filename=tempfile())
+			for (r in 1:nrow(e2)) {
+				raster <- setValues(raster, callGeneric( .getRowValues(e2, r), e1) , r)
+				raster <- writeRaster(raster)
+			}
+			return(raster)
+		}		
 	}
 )
 
 
+
 setMethod("max", signature(x='Raster'),
 	function(x, ..., na.rm=FALSE){
 		obs <- list(...)



More information about the Raster-commits mailing list