view sequences/@seq/meandata.m @ 61:eff6bddf82e3 tip

Finally implemented perceptual brightness thing.

function m=meandata(x,dim,varargin)
% MEANDATA - Mean of a sequence of arrays
%
% meandata ::
%    seq [D]      ~ sequence of E dimensional arrays of size D, ie D::[[E]]
%    I:[[N]->[E]] ~ array of dimension numbers (between 1 and E) to do mean over
%    options {
%       check_nans :: {0,1}/0 ~'whether to check for (and ignore) nans'
%    }
% -> [R]          ~ array of size R, where R=set(D,I,1).
%
% if X is a sequence of 5 x 12 x 50 arrays, then
%    M=meandata(X,[1,3])
% returns a 1 x 12 array.
%
% Options: See ITERATE for applicable options

	opts=options('check_nans',0,varargin{:});

	S.total=0;
	S.count=0;

	if opts.check_nans, % world of pain
		if isempty(dim)
			S=foldl(@mean0_nan_st,S,x,varargin{:});
		elseif isscalar(dim)
			S=foldl(@mean_nan_st,S,x,varargin{:});
		else
			S=foldl(@meandims_nan_st,S,x,varargin{:});
		end
		m=S.total./S.count;
	else
		if isempty(dim)
			S=foldl(@mean0_st,S,x,varargin{:});
		elseif isscalar(dim)
			S=foldl(@mean_st,S,x,varargin{:});
		else
			S=foldl(@meandims_st,S,x,varargin{:});
		end
		m=S.total/S.count;
	end


	% non-nan checking accumulators
	function S=meandims_st(S,X)
		S.total=S.total+sumdims(X,dim);
		S.count=S.count+prod(sizedims(X,dim));
	end

	function S=mean_st(S,X)
		S.total=S.total+sum(X,dim);
		S.count=S.count+size(X,dim);
	end
	
	function S=mean0_st(S,X)
		S.total=S.total+X;
		S.count=S.count+1;
	end

	% nan checking accumulators
	function S=meandims_nan_st(S,X)
		I=isnan(X); X(I)=0;
		S.total=S.total+sumdims(X,dim);
		S.count=S.count+sumdims(~I,dim);
	end

	function S=mean_nan_st(S,X)
		I=isnan(X); X(I)=0;
		S.total=S.total+sum(X,dim);
		S.count=S.count+sum(~I,dim);
	end
	
	function S=mean0_nan_st(S,X)
		I=isnan(X); X(I)=0;
		S.total=S.total+X;
		S.count=S.count+(~I);
	end
end