海洋緯度經度距離海岸的距離

Question

我開始了一個「免費」開源項目，爲海洋pH值創建一個新的數據集。

我從NOAA的設置數據開放開始，創造了2.45億行數據集與列：

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

方法文檔HERE。

數據集HERE。

我現在的目標是「限定」每一行（2.45米）......爲此，我需要計算從經緯度的每個點到最近的岸的距離。

所以我要尋找，將採取 在一個方法：緯度/龍 出：距離（離岸邊公里）

有了這個，我可以出線，如果數據點可以從岸上污染的影響，比如附近的城市流量。

我有搜索一個方法來做到這一點，但似乎都需要我沒有的軟件包/軟件。

如果有人願意幫忙，我將不勝感激。 或者，如果你知道一個簡單的（免費）的方法來做到這一點，請讓我知道...

我可以一個R編程，Shell腳本的東西的工作，而不是那些專家....

Answer 1

所以有幾件事情在這裏進行。首先，你的數據集似乎有pH值與深度。因此，雖然有〜2.5MM行，但只有~20萬行，深度= 0 - 仍然很多。

其次，要獲得距最近的海岸的距離，您需要海岸線的shapefile。幸運的是，這是可用here，在優秀的Natural Earth website。第三，您的數據是long/lat（so，units = degrees），但是您希望以km爲單位的距離，所以您需要轉換您的數據（上面的海岸線數據也是long/lat，並且還需要被改變）。轉換的一個問題是您的數據顯然是全球性的，任何全球轉型必然是非平面的。所以準確度將取決於實際的位置。做到這一點的正確方法是對數據進行網格劃分，然後使用一組適合於你的點的網格的平面變換。然而，這超出了這個問題的範圍，因此我們將使用全局變換（mollweide）只給你它是如何在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") 

因此，這將讀取的數據集，做了一個想法提取行，其中Depth==0，並將其轉換爲SpatialPoints對象。然後我們將從上面鏈接下載的海岸線數據庫讀入SpatialLines對象。然後我們使用​​將兩者轉換爲Mollweide投影，然後我們在rgeos包中使用gDistance(...)來計算每個點與最近的海岸之間的最小距離。

同樣重要的是要記住，儘管所有的小數位，這些距離只是大約。

一個非常大的問題是速度：這個過程需要大約1000分鐘的距離（在我的系統上），所以要運行所有200,000個距離大約需要6.7個小時。理論上，一種選擇是找到分辨率較低的海岸線數據庫。

下面的代碼將計算所有201,000個距離。

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

編輯：OP的有關核評論讓我想着，這可能是這樣的情況：從並行的改善可能是值得的。所以這裏是你如何運行這個（在Windows上）使用並行處理。

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 

使用4個芯由3.所以約爲係數提高了處理速度，因爲1000米的距離需要大約一分鐘，100000應該大於2小時少一點。

請注意，使用times=1實際上是濫用microbenchmark(...)，因爲整個過程是多次運行該過程並對結果取平均值，但我只是沒有耐心。