MATLAB: Nlinfit with modelfun as an integral2 and variable integral limits

fminsearchintegral2nlinfitpass over integral limits

I have several volumes over specific areas and want to know which 3D bell shape fits best for those partial volumes. I try to estimate the peak value and sigma of the 3D bell shape by a double integral (‘integral2’) within ‘nlinfit’. I have seen a solution for a single integral ("nlinfit with modelfun as an integral"). and tried to modify it.
But I do not know how to pass over the integral limits.
%Given e.g. 3 volumes
%Unknown: peak and sigma for bell shape that fits best those volume parts
Vin=[0.1503 0.0949 0.0238];%Volume at xin,yin
xin=[0 1 2];%x-position of volume part under the 3D bell shape
yin=[0 0 0];%y-position of volume part under the 3D bell shape
parIn=[0.1632 1];%initial guess for peak and sigma of bell shape
%In this case it is also the result.
%Formula for bell volume, double integral. par(1)=peak, par(2)=sigma:
Integrand = @(par,x,y) par(1)*exp(-(x.^2+y.^2)/(2*par(2)^2));
Integralfun=@(par,x,y,xin,yin) integral2(@(x,y) Integrand(par,x,y),xin-0.5,xin+0.5,yin-0.5,yin+0.5);
parOut= nlinfit(xin,Vin,Integralfun,parIn);

Best Answer

I could not figure out how to solve the problem with 'nIinfit' , therefore, I solved the problem by using fminsearch’.
%Given e.g. 3 volumes
%Unknown: peak and sigma for bell shape that fits best those volume parts
Vin=[0.1503 0.0949 0.0238];%Volume at xin,yin
xin=[0 1 2];%x-position of volume part under the 3D bell shape
yin=[0 0 0];%y-position of volume part under the 3D bell shape
% parIn=[0.1632 1];%initial guess for peak and sigma of bell shape
parIn=[1 2];%initial guess for peak and sigma of bell shape
%Formula for bell volume, double integral. par(1)=peak, par(2)=sigma:
Integrand = @(par,x,y) exp(-(x.^2+y.^2)/(2*par(2)^2));
%Build the function for the sum of absolute differences:
funStr='@(par) ';
for i=1:numel(Vin)
is=num2str(i);%i as string
funStr=[funStr 'abs(par(1)*integral2(@(x,y) Integrand(par,x,y),'...
'xin(' is ')-0.5,xin(' is ')+0.5,yin(' is ')-0.5,yin(' is ')+0.5)-Vin(' is '))+'];
end
funStr=funStr(1:end-1);%remove last '+'
mystr2func=@(Str)evalin('base',Str);%dirty trick to include workspace variables
fun=mystr2func(funStr);%convert string to function
[parOut,fval,exitflag,output] = fminsearch(fun,parIn)
%parOut expected: [0.1632 1]