Problem with fitcurve in Igor
pjcau
only recently I started programming with Igor and I would like a little help. I open with Waveload with a file. Dat in 2D. Intact across the surface and I get a spectrum. So far so good. Now I wish I did the fit of this spectrum, with double exponetial between two intervals, the highest point and the end of the curve in profile.
How could we do?
Thanks to all
(1) you have data that is exponential for small and for large values of x, but not in between
(2) you suspect that the exponential decay is described by different decay constants tau1 and tau2
(3) you want to fit your data using two non-contiguous subranges, one at small x and one at large x
If that's correct, I think there are two ways to go about it:
(1) the easy way is to simply fit your data twice with a single exponential function, using a subrange
(2) the hard way is to fit your data using a mask wave (recommended) and also, possibly, weighting
Method (1) is appropriate if you are interested in the decay constants only (and if you need them both). If you intend to fit the data in order to subtract a background, you may need method (2). This method is also appropriate if you really wanted to fit a single tau parameter, but for different ranges of your data.
ad (1): To learn more about fitting to a subrange, execute
DisplayHelpTopic "Fitting a Subset of the Data"on the command line. Briefly, in the curve fitting dialog, on the "data options" tab, you can specify the point range manually, or you can click the "cursors" button. This button is active only when you have the cursors displayed in your graph (bring the graph to the front and then press Ctrl+I ...). Place the cursors to delimit the range of interest at the beginning (small x), then do the curve fit with a SINGLE exponential function. Adjust the cursors and repeat the curve fitting. For the second fit, you might need to hold some of your coefficients (especially the constant offset), the exp_Xoffset function will probably give wrong results. Keep in mind that fitting an "exponential tail" can be very tricky and give completely wrong results if you have noise or an unknown vertical offset in your data (see below for an example).
ad (2): If you abolutely need to fit everything at once, you can use a mask wave (recommended) and possibly add weighting to your data. You can look up these topics by executing
DisplayHelpTopic "weighting"orDisplayHelpTopic "Using a Mask Wave".I suspect that you might need this from your message (clue "the first point" - this makes sense only if you want a single exponential). Here's a short step-by-step procedure, anyway. Bring your graph to the front. Place the graph cursors so that they delimit the "end of the curve or profile". Then, execute on the command line
mask=0; mask[0]=1; mask[pcsr(A),pcsr(b)]=1
To visualize your mask, you can execute
ModifyGraph mode(mask)=7,hbFill(mask)=4,lsize(mask)=0;DelayUpdate
ModifyGraph usePlusRGB(mask)=1,plusRGB(mask)=(0,0,65280)
To fit your data with a single exponential using the mask wave, use
CurveFit/L=(numpnts(CsrWaveRef(A)))/X=1/H="100"/TBOX=768 exp CsrWaveRef(A) /M=mask /DTo subtract the fit from the data, use
$CsrWave(A) -= $("fit_"+CsrWave(A))To undo the subtraction, use
$CsrWave(A) += $("fit_"+CsrWave(A))If this is not satisfactory, you can add further regions to your mask wave by adjusting the cursors and again using
mask[pcsr(A),pcsr(b)]=1(remove regions by sayingmask[pcsr(A),pcsr(b)]=0).If the fit is still bad (overemphasizing the initial region, as exponential fits are bound to do), you can add weighting to your data. One possibility (quick and dirty) is to use
CurveFit/L=(numpnts(CsrWaveRef(A)))/X=1/H="100"/TBOX=768 exp CsrWaveRef(A) /W=weight /I=1 /M=mask /D
Best of luck,
Wolfgang Harneit
P.S.: A short example of fitting pitfalls with exponential tails
[Preparation]
On the command line, execute
Display dexp; ModifyGraph log(left)=1; ShowInfo
This puts up a nice double exponential curve with distinguishable decay constants tau1 and tau2, an "interesting" gaussian peak in the middle, plus some noise and a vertical offset. I use a semi-logarithmic plot for better inspection of the data. Note that there is significant "noise rectification", i.e. the curve decays rather slowly at the end (to force this, I included the small constant vertical offset). We will see this better in a bit.
[Fitting the initial part]
To fit the initial decay, place the cursors on the first point in the graph (point 0) and on the fourth one (point 3), or execute
Cursor/P A, dexp, 0; Cursor/P B, dexp, 3Choose "exp" from the "QuickFit" submenu in the "Analysis" menu. Usually, this also puts up the fitting coefficients in an annotation box in your graph; if not, you can look at the history in the "command" window. The "invTau" parameter is fairly accurate (0.5 +/- 0.02), even though we fitted only 4 data points with 3 independent parameters!
[Fitting the tail]
To fit the final decay, it is "natural" to place the cursors somewhere after the gaussian peak. Choosing again "exp" from the "QuickFit" menu, we however get an "invTau" parameter that is wide off the mark (it should be 0.05), even when we fit many data points. You may have to vary the noise a bit to see this effect.
[Repeat with different noise]
The easiest way to prepare the noisy data and do the fit repeatedly is to do it once and then use the command/history window. I now assume that you have the graph window visible, have placed the cursors somewhere on the tail as mentioned before, and have done the quickfit. Now, bring up the command window again by pressing Ctrl-J, then press the "up" button of your keyboard as many times as needed to bring up the line where it says
Make/D/O/N=101 dexp=0.02+enoise(0.02)+0.4*exp(-x/20)+exp(-x/2)+3*gauss(x,30,3)Then, simply press "return" twice. The first press copies the line from the history to the command line, the second press executes the line.
Repeat the "up-up...return-return" scheme searching for the line saying something like
CurveFit/M=2/W=0/TBOX=(0x300) exp, dexp[55,91]/DThis line was placed in the command window by the quickfit routine. Your /TBOX parameter may be different, and the data range (in [] brackets) probably will be. For further repetitions, you just need "up-up-return-return" if you don't do anything else in between. Once you get the rhythm, you can keep your eye on the fitting parameter in the annotation box and do the keyboard manipulation blindly.
[How to get better fit parameters]
If you know the vertical offset in your data a priori, you can hold the corresponding fit parameter. In the ordinary (not quick-) curve fitting dialog (with the "exp" function selected), the "coefficients" tab tells you that it is the first parameter, called "y0". You can click on the little "Hold" box and type in the a-prioi offset. Then, execute the fit. Or you can type on the command line
K0=0;CurveFit/Hold="100"/M=2/W=0/TBOX=(0x300) exp, dexp[55,91]/DNote that Igor's internal name for the first parameter is "K0" no matter what the fitting function is; also note the "/Hold=pattern" parameter. Even using the incorrect guess "there is no appreciable offset in my data" (K0=0) gives slightly better results than before.
[The best way...]
...is probably to use the more slowly decaying part before the gaussian peak, even if this region is small:
K0=0;CurveFit/Hold="100"/M=2/W=0/TBOX=(0x300) exp, dexp[12,20]/Dor at least, to avoid the region of high noise:
K0=0;CurveFit/Hold="100"/M=2/W=0/TBOX=(0x300) exp, dexp[39,54]/DBut maybe you have your own experiences and opinions here. At last, how could you have told that this is the best way to do it without knowing the proper function in advance?
January 14, 2010 at 01:25 pm - Permalink
but what I do in my opinion is much simpler. I have a file IPF, where I wrote a routine, but I would like to upload more than that now an array of data and extract a profile, let me note that I fitcurve that profile. As you can see from the figure that I have attached the details of which I would fit were the highest point of the profile and the last point of the profile. How can I write without using the mask but with pcsr ().
Sorry if I was unclear.
Thanks again
January 15, 2010 at 06:16 am - Permalink
Are you asking how to do this fit using Igor code (rather than manually)? If you can place cursors manually, you can put the A (round) cursor at the peak, and the B (square) cursor on the last point, then do a fit between cursors.
To do the equivalent thing from Igor code, the hard part is finding the peak. If that point is always the largest Y value in the wave, then you can use WaveStats/P to find the point number of the maximum point. If there is noise, that might not be completely reliable. In that case you may get better results using FindPeak.
Then use CurveFit on a subrange:
CurveFit/N/W=0 dblExp_XOffset, ywave[V_peakLoc, numpnts(ywave)-1] /X=...
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
January 15, 2010 at 02:26 pm - Permalink