Multivariate fitting with structure fit function
billyblack
Does anyone know if this can be done, or where I may be going wrong? The basic code I'm using follows. Thanks!
Structure mdfitstruct
wave coef
wave yw
wave xw1
wave xw2
wave theorycalc
endstructure
Function mdfitfunc(s) : FitFunc
struct mdfitstruct&s
s.yw[][]=... //function of s.theorycalc, s.xw1, s.xw2
end
Function mdfitdriver(coefwave,datamat,theorymat)
wave coefwave
wave datamat
wave theorymat
struct mdfitstruct fs
wave fs.theorycalc=theorymat
// The following lines are to make waves with appropriate x-values, if necessary
// make/o/n=(dimsize(datamat,0)) timescale
// timescale[]=p*dimdelta(datamat,0)
// make/o/n=(dimsize(datamat,1)) numlines
// numlines[]=p
FuncFitMD mdfitfunc, coefwave, datamat /D /STRC=fs
end
wave coef
wave yw
wave xw1
wave xw2
wave theorycalc
endstructure
Function mdfitfunc(s) : FitFunc
struct mdfitstruct&s
s.yw[][]=... //function of s.theorycalc, s.xw1, s.xw2
end
Function mdfitdriver(coefwave,datamat,theorymat)
wave coefwave
wave datamat
wave theorymat
struct mdfitstruct fs
wave fs.theorycalc=theorymat
// The following lines are to make waves with appropriate x-values, if necessary
// make/o/n=(dimsize(datamat,0)) timescale
// timescale[]=p*dimdelta(datamat,0)
// make/o/n=(dimsize(datamat,1)) numlines
// numlines[]=p
FuncFitMD mdfitfunc, coefwave, datamat /D /STRC=fs
end
When you post the experiment file, please include a note with a guide to what's in it. Be sure to include all the procedure code needed to run it. If some code is in external procedure files, you can "adopt" them by holding down the shift key while pulling down the File menu. Select Adopt All and make sure that "Adopt User Procedure Files" is checked, and that "Adopt WaveMetrics Procedures" is not checked.
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
October 30, 2013 at 09:01 am - Permalink
Hi John, thanks for your prompt response. I've cut down the experiment to mostly essential stuff. You can look through the history and main procedure file as well to see a bit of how I got to where I am. Hopefully this is set up in an understandable way - I'm new to this forum, so please let me know if you require other information. Below is a short guide to the experiment, which is also included as a notebook in the pxp file.
In this first attempt, I am simply trying to fit my experimental data to some theory calculations. The theory is necessarily calculated outside of Igor, since the calculation for a single set of parameters can take several weeks. In this simplest case, the only remaining adjustable fit parameter is the time-scaling (x-axis) of the theoretical calculations.
So I have a matrix of experimental data called “w2expmat” (NxD, where N is the number of points in time, and D is the number of detector positions in my experiment – here I have 3), and I also have a matrix of theoretical calculations called “w2calcmat” (MxD, where M is just the number of points for which I have calculated the theory, and D is as above). The coefficient wave, called “td” has only one element. The function that drives the fitting is called “DoW2Fit”, and I think the necessary inputs are easily understood from what I’ve said here.
I would like to fit all D columns of the data with the same fit coefficient, which scales all of the theoretical data in the same way. Perhaps there is a better way to do this than I am trying so far, but the structure fit function seemed to be the most appropriate for reading in the theoretical calculations to be interpolated (this will all eventually be a part of a much bigger procedure). I would eventually also like to include weighting by another matrix, and even a global fit with two different data/theory sets with some joint parameters.
I’ve also included a multi-dimensional data viewer GUI that I built, just because I find it convenient to use to view our data (I have an Image Plot version as well; I can provide more info on these or share code if anyone is interested).
October 30, 2013 at 12:04 pm - Permalink
It may seem strange, but when you fit to data in a matrix using FuncFitMD, you don't get a matrix as the yw wave in your fit function. Igor's curve fitting internals are geared to the most general case, that of points randomly scattered within the N-dimensional space represented by your N independent variables (FuncFit and FuncFitMD use the same fitting engine). So your matrix w2expmat gets unfolded into a 1D wave with rows*columns points in it. The X waves then contain the corresponding scaled X and Y values from w2expmat wave scaling, also unfolded. So when you write the assignment statement for s.yw, you need to computationally fold the p values on the RHS of the wave assignment into p and q values appropriate for s.inputw.
Second, for a multivariate fit to a structure fit function, you pass an array of x waves, not a list of x waves. So instead of this:
wave coef
wave yw
wave xw1
wave xw2
wave theorycalc
endstructure
you need this:
wave coef
wave yw
wave xw[2]
wave theorycalc
endstructure
That gives you an X wave for each independent variable. You can access then within the fit function like
Wave xw2 = s.xw[1]
I'm not sure the () subrange notation for X and [] for Y works with FuncFitMD. I recommend computing the row numbers corresponding to your time range and using [] subranges.
Be aware that since you are expecting to fit exactly 3 columns, you will need to avoid the /D flag for the fit autodestination wave.
I have attached a copy of your experiment file in which I attempted to fix most of this. I wasn't able to untangle the relationships of the X scaling if the exp wave compared to the computed wave; instead as a placeholder I simply used a "dilation factor"- that is, a factor used to go from rows in the exp wave to rows in the computed wave. You will need to figure out the appropriate computation yourself.
You must be careful not to try to access points outside the size of the computed wave. You will see that I used the min() function to avoid going off the end of the computed wave.
I found that with my dilation factor I needed an epsilon wave. The actual value isn't critical; I found that 1e-6 works well. I can tell you all about epsilon waves if you want...
Hope you can figure out what I did!
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
October 30, 2013 at 03:38 pm - Permalink
variable pdense=8
setscale/p x 0,s.coef[0]/pdense,"s" s.theoryinw2
wave xw0=s.xw[0]
wave xw1=s.xw[1]
s.yw[]=s.theoryinw2(xw0[p])[xw1[p]]
The () subrange notation worked just fine for FuncFitMD as far as I can tell, I used the following line without any problems:
I didn't encounter any of the problems you did with out-of-range indices or epsilon waves, but I'll keep what you have here as a reference if things start to go awry.
Now, is it possible to have several functions like this and fit several data sets with joint parameters using the Global Fit package? I can't seem to access it from the fit function menu in the Global Analysis dialog, but admittedly I haven't dug much deeper than that yet...
October 30, 2013 at 05:51 pm - Permalink
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
October 31, 2013 at 09:06 am - Permalink
Regards,
-Scott
#pragma rtGlobals=3 // Use modern global access method and strict wave access.
Function Gauss2Dtest2(pw,yw,xw) : FitFunc
wave pw,yw,xw
//pw[0]=Amplitude
//pw[1]=center frequency t
//pw[2]=sigma t standard deviation
//pw[3]=centerfrequency u
//pw[4]=sigma u standard deviation
yw=pw[0]*exp(-(xw[p][0]-pw[1])^2/2/pw[2]^2)*exp(-(xw[q][1]-pw[3])^2/2/pw[4]^2)
end
Function Gauss2Dtest2Driver(pw,xw,yw)
wave pw,yw,xw
FuncfitMD Gauss2Dtest2, pw,yw
end
April 20, 2015 at 06:23 pm - Permalink
I admit that the info on all-at-once multivariate functions is not as detailed as it might be...
You need one x wave for each independent variable:
wave pw,yw,xw1, xw2
//pw[0]=Amplitude
//pw[1]=center frequency t
//pw[2]=sigma t standard deviation
//pw[3]=centerfrequency u
//pw[4]=sigma u standard deviation
yw=pw[0]*exp(-(xw1[p]-pw[1])^2/2/pw[2]^2)*exp(-(xw2[q]-pw[3])^2/2/pw[4]^2)
end
When you post Igor code, it needs to be enclosed in the igor tags: http://www.igorexchange.com/node/25
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
April 21, 2015 at 09:15 am - Permalink