I have been trying to load a NetCDF with the rast
function from terra
but I get the following warning message:
> SM_rast <- rast("Soil_moisture_v3.nc", "soil_moisture")
Warning message:
[rast] unknown extent
The corresponding SpatRaster has no crs, wrong resolution and wrong coordinates:
> SM_rast
class : SpatRaster
dimensions : 2160, 4320, 172 (nrow, ncol, nlyr)
resolution : 1, 1 (x, y)
extent : 0, 4320, 0, 2160 (xmin, xmax, ymin, ymax)
coord. ref. :
source : Soil_moisture_v3.nc
names : Soil_~_v3_1, Soil_~_v3_2, Soil_~_v3_3, Soil_~_v3_4, Soil_~_v3_5, Soil_~_v3_6, ...
You can download the NetCDF file here.
Following the answer to this question I tried adding drivers="NETCDF"
to the rast()
function, but it couldn't even read the file:
> SM_rast <- rast("Soil_moisture_v3.nc", "soil_moisture", drivers="NETCDF")
Error: [rast] cannot open this file as a SpatRaster: Soil_moisture_v3.nc
In addition: Warning message:
`Soil_moisture_v3.nc' not recognized as a supported file format. (GDAL error 4)
When converting the SpatRaster to a dataframe, all cells containing NAs are missing and the coordinate columns are wrong:
> df_rast <- terra::as.data.frame(SM_rast[[1]], xy = TRUE)
> summary(df_rast)
x y Soil_moisture_v3_1
Min. : 2.5 Min. : 414.5 Min. :0.02000
1st Qu.:1525.5 1st Qu.: 937.5 1st Qu.:0.07621
Median :2403.5 Median :1240.5 Median :0.17518
Mean :2340.6 Mean :1197.1 Mean :0.20558
3rd Qu.:2841.5 3rd Qu.:1440.5 3rd Qu.:0.30225
Max. :4317.5 Max. :1960.5 Max. :0.91196
If I do the same using brick()
from the raster
package, everything works fine:
> SM_brick <- brick("data-raw/TROPOMI/Soil_moisture_v3.nc", varname = "soil_moisture")
> SM_brick
class : RasterBrick
dimensions : 2160, 4320, 9331200, 172 (nrow, ncol, ncell, nlayers)
resolution : 0.0833333, 0.0833333 (x, y)
extent : -180, 179.9998, -90, 89.99993 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs
source : Soil_moisture_v3.nc
names : X2018.02.02, X2018.02.10, X2018.02.18, X2018.02.26, X2018.03.06, X2018.03.14, X2018.03.22, X2018.03.30, X2018.04.07, X2018.04.15, X2018.04.23, X2018.05.01, X2018.05.09, X2018.05.17, X2018.05.25, ...
Date : 2018-02-02, 2021-10-24 (min, max)
varname : soil_moisture
> df_brick <- terra::as.data.frame(SM_brick[[1]], xy = TRUE)
> summary(df_brick)
x y X2018.02.02
Min. :-179.95833 Min. :-89.95834 Min. :0
1st Qu.: -89.97920 1st Qu.:-44.97919 1st Qu.:0
Median : -0.00008 Median : -0.00004 Median :0
Mean : -0.00008 Mean : -0.00004 Mean :0
3rd Qu.: 89.97905 3rd Qu.: 44.97911 3rd Qu.:0
Max. : 179.95818 Max. : 89.95826 Max. :1
NA's :8181031
I still need to use terra
because it's faster and has more useful functions not present in raster
.
It all looks good to me when I do:
library(terra)
# terra 1.7.20
rast("Soil_moisture_v3.nc")
#class : SpatRaster
#dimensions : 2160, 4320, 172 (nrow, ncol, nlyr)
#resolution : 0.0833333, 0.0833333 (x, y)
#extent : -180, 179.9998, -90, 89.99993 (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 (EPSG:4326)
#source : Soil_moisture_v3.nc
#varname : soil_moisture
#names : soil_~ure_1, soil_~ure_2, soil_~ure_3, soil_~ure_4, soil_~ure_5, soil_~ure_6, ...
#unit : m3/m3, m3/m3, m3/m3, m3/m3, m3/m3, m3/m3, ...
#time (days) : 2018-02-02 to 2021-10-24
gdal()
#[1] "3.6.2"
g <- gdal(drivers=T)
g[g$name=="netCDF", ]
# name type can vsi long.name
#88 netCDF raster read/write FALSE Network Common Data Format
If you do not see this, that may be because you have another version of terra, and/or GDAL, and/or because of your operating system and how you installed "terra".
From the contents of your zipfile I can see that you are using a Mac. The problem may be related to the GDAL that comes with the binary version on CRAN. You may try to install from R-Universe instead with
install.packages('terra', repos='https://rspatial.r-universe.dev')
Or, if that does not work, install from source code. See the instructions.
There is also this work-around:
library(raster)
# add zero to force the values out of the file.
# that takes a while
b <- brick("Soil_moisture_v3.nc") + 0
r <- rast(b)