Ozean Breite Länge Punkt Entfernung vom Ufer

Question

Ich begann ein "freies" Open-Source-Projekt, um einen neuen Datensatz für den pH-Wert der Weltmeere zu erstellen.

ich aus dem offenen Datensatzes von NOAA gestartet und erstellt eine 2,45 Millionen Zeilen-Daten-Set mit diesen Spalten:

colnames(NOAA_NODC_OSD_SUR_pH_7to9) 
[1] "Year" "Month" "Day" "Hour" "Lat" "Long" "Depth" "pH" 

Methode Dokument HERE.

Datensatz HERE.

Mein Ziel ist es nun, jede Reihe (2,45m) zu "qualifizieren" ... um dies zu tun, muss ich den Abstand von jedem Punkt von Lat/Long zum nächsten Ufer berechnen.

Also für eine Methode, die ich bin auf der Suche, die In dauern würde: Breite/Länge Out: Entfernung (km von der Küste)

Damit kann ich, wenn der Datenpunkt qualifizieren kann vom Ufer Kontamination betroffen sein, wie zum Beispiel der nahegelegene Stadtausfluss.

Ich habe nach einer Methode gesucht, um dies zu tun, aber alles scheint Pakete/Software zu brauchen, die ich nicht habe.

Wenn jemand bereit wäre zu helfen, würde ich mich freuen. Oder wenn Sie von einem einfachen (kostenlos) Methode wissen dies zu erreichen, lassen Sie es mich wissen ...

I in R Programmierung arbeiten können, Skripte Shell Material, aber kein Experte von denen ....

Answer 1

So gibt es hier einige Dinge. Erstens scheint Ihr Datensatz einen pH-Wert in Abhängigkeit von der Tiefe zu haben. Während also ~ 2.5MM Reihen sind, gibt es nur ~ 200.000 Reihen mit Tiefe = 0 - immer noch eine Menge.

Zweitens, um die Entfernung zur nächsten Küste zu erhalten, benötigen Sie ein Shapefile von Küstenlinien. Zum Glück ist dies here, in der ausgezeichneten Natural Earth website.

Drittens sind Ihre Daten in long/lat (also Einheiten = Grad), aber Sie wollen die Entfernung in km, also müssen Sie Ihre Daten transformieren (die Küsten Daten oben sind auch in long/lat und müssen auch umgewandelt werden). Ein Problem bei Transformationen besteht darin, dass Ihre Daten offensichtlich global sind und jede globale Transformation notwendigerweise nicht-planar ist. Die Genauigkeit hängt also vom tatsächlichen Standort ab. Der richtige Weg dazu besteht darin, Ihre Daten zu rasterieren und dann eine Reihe planarer Transformationen zu verwenden, die zu dem Raster passen, in dem sich Ihre Punkte befinden. Dies ist jedoch außerhalb des Rahmens dieser Frage. Daher verwenden wir eine globale Transformation (Mollweide). nur, um Ihnen eine Vorstellung davon, wie es in R.

library(rgdal) # for readOGR(...); loads package sp as well 
library(rgeos) # for gDistance(...) 

setwd(" < directory with all your files > ") 
# WGS84 long/lat 
wgs.84 <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0" 
# ESRI:54009 world mollweide projection, units = meters 
# see http://www.spatialreference.org/ref/esri/54009/ 
mollweide <- "+proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs" 
df  <- read.csv("OSD_All.csv") 
sp.points <- SpatialPoints(df[df$Depth==0,c("Long","Lat")], proj4string=CRS(wgs.84)) 

coast <- readOGR(dsn=".",layer="ne_10m_coastline",p4s=wgs.84) 
coast.moll <- spTransform(coast,CRS(mollweide)) 
point.moll <- spTransform(sp.points,CRS(mollweide)) 

set.seed(1) # for reproducible example 
test <- sample(1:length(sp.points),10) # random sample of ten points 
result <- sapply(test,function(i)gDistance(point.moll[i],coast.moll)) 
result/1000 # distance in km 
# [1] 0.2185196 5.7132447 0.5302977 28.3381043 243.5410571 169.8712255 0.4182755 57.1516195 266.0498881 360.6789699 

plot(coast) 
points(sp.points[test],pch=20,col="red") 

So dies liest Dataset getan, extrahiert Zeilen, in denen Depth==0 und wandelt die zu einem SpatialPoints Objekt. Dann lesen wir die Küstenlinien-Datenbank, die von dem obigen Link heruntergeladen wurde, in ein SpatialLines-Objekt. Dann transformieren wir beide zu der Mollweide-Projektion unter Verwendung von spTransform(...), dann verwenden wir gDistance(...) in dem rgeos-Paket, um den minimalen Abstand zwischen jedem Punkt und der nächsten Küste zu berechnen.

Auch hier ist es wichtig zu bedenken, dass trotz aller Dezimalstellen diese Abstände nur ungefähr sind.

Ein sehr großes Problem ist die Geschwindigkeit: Dieser Vorgang dauert ~ 2 min für 1000 Entfernungen (auf meinem System), so dass alle 200.000 Entfernungen etwa 6,7 ​​Stunden dauern würde. Eine Option wäre theoretisch, eine Küstendatenbank mit einer niedrigeren Auflösung zu finden.

Der folgende Code berechnet alle 201.000 Entfernungen.

## not run 
## estimated run time ~ 7 hours 
result <- sapply(1:length(sp.points), function(i)gDistance(sp.points[i],coast)) 

EDIT: OP Kommentar über den Kern hat mich zu denken, dass dies ein Fall sein könnte, wo die Verbesserung von Parallelisierung der Mühe wert sein könnte. So, hier ist, wie Sie dies (unter Windows) mit der parallelen Verarbeitung ausführen würden.

library(foreach) # for foreach(...) 
library(snow)  # for makeCluster(...) 
library(doSNOW) # for resisterDoSNOW(...) 

cl <- makeCluster(4,type="SOCK") # create a 4-processor cluster 
registerDoSNOW(cl)    # register the cluster 

get.dist.parallel <- function(n) { 
    foreach(i=1:n, .combine=c, .packages="rgeos", .inorder=TRUE, 
      .export=c("point.moll","coast.moll")) %dopar% gDistance(point.moll[i],coast.moll) 
} 
get.dist.seq <- function(n) sapply(1:n,function(i)gDistance(point.moll[i],coast.moll)) 

identical(get.dist.seq(10),get.dist.parallel(10)) # same result? 
# [1] TRUE 
library(microbenchmark) # run "benchmark" 
microbenchmark(get.dist.seq(1000),get.dist.parallel(1000),times=1) 
# Unit: seconds 
#      expr  min  lq  mean median  uq  max neval 
#  get.dist.seq(1000) 140.19895 140.19895 140.19895 140.19895 140.19895 140.19895  1 
# get.dist.parallel(1000) 50.71218 50.71218 50.71218 50.71218 50.71218 50.71218  1 

Unter Verwendung von 4 Kerne Verarbeitungsgeschwindigkeit um etwa einen Faktor 3 verbessert also seit 1000 Entfernungen etwa eine Minute dauert, sollte 100.000 nehmen etwas weniger als 2 Stunden.

Beachten Sie, dass die Verwendung von ist ein Missbrauch von microbenchmark(...) wirklich, wie der ganze Punkt ist, den Prozess mehrmals auszuführen und die Ergebnisse durchschnittlich, aber ich hatte einfach nicht die Geduld.