[R-sig-Geo] while loop help: SOLVED

[Takatsugu Kobayashi]

I figured it out. Thanks for taking a look. thanks

Takatsugu Kobayashi wrote:
> 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
>       }
> }
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>