Understanding Mathematica contour lines extraction

Question

My Professor sent me a Mathematica Notebook, where he plots a 2d contour and then extracts the contour as lines (i.e. a list of lists of coordinate-tuples). The specific code segment looks like this:

xdiv = 0.05
zt = 1.
lis = Table[SomeFunction[x, y, zt], {y, -xmax, xmax, xdiv}, {x, -xmax, xmax, xdiv}];
plot = ListContourPlot[lis, PlotLegends -> Automatic, Contours -> {0}, ContourShading -> None];
lines = Cases[Normal@plot, Line[pts_, ___] :> xdiv*pts, Infinity];

I do not fully understand how exactly the code does what it does and an I'm seeking help for an explanation. I want to write a similar code in python and understand how. Also I want to know if I can extract the lines without using the ContourPlot function. I am not particularly interested in plotting the contours, I only need a list of its lines.

Edit: Rephrased the question. Also: matplotlib - extracting data from contour lines seems to explain how to achieve my goal in python. However, I don't really understand how it does what it does and judge if there is a better way, in terms of performance, to achieve it (I am dealing with very large lists, so this contour extraction seems like quite the bottleneck). I was looking for something more of an explanation to understand what exactly happens. The provided answer does really well at explaining the Mathematica code. I will read more into the matplotlib methods.

Well yes, but none of us are computer scientists so he couldn't explain how exactly the last line, using "Cases", does what it does. — Banana, Dec 20 '17 at 05:20
@Banana `Cases` is both a pattern matching function and a destructuring function. The easiest thing to do will be to approach the problem in a pythonic way, rather than trying to port Mathematica code over. Alternately, since it seems you're willing to use `subprocess` to get the job done, look into [Wolfram Script](http://reference.wolfram.com/language/ref/program/wolframscript.html). Then dump the lines you find into a more python-readable format on the Mathematica side. — b3m2a1, Dec 20 '17 at 05:26
@b3m2a1 Do you have an idea for a pythonic way to achieve this? What do you mean by "dumping into a more python-readable format"? Wolfram Script will only run the script, won't it? — Banana, Dec 20 '17 at 05:34
@Banana yeah look at how matplotlib (or similar) generates its contour plots. That will give you a sense for how to do this (or even just extract it from a matplotlib plot). Wolfram Script can run a script that can export to file or to stdout, both of which python could read. By more python-readable I mean convert all of the points `List`s in Mathematica to strings and then `StringReplace` all of the `{` and `}` with `[` and `]` so you can import it directly into python. — b3m2a1, Dec 20 '17 at 05:42
@b3m2a1 I actually already have the list exported to a file in the correct python readable format, the import as a threefold nested list of tuples is where I get stuck. I will look into matplotlib contour plots, as it seems it is the best option. I already had found an option to extract the contour points (here: https://stackoverflow.com/questions/5666056/matplotlib-extracting-data-from-contour-lines), but the problem is: it just returns all points in one array. It doesn't distinguish different lines. — Banana, Dec 20 '17 at 06:14
@Banana what do you mean `"the import as a threefold nested list of tuples is where I get stuck"`? Surely if you formatted it right `eval` will do the job. — b3m2a1, Dec 20 '17 at 06:20
@b3m2a1 I searched and "ast.literal_eval" did the job. However, everybody advises not to use it and as it is more convenient on the long run I will further try to figure out how to extract the paths directly from matplotlib. It shouldn't be too hard and I'm sure someone already wrote down the answer somewhere, I'm just too stupid to find it. — Banana, Dec 20 '17 at 06:36
this question really has nothing to do with mathematica. Your python version will not use anything remotely similar to `Cases` . You should delete the code from the question and just describe in words what you want to do. (at least fix the title ) — agentp, Dec 20 '17 at 14:35
https://stackoverflow.com/questions/5666056/matplotlib-extracting-data-from-contour-lines hth — Alan, Dec 20 '17 at 16:12

Chris Degnen · Accepted Answer · 2017-12-20T18:29:14.597

3

Here is some explanation of what is happening, with an example function:

SomeFunction[x_, y_, zt_] := Sin[3 x + y^2]
xmax = 4;

xdiv = 0.05;
zt = 1.;
lis = Table[SomeFunction[x, y, zt], {y, -xmax, xmax, xdiv}, {x, -xmax, xmax, xdiv}];
plot = ListContourPlot[lis, PlotLegends -> Automatic, Contours -> {0},
  ContourShading -> None];
normalplot = Normal@plot

cases = Cases[normalplot, Line[pts_, ___], Infinity];
Length[cases]
First[cases]

17
Line[{{16.2216, 1.}, {17., 1.29614}, ... , {16.7818, 160.782}, {16.2216, 161.}}]

The Cases statement extracts each line from the normalized plot. (Normal simplifies the graphics form to one amenable to Cases.) There are 17 separate lines of the form Line[List[{x, y}, {x, y}, ...]].

For each line List[{x, y}, {x, y}, ...] is represented by pts, so

lines = Cases[normalplot, Line[pts_, ___] :> xdiv*pts, Infinity]

multiplies each list of pts by xdiv.

0.05 * {{16.2216, 1.}, {17., 1.29614}, ... , {16.7818, 160.782}, {16.2216, 161.}}
= {{0.81108, 0.05}, {0.85, 0.0648069}, ... {0.839088, 8.03909}, {0.81108, 8.05}}

The lines can be plotted.

ListLinePlot[lines, AspectRatio -> 1]

edited Dec 20 '17 at 18:29

answered Dec 20 '17 at 16:05

Chris Degnen

8,443
2
23
40

Thank you @chris this really helps. Can you explain what the role of the Infinity keyword is? I can't judge from the documentation. Do you know how the Contourplot "finds" the lines and if there is a faster way to get a list of the lines directly? I don't necessarily need the plot – Banana Dec 21 '17 at 02:41
`Infinity` is the *levelspec* value. It specifies how many levels down into the expression the `Line[pts_, ___]` pattern should be matched. See the syntax for [Cases](http://reference.wolfram.com/language/ref/Cases.html). For obtaining contour lines I would say `ListContourPlot` is the ideal method. The plot is incidental - just a container for the calculated contours - and it needn't be displayed. – Chris Degnen Dec 21 '17 at 09:23
Thank you. I was confused because in the docs it seems like Infinity was the default value anyways, making this Argument redundant. – Banana Dec 21 '17 at 09:49
@Banana I don't think there is anything batter/faster than using the contour plot routine. you should go to mathematica.stackexchange.com for further questions on that. ( search first ). – agentp Dec 21 '17 at 17:33

Understanding Mathematica contour lines extraction

1 Answers1