Distância da superfície no DEM entre duas coordenadas

11

Dadas duas coordenadas em um modelo de elevação digital (DEM), como computo a distância real percorrida entre os dois locais, assumindo uma rota em linha reta? Como a alteração da resolução do meu DEM afetará os resultados?

scw
fonte

Respostas:

10

Se você não possui esse recurso incorporado ao seu GIS, mas pode executar algumas operações básicas de grade ("álgebra de mapa"), ainda existe uma solução.

O cálculo se resume a encontrar a inclinação da rota em todos os pontos. Se você soubesse exatamente , sem nenhum erro de discretização, integrava a secante da inclinação. Em uma grade, a integral é estimada obtendo a média do secante para as células interceptadas pela rota e multiplicando a média pelo comprimento da rota. (No mapa de álgebra, isso seria uma "média zonal" multiplicada pelo comprimento da rota).

A inclinação da rota não é igual à inclinação do DEM! Depende exatamente de como a rota corta a superfície. Portanto, você precisa de informações completas sobre a "direção" da superfície, que pode ser descrita em termos de batida e inclinação, inclinação e aspecto, ou por meio de um vetor normal unitário ( ou seja , um campo vetorial 3D perpendicular à superfície). A maneira mais confiável é reduzir o problema para um onde você conheça o campo vetorial normal. Isso significa que você tem um triplo de números em cada célula - representados como três grades separadas, é claro - que eu chamarei (Nx, Ny, Nz). A direção da rota (no plano) pode ser representada como um vetor de unidade (x, y, t) onde (x, y) indica sua direção no mapa. O valor de t é o "aumento" na direção vertical:a taxa na qual a rota deve subir para permanecer na superfície . Assim, como a velocidade 2D da rota - sua "corrida" - é igual a Sqrt (x ^ 2 + y ^ 2), a inclinação é dada por

(1) tan (declive) = subida / corrida = t / Sqrt (x ^ 2 + y ^ 2) .

Nos cálculos t será uma grade, mas o denominador, Sqrt (x ^ 2 + y ^ 2), é apenas um número. Se você calculá-lo usando a fórmula (4) abaixo, será igual a 1, para que você possa esquecer: t será a tangente da grade de inclinação da rota e sec (inclinação) = sqrt (1 + t ^ 2) será a grade cuja média zonal você calcula.

É fácil encontrar t. Por definição, o vetor de direção (x, y, t) é perpendicular ao vetor normal. Isso significa

0 = x * Nx + y * Ny + t * Nz, então

(2) t = - (x * Nx + y * Ny) / Nz .

No cálculo, Nx, Ny e Nz são grades, mas x e y são números. Portanto t é uma grade, como pretendido. (Não haverá nenhum problema com a divisão, porque não é possível para Nz = 0: seria um penhasco perfeitamente vertical, que não pode ser representado em um DEM.)

Então: como você encontra o vetor normal (Nx, Ny, Nz) e o vetor de direção (x, y)? Normalmente, um GIS computa as grades de inclinação (s) e de aspecto (a) do DEM. Expresse cada um como um ângulo. Estas são basicamente coordenadas esféricas para o vetor normal da unidade. Para aspectos a leste do norte, a unidade normal é obtida pela conversão de coordenadas esférica-cartesiana usual,

(3) (Nx, Ny, Nz) = (pecado (s) * pecado (a), pecado (s) * cos (a), cos (s)) .

Nesse cálculo, s e a são grades , portanto, ele descreve três expressões de álgebra de mapa separadas para criar três grades Nx, Ny e Nz.

Como verificação, observe que quando a inclinação é zero (s = 0), o vetor normal é (0,0,1), apontando para cima, como deveria. Quando o aspecto é zero, o vetor normal é (0, sen (es), cos (s)) que evidentemente aponta para o norte (direção y) e se inclina da vertical por um ângulo de s, o que implica que a superfície se inclina de a horizontal por um ângulo de s: essa é realmente a sua inclinação.

Finalmente, deixe o rumo da rota ser b (um ângulo constante, a leste do norte). Seu vetor de direção é

(4) rolamento = (x, y) = (sin (b), cos (b)).

Observe que o rumo é um par de números , não um par de grades, porque descreve a direção da rota.


À medida que a resolução do DEM aumenta, é possível observar mais variações locais nas encostas, fazendo com que a encosta estimada aumente, como observa @johanvdw. Estudei esse fenômeno, aumentando sucessivamente os DEMs de alta resolução e comparando os DEMs de uma área obtida de fontes diferentes. Descobri que em áreas de alta inclinação as diferenças nas estimativas de inclinação podem ser substanciais . Isso se traduzirá em diferenças substanciais nas estimativas de comprimento de rota terrestre. Caso contrário, em áreas com declives uniformes, as diferenças podem não ter conseqüências.

Uma maneira de avaliar o efeito da resolução para o seu DEM é realizar um estudo semelhante. É preciso pouco esforço. Por exemplo, estime o comprimento terrestre de uma rota usando o DEM e, em seguida, reestime o comprimento após agregar esse DEM em 2 x 2 blocos (aumento por um fator de 2). Se houver uma diferença inconseqüente entre as duas estimativas, você deve estar bem; se a diferença for importante, pode valer a pena obter um DEM de resolução mais precisa para o seu trabalho. (Existem métodos mais sofisticados para melhorar as estimativas de inclinação e comprimento explorando o DEM que você possui, mas levaria muito tempo para descrevê-las aqui.)

whuber
fonte
7

O SAGA GIS possui um módulo para isso: Perfil interativo

http://www.saga-gis.org/saga_modules_doc/ta_profiles/index.html

Os pontos resultantes conterão a distância e a distância terrestre. Se o DEM tiver uma resolução mais grossa, sua distância por terra será sempre um pouco menor (a menos que você tenha condições de fronteira estranhas), mas, na realidade, essa diferença provavelmente não é importante. Se a área for bastante plana, até a distância terrestre e a distância normal serão quase as mesmas: se a inclinação entre dois pontos ao longo da linha for de 20%, a distância terrestre será apenas 2% maior que a distância normal (sqrt ( 1 ^ 2 + 0,2 ^ 2) = 1,019).

johanvdw
fonte
O link acima parece estar quebrado, mas aqui estão as instruções sobre como criar um perfil de terreno interativo no SAGA.
perfil completo de cengel