ラスターデータのK-meansクラスタリング

試行錯誤の上ようやくできました。
GRASS+Rです。
参考：R-Forge: Exiting with error

library(raster)
library(spgrass6)
library(cluster)

#import from GRASS
d=readRAST6(c(“A2002161_NDVI”),plugin=FALSE)
#convert to raster
rast=raster(d)
#convert ot rasterstack (これをやらないとgetvaluesができなかった)
rst.stack=stack(rast)
# get values from raster stack
v=getValues(rst.stack)
# NA
i <- which(! is.na(v[,1]))
#clustering (k=3)
clus=clara(v[i,],3)
#write the result to raster
#空のラスターを作成
r=raster(rast)
#上の結果をrに代入
r[i]=clus$clustering
#plot
plot(r)

#write to GRASS
spr=as(r,”SpatialGridDataFrame”)
writeRAST6(spr,”KA2002169_NDVI”,zcol=1,useGDAL=TRUE)