+ Näytä koodi- Piilota koodi
library(OpasnetUtils)
library(OpasnetUtilsExt)
library(xtable)
library(ggplot2)
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(raster)
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)
}
data <- tidy(op_baseGetData("opasnet_base", "Op_fi3129"))
for(i in 1:nrow(data)) {
temp <- as.character(data$Result[i])
if(!data$Unit[i] %in% c("-", "Teksti", "Lista")) {temp <- as.numeric(temp)}
if(data$Unit[i] == "-") {temp <- strsplit(temp, ",")[[1]]}
assign(as.character(data$Lyhenne[i]), temp)
}
cat("Louhintaprosessi ", louhinta, " (ton per a) synnyttää pienhiukkaspääston kohdassa LO ", LO, ", LA ", LA, "\n")
M <- louhinta / 365
A <- 1
objects.get("8EV2qsAPaXuk7ELI")
MetMalMurskPol <- EvalOutput(MetMalMurskPol)
cat("Hiukkaspäästo:\n")
print(xtable(MetMalMurskPol@output), type = 'html')
riippuvuudet <- data.frame(Name = "MetMalMurskPol")
kaava <- function(dependencies, ...) {
temp <- MetMalMurskPol@output[MetMalMurskPol@output$Hiukkaskoko == "PM10" , ]
temp <- as.data.frame(as.table(tapply(temp["MetMalMurskPolResult"], temp["MetMalMurskPolSource"], sum)))
colnames(temp)[colnames(temp) == "Freq"] <- "PäästöPM2.5Result"
return(temp)
}
Päästö.PM2.5 <- new(
"ovariable",
name = "PäästöPM2.5",
dependencies = riippuvuudet,
formula = kaava
)
# Ovariablet
## Pitoisuudet
riippuvuudet1 <- data.frame(
Name = c("Päästö.PM2.5", "LO", "LA")
)
funktio1 <- function(dependencies, ...) {
ComputeDependencies(dependencies, ...)
temp <- GIS.Concentration.matrix(Päästö.PM2.5, LO, LA, ...)
return(temp)
}
Pitoisuus <- new(
"ovariable",
name = "Pitoisuus",
formula = funktio1,
dependencies = riippuvuudet1
)
temp <- EvalOutput(Pitoisuus, N = 1)
temp2 <- gsub("\\(", "", temp@output$LObin)
temp2 <- gsub("\\]", "", temp2)
temp2 <- strsplit(temp2, ",")
temp3 <- c()
for(i in 1:length(temp2)){
a <- as.numeric(temp2[[i]][1])
b <- as.numeric(temp2[[i]][2])
temp3[i] <- ((a + b) / 2)
}
temp@output$LO <- temp3
temp2 <- gsub("\\(", "", temp@output$LAbin)
temp2 <- gsub("\\]", "", temp2)
temp2 <- strsplit(temp2, ",")
temp3 <- c()
for(i in 1:length(temp2)){
a <- as.numeric(temp2[[i]][1])
b <- as.numeric(temp2[[i]][2])
temp3[i] <- ((a + b) / 2)
}
temp@output$LA <- temp3
#print(temp)
data <- data.frame(LO=temp@output$LO, LA=temp@output$LA, concentration=temp@output$PitoisuusResult)
truncated = c()
for(i in 1:nrow(data)){
if (data$concentration[i] > 1)
{
truncated[i] <- 1
}
else
{
truncated[i] <- data$concentration[i]
}
}
data$truncated_concentration <- truncated
# Plot the data
coordinates(data)=c("LO","LA")
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) <- 41
nrow(rast) <- 41
#cat("<span style='font-size: 1.2em;font-weight:bold;'>PM2.5</span>\n")
#Rasterize point data
rast2<-rasterize(shp, rast, shp$truncated_concentration, fun=mean)
steps <- c(0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1)
colors <- rev(rainbow(length(steps), start=0, end=0.50))
colorstrip(colors, steps)
#Plot data
google.show_raster_on_maps(rast2, col=colors, style="height:500px;")
#Plot a ggplot graph
ggplot(temp@output, aes(LO, LA, fill = PitoisuusResult)) + geom_tile() +
xlab("Longitude (degrees)") + ylab("Latitude (degrees)") +
opts(title = "Outdoor PM concentration due to the emission") +
scale_fill_gradient(limits = c(0, 1), low = "cyan", high = "red") +
scale_x_continuous(expand = c(0,0)) +
scale_y_continuous(expand = c(0,0))
##############################################33
funktio = function(dependencies, ...){
ComputeDependencies(dependencies, ...)
input <- tidy(op_baseGetData("opasnet_base", "Op_fi1234"), objname = "Arviointi") # Vaihda tekstin Op_fi1234 paikalle oman arviointisi sivutunniste.
# Tähän muu tarpeellinen laskenta.
return(out@output)
}
Arviointi <- new("ovariable",
name = "Arviointi",
dependencies = riippuvuudet,
formula = funktio
)
temp <- EvalOutput(Arviointi)
print(xtable(temp@output),type = "html")
objects.put(Arviointi) # Poista tämä ja seuraava rivi, jos et halua käyttää tuloksia jatkolaskennassa.
cat("Muuttuja alustettu. Kopioi sivun osoitteen avain talteen käyttöä varten.\n")
| |