# Fitting complex functions

mtaylor

Properties that should make this possible:

1. Complex functions take in the same parameters and independent variable for the real and imaginary parts. This simplifies the situation somewhat.

2. The independent variable is the same wave for the real and imaginary parts.

With this in mind, all I need to do is format my input data appropriately and construct my fit function such that

1. It is an all-at-once function. Input waves : parameter wave, frequency wave. Output wave: calculated value wave.

2. The calculated wave is constructed like this: {r0,r1,r2,r3...rn,i0,i1,i2,i3...in} where ri is the real value for point i and ii is the imaginary value for point i.

I have had success with all-at-once functions when I don't have to take apart the input waves or assemble the output wave. Here is an example:

Function ZrealQandles2U(w, realout, frequency) : FitFunc

Wave w, realout, frequency

//CurveFitDialog/ Independent Variables 1

//CurveFitDialog/ frequency

//CurveFitDialog/ Coefficients 7

//CurveFitDialog/ w[0] = RP

//CurveFitDialog/ w[1] = Qdl

//CurveFitDialog/ w[2] = ndl

//CurveFitDialog/ w[3] = RS

//CurveFitDialog/ w[4] = Qbl

//CurveFitDialog/ w[5] = nbl

//CurveFitDialog/ w[6] = Rbl

Variable RP,Qdl,RS,ndl,Qbl,Rbl,nbl

RP=w[0]

Qdl=w[1]

ndl=w[2]

RS=w[3]

Qbl=w[4]

nbl=w[5]

Rbl=w[6]

realout = RS+ real(1/(1/RP +Qdl*(2*Pi*10^frequency*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^frequency*sqrt(-1))^nbl))

End

Wave w, realout, frequency

//CurveFitDialog/ Independent Variables 1

//CurveFitDialog/ frequency

//CurveFitDialog/ Coefficients 7

//CurveFitDialog/ w[0] = RP

//CurveFitDialog/ w[1] = Qdl

//CurveFitDialog/ w[2] = ndl

//CurveFitDialog/ w[3] = RS

//CurveFitDialog/ w[4] = Qbl

//CurveFitDialog/ w[5] = nbl

//CurveFitDialog/ w[6] = Rbl

Variable RP,Qdl,RS,ndl,Qbl,Rbl,nbl

RP=w[0]

Qdl=w[1]

ndl=w[2]

RS=w[3]

Qbl=w[4]

nbl=w[5]

Rbl=w[6]

realout = RS+ real(1/(1/RP +Qdl*(2*Pi*10^frequency*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^frequency*sqrt(-1))^nbl))

End

here is an example of what I'd like to do, but it doesn't seem to be working:

Function ZrealQandles2U(w, compout, frequency) : FitFunc

Wave w, compout, frequency

Wave temp1,temp2,realout,imagout

variable halfpoint,endpoint

//CurveFitDialog/ Independent Variables 1

//CurveFitDialog/ frequency

//CurveFitDialog/ Coefficients 7

//CurveFitDialog/ w[0] = RP

//CurveFitDialog/ w[1] = Qdl

//CurveFitDialog/ w[2] = ndl

//CurveFitDialog/ w[3] = RS

//CurveFitDialog/ w[4] = Qbl

//CurveFitDialog/ w[5] = nbl

//CurveFitDialog/ w[6] = Rbl

Variable RP,Qdl,RS,ndl,Qbl,Rbl,nbl

RP=w[0] //assign easy to code variable names

Qdl=w[1]

ndl=w[2]

RS=w[3]

Qbl=w[4]

nbl=w[5]

Rbl=w[6]

//calculate some things about the input wave structure

halfpoint=numpnts(frequency)/2 //this is always an integer, because of the input structure.

endpoint=numpnts(frequency)

//make two temporary frequency waves

duplicate frequency temp1

duplicate frequency temp2

//make the output wave

duplicate frequency compout

//trim the temporary frequency waves

deletepoints halfpoint, halfpoint, temp1

deletepoints 0, (halfpoint-1), temp2

//compute the real part

realout = RS+ real(1/(1/RP +Qdl*(2*Pi*10^temp1*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp1*sqrt(-1))^nbl))

//compute the imagimary part

imagout = imag(1/(1/RP +Qdl*(2*Pi*10^temp2*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp2*sqrt(-1))^nbl))

//combine them into the output wave

concatenate/np/o {realout,imagout},compout

//clean up temporary waves

killwaves realout, imagout, temp1, temp2

End

Wave w, compout, frequency

Wave temp1,temp2,realout,imagout

variable halfpoint,endpoint

//CurveFitDialog/ Independent Variables 1

//CurveFitDialog/ frequency

//CurveFitDialog/ Coefficients 7

//CurveFitDialog/ w[0] = RP

//CurveFitDialog/ w[1] = Qdl

//CurveFitDialog/ w[2] = ndl

//CurveFitDialog/ w[3] = RS

//CurveFitDialog/ w[4] = Qbl

//CurveFitDialog/ w[5] = nbl

//CurveFitDialog/ w[6] = Rbl

Variable RP,Qdl,RS,ndl,Qbl,Rbl,nbl

RP=w[0] //assign easy to code variable names

Qdl=w[1]

ndl=w[2]

RS=w[3]

Qbl=w[4]

nbl=w[5]

Rbl=w[6]

//calculate some things about the input wave structure

halfpoint=numpnts(frequency)/2 //this is always an integer, because of the input structure.

endpoint=numpnts(frequency)

//make two temporary frequency waves

duplicate frequency temp1

duplicate frequency temp2

//make the output wave

duplicate frequency compout

//trim the temporary frequency waves

deletepoints halfpoint, halfpoint, temp1

deletepoints 0, (halfpoint-1), temp2

//compute the real part

realout = RS+ real(1/(1/RP +Qdl*(2*Pi*10^temp1*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp1*sqrt(-1))^nbl))

//compute the imagimary part

imagout = imag(1/(1/RP +Qdl*(2*Pi*10^temp2*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp2*sqrt(-1))^nbl))

//combine them into the output wave

concatenate/np/o {realout,imagout},compout

//clean up temporary waves

killwaves realout, imagout, temp1, temp2

End

Unfortunately, it doesn't work. The fit proceeds as if there are no problems, but the parameters do not vary during the fit so there is obviously something wrong here.

Is there something fundamentally flawed with my method here? Does anyone have any experience fitting complex data without the aid of the global fit panel?

Wave w, compout, frequency

variable halfpoint, temp1, temp2, ii

//CurveFitDialog/ Independent Variables 1

//CurveFitDialog/ frequency

//CurveFitDialog/ Coefficients 7

//CurveFitDialog/ w[0] = RP

//CurveFitDialog/ w[1] = Qdl

//CurveFitDialog/ w[2] = ndl

//CurveFitDialog/ w[3] = RS

//CurveFitDialog/ w[4] = Qbl

//CurveFitDialog/ w[5] = nbl

//CurveFitDialog/ w[6] = Rbl

Variable RP,Qdl,RS,ndl,Qbl,Rbl,nbl

RP=w[0]

Qdl=w[1]

ndl=w[2]

RS=w[3]

Qbl=w[4]

nbl=w[5]

Rbl=w[6]

halfpoint=numpnts(frequency)/2 //this is always an integer, because of the input structure.

for(ii = 0 ; ii < halfpoint ; ii += 1)

temp1=frequency[ii]

compout[ii]=RS+ real(1/(1/RP +Qdl*(2*Pi*10^temp1*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp1*sqrt(-1))^nbl))

endfor

for(ii = 0 ; ii < halfpoint ; ii += 1)

temp2=frequency[ii+halfpoint]

compout[ii+halfpoint] = imag(1/(1/RP +Qdl*(2*Pi*10^temp2*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^temp2*sqrt(-1))^nbl))

endfor

End

But I got the same behavior.

February 13, 2012 at 08:34 am - Permalink

Your second attempt is basically what you should be doing. You can make it more efficient using wave assignments:

compout[halfpoint, ] = imag(1/(1/RP +Qdl*(2*Pi*10^ frequency*sqrt(-1))^ndl)+ 1/(1/Rbl +Qbl*(2*Pi*10^ frequency*sqrt(-1))^nbl))

But that won't change the result. I find it helpful with things like this to run my fit function outside of FuncFit to make sure that it's producing reasonable output. With an all-at-once function it's a bit tricky: you have to create appropriate wave for compout and frequency in order to run it. But it's easier to debug a fit function on its own than in the context of an actual fit.

John Weeks

WaveMetrics, Inc.

support@wavemetrics.com

February 13, 2012 at 11:14 am - Permalink

Thank you for your assistance. I used the wave assignment method you suggested, which shortened my code considerably. I did some testing and It appears that my function is working as expected outside of any attempts to fit it, so it must be an issue with the function calling the fit. I will look into this and see if I can figure it out.

February 13, 2012 at 01:29 pm - Permalink

`sqrt(-1)`

. For example, in your function define once, and then useThis saves frequent recalculation of a complex square-root in several places. It is a pity that Igor's method of Constant definition does not allow complex-number assignments.

February 14, 2012 at 08:47 am - Permalink

February 14, 2012 at 03:31 pm - Permalink

1. format your data such that the complex part has been seperated into two real waves. concatenate those waves using the /np flag to make a double-length wave containing your data.

2. Do the same for your independent variable, this is usually frequency, making a double-length wave of your frequency spectrum.

3. Format your fit function as below:

wave complexcat, w, logfreq //declare wave names. Don't edit this

variable halfpoint=numpnts(logfreq)/2 //finding the break between the real and imaginary data. don't edit this

variable/c ki = cmplx(0,1) //declare complex constant. Don't edit this. Do use "ki" as a replacement for sqrt(-1)

variable , Rs, Rp, Cdl //declare variable names; you can edit this. This is optional if you choose to use w[n] instead of easy to read variable names in your function.

Rs=w[0] //assign your variable names to the solution wave. Edit this.

Cdl=w[1]

Rp=w[2]

//put your function in here. In my case, I've got a simple function, so it just writes complexcat on the next lines.

complexcat[0,halfpoint-1]=Rs+real(1/(1/Rp+1/(ki/-(Cdl*2*Pi*10^logfreq)))) //use real() to get the real part for the first half

complexcat[halfpoint, ]=imag(1/(1/Rp+1/(-ki/(Cdl*2*Pi*10^logfreq)))) //use imag() to get the imaginary part for the second half

//note that no return command is needed for an all-at-once function.

End

Now use it as you normally would with the global fit dialog, the curvefit dialog, or whatever you use to fit curves these days (I use gencurvefit, a project on igorexchange.)

Thanks to everyone for helping me put this together as I stumbled through the early days of my igor eduction.

I sincerely hope that this saves future researchers from a major headache.

February 16, 2012 at 09:11 am - Permalink

Function Execution Error

While executing a wave read, the following error occured: index out of range.

At andyfaff's suggestion, I tried turning on the procedure debugger, but don't seem to be getting any breaks when running the function, so I started commenting out lines to see where the problem was.

1. the function does not give the error when run independently from gencurvefit

2. if I comment out the line marked below, where the imaginary part is assigned, the error goes away (of course, the curve fit is then bogus)

3. however, if I similarly comment out the line where the real part is assigned, the error stays put.

4. if I replace the last 3 lines with the two commented lines at the end, there is no error.

5. generating the complex wave during the function execution instead of referencing a pre-created wave does not fix the issue.

So I think that I am making a mistake in my wave assignment on that line, specifically "complexer[p-midpointt]", but I can't seem to sort it. I've obtained the same error message with various index shifts and using my best knowledge of igor syntax, it should be correct. The function's code works in its current form, giving a correct result as far as I can see, but kicks out the error message when the fit is complete. All I can figure is something different happens when fitting an allatonce function that I didn't account for. Even more puzzling is the fact that the curve fit gives a perfect result despite the error message. Can anyone shed some light on what may be happening here? Thanks in advance.

wave complexcat, w, freqw //declare wave names. Don't edit this

wave/d/c complexer //pull in temporary complex wave

variable/g midpointt //pull in midpoint

variable/g/c ki, kii //pull in complex constant

variable Rs, Rp, Qdl,n,sigmaW //declare variable names; you can edit this

Rs=w[0] //assign your variable names to the solution wave. Edit this.

Qdl=w[1]

Rp=w[2]

n=w[3]

sigmaW=w[4]

//computing complex and then assigning

complexer[0,midpointt-1] = Rs+1/(1/(Rp+sigmaW*kii/sqrt(freqw))+(Qdl*(ki*freqw)^n))

complexcat[0,midpointt-1]=real(complexer[p]) //if I comment out this line, I still get the error.

complexcat[midpointt, ]=imag(complexer[p-midpointt]) //if i comment out this line, I no longer get the error.

//computing seperately

//complexcat[0,midpointt-1]=real(Rs+1/(1/(Rp+sigmaW*kii/sqrt(freqw))+(Qdl*(ki*freqw)^n)))

//complexcat[midpointt, ]=imag(Rs+1/(1/(Rp+sigmaW*kii/sqrt(freqw))+(Qdl*(ki*freqw)^n)))

End

April 10, 2012 at 12:31 pm - Permalink

Be sure that in addition to enabling the debugger, go back to the procedure window right-click menu and enable Debug on Error. It can also help lots of bugs (probably not this one) to enable NVAR, SVAR Wave Checking.

John Weeks

WaveMetrics, Inc.

support@wavemetrics.com

April 10, 2012 at 02:14 pm - Permalink

April 10, 2012 at 04:45 pm - Permalink

You should never resize/kill/delete/change the type of/etc the waves supplied to the fit function. This will cause an error. The fit waves are not only calculated at the end, they can also be calculated during the fitting process. The fits are (by default) 200 points long and span the range of the input data. The error you were getting is not a gencurvefit problem.

From the help file:

April 10, 2012 at 05:55 pm - Permalink

during normal function evaluations, the value printed was correct. Periodically (I assume at the end of a generation), it would output as the value of L (default 200). So my code was measuring the length of the wave to be different from what was expected.... and this is what was causing the error message as evidenced by setting L to the same length as complexcat (wy) eliminating the error message, because the wave never changed length. Another way around this would be to calculate the length of the wave within the function, but that wastes a few cpu cycles, so I decided to calculate it beforehand and just reference it.

So yes, you are correct that it is not a bug in gencurvefit, but it is a good thing to know about if you are using a non-standard problem formulation.

April 10, 2012 at 06:37 pm - Permalink

The bottom line is that you should not be depending on a given numpnts for the dependent and independent waves, they are variable in size (for gencurvefit and funcfit). You should work out the wave length in the fit function.

April 11, 2012 at 12:10 am - Permalink

August 20, 2015 at 04:40 pm - Permalink