Difference between revisions of "Opasnet map"

From Testiwiki
Jump to: navigation, search
m (Rationale)
(See also)
 
(10 intermediate revisions by 2 users not shown)
Line 1: Line 1:
 
[[Category:Opasnet]]
 
[[Category:Opasnet]]
 
[[Category:Code under inspection]]
 
[[Category:Code under inspection]]
{{variable|moderator=Jouni|stub=Yes}}
+
[[Category:R tool]]
 +
[[Category:Modelling]]
 +
{{method|moderator=Jouni|stub=Yes}}
  
 
== Question ==
 
== Question ==
Line 24: Line 26:
  
 
Key guidance:
 
Key guidance:
 +
* [[OpasnetUtils/Drafts]] contains many functions needed.
 +
* [[Help:Drawing graphs]] contains guidance for making non-map graphs.
 
* [http://cran.r-project.org/web/packages/rgdal/rgdal.pdf rgdal package for R]
 
* [http://cran.r-project.org/web/packages/rgdal/rgdal.pdf rgdal package for R]
 
** [http://www.gdal.org/ GDAL - Geospatial Data Abstraction Library]
 
** [http://www.gdal.org/ GDAL - Geospatial Data Abstraction Library]
Line 31: Line 35:
  
 
== Rationale ==
 
== Rationale ==
 +
 +
===Wiki maps===
 +
 +
{{#display_map: 62.91826544557644,28.69525909423828~ ~Valuma-alue; 62.92647012476198,28.637065887451172
 +
|rectangles=62.94068612542196,28.726844787597656:62.92373548704899,28.702983856201172~Luikonlahden kaivos~ ~ ~ ~ ~#F2E1E1
 +
|zoom=10
 +
|centre=62.91826544557644,28.69525909423828
 +
}}
 +
 +
 +
 +
{{#display_map:
 +
62.928986, 28.706245;
 +
62.937734, 28.712769
 +
}}
 +
 +
[http://semantic-mediawiki.org/wiki/Special:MapEditor Map editor and instructions]
  
 
=== Calculations ===
 
=== Calculations ===
  
 
==== Kuopio buildings with Google pin map ====
 
==== Kuopio buildings with Google pin map ====
 +
 +
{{attack|#|In function google.point_kml, parameter name is not used for anything!|--[[User:Jouni|Jouni]] ([[User talk:Jouni|talk]]) 08:42, 26 December 2013 (EET)}}
 +
 +
{{defend|#|What we need is a function that eats SpatialPointDataFrames and gives KML files with the data points in different formats (pins, coloured marks etc.|--[[User:Jouni|Jouni]] ([[User talk:Jouni|talk]]) 08:42, 26 December 2013 (EET)}}
  
 
<rcode name='kuorakonmaps'>
 
<rcode name='kuorakonmaps'>
Line 45: Line 70:
 
library(OpasnetUtilsExt)
 
library(OpasnetUtilsExt)
  
objects.latest("Op_en6007", code_name = "answer") # [[OpasnetUtils/Drafts]]; we need MyPointsKML and ova2spat.
+
shp<-readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house') # Read building data
 
 
shp<-readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house')
 
  
 
shp <- shp[1:100 , ] # Save only 100 buildings for demonstration
 
shp <- shp[1:100 , ] # Save only 100 buildings for demonstration
Line 53: Line 76:
 
plotvar <- shp@data$ika
 
plotvar <- shp@data$ika
 
nclr <- 8
 
nclr <- 8
plotclr <- brewer.pal(nclr,"BuPu")
+
plotclr <- brewer.pal(nclr, "BuPu")
class <- classIntervals(plotvar,nclr,style="quantile")
+
class <- classIntervals(plotvar, nclr, style = "quantile")
colcode <- findColours(class,plotclr)
 
  
 
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
 
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
Line 61: Line 83:
 
shp2 <- spTransform(shp, epsg4326String)
 
shp2 <- spTransform(shp, epsg4326String)
  
kmlname <- "Kuopio house data"
+
dat <- google.point_kml(
kmldescription <- "Random stuff about here"
+
shp2,
icon <- "http://maps.google.com/mapfiles/kml/pal2/icon18.png"
+
kmlname = "Kuopio house data",
name <- paste("ika value: ", shp2$ika)
+
kmldescription = "Random stuff about here",
description <- paste("<b>Value:</b>",shp2$ika,"<br><b>Description:</b>",shp2$kayttotark)
+
name = paste("ika value: ", shp2$ika),
 
+
description = paste("Age:", shp2$ika, "<br>Description:", shp2$kayttotark),
dat <- MyPointKML(shp2, kmlname, kmldescription, name, description, icon, colcode)
+
icon = "http://maps.google.com/mapfiles/kml/pal2/icon18.png",
 +
col = findColours(class,plotclr)
 +
)
  
 
google.show_kml_data_on_maps(dat)
 
google.show_kml_data_on_maps(dat)
  
 
</rcode>
 
</rcode>
 
  
 
==== GoogleMaps Sorvi MML ====
 
==== GoogleMaps Sorvi MML ====
 
{{attack|#|Something wrong with sorvi.|--[[User:Jouni|Jouni]] ([[User talk:Jouni|talk]]) 12:54, 23 December 2013 (EET)}}
 
  
 
<rcode name='gmapspsqltest3'>
 
<rcode name='gmapspsqltest3'>
Line 84: Line 105:
 
library(sorvi)
 
library(sorvi)
 
library(rgdal)
 
library(rgdal)
 +
library(maptools)
  
data(MML)
+
shp <- LoadData("kuntarajat.maa.shp")
 
 
shp <- MML[["1_milj_Shape_etrs_shape"]][["kunta1_p"]]
 
  
 
#epsg3857String <- CRS("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
 
#epsg3857String <- CRS("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
Line 94: Line 114:
 
shp2<-spTransform(shp,epsg4326String)
 
shp2<-spTransform(shp,epsg4326String)
  
out<-sapply(slot(shp2,"polygons"),function(x
+
out<-sapply(slot(shp2,"polygons"),function(x)
 
{kmlPolygon(
 
{kmlPolygon(
 
x,
 
x,
Line 202: Line 222:
  
 
google.circles(data)
 
google.circles(data)
 +
 +
</rcode>
 +
 +
====Rasters on Google maps====
 +
 +
<rcode graphics=1>
 +
library(OpasnetUtils)
 +
library(OpasnetUtilsExt)
 +
library(rgdal)
 +
library(maptools)
 +
library(RColorBrewer)
 +
library(classInt)
 +
library(raster)
 +
 +
altiste <- "Ni"
 +
 +
pitoisuus <- new("ovariable",
 +
name = "pitoisuus",
 +
data = tidy(op_baseGetData("opasnet_base", "Op_fi3280", series_id = 5509), objname= "pitoisuus")
 +
)
 +
 +
# Remove empty rows.
 +
pitoisuus@data <- pitoisuus@data[pitoisuus@data$pitoisuusResult != "" , ]
 +
 +
# Convert values that are below detection limit to ranges 0 - detection limit.
 +
pitoisuus@data$pitoisuusResult <- gsub("<", "0-", pitoisuus@data$pitoisuusResult)
 +
 +
# Convert different time formats used to POSIXct time format.
 +
temp <- as.POSIXct(pitoisuus@data$Aika, format = "%d.%m.%Y %H:%M")
 +
temp <- ifelse(is.na(temp), as.POSIXct(pitoisuus@data$Aika, format = "%d.%m.%Y %H.%M"), temp)
 +
temp <- ifelse(is.na(temp), as.POSIXct(pitoisuus@data$Aika, format = "%d.%m.%Y"), temp)
 +
temp <- as.POSIXct(temp, origin = "1970-01-01 02:00.00 UTC")
 +
pitoisuus@data$Aika <- temp
 +
 +
pitoisuus <- EvalOutput(pitoisuus, N = 1)
 +
 +
data <- pitoisuus@output[pitoisuus@output$Altiste == altiste, ]
 +
data.time <- data
 +
 +
########### Concentrations on map
 +
 +
locs <- data.frame(Paikka = c(
 +
"Lumijärvi",
 +
"Ylä-Lumijärvi",
 +
"Kivijärvi",
 +
"Lumijoki, Kivijärveä lähin piste",
 +
"Lumijoki, keskimmäinen piste",
 +
"Lumijoki, Lumijärveä lähin piste",
 +
"Lampi kaivosalueen reunalla 1",
 +
"Lampi kaivosalueen reunalla 2",
 +
"Puro Kalliojärven ja Salmisen välillä",
 +
"Tuhkajoki",
 +
"Kortelampi",
 +
"Kärsälampi",
 +
"Rengaskaivo, Pappila",
 +
"Kalliojärvi",
 +
"Kalliojoki",
 +
"Lumijoki",
 +
"Salmisenpuro"
 +
),
 +
LA = c(
 +
63.942533,
 +
63.945662,
 +
63.931893,
 +
63.925618,
 +
63.936103,
 +
63.945605,
 +
63.986605,
 +
63.992419,
 +
64.004996,
 +
64.044465,
 +
63.950091,
 +
63.992843,
 +
63.968536,
 +
64.012646,
 +
64.018907,
 +
63.935284,
 +
64.003882
 +
),
 +
LO = c(
 +
27.936602,
 +
27.971964,
 +
27.907677,
 +
27.929521,
 +
27.926474,
 +
27.953424,
 +
27.971706,
 +
27.983251,
 +
28.002756,
 +
28.104916,
 +
27.984538,
 +
27.984581,
 +
27.953897,
 +
28.001447,
 +
28.022261,
 +
27.937717,
 +
27.999601
 +
))
 +
 +
oprint(locs)
 +
 +
data <- merge(data, locs)
 +
result <- "pitoisuusResult"
 +
 +
title <- paste(altiste, "Talvivaaran läheisyydessä (", data$pitoisuusUnit[1], ")", sep = "")
 +
coordinates(data)=c("LO","LA")
 +
 +
# Plot the data
 +
proj4string(data)<-("+init=epsg:4326")
 +
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
 +
shp<-spTransform(data,epsg4326String)
 +
 +
#Create blank raster
 +
rast<-raster()
 +
 +
#Set raster extent to that of point data
 +
extent(rast)<-extent(shp)
 +
 +
#Choose number of columns and rows
 +
ncol(rast) <- 32
 +
nrow(rast) <- 32
 +
 +
#Rasterize point data
 +
rast2<-rasterize(shp, rast, shp[[result]], fun=mean)
 +
 +
start <- 0 # min(shp[[result]])
 +
end <- max(shp[[result]])
 +
 +
steps <- approx(c(start,end),n=6)$y
 +
colors <- rev(rainbow(length(steps), start=0, end=0.50))
 +
 +
par(mfrow=c(6,1), mar=c(3,1,0,1), cex=1.5)
 +
 +
colorstrip <- function(colors, labels)
 +
{
 +
count <- length(colors)
 +
m <- matrix(1:count, count, 1)
 +
image(m, col=colors, ylab="", axes=FALSE)
 +
axis(1,approx(c(0, 1), n=length(labels))$y, labels)
 +
}
 +
 +
cat("<span style='font-size: 1.2em;font-weight:bold;'>", title, "</span>\n")
 +
 +
colorstrip(colors, steps)
 +
 +
#Plot data
 +
 +
google.show_raster_on_maps(rast2, col=colors, style="height:500px;")
 +
 +
ggplot(data.time, aes(x = Paikka, y = pitoisuusResult, colour = pitoisuusResult)) +
 +
geom_point(size = 5) + # point size
 +
theme_grey(base_size = 24) + # font etc. http://sape.inf.usi.ch/quick-reference/ggplot2/themes
 +
labs(
 +
title = paste(altiste, "Talvivaaran läheisyydessä"),
 +
y = paste("Pitoisuus,", data$pitoisuusUnit[1]),
 +
x = "Mittauspaikka"
 +
)
 +
 +
########### Concentrations in time
 +
 +
ggplot(data.time, aes(x = Aika, y = pitoisuusResult, colour = Paikka)) +
 +
geom_point(size = 5) +
 +
theme_grey(base_size = 24) +
 +
labs(
 +
title = paste(altiste, "Talvivaaran läheisyydessä"),
 +
y = paste("Pitoisuus,", data$pitoisuusUnit[1]),
 +
x = "Aika"
 +
)
  
 
</rcode>
 
</rcode>
Line 234: Line 422:
  
 
* [https://gist.github.com/2292662 Plotting sp-objects to Google Maps with different CRSs]
 
* [https://gist.github.com/2292662 Plotting sp-objects to Google Maps with different CRSs]
 +
* [http://semantic-mediawiki.org/wiki/Special:MapEditor Map editor]
 +
* [http://www.findlatitudeandlongitude.com/batch-geocode/ Batch geocoding]: give addresses, the website will provide geolocations. Works for a few hundred addresses.
  
 
==Keywords==
 
==Keywords==

Latest revision as of 11:29, 1 November 2016



Question

How should GIS data be handled and visualised in Opasnet?

Answer

  • Original GIS data is stored in a PostgreSQL database.
  • Data is accessed by R using the RGDAL package and ?? connection.
  • Data in manipulated in R.
  • Ovariables are converted to SpatialPointsDataFrame objects using ova2spat function.
  • Data is displayed by producing a KML file with ?? package.
    • The KML file is created with MyPointKML function, if pins are shown.
    • The KML file is saved at the R-tools server.
    • Google Maps is used to show the KML file on a web page.

Currently, ova2spat and MyPointKML functions are located in OpasnetUtils/Drafts. Therefore, you need this command:

objects.latest("Op_en6007", code_name = "answer") # OpasnetUtils/Drafts. We need MyPointsKML and ova2spat.

This is the projection that the National Land Survey Finland uses: [1].

Key guidance:


Rationale

Wiki maps

Loading map...


Loading map...

Map editor and instructions

Calculations

Kuopio buildings with Google pin map

#: In function google.point_kml, parameter name is not used for anything! --Jouni (talk) 08:42, 26 December 2013 (EET)

#: What we need is a function that eats SpatialPointDataFrames and gives KML files with the data points in different formats (pins, coloured marks etc. --Jouni (talk) 08:42, 26 December 2013 (EET)

+ Show code

GoogleMaps Sorvi MML

+ Show code

Google with shapefiles

+ Show code

Google show data from url on map

+ Show code

Google circles

+ Show code

Rasters on Google maps

+ Show code

Static GoogleMaps

+ Show code

See also

Keywords

References


Related files

<mfanonymousfilelist></mfanonymousfilelist>