I am trying to extract pixel values from a time series of raster images to xy points.
The object tile.evib is a raster stack of Modis EVI data for 2 decades, while tairb is a raster brick time series of air temperature values for 2 decades as well. I would like to extract pixel values to a set of xy points.
library(raster)
> tile.evib
class : RasterStack
dimensions : 2307, 1825, 4210275, 457 (nrow, ncol, ncell, nlayers)
resolution : 0.01168239, 0.008333333 (x, y)
extent : 11.73743, 33.05778, -34.83383, -15.60883 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
names : MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, MOD13A2.A//EVI_masked, ...
min values : -1999, -1999, -2000, -1999, -2000, -2000, -1999, -1999, -1998, -2000, -1998, -1998, -1985, -1997, -1997, ...
> tairb
class : RasterBrick
dimensions : 1801, 3600, 6483600, 240 (nrow, ncol, ncell, nlayers)
resolution : 0.1, 0.1 (x, y)
extent : -0.05, 359.95, -90.05, 90.05 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : /home/fbiome/Attribution/TTRforceTime.nc
names : X2000.01.01.00.00.00, X2000.02.01.00.00.00, X2000.03.01.00.00.00, X2000.04.01.01.00.00, X2000.05.01.01.00.00, X2000.06.01.01.00.00, X2000.07.01.01.00.00, X2000.08.01.01.00.00, X2000.09.01.01.00.00, X2000.10.01.01.00.00, X2000.11.01.00.00.00, X2000.12.01.00.00.00, X2001.01.01.00.00.00, X2001.02.01.00.00.00, X2001.03.01.00.00.00, ...
Date/time : 2000-01-01 00:00:00, 2019-12-01 00:00:00 (min, max)
varname : t2m
I then did the following:
tile.xyb <-coordinates(tile.evib) #get the lon lat coords of the raster stack
tile.xyb<-SpatialPoints(tile.xyb,proj4string=CRS(projection(tile.evib))) #make them a SP data frame
tile.xyb<-spTransform(tile.xyb, CRS(projection(tile.evib))) #project them into same projection as stack
> tile.xyb
class : SpatialPoints
features : 4210275
extent : 11.74327, 33.05194, -34.82966, -15.613 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
#Extract
ts.evib<-extract(tile.evib,tile.xyb)
ts.tairb<-extract(tairb,tile.xyb)
> head(ts.evib)
MOD13A2.A2000049.1_km_16_days_EVI_masked MOD13A2.A2000065.1_km_16_days_EVI_masked
[1,] NA NA
[2,] NA NA
MOD13A2.A2000081.1_km_16_days_EVI_masked MOD13A2.A2000097.1_km_16_days_EVI_masked
[1,] NA NA
[2,] NA NA
MOD13A2.A2000113.1_km_16_days_EVI_masked MOD13A2.A2000129.1_km_16_days_EVI_masked
[1,] NA NA
[2,] NA NA
MOD13A2.A2000145.1_km_16_days_EVI_masked MOD13A2.A2000161.1_km_16_days_EVI_masked
[1,] NA NA
> head(ts.tairb)
X2000.01.01.00.00.00 X2000.02.01.00.00.00 X2000.03.01.00.00.00
[1,] NA NA NA
[2,] NA NA NA
[3,] NA NA NA
[4,] NA NA NA
X2000.04.01.01.00.00 X2000.05.01.01.00.00 X2000.06.01.01.00.00
[1,] NA NA NA
[2,] NA NA NA
[3,] NA NA NA
[4,] NA NA NA
X2000.07.01.01.00.00 X2000.08.01.01.00.00 X2000.09.01.01.00.00
[1,] NA NA NA
I do not understand why I am getting NAs. Any assistance?
Regards
Edward