2014-02-11 166 views
0

作为a question I posted yesterday,的跟进,是否有R中的包/功能可以给我国(和洲),其中由经度和纬度定义的点位于?沿着什么样的线东西做here in MatLab.从R经度和纬度点获取国家(和大陆)

数据框看起来像这样...

Point_Name    Longitude Latitude 
University of Arkansas 36.067832 -94.173655 
Lehigh University  40.601458 -75.360063 
Harvard University  42.379393 -71.115897 

而且我想与国家和大洲列添加到它的输出上面的数据帧。 作为的额外奖励,美国专栏针对美国(而“其他”针对美国以外的专栏)?

+1

看一看[**这个问题**](http://stackoverflow.com/questions/14334970/convert-latitude-and-longitude-coordinates-to-country-name-in-r ?rq = 1),从本页右边空白处的“相关”列表中挑选出来。 – Henrik

+0

[这个答案](http://stackoverflow.com/a/8751965/980833)会让你的状态。 –

+0

谢谢,我不知道我是如何错过的(我曾看过相当多的数字,但没有找到那些) - 只要我根据旧的答案得到脚本,我将删除我的Q. – Ell

回答

1

因此,这是在Google上使用reverse geocoding API的替代方案。该代码部分基于this reference

呼唤你的数据帧以上df

reverseGeoCode <- function(latlng) { 
    require("XML") 
    require("httr") 
    latlng <- as.numeric(latlng) 
    latlngStr <- gsub(' ','%20', paste(round(latlng,2), collapse=",")) 
    url <- "http://maps.google.com" 
    path <- "/maps/api/geocode/xml" 
    query <- list(sensor="false",latlng=latlngStr) 
    response <- GET(url, path=path, query=query) 
    if (response$status !=200) { 
    print(paste("HTTP Error:",response$status),quote=F) 
    return(c(NA,NA)) 
    } 
    xml <- xmlInternalTreeParse(content(response,type="text")) 
    status <- xmlValue(getNodeSet(xml,"//status")[[1]]) 
    if (status != "OK"){ 
    print(paste("Query Failed:",status),quote=F) 
    return(c(NA,NA)) 
    } 
    xPath <- '//result[1]/address_component[type="country"]/long_name[1]' 
    country <- xmlValue(getNodeSet(xml,xPath)[[1]]) 
    xPath <- '//result[1]/address_component[type="administrative_area_level_1"]/long_name[1]' 
    state <- xmlValue(getNodeSet(xml,xPath)[[1]]) 
    return(c(state=state,country=country)) 
} 
st.cntry <- t(apply(df,1,function(x)reverseGeoCode(x[2:3]))) 
result <- cbind(df,st.cntry) 
result 
#    Point_Name Longitude Latitude   state  country 
# 1 University of Arkansas 36.06783 -94.17365  Arkansas United States 
# 2  Lehigh University 40.60146 -75.36006 Pennsylvania United States 
# 3  Harvard University 42.37939 -71.11590 Massachusetts United States 

在API定义中, “administrative_area_level_1” 低于国家的最高行政区域。在美国,这些是国家。在其他国家,定义各不相同(例如,可能是省份)。

顺便说一句,我相当肯定你有经度和纬度逆转。

4

要得到大陆,您可以修改使用rworldmap从这个answer创建coords2continent功能,如下图所示的coords2country函数的最后一行。选择是否需要6或7洲模式。我会考虑将此代码放入rworldmap

library(sp) 
library(rworldmap) 

# The single argument to this function, points, is a data.frame in which: 
# - column 1 contains the longitude in degrees 
# - column 2 contains the latitude in degrees 
coords2continent = function(points) 
{ 
    countriesSP <- getMap(resolution='low') 
    #countriesSP <- getMap(resolution='high') #you could use high res map from rworldxtra if you were concerned about detail 

    # converting points to a SpatialPoints object 
    # setting CRS directly to that from rworldmap 
    pointsSP = SpatialPoints(points, proj4string=CRS(proj4string(countriesSP))) 


    # use 'over' to get indices of the Polygons object containing each point 
    indices = over(pointsSP, countriesSP) 

    #indices$continent # returns the continent (6 continent model) 
    indices$REGION # returns the continent (7 continent model) 
    #indices$ADMIN #returns country name 
    #indices$ISO3 # returns the ISO3 code 
} 

这是一个测试。

points = data.frame(lon=c(0, 90, -45, -100, 130), lat=c(52, 40, -10, 45, -30)) 

coords2continent(points) 
#[1] Europe  Asia   South America North America Australia 
coords2country(points) 
#[1] United Kingdom China Brazil United States of America Australia