Obtenção de Coordenadas Por Meio da Monorrestituição
Descrição do Produto
0
UNIVERSIDADE FEDERAL DO PARANÁ ENGENHARIA CARTOGRÁFICA E DE AGRIMENSURA
KAUÊ DE MORAES VESTENA
DETERMINAÇÃO DE COORDENADAS TRIDIMENSIONAIS DE PONTOS EM UMA IMAGEM EMPREGANDO A TÉCNICA DA MONORRESTITUIÇÃO.
CURITIBA 2015
1
KAUÊ DE MORAES VESTENA
DETERMINAÇÃO DE COORDENADAS TRIDIMENSIONAIS DE PONTOS EM UMA IMAGEM EMPREGANDO A TÉCNICA DA MONORRESTITUIÇÃO.
Trabalho realizado para a disciplina de Fotogrametria III, do curso de Engenharia Cartográfica e de Agrimensura, sob a orientação do Profº Dr. Edson Aparecido Mitshita,
como
requisito
avaliação semestral.
CURITIBA 2015
parcial para
2
LISTA DE FIGURAS
FIGURA 01: IMAGEM AÉREA UTILIZADA NO TRABALHO ............................................................... 11 FIGURA 02: MDT DA REGIÃO ............................................................................................................. 13 FIGURA 03: EMPREGO DO QGIS NA COLETA DAS COORDENADAS NO REFERENCIAL DIGITAL ............................................................................................................................................................... 14 FIGURA 04: ILUSTRAÇÃO DO PROCESSO DA MONORRESTITUIÇÃO.......................................... 24
3
LISTA DE TABELAS
TABELA 01: PARÂMETROS DE ORIENTAÇÃO INTERIOR E DIMENSÕES FORNECIDAS PARA O SENSOR ............................................................................................................................................... 10 TABELA 02: COORDENADAS DOS PONTOS DE APOIO ................................................................ 12 TABELA 03: COORDENADAS DOS PONTOS DE APOIO NA IMAGEM .......................................... 15 TABELA 04: FOTOCOORDENADAS EM REFERENCIAIS MÉTRICOS ........................................... 27 TABELA 05: FOTOCOORDENADAS EM REFERENCIAIS MÉTRICOS ........................................... 27 TABELA 06: MVC DOS POE AJUSTADOS ........................................................................................ 28 TABELA 07: FOTOCOORDENADAS AJUSTADAS E OS RESPECTIVOS RESÍDUOS .................. 28 TABELA 08: TESTES ESTATÍSTICOS PARA O AJUSTAMENTO DOS POE .................................. 29 TABELA 09: VALORES DAS COORDENADAS TRIDIMENSIONAIS DOS PONTOS DE APOIO DETERMINADAS PELA MONORRESTITUIÇÃO ............................................................................... 30 TABELA 10: DISCREPÂNCIAS ENTRE AS COORDENADAS DETERMINADAS E AS FORNECIDAS ....................................................................................................................................... 30 TABELA 11: ANÁLISE DAS DISCREPÂNCIAS ................................................................................. 31
4
LISTA DE ACRÔNIMOS
CP
Centro Perspectivo da câmera: abstração pontual do sistema de lentes
f
Distância focal calibrada da câmera;
ECI
Equações de Colinearidade Inversa
MVC
Matriz de Variância-Covariância;
MDT
Modelo Digital do Terreno
POI
Parâmetos de Orientação Interior;
POE
Parâmetros de Orientação Exterior;
pp
ponto principal
x
Coordenada do eixo das abcissas no referencial fotogramétrico
x0
Coordenada do eixo das abcissas do ponto principal, no referencial com origem no centro da fotografia.
X
Coordenada do eixo das abcissas no referencial geodésico local;
X0
Coordenada do eixo das abcissas do CP no referencial geodésico local;
y
Coordenada do eixo das ordenadas no referencial fotogramétrico;
y0
Coordenada do eixo das ordenadas do ponto principal, no referencial com origem no centro da fotografia
Y
Coordenada do eixo das ordenada no referencial geodésico local;
Y0
Coordenada do eixo das ordenadas do CP no referencial geodésico local;
Z
Coordenada do eixo da Normal, no referencial geodésico local;
Z0
Coordenada do eixo da normal do CP no referencial geodésico local;
ω
“Ômega”, ângulo de rotação ao longo do eixo das abcissas;
φ
“Phi”, ângulo de rotação ao longo do eixo das ordenadas;
κ
“Kappa”, ângulo de rotação ao longo do eixo da Normal;
X²
Qui-quadrado.
5
SUMÁRIO
1 INTRODUÇÃO ........................................................................................................ 4 2 OBJETIVOS ........................................................................................................... 5 2.1 OBJETIVOS GERAIS ........................................................................................... 5 2.2 OBJETIVOS ESPECÍFICOS ................................................................................. 5 3 REVISÃO DE LITERATURA .................................................................................. 6 3.1 SISTEMA REFERENCIAL DIGITAL ..................................................................... 6 3.2 SISTEMA REFERENCIAL FOTOGRAMÉTRICO ................................................ 6 3.5 ORIENTAÇÂO EXTERIOR PELA RESSEÇÃO ESPACIAL .................................. 8 3.6 MODELO DIGITAL DE TERRENO ...................................................................... 9 4 MATERIAIS E MÉTODOS .................................................................................... 10 4.1 DADOS INICIAIS ................................................................................................ 10 4.1 COLETA DOS PONTOS NO REFERENCIAL DIGITAL .................................... 13 4.2 CONVERSÃO DO REFERENCIAL DIGITAL PARA O FOTOGRAMÉTRICO.... 15 4.3 A RESSEÇÃO ESPACIAL.................................................................................. 17 4.4 A MONORRESTITUIÇÃO ................................................................................... 22 4.5 A IMPLEMENTAÇÃO ......................................................................................... 24 5 RESULTADOS E ANÁLISE................................................................................... 27 7 CONCLUSÃO........................................................................................................ 32
4
1 INTRODUÇÃO
A capacidade de se localizar, identificando tudo que está no seu entorno, é uma capacidade única do ser humano. Com o passar do tempo, várias técnicas foram sendo desenvolvidas para realizar o posicionamento relativo a um sistema de coordenadas. Com o passar do tempo, diversas técnicas de mensuração foram surgindo, inicialmente técnicas de posicionamento astronômicas e geodésicas. Com o desenvolvimento. Com o avanço da civilização e advento dos sistemas de imageamento, o ser humano tornou-se capaz de registrar abstrações de “pedaços” do mundo em que vive, armazenadas inicialmente em papel, e muitos anos depois armazenadas em meio digital. Tais abstrações são denominadas fotografias, onde encontra-se o registro da energia eletromagnética (de uma determinada faixa do espectro eletromagnético), que é conduzida por meio do sistema de lentes até o sensor. A percepção, da aplicabilidade das fotografias para uso profissional para estudos em diversas áreas do conhecimento, veio desde o registro das primeiras cenas que se tem relato. De fato percebeu-se, que poderiam ser utilizadas para realizar medidas relativas a diferentes entidades ou elementos contidos na superfície física da terra. Assim surgiu a fotogrametria, que basicamente, estuda a obtenção de informações métricas à partir de fotografias. Grande parte das aplicações da fotogrametria baseia-se no emprego de duas ou mais imagens, assim, para obter medidas de um objeto, geralmente o mesmo deve ser fotografado sob dois pontos de vistas diferentes. Entretanto, há técnicas que possibilitam a obtenção de informação métrica a partir de apenas uma imagem, a grande maioria de caráter mais aproximado, geralmente voltadas à estudos expeditos. No presente trabalho será abordado um exemplo que visa a determinação de coordenadas tridimensionais, de pontos de interesse em uma fotografia, por meio da técnica da monorrestituição, que se destaca entre as técnicas “monoscópicas” como a mais rigorosa de todas, cujo rigor em um caso específico será alvo de estudo no presente trabalho.
5
2 OBJETIVOS
2.1 OBJETIVOS GERAIS
-Determinar coordenadas tridimensionais de pontos de interesse em uma fotografia (orientada, tomada com uma câmera cujos parâmetros de orientação interior são conhecidos, com um modelo digital de terreno disponível), por meio da técnica da monorrestituição;
2.2 OBJETIVOS ESPECÍFICOS
Coletar as coordenadas dos pontos de apoio no referencial digital, por meio de croqui previamente disponibilizado;
Conversão destas para o referencial fotogramétrico, aplicando
devidas
correções;
Selecionar entre os pontos de apoio, cinco para determinar os POE (Parâmetros de Orientação Exterior) por meio da Resseção Espacial. Encontrar uma única estimativa, bem como a estimativa da qualidade da solução, utilizando ajustamento pelo método dos mínimos quadrados;
Realizar o procedimento da monorrestituição para todos os pontos de apoio;
Comparar as coordenadas obtidas pela técnica explicitada no item anterior, com as fornecidas e a partir disso, proceder com uma análise estatística dos resultados.
6
3 REVISÃO DE LITERATURA
3.1 SISTEMA REFERENCIAL DIGITAL
O sistema de coordenadas digital consiste em uma matriz onde cada linha e coluna é atribuído de um valor de nível de cinza para determinada banda, portanto cada pixel possui um par de coordenadas linha e coluna com um valor de cinza. A origem deste sistema se dá na extremidade superior esquerda da imagem, eixo das colunas com inicio na origem em direção a extremidade superior direita, e eixo das linhas com inicio na origem com sentido positivo na direção inferior esquerda. De acordo com Scuri (2002), essa localização preferencial para o canto superior esquerdo, com eixo y invertido em relação ao padrão da matemática, é herdada dos próprios dispositivos de visualização, onde a primeira linha a ser redesenhada na tela é a do topo da tela, que segue o padrão de escrita ocidental, da esquerda para a direita e de cima para baixo.
3.2 SISTEMA REFERENCIAL FOTOGRAMÉTRICO
De acordo com SANTOS (2008) o referencial fotogramétrico, é um sistema cartesiano tridimensional, cuja orientação dos eixos CPx e CPy é paralela, respectivamente, aos eixos do referencial com origem no centro da imagem, O eixo CPz é coincidente com o eixo óptico da câmera, cujo sentido é tal que torna o sistema dextrógiro. A origem do sistema fotogramétrico é o centro de projeção do sistema de lentes, ou seja, o ponto nodal anterior, caso se esteja trabalhando com o diapositivo ou o ponto nodal posterior se for utilizado o negativo. Como a coordenada z é constante e igual à distância focal da câmera, utiliza-se um referencial plano cuja origem é o ponto principal (pp).
7
3.3 PARÂMETROS DE ORIENTAÇÂO INTERIOR
Segundo SANTOS (2008) “uma câmera possui um conjunto de parâmetros de orientação interior, sendo eles, a distância focal da câmera (f), as coordenadas do ponto principal (denominados xo e yo), os coficientes de distorção radial simétrica (k1, k2, k3) e os coeficientes de distorção descentrada (P1, P2), cujas componentes são divididas em tangencial e radial simétrica. Os coeficientes de distorção das lentes são obtidos devido às imperfeições na construção de lentes perfeitamente concôvas e convexas e a imprecisão no alinhamento dos eixos ópticos no momento da formação do sistema de lentes.” (adaptado) De acordo com Merchant (1979) o ponto principal pode ser definido como ponto principal de autocolimação (ppa) e ponto principal calibrado. O primeiro é definido como a intersecção de um raio de luz com o plano da fotografia, sendo perpendicular a este plano (da fotografia) antes de passar pelo sistema de lentes. Já o ponto principal calibrado é aquele ponto de melhor simetria das curvas de distorção.
3.4
AJUSTAMENTO
DE
OBSERVAÇÕES
POR
MÍNIMOS
QUADRADOS
UTILIZANDO O MÉTODO PARAMÉTRICO
Quando realizamos uma medida mais de uma vez, ela jamais repetirá o valor estabelecido. Por exemplo, ao medir uma dimensão linear de um objeto, mesmo que a medida seja efetuada por um instrumento com a graduação na ordem de 10³ unidades e o valor se repita todas as vezes, de fato, significa que as incertezas estão em todas as ordens de grandeza inferior, onde os valores jamais se repetirão.Segundo GEMAEL (1994) “As observações possuem uma propriedade inerente a elas, conhecida como flutuações probabilísticas, pois, quando se repete "n" vezes a medida de uma grandeza, os n valores não são idênticos, mas estão dispersos numa certa região ou intervalo. Essas flutuações, tradicionalmente, são classificadas como erros aleatórios e resultam na inconsistência do sistema de equações”.
8
Segundo SANTOS (2008) “O Método de estimação por Mínimos Quadrados (MMQ) tem como objetivo encontrar solução única para os parâmetros a serem estimados através da minimização da soma dos quadrados dos resíduos”. Este é o método de ajustamente geralmente empregado em ciências geodésicas. Segundo AMORIN (2005) “no ajustamento de observações denominam-se parâmetros (ou observações indiretas) as grandezas que não são obtidas diretamente, ou seja, aquelas que são calculadas em função de grandezas medidas diretamente. O método dos mínimos quadrados estima o valor dos parâmetros, dentre outros, usando o método de ajustamento paramétrico, a aplicação deste método exige a formulação das equações de observações que correspondem aos modelos
matemáticos
que
relacionam os
parâmetros
e as
observações”
(ADAPTADO). Assim, é por meio do modelo matemático, onde as observações são funções explícitas dos parâmetros, que todo o modelo paramétrico é baseado. Em modelos não lineares, a solução é repetidamente refinada por meio de processos iterativos, uma vez que não há solução para tais sistemas: daí vem a necessidade de linearizar o modelo matemático.
3.5 ORIENTAÇÂO EXTERIOR PELA RESSEÇÃO ESPACIAL
Segundo BRITO & COELHO (2002), há seis parâmetros que localizam a imagem no espaço: X0, Y0, Z0, ω,φ,κ. “X0“ representa a posição do centro de perspectiva no eixo X do sistema de coordenadas do espaço objeto. “Y0” representa a posição do centro de perspectiva no eixo Y do sistema de coordenadas do espaço objeto. “Z0” representa a posição do centro de perspectiva no eixo Z do sistema de coordenadas do espaço objeto. Os ângulos de Euler (ω,φ,κ) representam rotações sofridas pelo sistema local de coordenadas fotogramétrico xf, yf e zf (de cada câmara) em relação ao sistema global do terreno X, Y e Z. Rotacionando-se x, y, z de - ω, -,φ. e -k, pode-se torná-lo paralelo a X, Y, Z. Os mesmos autores indicam que A partir das equações de colinearidade, pode-se determinar os seis elementos de orientação exterior de uma fotografia, partir de, no mínimo três pontos de controle. Como os pontos de controle foram identificados na imagem, são conhecidas as suas “coordenadas pixel”. A partir dos parâmetros de orientação interior, chega-se às suas coordenadas no sistema
9
fotogramétrico. Como dados, também encontram-se disponíveis suas coordenadas tridimensionais (pois são pontos de controle, ou de campo) X1, Y1, Z1, X2, Y2, Z2, X3, Y3, Z3. A distância focal já é conhecida, uma vez que está contida no certificado de calibração da câmara. As coordenadas do ponto principal x0 e y0 (no espaço imagem) também encontram-se no mesmo, porém. Com todos estes valores conhecidos, resta apenas determinar as incógnitas, empregando-se um ajustamento pelo método paramétrico, esse processo é conhecido como resseção espacial.
3.6 MODELO DIGITAL DE TERRENO
Modelo digital de terreno ou D.T.M. ("digital terrain model") é o termo genérico empregado para referir-se ao modelamento matemático de superfícies. Podese definir modelo digital de terreno como sendo um conjunto de pontos amostrados da superfície real, com coordenadas espaciais (X,Y,Z) determinadas num dado referencial e algoritmos quê possibilitem construir um modelo matemático que reproduza da melhor maneira possível o comportamento altimétrico da superficie real (OSTMAN, 1987, p. 04).
10
4 MATERIAIS E MÉTODOS
Neste capítulo serão apresentadas as técnicas empregadas para se chegar aos objetivos propostos. Foi feita a seguinte subdivisão: apresentação dos dados fornecidos inicialmente; coleta dos pontos no referencial digital; conversão destes últimos para o referencial fotogramétrico; determinação dos POE pela resseção espacial; determinação das coordenadas no espaço-objeto pela monorrestituição; e implementação computacional.
4.1 DADOS INICIAIS
O presente trabalho, foi desenvolvido inteiro em ambiente computacional, sem emprego de qualquer procedimento de campo, assim sendo, o material necessário limita-se ao software empregado. Os dados são de uma fotografia que registra uma parte do bairro Nova Rússia de Ponta Grossa-PR. Para realizar os procedimentos demandados no capítulo 2, foi fornecida uma fotografia tomada por câmera digital não métrica, embarcada em uma aeronave. Foram fornecidos os parâmetros de orientação interior de tal sensor, assim como suas dimensões, informações, disponíveis na tabela 01. Tabela 01: Parâmetros de Orientação Interior e dimensões fornecidas para o sensor POI
valor
desvio-padrão
Distância Focal Calibrada
34.1450 mm
0,005 mm
Coordenada x0 do ponto principal
0,217mm
0,006 mm
Coordenada y0 do ponto principal
-0,089mm
0,005 mm
coeficiente K1 da distorção radial simétrica
-9,3625961e-05 mm-2
3,172455E-14 mm-6
coeficiente K2 da distorção radial simétrica
1,0931535e-07 mm-4
4,1924465E-14 mm-8
DIMENSÕES DO SENSOR Linhas
4500
Colunas
3000
Tamanho de cada célula (pixel)
0,008 mm
FONTE: fornecido pelo orientador da disciplina
11
Foi fornecida a altitude (contada do solo até o sensor) média de vôo: 1000 metros. A imagem fornecida, utilizada para a determinação das coordenadas fotogramétricas dos pontos no espaço-imagem, (empregadas em todos os procedimentos), tem sua miniatura disponibilizada na figura 01. Foi considerado que a precisão (desvio-padrão) de qualquer observação (fotocoordenada) é de meio pixel (0,004mm). Figura 01: Imagem aérea utilizada no trabalho
FONTE: fornecido pelo orientador da disciplina
12
Assim como também foram fornecidas as coordenadas dos pontos de apoio. Para as coordenadas “X” e “Y” em um “referencial terrestre”. de fato são coordenadas planas do sistema de projeção UTM, no fuso 22S associadas a uma altitude elipsoidal. As 3 componentes (E,N,h) estão referidas ao SIRGAS2000. Tais coordenadas são apresentadas na tabela 02. No anexo B está disponibilizado o croqui utilizado para localização dos pontos na imagem. Tabela 02: Coordenadas dos pontos de apoio Ponto
E (m)
N (m)
h (m)
1
583617.3780
7225141.4210
949.9260
2
584053.4540
7225220.5930
961.4010
3
584220.8590
7225799.7120
957.8920
4
583527.2300
7224979.6700
940.7000
5
583492.0540
7225529.0380
936.5200
6
583715.0300
7225763.4360
942.0150
7
583499.1680
7225362.3280
936.4570
8
583629.8380
7225202.2730
947.9810
9
583881.5810
7225582.5940
951.9500
10
583441.5770
7225098.1700
943.3080
11
583605.7010
7225599.1900
943.9770
12 583977.1040 7225142.0580 958.7370 FONTE: fornecido pelo orientador da disciplina
Foi fornecida também uma nuvem com 376853 pontos, que caracteriza a superfície da região onde a fotografia foi tomada. Todos os pontos tem suas coordenadas graduadas na casa do centímetro . Caso renderizada graficamente sob uma estrutura geralmente regular ou irregular (como uma TIN1), seria um Modelo Digital do Terreno. Uma vez que são muitos pontos, não será exibida uma tabela com estes. Com auxílio do software SURFER 12, foi renderizado o MDT, empregando a técnica da krigagem para interpolação da superfície, tal MDT encontra-se apresentado na figura 02.
1
Sigla em inglês para uma rede irregular de triângulos, cada qual contido em um plano.
13
Figura 02: MDT da região
FONTE: autoria própria (2015)
4.1 COLETA DOS PONTOS NO REFERENCIAL DIGITAL
Para se coletar as coordenadas no referencial digital é necessário algum software que seja capaz de ler o arquivo de imagem e retornar o índice do pixel selecionado (ou suas coordenadas Coluna e Linha). Há softwares que possibilitam apenas uma coleta manual, isto é, é selecionado (“clicado”) um pixel e se pode visualizar na tela seus valores, então procede-se com a anotação de tais valores. Tal
14
procedimento abre precedentes para erros de digitação, erros grosseiros que podem ser de difícil detecção. Assim, optou-se por um processo de coleta manual com anotação automatizada, para tanto, foi escolhido o software Quantum GIS, pelo fato de
ser
livre
e
de
fácil
manipulação.
Foi
simulado
um
processo
de
georreferenciamento de imagem, então foi gerado o arquivo contendo as coordenadas dos “pontos de controle” a partir do qual se extraiu a informação de interesse. Na figura 03 está o pedaço de uma captura de tela exemplificando tal uso, repare que foram inseridas falsas coordenadas para os “pontos de controle”.
Figura 03: emprego do QGIS na coleta das coordenadas no referencial digital
FONTE: Autoria Própria (2015)
Pode parecer estranho o emprego de um software de SIG para tal aplicação. Foram realizados testes que atestaram a qualificação do software para a atividade descrita. Assim, foram coletados os pontos no referencial digital. No software são fornecidos valores no subpixel (decimais) em um sistema dextrógiro, uma vez que o
15
interesse está no centro do pixel (a origem do sistema está no centro do primeiro pixel), foi feita a conversão para o referencial “coluna-e-linha” padrão, levógiro e apenas com números naturais. Ambos estão apresentados na tabela 03. Tabela 03: Coordenadas dos pontos de apoio na imagem Ref Digital Coletado
Ref Digital "Clássico"
Coluna
Linha
Coluna
Linha
2285.679
-3249.364
2285
3249
2487.450
-1444.497
2487
1444
419.493
-164.501
419
164
2796.514
-3766.461
2796
3766
651.788
-3269.156
651
3269
6.505
-2157.491
6
2157
1294.480
-3435.565
1294
3435
2062.497
-3130.604
2062
3130
867.487
-1712.520
867
1712
2242.382
-3963.643
2242
3963
501.478
-2761.528
501
2761
2708.439 -1838.612 2708 FONTE: Autoria Própria (2015)
1838
4.2 CONVERSÃO DO REFERENCIAL DIGITAL PARA O FOTOGRAMÉTRICO
Para ser possível empregar as equações de colinearidade, para os processos subsequentes da resseção espacial e da monorrestituição, é necessário que os pontos estejam em um referencial em unidades métricas, além de terem sido corrigidos de erros sistemáticos. Assim, foram empregados modelos matemáticos, para sequencialmente realizar as seguintes transformações: para o referencial com origem no centro da imagem; para o referencial fotogramétrico distorcido; e para o referencial fotogramétrico corrigido das distorções de lentes. Para realizar a conversão do referencial digital para o referencial com origem no centro da imagem, são empregadas as equações 01 e 02. As equações são apresentadas por SANTOS (2008), com alguns termos adaptados.
16
(
)
(01)
(
)
(02)
Onde: Xmm:
coordenada abcissa do ponto-imagem no novo referencial
Xmm:
coordenada ordenada do ponto-imagem no novo referencial
tp:
Tamanho do pixel em unidade métrica, que representa um
C:
Coluna do ponto em questão;
L
Linha do ponto em questão;
TC
Total de colunas da imagem;
TL
Total de linhas da imagem;
Sequencialmente, é necessário transladar a origem, do centro geométrico da imagem para o ponto principal desta, materializando o referencial fotogramétrico (ainda distorcido, junto com a consideração de “-f” como coordenada “z” fixa). Assim, basta apenas aplicar translações, como mostrado nas equações 03 e 04, disponíveis em SANTOS (2008).
x’ = xmm – x0
(03)
y’ = ymm – y0
(04)
Onde:
x’
coordenada x no espaço imagem no referencial fotogramétrico com distorções
y’
coordenada y no espaço imagem no referencial fotogramétrico com distorções
xmm
coordenada x no espaço imagem no referencial com origem no centro da imagem
ymm
coordenada y no espaço imagem no referencial origem no centro da imagem
17
O próximo passo é aplicar a correção devida a distorção radial simétrica, apresentadas de acordo com as equações 05 e 06.
xf = x’ – x’(k1r²+k2r4+k3r6)
(05)
yf = y’ – y’(k1r²+k2r4+k3r6)
(06)
Onde: r = (x’²+y’²)0,5
xf
coordenada x no referencial fotogramétrico.
yf
coordenada y no referencial fotogramétrico
x’
coordenada x no espaço imagem no referencial fotogramétrico com distorções
y’
coordenada y no espaço imagem no referencial fotogramétrico com distorções
k1, k2, k3,
Coeficientes da distorção radial simétrica
Para uma materialização mais rigorosa do referencial fotogramétrico, ainda haveriam correções que à serem feitas como a refração atmosférica e os parâmetros de afinidade do sensor digital, que para fins práticos serão desconsiderados. Vale lembrar que os valores fornecidos para a distorção descentrada, para este caso, são iguais a zero, logo não há porque, no contexto do presente trabalho realizar tal procedimento.
4.3 A RESSEÇÃO ESPACIAL
De posse das coordenadas da “imagem” dos pontos de apoio, já no referencial fotogramétrico corrigido dos erros sistemáticos tratados, juntamente com dados das tabelas 01 e 02, foi possível determinar os POE da imagem por meio da resseção espacial.
18
Foi determinado que dos 12 pontos de apoio disponíveis, fossem empregados apenas 5, então estes foram escolhidos de maneira semi-automática, considerando os 4 pontos com maiores e menores componentes x e y, e o ponto com menor variação nas distâncias para os quatro pontos. Para solução da resseção espacial, é necessária a resolução de um sistema de equações que tem as equações de colinearidade direta como modelo matemático. Tal equação baseia-se na condição de colinearidade, em que um vetorimagem2 é colinear a um vetor-objeto. Consiste em uma transformação entre dois espaços cartesianos tridimensionais e triortogonais ambos dextrógiros ou ambos levógiros, é feita por meio da aplicação de três rotações, três translações e um fator de escala (este último é substituido por termos equivalentes após uma manipulação algébrica). Suas formulações estão apresentadas nas equações 07 e 08.
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(
)
(07)
(08)
Onde: mij são elementos da matriz de rotação ou dos cossenos-diretores, apresentada na sequência:
Para resolver o sistema com solução única, são necessárias seis equações. Uma vez que cada ponto fornece duas equações, haverá um sistema com dez equações. Assim o sistema será redundante, com uma solução para cada trio de duplas possível, cada uma inconsistente com as demais devido aos erros aleatórios e a prevalência de erros sistemáticos residuais e mesmo não-modelados. Para determinar uma solução única, assim como a estimativa da qualidade de tal estimativa, é necessário um ajustamento de observações, que será realizado pelo usual método dos mínimos quadrados, que minimiza a soma do quadrado dos 2
Conforme abordado por SANTOS (2008)
19
resíduos. Foi adotada a simplificação de que os pontos de apoio de campo possuem desvio-padrão igual a zero, assim, a situação apresentada é propícia para a aplicação do modelo paramétrico, cujo modelo funcional, tem as observações ajustadas (La) como funções explícitas dos parâmetros ajustados (Xa), como demonstrado na equação 09
(
)
(09)
O ajustamento, é facilmente implementado computacionalmente, seguindo alguns passos (dos quais não são todos inerentes ao método paramétrico propriamente dito, mas necessário para os procedimentos do presente trabalho) apresentados na sequência. I.
É efetuada a leitura dos dados referentes a observações (fotocoordenadas), assim como dos demais dados: coordenadas dos pontos de apoio, distância focal calibrada... Todas os valores lineares são convertidos para a mesma unidade (metros);
II.
Em um vetor-coluna, geralmente denominado “Lb”, todas as observações (fotocoordenadas) são organizadas, na ordem: xi,yi, xi+1, yi+1...;
III.
É declarada a matriz dos pesos das observações geralmente denominada “P”;. Correspondente ao inverso da matriz de variância (de cada observação) e covariâncias (entre cada par de observações)- MVC. Seguindo ordem coerente com o vetor “Lb”, aqui todos os pontos possuíam mesmo peso, isto é 62500 (inverso do quadrado da precisão adotada);
IV.
É definido o modelo matemático a ser empregado (equações 07 e 08), nele, as observações, devem ser funções matemáticas explícitas, dos parâmetros;
V.
É
declarado
o
vetor
dos
parâmetros
aproximados,
genericamente
denominado “Xo”. As colunas 1 a 6 são preenchidas com os POE aproximados. Para X0, Y0 e Z0 foi adotada a média da componente correspondente dos pontos de apoio, para Z0 foi somada a altura de vôo fornecida; para ω e φ foram considerados nulos, já que se trata de um fotografia quase-nadiral; e κ foi obtido pela diferença de azimutes no espaçoimagem e no espaço-objeto; VI.
É estabelecida uma estrutura de repetição, uma vez que a partir do passo seguinte, será dado início ao processo iterativo, demandado devido a
20
natureza não-linear do modelo matemático (ver VIII). Antes de inicializar tal estrutura é estabelecido um contador de iterações e um critério de parada. VII.
São calculadas aproximações iniciais para as observações, utilizando o modelo matemático e os valores aproximados dos parâmetros obtidos na iteração anterior. Tais valores são alocados em um vetor-coluna geralmente denominado “Lo”. É a primeira etapa da linearização por série de Taylor 3;
VIII.
São calculadas as diferenças entre as observações (constantes) e suas aproximações da iteração corrente. Tais valores são alocados em um vetorcoluna geralmente denominado “L”;
IX.
São determinadas as derivadas parciais do modelo matemático, com relação aos parâmetros. Com as equações das mesmas, para cada observação, os valores numéricos de tais derivadas são calculados com os valores dos parâmetros da iteração anterior. Tais valores são alocados em uma matriz geralmente denominada “A”, com as linhas em ordem coerente com o vetor “Lb” e as colunas em ordem coerente com o vetor “Xo”. É a segunda etapa da linearização por série de Taylor;
X.
São resolvidas as equações normais do modelo paramétrico, geralmente denominadas “N” e “U”, demonstradas nas equações 11 e 12;
XI.
N= AT P A
(11)
U=AT P L
(12)
São calculadas as correções aos parâmetros ajustados da iteração anterior. São alocados em um vetor-coluna geralmente denominad0 “X”, apresentado na equação 13.
X = - (N)-1 U XII.
(13)
São calculados os parâmetros ajustados para a iteração corrente, aplicando a correção aos parâmetros da iteração anterior, conforme demonstrado na equação 14.
Xa=Xo+X XIII.
(14)
O vetor “Xo” assume o valor de “Xa” da iteração corrente, para servir como aproximação na iteração posterior, que será necessária caso o critério de
3
DALMOLIN (2002), p 94. Por não haver solução para sistemas de equações não-lineares, é necessário realizar uma aproximação linear, refinada a cada iteração.
21
parada não seja atingido, de modo que os procedimentos VI até este (XIII) são repetidos até que o critério seja atingido. O critério de parada aqui estabelecido, é que o maior valor (em módulo) do vetor das correções seja numericamente inferior a 10-10 ou até a ocorrência de 1000 iterações; XIV.
Com o fim do processo Iterativo, inicia-se a parte estocástica, com cálculo dos resíduos, alocados em um vetor geralmente denominado “V”, de mesma ordem de “Lb”, demonstrado equação 15;
V = AX + L XV.
(15)
Podem ser calculadas as observações ajustadas, alocadas em um vetor denominado “La”, somando-se os resíduos às observações, de acordo com a equação 16
La = Lb + V XVI.
(16)
É calculada estimativa da variância a posteriori, geralmente denotada por “vp” que é um escalar calculado pela equação 17
vp=VT P V /gl
(17)
Onde “gl”, são os graus de liberdade, valor correspondente ao número de observações redundantes, que pode ser calculado pela dimensão de diferença em módulo das dimensões da matriz A; Multiplicando-se “vp” por “gl”, obtem-se o X² calculado. XVII.
É calculada a estimativa da MVC dos parâmetros ajustados, alocadas em uma matriz geralmente denominada “MVCXa”, calculada de acordo com a equação 18;
MVCXa = vp(N)-1
(18)
Pode também ser estimada a MVC das observações ajustadas, alocadas em uma matriz geralmente denominada “MVCLa”, calculada de acordo com a equação 19;
MVCLa = A MVCXa A’ XVIII.
(19)
Da MVC dos parâmetros ajustados, são extraídos os desvios-padrões dos parâmetros ajustados, calculando raiz-quadrada dos elementos da diagonal principal de “MVCXa”, procedimento análogo pode ser feito pára a MVC das observações ajustadas; e
XIX.
Os resultados são convertidos para as unidades de interesse.
22
Assim, são determinados os parâmetros de orientação exterior da imagem.
4.4 A MONORRESTITUIÇÃO
De posse dos POE e dos pontos de interesse no referencial fotogramétrico, além da nuvem de pontos, é possível, determinar as coordenadas de feições pontuais presentes na imagem empregando o procedimento da monorrestituição. O modelo matemático fundamental empregado para a principal etapa de tal procedimento, são as Equações de Colinearidade Inversa, que com o conhecimento dos POE possibilitam a obtenção de coordenadas do espaço-objeto, desde que seja fornecido um valor para a coordenada Z. Tais equações, também estão baseadas na condição de colinearidade assim, vale reforçar que se faz necessária a correção das distorções, que afastam o raio de luz da condição de colinearidade, etapa já abordada em 4.2. Suas formulações, válidas para o diapositivo, como apresentado por SHENG (2005) estão apresentadas como as equações 20 e 21.
(
)
(
)
(
)
(
)
(
)
(
)
(20)
(21)
A monorrestituição, faz uso de um MDT para a obtenção da referida coordenada Z. Uma vez que o valor de X e Y consultado provavelmente não terá um correspondente no modelo digital, será necessária a aplicação de algum algoritmo de interpolação. Há diversos algoritmos, um dos mais empregados é o que faz o uso de uma TIN4. Para o presente trabalho, foi empregado a chamada interpolação pelo vizinho mais próximo, que consiste em, sem gerar novos valores para a variável a ser interpolada, consultar a entidade mais próxima para copiar dela tal valor. Para uma estrutura vetorial, com entidades do tipo ponto (como a nuvem disponibilizada), 4
Vide MITSHITA 1997.
23
o algoritmo consiste em encontrar o ponto menos distante do solicitado, por meio do cálculo de distâncias no plano XY daquele ponto para todos os disponíveis. Sua implementação em pseudocódigo está apresentada na sequência:
para i de 1 até número (n) de pontos da nuvem, armazenados em uma matriz n x 3 {
-Calcular distância plana entre oponto na posição i e o ponto de interesse, no plano XY; -Armazenar a distância na posição i de um vetor n x 1
}
encontrar em qual posição do vetor está contida a menor distância solicitar o ponto XYZ em tal posição guardar em uma variável o valor de Z de tal ponto
A monorrestituição, consiste na aplicação sequêncial das ECI seguidas do algoritmo de interpolação, por meio de um processo iterativo, parado quando um critério preestabelecido é atingido (foi adotado como menos de 10cm de diferença entre o ponto da iteração vigente com o da anterior). Na primeira iteração, é atribuído o valor de zero como primeira estimativa para Z. A monorrestituição é apresentada em pseudocógido na sequência:
Z=0; X=0; Y=0; declaração do valor numérico critério de parada para um ponto de interesse com coordenadas fotogramétricas conhecidas, enquanto critério de parada for falso {
-aplicação das ECI para obtenção de nova estimativa para X e Y; -aplicação do algoritmo de interpolação para obtenção de nova estimativa para Z; -teste da condição de parada {se verdadeiro terminar a estrutura de repetição}
}
Assim, repetindo o procedimento para uma lista de pontos de interesse é possível obter suas coordenadas tridimensionais no referencial terrestre, por meio do conhecimento de suas coordenadas no referencial fotogramétrico.
24
O procedimento da Monorrestituição está ilustrado na figura 04. Por meio da figura pode-se observar que a solução vai sendo refinada ao longo das iterações percorrendo um caminho de o que seria uma “espiral quadrada”. Figura 04: ilustração do processo da monorrestituição
Fonte: SHENG (2005)
4.5 A IMPLEMENTAÇÃO
Uma vez definidos os procedimentos, para realizá-los na prática, foi necessário implementá-los empregando uma linguagem de programação, uma vez que devido ao grande volume de operações e dados, torna-se inviável e completamente
sem
sentido
tentar
proceder
manualmente.
Dentre
várias
linguagens, foi optado pela utilização da linguagem matlab, uma linguagem rodada à partir de um interpretador homônimo. Algumas das vantagens da linguagem:
A sintaxe simplificada sob vários aspectos, sem por exemplo existir a necessidade de declarar o tipo de uma variável antes de declará-la;
A facilidade de efetuar cálculos matriciais, como por exemplo, realizar a inversão ou transposição de uma matriz;
A facilidade de importar e exportar dados;
A grande biblioteca de funções prontas, como size e length;
25
A simplicidade na hora de programar funções customizadas, podendo criar uma a partir da composição de outras, muito facilmente;
A existência de uma janela de comando, que permite executar apenas partes do código, facilitando os testes e depuração de erros;
Para maior organização, os objetivos abordados nas seções 4.2, 4.3 e 4.4 foram implementados,
cada
qual
em
um
script
diferente,
estando
disponíveis
respectivamente nos apêndices B, C e D. No apêndice A encontram-se todas as funções criadas para facilitar a solução de problemas repetidos e tornar o código mais elegante. Todos os códigos possuem comentários ao longo de seu corpo, de modo a facilitar o entendimento por parte do leitor. Antes de prosseguir, vale a pena comentar a respeito de artifícios interessantes empregados na implementação, visando otimizar o processo de solução dos problemas estipulados: a) Leitura dos dados a partir de um arquivo externo. Como na linha Possibilitando que sejam adicionados pontos de controle, caso necessário. Sendo possível manipular os dados com software de manipulação de planilhas, externo ao matlab; b) Leitura das dimensões das matrizes. Como em: “[m n1]=size(D)” e “u=length(POEi)” Assim, todo o código torna-se adaptável. De modo que se forem adicionados, ou excluídos pontos, nenhuma adaptação será necessária ao código; c) Escrita automatizada de matrizes, utilizando-se de for com duplo-contador, sendo um com incremento unitário e o outro não, assim, um percorre uma linha e outro mais de uma, em diferentes matrizes. Como no trecho que realiza a organização das observações no vetor Lb: “j=1; fori=1:2:(m*2) Lb(i,1)=xy(j,1); Lb(i+1,1)=xy(j,2); j=j+1; end” De modo a economizar laborioso trabalho manual, reduzindo a chance de ocorrerem erros de digitação. d) A utilização de funções anonymous. Como na linha: “eq=@(ma,mb,mc,X,X0,Y,Y0,Z,Z0) ma*(X-X0)+mb*(Y-Y0)+mc*(Z-Z0)”
26
Sendo próprias da linguagem, facilitam a implementação de funções e a sintaxe posterior. Possibilitam, incluir no corpo da função, variáveis declaradas anteriormente, coisa que não é possível em funções criadas em arquivos externos. e) Exportação dos resultados, com justaposição de matrizes/submatrizes, e operações realizadas diretamente. O que facilita o “manuseio” dos resultados obtidos.
27
5 RESULTADOS E ANÁLISE
Após a aplicação dos procedimentos relatados no capítulo anterior, foram obtidos os resultados que serão demonstrados neste capítulo. As coordenadas da tabela 03 transformadas para o referencial fotogramétrico, já corrigido dos erros sistemáticos abordados, que representam as fotocoordenadas dos pontos de apoio na imagem, estão apresentadas na tabela 04, juntamente com as coordenadas no referencial com origem no centro da imagem e no referencial fotogramétrico com distorções. Tabela 04: fotocoordenadas em referenciais métricos PONTO
ref Centro da Imagem
ref. Fotogramétrico dist.
ref Fotogrametrico
x (mm)
y (mm)
x' (mm)
y' (mm)
xf (mm)
yf (mm)
1
6.2840
-7.9960
6.0670
-7.9070
6.11687833
-7.97200543
2
7.9000
6.4440
7.6830
6.5330
7.74747371
6.58782321
3
-8.6440
16.6840
-8.8610
16.7730
-9.03410716
17.10067480
4
10.3720
-12.1320
10.1550
-12.0430
10.32257908
-12.24173510
5
-6.7880
-8.1560
-7.0050
-8.0670
-7.06988564
-8.14172269
6
-11.9480
0.7400
-12.1650
0.8290
-12.30493975
0.83853638
7
-1.6440
-9.4840
-1.8610
-9.3950
-1.87527096
-9.46704498
8
4.5000
-7.0440
4.2830
-6.9550
4.30766919
-6.99505936
9
-5.0600
4.3000
-5.2770
4.3890
-5.29899514
4.40729385
10
5.9400
-13.7080
5.7230
-13.6190
5.81013793
-13.82636177
11
-7.9880
-4.0920
-8.2050
-4.0030
-8.26279591
-4.03119708
12
9.6680
3.2920
9.4510
3.3810
9.52966433
3.40914137
FONTE: Autoria própria (2015)
Pode-se observar que as diferenças entre o referencial milimétrico original e o fotogramétrico corrigido podem estar na casa de vários pixels. Os parâmetros de orientação exterior calculados pelo ajustamento, assim como seus respectivos desvios-padrões estão apresentados na tabela 05. Tabela 05: fotocoordenadas em referenciais métricos POE
valor
X0
583809.6
desviopadrão 0.771539593
Y0
7225362
0.762449964
m
Z0
1993.221
0.123931955
m
ω
0.025623
0.00068265
rad
unidade m
28
0.01114
desviopadrão 0.000704398
-1.28527
0.000111809
rad
1.468076
0.039112975
graus
0.638294
0.040359054
graus
-73.6406 0.006406179 FONTE: autoria própria (2015)
graus
POE
valor
φ κ ω φ κ
unidade rad
Os valores para ω e φ foram determinados dentro do intervalo esperado, assim como todos os desvios-padrões estiveram na casa da subunidade. O ajustamento convergiu após 6 iterações. A MVC dos POE ajustados está apresentada na tabela 06 Tabela 06: MVC dos POE ajustados parâmetro
X0
Y0
Z0
ω
φ
κ
X0
0.5952733440
-0.3551055307
-0.0177146952
0.0003356802
0.0005418573
-0.0000075163
Y0
-0.3551055307
0.5813299482
0.0322169207
-0.0005188881
-0.0003417113
0.0000042160
Z0
-0.0177146952
0.0322169207
0.0153591296
-0.0000288102
-0.0000173259
0.0000002144
ω
0.0003356802
-0.0005188881
-0.0000288102
0.0000004660
0.0000003217
-0.0000000040
φ
0.0005418573
-0.0003417113
-0.0000173259
0.0000003217
0.0000004962
-0.0000000066
κ
-0.0000075163
0.0000042160
0.0000002144
-0.0000000040
-0.0000000066
0.0000000125
FONTE: autoria própria (2015)
As fotocoordenadas ajustadas, assim como os resíduos que as afastam das observações estão apresentadas na tabela 07. Tabela 07: fotocoordenadas ajustadas e os respectivos resíduos fotocoordenada ajustada
valor (mm)
resíduo (mm)
x1
0.7748022598
0.0000548890
y1
0.6587226879
-0.0000596332
x2
-0.9034496270
-0.0000389112
y2
1.7099611566
-0.0001063230
x3
1.0322083314
-0.0000495769
y3
-1.2246176509
-0.0004441412
x4
-1.2306666631
-0.0001726885
y4
0.0839711863
0.0001175485
x5
0.5812166923
0.0002028995
y5
-1.3821395429
0.0004966346
soma dos resíduos FONTE: autoria própria (2015)
0.0000006974
29
Em todas as fotocoordenadas, o resíduo ficou muito inferior à 1 pixel, de fato o valor de maior magnitude (para a coordenada y5) apresenta o valor de um vigésimo de pixel, indicando a ausência de erros grosseiros nas observações. O somatório dos resíduos é muito próximo a zero, o que indica a observância da distribuição normal nas observações, o que é um preceito que deve ser válido para o emprego do método dos mínimos quadrados seja a opção coerente. Os indicadores estatísticos variância a posteriori e o X² calculado, assim como os valores calculados para a quantidade de graus de liberdade a alguns intervalos de confiança, estão apresentados na tabela 08. Tabela 08: testes estatísticos para o ajustamento dos POE valor variância posteriori
0.860249952
Graus de liberdade
4
Ki-quadrado calculado
3.440999806
Ki-quadrado tab. 97.5%
11.1433
Ki-quadrado tab. 2.5%
0.4844
Ki-quadrado tab. 95%
9.4877
Ki-quadrado tab. 5%
0.7107
Ki-quadrado tab. 70%
4.8784
Ki-quadrado tab. 30% 2.1947 FONTE: autoria própria (2015)
O valor da variância a posteriori indica sua consistência como estimador imparcial da variância à priori. O teste estatístico, confirma a aceitação da hipótese básica com graus de significância de até 60%. Isto na prática indica a validade do modelo matemático empregado, assim como a ausência de erros grosseiros, que os pesos estão adequados e o sistema está bem condicionado5. Com a monorrestituição, foram determinadas as coordenadas de todos dos 12 pontos com valores nas últimas colunas da tabela 04. Os valores calculados para todos os pontos, estão apresentados na tabela 09, assim como a distância plana para o vizinho mais próximo e o número de iterações necessário para atingir o critério de parada.
5
segundo GEMAEL (1994)
30
Tabela 09: valores das coordenadas tridimensionais dos pontos de apoio determinadas pela monorrestituição PONTO
E (m)
N (m)
h (m)
dist. vizinho mais prox.
nº de iterações
1
583617.6447
7225141.3288
949.88
1.020
3
2
584053.4692
7225220.6123
961.39
1.549
3
3
584220.8941
7225799.7104
957.89
1.103
3
4
583527.3462
7224979.6687
940.63
2.163
4
5
583492.3409
7225528.7238
936.76
0.426
4
6
583715.0187
7225763.3374
942.11
0.117
4
7
583499.4807
7225361.8520
936.42
0.260
4
8
583630.0859
7225202.0175
947.92
0.286
3
9
583881.7632
7225582.1620
952.01
0.420
4
10
583441.6237
7225098.3373
943.91
28.576
3
11
583606.0473
7225598.9592
944.02
0.389
4
583977.2092 12 FONTE: autoria própria (2015)
7225142.0954
958.81
1.351
4
Aqui pode-se observar a rápida convergência do método, e o efeito do interpolador vizinho mais próximo, que não insere novos valores. O tempo total para processamento dos 12 pontos foi de 1.1324 segundo a função tic&toc do Matlab. As discrepâncias calculadas, entre as coordenadas fornecidas (tabela 02) e as calculadas (tabela 09), estão apresentadas na tabela 10. Tabela 10: discrepâncias entre as coordenadas determinadas e as fornecidas PONTO
ΔE (m)
ΔN (m)
Δh (m)
1
0.266738 0.015229 0.03508 0.116218 0.286868 -0.01128 0.312738 0.247912 0.182229 0.046723 0.346315 0.105178
-0.09215 0.019293 -0.00158 -0.00131 -0.31421 -0.09865 -0.47597 -0.25545 -0.43195 0.167274 -0.23082 0.037417
-0.046 -0.011 -0.002 -0.07 0.24 0.095 -0.037 -0.061 0.06 0.602 0.043 0.073
2 3 4 5 6 7 8 9 10 11
12 FONTE: autoria própria (2015)
31
Trata-se de valores pequenos, dentro do esperado. Para análise de tais discrepâncias, foram gerados alguns indicadores estatísticos simples, assim como foi calculado o erro médio quadrático (RMS) para cada uma das componentes. Tais valores estão apresentados na tabela 11. Tabela 11: análise das discrepâncias
em E (m) Mínimo -0.01128086 Máximo 0.34631483 Média 0.16249564 Mediana 0.14922373 somatório 1.94994770 RMS 0.20293501
em N (m) -0.47596868 0.16727391 -0.13984298 -0.09540017 -1.67811570 0.23776238
em h (m) -0.07000000 0.60200000 0.07383333 0.02050000 0.88600000 0.19408718
FONTE: autoria própria (2015)
É de fato pequena a amplitude das discrepâncias, havendo apenas um ponto onde o erro esteve acima do esperado. O somatório apresenta que as discrepâncias não seguem uma distribuição normal. O erro médio quadrático é um indicativo do grau de acurácia do método nas condições em que foi empregado, estando dentro do esperado. Assim sendo, o erro mínimo que se pode esperar em qualquer produto gerado por tal técnica nas condições deste trabalho, é de cerca de vinte centímetros.
32
7 CONCLUSÃO
Após a apresentação dos resultados obtidos e sua posterior análise, pode-se dizer que o trabalho proposto atingiu com êxito os objetivos estabelecidos. Uma vez que a resseção espacial ficou com qualidade além da esperada; na monorrestituição os valores ficaram dentro do esperado. A monorrestituição destaca-se como a técnica fotogramétrica mais acurada envolvendo apenas uma fotografia, sendo empregada em procedimentos como a ortorretificação direta de imagens. Mas é de fato viável para o mapeamento planimétrico e compilação de cartas e produtos cartográficos como os voltados ao cadastro. Assim,
há
no
apêndice
E,
um
exemplo
prático
do
emprego
da
monorrestituição para construção de mapas cadastrais. Uma vez que os dados já estão em um sistema de coordenadas UTM associados a alturas elipsoidais, os dados de saída já estarão nesse sistema. Foi produzida uma carta simplificada, mas que com a inserção de algumas informações faltantes, pode servir para diversas aplicações métricas. O processo é altamente automatizável, sendo possível, sim a conversão de grandes conjuntos de dados, empregando a monorrestituição para obter uma quantia densa de pontos.
33
REFERÊNCIAS
AMORIN, G.P. Curso de Formação Continuada em GeorreferenciamentoAplicado ao Cadastro Rural, Módulo 7: Ajustamento de Observações. CEFETES, coordenadoria de geomática, 2005. Disponível em: . Acesso em 25/04/2014
DALMOLIN, Q. Ajustamento por Mínimos Quadrados. Universidade Federal do Paraná, Curso de pós Graduação em Ciências Geodésicas. 2ª ed. Curitiba, 2002.
GEMAEL, C. Introdução ao ajustamento de observações: aplicações geodésicas. Curitiba, Editora da UFPR, 1994.
MITISHITA, E. A., Monorrestituição digital de fotos associada com sistema de computação gráfica CAD, para fins de mapeamento na área florestal. Curitiba, 1997. Tese de doutorado, Curso de Pós-Graduação em Engenharia Florestal, UFPR.
OSTMAN, A. A Quality Control and Accuracy Estimation of Digital Elevation Models. Papers on Digital Elevation, Image Correlation and Map Maintenance. Stockholm, Sweden, 1987. The Royal Institute of TechnologyDepartment of Photogrammetry SANTOS D. R. – Apostila de fotogrametria II– Departamento d e Geomática – Universidade Federal do Paraná – 1ª edição – Curitiba, 2008. SCURI A. E. –Fundamentos da Imagem Digital – Pontifícia Universidade Católica do Rio de Janeiro – Rio de Janeiro, 2002.
SHENG, Y. Theoretical Analysis of the Iterative Photogrammetric Method to Determining Ground Coordinates from Photo Coordinates and a DEM. American Society for Photogrammetry and Remote Sensing. Disponível em < http://www.asprs.org/a/publications/pers/2005journal/july/2005_july_863-871.pdf>. Acesso em 29/04/2015.
34
APÊNDICES
APÊNDICE A: FUNÇÕES IMPLEMENTADAS ........................................................ 35 APÊNDICE B: SCRIPT PARA CONVERSÃO PARA REF FOTOGRAMÉTRICO ... 37 APÊNDICE C: SCRIPT PARA A RESSEÇAO ESPACIAL ...................................... 39 APÊNDICE D: SCRIPT PARA REALIZAÇÃO DA MONORRESTITUIÇÃO ............ 43 APÊNDICE E: EXEMPLO PRÁTICO DE EMPREGO DA MONORRESTITUIÇÃO . 45
35
APÊNDICE A: FUNÇÕES IMPLEMENTADAS
Aqui serão apresentadas todas as funções desenvolvidas, que são empregadas nos módulos implementados. Para determinar a matriz de rotação a partir dos ângulos de Euler, foi criada a matriz função m_rot. Apresentada na sequência: %função para cálculo da matriz de rotação %na ordem ômega, phi e kappa function [M]=m_rot(o,p,k) m11= cos(p)*cos(k); m12= cos(o)*sin(k)+sin(o)*sin(p)*cos(k); m13= sin(o)*sin(k)-cos(o)*sin(p)*cos(k); m21= -cos(p)*sin(k); m22= cos(o)*cos(k)-sin(o)*sin(p)*sin(k); m23= sin(o)*cos(k)+cos(o)*sin(p)*sin(k); m31= sin(p); m32= -sin(o)*cos(p); m33= cos(o)*cos(p); M=[m11,m12,m13;m21,m22,m23;m31,m32,m33];
Para determinar o azimute fornecendo as coordenadas planas de dois segmentos de reta, foi criada a função Az_Ptos, Apresentada na sequência: function[Az]=az_ptos(x1,x2,y1,y2) vx=(x2-x1); vy=(y2-y1); %normalização do vetor: vx=vx/sqrt(vx^2+vy^2); vy=vy/sqrt(vx^2+vy^2); alfa=atand( abs(vx) / abs(vy) ); if vx > 0 & vy > 0 Az=alfa; elseif vx > 0 & vy < 0 Az=180-alfa; elseif vx < 0 & vy < 0 Az=180+alfa; elseif vx < 0 & vy > 0 Az=360-alfa; end
Para determinar os valores numéricos a partir da aplicação das derivadas parciais das equações de colinearidade direta, foi criada a função m_jacob. Apresentada na sequência: %função para cálculo da matriz jacobiana do modelo de colinearidade direto function [J]=mat_jacob(X0,Y0,Z0,o,p,k,X,Y,Z,f)
36
%elementos da matriz de rotação m11= cos(p)*cos(k); m12= cos(o)*sin(k)+sin(o)*sin(p)*cos(k); m13= sin(o)*sin(k)-cos(o)*sin(p)*cos(k); m21= -cos(p)*sin(k); m22= cos(o)*cos(k)-sin(o)*sin(p)*sin(k); m23= sin(o)*cos(k)+cos(o)*sin(p)*sin(k); m31= sin(p); m32= -sin(o)*cos(p); m33= cos(o)*cos(p); %variáveis auxiliares den= m31*(X-X0)+m32*(Y-Y0)+m33*(Z-Z0); num1= m11*(X-X0)+m12*(Y-Y0)+m13*(Z-Z0); num2= m21*(X-X0)+m22*(Y-Y0)+m23*(Z-Z0); %derivadas com relação à projeção do CP (para x) dxdX0=-f*(((den*-m11)-(num1*-m31))/(den^2)); dxdY0=-f*(((den*-m12)-(num1*-m32))/(den^2)); dxdZ0=-f*(((den*-m13)-(num1*-m33))/(den^2)); %derivadas com relação aos ângulos de euler (para x) dxdo=-f*((den*(((Y-Y0)*-m13)+((Z-Z0)*m12)))-(num1*(((Y-Y0)*-m33)+((Z-Z0)*m32)))) /(den^2); t1=((X-X0)*-sin(p)*cos(k))+((Y-Y0)*sin(o)*cos(p)*cos(k))+((Z-Z0)*-cos(o)*cos(p)*cos(k)); t2=((X-X0)*cos(p))+((Y-Y0)*(sin(o)*sin(p)))+((Z-Z0)*-cos(o)*sin(p)); dxdp=(-f*((den*t1)/(den^2)))+(f*((num1*t2)/(den^2))); dxdk=-f*((den*(((X-X0)*-cos(p)*sin(k))+((Y-Y0)*m22)+((Z-Z0)*m23)))/(den^2)); %dxdX=-dxdX0 ,dxdY=-dxdY0 ,dxdZ=-dxdZ0 %derivadas com relação à projeção do CP (para y) dydX0= -f*(((den*-m21)-(num2*-m31))/(den^2)); dydY0= -f*(((den*-m22)-(num2*-m32))/(den^2)); dydZ0= -f*(((den*-m23)-(num2*-m33))/(den^2)); %derivadas com relação aos ângulos de euler (para x) dydo=-f*((den*(((Y-Y0)*-m23)+((Z-Z0)*m22)))-(num2*(((Y-Y0)*-m33)+((Z-Z0)*m32)))) /(den^2); t3= ((X-X0)*sin(p)*sin(k))+((Y-Y0)*-sin(o)*cos(p)*sin(k))+((Z-Z0)*cos(o)*cos(p)*sin(k)); t4= ((X-X0)*cos(p))+((Y-Y0)*(sin(o)*sin(p)))+((Z-Z0)*-cos(o)*sin(p)); dydp= (-f*((den*t3)/(den^2))) + (f*((num2*t4)/(den^2))); dydk= -f*((den*(((X-X0)*-cos(p)*cos(k))+((Y-Y0)*-m12)+((Z-Z0)*-m13)))/(den^2)); %dydX=-dydX0 ,dydY=-dydY0 ,dydZ=-dydZ0 J=[ dxdX0,dxdY0,dxdZ0,dxdo,dxdp,dxdk; dydX0,dydY0,dydZ0,dydo,dydp,dydk];
37
APÊNDICE B: SCRIPT PARA CONVERSÃO PARA REF FOTOGRAMÉTRICO
Aqui será apresentado o código utilizado para implementar o descrito no capítulo 4.2. clear clc format long %script para geração das fotocoordenadas corrigidas para um conjunto de %pontos na imagem, todas as unidades métricas estão em milímetros %leitura da fotografia aérea F=imread('054.tif'); [totLin totCol totBan]=size(F) %leitura do arquivo contendo as coordenadas coletadas no QGIS P=csvread('fotocoord.csv') [m,n]=size(P) %modificando para o referencial Coluna e Linha (C,L) padrão P1=P; P=floor(abs(P)) %definição do tamanho do pixel tp=0.008 %convertendo para o referencial com origem no centro da imagem xy0=zeros(m,2); for i=1:m xy0(i,1)=tp * (P(i,2)-((totCol-1)/2)); xy0(i,2)=-tp * (P(i,3)-((totLin-1)/2)); end xy2=xy0 %transladando a origem para o ponto principal x0=0.217; y0=-0.089; xy0(:,1)=xy0(:,1)-x0; xy0(:,2)=xy0(:,2)-y0; %aplicando a correção para o deslocamento devido a distorção radial simétrica k1=-9.3625961e-05 k2=1.0931535e-07 k3=0 xy=zeros(m,2); for i=1:m r=sqrt(xy0(i,1)^2+xy0(i,2)^2); xy(i,1)=xy0(i,1)-xy0(i,1)*(k1*r^2+k2*r^4+k3*r^6); xy(i,2)=xy0(i,2)-xy0(i,2)*(k1*r^2+k2*r^4+k3*r^6); end
38
%exportando os resultados dlmwrite('fotocoord5.csv',xy,'precision',50); %para o matlab dlmwrite('fotocoord5b.csv',[xy2,xy0,xy],'precision',50) dlmwrite('refDigital.csv',[P1,P],'precision',50); %csvwrite('fotocoord2.csv',xy) %para o freemat
39
APÊNDICE C: SCRIPT PARA A RESSEÇAO ESPACIAL
Aqui será apresentado o código utilizado para implementar o descrito no capítulo 4.3. clear clc format long %orientação exterior pela resseção espacial %leitura dos dados (todos os pontos): D1=csvread('fotocoord5.csv') [m1 n1]=size(D1) %leitura do tamanho da matriz %seleção semi-automatizada de 5 pontos bem distribuídos pts=zeros(1,5); %os 4 primeiros pontos serão os pontos extremos da imagem Xmax=max(D1(:,2));pts(1)=find(D1(:,2)==Xmax) Xmin=min(D1(:,2));pts(2)=find(D1(:,2)==Xmin) Ymax=max(D1(:,3));pts(3)=find(D1(:,3)==Ymax) Ymin=min(D1(:,3));pts(4)=find(D1(:,3)==Ymin) %será matematicamente definido o ponto mais ao centro distancias1=zeros(m1,m1); for i=1:m1 x=D1(i,2); y=D1(i,3); for j=1:m1 distancias1(i,j)=sqrt( (D1(j,2)-x)^2+(D1(j,3)-y)^2 ) end end %seprarado apenas as distancias para os pontos de interesse distancias=[distancias1(:,pts(1)),distancias1(:,pts(2)),distancias1(:,pts(3)),distancias1(:,pts(4))] distancias=distancias' %o quinto elemento, será o ponto cujas distâncias para os 4 pontos de %interesse sejam mais próximas entre si --> menor desvio padrão pts(5)=find(std(distancias)==min(std(distancias))); %organizando apenas para conservar a ordem original pts=sort(pts) %separando de todos os pontos, os 5 escolhidos for i=1:length(pts) D(i,:)=D1(pts(i),:); end [m n2]=size(D) %dlmwrite('pts_escol.csv',D,'precision',50) %separação e conversão dos dados %coordenadas no espaço-imagem: xy=D(:,2:3) %(observações) xy=xy./1000 %em metros plot(xy(:,1),xy(:,2),'*') axis equal %coordenadas no espaço-objeto (já em metros): XYZ=D(:,4:6)
40
%declaração de constantes f=34.1450 %distância focal da câmera em milímetros f=f/1000 %em metros % %precisões %desvio padrão de todas as coord. xy: (na mesma unidade) dpxy=0.004/1000 % %PESOS %peso de todas as coord. xy: pxy=1/(dpxy^2) %matriz dos pesos das observações P=eye(m*2)*pxy % %POE iniciais: POSi=mean(XYZ) h0=1000 %kappa inicial: %azimutes do primeiro ponto para um arbitrário, como o quarto: azp1pnA=az_ptos(xy(4,1),xy(1,1),xy(4,2),xy(1,2)); azp1pnB=az_ptos(XYZ(4,1),XYZ(1,1),XYZ(4,2),XYZ(1,2));
%ordem empregada: X0,Y0,Z0,omega,phi,kappa Xo=[POSi(1),POSi(2),POSi(3)+h0,0,0,(azp1pnA-azp1pnB)*(pi/180)]' u=length(Xo) %coordenadas no espaço-objeto organizadas em um vetor-linha XYZv=zeros(m*3,1) j=1; for i=1:3:m*3 XYZv(i,1)=XYZ(j,1); XYZv(i+1,1)=XYZ(j,2); XYZv(i+2,1)=XYZ(j,3); j=j+1; end %organização das observações (xy) no vetor Lb: j=1; for i=1:2:(m*2) Lb(i,1)=xy(j,1); Lb(i+1,1)=xy(j,2); j=j+1; end n=length(Lb) %número de observações %modelo matemáticos empregado: %as equações de colinearidade eq=@(ma,mb,mc,X,X0,Y,Y0,Z,Z0) ma*(X-X0)+mb*(Y-Y0)+mc*(Z-Z0) x=@(den,num1) -f*(num1/den) y=@(den,num2) -f*(num2/den) %contador de iterações e auxiliar de condição de parada itera=-1 fim=false %predeclaração das matrizes a serem utilizadas Lo=zeros(n,1); A=zeros(n,u);
41
%início do processo iterativo while (~fim) itera = itera+1 %matriz de rotação para a iteração corrente M=m_rot(Xo(4),Xo(5),Xo(6)); %observações recalculadas, com os parametros da iteração anterior j=1; for i=1:2:n num1=eq(M(1,1),M(1,2),M(1,3),XYZv(j),Xo(1),XYZv(j+1),Xo(2),XYZv(j+2),Xo(3)); num2=eq(M(2,1),M(2,2),M(2,3),XYZv(j),Xo(1),XYZv(j+1),Xo(2),XYZv(j+2),Xo(3)); den= eq(M(3,1),M(3,2),M(3,3),XYZv(j),Xo(1),XYZv(j+1),Xo(2),XYZv(j+2),Xo(3));
Lo(i,1) =x(den,num1); Lo(i+1,1)=y(den,num2);
j=j+3; end %vetor contendo as diferenças entre Lo e Lb L=Lo-Lb; %matriz das derivadas parciais (linearização do modelo não-linear, por série de taylor) com relação aos parâmetros j=1; for i=1:2:n J=mat_jacob(Xo(1),Xo(2),Xo(3),Xo(4),Xo(5),Xo(6),XYZv(j),XYZv(j+1),XYZv(j+2),f); A(i,1)=J(1,1); A(i,2)=J(1,2); A(i,3)=J(1,3); A(i,4)=J(1,4); A(i,5)=J(1,5); A(i,6)=J(1,6);
A(i+1,1)=J(2,1); A(i+1,2)=J(2,2); A(i+1,3)=J(2,3); A(i+1,4)=J(2,4); A(i+1,5)=J(2,5); A(i+1,6)=J(2,6); j=j+3; end
%As equações normais do modelo paramétrico: N=A'*P*A; U=A'*P*L; %cálculo do vetor correção para a iteração corrente X=-inv(N)*U; %cálculo dos parâmetros ajustados para a iteração corrente
42
Xa=Xo+X; %para início da próxima iteração Xo=Xa; %critérios de parada fim= (max(abs(X))) 1000) end %vetor dos resíduos: V=A*X+L %observações ajustadas: La=Lb+V %cálculo da variância a posteriori gl=length(Lb)-u vp=(V'*P*V)/gl %teste estatístico de qui-quadrado kquad=vp*gl % %MVC dos parâmetros ajustados MVCXa=vp*inv(N) %MVC das observações ajustadas MVCLa=A*MVCXa*A' %parâmetros de orientação exterior ajustados: POE=Xa; %com o,p, e k em graus decimais POEg=POE; POEg(4:6,1)=POEg(4:6,1)*180/pi %seus respectivos desvios padrões dpPOE=sqrt(diag(MVCXa)) %com dpo,dpp e dpk em graus dpPOEg=dpPOE; dpPOEg(4:6,1)=dpPOEg(4:6,1)*180/pi %estimativa dos pesos adequados P2=P*(1/vp); dpxy2=sqrt(1/P2(1,1))*1000; % %exportação dos resultados dlmwrite('POE.csv',POE,'precision',50) dlmwrite('MVCXa.csv',MVCXa,'precision',50) dlmwrite('poe_dppoe.csv',[POE dpPOE POEg dpPOEg],'precision',50) dlmwrite('la_v.csv',[La*100 V*100],'precision',50) dlmwrite('MVCpoe.csv',MVCXa(1:6,1:6),'precision',50) dlmwrite('kiquadrado.csv',[vp;gl;kquad],'precision',50) dlmwrite('MVCLa.csv',MVCLa,'precision',50)
43
APÊNDICE D: SCRIPT PARA REALIZAÇÃO DA MONORRESTITUIÇÃO
Aqui será apresentado o código utilizado para implementar o descrito no capítulo 4.4. clear clc format long % %MONORRESTITUIÇÂO % %leitura de dados: %Parâmetros de Orientação Exterior da imagem, na ordem: X0,Y0,Z0,omega,phi,kappa POE=csvread('POE.csv'); M=m_rot(POE(4),POE(5),POE(6))%matriz de rotação %nuvem de pontos da região: N=csvread('nuvem.csv'); [m n]=size(N) %Coordenadas dos pontos cujas coordenadas no ref. terrestre serão calculadas % e comparadas com as determinadas por técnicas geodésicas P=csvread('fotocoord5.csv') [m2 n2]=size(P) x=P(:,2)./1000;%separação das fotocoordenadas x e conversão para metros y=P(:,3)./1000;%separação das fotocoordenadas y e conversão para metros % %declaração de constantes f=34.1450 %distância focal da câmera em milímetros f=f/1000 %em metros %modelo matemático: as equações de colinearidade inversa eq=@(ma,mb,mc,x,y) ma*x+mb*y+mc*(-f) X=@(Z,X0,Z0,den,num1) X0+(Z-Z0)*(num1/den) Y=@(Z,Y0,Z0,den,num2) Y0+(Z-Z0)*(num2/den) %inicialização do processo iterativo Zit=0;Xit=0;Yit=0; Xit0=0;Yit0=0; %criação da matriz para guardar os valores calculados Pmono=zeros(m2,5); tic %para cômputo do tempo de execução %a monorrestituição dos 12 pontos for i=1:m2 Zit=0 itera=-1; dif=1000; while(dif>0.01) % while(itera
Lihat lebih banyak...
Comentários
