Nitin
Nitin

Reputation: 2874

Plotting points converted to grids on map in R

I am new to R. I have a set of lat,lon coordinates and I have generated a map with those points on it. enter image description here

EDIT: Here is the data for the points,

lat <- c(44.3672, 39.3421, 33.978, 36.3901, 32.8388, 38.7519, +
         37.1863, 28.2408, 41.7098, 30.4127)
lon <- c(-83.5215, -86.0034, -84.579, -83.8163, -83.6077, -85.4381,+
         -85.2338, -82.7283, -85.0308, -88.4076)
df <- data.frame(lon, lat)

This is the code I used to plot the points, (obtained from this post)

library(maps)
library(ggplot2)
us_states<-map_data('state')
p <- ggplot(legend=FALSE) +
  geom_polygon( data=us_states, aes(x=long, y=lat,group=group)) +
  theme(panel.background = theme_blank()) +
  theme(panel.grid.major = theme_blank()) +
  theme(panel.grid.minor = theme_blank()) +
  theme(axis.text.x = theme_blank(),axis.text.y = theme_blank()) +
  theme(axis.ticks = theme_blank()) +
  xlab("") + ylab("")
# add a single point
p <- p + geom_point(data=df,aes(lon,lat),colour="blue",size=3)
p

Now I want to take each point, compute a 20km by 20km cell centered around the point and then plot the generated cells on the map. How do I then build a grid over the map using such cell sizes? Thanks.

Upvotes: 0

Views: 997

Answers (1)

IRTFM
IRTFM

Reputation: 263301

poly.df <- as.data.frame( t( apply(df, 1, 
              function(x) c(                        # for longitude, 
               xmin=x[1]-cos(x[2]*2*pi/360)*20/111.2,  # km/degr varies
               xmax=x[1]+cos(x[2]*2*pi/360)*20/111.2,  # with latitude
               ymin=x[2]-20/111.2,            # but not for lattitude
               ymax=x[2]+20/111.2) ) ) )

To be frank, I wasn't sure it's is a simple cos(latitude) adjustment. I checked with the table in Wikipedia and it's pretty close:

 round( cos(seq(90, 0, by=-15)*2*pi/360)*111.2 , 1)
[1]   0.0  28.8  55.6  78.6  96.3 107.4 111.2

But this does get rectangles that look roughly accurate. I'm thinking any off-centered-ness comes from the positioning of the dots.

p+geom_rect(data=poly.df, color="red", fill="transparent", 
       aes(xmin=xmin.lon,xmax=xmax.lon, ymin=ymin.lat,ymax=ymax.lat) )

To construct a grid:

    p+geom_vline(xintercept= seq(-50, -130, by=-1), color="red")+
      geom_hline(yintercept= seq(25,50, by=1), color="red")

Upvotes: 1

Related Questions