quarta-feira, 26 de agosto de 2009

Sql Server 2008 e o sistema de uma projeção só.

Boa noite pessoal,
Estava trabalhando em alguns assuntos no Sql Server 2008 e até aí tudo bem. Tinha alguns dados em outro banco de dados, representativos de pontos. Tudo UTM SAD69. O banco do GIS em Lat/Long SAD69.


As tabelas, eram mais ou menos assim:
ID
X
Y
Z
Sistema_coordenada


Sempre com diversos pontos e diversos sistemas de coordenadas "misturados". Até aí tudo bem.


Bem, a idéia era montar uma visão no banco de dados, registar a parada com o ArcSDE, e disponibilizar os dados que eu queria em tempo real para os usuários do GIS, inclusive eu.


Liguei em uma certa empresa que comercializa um certo software, no suporte. Expliquei a situação, perguntei se o rapaz ou alguém lá saberia de uma solução para realizar o idealizado. De forma alguma (eu) estava viajando na maionese, pela minha experiência (que não é muita) com o PostGIS, sabia que era fazível.


O certo indivíduo que me atendeu, deve ter feito uma cara do outro lado de "What the hell?" sugeriu que eu jogasse minha tabela no ArcMap e rodasse um XY para plotar os dados. Valeu pela super dica amigão. Se você ler isto e lembrar de mim sabe que existe uma forma!


Mas vamos lá. O PostGIS, que é um tanto mais avançado que o Sql Server 2008 nesses quesitos, bastaria o seguinte para criar uma visão espacial (com um ponto de verdade, e não um par de coordenadas) e transformar a parada para o datum que eu queria:

CREATE OR REPLACE VIEW view_foo AS
SELECT ST_Transform(ST_Point(X,Y),4618), coluna1, coluna2 from foo;


Bem, simples né? Acontece que o SQL Server 2008 NÃO POSSUI funções nativas para conversão de coordenadas/projeção de um sistema para outro.


Um camarada do trabalho, Danilo, me falou que já havia visto no CodePlex, um site com um repositório de código da Microsoft, um projeto que realizava o serviço. Bem, o projeto, Sql Spatial Tools, tem muita utilidade e outras funções, como agregador de geometrias, mas as funções de conversão de coordenadas estava furada. Pelo menos a que nós precisávamos usar.


Lembrei que o PostGIS, que também não tem esta capacidade nativa, mas é fornecida pela Proj4, e talvez poderíamos utilizar a Proj4. Mais uma busca no Codeplex e achamos a ProjNET, uma adaptação da Proj4 (acho que é escrita em C) para .NET, em C#.


Utilizamos o modelo da Sql Spatial Tools para registrar bibliotecas externas (muito top, parabéns Sql Server, muito tetinha de fazer) e para expor uma única função que criamos para o Sql Server. A ProjNET, ao contrário da Sql Spatial Tools, realmente converte as coordenadas.


Dicas
  1. Sempre olhe no Codeplex.
  2. Para registrar uma dll externa (dentro do Sql Server) utilize (isto é código SQL): CREATE ASSEMBLY ProjNet 'caminhoParaDll'
  3. Para registrar uma função: CREATE FUNCTION nomeFuncao (parametros) returns tipoRetorno AS external ProjNet.Namespace.Funcao
  4. É só.
No final nossa idéia ficou similar a do PostGIS, and it works. Resultado: utilizamos os dados que temos, sem a necessidade de administrar/atualizar FeatureClasses e realizar tarefas tediosas de Add XY data-Project Data e propensas à erros. E o melhor de tudo, é tudo em tempo real. Se o outro sistema for atualizado às 10:01:01:0001, nós iremos ver os dados atualizados às 10:01:01:0002.


Lições aprendidas
  • Novamente, sempre olhe no CodePlex.
  • É muito fácil desenvolver/linkar APIs e outros frameworks ao Sql Server 2008.
  • Nunca, nunca confie no suporte técnico oferecido.


Thank you notes/Agradecimentos
Danilão! Pelo tempo gasto nisso e paciência.
Onald McDonald e Ogerio McDonald, pela torcida.
Ed Katibah, lead developer at Microsoft, who exchanged several emails with me and noticed that their implementation of Sql Spatial Tools is not correct (for UTM projections), and opened a ticket for fixing. Soon Sql Spatial Tools will be ready to handle all conversions.

domingo, 23 de agosto de 2009

Geoprocessamento e acidentes de trânsito

Olá novamente pessoal,
Venho falar com vocês de um assunto sério: acidentes de trânsito. Matam milhares de pessoas, ferem outras milhares e custam aos cofres públicos milhares de reais anualmente. Claro que o dinheiro neste caso não é o denominador comum, mas imaginem:
De acordo com o IPEA (Instituto de Pesquisas Econômicas Aplicadas) cada acidente de trânsito custa, em média, por vítima envolvida:
  • vítima ilesa: R$1.040,00
  • vítima ferida: R$36.305,00
  • vítima fatal: R$270.165,00
São valores absurdos. Isso, de acordo com o IPEA, envolve não só custos materiais e despesas hospitalares, mas também perda de produção, ou seja, dias parados no trabalho.
Bem, e o que Geoprocessamento ou SIG tem à ver com isso? Bem, um acidente é um evento extremamente complexo, com diversas variáveis importantes para uma macro-análise do trânsito e da segurança viária. Então, como conhecer um acidente de trânsito e suas variáveis? Como analisar a frequência de acidentes e as condições de tráfego de cada via? Podemos correlacionar os acidentes com mau tempo?
Um sistema que cadastre e analise acidentes é o ideal. Mas como localizar cada acidente? Temos à nossa disposição em quase todos os softwares de SIG/GIS ferramentas que executam a geocodificação. A geocodificação é basicamente o processo de se assinalar uma localização pontual (um ponto) à uma descrição textual, geralmente um endereço. Podemos utilizar ferramentas de geocodificação para localizar acidentes? Afirmativo. E focos de dengue? Afirmativo. E...qualquer coisa que possua um endereço em seu cadastro.
Um sistema deste tipo pode ajudar a salvar vidas e a diminuir o gasto público com obras ou medidas de segurança desnecessárias. Um sistema que realize este trabalho pode dar aos responsáveis pelo planejamento viário, subsídio técnico-científico para a proposição de medidas concretas.
Durante minha monografia, tentei desenvolver um sistema que fizesse este tipo de trabalho. Utilizei para o trabalho o banco de dados PostgreSQL e PostGIS e um algoritmo de geocodificação desenvolvido por David Bitner (não tenho o endereço do algoritmo, mas se quiserem dar uma olhada só me pedir por email). Funcionou?
Funcionou, mas existem alguns percalços no caminho: a polícia militar e outros órgãos que registram os acidentes de trânsito em boletins de ocorrência, algumas vezes registram somente a via em que o acidente ocorreu e as vias que interseccionam aquele trecho (Rua X, entre rua A e Avenida B) ou somente o cruzamento em que o mesmo ocorreu. Ferramentas de geocodificação tradicionais não suportam este tipo de endereçamento, portanto utilizei um pouquinho de código pl/pgsql para adicionar esta funcionalidade. Seria um sistema muito ruim se mais da metade dos acidentes cadastrados não fossem possíveis de serem geocodificados. Qualquer análise espacial realizada seria furada, pois os dados não estão completos.
Certo, mas o que eu preciso para montar algo similar? Ok, vamos lá.
  1. Uma base de referência de logradouros, com os seguintes atributos (ou similares): ID_Trecho, ID_Logradouro, Tipo_Logradouro, Nome_Logradouro, Numero_Inicial_Esquerdo, Numero_Inicial_Direito, Numero_Final_Esquerdo e Numero_Final_Direito, Intersecao_Anterior e Intersecao_Posterior. A base deve conter os trechos de logradouros devidamente preenchidos. Você pode usar um Guia Sei ou IR A CAMPO, e anotar estas numerações. Claro que para um projeto maior isto não funcionaria, mas é um começo. Nos campos Numeracao*, devem constar as numerações iniciais/finais de cada lado da via e em Intersecao* os logradouros que interseccionam aquele trecho.
  2. O algortimo geocodificador e um serviço geocodificador. Tudo isto vem com o algoritmo. O algoritmo é um pedaço de código que é orientado pelo serviço. O serviço é criado preenchendo-se a tabela correta, onde indicaremos qual é a base de referência e quais são os campos de interesse.
  3. Não tem três. É só isso.
 Bem cada acidente tem suas características, um número x de veículos envolvidos e outras informações relativas à cada motorista e gravidade de ferimentos (por veículo). A maioria dos algoritmos geocodificadores possuem ferramentas para geocodificação ad-hoc (propósito em si) e podemos desenvolver maneiras para geocodificador dados em batch (ou fornadas, levas, etc.).
No caso deste trabalho, foi desenvolvido um trigger ou gatilho, disparado assim que um usuário inserir alguma coisa na tabela que registra os acidentes, realizando uma busca pela representação pontual daquele endereço. Este trigger também avalia o tipo de endereço inserido pelo usuário, se ele representava um endereço completo (Rua A, 768) ou se representava um cruzamento (Rua A com Rua B) ou se apenas continha os valores das interseções (Rua A, entre Rua B e Rua C).
As funções complementares, que geocodificam entre intersecções ou em cruzamentos são apenas o uso simples de SQLs comuns e da função nativa do PostGIS ST_Line_Interpolate_Point, que recebe como argumentos a linha que se deseja interpolar a uma porcentagem.
Certo, como localizar um cruzamento? Bem, lembra que nossa base de referência possui o atributo Sufixo_Direcao e que este representa para qual linha "estamos indo"? Então, para selecionar a linha e mandar ela para a função, foi utilizando o seguinte:
SELECT * FROM logradouros WHERE Tipo_Logradouro = 'TipoEscolhido AND Nome_Logradouro = 'NomeEscolhido' AND Sufixo_Direcao = 'IntersecaoPosteriorEscolhida'.
 
À partir daí é só o caso de jogar na função: ST_Line_Interpolate_Point(ResultadoQuery,.99). A função já lhe retorna o ponto desejado.
E como usar esta função para localizar um acidente que ocorreu entre intersecções? Bem, como não temos idéia de onde o acidente ocorreu na via, utilizamos o ponto médio do logradouro para localizá-lo. A função utilizada é a mesma, a única coisa que muda é o SQL, onde especificamos também um Prefixo_Direcao e ao invés de .99 na porcetagem da função, utilizamos .5, nos retornando o exato ponto médio daquele trecho.
Bem isso foi apenas uma pequena introdução à geocodificação e ao pequeno trabalho que realizei durante o ano de 2008. O sistema está pronto e funcioona! Mas ninguém usa :P. Já existe um sistema para cadastro de acidentes em Uberlândia, mas é um sistema proprietário, e quem o desenvolvou não tem nem idéia do que Geoprocessamento significa. Ainda vai levar um bom tempo para os órgãos públicos descobrirem o tanto que as novas tecnologias podem ajudar, a poupar e neste caso, até a salvar vidas.

George

#Python: parte 2 e 1/2

Bom dia pessoal,

No post #python parte 2, eu deixei como exercício para vocês uma pequena tarefa. Como podemos perguntar ao usuário qual é pasta que ele quer projetar e qual a projeção. Existem duas maneiras, uma delas é criar o script para rodar somente dentro do ArcGIS ou utilizar algum tipo de manipulação dentro do script, que roda em DOS/Python.
Vou mostrar a segunda maneira. Os scripts criados pelo ArcGIS são bastante simples, mas é preciso um pouco mais de base teórica, então fica pra depois.

Bem, para isso precisamos utilizar umas manhazinhas. Esta solução vai em duas partes, e poderia ser muito mais complexa, mas o Python vai nos ajudar.

Primeiramente, como perguntar qualquer coisa ao usuário? O Python tem uma pequena função chamada raw_input.

A função raw_input tem como único argumento opcional uma mensagem à ser passada ao seu usuário, mostrando a ele o que ele tem de digitar:

pasta = raw_input("Digite a pasta onde se encontram os shapefiles. ")
#funciona
pasta2 = raw_input()
#funciona também, mas sem mensagem

Bem, antes do script rodar e começar a tentar processar as pastas, é bom testar se o valor que o usuário digitou é válido. Existem diversas validações que poderiam ser executadas, mas vamos usar uma em especial.

#lembre-se que cada - equivale à um espaço e um conjunto com ---- equivale à uma tabulacao.

import sys, os, arcgisscripting

pasta = raw_input("Digite a pasta que deseja projetar. Lembre-se de duplicar as barras \\ senão não acho o caminho. ")

if os.path.isdir(pasta) == true:
----#é uma pasta válida e existe.
----#processe
else:
----#não é uma pasta válida.
----#imprima uma mensagem de erro e saia do script.
----print "Erro, a pasta indicada não é válida."
----exit

Bem pessoal, esse é o começo de tudo. Não repetirei os processamentos com o arcgisscripting, e isso é pra vocês estudarem e tentarem por aí.

Não se esqueçam, o F1 do Python traz a ajuda da versão utilizada, que contém as assinaturas de funções (quais são os paramêtros de entrada e qual é a saída), notas importantes e exemplos de como usá-las.

Agora, como permitir o usuário escolher o datum? Essa é um pouco mais difícil. Teremos de usar um outro tipo de dados do Python, chamado dicionário, que é tão poderoso quanto as listas, mas existem algumas diferenças sobre ele. O dicionário não suporta aquela maravilha de "for x in lista:", pois ele é um...dicionário, e não uma lista.

Os dicionários tem alguns métodos que nos ajudam a achar o que queremos:

a = {'SAD6922SUL':'\\Coordinate Systems\\Projected Coordinate Systems\\Utm\\Other GCS\\South American 1969 UTM Zone 22S.prj'}

Este é um dicionário com somente uma entrada. Enquanto as listas utilizam a posição de cada item (de 0 à n), os dicionários utlizam as keys. No caso supracitado, a key para este endereço monstro é SAD6922SUL, ou seja, ela armazena um valor, que é o diretório para este datum.

Métodos:
len(dicionario) retorna o tamanho do dicionario.
k in dicionario. Verdadeiro se a chave k está presente no dicionario dicionario.
k not in dicionario. Verdadeiro se a chave k NÃO está presente em dicionario.
dicionario.has_key(k). Mesmo que k in dicionario. Versão nova da função.
dicionario.keys(). Copia as keys para uma lista.
dicionario.values(). Copia os valores para uma lista.
dicionario[k] = v. Adiciona uma chave k com o valor de v no dicionario.

Existem outros métodos, mas estes são os mais básicos e necessários para começar.

Imagine o seguinte dicionario:

dicionarioDatums = {'sad69utm22s':\\Coordinate Systems\\Projected Coordinate Systems\\Utm\\Other GCS\\South American 1969 UTM Zone 22S.prj',
'sad69utm23s':'\\Coordinate Systems\\Projected Coordinate Systems\\Utm\\Other GCS\\South American 1969 UTM Zone 23S.prj',
'sad69utm24s':'\\Coordinate Systems\\Projected Coordinate Systems\\Utm\\Other GCS\\South American 1969 UTM Zone 24S.prj'}

e o código a seguir:

chaveDatums = dicionarioDatums.keys()

for chave in chaveDatums:
----print chave

datumEscolhido = raw_input("Digite o datum de sua escolha. ")

if datumDicionario.has_key(datumEscolhido):
----#processo
else:
----print "Você não escolheu um datum válido. Reinicie e tente novamente."

Bem, veja que hoje apenas estruturei a solução. Quero ver algum código aí nos comentários. Alguém tem uma outra idéia para fazer isso?

É bom notar que existem melhoras a serem feitas neste script, por exemplo, forçar o usuário a digitar algo correto no último raw_input, e software não sair daí enquanto ele não o fizer. Existe também a possibilidade de se mapear a pasta Coordinate System automaticamente, sem a necessidade de sempre que precisarmos de um datum novo, que não está na lista, reescrever o código. Deem suas sugestões e espero ter ajudado.

Abraço

George

segunda-feira, 17 de agosto de 2009

#Python, parte 3

Boa noite pessoal,

Como disse anteriormente o Python é a linguagem "preferida" pela ESRI para se construir ferramentas e modelos de geoprocessamento.

Na última vez, usei um exemplo pequeno, mas até útil.

Hoje vou mostrar como se pode fazer múltiplas operações, com uma só Feature Class.

O Python é ótimo em executar tarefas tediosas e repetitivas, daquelas que gastaríamos anos para terminar se feito na mão.

Imagine que você precise converter uma série de SRTMs para pontos, juntá-los, e depois interpolar tudo isso (espero que seu computador seja bom :P)

Vamos ao exemplo:

#lembrem-se que o python faz questão de código identado, se ele não tiver identado, não funcionará.
#aqui ao invez de espaços usaremos -

#vamos importar nossas coisas
import sys, os, arcgisscripting, string

gp = arcgisscripting.create(9.3)

diretorio = raw_input("Digite o diretório que contém os rasters. Não se esqueça de duplicar as barras \\, senão não conseguirei achar o caminho...")
#aqui você pode indicar qualquer coisa, um diretório, um banco de dados SDE e até um banco de dados file/personal.
if diretorio:
----pass
else:
----print "Workspace/Espaço de trabalho inválido. Reinicie o script e tente denovo."
----exit()

gp.Workspace = diretorio

ListaRasters = gp.ListDatasets("*","Raster")
#se tivessemos escolhido "SRTM*" ao invés de "*" o Python listaria todos os rasters que começam com as letras SRTM

gp.AddToolbox("conversion") #adiciona a toolbox conversion
gp.AddToolbox("sa") #adiciona a toolbox do spatial analyst - só funciona para quem TEM a licença eim pessoal.
gp.AddToolbox("management") #adiciona management

#já temos nossas ferramentas, mãos à obra.

for Raster in ListaRasters:
----OutFeature = diretorio + "\\" + "ponto_convertido_" + Raster
----gp.RasterToPoint_conversion(Raster,OutFeature)
----print Raster + " convertido para ponto."

#vamos listar nossas feature classes de ponto. lembra do "*"?
ListaFeatureClasses = gp.ListFeatureClasses("ponto_convertido_*","POINT")
#assim garantimos que só vamos mosaicar as featureclasses convertidas e não outras perdidas no mesmo workspace, e claro só do tipo ponto :D

#a ferramenta merge exige uma lista das featureclasses a serem juntadas separadas por ponto e virgula, entao vamos dar a ela o que ela quer!

OutputMerge = gp.Workspace+"\\" + "MergeFinal"

gp.Merge_management(string.join(ListaFeatureClasses,";"),OutputMerge)
#a funcao join concatena uma lista utilizando um caractere (ou caracteres separadores) ou seja: ['shape1','shape2',shape3'] viram shape1;shape2;shape3

print "Já está tudo junto, vamos interpolar?"
print "interpolando. deve demorar, vá tomar um café..."

gp.Idw_sa(OutputMerge,"VALUE",gp.workspace+"\\"+"raster_interpolado",30,1)

print "interpolado. parabens!"

Bem pessoal, não testei o código, mas deve funcionar de acordo. O python é muito poderoso e suas funções internas se encaixam muito bem com a API da ESRI. Vejam na parte do Merge. Em outras linguagens teriamos de escrever loops e loops para concatenar as palavras daquela maneira, mas em python tudo é feito rapidinho!

Espero que gostem e que comecem a utilizar Python no seu dia a dia. É de fato um negócio muito útil. 1 hora para escrever um script, mas te salva 80mil horas de trabalho chato, cansativo e passível de erro humano.

Ahn! Não se esqueçam de consultar a documentação das funções utilizadas (join do Python, Merge_management, RasterToPoint_conversion e IDW_sa). Todas as versões do Help do ArcGIS possuem essa documentação. É só apertar F1 e correr para o help de cada ferramentinha, na última parte delas, Scripting. Se acharem algum erro no código postados, favor dar um toque que eu corrijo.

Python Rocks!

quarta-feira, 12 de agosto de 2009

#Python parte 2 - utilizando objetos de geoprocessamento

Buenas pessoal,

Bem, Python, como expliquei no post passado é uma linguagem de alto nível e muito recomendada pela ESRI para se construir modelos de geoprocessamento (semelhantes aos construídos no ModelBuilder) e scripts para processamento de informações.

Primeiro, para quem deseja aprender Python:

python.org e python.org.br

Aqui encontramos algumas coisas fundamentais, como documentação e links para trocar idéias com a comunidade que desenvolve em Python.

Segundo: irc.freenode.net nos canais #python e #python-br . Nestes locais é possível conversar com desenvolvedores brasileiros e internacionais e resolver possíveis dúvidas.

Terceiro: este aplicado a quem quer integrar esta poderosa ferramenta aos seus métodos de trabalho - WebHelp ESRI e WebHelp Esri 2.

Não tenha medo de errar e não tenha medo de buscar ajuda.

Então vamos lá.

Para qualquer coisa que você queira fazer em python, é necessário importar o módulo que cuida dos objetos geográficos da ESRI, o arcgisscripting.

#importando o módulo arcgisscripting
import arcgisscripting

Quando um script .py contém esta instrução (ela deve ser feita no ínicio do seu script), o python deixa a sua disposição os objetos da ESRI. No webhelp você pode encontrar uma descrição ampla de TODOS eles, bem como exemplos.

Certo, mas depois que este módulo está disponível, o que fazemos para começar a mexer com o ArcGIS?

É necessário criar um objeto do tipo geoprocessing, fornecido pelo módulo arcgisscripting. Olhem o exemplo

gp = arcgisscripting.create()
#este método, sem parâmetros serve tanto para 9.2 quanto para 9.3
#para 9.3, é melhor especificar o paramêtro 9.3 e.g.: gp = arcgisscripting.create(9.3)

O que podemos fazer com o Python? Bem, em teoria, tudo o que estiver disponível dentro de uma Toolbox do ArctoolBox também está disponível para o Python. O Python na verdade, consegue fazer qualquer coisa que esteja em uma toolbox e ainda mais!

Lembre-se que o Python contém muitos módulos, como acesso à sítios remotos, módulos para gerenciar emails, matemáticos, FTP, conversores para PDF, acesso à banco de dados (MSSQL, PostgreSQL, SQLite, entre outros), manipuladores XML...ou seja, é muito flexível e extensível.

Para realizarmos operações dentro do ArcGIS, temos que utilizar o módulo sys e geralmente o módulo os, ambos nativos no Python. Hoje não falarei de como construir scripts que rodam à partir do ArcGIS, apenas scripts que rodam diretamente da shell python (alguém conhece alguma utilidade pra isso? acessar objetos geográficos sem abrir o ArcGIS?)

Bem, imaginemos que queremos Projetar ou Definir projeção de diversos itens, ao mesmo tempo:

#projetar diversos shapefiles para uma única projeção
import sys, os, arcgisscripting

Caminho = "C:\\Pasta\\Shapefiles\\"
#no Python, temos de duplicar as barras, pois uma barra \ é interpretada como alguns
#caracteres especiais, como nova linha, tabulação, entre outros

gp = arcgisscripting.create(9.3)

#vamos projetar shapefiles/featureclasses. temos de utilizar
#as toolboxes correspondentes, entao vamos adicioná-las

gp.AddToolbox("C:\\CaminhoAtéAToolbox\\Nome da Toolbox.tbx")

#podemos adicionar tantas toolboxes quanto precisarmos
#só devemos lembrar que mesmo através de scripts Python, as licensas de extensões
#SÃO CHECADAS. se você não tiver o spatial analyst, o Python irá dar um chilique

gp.Workspace = Caminho
#vamos definir o workspace. No ArcGIS + Python, muita coisa depende do workspace.
#ele é basicamente um apontamento para o diretório correto.

ListaShapefiles = gp.ListDatasets("*","ALL")
#vamos criar um objeto do tipo lista, e inserir nela todos os objetos do arcgis.
#os parametros "*" e "ALL" está dizendo para o ArcGIS inserir na lista todos
#os objetos, sem especificação de tipo (que podem ser Feature, TIN, Raster e CAD)

print ListaShapefiles
#me mostre o que será projetado

projecao = "C:\\CaminhoAteAProjecao\\Projecao.prj"

for Shapefile in ListaShapefiles:
----gp.Project_management(Shapefile,Shapefile+"_projetado",projecao)
----print "Projetei " + Shapefile

#considere os quatro - como uma identação, um toque na tecla tab (o blog está engolindo os espaços)
#fim

Este script é bem simples, deve ser salvo como .py, e ser rodado através da shell python ou com dois cliques no mesmo (desde que o Python esteja devidamente e corretamente instalado).

Algumas notas no entanto: note a forma que duplicamos as barras. Se você utilizar somente uma barra o Python não conseguirá achar os caminhos corretamente.

Os nomes de saída da função project são Shapefile_projetado.

A identação no trecho "gp.Project_management..." é PROPOSITAL E DEVE ESTAR LÁ. SEM ELA O SCRIPT NÃO FUNCIONARÁ. O Python usa a identação obrigatória como forma de controle de código.

Não se esqueça, tudo o que vem depois de # é considerado pelo Python um comentário, não influenciando no programa.

Note que para projetar uma pasta/set de shapefiles diferentes, temos que alterar a variável Caminho, lá em cima. O que podemos fazer para dinamizar esta questão? Como podemos perguntar ao usuário qual é a pasta que precisamos projetar? E a projeção, como definir ela dinamicamente? Essa fica de dever de casa :P

E agora senhores, o que podemos fazer com Python? Podemos utilizar todas as Toolboxes do ArcGIS para realizar milhares de operações com diversos shapes/featureclasses, e até mesmo cosntruir modelos complexos, do tipo que vemos em ModelBuilder, simulando eventos.

Vamos aprender Python?

segunda-feira, 10 de agosto de 2009

#Python: o que é?

Olá pessoal,

Primeiro gostaria de falar sobre a Geo Rede. GeoRede é um projeto do site Geoprocessamento.net, um agregador de feeds de diversos blogs sobre Geoprocessamento/Sensoriamento Remoto em geral.

Estou lá também! Depois passem lá e deêm uma olhada. Se você não conhece o Geo.NET é mais uma coisa bacana. Fórum, seção de downloads e muito mais.

Bem, hoje falo para os usuários de ArcGIS. Vocês já devem ter notado, que durante a instalação do ArcGIS ele te pergunta/pede para instalar uma linguagem de programação chamada Python. Bem, se algum de vocês já abriram a aba dela no menu Iniciar do windows devem ter se deparado com algo bastante simples, a IDLE. A IDLE é um shell que intepreta comandos em Python. Algo como o DOS.

Bem, não dá pra fazer muita coisa pela IDLE, já que estamos limitados a "programas" de uma linha. Mas o porque estou falando de Python? O ArcGIS, felizmente, tem um ambiente customizável e permite ao usuário mais avançado a desenvolver suas próprias ferramentinhas, para executar uma infinidade de ações. Desde tarefas de "geoprocessamento" (geoprocessing) à customizar suas layers e representações.

Existem duas maneiras de se estender o ArcGIS. Uma delas é através de uma linguagem compilada .NET (ou VBA) e através do Python.

De acordo com a própria ESRI, a linguagem recomendada para se trabalhar com geoprocessing é Python. O .NET é mais fácil caso você precise estender a GUI (Graphic User Interface) e montar coisinhas mais complexas. O Python é recomendado para geoprocessing e a criação de ferramentas mais simples justamente pelo Python ser simples. A linguagem Python, criada em 1991, é considerada VHLL (Very High Level Language) e traz enormes facilidades para desenvolvedores iniciantes.

Primeiramente, ela é dinamicamente tipada (what the hell?!). Bem isso significa que não é necessário declarar com QUE tipos de variáveis queremos trabalhar antecipadamente. Você coloca alguma coisa na variável e o Python identifica o que ela é pra você. Ou seja, não existe muita confusão no momento de falar que A = 3 ou A = "spam". (neste exemplo, A = 3, representa automaticamente que A é uma variável do tipo inteiro e que em segundo momento, ela é do tipo string)

Segundamente, Python não é uma linguagem compilada. Não é necessário desenvolver um arquivo executável ou algo do tipo, ela é interpretada conforme o Python lê o código digitado. Isso traz uma perda de velocidade (nada muito importante) mas o ganho em agilidade no desenvolvimento compensa.

Python tem tipos muito flexíveis, como listas e dicionários. Não entrarei em detalhes, mas você pode com uma única expressão, criar uma lista com valores que você esteja trabalhando e de forma mais fácil ainda, ler os mesmos.

Um exemplo de lista, que sempre estão entre COLCHETES:

ListaTeste = ['spam','geoprocessamento',2,['arcigs','envi','arcsde']]

Temos na lista acima os items 'spam', 'geoprocessamento',2 e uma outra lista! aninhada, com os objetos 'arcgis','envi' e 'arcsde'. Listas são o carro-chefe do Python e devem ser aprendidas o quanto antes. Através delas podemos fazer o capeta com muito pouco código.

Agora, como integro ArcGIS com Python? A ESRI disponibiliza uma "API" própria para manipulação (de seus) objetos geográficos, como Feature Classes, Layers, Datasets, Geodatabases entres outros.

tudo o que tiver após # são considerados comentários, e não são interpretados pelo Python.

Um exemplinho simples de um script em Python (embora não muito útil) manipulando dados do ArcGIS seria:

import os, sys, arcgisscripting

gp = arcgisscripting.create(9.3) #vamos criar um objeto do tipo geoprocessing, versão 9.3

ListaCamadas = gp.ListDatasets(r"C:\teste python\geodatabase.gdb")
#pega uma lista dos featuredatasets no geodatabase.gdb

for camada in ListaCamadas: #para cada camada em ListaCamadas
print camada #imprima camada

#a resposta do Python vai ser algo assim.

>>>>['Hidrografia','Altimetria','Geologia','Transportes']

Bem, vejam como é fácil andar pelas listas e realizar operações em série. Os loops em Python são descomplicadíssimos.

Sugiro que agora, corram ao site da ESRI, e deem uma olhada nos geoprocessing Objetcs. Eles podem facilitar muito a vida de um utilizador de SIG.

Um abraço

sábado, 1 de agosto de 2009

Leituras...

Sábado feio de São Paulo...

Ontem comprei dois livros: Aprendendo Python da O'Reilly e o SQL Server 2008 para Desenvolvedores.

Tem que estudaaaar. Mas agora é hora de Antartica.

Abraços