Just a note that Ive uploaded the initial version of vfuncs to google code. Ive released under a BSD license so you can use it in your commercial and noncommercial code easily.
Download from here [I’ll import to SVN sometime soon]. See my previous post for a description of vfuncs.
This version contains an example of a digital filter. This can be used to smooth the series data, or apply other signal processing operations. If your familiar with applying a blur filter in photoshop or gimp, using a gaussian filter kernel, this is exactly the same idea (except in one dimension). Gaussian filter is basically just a moving average of the data.
Think of the algorithm as applying a sliding window across the data – the sliding window contains the filter weights, and at each position you apply the weighted average [dot product] of the filter weights against each data point in the window.
If the filter contains a single element of weight 1.0, then the result is just the input (the filter is just the Dirac delta function in that case). If the filter contains [0.25 0.50 0.25] its going to mix each element with its neigbours and take a weighted average, thus smoothing the data.
To make wider gaussian filters, I use a primitive to generate pascals triangle or order n [eg 1 4 6 4 1 for order 4] and then normalize so it sums to 1.0. Normalizing to 1.0 means the average is unchanged, ignoring handling of the ends for the moment. See the main.c for comments.
This naive implementation can smooth series at more than 6million items per second on a current laptop. My initial impl just allocates a single large array, and its easy to make copies of that [vnewXXX] – so my approach has been to write readable code, get it working then use valgrind and printf debug output to check allocations. In production code you’d want to have a demand load system… and probably also reference count or snapshots to avoid many explicit calls to vdel. Currently, for each V a=vnewXXX() youll need to call vdel(a).
As an aside… filters are a subset of convolution, a simple and important mathematical technique. Also notice that as well as smoothing – which is a low pass filter, as it removes peaks [like integration operator] – you can also accentuate the peaks in some useful way. For example, the filter [-1.0, 1.0] pulls out the volatility between adjacent items while bringing the average down to zero. This is the simplest hi-pass filter [akin to differentiation operator] … think edge detection in 2D images. Actually you can approximate any of the nth derivatives with the right filter… but thats enough for todays rant.