Difference between revisions of "Opasnet map"

library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(OpasnetUtils)
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')

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")
colcode	<- findColours(class,plotclr)

epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(shp)<-("+init=epsg:3067")
shp2	<- spTransform(shp, epsg4326String)

kmlname	<- "Kuopio house data"
kmldescription <- "Random stuff about here"
icon	<- "http://maps.google.com/mapfiles/kml/pal2/icon18.png"
name	<- paste("ika value: ", shp2$ika)
description <- paste("<b>Value:</b>",shp2$ika,"<br><b>Description:</b>",shp2$kayttotark)

dat <- MyPointKML(shp2, kmlname, kmldescription, name, description, icon, colcode)

google.show_kml_data_on_maps(dat)

library(OpasnetUtils)
library(OpasnetUtilsExt)
library(sorvi)
library(rgdal)

data(MML)

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

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

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

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)

library(RgoogleMaps)
library(rgdal)
library(maptools)
library(RColorBrewer)
library(classInt)
library(OpasnetUtils)
library(OpasnetUtilsExt)

shp	<- readOGR('PG:host=localhost user=postgres dbname=spatial_db','kuopio_house')
#oprint(shp@data[1:100 , ])

plotvar	<- shp@data$ika
nclr	<- 8

epsg4326String <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(shp) <- ("+init=epsg:3067")
shp2	<- spTransform(shp,epsg4326String)

# MyMap function without the "file destination" parameter
# Requires RgoogleMaps package

MyRmap <- function (
	shp2, # a spatial data object
#	lonR, # longitude of data points
#	latR, # latitude of data points
	center, # center of the map
	size = c(640, 480), # size of the map. This produces the right dimensions in Opasnet.
	MINIMUMSIZE = FALSE, 
	RETURNIMAGE = TRUE, 
	GRAYSCALE = FALSE, 
	NEWMAP = TRUE, 
	zoom, 
	verbose = 1, 
	pch = 19, 
	cex = 0.3,
	legend_title = "",
	...
) {
	plotclr	<- brewer.pal(nclr,"Spectral")
	class	<- classIntervals(plotvar,nclr,style="quantile")
	colcode	<- findColours(class,plotclr)

latR <- shp2@coords[ , 2]
lonR <- shp2@coords[ , 1]

	mymarkers <- data.frame(
		lat = c(shp2@coords[,2]),
		lon = c(shp2@coords[,1]),
		color = colcode
	)

	#get the bounding box:

	bb <- qbbox(lat = mymarkers[,"lat"], lon = mymarkers[,"lon"])

    if (missing(zoom)) 
        zoom <- min(MaxZoom(latR, lonR, size))
    if (missing(center)) {
        lat.center <- mean(latR)
        lon.center <- mean(lonR)
    }
    else {
        lat.center <- center[1]
        lon.center <- center[2]
    }
    if (MINIMUMSIZE) {
        ll <- LatLon2XY(latR[1], lonR[1], zoom)
        ur <- LatLon2XY(latR[2], lonR[2], zoom)
        cr <- LatLon2XY(lat.center, lon.center, zoom)
        ll.Rcoords <- Tile2R(ll, cr)
        ur.Rcoords <- Tile2R(ur, cr)
        if (verbose > 1) {
            cat("ll:")
            print(ll)
            print(ll.Rcoords)
            cat("ur:")
            print(ur)
            print(ur.Rcoords)
            cat("cr:")
            print(cr)
        }
        size[1] <- 2 * max(c(ceiling(abs(ll.Rcoords$X)), ceiling(abs(ur.Rcoords$X)))) + 
            1
        size[2] <- 2 * max(c(ceiling(abs(ll.Rcoords$Y)), ceiling(abs(ur.Rcoords$Y)))) + 
            1
        if (verbose) 
            cat("new size: ", size, "\n")
    }

	MyMap <- GetMap(
		center = c(lat.center, lon.center),
		zoom = zoom, 
		size = size, 
		RETURNIMAGE = RETURNIMAGE, 
		GRAYSCALE = GRAYSCALE, 
		verbose = verbose, 
		...
	)

	PlotOnStaticMap(MyMap)

	PlotOnStaticMap(
		MyMap,
		lat = mymarkers[,"lat"],
		lon = mymarkers[,"lon"],
		pch = pch,
		cex = cex,
		col = colcode,
		add = T
	)

	legend(
		"topleft", 
		legend	= names(attr(colcode, "table")),
		title	= legend_title, 
		fill	= attr(colcode, "palette"), 
		cex	= 1.0, 
		bty	= "y",
		bg	= "white"
	)
}

temp <- MyRmap(shp2, maptype = "mobile", legend_title = "Ikä") # bb$lonR, bb$latR