function [one,two] = twoCircles(x1,y1,r1,x2,y2,r2)
% [one,two] = twoCircles(x1,y1,r1,x2,y2,r2)
% Finds intersection points of two circles
one = [(-(power(r1,2)*power(x1,2)) + ...
        power(r2,2)*power(x1,2) + power(x1,4) + ...
        2*power(r1,2)*x1*x2 - 2*power(r2,2)*x1*x2 - ...
        2*power(x1,3)*x2 - power(r1,2)*power(x2,2) + ...
        power(r2,2)*power(x2,2) + 2*x1*power(x2,3) - ...
        power(x2,4) + power(x1,2)*power(y1,2) - ...
        power(x2,2)*power(y1,2) - 2*power(x1,2)*y1*y2 + ...
        2*power(x2,2)*y1*y2 + power(x1,2)*power(y2,2) - ...
        power(x2,2)*power(y2,2) + ...
        y1*sqrt(-(power(x1 - x2,2)*...
             (power(r1,4) + ...
               power(-power(r2,2) + power(x1,2) - ...
                 2*x1*x2 + power(x2,2) + power(y1,2) - ...
                 2*y1*y2 + power(y2,2),2) - ...
               2*power(r1,2)*...
                (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                  power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                  power(y2,2))))) - ...
        y2*sqrt(-(power(x1 - x2,2)*...
             (power(r1,4) + ...
               power(-power(r2,2) + power(x1,2) - ...
                 2*x1*x2 + power(x2,2) + power(y1,2) - ...
                 2*y1*y2 + power(y2,2),2) - ...
               2*power(r1,2)*...
                (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                  power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                  power(y2,2))))))/...
      (2.*(x1 - x2)*(power(x1,2) - 2*x1*x2 + ...
          power(x2,2) + power(y1 - y2,2))),...
  (-(power(r1,2)*y1) + power(r2,2)*y1 + ...
        power(x1,2)*y1 - 2*x1*x2*y1 + power(x2,2)*y1 + ...
        power(y1,3) + power(r1,2)*y2 - power(r2,2)*y2 + ...
        power(x1,2)*y2 - 2*x1*x2*y2 + power(x2,2)*y2 - ...
        power(y1,2)*y2 - y1*power(y2,2) + power(y2,3) - ...
        sqrt(-(power(x1 - x2,2)*...
            (power(r1,4) + ...
              power(-power(r2,2) + power(x1,2) - ...
                2*x1*x2 + power(x2,2) + power(y1,2) - ...
                2*y1*y2 + power(y2,2),2) - ...
              2*power(r1,2)*...
               (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                 power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                 power(y2,2))))))/...
      (2.*(power(x1,2) - 2*x1*x2 + power(x2,2) + ...
          power(y1 - y2,2)))];
   two = [(-(power(r1,2)*power(x1,2)) + ...
        power(r2,2)*power(x1,2) + power(x1,4) + ...
        2*power(r1,2)*x1*x2 - 2*power(r2,2)*x1*x2 - ...
        2*power(x1,3)*x2 - power(r1,2)*power(x2,2) + ...
        power(r2,2)*power(x2,2) + 2*x1*power(x2,3) - ...
        power(x2,4) + power(x1,2)*power(y1,2) - ...
        power(x2,2)*power(y1,2) - 2*power(x1,2)*y1*y2 + ...
        2*power(x2,2)*y1*y2 + power(x1,2)*power(y2,2) - ...
        power(x2,2)*power(y2,2) - ...
        y1*sqrt(-(power(x1 - x2,2)*...
             (power(r1,4) + ...
               power(-power(r2,2) + power(x1,2) - ...
                 2*x1*x2 + power(x2,2) + power(y1,2) - ...
                 2*y1*y2 + power(y2,2),2) - ...
               2*power(r1,2)*...
                (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                  power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                  power(y2,2))))) + ...
        y2*sqrt(-(power(x1 - x2,2)*...
             (power(r1,4) + ...
               power(-power(r2,2) + power(x1,2) - ...
                 2*x1*x2 + power(x2,2) + power(y1,2) - ...
                 2*y1*y2 + power(y2,2),2) - ...
               2*power(r1,2)*...
                (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                  power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                  power(y2,2))))))/...
      (2.*(x1 - x2)*(power(x1,2) - 2*x1*x2 + ...
          power(x2,2) + power(y1 - y2,2))),...
   (-(power(r1,2)*y1) + power(r2,2)*y1 + ...
        power(x1,2)*y1 - 2*x1*x2*y1 + power(x2,2)*y1 + ...
        power(y1,3) + power(r1,2)*y2 - power(r2,2)*y2 + ...
        power(x1,2)*y2 - 2*x1*x2*y2 + power(x2,2)*y2 - ...
        power(y1,2)*y2 - y1*power(y2,2) + power(y2,3) + ...
        sqrt(-(power(x1 - x2,2)*...
            (power(r1,4) + ...
              power(-power(r2,2) + power(x1,2) - ...
                2*x1*x2 + power(x2,2) + power(y1,2) - ...
                2*y1*y2 + power(y2,2),2) - ...
              2*power(r1,2)*...
               (power(r2,2) + power(x1,2) - 2*x1*x2 + ...
                 power(x2,2) + power(y1,2) - 2*y1*y2 + ...
                 power(y2,2))))))/...
      (2.*(power(x1,2) - 2*x1*x2 + power(x2,2) + ...
          power(y1 - y2,2)))];