% seq - Base class for sequences
classdef seq 
%   properties (Constant, GetAccess=public)
%	   NIL = seq.nil
%	end

	methods (Abstract)
		x=head(a)     % seq(A) -> A.
		b=next(a)     % seq(A) -> seq(A).
		s=elsize(a)   % seq([E->T]) -> E:[[1,D]].
		s=tostring(a) % seq(A) -> string. 
	end

	methods (Sealed=true)
		function varargout=size(a,n)
			s=elsize(a);
			if nargin>1, varargout={s(n)}; 
			elseif nargout<=1, varargout={s};
			else 
				varargout=[num2cell(s(1:nargout-1)),{prod(s(nargout:end))}];
			end
			if any(isnan(cell2mat(varargout)))
				error('seq:Element size indeterminate');
			end
		end

		function l=length(a), l=max(size(a)); end
		function n=numel(a, varargin), n=prod(size(a)); end
		function n=end(a,k,m)
			if m==length(elsize(a)), n=size(a,k);
			else n=prod(size(a)); end
		end

		function n=count(a), 
			n=0; while ~isempty(a), n=n+1; a=next(a); end
		end

		function display(o)
			display(sprintf(' %s :: seq([[%s]])\n', tostring(o), ...
						catsep(',',cellmap(@fmt_dim,num2cell(elsize(o))))));
			function s=fmt_dim(x), if isnan(x), s='_'; else s=num2str(x); end; end
		end

		function y=subsref(a,S)
			switch S(1).type
			case '()', y=paren(a,S(1));
			case '.', fn=S(1).subs; y=map(@(z)getfield(z,fn),a);
			end
			if length(S)>1, y=subsref(y,S(2:end)); end
		end
	end

	methods
		function [x,b]=decons(a)
			% decons - head and tail of sequence
			% decons :: seq(A) -> A, seq(A).
			x=head(a);
			b=next(a);
		end

		function z=last(y)
			% last - get the last array in a sequence (there must be at least one element)
			% last :: seq(A) -> A.
			while ~isempty(y), x=y; y=next(x); end
			z=head(x);
		end

		function f=headfn(o), f=@()head(o); end


		% ------------- WRAPPERS FOR OTHER SEQ CLASSES ---------------------------
		function y=take(n,x), y=seq.taker.make(n,x); end; 
		function y=takewhile(f,x), y=seq.takewhile.make(f,x); end 

		function s=drop(n,s), 
			% drop - Drop the first n alements of sequence 
			% drop :: natural, seq(A) -> seq(A)
			for i=1:n
				if ~isempty(s), s=next(s); end; 
			end 
		end

		function s=dropwhile(f,s)
			% dropwhile - Drop elements of sequence that satisfy condition 
			% drop :: (A->bool), seq(A) -> seq(A)
			while ~isempty(s) && f(head(s)), s=next(s); end
		end

		function [x,tail]=nth1(n,X)
			% nth1 - Get nth element of a sequence starting at 1
			% nth1 :: natural, seq(X) -> X, seq(X) ~'tail of sequence'.
			z=drop(n-1,X);
			x=head(z);
			if nargout>=2, tail=next(z); end
		end	

		function Y=once(X), Y=cons(head(X),[]); end
		% once - equivalent to take(1,...)
		% once :: seq(A) -> seq(A).

		function x=map(fn,y) 
		% map - seq where values is a stateless function of another sequence
		%
		% map :: (A->B), seq(A) -> seq(B).
			if isempty(y), x=nil; else x=seq.mapper(fn,y); end; 
		end

		function x=select(fn,y)
			% select - Use boolean function to keep only matching elements from sequence
			%
			% select :: (A->bool), seq(A) -> seq(A).
			x=seq.selector.make(fn,y); 
		end 

		function y=cache(x), if isempty(x), y=nil; else y=seq.cacher(x); end; end
		function y=buffer(x,varargin), y=seq.bufferer.make(x,varargin{:}); end

		function x=mapaccum(fn,s0,y)
			if iscell(y), x=cmapaccum(fn,s0,y);
			elseif isempty(y), x=nil;
			else x=seq.stmapper(fn,s0,y); 
			end
		end

		function x=scanl(fn,e,y,varargin)
			% scanl - scanl combinator for sequences
			%
			% This function applies an associative operator to a list of arguments,
			% starting from the left using the given starting element.
			%
			% scanl :: 
			%    (X,Y->X)	~'associative binary operator',
			%    X         ~'initial element',
			%    seq(Y)    ~'a sequence'
			% -> seq(X).
			if iscell(y), x=cscanl(fn,e,y,varargin{:});
			elseif isempty(y), x=nil;
			else x=seq.scanner(fn,e,y,varargin{:}); 
			end
		end

		% concat - Concatenate sequences
		%
		% concat :: seq(seq(A)) -> seq(A).
		function y=concat(x), y=seq.concatenator.make(x); end 

		% append - Append one or more sequences
		%
		% append :: seq(A), seq(A), ... -> seq(A).
		function x=append(varargin), x=concat(cellseq(varargin)); end 

		function z=and(x,y), 
			if isempty(x)
				if isempty(y), z=nil; else z=y; end
			else
				if isempty(y), z=x; else z=seq.concatenator(cons(x,cons(y,nil))); end 
			end
		end

		function y=cycle(x)
			% cycle - cycles through input sequence repeatedly
			% 
			% cycle :: seq(A) -> seq(A).
			if isempty(x), error('seq: cannot cycle empty sequence'); end; 
			y=seq.cycler(x); 
		end

		function y=bindcat(x,f) 
			% bindcat - sort of monadic bind for sequences.
			%
			% bindcat ::
			%    seq(A)				~ 'the first sequence',
			%	  (A->seq(A))     ~ 'function to return second sequence given last element of first'
			% -> seq(A)          ~ 'resultant sequence'.
			%
			% The resulting sequence consists of the entire sequence represented by the
			% first parameter, followed by the sequence obtained by applying the second
			% parameter to the last element of the first sequence. 
			%
			% Example:
			%
			%    gather(2,bindcat(cellseq({1,2,3,4}),@(x)take(head(x),0)))
			%
			% ans = 1  2  3  4  0  0  0  0
			y=seq.binder(x,f); 
		end 

		function y=subsample(n,x), if isempty(x), y=nil; else y=seq.subsampler(n,x); end; end

		function x=zip(varargin)
			% zip - combine several sequences into one
			%
			% zip :: seq(A), seq(B), ... -> seq(cell {A,B,...}).
			x=seq.zipper(@tuple,varargin);
			function z=rtuple(varargin), z=varargin; end
		end

		function x=zipwith(fn,varargin), 
			% zipwith - apply function to several sequences
			%
			% zipwith :: 
			%    (A,B,...->X),
			%    seq(A), seq(B), ...
			% -> seq(X).
			if any(cell2mat(cellmap(@isempty,varargin))), x=nil;
			else x=seq.zipper(fn,varargin); 
			end
		end

		function x=zipaccum(fn,s0,varargin), 
			% zipaccum - apply stateful function to several sequences
			%
			% zipaccum :: 
			%    (A,B,...,S->X,S),
			%    S,
			%    seq(A), seq(B), ...
			% -> seq(X).
			if isseq(s0) && any(map(@iscell,varargin)) % catch dispatching errors
				z=czipcaccum(fn,s0,varargin{:});
			elseif any(cell2mat(cellmap(@isempty,varargin))), x=nil;
			else x=seq.stzipper(fn,s0,varargin); 
			end
		end

		function xx=unzip(y)
			% unzip - Separate sequence of tuples into several sequences
			%
			% unzip :: 
			%    seq({I:[D]->A(I)}).
			% -> {I:[D]->seq(A(I))}.
			%
			% Note: input MUST be a non-empty sequence of constant size cell arrays.
			% Output is a cell array of sequences of the same size and shape.

			xx=cell(size(y));
			for i=1:numel(xx)
				xx{i}=map(@(a)a{i},y);
			end
		end

		function y=merge(f,varargin), 
			if all(cell2mat(cellmap(@isempty,varargin))), y=nil;
			else y=seq.merger(f,varargin); end
		end


		% ------------------ MAPPERS ------------------------------------

		% BINOP - Binary operation
		%
		% binop :: seq(A), seq(B), (A,B->C), string -> seq(C).
		%
		% Three cases
		%	A is seq, B is an array
		%	A is array, B is seq
		function o=binop(A,B,fn,opstr)
			o=binfun(A,B,fn,@(a,b)sprintf('(%s%s%s)',a,b,opstr));
		end

		function o=binfun(A,B,fn,fmtstr)
			if isseq(A),
				if isseq(B), o=zipwith(fn,A,B); %@(o)fmtstr(tostring(A),tostring(B)));
				else
					o=map(@(x)fn(x,B),A); %@(o)fmtstr('.',tostring(B)));
				end
			else
				o=map(@(y)fn(A,y),B); %,@(o)fmtstr(tostring(A),'.'));
			end
		end

		% vecop - apply binary function to different sized array sequences
		%
		% vecop ::
		%    ([[D]],[[D]]->[[D]])  ~'some function requiring equal size args',
		%    seq [[DX]]            ~'first arg of size DX',
		%    seq [[DY]]            ~'second arg of size DY'
		% -> seq [[DZ]]            ~'result of size DZ' :- DZ=max(DX,DY).
		%
		% The input sequences must maintain the same size throughout.
		function Z=vecop(F,X,Y)
			DX=size(X); DY=size(Y); 
			E=max(length(DX),length(DY));
			EDX=pad1s(E,DX);
			EDY=pad1s(E,DY);
			if all(EDX>=EDY)
				S=EDX./EDY;
				Z=binop(X,Y,@(x,y)F(x,repmat(y,S)),['<' tostring(F) '>']);
			elseif all(EDY>=EDX)
				S=EDY./EDX;
				Z=binop(X,Y,@(x,y)F(repmat(x,S),y),['<' tostring(F) '>']);
			else
				DZ=max(EDX,EDY);
				Z=binop(X,Y,@(x,y)F(repmat(x,DZ./EDX),repmat(y,DZ./EDY)),['<' tostring(F) '>']);
			end
		end


		function h=plot(A,varargin), h=plotseq(@(x)plot(x,varargin{:}),A); end
		function h=imagesc(A,varargin), h=plotseq(@(x)imagesc(x,varargin{:}),A); end

		function o=isfinite(A),  o=map(@isfinite, A); end
		function o=isinf(A),     o=map(@isinf, A); end
		function o=isnan(A),     o=map(@isnan, A); end
		function y=powspec(x),   y=map(@powspec,x); end
		function y=magspec(x),   y=map(@magspec,x); end
		function y=phasespec(x), y=map(@phasespec,x); end
		function o=uminus(A),    o=map(@uminus, A); end
		function y=exp(x),       y=map(@exp,x); end
		function y=cos(x),       y=map(@cos,x); end
		function y=sin(x),       y=map(@sin,x); end
		function y=abs(x),       y=map(@abs,x); end
		function y=sqrt(x),      y=map(@sqrt,x); end
		function y=tanh(x),      y=map(@tanh,x); end
		function y=log(x),       y=map(@log,x); end
		function y=log10(x),     y=map(@log10,x); end
		function y=log2(x),      y=map(@log2,x); end
		function o=ctranspse(A), o=map(@ctranspose,A); end
		function o=transpse(A),  o=map(@transpose,A); end
		function y=fft(x),       y=map(@fft,x); end
		function y=ifft(x),      y=map(@ifft,x); end

		function o=reshape(source,varargin)
			% reshape - Map reshape over elements of sequence
			%
			% reshape :: seq [Size->A], ... - > seq [Size1->A].
			% Works exactly like the usual reshape function but when applied
			% to a sequence object, returns a new sequence.

			sz=tosize(varargin{:});
			o=map(@(x)reshape(x,varargin{:}),source);

			function s=charfn(sz,o)
				s=sprintf('%s >> reshape[%s]',tostring(source(o)),tostring(sz));
			end
		end

		function y=paren(a,S)
			% paren - Map application of subsref with parentheses to sequence
			%
			% paren :: seq(A), subs -> seq(B) :- (subsref :: A, subs -> B).

			% NOTE TO SELF: it would be good to work out the size of the
			% array that will result when the function is evaluated, to
			% save map evaluating it once on construction.
			y=map(@(z)subsref(z,S),a); % 'charfn',@(o)charfn(tostring(S.subs{:}),o));
						
			function s=charfn(argstr,o)
				s=sprintf('%s >> (%s)',tostring(source(o)),argstr);
			end
		end
			
		function o=plus(A,B),    o=binop(A,B,@plus,'+'); end
		function o=power(A,B),   o=binop(A,B,@power,'.^'); end
		function o=eq(A,B),      o=binop(A,B,@eq,'=='); end
		function o=ge(A,B),      o=binop(A,B,@ge,'>='); end
		function o=gt(A,B),      o=binop(A,B,@gt,'>'); end
		function o=ldivide(A,B), o=binop(A,B,@ldivide,'.\'); end
		function o=le(A,B),      o=binop(A,B,@le,'<='); end
		function o=lt(A,B),      o=binop(A,B,@lt,'<'); end
		function o=times(A,B),   o=binop(A,B,@times,'.*'); end
		function o=minus(A,B),   o=binop(A,B,@minus,'-'); end
		function o=mldivide(A,B),o=binop(A,B,@mldivide,'\'); end
		function o=mod(A,B),     o=binop(A,B,@mod,'mod'); end
		function o=mrdivide(A,B),o=binop(A,B,@mrdivide,'/'); end
		function o=rdivide(A,B), o=binop(A,B,@rdivide,'./'); end
		function o=mtimes(A,B),  o=binop(A,B,@mtimes,'*'); end


		% max - max mapped over sequence (ie NOT aggregate)
		function o=cat(dim,varargin)
			if length(varargin)==2
				o=binfun(varargin{1},varargin{2},@(a,b)cat(dim,a,b),@(a,b)sprintf('cat(%d,%s,%s)',dim,a,b));
			else
				o=zipwith(@catdim,varargin{:});
			end
			function x=catdim(varargin), x=cat(dim,varargin); end
		end

		function o=vertcat(varargin)
			if length(varargin)==2
				o=binfun(varargin{1},varargin{2},@vertcat,@(a,b)sprintf('[%s;%s]',a,b));
			else
				o=zipwith(@vertcat,varargin{:});
			end
		end

		function o=horzcat(varargin)
			if length(varargin)==2
				o=binfun(varargin{1},varargin{2},@horzcat,@(a,b)sprintf('[%s,%s]',a,b));
			else
				o=zipwith(@horzcat,varargin{:});
			end
		end
	end
end
