    function [c,s] = Rotate(x,y)
  % x and y are column n-vectors.
  % c and s are a cosine-sine pair with the property that the columns
  % of [x y]*[c s; -s c] are orthogonal to each other.
  
  % If U'*A*V = S is the svd, then from  A*V = U*S it follows that the 
  % columns of A*V are orthogonal.
    [U,S,V] = svd([x y],0);
    
  % A 2-by-2 orthogonal matrix is either a rotation [c s;-s c] or a 
  % reflection [c s ; s -c].
  
  % If V is a reflection, make it a rotation...
    if V(1,2)*V(2,1) >= 0
       V(:,1) = -V(:,1);
    end
  % Read off c and s..
    c =  V(1,1);
    s = -V(2,1);
    
  % However, since (-c,-s) works if (c,s) works, the V update
  % to make it a rotation is irrelevant.
  
  % Also note that if (c,s) works so does (-s,c)