alp0=-0.1896;
alp1= 0.0919;
bet0=-0.1757;
bet1= 0.0882;
gam0=-0.1788;
gam1= 0.0882;
gam2= 0.0919;
%
rho1=(-alp1+sqrt(alp0*alp0-4*alp1*alp1))/(2*alp1);
if abs(rho1)>1 ;
rho1=(-alp1-sqrt(alp0*alp0-4*alp1*alp1))/(2*alp1);
if abs(rho1)>1 break; end;
end;
rho0=-alp1/rho1;
%
a2=alp1*bet1-gam2*gam1;
a1=alp0*bet1+alp1*bet0-gam0*(gam1+gam2);
a0=alp0*bet0+2*alp1*bet1-gam2*gam2-gam1*gam1-gam0*gam0;
%
if a2~=0;
  kaph1=(-a1/a2+sqrt(a1*a1/a2/a2-4*(a0/a2-2)))/2;
  kaph2=(-a1/a2-sqrt(a1*a1/a2/a2-4*(a0/a2-2)))/2;
  kap1=(kaph1+sqrt(kap1*kap1-4))/2;
  if abs(kaph1)>1;
    kap1=(kaph1-sqrt(kap1*kap1-4))/2;
    if abs(kaph1)>1 break; end;
  end;
%
  kap2=(kaph2+sqrt(kap2*kap2-4))/2;
  if abs(kaph2)>1 ;
    kap2=(kaph2-sqrt(kap2*kap2-4))/2;
    if abs(kaph1)>1 break; end;
  end;
  kap0=a2/kap1/kap2;
%
else;
  kap1=0;
  kap2=(-a0+sqrt(a0*a0-4*a1*a1))/2/a1;
  if kap2>1;
    kap2=(-a0-sqrt(a0*a0-4*a1*a1))/2/a1;
    if kap2>1; break; end;
  end;
  kap0=-a1/kap2;
end;
%
xi0 =gam2+gam0*rho1+gam1*rho1*rho1;
xi1 =sqrt(kap0)*(1-kap1*rho1)*(1-kap2*rho1);
xi2 =sqrt(xi0*xi0+xi1*xi1);
psi0=xi0/xi2;
psi1=xi1/xi2;
kap3=-(kap1+kap2);
%
ome0=-rho1*psi0*sqrt(rho0);
ome1= psi0*sqrt(rho0);
zet0=-psi1*sqrt(rho0);
zet1= rho1*psi1*sqrt(rho0);
eta0=(psi0*gam2+psi1*sqrt(kap0))/sqrt(rho0);
eta1=(psi0*gam0+rho1*psi0*gam2+rho1*psi1*sqrt(kap0)+psi1*sqrt(kap0)*kap3)/...
      sqrt(rho0);
pi0 =(psi1*gam2-psi0*sqrt(kap0))/sqrt(kap0)/rho1;
pi1 =(psi1*gam0+psi1*gam2/rho1-psi0*sqrt(kap0)*kap3-psi0*sqrt(kap0)/rho1)/...
      sqrt(kap0)/rho1;
%
[ome0 ome1]
[eta0 eta1]
[zet0 zet1]
[pi0  pi1]

