function [xfeas,ubout] = usmuwcp2(m,blk,bidx,lbnd,beta,ftdidx,opt)

% Copyright 2003-2004 The MathWorks, Inc.

% Use in conjunction with UBALL, where the ordering of the variables
% is as done below.  typically, BETA should be a bound that easily
% passes as a worst-case performance bound
if nargin==5
   opt = 'cuw';
end

mu3blk = [];
nblk = size(blk,1);
for i=1:nblk
   if blk(i,1)==1 & blk(i,2)==1
      % repeated real or complex scalar
      mu3blk = [mu3blk;blk(i,3) 0];
   elseif blk(i,3)>0
      % full block, make it not repeated
      mu3blk = [mu3blk;blk(i,3)*blk(i,[1 2])];
   end
end  
DGbeta = max([ftdidx.grealidx ftdidx.gimagidx ftdidx.drealidx ftdidx.dimagidx]) + 1;
        
dims = ndims(m);
szm = size(m);
ny = szm(1);
nu = szm(2);
ne = szm(1) - bidx.rdimm;
nd = szm(2) - bidx.cdimm;
mublk = [mu3blk;nd ne];
cnt = 0;
fnd = 0;
CNTMAX = 2;
lower = lbnd;
upper = 2*beta - lower;
go = 1;
xydata = [];

m11 = m(1:bidx.rdimm,1:bidx.cdimm);
optnowarning = opt(find(opt=='d'));
bnds11 = mussv(m11,mu3blk,optnowarning);
if bnds11(1) < 1
       % OK
else
       % bound will never work
       if bnds11(2)>=1
           % really have a problem
           % XXXMAJOR - will get INF later
       end
end

while cnt<CNTMAX & go==1
   if cnt>=3
      % CNTMAX needs to be at least 4 for this to activate 
      betaval = LOCALestibeta(xydata,lbnd,beta);
   else
      betaval = (lower+upper)/2;   
   end
   fac = 1/sqrt(betaval);
   
   sclm = m;
   sclm(bidx.rdimm+1:bidx.rdimm+ne,:) = fac*sclm(bidx.rdimm+1:bidx.rdimm+ne,:);
   sclm(:,bidx.cdimm+1:bidx.cdimm+nd) = fac*sclm(:,bidx.cdimm+1:bidx.cdimm+nd);
   
   [bnds,info] = mussv(sclm,mublk,opt);
   xydata = [xydata;betaval bnds(1)];

   bpass = min(xydata(find(xydata(:,2)<1),1));
   loc = find(xydata(:,1)==bpass);
   bfail = max(xydata(find(xydata(:,2)>1),1));

   if bnds(1)<1
      fnd = 1;
      gbnds = bnds;
      gdvec = info.dvec;
      ggvec = info.gvec;
      gfac = fac;
      gbetaval = betaval;
      upper = betaval;
      if bnds(1)>0.9995
         go = 0;
      end
   else
      if cnt==0
         xfeas = [];
         ubout = inf;
         go = 0;
      end
      lower = betaval;
   end
   cnt = cnt + 1;   
end

if fnd ==1
   ubout = gbetaval;
   xfeas = zeros(DGbeta,1);
   [Dr,Dc,Grc,Gcr] = ynftdam2(gdvec,ggvec,mublk,gbnds(1));
   %tmp = sclm'*Dr*sclm - (1.0000001*bnds(1))^2*Dc + sqrt(-1)*(Gcr*sclm-sclm'*Grc);
   %eig(tmp)
   % Dr is NY x NY, Dc is NU x NU, Gcr is NU x NY, Grc is NY x NU
   alpha = Dr(ny,ny);
   fix = gbetaval/alpha;
   xfeas(ftdidx.grealidx) = fix*real(Grc(ftdidx.grealselect));
   xfeas(ftdidx.gimagidx) = fix*imag(Grc(ftdidx.gimagselect));
   xfeas(ftdidx.drealidx) = fix*real(Dr(ftdidx.drealselect));
   xfeas(ftdidx.dimagidx) = fix*imag(Dr(ftdidx.dimagselect));
   xfeas(end) = gbetaval^2;
   %Dr = fix*Dr;
   %Grc = fix*Grc;
   %eig(m'*DD*m - DDD + sqrt(-1)*(Grc*m-m'*Grc))
end

function [bguess] = LOCALestibeta(bfdata,blow,bhigh)

szd = size(bfdata);
targetvalue = .999999;

if szd(1)==2
   mI = [bfdata(:,1) [1;1]]\bfdata(:,2);
   bguess = max([eps (targetvalue-mI(2))/mI(1)]);
else
   [dum,idx] = sort(bfdata(:,2));
   b = bfdata(idx,1);
   f = bfdata(idx,2);
   gt = find(f>1);
   lt = find(f<1);
   if length(gt)>0
      newdata = [b(lt(end)) f(lt(end));b(gt(1)) f(gt(1))];
      bguess = LOCALestibeta(newdata,blow,bhigh);
   else
      newdata = [b(lt(end-1:end)) f(lt(end-1:end))];
      bguess = LOCALestibeta(newdata,blow,bhigh);
   end
end