I still have to make an addition:

When data of more realistic structures are used, e.g.:

the method above simply does not work! One reason is the possible big difference in the length of neighboring "wires". This graphic shows as an example an imprint (i.e. negative form) of an human ear - with some typical "complications" added like show by these slices:

The problem here is that those slices contain more than just a single line/wire - and they have different meaning. What I want to end up with is of course:

Hence a completely different approach was needed:

• Basically I create from a stack of those images a 3D array;
• from this an interpolated function is defined;
• using RegionPlot3D a graphic is generated;
• finally this graphic is converted into a mesh.

The function mkImage does the relevant job of creating the correct images. Here I had to use high level functions like FillingTransform or SelectComponents. My full (but short) code now reads:

ClearAll["Global`*"]

mkImage[wire2D_, xyRange_] := Module[{imgGraphics, img0, img1},
  imgGraphics = 
   Binarize[
    ColorNegate[(ColorConvert[#, "Grayscale"] &) /@ 
      Image[Graphics[{White, EdgeForm[Black], Polygon[wire2D]}, 
        PlotRange -> xyRange]]]];
  img0 = FillingTransform[imgGraphics];
  img1 = FillingTransform@
    SelectComponents[imgGraphics, #EnclosingComponentCount == 1 &];
  ImageSubtract[img0, img1]
  ]

SetDirectory[NotebookDirectory[]]

sliceData = << "Bolus3D_complex.txt";

planarAssoc = 
  GroupBy[sliceData, #[[1, 3]] &] /. {x_, y_, z_} :> {x, y};

border = 2;
range3D = {{-border, border}, {-border, border}, {0, 0}} + 
  CoordinateBounds@sliceData

imgs = mkImage[#, Most[range3D]] & /@ Values[planarAssoc];

data3D = Transpose[ImageData /@ imgs, {3, 2, 1}];

func3D = ListInterpolation[data3D, {#1, Reverse[#2], #3} & @@ range3D,
   InterpolationOrder -> 3]

gr = RegionPlot3D[func3D[\[FormalX], \[FormalY], \[FormalZ]] > .5, 
  Evaluate[
   Sequence @@ 
    MapThread[
     Prepend, {range3D, {\[FormalX], \[FormalY], \[FormalZ]}}]],
  BoxRatios -> Automatic, PlotPoints -> {100, 100, 100}, 
  ImageSize -> 700]

mesh = DiscretizeGraphics[gr, ImageSize -> 700]

Export["bolus.stl", mesh, "STL"]

And this now seems to work! I end up with:

Many thanks again to all who have helped!

Regards -- Henrik

Attachments:

Great idea, @Henrik Schachner and a good choice for a function to minimize.

BTW @Henrik Schachner, where are this data from? What project are you working on and what is your final goal? It is always interesting to know the story behind the data.

Happy to help, @Henrik Schachner ! And there is no "lengthy explanations!", just the detailed ones - and that's a good thing. We added the data description in the top post, - this is what makes this case interesting and potentially useful to other professionals. Thank you for sharing!

Very neat idea!

If OP wants to print the solid head, he can fill the holes using RepairMesh.

plot = ParametricPlot3D[surfBSF[u, v], {u, 0, 1}, {v, 0, 1}, PlotRange -> All, PlotPoints -> 100];
RepairMesh[DiscretizeGraphics[plot], PlotTheme -> "SmoothShading"]

RepairMesh[DiscretizeGraphics[plot], PlotTheme -> "SmoothShading"]

@Chip Hurst, thank you again for this nice hint! I think the option PlotTheme -> "SmoothShading" is important, but I cannot find it anywhere in the documentation. It would be a pity if it were just undocumented - or am I again missing something?

Oooops .... thanks again! Henrik

The second idea was to transform everything into slices and then use Image3D.

sliceData = ToExpression[Import["/Users/thiel/Desktop/SliceData.txt"]];
range = {MinMax[#[[All, 1]]], MinMax[#[[All, 2]]]} &@Flatten[sliceData[[All, All, {1, 2}]], 1];
imgs = Table[Image@Graphics[Polygon[sliceData[[k, All, {1, 2}]]], PlotRange -> range], {k, 1, Length[sliceData]}];

This gives slices like these:

Partition[imgs[[1 ;; ;; 10]], 3] // TableForm

Then something like this should work:

Image3D[imgs, "Bit", BoxRatios -> {1, 1, 3}, Background -> White, 
 ColorFunction -> "Grayscale"]

But it doesn't. So I tried something more like this:

graphics = Table[Graphics3D[Polygon[sliceData[[k, All, {1, 2, 3}]]]], {k, 1, Length[sliceData]}];
Show[graphics]

The thing is that this is not a "real" 3D object but just a bunch of 2D planes in 3D space. This here would "fill the gaps".

Show[Table[
  Graphics3D[
   Polygon[Riffle[
     Append[#, sliceData[[k, 1, -1]] + 1.5] & /@ 
      sliceData[[k, All, {1, 2}]], 
     Append[#, sliceData[[k, 1, -1]] - 1.5] & /@ 
      sliceData[[k, All, {1, 2}]]]]], {k, 1, Length[sliceData]}]]

Cheers,

Marco

Hi Henrik,

there is this:

ListSurfacePlot3D[Flatten[sliceData, 1], MaxPlotPoints -> 75]

but there are obvious artefacts. I have another approach which I will implement next.

Cheers,

Marco