oceano latitudine longitudine punto di distanza dalla costa

Question

Ho avviato un progetto "libero" open source per creare un nuovo set di dati per il pH degli oceani terrestri.

ho iniziato dal aperto set di dati dal NOAA e ha creato un 2,45 milioni le righe di dati-set con quelle colonne:

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

documento Metodo HERE.

Data set HERE.

Il mio obiettivo ora è "qualificare" ogni riga (2,45 m) ... per farlo, ho bisogno di calcolare la distanza da ogni punto di Lat/Long alla riva più vicina.

, quindi sono alla ricerca di un metodo che avrebbe preso A: Lat/Long Out: Distanza (km dalla costa)

Con questo, posso qualificare se il punto di dati può essere influenzata dalla contaminazione riva, come l'effluenza della città vicina, per esempio.

Ho una ricerca per un metodo per farlo, ma tutto sembra aver bisogno di pacchetti/software che non ho.

Se qualcuno fosse disposto a dare una mano, apprezzerei. Oppure, se siete a conoscenza di un metodo semplice (gratuito) per raggiungere questo obiettivo, per favore fatemelo sapere ...

posso lavorare nella programmazione R, Shell script roba, ma non un esperto di quelli ....

Answer 1

Quindi ci sono diverse cose qui. Innanzitutto, il set di dati sembra avere pH vs. profondità. Quindi, mentre ci sono ~ 2.5 MM di file, ci sono solo ~ 200.000 righe con profondità = 0 - ancora molto.

In secondo luogo, per ottenere la distanza dalla costa più vicina è necessario uno shapefile delle linee costiere. Fortunatamente questo è disponibile here, all'eccellente Natural Earth website.

In terzo luogo, i dati sono in long/lat (quindi, unità = gradi), ma si desidera la distanza in km, quindi è necessario trasformare i dati (i dati di costa sopra sono anche in long/lat e devono anche essere trasformato). Un problema con le trasformazioni è che i tuoi dati sono evidentemente globali e qualsiasi trasformazione globale sarà necessariamente non planare. Quindi la precisione dipenderà dalla posizione attuale. Il modo giusto per farlo è quello di grigliare i dati e quindi utilizzare una serie di trasformazioni planari appropriate a qualsiasi griglia in cui si trovano i punti. Questo va oltre lo scopo di questa domanda, quindi, useremo una trasformazione globale (mollweide) solo per darvi un'idea di come è fatto 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") 

Quindi questa legge il set di dati, estrae le righe in cui Depth==0, e converte che a uno SpatialPoints oggetto. Quindi leggiamo il database delle linee costiere scaricato dal link sopra in un oggetto SpatialLines. Quindi trasformiamo entrambi alla proiezione di Mollweide utilizzando spTransform(...), quindi utilizziamo gDistance(...) nel pacchetto rgeos per calcolare la distanza minima tra ciascun punto e la costa più vicina.

Ancora una volta, è importante ricordare che nonostante tutte le cifre decimali, queste distanze sono solo approssimative.

Un problema molto grande è la velocità: questo processo richiede ~ 2 minuti per 1000 distanze (sul mio sistema), quindi per eseguire tutte le 200.000 distanze occorrerebbero circa 6,7 ​​ore. Un'opzione, in teoria, sarebbe quella di trovare un database di costa con una risoluzione più bassa.

Il codice seguente calcolerà tutte le 201.000 distanze.

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

EDIT: il commento di OP sui nuclei mi ha fatto pensare che questo potrebbe essere un caso in cui il miglioramento dalla parallelizzazione potrebbe essere valsa la pena. Ecco come eseguiresti questo (su Windows) usando l'elaborazione parallela.

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 

Utilizzando 4 core migliora la velocità di elaborazione di circa un fattore di 3. Quindi, dal momento che 1000 distanze richiede circa un minuto, 100.000 dovrebbe prendere un po 'meno di 2 ore.

Si noti che l'utilizzo di times=1 è un abuso di microbenchmark(...) in realtà, poiché l'intero punto è eseguire il processo più volte e mediare i risultati, ma non ho avuto la pazienza.