Difference between revisions of "OpasnetUtils/GIS"

Latest revision as of 09:54, 26 August 2013

library(OpasnetUtils)

concentration <- GIS.Concentration.matrix(Emission = emissions, LA = LA, LO = LO)

oprint(head(concentration@output))

exposure <- GIS.Exposure(concentration, LA = LA, LO = LO)

oprint(head(exposure@output))

Description

This section contains a couple of general GIS related functions

Code

https://www.opasnet.org/svn/opasnet_utils/trunk/R/GIS%20methods.r

Kehitysehdotus

Tarkastellaan GIS.Exposure-funktion tätä osaa:

temp <- Population * Concentration.matrix
	
	if(dbug) cat(colnames(temp@output), "\n")
	
	#temp <- oapply(temp, cols = c("LObin", "LAbin"), sum)
	# Sum over spatial data.
	out <- tapply(
		temp@output$Result, 
		temp@output[,colnames(temp@output)[temp@marginal & !colnames(temp@output) %in% c("LObin", "LAbin")]], 
		sum,
		na.rm = TRUE
	)

Väestö * pitoisuus -sarakkeen lisäksi lasketaan

• Ryhmittely quantiles kvantiiliin pitoisuuden perusteella. Oletusarvo on 1 eli nykyinen käytäntö. Jos quantiles = 0, ei tapplyta lainkaan vaan lasketaan alkuperäisillä riveillä. Joka tapauksessa kuitenkin LO ja LA poistetaan, koska paikkatietoa ei haluta antaa funktion läpi. --# : Koska kvantiilit voivat muiden indeksien takia olla hyvinkin erisuuruisia, voi kuitenkin olla järkevämpää laskea yksinkertaisesti cut:lla pilkottuja tasalevyisiä pitoisuusluokkia. --Jouni 10:46, 17 March 2013 (EET)
• tapplytaan yllä olevalla ehdolla mutta lisättynä juuri luodulla kvantiilisarakkeella:
  • väestö (funktio: sum)
  • pitoisuus (funktio: mean)

Tällä tavalla saadaan output, jossa on tarvittavat tiedot HIA-laskentaan, esim. Exposures in Finland-tyyppiseen tauluun Exposure, ja väestön koko.

Functions

GIS.Exposure

Exposure computes exposure using a given concentration matrix and protected population data from Heande.

Inputs:

• Concentration.matrix - A matrix containing spatially dependent concentratio data;
• LO & LA - coordinates of the center of the concentration matrix;
• distx & disty - maximum displacement from center of concentration matrix, assumed symmetrical, defaults to 10.5 km (Piltti source-receptor matrix);
• resolution - resolution of concentration matrix, length of side of grid element which are assumed squares, defaults to 1 km (Piltti source-receptor matrix)

Output:

• An ovariable containing result of Population * Concentration. Output and marginal slots are defined. Spatial information lost in summation (though there is an easy way around it).

GIS.Concentration.matrix

Computes a concentration matrix from given emission and coordinates, based on random sampling Piltti source-receptor matrix.

Inputs:

• Emission - emission of substance in Mga^-1;
• LO & LA - coordinates where emission occurs;
• distx & disty - maximum displacement in kilometers from center of desired matrix, assumed symmetrical, defaults to 10.5 km (Piltti source-receptor matrix);
• resolution - resolution of desired matrix (length in kilometers of side of grid element which are assumed squares), defaults to 1 km (Piltti source-receptor matrix);
• N - number of iterations to be run

Output:

• An ovariable containing spatially dependent concentration data. Output and marginal slots are defined.

See also

• OpasnetBaseUtils
• Object-oriented programming in Opasnet
• Opasnet (R library)
• Piltti source-receptor matrix

Calculations

This code creates a function that is used to plot a concentration field on a Google map.

library(OpasnetUtils)

# Function concmap draws a concentration field on a dynamic Google map.
concmap <- function(Pitoisuus, ncol_rast = 41, nrow_rast = 41, ...) {

# Required packages:
#library(OpasnetUtils)
#library(ggplot2)
#library(xtable)
#library(OpasnetUtilsExt)
#library(rgdal)
#library(maptools)
#library(RColorBrewer)
#library(classInt)
#library(raster)
#par(mfrow=c(6,1), mar=c(3,1,0,1), cex=1.5) # Graphics settings


	cat("Concentration field.\n")

	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)
	}

	# x is a factor with ranges. Function range2sum averages the ranges and gives averages as the output.
	range2num <- function(x) {
		temp <- gsub("[\\(\\)\\[ ]|\\]", "", levels(x)) # Remove characters "()[] "
		temp <- strsplit(temp, ",") # Split into two from the comma (assuming exactly one comma)

		temp2 <- c()

		for(i in 1:length(temp)) {
			a <- as.numeric(temp[[i]][1])
			b <- as.numeric(temp[[i]][2])
			temp2[i] <- ((a + b) / 2) # Average of the lower and upper part of range.
		}
		levels(x) <- temp2
		x <- as.numeric(as.character(x))
		return(x)
	}

	data <- data.frame( # Create the data.frame that is used as the input for the concentration field.
		LO = range2num(Pitoisuus@output$LObin), 
		LA = range2num(Pitoisuus@output$LAbin), 
		log_concentration = (log10(Pitoisuus@output$PitoisuusResult) - log10(0.001)) / 4 # Miksi jaetaan neljällä?
	)

	# 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) <- ncol_rast
	nrow(rast) <- nrow_rast
	
	#Rasterize point data
	rast2 <- rasterize(shp, rast, shp$log_concentration, fun=mean)
	
	steps <- c(0.001,0.003,0.01,0.03,0.1,0.3,1,3,10,30)
	steps <- c(steps[steps < max(data$log_concentration)], max(data$log_concentration)) # Limit upper part to max value.

	colors <- rev(rainbow(length(steps), start=0, end=0.50))

	cat(colorstrip(colors, steps)) # Plot the colorstrip (legend).

	#Plot data
	return(google.show_raster_on_maps(rast2, col=colors, style="height:500px;"))
}

objects.store(concmap)
cat("Object concmap stored. It plots a concentration field on a Google map.\n")