use_reduceat.py

from __future__ import division

import numpy as np
import Ska.Numpy

times = np.arange(0, 86400 * 365, 10)  # 10 second sampling for a year
vals = np.random.uniform(size=len(times))
on_mask = vals < 0.5  # 50% duty cycle
on_vals = vals[on_mask]


def on_freq1():
    """
    For N samples and M days this runs in time O(N * M)
    """
    days = times // 86400
    on_days = days[on_mask]
    out = np.array([np.sum(on_vals[on_days == day]) for day in np.arange(365)])
    return out


def on_freq2():
    """
    For N samples and M days this runs in time O(N)
    """
    on_times = times[on_mask]
    day_times = np.arange(365) * 86400
    indices = Ska.Numpy.search_both_sorted(on_times, day_times)
    out = np.add.reduceat(on_vals, indices)
    # To get mean, do:
    # nvals = np.add.reduceat(np.ones_like(on_vals), indices)
    # mean = out / n_vals
    return out

# In [10]: timeit out1 = on_freq1()
# 1 loops, best of 3: 5.01 s per loop
#
# In [11]: timeit out2 = on_freq2()
# 10 loops, best of 3: 40.5 ms per loop