Thanks Antoine, it does work after I restart Scilab.
However, my actual code consists of more than 1 data set, and I would have
to normalise them according the their respective areas under the curve. The
code only works without the normalisation part (replacing ri with r), but
error 21 arises when the part is included.
funcprot(0)
workbook = readxls('C:\Users\Kying\Desktop\Stray_light\Data.xls')
sheet = workbook(1) //old
k=6.63e-34*3e8/1.38e-23
x=sheet(:,2)
ri=sheet(:,3)
gi=sheet(:,4)
bi=sheet(:,5)
uvi=sheet(:,6)
//normalization
norm=max([inttrap(x,ri),inttrap(x,gi),inttrap(x,bi),inttrap(x,uvi)])
r=norm/inttrap(x,ri)*ri
g=norm/inttrap(x,gi)*gi
b=norm/inttrap(x,bi)*bi
uv=norm/inttrap(x,uvi)*uvi
//end of normalization
dims=1
ur=min(find(r>0.71*max(r))) //uppler limit
lr=max(find(find(r<10)<find(r==max(r)))) //lower limit
wr=matrix((cat(dims,zeros(lr-1,1),ones(ur-lr+1,1),zeros(1024-ur,1))),1024,1)
//weighting
//exponential function
function y=yth(x,a)
y=a(1)*exp(-k./x/a(2))
endfunction
//initial parameters
a0=[1.0;300]
function e=myfun(a,x,r,wr)
e=wr.*(yth(x,a)-r)
endfunction
//call leastsq
[f,xopt,gopt]=leastsq(list(myfun,x,r,wr),a0)
scf();
plot2d(x,r,logflag="nl",style=[color("red")])
//plot2d(x,yth(x,xopt))
plot2d(x,g,logflag="nl",style=[color("green")])
plot2d(x,b,logflag="nl",style=[color("blue")])
plot2d(x,uv,logflag="nl",style=[color("purple")])
Error message:
-->exec('C:\Users\Kying\Desktop\Stray_light\testing2.sce', -1)
!--error 21
Invalid index.
at line 118 of function numderivative called by :
at line 2 of function Dfun called by :
at line 2 of function %opt called by :
at line 92 of function leastsq called by :
[f,xopt,gopt]=leastsq(list(myfun,x,r,wr),a0)
at line 38 of exec file called by :
exec('C:\Users\Kying\Desktop\Stray_light\testing2.sce', -1)
Do you have any idea why that is?
Thanks
--
View this message in context:
http://mailinglists.scilab.org/Scilab-leastsq-exponential-fitting-tp4030949p4030973.html
Sent from the Scilab users - Mailing Lists Archives mailing list archive at
Nabble.com.
_______________________________________________
users mailing list
[email protected]
http://lists.scilab.org/mailman/listinfo/users
