从Python运行R脚本

Question

我希望能够从Python脚本中运行一个R脚本。这个R脚本用于在不同的坐标系统中投影lat lon坐标。我考虑了两种方式来实现这个目标。第一种方式是将lat和lon坐标解析到下面的R脚本中，最后让R脚本将x和y返回给Python脚本，但我不知道该如何实现。

project<-function(lat,lon){

library(sp)
library(rgdal)

xy <- cbind(x = lon, y = lat)
S <- SpatialPoints(xy)
proj4string(S) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
Snew <- spTransform(S, CRS("+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs"))
x <- coordinates(Snew)[1]
y <- coordinates(Snew)[2]

return(x,y)
}

对于我的第二个选择，我考虑使用底部已经有纬度和经度的R脚本。我尝试使用subprocess.Popen('Rscript project.r', shell=True)从python运行它。但是好像并没有起作用。它没有写入xy.txt文件。然而，如果我从命令行运行这个脚本，它可以正常工作。谁能帮我解决其中的一个选择？

library(sp)
library(rgdal)

lat <- 52.29999924
lon <- 4.76999998

xy <- cbind(x = lon, y = lat)
S <- SpatialPoints(xy)
proj4string(S) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
Snew <- spTransform(S, CRS("+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs"))
x <- coordinates(Snew)[1]
y <- coordinates(Snew)[2]

cat(x, file="xy.txt",sep="")
cat(y,file="xy.txt",append=TRUE)

- user2357183

---

只是出于好奇，您不使用Python GDAL绑定的原因是什么？（https://pypi.python.org/pypi/GDAL/）看起来会更容易。 - rici

如果这样做更容易的话，我肯定会尝试。然而，我不熟悉在Python中使用GDAL。是否有一个示例代码，可以将一个坐标系统中的项目坐标投影到另一个坐标系统中？ - user2357183

我认为它与R绑定没有太大区别。如果您只需要转换坐标，则只需要proj4绑定；就像这样简单：import pyproj; p1 = pyproj.Proj("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"); p2 = pyproj.Proj("+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs"); pyproj.transform(p1, p2, 52.3, 4.77) lat/lon参数可以是列表、元组、numpy数组或标量。请参见http://pyproj.googlecode.com/svn/trunk/docs/index.html - rici

Answer 1

看起来你已经快完成了，这里有一些注意事项：

• x <- coordinates(Snew)[1] 选择第一个项目，请使用 x <- coordinates(Snew)[,1] 来获取整个向量。这假设你传递的是 lat 和 lon 的向量，而不是单个值。在这里，我们没有理由不允许传递 lat lon 的向量。

• return(x,y) 是无效的，R 不支持返回多个对象。相反，只需使用 return(cbind(x,y))。

• 在函数内部，我会使用 require。

关于从 Python 运行代码，我会看看 Rpy。这使你可以从 Python 中调用 R 代码。Rpy 为你传递 R 对象来回做了很多艰难的工作，因此这是一个不错的选择。

我会在Rpy中采取的方法是将project放在一个.R文件中，source该文件，并调用project。

- Paul Hiemstra

---

library 不会重新加载已经加载的包。 - G. Grothendieck

Answer 2

针对您的第一个解决方案--让R将坐标输出到控制台(使用print、cat或其他在R中使用的命令)，然后在Python中捕获该输出(请参见此处: Running shell command from Python and capturing the output)

这将为您提供一个包含经纬度坐标的字符串，您只需要使用适当的Python函数解析它们即可。

- Colin Nichols

---

不需要手动解析R函数的输出，Rpy可以更轻松地为您完成。 - Paul Hiemstra

大家好，我想感谢你们提供的非常好的答案。我认为Rpy可以为我完成工作，但是，Rpy似乎无法与Python版本2.5.4和R版本2.14.2一起工作？ - user2357183