Leitura e tratamento das instâncias do Metrofor

Dados e Pacotes

Atilizaremos os pacotes: geosphere, RgoogleMaps e leaflet.

install.packages("geosphere")
install.packages("RgoogleMaps")
install.packages("leaflet")

library(geosphere)
library(RgoogleMaps)
library(leaflet)

A função distm() da biblioteca geosphere, recebe um vetor de coordenadas de origem, e um vetor de coordenadas de destino e retorna a matrix de distâncias entre os dois conjuntos de pontos. Por exemplo:

library(geosphere)

## Loading required package: sp

distm(c(-38.53453523, -3.72019331), c(-38.53494246, -3.72004400), fun = distHaversine)

##          [,1]
## [1,] 48.19396

Carregando o conjunto de paradas no objeto estacoes.

setwd("C:/Users/Lucas/Documents")

estacoes = read.csv("Estacoes Metrofor.txt", dec = ".", sep = ",", header = TRUE)

head(estacoes)

##                           id       lon       lat
## 1          Estação Parangaba -38.56367 -3.776055
## 2            Estação Montese -38.55284 -3.773497
## 3         Estação Vila União -38.53412 -3.769700
## 4     Estação Borges de Melo -38.52872 -3.760782
## 5 Estação São João do Tauape -38.51207 -3.760971
## 6      Estação Pontes Vieira -38.50021 -3.751978

Carregando o conjunto de rotas no objeto rota.

setwd("C:/Users/Lucas/Documents")

rota = read.csv("Trejeto metrofor.txt", dec = ".", sep = "," , header = TRUE)

head(rota)

##         lon       lat
## 1 -38.56411 -3.776212
## 2 -38.56214 -3.774653
## 3 -38.56198 -3.774554
## 4 -38.56178 -3.774458
## 5 -38.56158 -3.774388
## 6 -38.56138 -3.774335

Matriz de distâncias

Agora iremos criar uma matriz ‘D’ de distâncias das estações.

##          [,1]     [,2]     [,3]     [,4]     [,5]     [,6]
## [1,]    0.000 1236.512 3358.663 4238.166 5973.442 7542.114
## [2,] 1236.512    0.000 2122.549 3029.840 4739.079 6318.238
## [3,] 3358.663 2122.549    0.000 1159.337 2634.968 4251.848
## [4,] 4238.166 3029.840 1159.337    0.000 1850.505 3315.767
## [5,] 5973.442 4739.079 2634.968 1850.505    0.000 1654.462
## [6,] 7542.114 6318.238 4251.848 3315.767 1654.462    0.000

De maneira idêntica, definimos uma matriz ‘M’ de distâncias para as rotas.

##          [,1]      [,2]      [,3]      [,4]      [,5]      [,6]
## [1,]   0.0000 279.50030 300.07117 324.83019 347.39736 368.58161
## [2,] 279.5003   0.00000  20.67966  45.95494  69.36237  91.63072
## [3,] 300.0712  20.67966   0.00000  25.37194  48.90685  71.31661
## [4,] 324.8302  45.95494  25.37194   0.00000  23.58925  46.08158
## [5,] 347.3974  69.36237  48.90685  23.58925   0.00000  22.52207
## [6,] 368.5816  91.63072  71.31661  46.08158  22.52207   0.00000

Construir a matriz de distâncias usando a função distHaversine()

xlon = cbind(estacoes$lon)
xlat = cbind(estacoes$lat)

ylon = cbind(estacoes$lon)
ylat = cbind(estacoes$lat)

dh <- function(xlon, xlat, ylon, ylat){
  tmp = matrix(1,nrow(xlon),nrow(ylon))
  k = 1
  for (i in 1:nrow(xlon)) {
    for(j in 1:nrow(ylon)){
      tmp[i,j] = distHaversine(c(xlon[i],xlat[i]),c(ylon[j],ylat[j]))
      k = k+1
    }
  }
  return(tmp)
}

head(dh(xlon,xlat,ylon,ylat)[,1:6])

##          [,1]     [,2]     [,3]     [,4]     [,5]     [,6]
## [1,]    0.000 1236.512 3358.663 4238.166 5973.442 7542.114
## [2,] 1236.512    0.000 2122.549 3029.840 4739.079 6318.238
## [3,] 3358.663 2122.549    0.000 1159.337 2634.968 4251.848
## [4,] 4238.166 3029.840 1159.337    0.000 1850.505 3315.767
## [5,] 5973.442 4739.079 2634.968 1850.505    0.000 1654.462
## [6,] 7542.114 6318.238 4251.848 3315.767 1654.462    0.000

Funções

A função inrote() verifica se um conjunto de pontos pertence ou não a rota.

lon = cbind(rota$lon)
lat = cbind(rota$lat)

elon = cbind(estacoes$lon)
elat = cbind(estacoes$lat)

inrote<-function(lon,lat,elon,elat){
  x = matrix(c(lon,lat),nrow = length(lon), ncol = 2)
  y = matrix(c(elon,elat),nrow = length(elon),ncol = 2)
  DM = distm(x,y,fun = distHaversine)
  return(DM)
}

Plotando os mapas interativos

Agora plotamos o mapa:

library(leaflet)

leaflet(data = estacoes) %>% addTiles() %>%
  addMarkers(~lon, ~lat, popup = ~as.character(id), label = ~as.character(id))

Vizualização alternativa com marcadores pontos.

m <- leaflet(estacoes) %>%
addTiles() %>%
addCircles(~lon, ~lat, popup=NA, weight = 3, radius=100, 
           color="black", stroke = TRUE, fillOpacity = 0.8) %>% 
addLegend("bottomright", colors= "#ffa500", labels="", title="Metrofor")%>%
  setView (lat=mean(estacoes$lat), lng = mean(estacoes$lon), zoom=12)
m

Também podemos incluir uma legenda e alterar a cor de cada linha do metrofor.

m <- leaflet(estacoes) %>%
  addTiles() %>%
  addCircles(~lon[1:10], ~lat[1:10], popup=estacoes$id[1:10], weight = 3, radius=100, color="black", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[11:20], ~lat[11:20], popup=estacoes$id[11:20], weight = 3, radius=100, color="red", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[21:40], ~lat[21:40], popup=estacoes$id[21:40], weight = 3, radius=100, color="blue", stroke = TRUE, fillOpacity = 0.8) %>% 
  addLegend("bottomright", colors= c("blue","red","black"), labels= c("Linha Sul","Linha Leste","VLT"), title="Metrofor")%>%
  setView (lat=mean(estacoes$lat), lng = mean(estacoes$lon), zoom=12)
m

Através do endereço: http://leaflet-extras.github.io/leaflet-providers/preview/index.html Podemos obter uma enorme variedade de estilos para o mata, por exemplo, visualização no modo satélite.

m <- leaflet(estacoes) %>%
addProviderTiles("Esri.WorldImagery") %>%
 addCircles(~lon[1:10], ~lat[1:10], popup=estacoes$id[1:10], weight = 3, radius=100, color="black", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[11:20], ~lat[11:20], popup=estacoes$id[11:20], weight = 3, radius=100, color="red", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[21:40], ~lat[21:40], popup=estacoes$id[21:40], weight = 3, radius=100, color="blue", stroke = TRUE, fillOpacity = 0.8) %>% 
  addLegend("bottomright", colors= c("blue","red","black"), labels= c("Linha Sul","Linha Leste","VLT"), title="Metrofor")%>%
  setView (lat=mean(estacoes$lat), lng = mean(estacoes$lon), zoom=12)
m

Por fim, um mapa topogáfico:

m <- leaflet(estacoes) %>%
addProviderTiles("Esri.WorldTopoMap") %>%
 addCircles(~lon[1:10], ~lat[1:10],popup=estacoes$id[1:10], weight = 3, radius=100, color="black", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[11:20], ~lat[11:20], popup=estacoes$id[11:20], weight = 3, radius=100, color="red", stroke = TRUE, fillOpacity = 0.8) %>% 
  addCircles(~lon[21:40], ~lat[21:40], popup=estacoes$id[21:40], weight = 3, radius=100, color="blue", stroke = TRUE, fillOpacity = 0.8) %>% 
  addLegend("bottomright", colors= c("blue","red","black"), labels= c("Linha Sul","Linha Leste","VLT"), title="Metrofor")%>%
  setView (lat=mean(estacoes$lat), lng = mean(estacoes$lon), zoom=12)
m

Algo interessante també é poder filtar os marcadores por tipo, além de incluir a opção de exibir o rastro do caminho:

map <- leaflet(estacoes) %>%
  # baseGroups
  addTiles(group = "Visualização padrão") %>%
  addProviderTiles("Esri.WorldTopoMap", group = "Topográfico") %>%
  addProviderTiles("Esri.WorldImagery", group = "Satélite") %>%
  # overlayGroups
  addCircles(~lon[1:10], ~lat[1:10], 10,popup=estacoes$id, color="black", stroke = TRUE, 
             fillOpacity = 0.8, group = "VLT") %>%
  addCircles(~lon[11:20], ~lat[11:20], 10,popup=estacoes$id, color="black", stroke = TRUE, 
             fillOpacity = 0.8, group = "linha leste") %>%
  addCircles(~lon[21:40], ~lat[21:40], 10,popup=estacoes$id, color="black", stroke = TRUE, 
             fillOpacity = 0.8, group = "linha sul") %>%
  addPolygons(data = estacoes, lng = ~lon, lat = ~lat,
              fill = F, weight = 2, color = "black", group = "rotas") %>%
  # Layers control
  addLayersControl(
    baseGroups = c("Visualização padrão", "Topográfico", "Satélite"),
    overlayGroups = c("VLT", "linha leste", "linha sul","rotas"),
    options = layersControlOptions(collapsed = FALSE)
  )
map

Leitura e tratamento das instâncias do Metrofor

Lucas Araújo da Silva

5 de setembro de 2017

Dados e Pacotes

Matriz de distâncias

Funções

Plotando os mapas interativos