1
$\begingroup$

I'm loading a dataset {x,y,z,v} from a text file. It is about a triangulated surface and then I am using ListSurfacePlot3D to plot it. However, ListSurfacePlot3D is smoothing the surface, that is something I do not want.

There is an option to avoid this? Should I use a different method to plot?

Thanks.

Data:
https://www.dropbox.com/sh/ll9hb3oxxreufsg/AABvU7nufHhyHY3W4ausanJ6a?dl=0

Code:

data = Import["geometry.txt", "Table"];
ListSurfacePlot3D[data[[;; , 1 ;; 3]]]

enter image description here

This is not the intended shape, which can be discerned by doing a ListPointPlot3D

ListPointPlot3D[data[[All, ;; 3]], PlotStyle -> PointSize[.008]]

enter image description here

You can see that ListSurfacePlot3D omits a large number of points by overlaying the types of plot,

Show[ListSurfacePlot3D[#], ListPointPlot3D[#]] &@data[[All, ;; 3]]

enter image description here

ListSurfacePlot3D clearly does some interpolation, but there is no documentation on how to apply the Method option, which controls "the method to use for interpolation and data reduction".

$\endgroup$
  • $\begingroup$ Without the data and code, it's hard to say, but have you tried using the option InterpolationOrder -> 0 ? $\endgroup$ – Jason B. Nov 30 '15 at 8:48
  • $\begingroup$ ListSurfacePlot3D does not take InterpolationOrder $\endgroup$ – Fabio Nov 30 '15 at 8:51
  • $\begingroup$ I think the problem here is worse than simply smoothing, there are many points that are simply omitted from the surface plot altogether. $\endgroup$ – Jason B. Nov 30 '15 at 9:13
  • 1
    $\begingroup$ Try ListSurfacePlot3D[data[[;; , 1 ;; 3]], MaxPlotPoints -> 50, PerformanceGoal -> "Quality", Mesh -> False] $\endgroup$ – user9660 Nov 30 '15 at 11:28
  • $\begingroup$ Surely you can just get a nice plot from Import["geometry.stl"] right away? $\endgroup$ – Taiki Dec 7 '15 at 16:43
1
$\begingroup$

Going off of Lou's comments, it looks like setting the option MaxPlotPoints to a value of 20 or higher gives a better plot than the default. (I had tried PlotPoints but that isn't an option for ListSurfacePlot3D)

ListSurfacePlot3D[data[[All, ;; 3]], MaxPlotPoints -> #, 
    ImageSize -> 400] & /@ {6, 8}~Join~Range[10, 70, 10];
Grid@(%~Partition~4)

enter image description here

The option , PerformanceGoal -> "Quality" did not have an effect.

$\endgroup$
  • $\begingroup$ I guess this is the real issue. I've solved in another way, not using ListSurfacePlot3D, however solved. Thanks! $\endgroup$ – Fabio Dec 8 '15 at 19:05
  • $\begingroup$ @Fabio Maybe you should tell other people who find this thread how you solved it! $\endgroup$ – ZachMcDargh Mar 16 '16 at 20:18
  • 1
    $\begingroup$ @Zach: I did it below $\endgroup$ – Fabio Mar 29 '16 at 19:50
2
$\begingroup$

The way I solved is quite complicated and use gmsh http://gmsh.info/ to generate the stl file.

I have this package that enables loading .geo file and generation of stl file with gmsh

BeginPackage["Gmsh`"];

GmshCommandLine="~/Tools/gmsh-2.11.0-Linux/bin/gmsh";

RunGmsh[ifile_,ofile_]:=Run[GmshCommandLine<>" -2 "<>ifile<> " -o "<>ofile];

GmshArc3D[{a_, m_, b_}, n_:10, prim_: Line] := 
 Module[{\[Alpha], lab, axis, aarc, tm, alpha}, 
  lab = m + Norm[a - m]*Normalize[b - m];
  axis = (a - m)\[Cross](b - m);
  aarc = (VectorAngle[a - m, b - m]);
  tm = RotationMatrix[alpha, axis];
  prim@Table[m + tm.(a - m), {alpha, 0, aarc, aarc/n}]];

GmshParameters/:(var_/;Head@var===GmshGeo)[GmshParameters]:=var[[1]];
GmshPoints/:(var_/;Head@var===GmshGeo)[GmshPoints]:=var[[2]];
GmshLines/:(var_/;Head@var===GmshGeo)[GmshLines]:=var[[3]];
GmshCircles/:(var_/;Head@var===GmshGeo)[GmshCircles]:=var[[4]];

GmshPointFormat:="Point("~~n__~~")={"~~x__~~","~~y__~~","~~z__~~","~~l__~~"}";
GmshLineFormat:="Line("~~n__~~")={"~~p1__~~","~~p2__~~ "}";
GmshCircleFormat:="Circle("~~n__~~")={"~~s__~~","~~c__~~","~~e__~~"}";

GmshImportGeo[ifile_,evaluate_:False]:=Module[

  {temp,parameters,points,lines,circles},

  temp=ReadList[ifile,"String"];

  temp=StringDelete[temp,{" ",";"}];

  temp=Select[
    temp,
    StringMatchQ[
      #,
      {"LineLoop"~~___,"PlaneSurface"~~___,"RuledSurface"~~___}
    ]==False&
  ];

  parameters=Select[
    temp,
    StringMatchQ[
        #,
        {"Point"~~___,"Line"~~___,"Circle"~~___}
    ]==False&];

  parameters=ToExpression[
    StringReplace[
      parameters,
      {"="->"\[Rule]"}
    ]
  ];

  points=Select[
    temp,
    StringMatchQ[
      #,
      "Point"~~__
    ]==True&];

  points=ToExpression[
    StringReplace[
      points,
      {GmshPointFormat:>ToString[{n,{x, y, z},True}]}
    ]
  ];

  lines=Select[
    temp,
    StringMatchQ[
      #,
      "Line"~~__
    ]==True&
  ];

  lines=ToExpression[
    StringReplace[
      lines,
      {GmshLineFormat:>ToString[{n,{p1, p2},True}]}
    ]

just change the GmshCommandLine to point your gmsh binary.

Let's say you want draw a sphere, you have to create sphere.geo that contains commands for gmsh (see the help)

lc=0.1;
Point(1)={0,0,0,lc};
Point(2)={1,0,0,lc};
Point(3)={-1,0,0,lc};
Point(4)={0,1,0,lc};
Point(5)={0,-1,0,lc};
Point(6)={0,0,1,lc};
Point(7)={0,0,-1,lc};

Circle(1) = {2, 1, 4};
Circle(2) = {4, 1, 3};
Circle(3) = {3, 1, 5};
Circle(4) = {5, 1, 2};
Circle(5) = {7, 1, 3};
Circle(6) = {3, 1, 6};
Circle(7) = {6, 1, 2};
Circle(8) = {2, 1, 7};
Circle(9) = {7, 1, 4};
Circle(10) = {4, 1, 6};
Circle(11) = {6, 1, 5};
Circle(12) = {5, 1, 7};
Line Loop(13) = {7, 1, 10};
Ruled Surface(14) = {13};
Line Loop(15) = {1, -9, -8};
Ruled Surface(16) = {15};
Line Loop(17) = {9, 2, -5};
Ruled Surface(18) = {17};
Line Loop(19) = {2, 6, -10};
Ruled Surface(20) = {19};
Line Loop(21) = {4, -7, 11};
Ruled Surface(22) = {21};
Line Loop(23) = {4, 8, -12};
Ruled Surface(24) = {23};
Line Loop(25) = {12, 5, 3};
Ruled Surface(26) = {25};
Line Loop(27) = {3, -11, -6};
Ruled Surface(28) = {27};

and then run this mathematica script

<< Gmsh`
geo = GmshImportGeo["sphere.geo"];
RunGmsh["sphere.geo", "sphere.stl"];
mesh = Import["sphere.stl", "STL"];
geo = geo /. geo[GmshParameters];
Graphics3D[{GmshCompileGeo[geo, False], mesh[[1]]}]

here you have the results

enter image description here

looking at the package you can hide geometric features, i.e.

GmshHideLine[geo, {1,2}];

hides line 1,2 (see the .geo file)

That's all folk! F

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.