Hi, Robert is right, obeyRes is not required as the function res() does exactly what obeyRes was tought for... I wrote that function because I had to create rasters with a given resolution to fit as good as possible into a given extent and I was trying to do that with raster() that has no "res" argument + obeyRes. So obeyRes was calculation the available ncol/nrow argument in raster()! But rasterExtent + res() does the job perfectly! Thanks Matteo
>>> "Robert J. Hijmans" 01.03.13 19.05 Uhr >>> There is some confusion here: To get a google map in a sinusoidial projection you can do library(dismo) # extent e <- extent(c(-1.87710, -1.772096, 55.61399, 55.682171)) # google map g <- gmap(e, type='hybrid') #projection prj <- "+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs" # project s1 <- projectRaster(g, crs=prj) # done # for more control # create empty raster r <- raster(ext=e,crs="+proj=longlat +ellps=WGS84 +datum=WGS84") # project extent r <- projectExtent(r,crs=prj) # set resolution to desired value # your value of 1000 seemed too large here. # I do not think you should obeyRes here # (I am not sure what it is useful for) res(r) <- 100 # project s2 <- projectRaster(g, r) Hope this clarifies. Robert On Fri, Mar 1, 2013 at 5:56 AM, Ross Ahmed wrote: > Hi all > > I have this raster: > > myExtent <- extent(c(-1.87710, -1.772096, 55.61399, 55.682171)) > r <- raster(ext=myExtent,crs="+proj=longlat +ellps=WGS84 +datum=WGS84") > > Ive projected the raster to an equal area projection (Sinusoidal): > > myRaster <- projectExtent(r,crs="+proj=sinu +lon_0=0 +x_0=0 +y_0=0 > +a=6371007.181 +b=6371007.181 +units=m +no_defs") > > Using the obeyRes function kindly given to me by Matteo Mattiuzzi > ( > https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20130224/3247ffd4/atta > chment.pl), I made each cell in the raster 1km^2: > > res1km <- obeyRes(extent(myRaster),resolution=sqrt(1)*1000) > myRaster1km <- raster(ext=res1km$extent,ncol=res1km$cols,nrows=res1km$rows, > crs="+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 > +b=6371007.181 +units=m +no_defs") > > I've downloaded a google map with long/lat projection: > > g <- gmap(myExtent, type='hybrid', interpolate=T, lonlat=T) > > I'm now trying to project the google map using Sinusoidal (ie, same as > myRaster1km). I've two functions: > > gg <- projectRaster(g, CRS("+proj=sinu +lon_0=0 +x_0=0 +y_0=0 > +a=6371007.181 > +b=6371007.181 +units=m +no_defs")) > gg <- projectExtent(g,crs="+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 > +b=6371007.181 +units=m +no_defs") > > How can I plot the google map using Sinusoidal? Ultimately, my aim is 1. > plot gg 2. plot the myRaster1km over gg 3. draw polygons in each cell of > the > myRaster1km and have R recognise which cell I'm drawing in. > > Once again, many thanks for help in advance and apologies for any > ignorance. > > Ross > > > > > > [[alternative HTML version deleted]] > > _______________________________________________ > R-sig-Geo mailing list > R-sig-Geo@r-project.org > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > [[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo