Opasnet map: Difference between revisions

From Opasnet
Jump to navigation Jump to search
(first draft)
 
 
(20 intermediate revisions by 3 users not shown)
Line 1: Line 1:
[[Category:Opasnet]]
[[Category:Opasnet]]
{{variable|moderator=Jouni|stub=Yes}}
[[Category:Code under inspection]]
[[Category:R tool]]
[[Category:Modelling]]
{{method|moderator=Jouni|stub=Yes}}


== Question ==
== Question ==
Line 11: Line 14:
* Data is accessed by [[R]] using the RGDAL package and ?? connection.
* Data is accessed by [[R]] using the RGDAL package and ?? connection.
* Data in manipulated in [[R]].
* Data in manipulated in [[R]].
* Ovariables are converted to SpatialPointsDataFrame objects using ova2spat function.
* Data is displayed by producing a [[:en:Keyhole markup language|KML]] file with ?? package.
* Data is displayed by producing a [[:en:Keyhole markup language|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.
** The KML file is saved at the R-tools server.
** Google Maps is used to show the KML file on a web page.
** 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: [http://spatialreference.org/ref/epsg/3067/].
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://www.gdal.org/ GDAL - Geospatial Data Abstraction Library]
* [http://cran.r-project.org/web/packages/sp/vignettes/intro_sp.pdf sp package for R]: contains descriptions about the key spatial objects.
* [http://rwiki.sciviews.org/doku.php?id=tips:spatial-data:change_crs Changing the coordinate system of Spatial objects]


== Rationale ==
== Rationale ==


All pieces of the puzzle exist already.
===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 ===


==== Kuopio buildings with Google pin map ====


=== Dependencies ===
{{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)}}


=== Formula ===
{{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'>
 
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(OpasnetUtils)
library(OpasnetUtilsExt)
 
shp<-readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house') # Read building data
 
shp <- shp[1:100 , ] # Save only 100 buildings for demonstration
 
plotvar <- shp@data$ika
nclr <- 8
plotclr <- brewer.pal(nclr, "BuPu")
class <- classIntervals(plotvar, nclr, style = "quantile")
 
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(shp)<-("+init=epsg:3067")
shp2 <- spTransform(shp, epsg4326String)
 
dat <- google.point_kml(
shp2,
kmlname = "Kuopio house data",
kmldescription = "Random stuff about here",
name = paste("ika value: ", shp2$ika),
description = paste("Age:", shp2$ika, "<br>Description:", shp2$kayttotark),
icon = "http://maps.google.com/mapfiles/kml/pal2/icon18.png",
col = findColours(class,plotclr)
)
 
google.show_kml_data_on_maps(dat)
 
</rcode>
 
==== GoogleMaps Sorvi MML ====
 
<rcode name='gmapspsqltest3'>
 
library(OpasnetUtils)
library(OpasnetUtilsExt)
library(sorvi)
library(rgdal)
library(maptools)
 
shp <- LoadData("kuntarajat.maa.shp")
 
#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")
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(shp)<-("+init=epsg:3047")
shp2<-spTransform(shp,epsg4326String)
 
out<-sapply(slot(shp2,"polygons"),function(x)
{kmlPolygon(
x,
name="nimi",
col='#df0000aa',
lwd=1,
border='black',
description="selite"
)}
)
 
data <- paste(
paste(
kmlPolygon(
kmlname =" This will be layer name",
kmldescription = "<i>More info about layer here</i>"
)$header,
collapse="\n"
),
paste(unlist(out["style",]), collapse="\n"),
paste(unlist(out["content",]), collapse="\n"),
paste(kmlPolygon()$footer, collapse="\n"),
sep=''
)
 
google.show_kml_data_on_maps(data)
 
</rcode>
 
==== Google with shapefiles ====
 
<rcode name='gmapspsqltest2'>
 
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(OpasnetUtils)
library(OpasnetUtilsExt)
 
shp <- readOGR('PG:host=localhost user=postgres dbname=spatial_db','watson_wkt')
 
plotvar <- shp@data$value_inhalation
nclr <- 8
plotclr <- brewer.pal(nclr,"BuPu")
class <- classIntervals(plotvar,nclr,style="quantile")
colcode <- findColours(class,plotclr)
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(shp)<-("+init=epsg:3035")
shp2 <- spTransform(shp,epsg4326String)
 
 
out<-sapply(
slot(shp2,"polygons"),
function(x){
kmlPolygon(
x,
name = as(shp2,"data.frame")[slot(x,"ID"),"country_code"],
col = colcode[[((as.numeric(slot(x,"ID"))+1))]],
lwd = 1,
border = 'black',
description = paste("Value:",as(shp2,"data.frame")[slot(x,"ID"),"value_inhalation"])
)
}
)
 
data <- paste(
paste(
kmlPolygon(
kmlname = "This will be layer name",
kmldescription = "<i>More info about layer here</i>")$header,
collapse="\n"
),
paste(unlist(out["style",]), collapse="\n"),
paste(unlist(out["content",]), collapse="\n"),
paste(kmlPolygon()$footer, collapse="\n"),
sep=''
)
 
cat("PostgreSQL test with coloured data areas.\n")
 
google.show_kml_data_on_maps(data)
 
cat("PostgreSQL Test with default data\n")
 
google.show_data_on_maps()
 
</rcode>
 
==== Google show data from url on map ====
 
<rcode name='gutilstest'>
library(OpasnetUtils)
library(OpasnetUtilsExt)
google.show_kml_on_maps(url='http://api.flickr.com/services/feeds/geo/?g=322338@N20&lang=en-us&format=feed-georss')
</rcode>
 
==== Google circles ====
 
<rcode name='gmaptest'>
library(OpasnetUtilsExt)
 
data = list()
data[[1]] = c(62.8925, 27.678333, 9683, '#ff0000', 'kuopio')
data[[2]] = c(65.016667, 25.466667, 14398, '#00ff00', 'oulu')
data[[3]] = c(60.170833, 24.9375, 59623, '#0000ff', 'helsinki')
 
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>
 
==== Static GoogleMaps ====
 
<rcode name='staticmapstest' graphics=1 embed=1>
 
library(RgoogleMaps)
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(OpasnetUtils)
library(OpasnetUtilsExt)
 
objects.latest("Op_en6007", code_name = "answer") # [[OpasnetUtils/Drafts]]. We need MyRmap.
 
shp <- readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house')
proj4string(shp) <- ("+init=epsg:3067") # The map projection of NLS Finland.
 
#oprint(shp@data[1:100 , ])
 
epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
shp2 <- spTransform(shp,epsg4326String) # Convert to longitude-latitude projection.
 
MyRmap(shp2, plotvar = "ika", legend_title = "Ikä", numbins = 8, pch = 19, cex = 0.3) # Draw the map.
 
</rcode>


==See also==
==See also==
* [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

{{#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 }}

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) (type: truth; paradigms: science: attack)

←--#: . 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) (type: truth; paradigms: science: defence)

+ 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>