function D1 = iocat(dim,D1,D2)
% Concatenates LFT models along input (2) or output (1) dimension.

%   Copyright 1986-2010 The MathWorks, Inc.
[nw1,nz1] = iosize(D1.Blocks);
[nw2,nz2] = iosize(D2.Blocks);
[ny1,nu1] = iosize(D1);
[ny2,nu2] = iosize(D2);
NumFlag = isnumeric(D1.IC);
if dim==1
   % Output concatenation
   nuw2 = nu2+nw2;
   nyz1 = ny1+nz1;
   uperm2 = [1:nu2 nuw2+1:nuw2+nw1 nu2+1:nuw2];
   if isnumeric(D1.IC)
      IC1 = blkdiag(D1.IC,zeros(0,nw2));
      IC2 = blkdiag(D2.IC,zeros(0,nw1));
      IC = [IC1 ; IC2(1:ny2+nz2,uperm2)];
   else
      IC1 = appendGain(D1.IC,zeros(0,nw2));
      IC2 = appendGain(D2.IC,zeros(0,nw1));
      IC = iocat(1,IC1,ioperm(IC2,1:ny2+nz2,uperm2));
   end
   uperm = 1:nu1+nw1+nw2;
   yperm = [1:ny1 nyz1+1:nyz1+ny2 ny1+1:nyz1 nyz1+ny2+1:nyz1+ny2+nz2];
else
   % Input concatenation
   nuw1 = nu1+nw1;
   nyz2 = ny2+nz2;
   yperm2 = [1:ny2 nyz2+1:nyz2+nz1 ny2+1:nyz2];
   if isnumeric(D1.IC)
      IC1 = blkdiag(D1.IC,zeros(nz2,0));
      IC2 = blkdiag(D2.IC,zeros(nz1,0));
      IC = [IC1 , IC2(yperm2,1:nu2+nw2)];
   else
      IC1 = appendGain(D1.IC,zeros(nz2,0));
      IC2 = appendGain(D2.IC,zeros(nz1,0));
      IC = iocat(2,IC1,ioperm(IC2,yperm2,1:nu2+nw2));
   end
   uperm = [1:nu1 nuw1+1:nuw1+nu2 nu1+1:nuw1 nuw1+nu2+1:nuw1+nu2+nw2];
   yperm = 1:ny1+nz1+nz2;
end

% Combine blocks
if NumFlag
   D1.IC = IC(yperm,uperm);
else
   D1.IC = ioperm(IC,yperm,uperm);
end
D1.Blocks = append(D1.Blocks,D2.Blocks);

% Try simplifying uncertain blocks based on their AutoSimplify settings
if any(isUncertain(D1.Blocks))
   D1 = simplify(D1,-1,Inf);
end
