[R-sig-Geo] while loop help

[Takatsugu Kobayashi]

Hi, I must apologize both for this lengthy thread and for posting a 
thread based on the same example I have been using. 

I successfully created a code that re-arrange 10*10 grid cells of the 
same size into new grid cells of varying sizes.  Now I am attempting to 
create grid cells of the same size that are larger than the orignal 
10*10 grid cells, which means the total number of cells is a lot less 
than 100.  So If I would like to re-arrange a 10*10 grid cells by 
combining 4 cells, the number of new grid cells will be 25.

The function I created below stops with errors. It stops when a single 
cell is isolated after its cell's 1st order nearest neighbors are 
already taken.  I guess I could come up with a mathematical way of 
avoiding it, but I just wanted to loop it until no single cells are 
isolated, which means all new cells have the size equal to n.nb (the 
number of cells to be combined).

I would like to use the function below in a while loop or repeat loop, 
and when single cells are isolated I want to tell a code to start over 
the while loop from the beginning.  Because when I set the number of 
cells to be combined to 4, there are 25 sets of new grid cells (100/4), 
I know there is a solution rather than a code being trapped in an 
infinite looping.

So my question is how can I set up a code such that it starts over when 
errors occur?

Again, I am sorry for this lengthy thread. 

Thank you.

Taka

### Include necessary packages
library(sp)
library(maptools)
library(PBSmapping)


### Define parameters for a grid
start.point <- -0
cco <- c(start.point,start.point)
csize <- c(1,1)
n <- 10 #abs(cco[1]*2)+1
cnum <- c(n,n)
grd <- GridTopology(cco, csize, cnum)
r.grd <- as(grd, "SpatialPolygons")
r.grd.ps <- SpatialPolygons2PolySet(r.grd)


### Subset the grid into n*n single cells
### I use r.grd.ps over r.grd because I did not know how to clip out 
Slot "coords" from SpatialPolygons
### I tried lapply(slot(r.grd, "polygons"), function(x) slot(x, 
"coords")), but cannot obtain this slot

subPoly <- list()
coords <- list()
for (i in 1:n^2)
 {
  beg <- 5*i-4; end <- 5*i
  sP <- r.grd.ps[beg:end,]
  subPoly[[i]] <- sP
  cds <- calcCentroid(sP)
  coords[[i]] <- cbind(cds$X, cds$Y)
 }


### 1st order continuity neighbors
knn.nb <- list()
res <- vector("numeric",length=n^2)

for (ki in 1:n^2)
 for (ji in 1:n^2)
  {{
   ## Find the cells whose boundaries share with jth cell
   res[ji] <- 
length(which(subPoly[[ki]][1:4,4]%in%subPoly[[ji]][1:4,4]&subPoly[[ki]][1:4,5]%in%subPoly[[ji]][1:4,5]))  

   res[ki] <- 0
   ## Identify which cells are k-nearest neighbors excluding the 
reference zone
   knn.nb[[ki]] <- which(res>1)  
  }}


### Clip individual polygons and merge some of them based on adjacency 
weights
## Fixed number of cells to be combined excluding the jth cell
n.nb <- 3

maup <- function (knn.nb, n.nb)
{
## j for subPoly ID: input
j <- 1
## k for randomPoly ID: output
k <- 1
## A vector of subPoly IDs
sample.cell <- 1:n^2
## A list in which to store combined Polygons
randomPoly <- list()
## A vector of cell areas: for making sure all cells have the size of 
n.nb =4
areaPoly <-numeric(n^2/n.nb)

  while (all(areaPoly != n.nb))
        {
         ## Retrieve at most 5 contiguous neighbors
         j <- min(which(!is.nan(sample.cell)))
 
         ## Pick one nearest neighbor of jth cell randomly, including 
the jth cell
         d1.id <- which(!is.nan(knn.nb[[j]]))
         d1 <- knn.nb[[j]][d1.id]
         d1.rand <- round(runif(1,0.5,length(d1)))
         d1 <-  d1[d1.rand]
         nb.list <- append(j, d1)

         for (k.nb in 2:n.nb)
             {
              ## Pick one nearest neighbor of j+1th cell randomly
              ## Repeat till n.nb and append every cells chosen
              ## d as nearest neighbors
              d_id <- nb.list[k.nb]
              d_k_id <- which(!is.nan(knn.nb[[d_id]]))
              d_k <- knn.nb[[d_id]][d_k_id]
              non_overlap_id <- which(d_k%in%nb.list)
              d_k <- d_k[-non_overlap_id]
              dk.rand <- round(runif(1,0.5,length(d_k)))
              dk <- d_k[dk.rand]
              nb.list <- append(nb.list, dk)
             }

        ## Remove cells that are already taken
        sample.cell[nb.list] <- NaN
   
        for (ord in 1: n^2)
            {
             id.check <- which(knn.nb[[ord]]%in%nb.list)
             if (length(id.check)>0)
                {  
                knn.nb[[ord]][id.check] <- NaN
                }   
            }


        ## JoinPolys allows to join two of the same polysets
        joinedPoly <- joinPolys(subPoly[[j]],subPoly[[nb.list[2]]],"UNION")
   
        for (pl in 3:(n.nb+1)) # +1 means the jth cell
            {      
             joinedPoly <- 
joinPolys(joinedPoly,subPoly[[nb.list[pl]]],"UNION")       
            }

        ## PID should be single and unique
        randomPoly[[k]] <- combinePolys(joinedPoly)
        randomPoly[[k]]$PID <- k

      areaPoly[k] <- calcArea(randomPoly[[k]])[,3]
 
      k <- k+1
      }
}