mymoment(x, y, n)

function [mom] = mymoment(x, y, n)
%MYMOMENT: Computes the nth moment mean of a function (needs more
%implentation).
%
%	m = moment(x, y, n)
%
%		m	= nth moment of y
%		x	= domain
%		y	= data
%		n	= nth moment == 1 by default
%
%	m = int(x^n * y) / int(y)
%
%	where int is integral over all x
%	integral is calculated using trapezium rule
%
%	v1.0

if ~exist('n', 'var')
	n = 1;
end

p = ((x.^n)*ones(1, size(y,2))).*y;

% mom = sum(xi^n * yi * dxi) / sum(yi * dxi)
%	i  = pixel number
%	dx = pixel width
mom = trapz(x, p)./trapz(x, y);