# How to find the integrand by fitting data - spectrum reconstruction

## Problem statement

I have some experimental data of measurement values vs. temperature. Now, I want to find the spectrum which causes the measured temperature dependence.

The measurement values are the integrals of ## The forward calculation

I will now illustrate the problem with an example. There is an yet unknown sensitivity spectrum.

spectrum[lam_] := Exp[-((lam - 5)/2)^4];
Plot[spectrum[lam], {lam, 0, 10}, AxesLabel -> {"wavelength", "sensitivity"}] And there is a known spectral intensity

intensity[lam_, T_] := 1/(T^4 lam^5) 1/(Exp[1/(lam *T)] - 1)
Plot[Table[intensity[lam, T], {T, 0.03, 0.1, 0.01}] // Evaluate, {lam, 0, 10},
AxesLabel -> {"wavelength", "intensity"}] // Quiet The measured signals are the areas under these curves

Plot[Table[spectrum[lam]*intensity[lam, T], {T, 0.03, 0.1, 0.01}] // Evaluate, {lam, 0, 10},
AxesLabel -> {"wavelength", "intensity"}] // Quiet

data = Table[{T,NIntegrate[spectrum[lam]*intensity[lam, T], {lam, 1, 10}]}, {T, 0.03, 0.1, 0.01}];
ListPlot[data, AxesLabel -> {"temperature", "signal"}]


## What's the solution?

I thought of nonlinear curve fitting for the control points of an interpolation function as was proposed in this question by Hugh and Alexey Popkov. Using the sample data, the result should look like the spectrum[lam] from above.

But how to make this code work?

nOfControlPoints = 5;
controlPoints = Subdivide[1, 10, nOfControlPoints - 1]

model[y : {__Real}] := Interpolation[Transpose[{controlPoints, y}]]
fctn = NIntegrate[model[Array[y, nOfControlPoints]][lam]*intensity[lam, T], {lam, 1, 10}];

nlm = NonlinearModelFit[data, fctn, Array[y, nOfControlPoints], T]


Fixing obvious mistakes and slightly redesinging the code in order to avoid NIntegrate error messages, we get:

nOfControlPoints = 5;
controlPoints = Subdivide[1, 10, nOfControlPoints - 1];

Clear[fctn, model];
model[y__Real] := Interpolation[Transpose[{controlPoints, {y}}]]
fctn[{y__Real}][T_Real] :=
NIntegrate[model[y][lam]*intensity[lam, T], {lam, 1, 10},
Method -> {Automatic, "SymbolicProcessing" -> 0}];

nlm = NonlinearModelFit[data, fctn[Array[y, nOfControlPoints]][T],
Array[y, nOfControlPoints], T] Now we can plot the model along with the data:

Show[ListPlot[data], Plot[nlm["Function"][T], {T, 0.03, 0.1}], Frame -> True] We see that NonlinearModelFit has done its job well.

But the fitted model doesn't reproduce the original spectrum well:

Show[Plot[spectrum[lam], {lam, 0, 10},
AxesLabel -> {"wavelength", "sensitivity"}],
Plot[Evaluate[model @@ Values@nlm["BestFitParameters"]][lam], {lam, 1, 10},
PlotStyle -> Directive[{Dashed, Red}]], PlotRange -> All] So you can try to increase the number of sample points, try other interpolation method, try to add constraints, better starting values, etc…