Message Boards Message Boards

0
|
3219 Views
|
8 Replies
|
6 Total Likes
View groups...
Share
Share this post:

Fitting peaks with a Voigt function

Dear All,

I am trying to fit a peak by using a Voigt function. If I define the Voigt function with NIntegrate like this:

Voigt[x_, Area_, wL_, wG_, x0_, y0_] := 
y0 + Area (2 Log[2])/
Pi^1.5 wL/(wG)^2 NIntegrate[
Exp[-t^2]/((Sqrt[Log[2]] wL/wG)^2 + (Sqrt[4 Log[2]] (x - x0)/wG - 
t)^2) , {t, -\[Infinity], \[Infinity]}]

I get some error messages during the fit procedure:

In[12]:= P1600 = 
 FindFit[B1600, {Voigt[x, Area, wL, wG, x0, y0], 1 < wL < 10, 
   1 < wG < 10, 1580 < x0 < 1620, 1000 < y0 < 3000}, {Area, wL, wG, 
   x0, y0}, x]

During evaluation of In[12]:= NIntegrate::inumr: The integrand E^-t^2/((-t+(2 (x+Times[<<2>>]) Sqrt[Log[<<1>>]])/wG)^2+(wL^2 Log[2])/wG^2) has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

During evaluation of In[12]:= FindFit::fitd: First argument B1600 in FindFit is not a list or a rectangular array.

During evaluation of In[12]:= NIntegrate::inumr: The integrand E^-t^2/((-t+(2 (x+Times[<<2>>]) Sqrt[Log[<<1>>]])/wG)^2+(wL^2 Log[2])/wG^2) has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

Out[12]= FindFit[B1600, {y0 + (
   0.248961 Area wL NIntegrate[
     Exp[-t^2]/(((Sqrt[Log[2]] wL)/
       wG)^2 + ((Sqrt[4 Log[2]] (x - x0))/wG - 
        t)^2), {t, -\[Infinity], \[Infinity]}])/wG^2, 1 < wL < 10, 
  1 < wG < 10, 1580 < x0 < 1620, 1000 < y0 < 3000}, {Area, wL, wG, x0,
   y0}, x]

If I try to replace "NIntegrate" with "Integrate" in the definition of the Voigt function, I do not get any error message in the fitting, but the fitting never ends!

Has anyone had a similar problem or have suggestions?

PS: with a Lorenzin function the script works.

all the best,

roberto
POSTED BY: Robertino Pilot
8 Replies
Sort By:

Robertino,

Mathematica already does know the Voigt function (i.e. VoigtDistribution)! So - it is probably better to just use this function. Let data be your data, then:

ClearAll[ampl, xShift, yShift, \[Delta], \[Sigma]];
model = {ampl Re[PDF[VoigtDistribution[\[Delta], \[Sigma]], x - xShift]] + yShift};
params = FindFit[data, model, {{ampl, 100000}, {xShift, 1600}, yShift, \[Delta], \[Sigma]}, x]

enter image description here

FindFit seems to have some problems, but the result looks quite reasonable:

Show[Plot @@ {model /. params, {x, 1550, 1650}, PlotStyle -> {Thick, Black}, PlotRange -> {0, All}}, 
 ListLinePlot[data, PlotRange -> All, PlotStyle -> Red], 
 GridLines -> Automatic, ImageSize -> Large]

enter image description here

Does that help? Regards -- Henrik

As a side note: In radiation physics (my field of work) the Voigt function can be of importance, because the point spread function (PSF) from incoming fluence to dose is assumed to be a Lorentzian and the PSF from dose to measured dose (using an ionization chamber) is approximated by a Gaussian. So one needs the convolution of both functions - which is the Voigt function!
POSTED BY: Henrik Schachner
Jim Baldwin
Jim Baldwin
Posted 1 year ago

Henrick Schachner's code gets you estimates of the parameters but you might want to look at the residuals to help assess whether the fit is adequate. If assessing the fit is desired other than by "I'll know it when I see it" from a plot of the data, then using NonlinearModelFit is recommended.

nlm = NonlinearModelFit[data, model, {{ampl, 183259}, {xShift, 1601.4881455559296`}, {yShift, 
    1451.3252395718675`}, {\[Delta], 2.822192930059779`}, {\[Sigma], 1.3677109564952756`}}, 
    x, MaxIterations -> 1000];

This results in a warning:

NonlinearModel warning

Plotting the residuals by the predictor variable results in the following:

Residuals vs predictor variable

This does not look like a good fit to me especially if the fit for the values in the middle of the dataset are important. This suggests that a regression with a translated shape of a Voigt distribution is not an adequate fit to the data.
POSTED BY: Jim Baldwin
Jim Baldwin
Jim Baldwin
Posted 1 year ago

Without the data it's pretty hard to help.

Note that B1600 and Voigt are in blue text which means they haven't been defined yet.
POSTED BY: Jim Baldwin
Jim Baldwin
Jim Baldwin
Posted 1 year ago

Henrik: Thanks! That clears it up. We are in agreement about definition of the Voigt distribution as to how it comes about as a convolution of a Gaussian and Lorentzian.

It's that from my experience in this forum and the Mathematica.StackExchange forum, that some physicists seem to think that in a regression situation (i.e. fitting some curve to measured data) if the curve form is either proportional to a known probability function (or a translation as above), then somehow the probabilistic properties of the associated probability distribution are somehow involved when those properties are not involved. And the fitting is labeled as "fitting a distribution" when it is more accurate to state that one is fitting a curve with a similar form as a known probability distribution.

It even seems the desire to perform a regression occurs even when one has random samples from a probability distribution: the data are binned (which loses information) and a regression is performed to estimate the parameters rather than using standard statistical methods such as maximum likelihood or the method of moments.

While the above question doesn't describe the process that generates the data, having the differences in the predictor variables always equal 0.87 or 0.88 (with an occasional 0.86) is consistent with binning and leaves me suspicious. (Although maybe this is just due to sampling at approximately equal time intervals.)
POSTED BY: Jim Baldwin

Dear Jim,

thanks a lot for your reply! I have copied the script in a more complete form here and in attachment.

thanks in advance anyone who can help :)

all the best

roberto pilot

(DATA TO BE FITTED)

In[30]:= 
B1600 = {{1539.64`, 1383.06`}, {1540.52`, 1368.56`}, {1541.4`, 
    1372.31`}, {1542.29`, 1383.06`}, {1543.17`, 1361.06`}, {1544.05`, 
    1379.06`}, {1544.93`, 1396.06`}, {1545.81`, 1389.06`}, {1546.69`, 
    1391.81`}, {1547.58`, 1405.56`}, {1548.46`, 1390.81`}, {1549.34`, 
    1407.31`}, {1550.22`, 1385.06`}, {1551.1`, 1425.56`}, {1551.98`, 
    1402.81`}, {1552.86`, 1391.31`}, {1553.74`, 1414.56`}, {1554.62`, 
    1429.56`}, {1555.5`, 1433.06`}, {1556.38`, 1442.31`}, {1557.27`, 
    1458.31`}, {1558.15`, 1444.81`}, {1559.03`, 1463.31`}, {1559.9`, 
    1468.81`}, {1560.78`, 1455.56`}, {1561.66`, 1502.81`}, {1562.54`, 
    1508.31`}, {1563.42`, 1472.81`}, {1564.3`, 1516.56`}, {1565.18`, 
    1529.56`}, {1566.06`, 1541.56`}, {1566.94`, 1556.31`}, {1567.82`, 
    1604.06`}, {1568.7`, 1632.56`}, {1569.58`, 1655.06`}, {1570.46`, 
    1670.56`}, {1571.33`, 1748.31`}, {1572.21`, 1776.31`}, {1573.09`, 
    1776.31`}, {1573.97`, 1817.81`}, {1574.85`, 1854.06`}, {1575.72`, 
    1842.31`}, {1576.6`, 1865.31`}, {1577.48`, 1927.81`}, {1578.36`, 
    1937.06`}, {1579.23`, 1989.56`}, {1580.11`, 2063.06`}, {1580.99`, 
    2124.81`}, {1581.87`, 2140.31`}, {1582.74`, 2258.06`}, {1583.62`, 
    2274.31`}, {1584.5`, 2382.31`}, {1585.37`, 2397.81`}, {1586.25`, 
    2417.06`}, {1587.13`, 2429.06`}, {1588, 2519.06`}, {1588.88`, 
    2585.06`}, {1589.76`, 2805.56`}, {1590.63`, 3004.56`}, {1591.51`, 
    3252.31`}, {1592.38`, 3607.56`}, {1593.26`, 3949.06`}, {1594.14`, 
    4439.06`}, {1595.01`, 5137.31`}, {1595.89`, 6125.06`}, {1596.76`, 
    7607.06`}, {1597.64`, 9685.31`}, {1598.51`, 12128.3`}, {1599.39`, 
    14735.3`}, {1600.26`, 17079.6`}, {1601.14`, 18656.6`}, {1602.01`, 
    18742.3`}, {1602.89`, 17297.8`}, {1603.76`, 14656.6`}, {1604.63`, 
    11683.1`}, {1605.51`, 9005.31`}, {1606.38`, 6929.56`}, {1607.26`, 
    5405.31`}, {1608.13`, 4430.31`}, {1609, 3766.06`}, {1609.88`, 
    3351.56`}, {1610.75`, 3082.06`}, {1611.63`, 2842.56`}, {1612.5`, 
    2694.06`}, {1613.37`, 2566.81`}, {1614.24`, 2456.81`}, {1615.12`, 
    2418.56`}, {1615.99`, 2290.31`}, {1616.86`, 2243.06`}, {1617.73`, 
    2195.56`}, {1618.61`, 2135.56`}, {1619.48`, 2090.31`}, {1620.35`, 
    2027.56`}, {1621.22`, 2006.56`}, {1622.1`, 1939.56`}, {1622.97`, 
    1904.06`}, {1623.84`, 1862.31`}, {1624.71`, 1789.31`}, {1625.58`, 
    1748.06`}, {1626.46`, 1720.31`}, {1627.33`, 1695.06`}, {1628.2`, 
    1671.06`}, {1629.07`, 1642.56`}, {1629.94`, 1640.56`}, {1630.81`, 
    1615.06`}, {1631.68`, 1600.81`}, {1632.55`, 1580.31`}, {1633.42`, 
    1538.56`}, {1634.3`, 1576.56`}, {1635.17`, 1599.06`}, {1636.04`, 
    1556.31`}, {1636.91`, 1565.06`}, {1637.78`, 1547.06`}, {1638.65`, 
    1551.06`}, {1639.52`, 1575.31`}, {1640.39`, 1555.31`}, {1641.26`, 
    1524.81`}, {1642.13`, 1552.56`}, {1643, 1540.31`}, {1643.87`, 
    1541.06`}, {1644.73`, 1514.06`}, {1645.6`, 1515.81`}, {1646.47`, 
    1499.56`}, {1647.34`, 1503.06`}, {1648.21`, 1492.06`}, {1649.08`, 
    1444.31`}, {1649.95`, 1464.81`}, {1650.82`, 1459.31`}, {1651.69`, 
    1447.56`}, {1652.55`, 1413.31`}, {1653.42`, 1420.56`}, {1654.29`, 
    1427.81`}, {1655.16`, 1409.56`}, {1656.02`, 1391.56`}, {1656.89`, 
    1384.31`}, {1657.76`, 1403.31`}, {1658.63`, 1397.06`}, {1659.5`, 
    1387.06`}, {1660.36`, 1372.56`}, {1661.23`, 1373.81`}, {1662.1`, 
    1380.81`}, {1662.96`, 1373.81`}, {1663.83`, 1376.81`}, {1664.7`, 
    1366.56`}, {1665.57`, 1367.06`}, {1666.43`, 1368.56`}, {1667.3`, 
    1362.81`}, {1668.16`, 1363.06`}, {1669.03`, 1385.56`}, {1669.9`, 
    1363.06`}, {1670.76`, 1351.06`}};

In[18]:= (*DEFINITION OF FIT FUNCTIONS;
VOIGT1 AND VOIGT2 DIFFER FOR THE CALCULATION OF THE INTEGRAL: \
NINTEGRATE VS INTEGRATE*)

In[31]:= 
Lorentz[x_, Area_, \[CapitalGamma]_, x0_, y0_] := 
 y0 + (2 Area)/Pi \[CapitalGamma]/(4 (x - x0)^2 + \[CapitalGamma]^2)

In[32]:= 
Voigt1[x_, Area_, wL_, wG_, x0_, y0_] := 
 y0 + Area (2 Log[2])/
   Pi^1.5 wL/(wG)^2 Integrate[
    Exp[-t^2]/((Sqrt[Log[2]] wL/wG)^2 + (Sqrt[4 Log[2]] (x - x0)/wG - 
        t)^2) , {t, -\[Infinity], \[Infinity]}]

In[33]:= 
Voigt2[x_, Area_, wL_, wG_, x0_, y0_] := 
 y0 + Area (2 Log[2])/
   Pi^1.5 wL/(wG)^2 NIntegrate[
    Exp[-t^2]/((Sqrt[Log[2]] wL/wG)^2 + (Sqrt[4 Log[2]] (x - x0)/wG - 
        t)^2) , {t, -\[Infinity], \[Infinity]}]

In[34]:= (*FIT WITH A VOIGT PROFILE*)

In[23]:= FitVoigt1 = 
 FindFit[B1600, {Voigt1[x, Area, wL, wG, x0, y0], 1 < wL < 10, 
   1 < wL < 10, 1580 < x0 < 1620, 1000 < y0 < 3000}, {Area, wL, wG, 
   x0, y0}, x]

Out[23]= $Aborted

In[24]:= FitVoigt2 = 
 FindFit[B1600, {Voigt2[x, Area, wL, wG, x0, y0], 1 < wL < 10, 
   1 < wL < 10, 1580 < x0 < 1620, 1000 < y0 < 3000}, {Area, wL, wG, 
   x0, y0}, x]

During evaluation of In[24]:= NIntegrate::inumr: The integrand E^-t^2/((-t+(2 (x+Times[<<2>>]) Sqrt[Log[<<1>>]])/wG)^2+(wL^2 Log[2])/wG^2) has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

During evaluation of In[24]:= NIntegrate::inumr: The integrand E^-t^2/((-t+(2 (<<1>>) Sqrt[Log[<<1>>]])/wG)^2+(wL^2 Log[2])/wG^2) has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

During evaluation of In[24]:= NIntegrate::inumr: The integrand E^-t^2/((-t+(2 (<<1>>) Sqrt[Log[<<1>>]])/wG)^2+(wL^2 Log[2])/wG^2) has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

During evaluation of In[24]:= General::stop: Further output of NIntegrate::inumr will be suppressed during this calculation.

Out[24]= $Aborted

In[58]:= (*FIT WITH LORENTZ*)

In[35]:= FitLorentz = 
 FindFit[B1600, {Lorentz[x, Area, \[CapitalGamma], x0, y0], 
   1 < \[CapitalGamma] < 10, 1580 < x0 < 1620, 
   1000 < y0 < 3000}, {Area, \[CapitalGamma], x0, y0}, x]

Out[35]= {Area -> 193881., \[CapitalGamma] -> 6.83586, x0 -> 1601.5, 
 y0 -> 1371.79}
Attachments:
Prove 2.nb
POSTED BY: Robertino Pilot
Jim Baldwin
Jim Baldwin
Posted 1 year ago

Henrick: I have a question related to your last comment: How do you justify fitting a regression with something that might happen to have a translated (shifted and multiplied) shape of a Voigt distribution pdf function and random samples from a Voigt distribution. Those are two completely different beasts. There are no random samples from a Voight distribution in this question so I'm not seeing how the convolution of a normal and Lorentzian distribution is relevant here. What am I missing?
POSTED BY: Jim Baldwin

Dear @Jim Baldwin,

thank you for the reference to NonlinearModelFit[]! I admit that I do not have any experience with that function so far (while the warning messages are the same as from simple FindFit[]). I just wanted to answer the OPs question ("How to fit the data with Voigt?" and not "Does Voigt give a good fit?"). But of course you are right: "I'll know it when I see it" is not a good criterion! So, again, thank you!

Maybe my last remark was a bit unclear. It has nothing to do with the original question. The convolution of Gaussian and Lorentzian is simply the definition of Voigt. I regard the Voigt function as quite an exotic thing - so I just wanted to mention how I learned about it, hoping that I hear about its context here.

Regards -- Henrik
POSTED BY: Henrik Schachner

Hi Jim and Henrik!,

thank you very much for contrinbuting to my question :).

The data are part of a Raman spectrum therefore in the x-axis there is Raman shift and in y-axis there is counts. As far as I know, a Voigt profile is the most general way to fit Raman peaks since it accounts for the natural linewidth (Lorentzian shape) and for the inhomogeneous broadening (convolution with a Gaussian).

The fitting itself may not be perfect because there could be more peaks nearby, as it can be seen form the residuals. However, what I find a bit weird is the error message in performing the fitting: If I fit the same data with a Voigt profile with OriginLab (Microcal) I do not get any error and the values are a bit different. Also, the peak is a very tall and I would not expect many difficulties in catching it. Maybe some settings in the fitting procedure I guess...I am not expert in the technical procedures.

Very helpful also the function NonlinearModelFit and its additional features compared to FindFit :)

wish you all a nice day!
POSTED BY: Robertino Pilot
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract