0

Since a couple of weeks ive been trying to get a NetCDF file into a dataframe. Even though Im succesfull in extracting the variables/dimensions and plotting a slice from the ncdf file, all the data is husseled when cbinding it into a dataframe and then plotting it. The data is weather data from Copernicus, with data for each longitude and latitude point in the world. The ultimate goal here is to raster the dataframe to be able to categorize weather per raster over time.

The data can be retrieved from https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-pressure-levels?tab=overview.

My code looks like this:

library(ncdf4)
library(raster)
library(ggplot2)

##dataset, 5 augustus
###bron: https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-pressure-levels?tab=overview

weatherpress <- nc_open("C:/Users/heinj/OneDrive/Documenten/Universiteit/Master/Master Thesis/dataset/allecoors5augweer.nc")
{
  sink('WeatherPressure5aug.txt')
  print(weatherpress)
  sink()
}
print(weatherpress)


long<- ncvar_get(weatherpress, "longitude")
lat <- ncvar_get(weatherpress, "latitude")
tijd <- ncvar_get(weatherpress, "time")

rain <- ncvar_get(weatherpress,"crwc")
temperature <- ncvar_get(weatherpress, "t") #temperature
Uwind <- ncvar_get(weatherpress, "u") #u wind
Vwind <- ncvar_get(weatherpress, "v") # v wind

#removal NA's
fillvaluecrwc <- ncatt_get(weatherpress, "crwc", "_FillValue")
rain[rain == fillvaluecrwc$value] <- NA
rain <- na.omit(rain)


fillvaluet <- ncatt_get(weatherpress, "t", "_FillValue")
temperature[temperature == fillvaluet$value] <- NA
temperature <- na.omit(temperature)

fillvalueu <- ncatt_get(weatherpress, "u", "_FillValue")
Uwind[Uwind == fillvalueu$value] <- NA
Uwind <- na.omit(Uwind)


fillvaluev <- ncatt_get(weatherpress, "v", "_FillValue")
Vwind[Vwind == fillvaluev$value] <- NA
Vwind <- na.omit(Vwind)
min(temperature)

#correcting longitude

nc_close(weatherpress)
#plotje
temperature_slice <- temperature[, ,1]
r_temperature <- raster(t(temperature_slice), xmn=min(long),
                        xmx=max(long), ymn=min(lat), ymx=max(lat),
                        crs=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs+ towgs84=0,0,0"))
plot(r_temperature)

#binding into a dataframe
weather <- cbind(long, lat, rain, tijd, temperature, Uwind, Vwind)
weather <- as.data.frame(weather)
View(weather)

#subset temperature
only_temperature <- weather
only_temperature <- subset(only_temperature[,c(1:2, 5)])
head(only_temperature)
summary(only_temperature)

##ggplot###
ggplot(as.data.frame(poging1.df), aes(x = long, y = lat)) +
   geom_raster(aes(fill = temperature))


### rastering ###
only_temp <- only_temperature
r = raster(xmn = min(long), xmx=max(long),
           ymn=min(lat), ymx=max(lat), res = 10)
p = as(r@extent, 'SpatialPolygons')

cordi <- only_temp[c("long", "lat")]
coordinates(cordi) <- ~long + lat
only_temp <- SpatialPointsDataFrame(cordi, only_temp)
meanr <- rasterize(only_temp, r, "temperature", fun = mean)
plot(meanr)

where the ggplot and meanr plot look different from the r_plot (the correct one) - with the temperature belonging to Antartica going dimensional through the plot.

Does anybody know where my problem lies?

Thanks in advance!

HeinB
  • 1
  • 2

1 Answers1

1

If your goal is to make a raster, you can use the terra package like this:

library(terra)
f <- "allecoors5augweer.nc"
r <- rast(f)

Or, to get a single variable, t in this case:

r <- rast(f, "t")

Also see sds(f). Alternatively, you can use raster::brick(f). To inspect the data, you can do

plot(r)

Just to answer your question about making a data.frame from a ncdf file (which does not seem to be necessary for your goals). It may depend a bit on the file at hand (for example, are there sub-datasets in the file, and how do you want to treat those --- it would have been very useful if you had shared a file). You can use tidync for any nc file, but if it is spatial (gridded) data, it may be easiest to use raster or terra. With terra you can do

library(terra)
f <- "allecoors5augweer.nc"
r <- rast(f)
d <- as.data.frame(r)
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63