Here is a set of structured 3D points. Now we can form a BSpline
using these points as knots.
dat=Import["3DFoil.mat", "Data"]
fu=BSplineFunction[dat]
Here we can do a ParametricPlot3D
with these points.
pic=ParametricPlot3D[fu[u,v],{u, 0, 1}, {v, 0, 1}, Mesh -> All, AspectRatio ->
Automatic,PlotPoints->10,Boxed-> False,Axes-> False]
Question
If we carefully look at the 3D geometry coming out of the spline we can see that it is a hollow structure. This hole appears in both side of the symmetric profile. How can we perfectly (not visually!) fill up this hole and create a unified Graphics3D
object where holes in both sides are patched.
What I am able to get so far is the following. Holes are not fully patched.
I am asking too many questions recently and I am sorry for that. But if any of you get interested I hope you will help.
Update
Here is the problem with belisarius method. It generates triangles with almost negligible areas.
dat = Import[NotebookDirectory[] <> "/3DFoil.mat", "Data"];
(*With your points in "dat"*)
fd = First@Dimensions@dat;
check = ParametricPlot3D[{BSplineFunction[dat][u, v],
BSplineFunction[{dat[[1]], Reverse@dat[[1]]}][u, v],
BSplineFunction[{dat[[fd]], Reverse@dat[[fd]]}][u, v]}, {u, 0,
1}, {v, 0, 1}, Mesh -> All, AspectRatio -> Automatic,
PlotPoints -> 10, Boxed -> False, Axes -> False]
output is here
Export[NotebookDirectory[]<>"myres.obj",check];
cd=Import[NotebookDirectory[]<>"myres.obj"];
middle=
check[[1]][[2]][[1]][[1(* Here are the numbers of different Graphics group*)]][[2,1,1,1]];
sidePatch1=check[[1]][[2]][[1]][[2]][[2,1,1,1]];
sidePatch2=check[[1]][[2]][[1]][[3]][[2,1,1,1]];
There are three Graphics
groups rest are empty. Now lets see the area of the triangles in those groups.
polygonArea[pts_List?
(Length[#]==3&)]:=Norm[Cross[pts[[2]]-pts[[1]],pts[[3]]-pts[[1]]]]/2
TriangleMaker[{a_,b_,c_}]:={vertices[[a]],vertices[[b]],vertices[[c]]}
tring=Map[polygonArea[TriangleMaker[#]]&,middle];
tring//Min
For the middle large group output is
0.000228007
This is therefore a permissible triangulation. But for the side patches we get zero areas.
Map[polygonArea[TriangleMaker[#]] &, sidePatch1] // Min
Map[polygonArea[TriangleMaker[#]] &, sidePatch2] // Min
Any way out here belisarius ?
My partial solution
First download the package for simplifying complex polygon from Wolfram archive.
fu = BSplineFunction[dat];
pic =(*ParametricPlot3D[fu[u,v],{u,0,1},{v,0,1},Mesh->None,
AspectRatio->Automatic,PlotPoints->25,Boxed->False,Axes->False,
BoundaryStyle->Red]*)
ParametricPlot3D[fu[u, v], {u, 0, 1}, {v, 0, 1}, Mesh -> None,
AspectRatio -> Automatic, PlotPoints -> 10, Boxed -> False,
Axes -> False, BoundaryStyle -> Black];
bound = First@Cases[Normal[pic], Line[pts_] :> pts, Infinity];
corners = Flatten[Table[fu[u, v], {u, 0, 1}, {v, 0, 1}], 1];
nf = Nearest[bound -> Automatic]; {a1, a2} =
Union@Flatten@(nf /@ corners);
sets = {bound[[2 ;; a1]], bound[[a1 ;; a2]],bound[[a2 ;; a2 + a1]]};
CorrectOneNodeNumber = Polygon[sets[[{1, 3}]]][[1]][[1]] // Length;
CorrectOneNodes1 =
Polygon[sets[[{1, 3}]]][[1]][[1]]; CorrectOneNodes2 =
Take[Polygon[sets[[{1, 3}]]][[1]][[2]], CorrectOneNodeNumber];
<< PolygonTriangulation`SimplePolygonTriangulation`
ver1 = CorrectOneNodes1;
ver2 = CorrectOneNodes2;
triang1 = SimplePolygonTriangulation3D[ver1];
triang2 = SimplePolygonTriangulation3D[ver2];
Show[Graphics3D[{PointSize[Large], Point[CorrectOneNodes1]},Boxed -> False,
BoxRatios -> 1], Graphics3D[{PointSize[Large], Point[CorrectOneNodes2]},
Boxed -> False, BoxRatios -> 1],
Graphics3D[GraphicsComplex[ver1, Polygon[triang1]], Boxed -> False,
BoxRatios -> 1],
Graphics3D[GraphicsComplex[ver2, Polygon[triang2]], Boxed -> False,
BoxRatios -> 1]]
We get nice triangles here.
picfin=ParametricPlot3D[fu[u,v],{u,0,1}, {v,0,1},Mesh->All,AspectRatio->Automatic,PlotPoints->10,Boxed->False,Axes->False,BoundaryStyle->None];pic3D=Show[Graphics3D[GraphicsComplex[ver1,Polygon[triang1]]],picfin,Graphics3D[GraphicsComplex[ver2,Polygon[triang2]]],Boxed->False,Axes->False]
Now this has just one problem. Here irrespective of the PlotPoints
there are four triangles always appearing that just shares only one edge with any other neighboring triangle. But we expect all of the triangles to share at least two edges with other trangles. That happens if we use belisarius method. But it creates too small triangles that my panel solver rejects as tingles with zero area.
One can check here the problem of my method. Here we will use the method from the solution by Sjoerd.
Export[NotebookDirectory[]<>"myres.obj",pic3D];
cd=Import[NotebookDirectory[]<>"myres.obj"];
polygons=(cd[[1]][[2]]/.GraphicsComplex-> List)[[2]][[1]][[1,1]];
pt=(cd[[1]][[2]]/.GraphicsComplex-> List)[[1]];
vertices=pt;
(*Split every triangle in 3 edges,with nodes in each edge sorted*)
triangleEdges=(Sort/@Subsets[#,{2}])&/@polygons;
(*Generate a list of edges*)
singleEdges=Union[Flatten[triangleEdges,1]];
(*Define a function which,given an edge (node number list),returns the bordering*)
(*triangle numbers.It's done by working through each of the triangles' edges*)
ClearAll[edgesNeighbors]
edgesNeighbors[_]={};
MapIndexed[(edgesNeighbors[#1[[1]]]=Flatten[{edgesNeighbors[#1[[1]]],#2[[1]]}];
edgesNeighbors[#1[[2]]]=Flatten[{edgesNeighbors[#1[[2]]],#2[[1]]}];
edgesNeighbors[#1[[3]]]=Flatten[{edgesNeighbors[#1[[3]]],#2[[1]]}];)&,triangleEdges];
(*Build a triangle relation table.Each'1' indicates a triangle relation*)
relations=ConstantArray[0,{triangleEdges//Length,triangleEdges//Length}];
Scan[(n=edgesNeighbors[##];
If[Length[n]==2,{n1,n2}=n;
relations[[n1,n2]]=1;relations[[n2,n1]]=1];)&,singleEdges]
(*Build a neighborhood list*)
triangleNeigbours=Table[Flatten[Position[relations[[i]],1]],{i,triangleEdges//Length}];
trires=Table[Flatten[{polygons[[i]],triangleNeigbours[[i]]}],{i,1,Length@polygons}];
Cases[Cases[trires,x_:>Length[x]],4]
Output shows always there are four triangles that shares only one edges with others.
{4,4,4,4}
In case of belisarius method we don't see this happening but there we get triangles with numerically zero areas.
BR
Import the data and construct the BSpline function as before:
dat = Import["Downloads/3DFoil.mat", "Data"];
fu = BSplineFunction[dat]
Generate the surface, making sure to include (only) the boundary line, which will follow the edge of the surface. Make sure to set Mesh
to either All
or None
.
pic = ParametricPlot3D[fu[u, v], {u, 0, 1}, {v, 0, 1}, Mesh -> None,
AspectRatio -> Automatic, PlotPoints -> 10, Boxed -> False,
Axes -> False, BoundaryStyle -> Red]
Extract the points from the boundary line:
bound = First@Cases[Normal[pic], Line[pts_] :> pts, Infinity]
Find the "corners", based on your parameter space:
corners = Flatten[Table[fu[u, v], {u, 0, 1}, {v, 0, 1}], 1]
Find the edge points best corresponding to the corners, keeping in mind that ParametricPlot3D
doesn't use the limits exactly, so we can't just use Position
:
nf = Nearest[bound -> Automatic];
nf /@ corners
Figure our which range of points on the boundary correspond to the areas you need to fill up. This step involved some manual inspection.
sets = {bound[[2 ;; 22]], bound[[22 ;; 52]], bound[[52 ;; 72]],
bound[[72 ;;]]}
Construct new polygons corresponding to the holes:
Graphics3D[Polygon[sets[[{1, 3}]]], Boxed -> False, BoxRatios -> 1]
Show[pic, Graphics3D[Polygon[sets[[{1, 3}]]]]]
Note that there is probably still a hole that can't be seen where the edge runs between the holes you mentioned, and I haven't tried to fill it in, but you should have enough information to do that if needed.
Your data set looks like this:
Graphics3D[Point@Flatten[dat, 1]]
It consists of 22 sections of 50 points.
Adding a mid-line in each end section (which is actually the end section flattened):
dat2 = Append[Prepend[dat,
Table[(dat[[1, i]] + dat[[1, -i]])/2, {i, Length[dat[[1]]]}]
],
Table[(dat[[-1, i]] + dat[[-1, -i]])/2, {i, Length[dat[[-1]]]}]
];
Graphics3D[{Point@Flatten[dat, 1], Red, Point@dat2[[1]], Green, Point@dat2[[-1]]}]
Now add some weights to the wingtip rim:
sw = Table[1, {24}, {50}];
sw[[2]] = 1000 sw[[1]];
sw[[-2]] = 1000 sw[[1]];
fu = BSplineFunction[dat2, SplineWeights -> sw];
Show[
ParametricPlot3D[fu[u, v], {u, 0, 1}, {v, 0, 1}, Mesh -> All,
AspectRatio -> Automatic, PlotPoints -> 20, Boxed -> False,
Axes -> False, Lighting -> "Neutral"
],
Graphics3D[{PointSize -> 0.025, Green, Point@dat2[[-1]], Red,Point@dat2[[-2]]}]
]
Note that I increased the PlotPoints
value to 20.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With