function [y,phi]=block_sr(o,phi,cutoff,f)
% block_sr - Generate block of blit signal data
%
% block_sr :: 
%    blit       ~'blit object'
%    0--1       ~'initial phase (in cycles)'
%    N:natural  ~'number of samples to compute',
%    0--1       ~'normalised frequency (1=sampling freq)',
%    0--0.5     ~'cutoff harmonics above this frequency [0.5]',
% -> [[1,N]]    ~'band-limited impulse train',
%    0--1       ~'initial phase for next block'.

m = max(2,2*floor(max(cutoff./f))); % number of harmonics to keep
u = cumsum([phi,f],2);
y = m*f.*diric(2*pi*u(1:end-1),m);
phi = mod(u(end),2);


