Extracting points from the MeshFunctions of Plot3D

Question

I am trying to extract the points from the mesh lines of a 3D plot. (Not my actual function.)

   p = Plot3D[Sin[x]^2*Cos[y], {x, -1, 1}, {y, -1, 1},  MeshFunctions -> {#1^2+#2^2 &}]

I found a similar question here. So I used the code suggested:

p2 = p // Cases[#, GraphicsComplex[Line_, ___] :> Line] &;  
p3 = Flatten[p2, 1];  
ListPlot[p3[[All, 1 ;; 2]]]  

As you can see, I am getting a bunch of points that are not on the mesh lines. When I increase the PlotPoints, it gets even worse.

How can I extract only the points on the mesh lines?

Answer 1

Expanding @Guess comment:

p1 = Join @@ Cases[Normal@p, Line[x1__] :> x1, Infinity];
ListPlot[Most /@ p1]

---

p1 = Join @@ Cases[Normal@p, Line[x1__] :> 
                  {RGBColor @@ RandomReal[{0, 1}, 3], Line[Most /@ x1]}, Infinity];
Graphics[p1, AspectRatio -> 1/GoldenRatio]

Answer 2

Working with GraphicsComplex retains a degree of flexibility. For instance,

Graphics3D[GraphicsComplex[p[[1, 1]], Line[Rest@Cases[p, Line[z__] :> z, Infinity]]]]

gives the Mesh in 3D. (Rest@ deletes the perimeter of the surface.) If, instead, a plot of the points in 3D is desired, use

Graphics3D[GraphicsComplex[p[[1, 1]], 
  Point[Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]]]]]

The same plot in 2D is obtained by dropping the last coordinate of each point.

p1 = Graphics[GraphicsComplex[Most /@ p[[1, 1]], 
   Point[Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]]]]]

Addendum

A List of the points themselves can be obtained directly from p.

Most /@ (p[[1, 1, #]] & /@ Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]])

or from p1, the 2D plot,

p1[[1, 1, #]] & /@ Flatten[Cases[p1, Point[z__] :> z, Infinity]]

which is equivalent to

Cases[Normal[p1], Point[z__] :> z, Infinity]

suggested by Guesswhoitis and belisarius.

Note: It might be tempting to try simply, Most /@ (p[[1, 1, #]], but doing so also recovers hundreds of additional points used to construct the 3D surface.

Explanation

In response to a comment below, this is how (more or less) I obtained the results above. To begin with, run

p//InputForm

to see what is inside, but be prepared for a lot of data. Basically, you will see the function GraphicsComplex, followed by several options at the very end. So, for instance,

Graphics3D[p[[1]]]

plots only the contents of the GraphicsComplex, which in this case is just the surface and the mesh lines.

The first argument of GraphicsComplex (see documentation) is a large array of 3D point locations, from which almost everything in the plot is constructed. The next argument consists of several other functions. What we care about here are Line functions, one for the boundary of the surface and one for each line in the mesh. The argument of Line is a list of indices, referring to the 3D points; i.e., Line[{3 ,127, 52}] draws a smooth curve from point 3 to point 127 to point 52.

The command

Line[Rest@Cases[p, Line[z__] :> z, Infinity]]

extracts a list of the lists of indices, throws away the first list (which is the outline of the surface), and puts them inside another Line function. Additionally p[[1, 1]] obtains the list of 3D points Finally, both of these items are repackaged in ``GraphicsComplex` and plotted.

Graphics3D[GraphicsComplex[p[[1, 1]], Line[Rest@Cases[p, Line[z__] :> z, Infinity]]]]

giving the mesh lines in 3D shown in the first plot above.

If, instead, the Line indices are repackaged inside Point,

Point[Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]]]

and combined with p[[1, 1]] inside GraphicsComplex, then

Graphics3D[GraphicsComplex[p[[1, 1]], 
  Point[Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]]]]]

gives a 3D plot of mesh points instead, as in the second plot above.

To obtain a 2D plot instead, use Most to discard the third dimension.

Graphics[GraphicsComplex[Most /@ p[[1, 1]], 
   Point[Flatten[Rest@Cases[p, Line[z__] :> z, Infinity]]]]]

which gives the third plot above. I hope that this brief explanation is helpful.