# Joining points into 3D surface

This seems like a simple problem, but I can't find any questions like it here so I'm making a new one. Apologies if I missed one.

I have a list of points contained in a file which can be visualized like so:

data = Import["http://www.sharecsv.com/dl/782ae64cbc8003d67f4ab98d784015d6/data.csv"];
{a, b} = {data[[1 ;; 335]], data[[336 ;; -1]]};
ListPointPlot3D[{a, b}, PlotStyle -> {Red, Blue}]


The problem is connecting/interpolating these points into a surface, ListSurfacePlot3D[] yields a garbled mess for any MaxPlotPoints.

ListSurfacePlot3D[data]


How can I transform this set of points into a (smooth) surface?

Update:

Alright fellas I'm nearly there. The best way I have come up with is to use DelaunayMesh with RunnyKine's code to produce concave meshes here. First I split the data points into non-overlapping sections by a plane between them. This is a naff hack but what can you do. Then I construct a 2D mesh for the top and bottom halves, replacing the third coordinate with Dispatch.

alphaShapes2D[points_,crit_]:=Module[{alphacriteria,del=Quiet@DelaunayMesh@points,tetras,tetcoords,tetradii,selectExternalFaces},alphacriteria[tetrahedra_,radii_,rmax_]:=Pick[tetrahedra,UnitStep@Subtract[rmax,radii],1];
selectExternalFaces[facets_]:=MeshRegion[points,facets];
tetcoords=MeshPrimitives[del,2][[All,1]];

data=Standardize@Import["http://www.sharecsv.com/dl/782ae64cbc8003d67f4ab98d784015d6/data.csv"];
{a,b}={data[[1;;335]],data[[336;;-1]]};
{bot, top} = {True, False} /. GroupBy[data, Last@# < .25 First@#&];
bottri=MeshPrimitives[alphaShapes2D[bot[[All,{1,2}]],.1],2][[All,1]]/.Dispatch[#[[{1,2}]]->#&/@bot];
toptri=MeshPrimitives[alphaShapes2D[top[[All,{1,2}]],.1],2][[All,1]]/.Dispatch[#[[{1,2}]]->#&/@top];
out=Graphics3D[{Polygon/@bottri,Polygon/@toptri}]


Now I'm just missing a way to replace the missing polygons. It would have been much better to cut the data along vertices belonging to both sets. Is there a simple way you can envisage doing this?

• Thank you very much for providing the points at the outset; a lot of people don't do that. If I may, where'd these points come from? – J. M. will be back soon Jun 6 '16 at 7:47
• Some very messy and inelegant code. I'm looking at fixed points in random 2D vector fields. Here the data[[All,3]] axis is one coordinate of a particular fixed point and the data[[All,{1,2}]] axes map out the domain that fixed point exists. I think Mathematica is getting tripped up by the bifurcation resulting in a fold in the surface. – Crêpo Jun 6 '16 at 7:55
• do you have access to the algorithm that generated the data, or just the raw points? It looks like its got a regular grid structure to it and if you had that grid info it would be useful. – george2079 Jun 6 '16 at 12:15
• I generated the data by first creating the boundary points, and then sampling uniformly in the data[[All,{1,2}]] plane. In principle I can generate as many points in whatever arrangement I want, although it doesn't seem to improve the result. – Crêpo Jun 6 '16 at 12:34
• I added a figure from a different viewpoint to show this is a single simple smooth surface. – george2079 Jun 6 '16 at 14:15

I think your best bet is to generate a triangualted polygon surface. Here is an initial stab at it:

p2d = b[[All, 1 ;; 2]];
tri1[i_, j_] :=
If[EvenQ[j], {{i, j}, {i + 1, j}, {i + 1/2, j + 1}} .05, {{i + 1/2,
j}, {i + 1 + 1/2, j}, {i + 1, j + 1}} .05]
tri2[i_, j_] :=
If[EvenQ[j], {{i, j}, {i + 1, j}, {i + 1/2, j - 1}} .05, {{i + 1/2,
j}, {i + 1 + 1/2, j}, {i + 1, j - 1}} .05]
triangles =
Flatten[Table[ tri2[i, j], {i, -45, 35}, {j, -55, 20}], 1]~Join~
Flatten[Table[ tri1[i, j], {i, -45, 35}, {j, -55, 20}], 1];
near = Nearest[p2d];
good = Select[triangles,
Total[Norm[(near[#, 1][[1]] - #)] & /@ #] < .001 &];
p3d[x_] := SortBy[b, Norm[#[[1 ;; 2]] - x] &][[1]]
sp3d = (p3d /@ #) & /@ good;
Graphics3D[{EdgeForm[None],
Polygon /@ (Select[sp3d, Variance[#[[All, 3]]] < .0001 &])},
BoxRatios -> {1, 1, 1}]


At this point I havn't dealt with the boundaries or the overlapping portions of the surface, but It should get you started.

• Yes, this is the kick in the arse I needed. I can actually get my bulk triangles at runtime with a little jiggery-pokery, although it won't shed much light on how to perform this procedure more generally though. I've accepted this answer and will update the question with the full result when I get it. – Crêpo Jun 6 '16 at 17:11

Not sure how useful will this be. Usually Standardize-d data is giving better results.

We can also cut off parts of the surface which are distant from our data.

rf = Nearest[Standardize@data]

sr = ListSurfacePlot3D[Standardize @ data,
BoxRatios -> 1,
MaxPlotPoints -> 100,
RegionFunction -> Function[{x, y, z},
Norm[{x, y, z} - rf[{x, y, z}][[1]]] < .2]
]


m = Mean @ data
sd = StandardDeviation @ data

Show[
Graphics3D[
GeometricTransformation[
First@sr,
TranslationTransform[m]@*ScalingTransform[sd]
]
]
,
ListPointPlot3D[{a, b}, PlotStyle -> {Red, Blue}  ]
,
SphericalRegion -> True, BoxRatios -> 1
]


• I think the main problem for the sr in your answer is smoothness. – xyz Jun 6 '16 at 14:31
• @ShutaoTANG don't claim it is best approach, but it is 1 line and looks way better than original approach. – Kuba Jun 6 '16 at 14:32