pythonpython-3.xnumpygrib

Revert time averaged values to instant values with numpy


I am having a 3D field of an averaged quantity. Now I am looking for a pythonic way to revert the averaging to get instantaneous values for each time stamp. Note: The averaging is made from the beginning of the whole period. It is like rolling mean with the window size adapted to the index of the value to be averaged.

For better clarification I give an 1D example:

import numpy as np

input_array = np.array([
       [0.      ],
       [0.5     ],
       [1.      ],
       [1.5     ],
       [2.      ],
       [2.5     ],
       [3.      ],
       [3.25    ],
       [3.333333],
       [3.3     ],
       [3.181818],
       [3.      ],
       [2.769231]
])

exp_result = de_average(input_array)

The expected result exp_result should look like this:

exp_result= np.array([
       [0],
       [1],
       [2],
       [3],
       [4],
       [5],
       [6],
       [5],
       [4],
       [3],
       [2],
       [1],
       [0]])

Solution

  • Something like this should work. You're problem is very related to the reverse cumsum problem, so I'm using part of the solution given here.

    from itertools import tee 
    
    def pairwise(iterable):
        "s -> (s0,s1), (s1,s2), (s2, s3), ..."
        a, b = tee(iterable)
        next(b, None)
        return zip(a, b)
    
    def inverse_cumsum(cumulative):
        """ Inverse an array obtained by cumulative sum
           [1, 3, 6, 10, 15, 21, 28, 36, 45] -> [1, 2, 3, 4, 5, 6, 7, 8, 9] 
        """
        yield cumulative[0]
        for a, b in pairwise(cumulative):
            yield b-a
    
    def inverse_average(a, decimals=1):
       """ Inverse an array averaged (where each entry i is the average up to i) """ 
       deav = inverse_cumsum([a * (i + 1) for i, a in enumerate(a)])
       return np.array(list(deav)).round(decimals)
    
    inverse_average(input_array)