You can get an expression for the centre in terms of 3D positions of 3 points on the target.
The centre is defined in terms of 3 planes:
- The plane in which the circle lies
- Any two of the planes half way between points 1 & 2, 2 & 3, or 3 & 1.
Put the coordinates in a 3x3 matrix
M1 = Through /@ Array[{X, Y, Z}, {3}]
(* {{X[1], Y[1], Z[1]}, {X[2], Y[2], Z[2]}, {X[3], Y[3], Z[3]}} *)
We can define a test for (x,y,z)
co-planar in terms of projective coordinates
Append[M1, {x, y, z}]
(* {{X[1], Y[1], Z[1]}, {X[2], Y[2], Z[2]}, {X[3], Y[3], Z[3]}, {x, y,
z}} *)
M = Append[#, 1] & /@ %
(* {{X[1], Y[1], Z[1], 1}, {X[2], Y[2], Z[2], 1}, {X[3], Y[3], Z[3],
1}, {x, y, z, 1}} *)
planar = FullSimplify[Det[M] == 0]
(* z (X[3] (Y[1] - Y[2]) + X[1] (Y[2] - Y[3]) + X[2] (-Y[1] + Y[3])) +
X[3] Y[2] Z[1] + x Y[3] Z[1] + x Y[1] Z[2] +
X[1] Y[3] Z[2] + (X[2] Y[1] - X[1] Y[2] + x (-Y[1] + Y[2])) Z[3] +
y (X[3] (-Z[1] + Z[2]) + X[2] (Z[1] - Z[3]) +
X[1] (-Z[2] + Z[3])) ==
x Y[2] Z[1] + X[2] Y[3] Z[1] + X[3] Y[1] Z[2] + x Y[3] Z[2] *)
Now find the condition that (x,y,z) is equidistant from the 3 points
# - {x, y, z} & /@ M1
(* {{-x + X[1], -y + Y[1], -z + Z[1]}, {-x + X[2], -y + Y[2], -z +
Z[2]}, {-x + X[3], -y + Y[3], -z + Z[3]}} *)
Equal @@ (# . # & /@ %)
(* (-x + X[1])^2 + (-y + Y[1])^2 + (-z + Z[1])^2 == (-x + X[2])^2 + (-y +
Y[2])^2 + (-z + Z[2])^2 == (-x + X[3])^2 + (-y + Y[3])^2 + (-z +
Z[3])^2 *)
equidistant = Simplify[Expand /@ %]
(* 2 x X[1] + X[2]^2 + 2 y Y[1] + Y[2]^2 + 2 z Z[1] + Z[2]^2 ==
X[1]^2 + 2 x X[2] + Y[1]^2 + 2 y Y[2] + Z[1]^2 + 2 z Z[2] &&
2 x X[2] + X[3]^2 + 2 y Y[2] + Y[3]^2 + 2 z Z[2] + Z[3]^2 ==
X[2]^2 + 2 x X[3] + Y[2]^2 + 2 y Y[3] + Z[2]^2 + 2 z Z[3] *)
Mathematica has no difficulty in solving this (but the result is quite large).
Note that much of the complexity comes because symbolic solution of 3x3 matrices gives complicated expressions (in general)
solution1 = Solve[equidistant && planar, {x, y, z}];
To verify the solution, define a circle
centre = {3, 4, 5};
radius = 2;
normal = {3/13, 4/13, 12/13};
normal . normal == 1
(* True *)
Get Mathematica to generate 3 points on the circle
pts = Block[{v = {x, y, z}},
FindInstance[
normal . v == normal . centre && (v - centre) . (v - centre) ==
radius^2, v, Reals, 3]
];
Thread /@ Thread[M1 -> ({x, y, z} /. pts)] // Flatten
(* {X[1] -> 33/26, Y[1] -> (4295 - 39 Sqrt[535])/1040,
Z[1] -> (5605 + 13 Sqrt[535])/1040, X[2] -> 119/26,
Y[2] -> (4037 - 39 Sqrt[879])/1040,
Z[2] -> (4831 + 13 Sqrt[879])/1040, X[3] -> 5/2,
Y[3] -> 1/80 (323 - 3 Sqrt[2391]), Z[3] -> 1/80 (409 + Sqrt[2391])} *)
Verify that this locates the solution correctly.
(solution1 /. %) // Simplify
(* {{x -> 3, y -> 4, z -> 5}} *)