Carta de Ocupação e Uso do Solo (2018) - v3

Nas minhas mais recentes experiências e ensinamentos de vida no que toca a trabalhar para outrem (salvo seja), todas as críticas e erros a mim apontados são encarados não como uma forma de rebaixamento mas sim como uma forma de auto enriquecimento e desafio a fazer melhor.

Nesse sentido, com vista a otimizar a vida e testar formas de assegurar um trabalho o mais ágil e corretamente possível, no âmbito da atualização da COS 2018 para a versão 3, desenvolvi uma série de ferramentas, scripts e, já em fase de testes, tanto um plug-in do QGIS como um programa compilado; todos com o fim de acelerar trabalho, monitorização de modo automático e resultados de qualidade, conforme esperado pela Direção-Geral do Território (DGT).

A concurso foram 3 empresas, sendo por elas distribuídos 3 lotes de Portugal Continental:

Fonte: DGT, Caderno de Encargos de Concurso Público CP/270/2023, Anexo A.

Todas têm o mesmo objetivo: assegurar uma exatidão temática superior a 85% com vista a assegurar um pagamento de 100% das tranches pela DGT.

Peripécias à parte durante o período de formação - sem entrar muito pela parte em que o próprio manual de fotointerpretação diz uma coisa e na realidade se aplica outra, e vice-versa; ou à divergência de opiniões de como segmentar polígonos, baseando-se em meras opiniões ou em dados meramente declarativos - existe um caderno de encargos que nuns pontos é lato, noutros é bastante sucinto. É neste ponto que a atualização da COS 2018, "através da correção de eventuais erros temáticos e geométricos" como ponto de partida para futuras versões da COS se justifica, surgindo novas classes e regras a aplicar à geometria e à temática, nomeadamente:

  • Unidade Mínima Cartográfica (UMC) passar de 1 hectare para meio hectare, o que levará à segmentação de muitos mais polígonos até então não possíveis de individualizar;
  • Distância mínima entre linhas (DMEL) e, por conseguinte, a largura dos polígonos, passar a respeitar um mínimo de 20 metros entre nós e edges, aplicando-se igualmente aos polígonos a terminar em bico, sendo restrito essa representação:
  • Correção da temática dos polígonos face ao novo número e tipologia de classes;
  • entre outros.
Assim, face ao exigido e de modo a que os trabalhos corram da melhor maneira, como já anunciado no topo desta página, desenvolvi uma série de estratégias para facilitar a vida aos foto-intérpretes. Como todo o trabalho será feito em QGIS, todas as implementações foram feitas tendo em conta esse conjunto de softwares.

1. Script para alívio dos servidores e rapidez de trabalho

Com vista a acelerar o trabalho dos foto-intérpretes, o primeiro script foi estudado de modo a que os foto-intérpretes tivessem que estar o mínimo possível "pendurados" ao servidor onde se encontram as ortofotos, pesadas e em elevado número para cada série temporal*

*friso série, pois cada ano é composto por varrimentos feitos em diferentes momentos e em diferentes coberturas, sem qualquer correção à refletância, logo diferentes tons de cores/falsas cores em zonas de transição, o que por sua vez e a título de exemplo aquilo que é vermelho numa cobertura, pode estar castanho noutra, veja-se a imagem seguinte:

Aqui consegue-se constatar de facto a diferença de refletâncias e a sua não correção para diferentes momentos de levantamento da mesma zona. Onde teoricamente em IRG seria normal ver-se sobreiros mais cor-de-rosas do que pinheiros mansos, em 2015 mostra, por exemplo, os sobreiros estarem mais em castanho e os pinheiros mansos com cores mais vivas; noutro momento a situação inverte-se.

Este pode ser um problema mediante treino de algoritmo de classificação supervisionada, uma vez que à falta de correção da refletância aplicada, dever-se-á recorrer a dados com essa mesma correção já aplicada e, só depois, treinar para as áreas cobertas pelas ortofotos com diferentes séries temporais.

Assim, tendo disponibilidade de espaço em disco nos computadores, o script, antes de ser corrido, requer correr um modelo do QGIS que navega pela CROA dos diferentes anos e dimensões de enquadramentos, interseta as áreas de trabalho com esses enquadramentos, e devolve um ficheiro temporário com o enquadramentos face às ortofotos de 2018, onde em dois dos seus campos contém os nomes dos ficheiros das ortofotos tanto de 2018 como <= 2015, sendo estes últimos inalteráveis independentemente do ano.

Com isto, o script lê esses campos, procura no servidor pelas ortofotos, se encontrar copia para o disco do foto-intérprete, e, mediante processo feito pelo modelo, caso haja clusters de áreas de trabalho intersetadas pelas células das grelhas de enquadramento, mas entre essas células (com outras áreas) não haja continuidade, isola em grupos, de forma a tornar o processo de rendering mais rápido e focar a atenção do foto-intérprete só naquelas imagens, isto é, daquele cluster de polígonos. Com isto, estrutura-as no Projeto de modo a estarem organizadas por grupos, e dentro de cada grupo, sub-grupos com o tipo de cores (RGB/IRG), e dentro desses sub-grupos, outros subgrupos por ano (2007, 2010, 2012, 2015 e 2018).

Exemplo de grupos seccionados de ortofotos, depois de corrido o modelo, e onde o script Python se irá centrar para fazer a cópia dos ortofotos do servidor, seguido de importação para o projeto por grupos. A quadrimetria usada corresponde às ortos 2018, que por sua vez cada enquadramento corresponde à metade vertical de cada enquadramento das ortos <= 2015.

Com este procedimento, apenas se perde um longo tempo uma vez, pois o script deve ser corrido na consola do Python do QGIS - usando módulos do QGIS - e, enquanto o mesmo não terminar o processo, o QGIS não está responsivo. Por outro lado, liberta o foto-intérprete de trabalho manual de inserção das ortofotos com diferentes dimensões, de diferentes anos e, assegura o acesso imediato a estas a partir do seu disco, libertando o servidor de sobrecargas e atrasos na transmissão.

2. Expressões QGIS para evitar erros e acelerar o procedimento de classificação de polígonos

Com vista a providenciar uma forma ágil e rápida de classificar os polígonos sem cair em erro de classificação, conforme secção II do Anexo D do Caderno de Encargos do já referido Concurso Público, desenvolveu-se uma série de expressões e filtros que evitam erro humano e, por conseguinte, perda de exatidão temática na atualização da COS.

Assim, tendo em conta a linguagem das expressões do QGIS e sua especificidade, o seguinte pseudocódigo demonstra o processo de classificação automática da espécie principal, que quando conjugado com outro tipo de expressões para a espécie secundária, eis que o foto-intérprete reduziu drasticamente o número de cliques por polígono, bem como não tem hipótese de se enganar na introdução da espécie florestal principal.

Pseudocódigo das expressões para espécies florestais principais, baseando-se em strings e em campos específicos para auto-atribuição de valores àquele campo.

Estas expressões não evitam, no entanto, que o foto-intérprete cometa erros de segmentação, se esqueça de segmentar alguma porção ou deixe algum dos campos por preencher, campos completamente alheios à classificação, mas respeitantes à gestão interna e totalmente alheios àquilo que a DGT irá por parte das empresas receber.

3. Teste de classificação automática de polígonos

Para tal, como ajuda ao processo de decisão sobre a correta ou incorreta segmentação por parte dos foto-intérpretes, encontra-se em desenvolvimento um algoritmo de classificação supervisionada tendo como base a rasterização de uma nuvem de pontos de dados LiDAR do projeto áGIL do ICNF com um tamanho de pixel de 50 cm, a fim de corresponder com a realidade das ortofotos disponibilizadas pela DGT para os anos de 2007 a 2015. Sabendo que as ortofotos de 2018 têm uma resolução de 25 cm, bastou redimensionar os pixels do raster produzido através dos dados LiDAR para a dimensão dos ortos de 2018, e o mesmo aplicar às ortos de 2007, 2010, 2012 e 2015. A área-amostra dos dados LiDAR diz respeito à interseção de uma área de 4 x 4 quadrículas da CROA 2018, correspondendo aos restantes anos a uma área de 2 x 4 na respetiva CROA; usou-se a dimensão dessa área amostra para fazer o ponto de partida dos pixels, com snap à mesma resolução e dimensão.

Exemplo da nuvem de pontos a ser utilizada, sendo rasterizada para uma resolução de pixel de 50 cm.
 

Usou-se uma área da Zona Centro do país para testes, extraindo os valores de refletância daquele levantamento LiDAR - que se assume corretamente feito (até prova em contrário) -, recorrendo à combinação IRG para aplicar às ortofotos de 2007 a 2018. Este processo não é perfeito, mas permitiu gerar alguma homogeneidade. Sabendo-se que o levantamento dos dados LiDAR foi feito em 2021, usar-se-á aqueles valores para aplicar a correção radiométrica às imagens dos diferentes anos, mas usando como ponto de partida as ortofotos de 2021 da DGT, existentes apenas para a metade Norte do país.

Na fase seguinte, e depois de ter as ortos "corrigidas", treinou-se o modelo de modo a poder ser feita uma amostragem com uma margem de erro de 5%, requerendo para tal o devido número de pontos amostra e a sua classificação manual. O total da exatidão do produtor vs. utilizador foi de 86%. Assim, recorrendo aos dados auxiliares do ICNF e do IFAP, aplicou-se um algoritmo que recorre a bibliotecas Python (maioritariamente Geopandas e suas dependências) que cruza a versão atual e em voga (os seus polígonos) da COS com esses dados auxiliares e o resultado em raster da classificação supervisionada. Assim, contando pixel a pixel, cruzando os dados auxiliares com aqueles valores "já sabidos" de correspondência às espécies amostradas, poligoniza-se um resultado raster de prováveis classificações, pureza ou mistura de espécies e, ainda sem prévia correção, dos limites, desrespeitando a DMEL e UMC.

Apresenta-se como uma primeira aproximação, e como tal, encontra-se em desenvolvimento em paralelo um programa compilado para acautelar essas falhas.

4. Programa de correção geométrica e topológica

Tendo em conta os vários momentos de atualização (neste caso, reclassificação) dos polígonos da COS, chega um momento em que a qualidade do produzido deve ser assegurado. Muitos métodos e procedimentos podem ser seguidos sem precisar sair do QGIS, o que não deixa de ser um método user-friendly e eficaz para resolver quaisquer incongruências. Por outro lado, eficácia não é sinónimo de eficiência, e é nesse sentido que tratei de fugir à lentidão do QGIS e compilar um programa que tem por base a linguagem Python, que em si embute as bibliotecas e ambiente virtual Python necessários e, que embora ainda sem uma interface gráfica (não totalmente necessário), faz o processamento da informação vetorial em shapefile ou geopackage mediante o utilizador carregar para as pastas a informação necessária; a única vez que o programa interage com algo externo a si é no momento de criação de um backup do ficheiro geral da COS, isto é, caso tenha sido intervido "ontem" para controlo de qualidade, adiciona um sufixo "_DDMMAAAA" ao nome do ficheiro em base de dados, importada por meio de chave de acesso em PostgreSQL. Assim, todo este procedimento é feito no disco do utilizador, assegurando celeridade e absente de pôr em risco quaisquer dados no servidor.

Exemplo de um dos passos iniciais do programa antes de entrar em processamento dos dados.

Este programa está estruturado para varrer os dados, passo a passo, exportando para pasta específica os dados corrigidos e criando variadas geopackages com diferentes geometrias assinalando onde foram encontrados esses mesmos problemas geométricos e/ou topológicos, caraterizando o tipo de problema por geopackage. Esses ficheiros servem para mostrar ao utilizador o tipo de problemas e em que polígonos, para despiste de eventuais correções automáticas feitas erradamente.

Algoritmo implementado pelo programa, em todos os seus passos, para utilização independente pelo plug-in e importação de resultados finais no QGIS.

Este programa produz resultados finais independentemente de ser corridos pelo plug-in ainda em testes.

5. Plug-in de processamento de dados da COS

Este plug-in não faz verdadeiramente correções aos polígonos (como anunciado no final do fluxograma anterior), mas exportando uma geopackage de linhas assinalando os locais onde os seus princípios são violados. Assim, este plug-in apenas analisa a componente geométrica, pelo que as correções mencionadas no ponto anterior devem estar previamente garantidas, assinalando:

  • locais onde a DMEL é violada (estreitamentos de polígonos);
  • "cantos" nos polígonos que embiquem e, por sua vez, desrespeitem a DMEL;
  • pares de polígonos vizinhos com a mesma classificação, espécies principal e secundária.

Assim, percorrendo iterativamente cada polígono de um universo geral ou isolado de polígonos (por exemplo, aplicando um filtro à base de dados), o resultado final apresenta-se por uma série de linhas mais ou menos perpendiculares tanto a edges como a nós onde a distância entre estes elementos seja inferior a 20 metros, unindo ao lado oposto onde se assinala a violação dessa regra.

Por outro lado, nos cantos dos polígonos, aplica regras trigonométricas tendo como base as coordenadas dos nós, determinando todos os ângulos inferiores a 90 graus e, pela Lei dos Senos, determina a distância dos lados opostos (união dos vértices do triângulo fictício) ao do ângulo. Se confirmadas menores que a raiz quadrada de 200 (14.1421... metros), então aplicando o Teorema de Pitágoras teremos o threshold mínimo para aplicar a regra.

lei dos senos
Fórmula da Lei dos Senos

O esboço seguinte demonstra o explanado anteriormente:

esboço ângulos
Por fim, faz-se a verificação aos polígonos com igual classificação, assinalando no mapa os polígonos em questão, unindo-os por setas a partir dos seus centróides (usando o algoritmo de Moment Centroid (Deakin et al.)), sempre assinalados por pares.

Par único de cardinalidades - dados alterados propositadamente para demonstração

cardinalidades3
Correspondência de cardinalidades de 1 para n - dados alterados propositadamente para demonstração.

Este procedimento assume, assim, a missão de encontrar as correspondências por leitura dos campos em geodataframe, encontrando primeiramente as interseções em 1 ou n nós de igual coordenada. Traduzindo sucintamente, todos os nós com igual coordenada são lidos e postos numa array e, caso a si estejam atribuídos dois polígonos com iguais campos de classe, espécie florestal principal e secundária, calcula o centróide e une-os com uma linha, à qual aplica um estilo, e exporta os polígonos visados para uma geopackage, onde se incluem os polígonos repetidos, os centróides e as linhas. 

 

[ post abandonado por desligamento a este projeto ]