function f = midpoint(x_edge, y, dim)
	% calculate integral of y with respect to x
	% using the midpoint method over dimension dim
	% given x_edges and y
	
	if nargin < 3
		dim = 1;
	end
	
	% get dx and reshape it
	dx		= diff(x_edge);
	S		= ones(1,ndims(y));
	S(dim)= numel(dx);
	dx		= reshape(dx, S);
	% calculate integral
	dx_y	= bsxfun(@times, dx, y);
	f		= sum(dx_y, dim);
end