fre-cli/fre/app/generate_time_averages/cdoTimeAverager.py at 3b994053ce870f0ee1b62212f6df6902f1f29405 · NOAA-GFDL/fre-cli · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
''' class using (mostly) cdo functions for time-averages '''

import logging

from netCDF4 import Dataset
import numpy as np

import cdo
from cdo import Cdo

from fre import log_and_raise
from .timeAverager import timeAverager

fre_logger = logging.getLogger(__name__)

class cdoTimeAverager(timeAverager):
    '''
    class inheriting from abstract base class timeAverager
    generates time-averages using cdo (mostly, see weighted approach)
    '''

    def generate_timavg(self, infile = None, outfile = None):
        """
        use cdo package routines via python bindings

        :param self: This is an instance of the class cdoTimeAverager
        :param infile: path to history file, or list of paths, default is None
        :type infile: str, list
        :param outfile: path to where output file should be stored, default is None
        :type outfile: str
        :return: 1 if the instance variable self.avg_typ is unsupported, 0 if function has a clean exit
        :rtype: int
        """

        if self.avg_type not in ['all', 'seas', 'month']:
            log_and_raise(f'requested unknown avg_type {self.avg_type}')

        if self.var is not None:
            fre_logger.warning('WARNING: variable specification not twr supported for cdo time averaging. ignoring!')

        fre_logger.info('python-cdo version is %s', cdo.__version__)

        _cdo = Cdo()

        wgts_sum = 0
        if not self.unwgt: #weighted case, cdo ops alone don't support a weighted time-average.

            nc_fin = Dataset(infile, 'r')

            time_bnds = nc_fin['time_bnds'][:].copy()
            # Ensure float64 precision for consistent results across numpy versions
            # NumPy 2.0 changed type promotion rules (NEP 50), so explicit casting
            # is needed to avoid precision differences
            time_bnds = np.asarray(time_bnds, dtype=np.float64)
            # Transpose once to avoid redundant operations
            time_bnds_transposed = np.moveaxis(time_bnds, 0, -1)
            wgts = time_bnds_transposed[1] - time_bnds_transposed[0]
            # Use numpy.sum for consistent dtype handling across numpy versions
            wgts_sum = np.sum(wgts, dtype=np.float64)

            fre_logger.debug('wgts_sum = %s', wgts_sum)

        if self.avg_type == 'all':
            fre_logger.info('time average over all time requested.')
            if self.unwgt:
                _cdo.timmean(input = infile, output = outfile, returnCdf = True)
            else:
                _cdo.divc( str(wgts_sum), input = "-timsum -muldpm "+infile, output = outfile)
            fre_logger.info('done averaging over all time.')

        elif self.avg_type == 'seas':
            fre_logger.info('seasonal time-averages requested.')
            _cdo.yseasmean(input = infile, output = outfile, returnCdf = True)
            fre_logger.info('done averaging over seasons.')

        elif self.avg_type == 'month':
            fre_logger.info('monthly time-averages requested.')
            outfile_str = str(outfile)
            _cdo.ymonmean(input = infile, output = outfile_str, returnCdf = True)
            fre_logger.info('done averaging over months.')

            fre_logger.warning(" splitting by month")
            outfile_root = outfile_str.removesuffix(".nc") + '.'
            _cdo.splitmon(input = outfile_str, output = outfile_root)
            fre_logger.debug('Done with splitting by month, outfile_root = %s', outfile_root)

        fre_logger.info('done averaging')
        fre_logger.info('output file created: %s', outfile)
        return 0