## part2 is the input object ## initclustsize <- 80 ## clustsizelimit <- 500 ## cmethod="max" library(geosphere) distmat <- distm (part2[,c("Longitude","Latitude")], fun = distHaversine)/ 1609 N <- nrow(distmat) coverview <- data.frame(ID=1:N, cluster = 1:N) curclust <- unique(coverview$cluster) kcluster <- length(curclust) ## initial clusters for(k in 1:N){ if(k %in% curclust){ coverview[which(distmat[k,]clust1] for(clust2 in clust2candi){ disttmp2 <- disttmp[,coverview$cluster==clust2] steptmp <- data.frame(cluster1=clust1, cluster2=clust2, distmin=min(disttmp2), distmax=max(disttmp2)) stepdat <- rbind(stepdat,steptmp) } } ## FULL WHILE GROUP while(nrow(stepdat)>0){ if(cmethod == "max"){ stepdatmin <- stepdat[which.min(stepdat$distmax),] }else{ stepdatmin <- stepdat[which.min(stepdat$distmin),] } tclust <- stepdatmin$cluster1 sclust <- stepdatmin$cluster2 coverview[coverview[,2]==sclust,2] <- tclust stepdat <- stepdat[!(stepdat$cluster1 == sclust | stepdat$cluster2 == sclust),] updatestep <- (stepdat$cluster1==tclust | stepdat$cluster2 == tclust) for(v in which(updatestep)){ clust1 <- stepdat$cluster1[v] clust2 <- stepdat$cluster2[v] disttmp2 <- distmat[coverview$cluster==clust1,coverview$cluster==clust2] steptmp <- data.frame(cluster1=clust1, cluster2=clust2, distmin=min(disttmp2), distmax=max(disttmp2)) stepdat[v,] <- steptmp } stepdat <- stepdat[stepdat$distmax<=clustsizelimit,] print(length(unique(coverview$cluster))) }