2015-12-03 2 views
2

Мне было интересно, использовал ли кто-либо ggmap радиус круга вокруг точки lat/long? Например, я хотел бы построить конкретную точку, а затем нарисовать и затенять в радиусе 2500 футов вокруг этой точки. У меня есть идея в моей голове о том, как это сделать, используя формулу окружности большего круга, но полагал, что я сначала проверил бы это.Радиус круга вокруг точки lat/long на созданном ggmap графике

ответ

1

Я работал над созданием простого хака, а не круга, который я себе представлял, но он будет работать пока.

library(ggmap) 
##__________________________________________________________________ 
### earth.dist I found on r-bloggers. I believe it now belongs to the fossil package. 
earth.dist <- function (lat1,long1,lat2, long2) 
{ 
    rad <- pi/180 
    a1 <- lat1 * rad 
    a2 <- long1 * rad 
    b1 <- lat2 * rad 
    b2 <- long2 * rad 
    dlon <- b2 - a2 
    dlat <- b1 - a1 
    a <- (sin(dlat/2))^2 + cos(a1) * cos(b1) * (sin(dlon/2))^2 
    c <- 2 * atan2(sqrt(a), sqrt(1 - a)) 
    R <- 6378.145 
    d <- R * c 
    # this I changed to return feet instead of KM 
    return(d* 3280.8) 
} 

##__________________________________________________________________ 
## Function to output polygon to map 
BoxGon <- function(Lat,Long,feet){ 

    for(i in 1:1000){ 
    point = Long - i/1000 
    Dist <- earth.dist(Lat,Long, Lat,point) 
    if(Dist > feet){ 
        West <- cbind(Lat,point) 
        break} 
    } 

    for(i in 1:1000){ 
    point = Long + i/1000 
    Dist <- earth.dist(Lat,Long, Lat,point) 
    if(Dist > feet){ 
        East <- cbind(Lat,point) 
        break} 
    } 

    for(i in 1:1000){ 
    point = Lat + i/1000 
    Dist <- earth.dist(Lat,Long, point,Long) 
    if(Dist > feet){ 
        North <- cbind(point,Long) 
        break} 
    } 

    for(i in 1:1000){ 
    point = Lat - i/1000 
    Dist <- earth.dist(Lat,Long, point,Long) 
    if(Dist > feet){ 
        South <- cbind(point,Long) 
        break} 
    } 

    return(rbind(West,North,East,South,West)) 
} 

##__________________________________________________________________ 
df = BoxGon(37.295844, -121.898057,5000) 
df = as.data.frame(df) 
colnames(df) <- c('Latitude','Longitude') 

map <- get_map(location = 'san,jose', zoom = 12) 
map <- ggmap(map) 

# we select - 1 because once we map in pairs. IE once we have the last record there is nothing for that record to map to 
for(i in 1:nrow(df)-1){ 
latlon <- head(df,2) 
map <- map + geom_polygon(data=latlon,aes(x=Longitude,y=Latitude),alpha=0.1,size = 1,colour="green",fill="green") 
df = df[-1,] 
print(i) 
} 
map <- map + geom_polygon(aes(x=-121.898057,y=37.295844),alpha=0.1,size = 6,colour='Purple',fill='Purple') 
Смежные вопросы