Exemplo: combinando dados de municípios e mesoregiões

Script para construir um banco com as cidades e as mesoregioes do parana a partir de um shapefile e um csv (arquivos já incluídos no pacote aRT). São criados vários temas com diferentes opções de visualização/fatiamento.

require(aRT)
con=openConn()
con
 
if(any(showDbs(con)=="meso"))
  deleteDb(con, "meso", force=T)

criando o banco

db=createDb(con,db="meso")

importando o shapefile e informando sobre projecao

proj="+proj=poly" +lon_0=0.0 +lat_0=0.0 +x_0=0.0 +y_0=0.0 +a=6378160.0 +f=0.003353"

criando layer a partir do shapefile (diretamente s/ impotar p/ R)

l=createLayer(db, l="cidades", proj=proj)
getProj(l)
addShape(l, tab="tab", file=system.file("meso/parana_pol.shp",package="aRT"), id="CODIGO")

visualizando diretamente do banco

getPolygons(l)->lpol
aRTplot(lpol)

abre a tabela do shapefile e importa dados para o R

tab = openTable(l,"tab")
tdata <- getData(tab)

e calcula as meso-regioes

mesogeom <- list()
mesonomes <- unique(tdata$CODIREGIAO)
 
for(i in 1:length(mesonomes))
{
  selecionando municipios dentro da regiao i
  id=tdata$CODIGO[which(tdata$CODIREGIAO == mesonomes[i])]
  cat(paste("Processando mesoregiao",mesonomes[i], "com",length(id), "cidades\n"))
  unindo os poligonos dos municipios para formar pol. da regiao
  set=getSetOperation(l,"union",id=id)
  set@polygons[[1]]@ID = paste(mesonomes[i])
  adicionando ao objeto (list) das regioes
  mesogeom[i] = set@polygons
  aRTplot(set, col=terrain.colors(length(mesonomes)+1)[i+1], add=T, lwd=2)
}
 
aRTplot(lpol, add=T)
res = SpatialPolygons(mesogeom, 1:length(mesogeom))

insere os dados das meso-regioes, e cria temas para visualizacao

l2=createLayer(db,l="mesoreg", proj=proj)
addPolygons(l2, res)

importa a tabela estatica de um csv

mr=read.csv(system.file("meso/parana_regions.csv",package="aRT"), sep=";")
tabmeso=importTable(l2, "tabmeso", id="CODIGO", data=mr)

##########NOVO ################################ carrega os dados da tabela no R

data=getData(tabmeso)

cria um dataframe com duas colunas: o CODIGO e a soma das producoes de QPM entre 2000 e 2003

somaQPM=data.frame(CODIGO=data$CODIGO,
                   somaQPM=as.integer(data$QPM_2000)+
				           as.integer(data$QPM_2001)+
						   as.integer(data$QPM_2002)+
						   as.integer(data$QPM_2003))

as colunas do banco. Para cada coluna que nao existe ele cria antes de atualizar, e o codigo é utilizado para manter as referencias corretas

updateColumns(tabmeso, somaQPM)

#######FIM NOVO ################################

nova classe em desenvolvimento para gerar temas coloridos para TV

v = visualPolygons(color="blue", transp=80, ccolor="black")
t = createTheme(l2, "mesoregioes",view="view")
setVisual(t, v)

modificando a visualizacao das cores (vai ser encapsulado em funcoes)

v@ccolor = col2rgb("#F0F0F0")
v@transp=as.integer(50)
t = createTheme(l, "cidades",view="view")
setVisual(t, v)
 
v1=visualRaster()
aRTsilent(TRUE)

criando temas com as atributos dos minicipios

#colnames = colnames(tdata)
#for(i in 1:length(colnames))
#{
#	if(!(any(colnames[i] == c("CODIGO", "AREA_1", "PERIMETRO_", "NOME", "CODIREGIAO"))))
#	{
#		cat(paste("Creating theme", colnames[i], "\n"))
#		t=createTheme(l, th=colnames[i], view = "view")
#		setVisual(t, v1, colnames[i])
#	}
#}

#######NOVO ################################ cria os temas para os cinco atributos, quatro producoes e a soma, coloridos. Estes temas sao criados a partir do layer de mesoregioes

#for(i in c("QPM_2000", "QPM_2001","QPM_2002","QPM_2003","somaQPM"))
#{
 cat(paste("Creating theme meso", i, "\n", sep=""))
 t=createTheme(l2, th=paste("meso", i, sep=""), view = "view")
 setVisual(t, v1, att=i)
#}

##########FIM NOVO ################################

for(i in c("e","q","s"))
{
  cat(paste("Creating theme ", "QPM_2003", i, "\n", sep=""))
  t=createTheme(l, paste("cidadesQPM_2003", i, sep=""), view = "view")
    #mode abaixo define o criterio de criacao da legend (padrao TL)
  setVisual(t, v1, att="QPM_2003", mode=i)
 setVisible(t, FALSE)
}

temas com atributos das mesoregioes

#colnames = colnames(mr)
#for(i in 1:length(colnames))
#{
#	cat(paste("Creating theme", colnames[i], "\n"))
#	if(!(any(colnames[i] == c("NOME", "CODIGO"))))
#	{
#		t=createTheme(l2, th=colnames[i], view = "view")
#		setVisual(t, v1, colnames[i])
#	}
#}

operacoes em R

#munic <- getGeometry(l)
#is.list(munic)
#names(munic)
#is.list(munic$polygons)
#length(munic$polygons)
#plot(munic)