clear all;
d=load('alu.dat');
r=d(:,1);
n=d(:,2)*10;
d=r;
dd=0.01;

%korrektur:
k=[1:numel(n)];
gt(k)=1./(10 - n(k).*(1.6.*10.^(-004)));
n=n.*gt';

%standartabweichung
dn=sqrt(n);


%Fit
%   f= y0 + A exp(-µ1 d) + B exp(-µ2 d)
%   f= y(1) + y(2) exp(-y(3) d) + y(4) exp(-y(5) d)
     a=[0,500,16,500,2];
 	eps = 0.0001;     %max relativer fehler
	
	ndata = numel(n);
	ma = numel(a);
	[a,sd_a,var,covar,c] = mrqmin(d,n,dn,ndata,ma,a,eps,'funcs_alu');
     disp('f= a(1) + a(2) exp(-a(3) d) + a(4) exp(-a(5) d)');
     disp(a);
     disp(sd_a);
    
        dd = linspace(0,max(d)*1.1);
	f= a(1) + a(2).* exp(-a(3) .* dd) + a(4).* exp(-a(5) .* dd);
  	
	chi2=var(1)*(numel(d)-ma);
    	disp('chi^2');
	disp(chi2);
    	disp('Varianz');
    	disp(var);
    	sigma_v=sqrt(2/(numel(d)-4));
    	disp('Standartabweichung der Varianz');
    	disp(sigma_v);




hold on;
	p3=plot(dd,f);
	p1=errorbar(d,n,dn,'.');
	ax1=gca;
	set(ax1,'YLim',[0,max(f)*1.1]);
	set(ax1,'XLim',[-1,max(dd)*1.1]);
	set(ax1,'YScale','log');
	
	le1=legend([p1(2),p3],'n','Doppelexponentialfit',1);
	xl1=xlabel('Aluminium- Dicke d');
	yl1=ylabel('n');

	grid on;
hold off;

