function I=midpoint(f,a,b,n)
h=(b-a)/n;
x=a+h/2;
I=0;
for i=1:n
    I=I+f(x);
    x=x+h;
end
I=I/n;