+ Show code- Hide code
library(OpasnetUtils)
library(ggplot2)
library(xtable)
library(OpasnetUtilsExt)
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(raster)
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
))
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"
)
| |
