Data interpolation in python
Quite often in my day to day work I end up having to flag out bad data. In many cases this just creates a gap and normally this is fine just to throw away. There are occasions however, such as when creating the bandpass of reciever that its still quite useful to have some idea as to what is going on in the normally flagged out channels – especially when its large chunks of data. To solve this we use interpolation. This post isn’t about explaining data interpolation – for that I suggest you take a look at wiki.
When interpolating I like to use splines. In linear interpolation we use a linear function for each intervale,s with splines we use low-degree polynomials in each of the intervals chosing a polynomial that nicely fits the sections together. Generally works quite well.
Python has many different ways todo interpolation, be it spline or not. There is a huge variety of interpolation functions and I’d urge you to explore the scipy help for an exhaustive list. In this example I’m going to use the spline functions from scipy.signal
Firstly, import the functions required (asciidata is to read in the data and is not needed):
from scipy import *
from numpy import r_, sin
from scipy.signal import cspline1d, cspline1d_eval
I’ve defined a function to do the actual work:
lastzero = max(where(beamdata == 0)) #first 146 channels are 0
startchan = lastzero+1
x = r_[startchan:len(beamdata)]
dx = x-x
newx = r_[0:startchan:1] # notice outside the original domain
y = beamdata[startchan:]
cj = cspline1d(y,lamb=0)
newy = cspline1d_eval(cj, newx, dx=dx,x0=x)
aa = list(newy)
bb = list(y)
mod = aa + bb
The function is called like:
output = spline_lowchan(originaldata)
Where originaldata is the initial bandpass (just a sequence of numbers in an array or list)
In the above example of interpolating the function is assuming that the first set of data is all set to 0, hence the where is 0. We then interpolate the data over the range where there is no signal. We then add this to the original signal, hence allowing the region without any useful data to have an interpolated data set. Not perfect but gives a better idea of what might be happening. Of course this is not really how you should treat the edges but the cspline1d method very well when you say lose a few data points out of a bandpass.