% Usage: r = boole(f,a,b,n)
% Composite Boole's (aka Bode's) rule (5-point closed Newton-Cotes)
%
% Input:
% f - Matlab inline function 
% a,b - integration interval
% n - number of subintervals (panels)
%
% Output:
% r - computed value of the integral
%
% Examples:
% r=boole(@sin,0,1,10);
% r=boole(@myfunc,0,1,10);          here 'myfunc' is any user-defined function in M-file
% r=boole(inline('sin(x)'),0,1,10);
% r=boole(inline('sin(x)-cos(x)'),0,1,10);

function r = boole(f,a,b,n)
h = (b - a) / (n * 4);
r = 7 * f(a);
x = a + h;
for i = 1 : n-1
    r = r + 32 * f(x);
    x = x + h;
    r = r + 12 * f(x);
    x = x + h;
    r = r + 32 * f(x);
    x = x + h;
    r = r + 14 * f(x);
    x = x + h;
end;
r = r + 32 * f(x);
x = x + h;
r = r + 12 * f(x);
x = x + h;
r = r + 32 * f(x);
r = r + 7 * f(b);
r = r * h*2/45;