function[A,degA] = spoldiv(G,degG,C,degC)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%       [A,degA]=spoldiv(G,degG,C,degC)   %
%                                         %
%	MATLAB function  spoldiv extract  %
%	given  left  divisor  G from any  %
%	polynomial    matrix    C         %
%	    C = G . A                     %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%       S. Pejchova,  March 11th 1993     %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%	size of C                         %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[rG,cG]=polsize(G,degG);
[rC,cC]=polsize(C,degC);
[m,n]=size(C);
[k,l]=size(G);
if m ~=k
	if k<m
	     B = zeros(m-k,l);
	     G = [G;B];
	     [k,l]=size(G);
	else
	     error('spoldiv: Inconsistent dimensions of input matrices');
	end
end

rA=cG;
cA=cC;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%	degA                              %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

degA=degC;
if degA < 0
	A = zeros(rC,cC);
	return;
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%       matrix A creating                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

A =[];
degB = degG;

for i = 1 : rA

    if i~=1
	    B=zeros(1,(i-1)*(degB+1));
	    f=1;
	    for s=1:(degB+1)
		for t=1:(i-1)
		    B(1,f)=G(i,cG*(s-1)+t); 
		    f=f+1;

		end
	    end


	[D,degD] = polmul(B,degB,A,degA);

	D = -1*D;
	[E,degE]=poladd(C(i,:),degC,D,degD);
	if degE<degC
		H=zeros(1,cC*(degC-degE));
		E=[E H];
		degE=degC;
	end

   else
       E=C(1,:);
       degE = degC;

   end


	[a,dega]=polen(G,degG,i,i);

	if a==0
	   A = [A;E];
	else        
	  for j = 1:cA
		[p,degp]=polen(E,degE,1,j);
		[u,degu,v,degv]=poldiv(p,degp,a,dega);
		if v ~=0
		   error('spoldiv: Division with remainder');
		end

	      [F,degF]=polseten(E,degE,u,degu,1,j);
	      E=F;
	  end

          [raA,caA]=size(A);
          [raE,caE]=size(E);
          if degA<degE
             J=zeros(raA,caE-caA);
             A=[A J];
             degA=degE;
             caA=caE;
          else
             if degA>degE
             J=zeros(raE,caA-caE);
             E=[E J];
             end
          end


          A = [A;E];
        end
end        

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%       setting degA                      %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

          cA=caA/(degA+1);

          if degA>0
             while A(:,cA*degA+1:cA*(degA+1))==zeros(rA,cA)
                   A(:,cA*degA+1:cA*(degA+1))=[];
                   degA=degA-1;
                   if degA<1,return,end
             end
          end