From 9029147433d0de2a26ce270051f27cfc3022bda3 Mon Sep 17 00:00:00 2001
From: Georges Kunstler <Georges.Kunstler@gmail.com>
Date: Wed, 4 Sep 2013 01:16:46 +1000
Subject: [PATCH] close #21 fixed geographic coordinates system
---
merge.data.NZ.R | 20 +++++++++++++-------
1 file changed, 13 insertions(+), 7 deletions(-)
diff --git a/merge.data.NZ.R b/merge.data.NZ.R
index 2bf3e99..8c5ab34 100644
--- a/merge.data.NZ.R
+++ b/merge.data.NZ.R
@@ -44,25 +44,31 @@ data.nz$htot <- rep(NA, nrow(data.nz)) ## Max height is already available so ha
data.nz$tree.id <- rep(NA, nrow(data.nz))
#data.nz$tree.id <- gsub("__", "_", data.nz$tree.id)
#data.nz$tree.id <- gsub("_", ".", data.nz$tree.id) ## tree unique id
-data.nz$weights <- 1/(pi*(0.5*data.nz$D/100)^2)
+data.nz$weights <- 1/(20*20) ## need to check with david
########################################## CHANGE COORDINATE SYSTEM DON'T KNOW THE EPSG CODE HERE change coordinates
########################################## system of Easting Northing to be in lat long WGS84
-data.nz <- merge(data.plot[, c(1, 3, 4)], data.nz, by = "plid")
library(sp)
-# library(dismo); library(rgdal);
+library(dismo)
+library(rgdal)
data.sp <- data.nz[, c("tree.id", "Easting", "Northing")]
coordinates(data.sp) <- c("Easting", "Northing") # define x y
-proj4string(data.sp) <- CRS("+init=epsg:23030") # define projection system of our data ## EPSG CODE 23030 ED50 / UTM zone 30N
+proj4string(data.sp) <- CRS("+init=epsg:27200") # define projection system of our data ## EPSG CODE 27200 NZGD49 / New Zealand Map Grid
summary(data.sp)
detach(package:rgdal)
data.sp2 <- spTransform(data.sp, CRS("+init=epsg:4326")) ## change projection in WGS84 lat lon
data.nz$Lon <- coordinates(data.sp2)[, "Easting"]
data.nz$Lat <- coordinates(data.sp2)[, "Northing"]
-## ## plot on world map library(rworldmap) newmap <- getMap(resolution = 'coarse')
-## # different resolutions available plot(newmap)
-## points(data.sp2,cex=0.2,col='red')
+## ## plot on world map
+pdf("./figs/map.NZ.pdf")
+library(rworldmap)
+newmap <- getMap(resolution = 'coarse')
+## # different resolutions available
+plot(newmap,xlim=c(166,177),ylim=c(-34,-47))
+points(data.sp2,cex=0.2,col='red')
+dev.off()
+
rm(data.sp, data.sp2)
###################### ECOREGION merge greco to have no ecoregion with low number of observation
--
GitLab