I've uploaded my data onto github here: https://github.com/dschmidt42/Benguela_Chla_Danielle_Schmidt.
It's concerning chl a data off the eastern coast of Africa and spans 26 years. I am trying to split the data into each month and then compare the variance within months. This is my first time working with a netCDF file and I'm having a lot of issues. The main one at the moment is that my longitude is on the x axis and latitude is on the y axis. How do I swap these?
What I've done
nc_chla_file <- "cmems_mod_glo_bgc_my_0.25_P1M-m_1673881119174.nc"
nc <- nc_open(nc_chla_file)
Here are the dimensions for longitude, latitude, depth, and months respectively. I only need depth 1 as that is the closest to the surface.
dim(ncvar_get(nc,varid="chl"))
Then I've separated chl data for January and turned it into a raster.
chl_jan <- ncvar_get(nc,varid="chl")[,,1,seq(1,312,12)]
chl_jan_rast <- rast(chl_jan)
Lastly I've taken the variance and made a plot.
variance_chl_jan <- terra::app(chl_jan_rast, fun = "var")
plot(variance_chl_jan, main = "Variance in chl a in January")
As you can see the latitude and longitude are the wrong way round. How do I fix this?
You can open the file like this:
library(terra)
#terra 1.6.47
f <- "cmems_mod_glo_bgc_my_0.25_P1M-m_1673881119174.nc"
r <- rast(f)
Inspect
r
#class : SpatRaster
#dimensions : 82, 48, 936 (nrow, ncol, nlyr)
#resolution : 0.25, 0.25 (x, y)
#extent : 8.875, 20.875, -35.875, -15.375 (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84
#source : cmems_mod_glo_bgc_my_0.25_P1M-m_1673881119174.nc
#varname : chl (Total Chlorophyll)
#names : chl_d~001_1, chl_d~855_1, chl_d~819_1, chl_d~001_2, chl_d~855_2, chl_d~819_2, ...
#unit : mg m-3, mg m-3, mg m-3, mg m-3, mg m-3, mg m-3, ...
#time : 1994-01-16 12:00:00 to 2019-12-16 12:00:00 UTC
names(r)[1:3]
#[1] "chl_depth=0.50576001_1" "chl_depth=1.555855_1" "chl_depth=2.6676819_1"
You can get the layer numbers for the depth you want like this
i <- grep("chl_depth=0.5", names(r))
#or like this
i <- seq(1, nlyr(r), 3)
And subset SpatRaster r
r1 <- r[[i]]
To get the variance within each month, you can now do
vm <- tapp(r1, "month", var)
vm
#class : SpatRaster
#dimensions : 82, 48, 12 (nrow, ncol, nlyr)
#resolution : 0.25, 0.25 (x, y)
#extent : 8.875, 20.875, -35.875, -15.375 (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84
#source(s) : memory
#names : m_1, m_2, m_3, m_4, m_5, m_6, ...
#min values : 0.0001622432, 2.360285e-05, 1.136451e-05, 4.821852e-05, 0.0002333084, 0.0009578676, ...
#max values : 3.0582257258, 3.222419e+00, 3.880621e+00, 4.295268e+00, 4.2899640566, 3.3797369950, ...
#time (mnts) : Jan to Dec
plot(vm)
Doing things "by hand" like you do is more error-prone, but you could fix your result by doing
vcj <- flip(t(variance_chl_jan), "v")