Fit 3d data to a torus

Question

I'd like to find the best fit coefficients "a" and "c" for a toroidal surface, to a list of {x,y,z} points which are approximately on the surface of the torus. I believe I need LinearModelFit for this, the function for the torus: z^2 == a^2 - (c - (x^2 + y^2)^(1/2))^2 and the data:
data={{x1,y1,z1},{x2,y2,z2},{x3,y3,z3},{x4,y4,z4},{x5,y5,z5}};   (* these would be passed as actual xyz values not variables. My fit function is more complicated than below so the actual values I have are not appropriate *)

In which case the command should look something like this but I'm clearly not expressing the torus function as a model correctly for LinearModelFit.  Also I will need to include some constraints, but once I have the syntax and functional form for the model correct I should be able to manage the constraints and starting values:
LinearModelFit[data, z^2 == a^2 - (c - (x^2 + y^2)^(1/2))^2,{a,c,,x,y,z}]

This answer
How to fit a surface to 3D data in Mathematica? is helpful but I couldn't see how to express the model for a torus in an equivalent form.

Cesareo · Accepted Answer

Assuming that data is obtained from translated and/or rotated models

Clear[a, c]
p = {x, y, z};
P = {X, Y, Z};
Txyz = z^2 - a^2 + (c - (x^2 + y^2)^(1/2))^2 /. Thread[p -> RollPitchYawMatrix[{alpha, beta, gamma}].(P - {x0, y0, z0})];

now using the excellent script from @flinty for random data generation

SeedRandom[1];
Torus = ResourceFunction["Torus"];
testTorus = Torus[{4, -2, 6}, {19, 4}];
(*pts on a torus plus some noise*)

pts = RandomPoint[DiscretizeGraphics@testTorus, 300] + RandomVariate[NormalDistribution[0, .5], {300, 3}];
gr1 = Graphics3D[Point@pts, Axes -> True];

and then following with a minimization procedure

error = Sum[(Txyz /. Thread[P -> pts[[k]]])^2, {k, 1, Length[pts]}];
sol = NMinimize[{error, -Pi <= alpha <= Pi, -Pi <= beta <= Pi, -Pi <= gamma <= Pi}, {a, c, x0, y0, z0, alpha, beta, gamma}]

and the results

gr0 = ContourPlot3D[(Txyz /. sol[[2]]) == 0, {X, -20, 20}, {Y, -20, 20}, {Z, -20, 20}, ContourStyle -> {Yellow, Opacity[0.2]}, Mesh -> None, BoundaryStyle -> None];
Show[gr1, gr0]

flinty · Answer

Shape fitting is quite difficult and has robustness issues without taking special care, but you can get quite far without any fitting. If your torus data is not rotated and oriented directly up, then the BoundingRegion commands provide very good parameters:
SeedRandom[1];
Torus = ResourceFunction["Torus"];
testTorus = Torus[{4, -2, 6}, {19, 4}];
(* pts on a torus plus some noise *)
pts = RandomPoint[DiscretizeGraphics@testTorus, 300] + 
   RandomVariate[NormalDistribution[0, .5], {300, 3}];

minball = BoundingRegion[pts, "MinBall"];
minbox = BoundingRegion[pts, "MinCuboid"];
pos = minball[[1]];
radius = minball[[2]];
holeradius = Min[EuclideanDistance[pos, #] & /@ pts];
piperadius = Min[Abs[minbox[[1]] - minbox[[2]]]]/2;
Graphics3D[{Point@pts, Opacity[.5], 
  Torus[pos, {(radius + holeradius), piperadius}]}]

If you torus data is oriented, then you should use the KarhunenLoeveDecomposition to get it into a face-up and centered form first, then fit, then rotate and translate back. Let me know if this is desired too.

Fit 3d data to a torus

2 Answers

Add your own answers!

Ask a Question