Calculate average time between timespec

Question

I am working with an averaging function following the formula

new average = old average * (n-1) / n + (new value / n)

When passing in doubles this works great. My example code for a proof of concept is as follows.

double avg = 0;

    uint16_t i;
    for(i=1; i<10; i++) {
        int32_t new_value = i;
        avg = avg*(i-1);
        avg /= i;
        avg += new_value/i;
        printf("I %d New value %d Avg %f\n",i, new_value, avg);
    }

In my program I am keeping track of messages received. Each time I see a message its hit count is increased by 1, it is them timestamped using a timespec. My goal is to keep a moving average (like above) of the average time between messages of a certain type being received.

My initial attempt was to average the tv_nsec and tv_sec separately as follows

static int32_t calc_avg(const int32_t current_avg, const int32_t new_value, const uint64_t n) {
    int32_t new__average = current_avg;
    new__average = new__average*(n-1);
    new__average /= n;
    new__average += new_value/n;
    return new__average;
}

void average_timespec(struct timespec* average, const struct timespec new_sample, const uint64_t n) {
    if(n > 0) {
        average->tv_nsec = calc_avg(average->tv_nsec, new_sample.tv_nsec, n);
        average->tv_sec = calc_avg(average->tv_sec, new_sample.tv_sec, n);
    }
}

My issue is I am using integers, the values are always rounded down and my averages are way off. Is there a smarter/easier way to average the time between timespec readings?

Answer 1

Below is some code that I've used a lot [in production S/W] for years.

The main idea is that just because clock_gettime uses struct timespec does not mean this has to be "carried around" everywhere:

1. It's easier to convert to a long long or double and propagate those values as soon as they're gotten from clock_gettime.

2. All further math is simple add/subtract, etc.

3. The overhead of the clock_gettime call dwarfs the multiply/divide time in the conversion.

Whether I use the fixed nanosecond value or the fractional seconds value depends upon the exact application.

In your case, I'd probably use the double since you already have calculations that work for that.

Anyway, this is what I use:

#include <time.h>

typedef long long tsc_t;                    // timestamp in nanoseconds

#define TSCSEC      1000000000LL
#define TSCSECF     1e9

tsc_t tsczero;                              // initial start time
double tsczero_f;                           // initial start time

// tscget -- get number of nanoseconds
tsc_t
tscget(void)
{
    struct timespec ts;
    tsc_t tsc;

    clock_gettime(CLOCK_MONOTONIC,&ts);
    tsc = ts.tv_sec;
    tsc *= TSCSEC;
    tsc += ts.tv_nsec;

    tsc -= tsczero;

    return tsc;
}

// tscgetf -- get fractional number of seconds
double
tscgetf(void)
{
    struct timespec ts;
    double sec;

    clock_gettime(CLOCK_MONOTONIC,&ts);

    sec = ts.tv_nsec;
    sec /= TSCSECF;
    sec += ts.tv_sec;

    sec -= tsczero_f;

    return sec;
}

// tscsec -- convert tsc value to [fractional] seconds
double
tscsec(tsc_t tsc)
{
    double sec;

    sec = tsc;
    sec /= TSCSECF;

    return sec;
}

// tscinit -- initialize base time
void
tscinit(void)
{

    tsczero = tscget();
    tsczero_f = tscsec(tsczero);
}

Answer 2

1. Use better integer math.

• Use signed math if new_value < 0 possible, else int64_t cast not needed below.

• Form the sum first and then divide.

• Round.

Sample code:

// new__average = new__average*(n-1);
// new__average /= n;
// new__average += new_value/n;

//              v-------------------------------------v  Add first
new__average = (new__average*((int64_t)n-1) + new_value + n/2)/n;
//                            Add n/2 to effect rounding  ^-^ 

---

2. On review, the whole idea of doing averages in 2 parts is flawed. Instead use a 64-bit count of nanoseconds. Good until the year 2263.

Suggested code:

void average_timespec(int64_t* average, struct timespec new_sample, int64_t n) {
  if (n > 0) {
    int64_t t = new_sample.tv_sec + new_sample.tv_nsec*(int64_t)1000000000;
    *average = (*average*(n-1) + t + n/2)/n; 
  }
}  

If you must form a struct timespec from the average, easy to do when average >= 0.

int64_t average;
average_timespec(&average, new_sample, n);
struct timespec avg_ts = (struct timespec){.tm_sec = average/1000000000, 
    .tm_nsec = average%1000000000);

Answer 3

If you wanted to keep integers, you need to do math with fixed point numbers. That's a bit tricky.

In most cases, this is done by reserving N bits after the decimal point. Unfortunately, with the timespec structure, the decimal point would be a number from 0 to 999,999,999, which is not a power of 2. I still think you can achieve it.

Using a compiler that supports the __int128 type, you could convert he timespec to an __int128 first, then do the math from that number.

__int128 t = ts.tv_nsec + ts.tv_sec * 1000000000;

Now you can compute an average using t which has a precision of 9 digits after the decimal point. If you want a little more precision, say another 2 digits, you can use:

__int128 t = ts.tv_nsec * 100 + ts.tv_sec * 100000000000;

i.e. multiple the left and right by another 100.

_note: the very large numbers probably need casting; so the above probably needs to be written as to make sure the math works as expected:

const __int128 one_hundred = 100;
const __int128 one_hundred_billion = 100000000000;
const __int128 seconds = ts.tv_sec;
const __int128 nanoseconds = ts.tv_nsec;
__int128 t = nanoseconds * one_hundred
           + seconds * one_hundred_billion;

There no way to write a full __int128 literal number (that I know of) and you may need to use some math to use such.

Now t as a precision of 11 digits. Do your math with these fixed point values -- i.e. additions and subtractions are just fine, multiplications and divisions require additional work...

__int128 t3 = t1 + t2;  // t3 is defined as expected
__int128 t3 = t1 - t2;

__int128 t6 = t4 * t5;  // t6 is not as expected, it was multiplied by another billion

__int128 t7 = t4 * t5 / one_billion; // or `one_hundred_billion`
                        // t7 is valid, same number of decimal digits

__int128 t10 = t9 / t8; // t10 is missing decimal digits
__int128 t11 = t9 * one_billion / t8;
                        // t11 is valid

__int128 t12 = t9 / n;  // this is valid if n is a regular integer (no decimal digits)

etc.

You will still have similar precision errors, only much less than using plain integers. It's just that the math is harder to follow between plain integers and fixed point numbers. If using doubles is more than enough, then Craig's solution is the simplest.