function [T]=TransisitionOperator(mask,M)
%  TranferOperator -- Build Matrix Representation of the transfrer operator
%                     associated with a vector subdivision operator
%  Usage
%    [T] = TransferOperatorMatrix(mask[,M])
%  Inputs
%    mask      subdivision mask, cell array of length N+1
%    M         size of chosen invariant subspace (see Description), default = N
%  Outputs
%    T         Matrix Representation of the transfer operator
%              associated with input mask
% 
%  Description
%    Given a finitely supported subdivision mask A=[A0,A1,...,AN]
%    of a subdivision operator (Ai rxr matrices):
%     S v(a) = \sum_{b} A(a-2b) v(b),
%    the transition operator associated with S is defined by:
%     T w(a) = \sum_{b,g} A(2a-b)^* w(b+g) A(g).
%    A careful inspection shows that W:=(l([-N,N]))^{rxr} (and, more
%    generally, any W:=(l([-M,M]))^r, M>=N) is
%    invariant under T.
%    This function provides a matrix representation of T|_W.
%    Notice that the dimension of this square
%    matrix is (2M+1)xr^2.
%

for i=1:length(mask), mask{i}=mask{i}'; end
% I was using the definition
% T w(a) = \sum_{b,g} A(2a-b) w(b+g) A(g)^*,
% hence the transpose is necessary.

N = length(mask)-1;
if nargin<2, M=N; end
r = size(mask{1},1);
isdouble=0;
if strcmp(class(mask{1}),'double')
   isdouble=1;
end

T=[]; w = zeros(r);
for th=(-M):M,
   for j=1:r, for i=1:r,
   	w(i,j) = 1;
      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   	COL=[];
   	for al=(-M):M,
    	  F_al = szeros(r,isdouble);
    	  for g=0:N,
      	   be = th-g;
         	if (2*al-be)>=0 & (2*al-be)<=N,
            	F_al = F_al + mask{2*al-be +1}*w*mask{g +1}.'/2;
         	end
      	end
      	COL = [COL; F_al(:)];
      end
      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      w(i,j)=0;
      T=[T COL];
   end, end
end

function Z=szeros(r,isdouble)
Z=zeros(r);
if ~isdouble, Z=sym(Z); end