Trying to fit a cosine function to measured data...

[psmorris]

Hi all,

I'm having a little trouble... I'd like to use labview to fit a cosine function of the form:

y = a0*cos(a1*x + a2) + a3,

to some measured data. I've been trying to use the lev-mar fit but it keeps throwing up errors (although it does still work...) it's also extremely sensitive to the initial guesses, which is a tad annoying...

I've attached my current test code, it generates a cosine based on the function above, and then adds noise and then runs the fit... (its directly based on an example found in the ni examples list)

Any suggestions for improvements or an alternative approach would be greatly appreciated!

Thanks

Paul

[johnsold]

Paul,

One thing I have learned is that fitting to sine or cosine functions does not work well if you have multiple cycles of the waveform in the data. It appears that the methods used do not change frequencies very well. Even if the initial estimate of frequency is very close, sometimes it does not work. It does seem to work fairly well for fractions of a cycle as shown in your example.

This is something I have been working on (at a low level of effort) for a couple of years. Unfortunately, I cannot offer any good suggestions yet.

Lynn

[altenbach]

Periodic functions can be very hard to fit if your guesses are wrong, because you can equally well fit to an alias frequency or a phase shifts corresponding to an integer multiple of the wavelength.

Anyway, If your initial estimates are reasonable, you should have no problem. In your case, just graph also the function calculated with the initial estimate and you can see that it is almost inverted (A1=-01, vs. A1=0.1, different sign!). Here's it can roll either way or go off into the woods. You can even interactively adjust the guesses until things look similar before fitting.

If your data typilcally looks similar (less than one full wave, one frequency, etc.) you should calculate the initial estimates from features in the data (average Y, max(y), min(y), x(max Y), x(min(y), etc.) Now things will be much more robust!

To make sure we're not dealing with an alias frequency, I also typically graph the final function calculated with much denser x-values and use an xy graph.

You might also use wires for the function as follows. Looks nice and is easier to debug.

Message Edited by altenbach on 04-17-2008 10:16 AM

[DSPGuy]

A reasonable starting point can be usually be found using the Extract Single Tone information.vi, and Basic Averaged DC-RMS.vi.  Also, you might consider using the constrained version of the Nonlinear Fit to control run-away parameter estimates.  Attached is a modified version of your VI that shows these ideas.

-Jim

[altenbach]

---

DSPGuy wrote:
A reasonable starting point can be usually be found using the Extract Single Tone information.vi, and Basic Averaged DC-RMS.vi. 

---

I was hesitant to suggest that here. How well do these tools work if we have significantly less than a full wave period in the data as in this case? (To test, set the number of points to e.g. 39 or less and it seems fail 100% even with zero noise).

OTOH, if we have many periods of data, we probably don't even need to fit and can get everything from the FFT based analysis directly. 🙂

[psmorris]

Thanks for all the replies!

I'll explain a little more what I'm doing... basically, I'm measuring the back reflection from an optical interferometer as a function of the light wavelength. As the interferometer is low finesse, the reflectivity varies as a cosine. I need to find the period of said cosine in order to calculate the thickness of the interferometer. It's highly possible that there will be less than one full period of cosine, or at least, there will likely be no more than one full period. Hence I stayed away from the "single tone information" approach

What I've done so far is to put the fitting code into a loop which displays the fit result and allows me to manually update the initial parameters until I get a fit... The data ought to be pretty "clean" in terms of noise so I'm hoping I can just give good enough initial values and progressively improve the fit by hand. I've set up the code to calculate A0 and A3 (the first and last fit parameters) based on the measured amplitude and offset of the cosine, so that should help (although if there's less than a full period, that may not work well... i'll just have to see I guess...

Anyways, thanks again for all the replies...

Paul

[johnsold]

Paul,

If you have two zero crossings you can get a pretty good estimate of the frequency from that. Similarly, the locations of the zero crossings in the dataset can give an estimate for the phase. Since you will be evaluating the fit manually, the time it would take the software to calculate these estimates would be insignificant.

Lynn

[altenbach]

---

johnsold wrote:
If you have two zero crossings you can get a pretty good estimate of the frequency from that.

---

The problem here is that he also fits for an offset in the data (a3). Unless this value is much smaller than the amplitude (a0), your idea will not work for a cosine with less than one period.

Still, I am pretty sure a good initial estimate can be made so the fit can converge successfully. (look for max, min, find places where the derivative in near zero, etc.). We also have an excellent estimate on the frequency, e.g. assuming about half a period. In the worst case, we could even do a grid search and try about 3 different phases.

[johnsold]

I forgot about the offset. As you said, with an offset my ideas won't work.

Lynn