[Raster-commits] r286 - pkg/raster/R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Feb 15 15:48:33 CET 2009
Author: jacobvanetten
Date: 2009-02-15 15:48:33 +0100 (Sun, 15 Feb 2009)
New Revision: 286
Modified:
pkg/raster/R/distance.R
Log:
distance disk-to-disk
syntax correct but not yet tested
Modified: pkg/raster/R/distance.R
===================================================================
--- pkg/raster/R/distance.R 2009-02-14 12:29:14 UTC (rev 285)
+++ pkg/raster/R/distance.R 2009-02-15 14:48:33 UTC (rev 286)
@@ -4,45 +4,116 @@
# Version 1.0
# Licence GPL v3
-setGeneric("distance", function(object, ...) standardGeneric("distance"))
+#setGeneric("distance", function(object, ...) standardGeneric("distance"))
-setMethod("distance", signature(object = "RasterLayer"), def = function(object) {
+#setMethod("distance", signature(object = "RasterLayer"), def =
+
+distance <- function(object, filename="") {
n <- ncell(object)
- fromCells <- which(!is.na(values(object)))
- toCells <- which(is.na(values(object)))
- accDist <- rep(0,times=n)
- accDist[toCells] <- Inf
- if (isLatLon(object)) {
- while(length(fromCells)>0)
- {
- adj <- adjacency(object,fromCells=fromCells,toCells=toCells,directions=8)
- coord <- cbind(xyFromCell(object,adj[,1]),xyFromCell(object,adj[,2]))
- distance <- apply(coord,1,function(x){distanceGreatcircle(x[1:2],x[3:4])})
-
- transitionValues <- accDist[adj[,1]] + distance
- transitionValues <- tapply(transitionValues,adj[,2],min)
- transitionValues <- transitionValues[transitionValues < Inf]
- index <- as.integer(names(transitionValues))
- fromCells <- index[transitionValues < accDist[index]]
- accDist[index] <- pmin(transitionValues,accDist[index])
+
+ if(dataContent(object)=='all'){
+ fromCells <- which(!is.na(values(object)))
+ toCells <- which(is.na(values(object)))
+ accDist <- rep(0,times=n)
+ accDist[toCells] <- Inf
+ if (isLatLon(object)) {
+ while(length(fromCells)>0)
+ {
+ adj <- adjacency(object,fromCells=fromCells,toCells=toCells,directions=8)
+ coord <- cbind(xyFromCell(object,adj[,1]),xyFromCell(object,adj[,2]))
+ distance <- apply(coord,1,function(x){distanceGreatcircle(x[1:2],x[3:4])})
+ #What follows is the same as for non-projected (below)
+ transitionValues <- accDist[adj[,1]] + distance
+ transitionValues <- tapply(transitionValues,adj[,2],min)
+ transitionValues <- transitionValues[transitionValues < Inf]
+ index <- as.integer(names(transitionValues))
+ fromCells <- index[transitionValues < accDist[index]]
+ accDist[index] <- pmin(transitionValues,accDist[index])
+ }
+ } else {
+ while(length(fromCells)>0) {
+ adj1 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions=4)
+ adj2 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions="Bishop")
+ distance <- c(rep(1,length=length(adj1[,1])),rep(sqrt(2),length=length(adj2[,1])))
+ adj <- rbind(adj1,adj2)
+ #What follows is the same as for LatLon
+ transitionValues <- accDist[adj[,1]] + distance
+ transitionValues <- tapply(transitionValues,adj[,2],min)
+ transitionValues <- transitionValues[transitionValues < Inf]
+ index <- as.integer(names(transitionValues))
+ fromCells <- index[transitionValues < accDist[index]]
+ accDist[index] <- pmin(transitionValues,accDist[index])
+ }
}
- } else {
- while(length(fromCells)>0) {
- adj1 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions=4)
- adj2 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions="Bishop")
- distance <- c(rep(1,length=length(adj1[,1])),rep(sqrt(2),length=length(adj2[,1])))
- adj <- rbind(adj1,adj2)
-
- transitionValues <- accDist[adj[,1]] + distance
- transitionValues <- tapply(transitionValues,adj[,2],min)
- transitionValues <- transitionValues[transitionValues < Inf]
- index <- as.integer(names(transitionValues))
- fromCells <- index[transitionValues < accDist[index]]
- accDist[index] <- pmin(transitionValues,accDist[index])
+ outRaster <- object
+ outRaster <- setValues(outRaster, accDist)
+ return(outRaster)
+ }
+ if(dataContent(object)=='nodata' & dataSource(object) =='disk'){ #to be tested
+ nrows <- nrow(object)
+ ncols <- ncol(object)
+ outRaster <- setRaster(object, filename)
+ for(r in 1:nrows)
+ {
+ rowValues <- values(readRow(object, rownr = r))
+ outRowValues <- rep(Inf,times=ncols)
+ outRowValues[is.na(rowValues)] <- 0
+ outRaster <- setValues(outRaster, outRowValues, r)
+ outRaster <- writeRaster(outRaster, overwrite=TRUE)
}
+ if(isLatLon(object)){
+ remainingCells <- TRUE
+ while(remainingCells){
+ remainingCells <- FALSE
+ oldRowValues <- integer(0)
+ rowWindow <- values(readRow(outRaster, rownr=1))
+ for(r in 1:nrows){
+ if(r<nrows){rowWindow <- c(rowWindow,values(readRow(outRaster, rownr=r+1)))}
+ adj <- adjacency(fromCells=(((max(1,r-1))*ncols)+1):(min(nrows,(r+2)*ncols)), toCells=((r-1)*ncols+1):(r*ncols),directions=8)
+ coord <- cbind(xyFromCell(object,adj[,1]),xyFromCell(object,adj[,2]))
+ distance <- apply(coord,1,function(x){distanceGreatcircle(x[1:2],x[3:4])})
+ adj <- adj-((r-1)*ncols+1)
+ transitionValues <- as.vector(rowWindow)[adj[,1]] + distance
+ transitionValues <- tapply(transitionValues,adj[,2],min)
+ transitionValues <- transitionValues[transitionValues < Inf]
+ index <- as.integer(names(transitionValues))
+ newValues <- pmin(transitionValues,rowWindow[index])
+ if(newValues != rowWindow[index]){remainingCells<-TRUE}
+ rowWindow[index]
+ outRaster <- setValues(outRaster, rowValues, startRow)
+ outRaster <- writeRaster(outRaster, overwrite=TRUE)
+ if(r>1){rowWindow <- rowWindow[-1:ncols]}
+ }
+ }
+ }
+ else{
+ remainingCells <- TRUE
+ while(remainingCells){
+ remainingCells <- FALSE
+ oldRowValues <- integer(0)
+ rowWindow <- values(readRow(outRaster, rownr=1))
+ for(r in 1:nrows){
+ if(r<nrows){rowWindow <- c(rowWindow,values(readRow(outRaster, rownr=r+1)))}
+ fromCells <- (((max(1,r-1))*ncols)+1):(min(nrows,(r+2)*ncols))
+ toCells <- ((r-1)*ncols+1):(r*ncols)
+ adj1 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions=4)
+ adj2 <- adjacency(object,fromCells=fromCells,toCells=toCells,directions="Bishop")
+ distance <- c(rep(1,length=length(adj1[,1])),rep(sqrt(2),length=length(adj2[,1])))
+ adj <- rbind(adj1,adj2)
+ adj <- adj-((r-1)*ncols+1)
+ transitionValues <- as.vector(rowWindow)[adj[,1]] + distance
+ transitionValues <- tapply(transitionValues,adj[,2],min)
+ transitionValues <- transitionValues[transitionValues < Inf]
+ index <- as.integer(names(transitionValues))
+ newValues <- pmin(transitionValues,rowWindow[index])
+ if(newValues != rowWindow[index]){remainingCells<-TRUE}
+ rowWindow[index]
+ outRaster <- setValues(outRaster, rowValues, startRow)
+ outRaster <- writeRaster(outRaster, overwrite=TRUE)
+ if(r>1){rowWindow <- rowWindow[-1:ncols]}
+ }
+ }
+ }
}
- result <- object
- result <- setValues(result, accDist)
- return(result)
}
-)
+#)
\ No newline at end of file
More information about the Raster-commits
mailing list