Question: Fitting cylinder to data

I was hoping to fit a cylinder to some data points in 3d space and find the equation of the long axis of the cylinder.

Sample data (actual data much larger matrix):

data := Matrix(1..20,1..3,{(1,1)=4.7,(1,2)=13,(1,3)=-28.2,(2,1)=9.5,(2,2)=11.7,(2,3)=-28.1,
(3,1)=14.2,(3,2)=10.2,(3,3)=-27.9,(4,1)=-0.2,(4,2)=12.5,(4,3)=-27.6,
(5,1)=-5.2,(5,2)=12,(5,3)=-26.9,(6,1)=-10.1,(6,2)=11.5,(6,3)=-26.2,(7,1)=16,
(7,2)=6,(7,3)=-26.2,(8,1)=-14,(8,2)=8.9,(8,3)=-24.7,(9,1)=17.6,(9,2)=1.6,
(9,3)=-24.4,(10,1)=-17.9,(10,2)=6.1,(10,3)=-23.1,(11,1)=15.7,(11,2)=-2.5,(11,3)=-22.4,
(12,1)=-0.9,(12,2)=13.1,(12,3)=-22.3,(13,1)=3.8,(13,2)=12.2,(13,3)=-22.3,(14,1)=8.3,
(14,2)=10.1,(14,3)=-21.8,(15,1)=-5.9,(15,2)=12.7,(15,3)=-21.6,(16,1)=-19.4,(16,2)=2.3,
(16,3)=-21.3,(17,1)=-10.8,(17,2)=12.3,(17,3)=-21,(18,1)=11.8,(18,2)=7.2,(18,3)=-20.8,
(19,1)=13.5,(19,2)=-6.5,(19,3)=-20.4,(20,1)=-15.1,(20,2)=10.1,(20,3)=-19.6});

I've managed to fit a cylinder to the data if I assume that the cylinder's axis is the z axis>>

n := LinearAlgebra[RowDimension](data);

sumSquareError := (r, X, Y) -> add(
 (data[i, 1]-X-r*(data[i, 1]-X)/
 sqrt((data[i, 1]-X)^2+(data[i, 2]-Y)^2))^2+
 (data[i, 2]-Y-r*(data[i, 2]-Y)/
 sqrt((data[i, 1]-X)^2+(data[i, 2]-Y)^2))^2, i = 1 .. n);
 

Optimization[Minimize](sumSquareError(radius, XOffset, YOffset));

dataPlot := pointplot3d(da,symbol=box);
display(PSurfPlot,cylinder([XOffset,1.YOffset, -100], radius, 500), orientation=[45, 75], scaling=constrained);
 

But I'm having lots of trouble trying to find out how to fit the cylinder for an arbitary axis. I've found this on another web site:
 

Suppose the axis goes through the point (X0,Y0,Z0) and lies in the
direction of the unit vector w; further suppose that u and v are unit
vectors that are (1) mutually perpendicular and (2) are perpendicular
to the axis. Parametric equations for the cylindrical surface are

x = X0 + cos(theta)u1 + sin(theta)v1 + t*w1
y = Y0 + cos(theta)u2 + sin(theta)v2 + t*w2
z = Z0 + cos(theta)u3 + sin(theta)v3 + t*w3

where theta varies over the interval 0 to 2pi and t ranges over the
the set of real numbers.

But I still can't figure out how to do it (ie. how to find X0,Y0,Z0, the unit vector w and the radius of the cylinder).
 

Sorry this is my second question in this forum in a short time.

Any help would be much appreciated.