nuvola=ReadList["/Users/falco/Desktop/mura aureliane/prova rid2.txt",Number,RecordLists->True];

nuvola[[1]]

nuvola1=Map[Drop[#,-6]&,nuvola];

Length[nuvola1]

maxr=Map[Max,Transpose[nuvola1]]

minr=Map[Min,Transpose[nuvola1]]

Graphics3D[Table[{PointSize[.002],RGBColor[Take[nuvola[[k]],{4,6}]],Point[Take[nuvola[[k]],3]]},{k,1,Length[nuvola]}],Axes->True,ViewPoint->Right]

Show[{Graphics3D[{PointSize[.0005],Point[nuvola1]}]},ImageSize->800,Axes->True]

diste[{y1_,z1_}]:=((Abs[b *y+c* z+ d])/(Sqrt[(b^2)+(c^2)]))

dist[{x_,y_,z_},{a_,b_,c_,d_}]:=((Abs[a*x+b *y+c* z+ d])/(Sqrt[(a^2)+(b^2)+(c^2)]))

dist2div[{x_,y_,z_},{a_,b_,c_,d_}]:=Expand[(a*x+b *y+c* z+ d)^2]

fun=Total[Map[dist2div[#,{a,b,c,d}]&,nuvola1]]/((a^2)+(b^2)+(c^2))

{min,sost}=FindMinimum[fun,{a,-2},{b,1},{c,-1},{d,0},MaxIterations->500]

Show[{ContourPlot3D[Evaluate[a*x+b*y+c*z+d==0/.sost],{x,minr[[1]],maxr[[1]]},{y,minr[[2]],maxr[[2]]},{z,minr[[3]],maxr[[3]]},Mesh->None
]},ImageSize->300,Axes->True]

Show[{%41,%9},PlotRange->All]

Manipulate[Show[{Graphics3D[{PointSize[.001],Point[Select[nuvola1,dist[#,{a,b,c,d+k}/.sost]<.001&]]}]},ImageSize->800,Axes->True],{k,-1,1}]

Show[{Graphics3D[{PointSize[.0005],Point[Select[nuvola1,dist[#,{a,b,c,d}/.sost]<.1&]]}],ContourPlot3D[Evaluate[a*x+b*y+c*z+d==0/.sost],{x,minr[[1]],maxr[[1]]},{y,minr[[2]],maxr[[2]]},{z,minr[[3]],maxr[[3]]},Mesh->None,ContourStyle->Directive[Opacity[0.5]]
]},ImageSize->800,Axes->True]

mrt=RotationMatrix[{{a,b,c}/.sost,{1,0,0}}]

nuvola1rt=Map[mrt.#&,nuvola1];

Graphics[{PointSize[.0005],Point[Map[Drop[#,{3}]&,nuvola1rt]]},Axes->True]

nuvola1rtm=Select[nuvola1rt,-56.43<#[[1]]<-56.33&];

minrtm=Map[Min,Transpose[nuvola1rtm]]

maxrtm=Map[Max,Transpose[nuvola1rtm]]

Graphics[{PointSize[.0005],Point[Map[Drop[#,{3}]&,nuvola1rtm]]}]

-56.342, -56.36

Manipulate[Graphics[{PointSize[.0005],Point[Map[Drop[#,{1}]&,Select[nuvola1rt,h<#[[1]]<k&]]]},ImageSize->500],{{h,minrtm[[1]]},minrtm[[1]],maxrtm[[1]]},{{k,maxrtm[[1]]},minrtm[[1]],maxrtm[[1]]}]

GraphicsRow[{Graphics[{PointSize[.0005],Point[Map[Drop[#,{2}]&,nuvola1rtm]]}],Graphics[{PointSize[.0005],Point[Map[Drop[#,{1}]&,nuvola1rtm]]}],Graphics[{PointSize[.0005],Point[Map[Drop[#,{3}]&,nuvola1rtm]]}]},ImageSize->1000]

ListSurfacePlot3D[nuvola1rtm,MaxPlotPoints->50,MeshFunctions->(#1&)]

figmuro1m=Graphics[Table[{PointSize[.001],Hue[10(nuvola1rtm[[k,1]]+56.43)],Point[Drop[nuvola1rtm[[k]],{1}]]},{k,1,Length[nuvola1rtm]}],ImageSize->600]

Export["figmuro1.jpg",figmuro1,ImageResolution->300]

Graphics[Table[{PointSize[.0005],Hue[10(nuvola1rtm[[k,1]]+56.43)],Point[Drop[nuvola1rtm[[k]],{3}]]},{k,1,Length[nuvola1rtm]}]]

minrtm

maxrtm

Map[Mean,Transpose[nuvola1rtm]]

nuvola1rtmsez=Select[nuvola1rtm,-11.63<#[[2]]<-11.61&];

Graphics3D[{PointSize[.0005],Point[nuvola1rtmsez]}]

puntipiano=Map[Drop[#,{2}]&,nuvola1rtmsez];

Length[puntipiano]

FindFit[puntipiano]

ListCurvePathPlot[puntipiano]

Graphics[{PointSize[.0005],Point[Map[Drop[#,{2}]&,nuvola1rtmsez]]}]

sez2=Select[nuvola1rtm,-56.36<#[[1]]<-56.357&];

Graphics3D[{PointSize[.0005],Point[sez2]}]

5./6Pi-2

3./5(3+Sqrt[5])

N[Sqrt[10]]

minrt=Map[Min,Transpose[Map[mrt.#&,nuvola1]]]

maxrt=Map[Max,Transpose[Map[mrt.#&,nuvola1]]]

Show[{Graphics3D[{PointSize[.0005],Point[nuvola1]}]},ImageSize->800]

Timing[{m,rt}=FindMinimum[funcsz,{a,1},{b,3},{c,-0.4},{d,5},{r,0.5},MaxIterations->500]]

nuvola1trrt=Map[mrt.#&,nuvola1tr];

nuvola1sez=Select[nuvola1trrt,-.2<#[[3]]<.8&];

nuvola1sezrt=Map[mrt.#&,nuvola1sez];

Length[nuvola1sez]

Timing[coefsz=Map[Sum[Apply[Times,nuvola1sez[[k]]^#],{k,1,Length[nuvola1sez]}]&,expo]]

funcsz=(Expand[Apply[Plus,Thread[coefsz*par]]])/(1+a^2+c^2)^2;

Timing[{m2sz,rt2sz}=FindMinimum[funcsz,{a,1},{b,3},{c,-0.4},{d,5},{r,0.5},MaxIterations->500]]

Show[{Graphics3D[{PointSize[.0005],Point[nuvola1sez]}],Graphics3D[{Opacity[.5],Cylinder[{{x,a*x+b,c*x+d}
/.x->-1,{x,a*x+b,c*x+d}/.x->1},r]}/.rt2sz]},PlotRange->{{-1,1},{-1,1},{-1,1}},ImageSize->800]

mrtsz=RotationMatrix[{{1,a*1+b,c*1+d}/.rt2sz,{0,0,1}}]