% azeromean -  additive normalisation arrow
%
% azeromean ::
%    options {
%       offset  :: real   /0  ~'initial offsets';
%       rate    :: nonneg/1e-7 ~'offset adaptation rate';
%       tension :: real   /0  ~'smoothing strength';
%    }
% -> arrow({[[N]]},{[[N]]},azeromean_state).

function o=azeromean(model,varargin)
	opts=options('offset',nan,'rate',1e-7,'tension',0,varargin{:});
	score=scorefn(model);

	rate=opts.rate;
	if opts.tension>0
		tension=opts.tension;
		o=loop(@update_t,@(s)repmat_to(opts.offset,[s(1),1]));
	else
		o=loop(@update,@(s)repmat_to(opts.offset,[s(1),1]));
	end

	function [y,offset]=update(x,offset)
		nans=isnan(offset);
		offset(nans)=x(nans);
		y = x-repmat(offset,1,size(x,2));
		delta = rate*sum(score(y),2);
		delta(~isfinite(delta))=0;
		offset = offset + delta;
	end

	function [y,offset]=update_t(x,offset)
		nans=isnan(offset);
		offset(nans)=x(nans);
		w = size(x,2);
		y = x;
		for i=1:w
			y(:,i)=y(:,i)-offset;
			delta = rate*(score(y(:,i))+reg(offset));
			delta(~isfinite(delta))=0;
			offset = offset + delta;
		end
	end

	function g=reg(h)
		zz=zeros(1,size(h,2));
		g=tension*[zz;diff(h,2,1);zz];
	end
end
