I want to stack 6 rasters in a list called allrasters but first must fix crs and extent inconsistencies. Here is my code attempt to set the second raster in list to the crs of the third raster in list:

projectRaster(allrasters[[2]], crs=crs(allrasters[[3]]))

However when I run this code and check, allrasters[[2]] is still proj.merc and nothing has changed...

Raster information:

crs(allrasters[[2]])
CRS arguments:
  proj=merc  a=6378137  b=6378137  lat_ts=0  lon_0=0
 x_0=0  y_0=0  k=1  units=m  nadgrids=@null  wktext
 no_defs

crs(allrasters[[3]])
CRS arguments:
  proj=aea  lat_0=0  lon_0=-120  lat_1=34  lat_2=40.5
 x_0=0  y_0=-4000000  datum=NAD83  units=m  no_defs 

CodePudding user response：

I assume that what you are after is:

allrasters[[2]] <- projectRaster(allrasters[[2]], crs=crs(allrasters[[3]]))

That is, you forgot to assign the output of projectRaster

CodePudding user response：

I think you need a few steps:

1. you need to get all of your rasters in the same projection
2. you need to find the full extent of all rasters as if they were mosaicked together
3. you need to resample your rasters so that they have the same extent, resolution, and projection
4. you will stack your rasters.

Here is an example I created with some fake data that should help you accomplish this:

##Loading necessary packages##
library(raster)
library(rgeos)
library(tmaptools)

#For reproducibility#
set.seed(52)

##Creating fake rasters with different extents and projections##
R1<-raster(nrow=100, ncol=100, xmn=44.52, xmx=45.1, ymn=-122.1, ymx=-121.2, crs=crs(" init=epsg:4267"))
R2<-raster(nrow=100, ncol=100, xmn=44.49, xmx=45.8, ymn=-122.0, ymx=-121.3, crs=crs(" init=epsg:4326"))
R3<-raster(nrow=100, ncol=100, xmn=44.48, xmx=45.1, ymn=-122.5, ymx=-121.5, crs=crs(" init=epsg:4979"))
R4<-raster(nrow=100, ncol=100, xmn=44.55, xmx=45.6, ymn=-122.2, ymx=-121.0, crs=crs(" init=epsg:4269"))

values(R1)<-rnorm(10000, 500, 10)
values(R2)<-rnorm(10000, 1000, 60)
values(R3)<-rnorm(10000, 300, 10)
values(R4)<-rnorm(10000, 2500, 70)

##Creating a list of the rasters##
tmp<-list(R1,R2,R3,R4)

##Looping to reproject the rasters all into the same projection##
allras<-list()
for (i in 1:length(tmp)){
  if(i==1){
    allras[[i]]<-tmp[[i]]
  }else{
    allras[[i]]<-projectRaster(tmp[[i]], crs=crs(tmp[[1]]))
  }
  
}

##Creating a function to make a polygon of each raster's extent##
fxn<-function(ras){
  bb<-bbox(ras)
  bbpoly<-bb_poly(bb)
  st_crs(bbpoly)<-crs(ras)
  return(as_Spatial(bbpoly))
}
ext<-lapply(allras, fxn)

##Aggregating and dissolving all extents to get the full extent of all rasters##
full.ext<-aggregate(do.call(bind, ext), dissolve=TRUE)

##Creating a blank raster with the full extent, the desired final projection, and the desired resolution##
blank<-raster(ext=extent(full.ext), nrow=allras[[1]]@nrows, ncol=allras[[1]]@ncols, crs=allras[[1]]@crs)

##Resampling all rasters in the list to the desired extent and resolution##
rastostack<-lapply(allras, resample, y=blank)

##Stacking the rasters##
Ras<-stack(rastostack)