从 netcdf 层获取两个变量的月平均值和总体平均值，多年来在 R 中

问题描述

提前道歉，我对使用 netcdf 文件很陌生，想象这很简单，但我似乎无法破解它。我对变量、维度和实际数据感到困惑。我想要尽可能简单的代码。

我有一个来自哥白尼的 netcdf 文件，当通过以下方式将其引入 R 时，其结构如下： chl <- nc_open(file) 使用 ncdf4 包。

所以我有 2 个变量 - 叶绿素和净初级生产力。 4个维度：时间、深度、纬度、经度。

 2 variables (excluding dimension variables):
    float chl[longitude,latitude,depth,time]   
        long_name: Total Chlorophyll
        standard_name: mass_concentration_of_chlorophyll_a_in_sea_water
        units: mg m-3
        unit_long: milligram of Chlorophyll per cubic meter
        _FillValue: 9.96920996838687e+36
        _ChunkSizes: 1
         _ChunkSizes: 15
         _ChunkSizes: 137
         _ChunkSizes: 288
    float nppv[longitude,time]   
        long_name: Total Primary Production of Phyto
        standard_name: net_primary_production_of_biomass_expressed_as_carbon_per_unit_volume_in_sea_water
        units: mg m-3 day-1
        unit_long: milligrams of Carbon per cubic meter per day
        _FillValue: 9.96920996838687e+36
        _ChunkSizes: 1
         _ChunkSizes: 15
         _ChunkSizes: 137
         _ChunkSizes: 288

 4 dimensions:
    time  Size:323
        long_name: Time (hours since 1950-01-01)
        standard_name: time
        calendar: gregorian
        units: hours since 1950-01-01 00:00:00
        axis: T
        _ChunkSizes: 1024
        _CoordinateAxisType: Time
        valid_min: 377316
        valid_max: 612503.9375
    depth  Size:1
        valid_min: 0.505760014057159
        valid_max: 0.505760014057159
        units: m
        positive: down
        unit_long: Meters
        long_name: Depth
        standard_name: depth
        axis: Z
        _ChunkSizes: 75
        _CoordinateAxisType: Height
        _CoordinateZisPositive: down
    latitude  Size:153
        valid_min: 34
        valid_max: 72
        step: 0.25
        units: degrees_north
        unit_long: degrees north
        long_name: Latitude
        standard_name: latitude
        axis: Y
        _ChunkSizes: 681
        _CoordinateAxisType: Lat
    longitude  Size:185
        valid_min: -27
        valid_max: 19
        step: 0.25
        units: degrees_east
        unit_long: degrees East
        long_name: Longitude
        standard_name: longitude
        axis: X
        _ChunkSizes: 1440
        _CoordinateAxisType: Lon

17 global attributes:
    product: GLOBAL_REANALYSIS_BIO_001_029
    producer: cmeMS - Global Monitoring and Forecasting Centre
    title: Monthly mean fields for product GLOBAL_REANALYSIS_BIO_001_029
    area: GLOBAL
    quality_information_document: http://marine.copernicus.eu/documents/QUID/cmeMS-GLO-QUID-001-029.pdf
    Conventions: CF-1.6
    credit: E.U. copernicus Marine Service information (cmeMS)
    contact: servicedesk.cmems@mercator-ocean.eu
    references: http://marine.copernicus.eu
    source: MERCATOR FREEBIORYS2V4
    licence: http://marine.copernicus.eu/services-portfolio/service-commitments-and-licence/
    dataset: global-reanalysis-bio-001-029-monthly
    institution: Mercator Ocean
    product_user_manual: http://marine.copernicus.eu/documents/PUM/cmeMS-GLO-PUM-001-029.pdf
    _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention
    comment: 
    history: Data extracted from dataset http://localhost:8080/thredds/dodsC/global-reanalysis-bio-001-029-monthly

我想获得一个简单的代码：

将此 netcdf 文件放入光栅堆栈，每个月一个用于 CHL，一个用于 PP。
将每个“月”（时间属性）拆分为年平均值。
给出整个数据集中所有月份的 Chl 的总体平均值和 PP 的平均值。
有一行代码要包含，它采用深度平均，用于上述 CHL 和 PP。目前有 1 个深度级别，深度大小：1，但我最多可以有 72 个，并且如果可能的话，我希望在我的代码中包括取深度平均值。**

到目前为止......我已经在代码中做到了这一点：

    # Opened the netcdf file
    nc <- nc_open(file) # Open file
    # Use attributes to get variable attributes
    attributes(nc)$names ```
    # Got the nc variable names
    attributes(nc$var)$names
    [1] "chl"  "nppv"
    
    #Get attributes of variables
    ncatt_get(nc,attributes(nc$var)$names[1])
    ncatt_get(nc,attributes(nc$var)$names[2])
    
    # Get  latitude and longitude values - dimensions
    nc_lat <- ncvar_get( nc,attributes(nc$dim)$names[3])
    nc_lon <- ncvar_get( nc,attributes(nc$dim)$names[4])
    nc_time <- ncvar_get( nc,attributes(nc$dim)$names[1])

但是从这里开始，我很迷茫如何将它变成可用的格式，例如一个光栅。任何帮助将不胜感激。

解决方法

应该没有必要像您一样直接使用 ncdf4 包（并且在本网站上看到几个类似的问题）。通常，您可以执行 b <- brick(file);但我在这里使用较新的 terra 包

 library(terra)
 chl <- rast("filename","chl")
 pp <- rast("filename","nppv")

从那里您可以使用，例如，tapp。但是，如果这可行，您也许可以编辑您的问题和 show(chl)，以便更容易看到您正在处理的内容，并且可能更容易描述您想要完成的任务。分享一个示例文件也会很有帮助。

netcdf r r raster spatial