# Exemplo de como criar um banco de dados
# no aRT, a partir do conjunto de dados
# parana, do pacote geoR

##
## Carregando pacotes/funcoes
##
require(geoR)
require(aRT)

##
## Exemplo 1: dados parana
##
data(parana)
points(parana, bor=borders)

## criando uma conexao com o SGBD (MySQL)
con <- openConn(u="pedro")
showDbs(con)
## opcoes para: host, port, user, password

if(any(showDbs(con)=="parana")) deleteDb(con, "parana", force=T)

## criando um novo BD
con
pr <- createDb(con, "parana")

##
## 1.1 criando um layer para armazenar os dados
##
l_dados <- createLayer(pr, "dados")

## carregando geometrias (pontos)
artcoords <- .gr2aRTpoints(parana)
addPoints(l_dados, artcoords)
t_dados = createTable(l_dados, "t_dados",gen=T)

plot(l_dados)

# inserindo atributos (so depois de fechar geometria)
artdata <- data.frame(.gr2aRTattributes(parana))

artdata[1:5,]

names(artdata)
createColumn(t_dados, "data")
updateColumns(t_dados, artdata)

## criando vistas e temas automaticas para o TV (opcional)
tema_dados <- createTheme(l_dados, "estacoes", view = "view")
setVisual(tema_dados, visualPoints(pch=22, color="black", size=5))

##
## 1.2 criando um layer para armazenar o contorno do estado
##

artpols <- .gr2aRTpolygons(parana)

l_pol<-createLayer(pr, "contorno")
addPolygons(l_pol, artpols)
createTable(l_pol, "t_pol")
plot(l_pol)

tema_pol<-createTheme(l_pol, "contorno", view="view")
setVisual(tema_pol, visualPolygons())

#
# Executando um algoritmo para análise espacial (krigagem)
#
ml <- likfit(parana, trend="1st", ini=c(1000, 100))
gx <- seq(50,800, by=10)
gy <- seq(-100,650, by=10)
loc0 <- expand.grid(gx,gy)
points(parana, bor=borders)
points(loc0, pch=".", col=2)
loc1 <- polygrid(loc0, bor=parana$bor)
points(loc1, pch="+", col=4)

KC <- krige.control(trend.d="1st", trend.l="1st", obj=ml)
kc <- krige.conv(parana, loc=loc0, krige= KC, bor=parana$borders)
save.image()

image(kc, col=terrain.colors(15))

# preparando objeto para aRT
georpred <- .prepare.graph.kriging(locations=loc0, borders=parana$borders, values=kc$pred)
artpred <- .gr2aRTraster(georpred)

l_pred <- new("aRTlayer", pr, layer="pred", create=T)

addRaster(l_pred, artpred)

# se o nome do tema for omitido, o seu nome sera o mesmo do layer
th<-createTheme(l_pred, "raster", view="view")

setVisual(th, visualRaster(), mode="r")

plot(l_pred, col=terrain.colors(15))
plot(l_dados, add=T)
plot(l_pol, add=T)