I wanted to implement running median (windowed). But google came up with all kinds of fancy algorithms.. I thought they were just too fancy. This is my first solution to update median, mean, min and max of the last N datapoints:
The main idea is simple.. keep a double buffer of the data where the entries are ordered. Then it is easy to remove and add samples while keeping the samples ordered.
This is simple test plot using the algorithm:
**broken link removed**
blue: raw data
red: 64 sample running average
green: 64 sample running median
Any ideas or interest how to improve this.. is this useful at all
?
I can post more code if somebody shows interest.
Test code:
The main idea is simple.. keep a double buffer of the data where the entries are ordered. Then it is easy to remove and add samples while keeping the samples ordered.
C:
struct rolling_window {
volatile uint16_t *data;
volatile uint8_t index;
volatile uint8_t size;
};
struct rolling_median {
struct rolling_window window;
uint16_t *ordered_data;
uint16_t median;
uint32_t accum;
};
/**
* Updates Rolling Median with one data entry.
*
* @param[in] object Pointer to rolling_median structure.
* @param[in] data New data.
*
*/
void
rolling_median_update(struct rolling_median *object, uint16_t data)
/* Updates Rolling Median with one data entry. */
{
uint16_t temp;
temp = rolling_window_put_and_get(&(object->window), data);
find_and_replace_ordered(object->ordered_data, temp, data, object->window.size);
object->accum += data;
object->accum -= temp;
object->median = (object->ordered_data)[(object->window.size)/2];
}
/******************************************************************************/
void find_and_replace_ordered(uint16_t *data, uint16_t old, uint16_t new, uint8_t size)
{
uint8_t i;
/* Best case */
if (old == new) { return; }
/* Find the (first) index place of the old entry */
/* This can be optimized using binary search etc. */
i=0;
while(old != data[i]) {
i++;
}
if (new > old) {
/* Move everything one place */
while (i<(size-1)) {
if (data[i+1]>=new) {
data[i] = new;
return;
}
data[i] = data[i+1];
i++;
}
}
else
{
/* Move everything one place */
while (i>(0)) {
if (data[i-1]<=new) {
data[i] = new;
return;
}
data[i] = data[i-1];
i--;
}
}
data[i]=new;
}
/******************************************************************************/
/**
* Writes one data entry in Rolling Window and returns the oldest entry.
*
* @param[in] object Pointer to rolling_window structure.
* @param[in] data Data to be written in buffer.
*
* Return uint16_t Oldest data entry that is removed from the window.
*
*/
uint16_t rolling_window_put_and_get(struct rolling_window *object, uint16_t data);
#define NEXT_INDEX(obj) \
(obj->index) ? (obj->index--) : ((obj->index) = ((obj->size)-1))
uint16_t
rolling_window_put_and_get(struct rolling_window *object, uint16_t data)
/* Writes one data entry in Rolling Window and returns the oldest entry. */
{
uint16_t temp;
temp = (object->data)[object->index];
(object->data)[object->index] = data;
NEXT_INDEX(object);
return temp;
}
/******************************************************************************/
This is simple test plot using the algorithm:
**broken link removed**
blue: raw data
red: 64 sample running average
green: 64 sample running median
Any ideas or interest how to improve this.. is this useful at all
I can post more code if somebody shows interest.
Test code:
C:
/* Initialize the sine table */
for (int i=0; i<256; i++)
{
sineTable[i] = (int8_t)(127.0 * sin(6.283*((float)i)/256.0));
}
initialize_uart(B38400);
sei();
/* Infinite loop */
for(;;)
{
fbi(PORTB, 0);
_delay_ms(50);
accumulator += increment;
acc2 += inc2;
temp = sineTable[accumulator>>8] + 128;
temp += sineTable[acc2>>8] + 128;
if(rand()>30000) {
temp+= (rand()>>4);
}
temp+= (rand()>>10);
rolling_median_update(&samples, (uint16_t)temp);
temp2 = (uint16_t)((samples.accum)/WINDOW_SIZE);
printf("%u\t", temp);
printf("%u\t", temp2);
printf("%u\n", (samples.median));
}
Last edited: