Find circle intersection points

circle_intersections.pde

PVector[] circleIntersections(PVector p1, float r1, PVector p2, float r2)
{
    float d, a, h;
    float dx, dy;

    dx = p2.x - p1.x;
    dy = p2.y - p1.y;

    PVector dvect = new PVector(p1.x - p2.x, p1.y - p2.y);
    d = dvect.mag();
    a = (r1 * r1 - r2 * r2 + d * d) / (2 * d);
    h = sqrt(r1 * r1 - a * a);

    // upper intersection
    PVector sp1 = new PVector();
    sp1.x = p1.x + (a/d) * dx - (h/d) * dy;
    sp1.y = p1.y + (a/d) * dy + (h/d) * dx;
    sp1.z = p1.z; // z-coordinate is not used

    // lower intersection
    PVector sp2 = new PVector();
    sp2.x = p1.x + (a/d) * dx + (h/d) * dy;
    sp2.y = p1.y + (a/d) * dy - (h/d) * dx;
    sp2.z = p1.z; // z-coordinate is not used

    // return both intersection points in an array
    return new PVector[] {sp1, sp2};
}