%calka = 0.03562*x^4 + 0.8672*x^3 + 0.572*x^2 + 7.24*x + 3.21

k = 1;
for i=-24:0.5:7
   x(k) = i;
   y(k) = f(i);
   k = k+1;
end
hold off;
plot(x, y);
hold on;

x0 = -24;
xn = 7;

wynik_ana = cal_nieozn(xn) - cal_nieozn(x0)

% Romberg
h = (xn-x0) * 2;
for i=1:1:4
    h = h / 2
    I1 = trapez(x0,xn,h);
    I2 = trapez(x0,xn,h/2);
    I = 4/3 * I2 - 1/3 * I1;
    eWzgl = (wynik_ana - I)/wynik_ana * 100;
end
eWzgl
I
% Gauss

h1 = (xn + x0)/2;
h2 = (xn - x0)/2;

GI1 = f(h1 + h2 * (-sqrt(3/5))) * h2;

GI0 = f(h1) * h2;

GI2 = f(h1 + h2 * sqrt(3/5)) * h2;
GI = 5/9 * GI1 + 8/9 * GI0 + 5/9 * GI2
eWzgl = (wynik_ana - GI)/wynik_ana * 100