Here is an attempt:
center1 = {1, 1, 3};
center2 = {2, 2, 1};
intersections = SolveValues[{h == 0, g == 0, i == 0}, {x, y, z}]
Graphics3D[{{Opacity[.5], Sphere[center1, 2],
Sphere[center2, Sqrt[2]],
InfinitePlane[{1, 1, 0}, {{0, 1, 0}, {0, 0, 1}}]},
Thickness[.01], PointSize[.03],
Point[intersections],
Line[{center1, intersections[[1]], center2}]}]
angle = VectorAngle[intersections[[2]] - center1, {0, 0, -1}]