function trendm = trendinterp(trendq) ;
%
% Interpolate a quarterly trendq to monthly trendm in such a way that
%  average of trendm(3*i-1:3*i+1) = trendq(i).
%
%

[T,k] = size(trendq) ;

% Use the Matlab diag routine to construct the matrix A quickly:

A = diag(repmat(7/9,T,1)) + diag(repmat(1/9,T-1,1),-1) + diag(repmat(1/9,T-1,1),1) ;
A(1,1:2) = [1 0] ;
A(T,T-1:T) = [0 1] ;

midpoints = A \ trendq ;

trendm(2:3*T-2,1) = ...
    reshape([1 2/3 1/3]'*midpoints(1:T-1)' + [0 1/3 2/3]'*midpoints(2:T)', 3*T-3,1) ;

trendm(1) = midpoints(1) - 1/3 * (midpoints(2)-midpoints(1)) ;
trendm(3*T-1) = midpoints(T) ;
trendm(3*T) = midpoints(T) + 1/3 * (midpoints(T)-midpoints(T-1)) ;
