radio/statistics.c

// --------------------------------------------------------------
#include <malloc.h>
#include <string.h>
#include "statistics.h"

// --------------------------------------------------------------
// Sliding Mean
// --------------------------------------------------------------
void SlidingMeanInit(sl_mean_t *pObj, int N)
{
	RBufInit(&pObj->buf, N);
	pObj->Nr = (radio_float_t)1.0/N;
	pObj->mean = 0.0;
	pObj->sum = 0.0;
	pObj->normalize_count = N;
}

void SlidingMeanFree(sl_mean_t *pObj)
{
	RBufFree(&pObj->buf);
}

radio_float_t SlidingMeanProcess(sl_mean_t *pObj, radio_float_t x)
{
	uint32_t i;
	radio_float_t last;

	last = RBufGetAfter(&pObj->buf, pObj->buf.size-1);
	RBufPut(&pObj->buf, x);
	pObj->sum += (x - last);
	pObj->mean = pObj->sum*pObj->Nr;

	if (--pObj->normalize_count == 0)
	{
		pObj->normalize_count = pObj->buf.size;
		pObj->sum = 0;
		for (i=0; i < pObj->buf.size; i++)
		{
			pObj->sum += pObj->buf.pData[i];
		}
	}
	return pObj->mean;
}

void SlidingMeanProcessV(sl_mean_t *pObj, radio_float_t *pX, uint32_t len)
{
	uint32_t i;

	for (i=0; i < len; i++)
		SlidingMeanProcess(pObj, *(pX++));
}

radio_float_t SlidingMeanGet(sl_mean_t *pObj)
{
	return pObj->mean;
}

// --------------------------------------------------------------
// Sliding Variance
// --------------------------------------------------------------
void SlidingVarInit(sl_var_t *pObj, int Npwr, int Nmean)
{
	RBufInit(&pObj->buf, Npwr);
	SlidingMeanInit(&pObj->sl_mean, Nmean);
	pObj->Nr = (radio_float_t)1.0/Npwr;
	pObj->var = 0.0;
	pObj->sum = 0.0;
	pObj->normalize_count = Npwr;
}

void SlidingVarFree(sl_var_t *pObj)
{
	RBufFree(&pObj->buf);
	SlidingMeanFree(&pObj->sl_mean);
}

radio_float_t SlidingVarProcess(sl_var_t *pObj, radio_float_t x)
{
	uint32_t i;
	radio_float_t last, xm, mean, var;

	last = RBufGetAfter(&pObj->buf, pObj->buf.size-1);
	mean = SlidingMeanProcess(&pObj->sl_mean, x);
	xm = x - mean;
	xm *= xm;
	RBufPut(&pObj->buf, xm);
	pObj->sum += (xm - last);

	var = pObj->sum*pObj->Nr;
	if (var >= 0)
		pObj->var = var;

	if (--pObj->normalize_count == 0)
	{
		pObj->normalize_count = pObj->buf.size;
		pObj->sum = 0;
		for (i=0; i < pObj->buf.size; i++)
		{
			pObj->sum += pObj->buf.pData[i];
		}
	}
	return pObj->var;
}

void SlidingVarProcessV(sl_var_t *pObj, radio_float_t *pX, uint32_t len)
{
	uint32_t i;

	for (i=0; i < len; i++)
		SlidingVarProcess(pObj, *(pX++));
}

radio_float_t SlidingVarGet(sl_var_t *pObj)
{
	return pObj->var;
}

// --------------------------------------------------------------
// Sliding Min/Max
// --------------------------------------------------------------
void SlidingMinMaxInit(sl_minmax_t *pObj, uint32_t L, radio_float_t P, radio_float_t mode)
{
	pObj->L = L;
	pObj->P = P;
	pObj->N = 0;
	pObj->max = 0;
	pObj->mode = mode;
	pObj->pP = (uint32_t*)malloc(L*sizeof(uint32_t));
	pObj->pU = (radio_float_t*)malloc(L*sizeof(radio_float_t));
	pObj->pP_last = (uint32_t*)malloc(L*sizeof(uint32_t));
	pObj->pU_last = (radio_float_t*)malloc(L*sizeof(radio_float_t));
	memset(pObj->pP, 0, L*sizeof(uint32_t));
	memset(pObj->pU, 0, L*sizeof(radio_float_t));
	memset(pObj->pP_last, 0, L*sizeof(uint32_t));
	memset(pObj->pU_last, 0, L*sizeof(radio_float_t));
	pObj->pU_last[1] = P;
}
	
void SlidingMinMaxFree(sl_minmax_t *pObj)
{
	if (pObj->pP)
		free(pObj->pP);

	if (pObj->pU)
		free(pObj->pU);

	if (pObj->pP_last)
		free(pObj->pP_last);

	if (pObj->pU_last)
		free(pObj->pU_last);

	pObj->pP = (uint32_t*)0;
	pObj->pU = (radio_float_t*)0;
	pObj->pP_last = (uint32_t*)0;
	pObj->pU_last = (radio_float_t*)0;
}

// Sliding Min/Max based on MAXLIST-Algorithm, based on:
// "AN EFFICIENT ALGORITHM FOR RUNNING MAX MIN CALCULATION", 1996, S.C. Douglas, University of Utah
// Dependeding on 'mode' this algorithm calculates either minimum (mode = -1) or maximum (mode = +1)
radio_float_t SlidingMinMaxProcess(sl_minmax_t *pObj, radio_float_t x)
{
	uint32_t m;
	uint32_t i;
	uint32_t j;

	x *= pObj->mode;

	m = 1;
	if (pObj->pP_last[pObj->N] == pObj->L)
	{
		m = 0;
		pObj->pU_last[pObj->N] = pObj->P;
	}

	i = 0;
	while (x >= pObj->pU_last[i+1])
	{
		i++;
	}
	pObj->N = pObj->N - i + m;

	for (j=1; j < pObj->N; j++)
	{
		pObj->pU[j+1] = pObj->pU_last[j+i];
		pObj->pP[j+1] = pObj->pP_last[j+i] + 1;
	}
	pObj->pU[pObj->N+1] = pObj->P;
	pObj->pU[1] = x;
	pObj->pP[1] = 1;
	pObj->max = pObj->pU[pObj->N];

	i=0;
	while(pObj->pU[i] != pObj->P)
	{
		pObj->pU_last[i] = pObj->pU[i];
		pObj->pP_last[i] = pObj->pP[i];
		i++;
	}
	pObj->pU_last[i] = pObj->pU[i];
	pObj->pP_last[i] = pObj->pP[i];

	return pObj->mode*pObj->max;
}

radio_float_t SlidingMinMaxGet(sl_minmax_t *pObj)
{
	return pObj->mode*pObj->max;
}

void SlidingMinMaxProcessV(sl_minmax_t *pObj, radio_float_t *pX, uint32_t len)
{
	uint32_t i;

	for (i=0; i < len; i++)
	{
		SlidingMinMaxProcess(pObj, pX[i]);
	}
}
	
// --------------------------------------------------------------
void TraceOpen(trace_t *pObj, const char *pFilename)
{
	pObj->pFile = fopen(pFilename, "w");
}

void TraceClose(trace_t *pObj)
{
	if (pObj->pFile)
		fclose(pObj->pFile);
}

void TracePrint(trace_t *pObj, const char *fmt, ... )
{

    va_list a_list;

	if (!pObj->pFile)
		return;

    va_start( a_list, fmt );

	vfprintf(pObj->pFile, fmt, a_list);
}