Removing a non linear noise trend from data
Rather than modelling the non-linearity analytically, I've directly attempted to remove it. The longer trend features of the spectrum have been extracted by smoothing out most of the noise and signal, then subtracting these from the original spectrum.
A convenient way to do this is with an interpolation function, Interpolation
, constructed from the smoothed spectrum.
freqs = dat[[All, 1]];
vals = Abs@dat[[All, 2]];
smoothedVals = MeanFilter[vals, 5000];
if = Interpolation[{freqs, smoothedVals}\[Transpose]];
trdat = {First@#, Last@# - if@First@#} & /@ ({freqs, vals}\[Transpose]);
I've chosen a value of 5000 to samples to smooth the spectrum with as this gives what looks like to me the 'nicest' non-linearity. (3000 samples produces a nicer spectrum, but the non-linearity looks less attractive, try the values for yourself ). Below is a plot of the interpolation function for a 5000 sample window:
ListPlot[{First@#, if@First@#} & /@ ({freqs, vals}\[Transpose]), AxesLabel -> {"Freq", "Amplitude"};
Subtracting the non-linearity from the spectrum gives us this the spectrum with the non-linearity removed:
ListPlot[trdat, PlotRange -> All, AxesLabel -> {"Freq", "Amplitude"}]
You can get a version of the spectrum with noise more centred on zero by reducing the number of samples to smooth by, but the non-linearity looks more disjointed and it disturbs the level of the response peaks.
If you drop the smoothing window size down to 100 samples you begin to see the level shift associated with too small a sample window and more noise, and the signal, appearing in the non-linearity:
and
Further low pass filtering of the spectrum shows what maybe some harmonic interaction:
ListPlot[{freqs, MeanFilter[trdat[[All, 2]], 50]}\[Transpose],
PlotRange -> All, AxesLabel -> {"Freq", "Amplitude"}]
It may be worth noting that noise amplitude seems to increase with frequency, so it may not be Gaussian.
Credit goes to Ondřej Grover from the Signal Processing stackexchange. He suggested using a High pass filter, as it clearly seems to be some low frequency signal. His solution involved a Python script (found here) but something similar can be done with Mathematica.
Using a forward-backward method and some experimentation with cutoff values, we can use Mathematica's Highpassfilter
to fix our data (which I originally defined as Sxx
:
Sxxfiltered =
HighpassFilter[HighpassFilter[Abs[Sxx], 0.0003], 0.0003];
If we now look at what the data looks like, we see that instead of the curvey data we had before, we instead get a nice, straight line.
As you can see on the y-axis the scale is no longer as before, but this can be corrected for without too much trouble.