samer@1: function s2=construct(sig)
samer@1: 
samer@1: 	fprintf('\nResampling %s to %g Hz...\n',tostring(sig.source),sig.rate);
samer@1: 	f1=rate(sig.source);
samer@1: 	f2=sig.rate;
samer@1: 	c=channels(sig);
samer@1: 
samer@1: 	[p,q]=rat(f2/f1,1e-12);
samer@1: 	sysobj=design(sig,p,q);
samer@1: 
samer@1: 
samer@1: 	delay=(length(sysobj.Numerator)-1)/2; 
samer@1: 	outdelay = ceil(delay/q);
samer@1: 	indelay  = ceil(delay/p);
samer@1: 	m=max(1,floor(sig.opts.bs/q)); % try to approximate requested input block size
samer@1: 	fprintf('Input/output delays are %d and %d.\n',indelay,outdelay);
samer@1: 	fprintf('Read block size is %d.\n\n',m*q);
samer@1: 
samer@1: 	s1=construct(sig.source & sigarray(zeros(channels(sig.source),indelay)));
samer@1: 	r1=s1.reader(m*q); % one and only input reader
samer@1: 	chunk=uint32(m*p); % filter output block size
samer@1: 	CHUNK=1:chunk;
samer@1: 	queue=[];
samer@1: 
samer@1: 	s2.start   = s1.start;
samer@1: 	s2.stop    = s1.stop;
samer@1: 	s2.dispose = @dispose;
samer@1: 	s2.reader  = @reader;
samer@1: 
samer@1: 	% drop some samples to account for filter delay
samer@1: 	s2.start();
samer@1: 	sigreadn(s2,outdelay);
samer@1: 	s2.stop();
samer@1: 
samer@1: 	function dispose, s1.dispose(); release(sysobj); end 
samer@1: 	function r2=reader(n)
samer@1: 		outbuf=zeros(c,n);
samer@1: 		r2=@next;
samer@1: 
samer@1: 		function [x,rem]=next 
samer@1: 			% transfer up to n queued samples to outbuf
samer@1: 			pos=uint32(size(queue,2));
samer@1: 			if pos>=n % enough samples already waiting
samer@1: 				outbuf=queue(:,1:n);
samer@1: 				queue=queue(:,n+1:end);
samer@1: 				rem=0;
samer@1: 			else
samer@1: 				if pos==0, toread=n;
samer@1: 				else % use up queue
samer@1: 					outbuf(:,1:pos)=queue; queue=[];
samer@1: 					toread=n-pos;
samer@1: 				end
samer@1: 
samer@1: 				% transfer complete chunks
samer@1: 				while toread>=chunk
samer@1: 					[outbuf(:,pos+CHUNK),rem]=filter_next;
samer@1: 					toread=toread-chunk;
samer@1: 					pos=pos+chunk;
samer@1: 					if rem>0 % we ran out of samples
samer@1: 						rem=rem+toread; % account for rest of missing samples (not just this chunk) 
samer@1: 						toread=0; % causes immediate exit after this
samer@1: 					end
samer@1: 				end
samer@1: 
samer@1: 				% transfer partial chunk if necessary
samer@1: 				if toread>0
samer@1: 					[y,rem]=filter_next;
samer@1: 					outbuf(:,pos+1:n)=y(:,1:toread);
samer@1: 					queue=y(:,toread+1:chunk-rem);
samer@1: 					rem=max(0,toread-(chunk-rem));
samer@1: 				end
samer@1: 			end
samer@1: 			x=outbuf;
samer@1: 		end
samer@1: 
samer@1: 		% returns the next block of m*p output samples 
samer@1: 		function [x,rem]=filter_next
samer@1: 			[y,rem1]=r1();
samer@1: 			if rem1>0,
samer@1: 				y(:,end-rem1+1:end)=zeros(c,rem1); % pad with zeros
samer@1: 				rem=uint32(ceil(rem1*p/q)); % round down the number of valid samples
samer@1: 			else
samer@1: 				rem=rem1;
samer@1: 			end
samer@1: 			x=step(sysobj,y')';
samer@1: 		end
samer@1: 	end
samer@1: end