function [ Q ] = MyQuadrature( x,y )
% returns quadrature of points x, f(x)
% based on different Newton-Cotes rules

n = length(x);

L = 0;
R = 0;
M = 0;
T = 0;

for i = 1:n-1
    
    L = L + y(i) * (x(i+1) - x(i));
    
    R = R + y(i+1) * (x(i+1) - x(i));
    
end

% Midpoint rule
% First point has half interval

M = y(1) * 0.5 * (x(2)-x(1));

for i = 2:n-1
    M = M + y(i) * 0.5 * (x(i+1) - x(i-1));
end

M = M + y(n) * 0.5 * (x(n)-x(n-1));
    

% Trapezoidal rule 

for i = 1 : n - 1
    T = T + 0.5 * (y(i)+y(i+1)) * (x(i+1)-x(i));
end

end