Programming for the evolution of spectrum with time

johnweeks wrote:

Perhaps what you want is MPF2_AutoMPFit(), which is part of the Multipeak Fit 2 package. It is somewhat difficult to use... If your spectra are in one wave per time step, you can feed all the waves into that function and it will do the whole batch. Creating the summary of all the data can be difficult, though.

In the Multipeak Fit 2 package control panel, click the Help button. Follow the link near the top of the help to the section on MPF2_AutoMPFit().

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

jjweimer wrote:

As I understand, you want to fit each spectrum, obtain the fitting parameters and their uncertainties, and analyze how those parameters change with time. Also, I understand that you have no programming experience (certainly not in Igor Pro) but that you have a strong interest to do the programming, are willing to devote the time needed, and have an urgent but not immediate need to get the results. In this case, I might recommend this approach to start.

* Create a function that you can call from the command line to fit just ONE of your spectra. For example, FitMySpectrum(...). The (...) input could be nothing or it could be the name of the spectrum wave. The return from this function would be a print of all of the fit parameters and their uncertainty values. You can used the Multipeak package but ... that is overkill. When you really want to learn how to automate this, you should start by learning how to write your own fitting function.
* Run FitMySpectrum(...) manually to fit each spectrum.
* Record the results in a spreadsheet-like manner either in Igor Pro or in a spreadsheet app in a way to import to Igor Pro. The results would be a set of columns time, K1, dK1, K2, dK2, ... where Kj is the fit constant and dKj is its uncertainty.
* Plot Kj versus time with error bars dKj in Igor Pro (or in your spreadsheet app).

Once you have convinced yourself that you have this manual process working, the next level of programming can begin. For this, you will be wise to make one fundamental change.

* Store your waves in Folders with names that do NOT have () and do NOT start with numbers. For example, instead of 1ps(2) the folder might be named t1ps_2. This naming convention avoids later issues about using "liberal" names (look up liberal names in the Igor help or manual).

After this, you can write a function to move through your folders starting at the root level and run FitMySpectrum(...) on the wave set(s) inside. Also, you can have the fitting results stored automatically in Igor Pro waves so that you do not have to print them, copy them over, and re-load them back in to Igor Pro by hand.

So, to recap ...

* Learn how to write the function required to fit ONE spectrum and use it to fit each spectrum in your data set.
--> print the results to the command line and graph or analyze them in Igor Pro or elsewhere
* Learn how to write a function that moves through different data folders.
* Learn how to store results inside Igor Pro.

Your best starting point when you have no foundation in coding but really do want to do this on your own is to stay firmly focused on the first task.

When you must have the results by yesterday ... well ... this becomes an entirely different problem.

Does this help focus on what to do next?

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH

Tika wrote:

Thank you Prof. Weimer for clear picture on how to proceed on. Since I am a beginner in programming, it took some time for me to go through the manual. I wrote a function 'FitMySpectrum' which does not show any error while compiling. I tried to call the function in the command line using 'Print FitMySpectrum()' but it returned NAN. In fact I do not know on how to use the function to fit the data 'col_40ps' because it does not pop up any box where I can enter the fitting parameters as in MultipeakFitting. Can you please suggest me further? Here is the function I wrote.
Thank you

#pragma rtGlobals=3     // Use modern global access method and strict wave access.
    Function FitMySpectrum()   
          Wave col_40ps
          Variable p1,p2,p3,p4 //fitting parameters
          Variable x //time
          //Variable p1 = xo           //Xoffset
         //Variable p2 = A            //amplitude
        //Variable p3 = s            // sigma parameter that determines the shape
        //Variable p4 = T           //Tau exponential decay constant
        Variable w1 = exp(0.5*(p3/p4)^2-(x-p1)/p4)
        Variable w2 = erfc(0.707*(p3/p4-(x-p1)/p3))
        Variable w3 = 0.5*p2/p4
        Variable w4 = w1*w2*w3     //exponentially modified gaussian function
   
    return w4
end

jtigor wrote:

Tika wrote:

Thank you Prof. Weimer for clear picture on how to proceed on. Since I am a beginner in programming, it took some time for me to go through the manual. I wrote a function 'FitMySpectrum' which does not show any error while compiling. I tried to call the function in the command line using 'Print FitMySpectrum()' but it returned NAN. In fact I do not know on how to use the function to fit the data 'col_40ps' because it does not pop up any box where I can enter the fitting parameters as in MultipeakFitting. Can you please suggest me further? Here is the function I wrote.
Thank you

#pragma rtGlobals=3     // Use modern global access method and strict wave access.
    Function FitMySpectrum()   
          Wave col_40ps
          Variable p1,p2,p3,p4 //fitting parameters
          Variable x //time
          //Variable p1 = xo           //Xoffset
         //Variable p2 = A            //amplitude
        //Variable p3 = s            // sigma parameter that determines the shape
        //Variable p4 = T           //Tau exponential decay constant
        Variable w1 = exp(0.5*(p3/p4)^2-(x-p1)/p4)
        Variable w2 = erfc(0.707*(p3/p4-(x-p1)/p3))
        Variable w3 = 0.5*p2/p4
        Variable w4 = w1*w2*w3     //exponentially modified gaussian function
   
    return w4
end

With the quick disclaimer that I have never used a custom fitting function I offer the following advice... a quick perusal of the manual (see "User-Defined Fitting Function" in the manual/help files) shows the function name should have " : FitFunc" appended in order for it to appear in the list of fit functions in the Curve Fitting dialog, as in:

        Function FitMySpectrum() : FitFunc

Additionally, your custom fit function must be set up to accept certain input parameters. It's best to read over the manual sections on user defined fitting functions for more precise info on this aspect.

Hope this helps move your project along!

Jeff

jjweimer wrote:

For two peaks, use the modification below. You must supply additional parameters in the coefficient wave to handle two peaks. Note: p, q, and x are reserved names in Igor Pro!

Function FitTwoPeaks(ww, xx) : FitFunc 
    WAVE ww     //fitting parameters
        variable xx

    //Variable ww[0] = xo           //Xoffset
    //Variable ww[1] = Ao            //amplitude
    //Variable ww[2] = so            // sigma parameter that determines the shape
    //Variable ww[3] = To           //exponential decay constant
    //Variable ww[4] = x1           //Xoffset
    //Variable ww[5] = A1            //amplitude
    //Variable ww[6] = s1            // sigma parameter that determines the shape
    //Variable ww[7] = T1           //exponential decay constant

    Variable xx                                      
    Variable w1=0, w2=0, w3=0, qq, ic                  // local variable
 
    // Body code
 
        for(ic=0;ic<2;ic+=1)
        w1 += exp(0.5*(ww[2+ic*4]/ww[3+ic*4])^2-(xx - ww[0+ic*4])/ww[3+ic*4])  
        w2 += erfc(0.707*(ww[2+ic*4]/ww[3+ic*4]-(xx - ww[0+ic*4])/ww[2+ic*4]))
        w3 += 0.5*ww[1+ic*4]/ww[3+ic*4]
        endfor
    qq = (w3*w1)*w2

    return qq
end

Assume your data is in a scaled wave called mySpectrum. Move to the folder that contains that wave. Here is a summary of the input to run this from the command line.

make/N=8/D/O wcoeff // this is the coefficient wave
wcoeff = {...} // provide initial guesses for coefficients ... read the manual to see how to do this
FuncFit FitTwoPeaks, wcoeff, mySpectrum /D

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH

johnweeks wrote:

In the experiment file you posted, I see two Multipeak Fit sets, and both give no problems when you click the Fit button. Please tell us how to reproduce the problem you are having with Multipeak Fit. If the EMG peak shape works for you, Multipeak Fit is probably the easiest way for you to proceed.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

johnweeks wrote:

My favorite technique for debugging such things is to put a test into the fit function and then break into the debugger when the test detects a bad number. So modify your function like this:

Function FitTwoPeaks(ww, xx) : FitFunc 
    WAVE ww     //fitting parameters
        variable xx
 
    //Variable ww[0] = xo           //Xoffset
    //Variable ww[1] = Ao            //amplitude
    //Variable ww[2] = so            // sigma parameter that determines the width
    //Variable ww[3] = To           //exponential decay constant
    //Variable ww[4] = x1           //Xoffset
    //Variable ww[5] = A1            //amplitude
    //Variable ww[6] = s1            // sigma parameter that determines the width
    //Variable ww[7] = T1           //exponential decay constant
 
    //Variable xx                                      
    Variable w1=0, w2=0, w3=0, qq, ic                  // local variable
 
    // Body code
 
        for(ic=0;ic<2;ic+=1)
        w1 += exp(0.5*(ww[2+ic*4]/ww[3+ic*4])^2-(xx - ww[0+ic*4])/ww[3+ic*4])  
        w2 += erfc(0.707*(ww[2+ic*4]/ww[3+ic*4]-(xx - ww[0+ic*4])/ww[2+ic*4]))
        w3 += 0.5*ww[1+ic*4]/ww[3+ic*4]
        endfor
    qq = (w3*w1)*w2
    if (numtype(qq))
        print "Oops"
    endif
 
    return qq
end

Now right-click in a procedure window and select Enable Debugger. Put a breakpoint on Print "Oops" and do the fit. When coefficients result in a return value that is NaN or Inf, you will break into the debugger. Now you can examine the coefficients and see what's happening.

When I did this, w1 was Inf. I didn't look at it in detail; my guess is that some coefficient got too big and the exp() function overflowed.

You will want to suppress the Curve Fit Progress window, as it may cause a lock-up condition with the debugger.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com