(*---------------------------start of LiveGraphics3D.m-------------------*) (*:Mathematica Version: 4.0 *) (*:Package Version: 1.15 *) (*:Title: LiveGraphics3D Toolbox *) (* :Context: LiveGraphics3D` *) (* :Author: Martin Kraus (mk) (Martin.Kraus@informatik.uni-stuttgart.de) some routines modified by Eric W. Weisstein (eww) (eww@wolfram.com) *) (* :Summary: This package implements several functions for three-dimensional graphics, which are especially useful in preparing graphics for the LiveGraphics3D applet. Master version available from: http://mathworld.wolfram.com/notebooks/LiveGraphics3D.m *) (*:History: v1.00 (1998-11-08): Written by Martin Kraus (mk) v1.01 (1998-11-11): Dir added (eww) v1.02 (1998-11-12): Minor mods to WriteLiveAnimationForm and reduced number of sides for circle and disk (eww) v1.03 (1998-11-15): improved NonIntersectingPolygons, reprogrammed DropInnerPolygons, removed old version of DropInnerPolygons and PlotCentersOfMasses (mk) v1.04 (1998-11-16): corrected ClipAtPlotRange (did not work in v1.03) (mk) v1.05 (1998-11-16): documentation updated slightly, tweaked circle, disk (eww) v1.06 (1998-11-25): new function UnifyPolygons and helpers (MergePolygons, PolygonInSamePlane) (mk) v1.07 (1998-11-28): NonIntersectingLines (mk) v1.08 (1998-12-05): documentation updated slightly (eww) v1.09 (1998-12-11): PolyhedronPolygon (eww+mk) v1.10 (1999-01-08): LinelessFaces (eww+mk) v1.11 (1999-01-17): OverlappingEdgesQ now more tolerant (mk) v1.12 (1999-03-30): Improved DropInnerPolygons to make the front face of the starting polygon facing away from the origin. (mk) v1.13 (1999-07-17): Added BoundingBox[{}]={} to stop occasional errors in working with duals of uniform polyhedra. At some point, I should figure out why this is happening and actually fix it. (eww) v1.14 (1999-08-18): Added workaround for writing only 5 significant numbers to files since NumberForm[InputForm[N[g]],5] does not work properly under Mathematica 4.0. v1.15 (2000-01-20): Updated emails and URL (eww) (c) 1998-2000 Martin Kraus *) (*:Keywords: Graphics3D, LiveGraphics3D *) (*:Sources: None. *) (* :Limitations: Several functions work only with flat, convex polygons. However, PolygonTriangulation routines can be used to break up concave polygons. Several functions deal only with polygons while they should also deal with lines. Many functions have a certain tolerance for numerical deviations; however, this tolerance assumes a certain scale (coordinates not much greater than 10 and not much smaller than 1/10). The output file consists of a single long line with no newlines, and therefore is difficult to view in Mathematica and some text editors. Nothing has been seriously tested. *) (* :Discussion: This package is still under construction. Martin Kraus's polygon triangulation routines available from http://library.wolfram.com/packages/polygontriangulation/ may be used in addition to or instead of some of the routines included here. Please send bug reports, comments and suggestions to Martin.Kraus@informatik.uni-stuttgart.de *) (*:References: *) BeginPackage["LiveGraphics3D`", {"Graphics`Animation`","Graphics`Graphics3D`","Graphics`Shapes`"}] CenterAtCenterOfMass::usage = "CenterAtCenterOfMass[graphics3d] resets the plot range such that the center of mass is in the center of the bounding box. (All polygons are assumed to have equal weights.)" ClipAtPlotRange::usage = "ClipAtPlotRange[graphics3d] clips all polygons at the bounding box."; Dir::usage = "Dir->directory is an option to WriteLiveForm and WriteLiveAnimationForm which specified in which directory the output file should be placed." DropInnerPolygons::usage = "DropInnerPolygons[graphics3d] tries to drop all polygons which are not part of the surface of the one solid in graphics3d. It also reverses the order of points in some polygons such that the same face (BACK or FRONT) will be on the surface of the solid. Therefore, DropInnerPolygons is also useful for solids without any inner polygons (like eww's obelisk.m example). The dropping of inner polygons should reduce the size of certain large graphics dramatically. DropInnerPolygons must often be called using the syntax DropInnerPolygons[NonIntersectingPolygons[graphics3D]], and often best results are obtained with applet parameters VISIBLE_FACES=BACK or VISIBLE_FACES=FRONT. DropInnerPolygons returns some progress information as it is running." Graphics2DTo3D::usage = "Graphics2DTo3D[graphics2D] gives a flat 3D object from a 2D object. Now that Mathematica 4.0 supports animiated GIF exporting, this function isn't needed so much." LinelessFaces::usage = "LinelessFaces[graphics3D] returns the solid with only exterior edges included." LiveAnimationForm::usage = "LiveAnimationForm[list[graphics3d]] converts the list of graphics3d objects to a form appropriate for an animation to be displayed by the LiveGraphics3D applet." LiveForm::usage = "LiveForm[graphics3d] converts g to a form appropriate for a rotatable 3-D solid to be displayed by the LiveGraphics3D applet." NonIntersectingPolygons::usage = "NonIntersectingPolygons[graphics3D] breaks up all intersecting polygons in such a manner that the resulting polygons do not intersect each other. NonIntersectingPolygons can take a long time to run (i.e., hours) for complicated solids such as complex uniform polyhedra." NonIntersectingLines::usage = "NonIntersectingLines[graphics3D] breaks up all lines which intersect polygons." PlotCentersOfMasses::usage = "PlotCentersOfMasses[graphics3d, p] shows a plot of points which shows the distance of the center of masses of all polygons to p. (The y coordinate of the points is random.)" PolyhedronPolygon::usage = "PolyhedronPolygon[polyhedron,polygon list] returns a Graphics3D object representing polyhedron with a line on the surface. The polygon list is a pure list of vertices of the intersecting plane, which is assumed to contain and be symmetric with respect to the origin." SplitLines::usage = "SplitLines[graphics3d,n] splits all lines into smaller lines. Each original straight line is broken into n pieces." Triangulate::usage = "Triangulate[g,n] takes a Graphics3D polyhedron and triangulates all faces with n sides." UnifyPolygons::usage = "UnifyPolygons[graphics3d] unifies neighbouring polygons which are in the same plane. Example: DropInnerPolygons[NonIntersectingPolygons[UnifyPolygons[Graphics3D[NumberedPolyhedron[21]]]]]." WriteLiveAnimationForm::usage = "WriteLiveAnimationForm[file,g,] takes a list of Graphics3D objects and writes an animation to file which can be displayed by the LiveGraphics3D applet. Any 2D graphics must first be converted using Graphics2DTo3D." WriteLiveForm::usage = "WriteLiveForm[file,g,] takes a Graphics3D object or list of such objects and writes an animation to file which can be displayed by the LiveGraphics3D applet. Any 2D graphics must first be converted using Graphics2DTo3D." (************* set your default directory for output files here *************) Options[LiveForm]={ Dir->"" }; Begin["`Private`"] addn[n_,k_,m_]:=Mod[n+k-1,m]+1; (* Under certain circumstances, e.g., NonIntersectingPolygons[Graphics3D[Dual[NumberedPolyhedron[18]]]]][[1]]], a couple BoundingBox[{}]'s seem to get thrown in *) BoundingBox[{}]:={} BoundingBox[Polygon[l_List]]:= Module[{xs,ys,zs},{xs,ys,zs}= Transpose[l];{{Min[xs],Max[xs]},{Min[ys],Max[ys]},{Min[zs],Max[zs]}}]; BoundingBox[Line[l_List]]:= Module[{xs,ys,zs},{xs,ys,zs}= Transpose[l];{{Min[xs],Max[xs]},{Min[ys],Max[ys]},{Min[zs],Max[zs]}}]; BreakLineWithPolygon[{Line[l1_List],bb1_List},{p_Polygon,bb_List}]:= Module[{i,breakpoint,broken={},part,rest}, i=1; rest=l1; While[i BreakLineWithPolygon[{Line[l],box},pl[[i]]]}; i=i+1];Flatten[broken/.{{Line[l_],_List}:>Line[l]}]]; BreakPolygonWithPolygon[{Polygon[l1_List],bb1_List},{Polygon[l2_List], bb2_List}]:=Module[{line,s1,s2,p1,p2,i1,i2,h,fulll1,result}, If[!OverlappingBoxesQ[bb1,bb2],Return[{{Polygon[l1],bb1}}]]; fulll1=Flatten[{l1,l1},1]; line=LineOfIntersection[Polygon[l1],Polygon[l2]]; If[line=={},Return[{{Polygon[l1],bb1}}]]; s1=ParametersOfClosestPoints[line,Polygon[l1]]; s2=ParametersOfClosestPoints[line,Polygon[l2]]; If[s1=={}||s2=={},Return[{{Polygon[l1],bb1}}]]; If[s1[[2,1]]i2,h=i1;i1=i2;i2=h;h=p1;p1=p2;p2=h]; result={DropDoubledPoints[Polygon[Flatten[ {{p1},Take[fulll1,{i1+1,i2}],{p2}} ,1]]],DropDoubledPoints[Polygon[Flatten[ {{p2},Take[fulll1,{i2+1,i1+Length[l1]}],{p1}} ,1]]]}; {{result[[1]],BoundingBox[result[[1]]]}, {result[[2]],BoundingBox[result[[2]]]}}]; BreakPolygonWithPolygons[{p_Polygon,bb_List},pl_List]:= Module[{broken={{p,bb}},i=1,l}, While[i<=Length[pl], broken=broken/.{{Polygon[l_],box_List}:> BreakPolygonWithPolygon[{Polygon[l],box},pl[[i]]]}; i=i+1];Flatten[broken/.{{Polygon[l_],_List}:>Polygon[l]}]]; CenterOfMass[Polygon[l_]]:=Plus@@N[l]/Length[l]; CenterOfMass[g_Graphics3D]:= Module[{com,p}, com=(Cases[g[[1]],_Polygon,Infinity]/.p_Polygon:>CenterOfMass[p]); Plus@@com/Length[com]]; CenterAtCenterOfMass[g_Graphics3D,opts___]:= Module[{com,pr,oldcenter,p,i}, com=(Cases[g[[1]],_Polygon,Infinity]/.p_Polygon:>CenterOfMass[p]); com=Plus@@com/Length[com];pr=FullOptions[g,PlotRange]; oldcenter=Apply[Plus,pr,2]/2.;i=1; While[i<=3, If[oldcenter[[i]]pr] ] centroid[l_List]:=(Plus@@l)/Length[l]; circle[{x_,y_},{a_,b_}]:=Module[{i,n=25}, Line[Table[{x+a Cos[2.Pi i/n],y+b Sin[2.Pi i/n]},{i,0,n}]] ] circle[{x_,y_},r_]:=circle[{x,y},{r,r}] ClipAtPlotRange[SurfaceGraphics[gl_,opts___]]:= ClipAtPlotRange[Graphics3D[SurfaceGraphics[gl,opts]]] ClipAtPlotRange[Graphics3D[gl_,opts___]]:=Module[ {g,l,xmin,xmax,ymin,ymax,zmin,zmax,polys,pr,polysbbs}, Off[Solve::svars]; Off[RowReduce::luc]; pr={{xmin,xmax},{ymin,ymax},{zmin,zmax}}= FullOptions[Graphics3D[gl,opts],PlotRange]; polys={Polygon[{{xmin,ymin,zmin},{xmin,ymin,zmax},{xmin,ymax,zmax},{xmin, ymax,zmin}}], Polygon[{{xmin,ymin,zmin},{xmin,ymin,zmax},{xmax,ymin,zmax},{xmax, ymin,zmin}}], Polygon[{{xmin,ymin,zmin},{xmin,ymax,zmin},{xmax,ymax,zmin},{xmax, ymin,zmin}}], Polygon[{{xmax,ymin,zmin},{xmax,ymin,zmax},{xmax,ymax,zmax},{xmax, ymax,zmin}}], Polygon[{{xmin,ymax,zmin},{xmin,ymax,zmax},{xmax,ymax,zmax},{xmax, ymax,zmin}}], Polygon[{{xmin,ymin,zmax},{xmin,ymax,zmax},{xmax,ymax,zmax},{xmax, ymin,zmax}}]}; polysbbs=polys/.{p_Polygon:>{p,BoundingBox[p]}}; g=gl/.{p_Polygon:>BreakPolygonWithPolygons[{p,BoundingBox[p]}, polysbbs]}; grownpr=Transpose[{ Transpose[pr][[1]]-0.001*(Transpose[pr][[2]]-Transpose[pr][[1]]), Transpose[pr][[2]]+0.001*(Transpose[pr][[2]]-Transpose[pr][[1]])}]; g=g/.{p_Polygon:>If[IsPointOfPolygonOutside[p,grownpr],{},p]}; On[Solve::svars]; On[RowReduce::luc]; Graphics3D[g,opts]]; CosAngleBetweenPolygons[normal1_,sign1_,dir1_,normal2_,dir2_]:= Module[{normal2p,into1,cosine}, If[dir1.dir2<0., normal2p=normal2, normal2p=-normal2]; into1=Cross[normal1,dir1]; cosine=N[Re[ArcCos[-normal1.normal2p]]]; If[normal2p.into1<0.,cosine=2.*Pi-cosine]; If[sign1<0,cosine=2.*Pi-cosine]; cosine]; disk[{x_,y_},{a_,b_}]:=Module[{i,n=25}, {EdgeForm[],Polygon[Table[{x+a Cos[2.Pi i/n],y+b Sin[2.Pi i/n]},{i,n}]]} ] disk[{x_,y_},r_]:=disk[{x,y},{r,r}] Distance[p1_List,p2_List]:=N[Sqrt[(p1-p2).(p1-p2)]]; DropDoubledPoints[Polygon[l_]]:= Module[{maxdist,newl={Last[l]},i=1}, maxdist=10.^-8* Max[Map[((#[[1]]-#[[2]]).(#[[1]]-#[[2]]))&, Transpose[{l,RotateRight[l]}]]]; While[i<=Length[l], If[(l[[i]]-Last[newl]).(l[[i]]-Last[newl])>maxdist, newl=Append[newl,l[[i]]]];i=i+1]; newl=Drop[newl,1]; If[Length[newl]>2,Polygon[newl],{}]]; DropInnerPolygons[g_Graphics3D]:= Module[{startpolys,xcoors,minx,current,polys,normals,edges,signs,mainpolygon, mainsign,surfaces,checkedsurfaces,i,gp},polys=Cases[g,_Polygon,Infinity]; normals=Map[NormalVector,polys];edges=Map[Edges,polys]; boxes=Map[BoundingBox,polys]; signs=Table[1,{Length[polys]}]; checkedsurfaces={}; Print["There are a total of ",Length[polys]," polygons."]; minx=Min[Transpose[Transpose[boxes][[1]]][[1]]]+10.^-4; i=1;startpolys={}; While[i2, startpolys= Append[startpolys,{Sort[xcoors][[2]],-Abs[normals[[i,1]]],i}]]]; i=i+1]; surfaces={Sort[startpolys][[1,3]]}; i=0; While[Length[surfaces]>0, current=First[surfaces]; i=i+1; If[i<4||Mod[i,10]==0, Print["processing polygon no.",current, " (",Length[checkedsurfaces], " processed, ",Length[surfaces]," in buffer)"], If[i==4,Print["Reporting only every 10th."]]]; {newsurfs,signs}= NeighbouringSurfaces[current,polys,normals,boxes,edges,signs,surfaces, checkedsurfaces]; checkedsurfaces=Append[checkedsurfaces,First[surfaces]]; surfaces=Join[Drop[surfaces,1],newsurfs]; (*Show[Graphics3D[{Extract[polys,Partition[checkedsurfaces,1]]}], Axes->True];*) ]; Print["Returning ",Length[checkedsurfaces]," of ",Length[polys], " polygons."]; mainpolygon=polys[[First[checkedsurfaces]]]; mainsign=Sign[CenterOfMass[mainpolygon].NormalVector[mainpolygon]]; i=1;gp=g; While[i<=Length[polys], pos=Position[gp,polys[[i]]]; If[!MemberQ[checkedsurfaces,i], gp=ReplacePart[gp,{},pos], If[signs[[i]]*mainsign<0, gp=ReplacePart[gp,Polygon[Reverse[polys[[i,1]]]],pos]]]; i=i+1]; gp]; Edges[p_Polygon]:= Transpose[{Transpose[{p[[1]],RotateLeft[p[[1]]]}], Map[Normalize,RotateLeft[p[[1]]]-p[[1]]]}]; Graphics2DTo3D[g_Graphics,r_:1,opts___]:=Module[{}, Off[StackGraphics::bdpts]; Show[StackGraphics[{g} /. {Circle->circle,Disk->disk}], opts,ViewPoint->{0,-3,0},Boxed->False,Axes->False,Ticks->None, BoxRatios->{1,1,r},PlotRange->All] ] Graphics2DTo3D[g_List,r_:1,opts___]:=Module[{i}, Table[Graphics2DTo3D[g[[i]],r,opts],{i,Length[g]}] ] Graphics2DTo3D[gl_GraphicsArray,r_:1,opts___]:=Module[{i,g=Flatten[gl[[1]]]}, Table[Graphics2DTo3D[g[[i]],r,opts],{i,Length[g]}] ] IsPointOutside[{x_,y_,z_},{{xmin_,xmax_},{ymin_,ymax_},{zmin_,zmax_}}]:= If[xxmax||yymax||zzmax,True,False]; IsPointOfPolygonOutside[Polygon[l_List],plotrange_List]:= Or@@Map[IsPointOutside[#,plotrange]&,l]; line[l_List]:=Line[Join[l,{l[[1]]}]] LinelessFaces[g_]:=Module[ {g1=Show[g,DisplayFunction->Identity], g2,g3,g4 }, Print["Step 1 of 4. Splitting intersecting and dropping inner polygons..."]; g2=DropInnerPolygons[NonIntersectingPolygons[g]]; Print["Step 2 of 4. Re-unifying polygons..."]; g3=UnifyPolygons[g2]; Print["Step 3 of 4. Determining Edges..."]; g4=Show[SplitLines[g3/.Polygon[l_]:>Line[Append[l,First[l]]],1], DisplayFunction->Identity]; Print["Step 4 of 4. Combining faces and edges..."]; Graphics3D[{EdgeForm[],g2[[1]],g4[[1]]}] ] LineOfIntersection[Polygon[l1_List],Polygon[l2_List]]:= Module[{n1=NormalVector[Polygon[l1]],n2=NormalVector[Polygon[l2]],n,p1,p2, p3},n=Cross[n1,n2];If[n.n<10.^-5,Return[{}]];n=n/Sqrt[n.n]; sol=Solve[{n1.({p1,p2,p3}-l1[[1]])==0,n2.({p1,p2,p3}-l2[[1]])==0, n.({p1,p2,p3}-l1[[1]])==0},{p1,p2,p3}]; If[Length[sol]=!=1||Length[sol[[1]]]=!=3, Return[{}]];{{p1,p2,p3}/.sol[[1]],n}]; LiveAnimationForm[gl_]:= HoldForm[ShowAnimation][LiveForm[gl]] LiveAnimationForm[gl_,opts___]:= HoldForm[ShowAnimation][LiveForm[gl],LiveForm[opts]] LiveForm[gr_]:=Module[ {g=Switch[gr, _List,gr, _Graphics3D,gr, _ContourGraphics,Graphics3D[SurfaceGraphics[gr]], _DensityGraphics,Graphics3D[SurfaceGraphics[gr]], _SurfaceGraphics,Graphics3D[gr], _,gr ]}, NumberForm[InputForm[N[g]],5] ] MergePolygons[p1_Polygon,e1_Integer,p2_Polygon,e2_Integer]:= Module[{newpol,dif,i,d,maxdist=10.^-8}, newpol = DropDoubledPoints[ Polygon[Flatten[{Take[p1[[1]],{1,e1}], Take[p2[[1]],{e2+1,Length[p2[[1]]]}],Take[p2[[1]],{1,e2}], Take[p1[[1]],{e1+1,Length[p1[[1]]]}]},1]]]; i=1; While[i <= Length[newpol[[1]]], dif=newpol[[1,Mod[i-1-1,Length[newpol[[1]]]]+1]]- newpol[[1,Mod[i+1-1,Length[newpol[[1]]]]+1]]; If[dif.dif0&&!MemberQ[surfaces,oedges[[1,2]]]&&! MemberQ[checkedsurfaces,oedges[[1,2]]], result=Append[result,oedges[[1,2]]]; newsigns[[oedges[[1,2]]]]=-signs[[pindex]]*Sign[oedges[[1,3]]]]; i=i+1]; {result,newsigns}]; NonIntersectingPolygons[Graphics3D[gl_]]:= NonIntersectingPolygons[Graphics3D[gl,{}]] NonIntersectingPolygons[Graphics3D[gl_,opts___]]:=Module[{g,l,polys}, Off[Solve::svars]; Off[RowReduce::luc]; polys=Cases[gl,_Polygon,Infinity]; polys=polys/.Polygon[l_]:>{Polygon[l],BoundingBox[Polygon[l]]}; g=Graphics3D[gl/.{Polygon[l_]:>BreakPolygonWithPolygons[ {Polygon[l],BoundingBox[Polygon[l]]},polys]}, opts]; On[Solve::svars]; On[RowReduce::luc]; g ]; NonIntersectingLines[Graphics3D[gl_,opts___]]:=Module[{g,l,polys}, Off[Solve::svars]; Off[RowReduce::luc]; polys=Cases[gl,_Polygon,Infinity]; polys=polys/.Polygon[l_]:>{Polygon[l],BoundingBox[Polygon[l]]}; g=Graphics3D[gl/.{Line[l_]:> BreakLineWithPolygons[{Line[l],BoundingBox[Line[l]]},polys]}, opts]; On[Solve::svars]; On[RowReduce::luc]; g ]; Normalize[v_List]:= Module[{length=N[Sqrt[v.v]]},If[length<10.^-20,v,v/length]]; NormalVector[Polygon[l_List]]:= Module[{i=1,k=1,v}, While[i10.^-5,v[k]=l[[i+1]]-l[[i]];k=k+1];i=i+1]; If[k!=3,{1,1,1},Normalize[Cross[v[1],v[2]]]]]; OverlappingBoxesQ[{{xmin1_,xmax1_},{ymin1_,ymax1_},{zmin1_,zmax1_}},{{xmin2_, xmax2_},{ymin2_,ymax2_},{zmin2_,zmax2_}},tol_:0]:= If[(xmin1>xmax2+tol||xmax1+tol < xmin2)||( ymin1>ymax2+tol||ymax1+tol < ymin2)||( zmin1>zmax2+tol||zmax1+tol < zmin2),False,True]; OverlappingEdgesPositions[edges_,{{p1a_,p1b_},dir1_}]:= Position[edges,{{p2a_List,p2b_List},dir2_List}/; OverlappingEdgesQ[{{p1a,p1b},dir1},{{p2a,p2b},dir2}]]; OverlappingEdgesQ[{{p1a_,p1b_},dir1_},{{p2a_,p2b_},dir2_}]:= Module[{dot}, If[Plus@@Abs[dir1-dir2]<10.^-4&& Plus@@Abs[(p2a-p1a)-((p2a-p1a).dir1)*dir1]<10.^-4, dot=(p1b-p1a).(p2a-p1a); If[dot>0, If[dot >(p1b-p1a).(p1b-p1a)-10.^-4, False, True], dot=(p1b-p1a).(p2b-p1a); If[dot>10.^-4, True, False]], If[Plus@@Abs[dir1+dir2]<10^-4&& Plus@@Abs[(p2a-p1a)-((p2a-p1a).dir1)*dir1]<10.^-4, dot=(p1b-p1a).(p2b-p1a); If[dot>0, If[dot>(p1b-p1a).(p1b-p1a)-10.^-4, False, True], dot=(p1b-p1a).(p2a-p1a); If[dot >10.^-4, True, False]], False]]]; ParameterOfClosestPoint[{p1_List,v1_List},{p2_List,v2_List}]:= Module[{},-(Cross[v1,v2].Cross[p1-p2,v2])/(Cross[v1,v2].Cross[v1,v2])]; ParameterOfClosestPointToPoint[{p1_List,v1_List},p_List]:=(p-p1).v1/Sqrt[v1.v1]; ParametersOfClosestPoints[{p_List,v_List},Polygon[l_List]]:= Module[{n=NormalVector[Polygon[l]],i,dir,fulll=Flatten[{l,l},1],oldsign, newsign,index,parameter},dir=Cross[v,n];i=2;k=1; oldsign=TolerantSign[(fulll[[i]]-p).dir]; While[i<=Length[fulll]&&k<3,newsign=TolerantSign[(fulll[[i]]-p).dir]; If[newsign!=0&&oldsign!=0&&newsign!=oldsign,index[k]=i-1; parameter[k]= ParameterOfClosestPoint[{p,v},{fulll[[i-1]], fulll[[i]]-fulll[[i-1]]}];k=k+1];If[newsign!=0,oldsign=newsign]; i=i+1];If[k!=3,{}, Sort[{{parameter[1],index[1]},{parameter[2],index[2]}}]]]; PlotCentersOfMasses[g_Graphics3D,origin_]:= Module[{pl,p}, pl=Cases[g[[1]],_Polygon,Infinity]/.p_Polygon:>CenterOfMass[p]; pl=Map[(Distance[origin,#])&,pl]; ListPlot[Transpose[{pl,Table[Random[],{Length[pl]}]}],Axes->{True,False}]]; PointInBoundingBoxQ[{x_,y_,z_}, {{xmin_,xmax_},{ymin_,ymax_},{zmin_,zmax_}},tol_:0]:= If[(xmin-tol<=x&&x<=xmax+tol)&&(ymin-tol<=y&&y<=ymax+tol)&&(zmin-tol<=z&&z<=zmax+tol),True,False]; PointOfIntersection[l1_,l2_,{poly_Polygon,bbox_}]:= Module[{l,p=poly[[1,1]],normal=NormalVector[poly],denom}, denom=((l2-l1).normal); If[Abs[denom]<0.000001,Return[{}]]; l=((p-l1).normal)/denom; point=l1+l*(l2-l1); If[l<=0.000001||l>=0.999999||!PointInBoundingBoxQ[point,bbox],{},point] ]; PolygonInSamePlane[pindex_,polys_,normals_,boxes_,edges_,checkedsurfaces_]:= Module[{result={},opolys,oedges,i}, opolys=Position[boxes,box_/;OverlappingBoxesQ[boxes[[pindex]],box,10.^-4], 1]; opolys=Intersection[opolys, Position[normals, normal_/;Abs[Abs[normal.normals[[pindex]]]-1]<10.^-4,1]]; If[Length[opolys]==0,Return[{}]]; i=1; oedges={}; While[i<=Length[edges[[pindex]]]&&Length[oedges]==0, currentedge=edges[[pindex,i]]; oedges= Map[({#[[1]], OverlappingEdgesPositions[edges[[#[[1]]]],currentedge]})&, opolys]; (* oedges has now elements of the form {index of polygon,{positions of overlapping edges}}*) oedges= Map[( If[#[[2]]=={}||#[[1]]==pindex|| MemberQ[checkedsurfaces,#[[1]]],{},{{ i,#[[1]],#[[2,1,1]]}}])&, oedges]; oedges=Flatten[oedges,1]; i=i+1]; If[Length[oedges]>0,oedges[[1]],{}]]; PolyhedronPolygon[solid_,polygon_List]:=Module[ {g=Show[Graphics3D[{solid,Polygon[polygon]}],DisplayFunction->$DisplayFunction]}, Graphics3D[ NonIntersectingPolygons[Show[{DropInnerPolygons[Graphics3D[g[[1,1]]]], AffineShape[Graphics3D[g[[1,2]]],{1.1,1.1,1.1}]}]][[1,1]] ] ] SplitPoly[p_]:=Module[{i,c=centroid[p],n=Length[p],f=1}, {Table[Polygon[{c+f(p[[i]]-c), c+f(p[[addn[i,1,n]]]-c),c}],{i,n}] (*,line[p] *)} ] SplitLines[g_,n_Integer]:= Module[{gl}, gl=g/.Line[l_List]:>Map[Line[#]&,Transpose[{Drop[l,-1],Drop[l,1]}]]; If[n>1,gl=gl/.Line[{a_,b_}]:>Line[Table[N[a+i*(b-a)],{i,0,1,1/n}]]; gl=gl/.Line[l_List]:>Map[Line[#]&,Transpose[{Drop[l,-1],Drop[l,1]}]];]; gl]; TolerantSign[x_]:=If[Abs[x]<10.^-4,0,Sign[x]]; Triangulate[g_,s_]:=Module[{poly=g[[1]] /.Polygon->Identity,base,sides}, {base=Select[poly,Length[#]==s&],sides=Select[poly,Length[#]!=s&]}; Graphics3D[{{EdgeForm[],SplitPoly/@base},Polygon/@sides}] ] UnifyPolygons[g_Graphics3D]:= Module[{current,polys,oldpolys,normals,edges,signs,surfaces,checkedsurfaces, mergedsurfaces,i,gp,mergethis},polys=Cases[g,_Polygon,Infinity]; oldpolys=polys;normals=Map[NormalVector,polys];edges=Map[Edges,polys]; boxes=Map[BoundingBox,polys]; checkedsurfaces={}; mergedsurfaces={}; Print["There are a total of ",Length[polys]," polygons."]; surfaces=Table[i,{i,Length[polys]}]; i=1; While[Length[surfaces]>0, current=First[surfaces]; (*If[i<4||Mod[i,10]==0, Print["processing polygon no.",current, " (",Length[checkedsurfaces], " processed, ",Length[surfaces]," in buffer)"], If[i==4,Print["Reporting only every 10th."]]];*) If[MemberQ[checkedsurfaces,current],mergethis={},mergethis= PolygonInSamePlane[current,polys,normals,boxes,edges, checkedsurfaces]]; If[Length[mergethis]==0, checkedsurfaces=Append[checkedsurfaces,current]; surfaces=Drop[surfaces,1]; i=i+1, {edgeindex,pindex,pedgeindex}=mergethis; Print["Merging polygon no.",current," (edge ",edgeindex,") with no.", pindex, " (edge ",pedgeindex,")"]; polys[[current]]= MergePolygons[polys[[current]],edgeindex,polys[[pindex]], pedgeindex]; edges[[current]]=Edges[polys[[current]]]; boxes[[current]]=BoundingBox[polys[[current]]]; checkedsurfaces=Append[checkedsurfaces,pindex]; mergedsurfaces=Append[mergedsurfaces,pindex]] (*Show[Graphics3D[{Extract[polys,Partition[checkedsurfaces,1]]}], Axes->True];*) ]; i=1;gp=g; While[i<=Length[polys], pos=Position[gp,oldpolys[[i]]]; If[MemberQ[mergedsurfaces,i], gp=Delete[gp,pos], gp=ReplacePart[gp,polys[[i]],pos]]; i=i+1]; gp]; WriteLiveAnimationForm[filename_,gl_,opts___]:=Module[ {outfile,dir=Dir/.{opts}/.Options[LiveForm]}, Unprotect[Real]; Format[x_Real, InputForm] := OutputForm[NumberForm[x, 5, NumberFormat -> (If[#3 == "", #1, SequenceForm[#1, "*^", #3]] &)]]; Protect[Real]; outfile=dir<>filename; WriteString[outfile, ToString[HoldForm[ShowAnimation][LiveForm[gl],LiveForm[opts]]]]; Close[outfile]; Unprotect[Real]; Format[x_Real, InputForm] =.; Protect[Real]; ] WriteLiveAnimationForm[filename_,gl_]:=Module[ {outfile,dir=Dir /. Options[LiveForm]}, Unprotect[Real]; Format[x_Real, InputForm] := OutputForm[NumberForm[x, 5, NumberFormat -> (If[#3 == "", #1, SequenceForm[#1, "*^", #3]] &)]]; Protect[Real]; outfile=dir<>filename; WriteString[outfile, ToString[HoldForm[ShowAnimation][LiveForm[gl]]]]; Close[outfile]; Unprotect[Real]; Format[x_Real, InputForm] =.; Protect[Real]; ] WriteLiveForm[filename_,g_,opts___]:=Module[ {outfile,dir=Dir/.{opts}/.Options[LiveForm]}, Unprotect[Real]; Format[x_Real, InputForm] := OutputForm[NumberForm[x, 5, NumberFormat -> (If[#3 == "", #1, SequenceForm[#1, "*^", #3]] &)]]; Protect[Real]; outfile=dir<>filename; WriteString[outfile,ToString[LiveForm[g]]]; Close[outfile]; Unprotect[Real]; Format[x_Real, InputForm] =.; Protect[Real]; ] End[] (* Protect[ ] *) EndPackage[] (*---------------------end of LiveGraphics3D.m-----------------------*)