[GIS] NA values in the raster after changing the resolution, extent, and origin

rraster

I'm working with two rasters that differ in their origin, extent, and resolution. I have a bathymetry raster, with a very high resolution (x=0.0008333333, y=0.0008333333) and a MUCH great spatial extent. I also have a sea surface temperature raster, which has a much coarser resolution (x=0.04166667, y=0.04166667). Both rasters have the same projection (longlat, datum=WGS84). Both rasters also have 'NA' values that correspond to landmasses.

Bathymetry data:

https://www.ngdc.noaa.gov/dem/squareCellGrid/download/4956

Sst data: June 25-July 2 (2016) 8 day composite, cropped to -131.8333, -127.6667, 50.125, 54.54166 (xmin, xmax, ymin, ymax)

https://oceancolor.gsfc.nasa.gov/cgi/l3?per=8D&prd=NSST_sst.nc&sen=A&res=9km&num=24&ctg=Standard&date=25Jun2016

I have manipulated the bathymetry raster (using R) to match the extent, origin, and resolution of the sea surface temperature raster. 

bathy2<-projectRaster(bathymetry, sst, method="bilinear")

I made a data frame so that I could compare the sea surface temperature values and the bathymetry values at each pixel. I then noticed that I have a problem. Some of the pixels have 'NA' values for bathymetry, but numeric values for SST. How can this be, given that the landmasses should be in the exact same place? The original bathymetry raster has a much finer resolution than the SST raster, so its outline of the landmasses should be accurate.

Is this a result of using projectRaster(), or an issue with the data itself? Anyway I can fix it?

Best Answer

Bathymetry raster that you give us as example isn't worldwide, such SST raster. So, if you try to project the raster to SST with projectRaster() you'll get most of NA values in resampled data, because function resample to new resolution/extent/CRS. First, crop SST raster and after that, project bathymetry raster:

library(raster)

bathy <- raster("~/Downloads/British_Columbia_DEM_5076/british_columbia_3sec.asc")

crs(bathy) <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs "

SST <- raster("~/Downloads/A20163372016344.L3m_8D_NSST_sst_9km.nc")

SST2 <- crop(SST, bathy)

plot(SST2, legend=F)

plot1

bathy2 <- projectRaster(bathy, SST2, method = "bilinear")

plot(bathy2, col=grey(1:100/100))
plot(SST2, add=T, alpha=0.6,legend=F)

plot2

NA data

From product description:

The digital coastline used in developing the British Columbia DEM was generated by merging vector coastlines from NOAA Electronic Navigational Charts (ENCs) and Natural Resources Canada (NRCAN) then edited based on ESRI world imagery layer. The final digital coastline was converted to xyz format with elevation set at zero and point spacing at 10 meters. The digital coastline was also converted to a polygon and ultimately a raster for masking topography and eliminating interpolated data.

df <- as.data.frame(stack(bathy2, SST2))

dim(df[complete.cases(df),])[1]
[1] 8070

length(df$british_columbia_3sec[!is.na(df$british_columbia_3sec)])
[1] 8710

length(df$Sea.Surface.Temperature[!is.na(df$Sea.Surface.Temperature)])
[1] 8264

test <- overlay(is.na(bathy2),is.na(SST2),fun=sum)

test[test!=2] <- NA ;test[test==2] <- 1

plot(test, col="yellow",legend=F)
plot(bathy2, col=scales::alpha("red",0.5),legend=F,add=T)
plot(SST2, col=scales::alpha("blue",0.5),legend=F,add=T)

plot3

So you can expect valid values in coastlines or inner terrain. Those NA values comes from mismatch between valid values of SST in coastline. If you want to complete values for SST, you need to merge bathymetry with a surface DEM.

Related Question