%% Nactení dat
data = readmatrix('data.txt');  % ocekava 3 sloupce: x y sigma
x = data(:,1);
y = data(:,2);
err = data(:,3);
N = length(x);

%% 1) Linearni regrese (analyticky postup)
A = [x ones(N,1)];               % design matrix
W = diag(1./err.^2);             % vahy podle sigma
cov_lin = inv(A'*W*A);           % kovariantni matice
params_lin = cov_lin*(A'*W*y);   % parametry

m_lin = full(params_lin(1));
c_lin = full(params_lin(2));
dm_lin = sqrt(cov_lin(1,1));
dc_lin = sqrt(cov_lin(2,2));

res_lin = y - (m_lin*x + c_lin);
chi2_lin = sum((res_lin./err).^2);
chi2r_lin = chi2_lin / (N-2);

fprintf('Linearni regrese (analyticky)\n');
fprintf('m = %.6g ± %.2g\n', m_lin, dm_lin);
fprintf('c = %.6g ± %.2g\n', c_lin, dc_lin);
fprintf('chi^2_red = %.3f\n\n', chi2r_lin);

%% 2) Nelinearni regrese (numericky fit s vahami)
xw = x;                       % nezávislá promenna bez vah
yw = y ./ err;                 % vážená závislá promenna
modelw = @(p,xx) (p(1)*xx + p(2)) ./ err;  % vazeny model

p0 = [1 0];                    % pocatecni odhad parametru
opts = optimoptions('lsqcurvefit','Display','off');

[popt,resnorm,residual,~,~,~,J] = lsqcurvefit(modelw, p0, xw, yw, [], [], opts);

% Explicitni doplneni nulami tam, kde jsou prazdna mista. Aby fungoval fprintf
popt = full(popt);
J = full(J);

m_nl = popt(1);
c_nl = popt(2);
cov_nl = inv(J'*J);
dm_nl = sqrt(cov_nl(1,1));
dc_nl = sqrt(cov_nl(2,2));

res_nl = y - (m_nl*x + c_nl);
chi2_nl = sum((res_nl./err).^2);
chi2r_nl = chi2_nl / (N-2);

fprintf('Nelinearni regrese (numericky)\n');
fprintf('m = %.6g ± %.2g\n', m_nl, dm_nl);
fprintf('c = %.6g ± %.2g\n', c_nl, dc_nl);
fprintf('chi^2_red = %.3f\n', chi2r_nl);

%% 3) Vykresleni
xx = linspace(min(x),max(x),300);

% Figure 1: linearni regrese
figure;
errorbar(x,y,err,'o'); hold on;
plot(xx, m_lin*xx + c_lin,'r','LineWidth',2);
xlabel('x'); ylabel('y');
title('Linearni regrese (analyticky)');
legend('data','fit','Location','northwest');
text(0.05,0.75,sprintf('m = %.3g ± %.2g\nc = %.3g ± %.2g\nchi^2_{red}=%.3f',...
    m_lin,dm_lin,c_lin,dc_lin,chi2r_lin),...
    'Units','normalized','BackgroundColor','white');

% Figure 2: nelinearni regrese
figure;
errorbar(x,y,err,'o'); hold on;
plot(xx, m_nl*xx + c_nl,'g','LineWidth',2);
xlabel('x'); ylabel('y');
title('Nelinearni regrese (numericky)');
legend('data','fit','Location','northwest');
text(0.05,0.75,sprintf('m = %.3g ± %.2g\nc = %.3g ± %.2g\nchi^2_{red}=%.3f',...
    m_nl,dm_nl,c_nl,dc_nl,chi2r_nl),...
    'Units','normalized','BackgroundColor','white');
