Como interpolar planos no MATLAB?

Confira nosso curso online de QGIS

Quando se fala em MATLAB, sempre nos lembramos dos seus belos gráficos. Nesta postagem, você irá aprender a interpolar dois planos utilizando o MATLAB e criar um belo gráfico.

Quando se fala em MATLAB, sempre nos lembramos daqueles gráficos em 3 dimensões que enchem os olhos de qualquer pessoa. Porém, executar tais gráficos no MATLAB não é fácil.

Estarei mostrando aqui uma abordagem para interpolar (linearmente) dois planos/superfícies no MATLAB, por exemplo, interpolar duas profundidades amostradas em um solo para avaliar a concentração de um determinado contaminante.

Visualização do gráfico e código desta postagem.
Visualização do gráfico e código desta postagem.

Nosso script irá iniciar com as funções clc e clear. Tais funções servem para limpar/apagar as variáveis existentes, além de limpar a janela de comando.

Posteriormente, será necessário carregar os nossos dados. Esses são originados da grade de pontos da nossa outra postagem. Caso você não tenha visualizado ela antes, dê uma olhada antes de prosseguir.

Aqui, converti o nosso arquivo .txt da postagem passada em .csv.

Para realizar isso, basta abrir seu programa de edição de planilhas (e.g. Excel ou Calc), e abrir o arquivo .txt. Um assistente de importação irá aparecer, clique em Avançar na primeira etapa, e na segunda etapa, nos delimitadores, marque a vírgula também (o quadro visualização de dados irá criar uma coluna para cada dado).

Na terceira etapa, caso o seu separador decimal seja um ponto, será necessário altera-lo clicando no botão Avançado e trocar o Separador Decimal por uma vírgula e o Separador de Milhar por um ponto. Depois é só salvar como um arquivo .csv.

O código abaixo irá carregar nossos dados.


filename = 'C:\Users\usuario\Desktop\idw_fishnet.csv';
dados = readtable(filename, 'Delimiter', ';');

Lembre-se de substituir o caminho “C:\Users\[….]” pelo caminho onde esta o seu arquivo. Em seguida, a função readtable() irá fazer a leitura do nosso arquivo csv.

Em seguida, vamos indicar o número de pontos que cada lado de nossa grade contém (npontos), além de readequar os dados às duas profundidades que estamos trabalhando e criar a grade de coordenadas.


npontos = 100;
C0 = reshape(dados.idw_pol50, [npontos, npontos]);
C1 = reshape(dados.idw_pol100, [npontos, npontos]);
[X, Y, Z] = meshgrid(1:1:npontos, 1:1:npontos, -50:-1:-100);

A função reshape() irá adequar nossos dados, que antes estavam dispostos em apenas 4 colunas (FID, id, concentração na profundidade 50 e 100cm), em uma grade com diversas colunas e linhas, onde cada célula é um ponto “amostrado”.

Perceba que o MATLAB faz referencias às nossas colunas colocando um ponto após o nome da nossa variável, seguida do nome da coluna.

Para descobrir o nome das colunas basta digitar na linha de comando dados.Properties.VariableNames.

A função meshgrid() é similar à anterior, porém cria diversas tabelas que são agrupadas formando um “cubo”.

Após criar a grade de pontos no MATLAB, é preciso criar a grade com os valores das concentração do contaminante que você mediu. Para tal, iremos interpolar os valores entre as nossas duas variáveis C0 e C1. Essa interpolação ocorrerá pela seguinte fórmula:

[Ci] = Valor Mínimo + ((Valor Máximo – Valor Mínimo) / Nº de Planos) * i)

Onde i será substituído por um valor a cada iteração/loop. No nosso caso, cada iteração é uma profundidade diferente, valor que inclusive acompanha a variável nplanos (número de planos).


nplanos = 51;
C = zeros(npontos, npontos, nplanos);
for m = 1:npontos
  for n = 1:npontos
    for i = 1:nplanos
      [Ci] = C0(m,n) + ((C1(m,n)-C0(m,n))/nplanos)*i;
      C(m,n,i) = Ci;
    end
  end
end

São três loops interligados, onde o primeiro irá passar por todos os pontos no eixo X, o segundo nos pontos no eixo Y e o terceiro nos pontos do eixo Z.

Dessa forma, podemos criar nosso gráfico em 3 dimensões com os comandos abaixo. As variáveis xslice, yslice e zslice  definem onde cortam os planos que você quer demonstrar, e aqui, vamos mostrar a terça parte do plano X (por isso dividimos por 3), não mostraremos o plano Y (por isso ele é igual a 0) e o plano Z vamos apresentar as profundidades -50 e -100cm.


xslice = npontos/3;
yslice = 0;
zslice = [-50, -100];
h = slice(X, Y, Z, C, xslice, yslice, zslice);
set(h, 'EdgeColor', 'none')
colormap jet

A função colormap indica as cores do seu gráfico, sendo possível utilizar: jet, parula, hsv, hot, entre outros. A partir disso, é só inserir alguns detalhes como rótulos nos eixos e adicionar um título.


zlabel('Profundidade (cm)');
ylabel('Eixo Y');
xlabel('Eixo X');
ylabel(colorbar, 'Concentração Poluente (ppm)');

Basta rodar o modelo clicando em Run e o nosso gráfico irá aparecer. Algumas alterações podem ser realizadas, como essas abaixo, para obter gráficos diferentes.

xslice = npontos/2;
yslice = npontos/2;
zslice = [-75];
Uma das modificações que podem ser realizadas no nosso gráfico.
Uma das modificações que podem ser realizadas no nosso gráfico.

Qualquer dúvida ou dificuldade, deixe seu comentário que estaremos respondendo.



Clique na figura abaixo e assine nossa lista de emails para receber nosso ebook "Como criar mapas de localização com ArcGIS 10.x".

Apostila Mapa de Localização Banner

Author: Fernando BS

Engenheiro Ambiental e de Segurança do Trabalho. Atua nas áreas de geoprocessamento, mineração e hidrologia. Busca soluções utilizando softwares como QGIS, R e Python.

Deixe um comentário

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *