function	[U2,U1,U0] = factor_SU2pow3_matrix(U)



%This function checks whether a matrix U is a

%tensor product of three SU(2) matrices,

%U=U2 otimes U1 otimes U0,

%If it is, the function finds U2, U1 and U0.





%Start by assuming that 

%U =

%[x*U1U0, y*U1U0]

%[-y'*U1U0, x'*U1U0]

%where U1U0 = U1 otimes U0.

%Will check validity of this assumption at the end.



x_sq = U(5:8,5)'*U(1:4,1);

y_sq = -U(5:8,1)'*U(1:4,5);

yh_x = U(1:4,5)'*U(1:4,1);

x = sqrt(x_sq);

y = sqrt(y_sq);



%sqrt ambiguous in its sign, 

%so the values just

%assigned to x and

%to y may have the wrong sign.

%Fix relative sign:

if(abs(y'*x-yh_x)>1e-9)

	x=-x;

end



U2=[x,y;-y',x'];

if(abs(x)>abs(y))

	U1U0 = U(1:4,1:4)/x;

else

	U1U0 = U(1:4,5:8)/y; 

end



[U1,U0] = factor_SU2pow2_matrix(U1U0);



%the moment of truth

kk = norm(U-kron(U2,kron(U1,U0)));





if(kk>1e-9)

	error("U doesnt factor into tensor product");

end