plot contour with levels from discrete data

[Rehamm]

I am using Python 2.7.3 and Matplotlib 1.2.0 in Linux.
I have three parameters x, y and z each having 61 values. I would like to create a contour plot with levels 0.5, -2.3, -4.61, -9.21. The problem is that Z is 1D where it should be 2D array and so I used griddata. The code that I used is
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.mlab import griddata
x=[4, 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8, 5.9, 6, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7, 6.8, 6.9, 7, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 8, 8.1, 8.2, 8.3, 8.4, 8.5, 8.6, 8.7, 8.8, 8.9, 9, 9.1, 9.2, 9.3, 9.4, 9.5, 9.6, 9.7, 9.8, 9.9, 10]
y=[-6.95138e-06, -9.09998e-07, 8.24384e-06, 3.4941e-06, 5.08276e-06, 7.82652e-06, -4.7378e-06, -1.40027e-05, -1.62638e-05, -3.97604e-06, 3.19294e-06, 2.50123e-06, -4.13063e-06, -6.60289e-06, -4.02982e-06, -2.32882e-06, -3.86464e-06, -1.09167e-05, -9.42387e-06, -3.48118e-07, 5.22e-06, 8.74445e-06, 1.35842e-05, 2.33632e-05, 2.71328e-05, 1.747e-05, 2.32177e-06, -7.25386e-06, -9.75881e-06, -2.99633e-06, 1.19281e-06, -4.24077e-06, -7.4252e-06, -4.54435e-07, 1.03078e-05, 1.14579e-05, 3.90613e-06, -4.77174e-06, -9.25321e-06, -8.36579e-06, -3.0257e-06, -1.69309e-06, -5.36534e-06, -4.01092e-06, 1.20577e-06, 5.13284e-06, 5.06792e-06, 4.81178e-06, 5.9607e-06, 6.70492e-06, 3.45118e-06, 2.51942e-06, 1.23012e-06, 2.09802e-06, 1.44658e-06, -8.93274e-08, -5.14753e-06, -9.93717e-06, -7.91692e-06, -4.12816e-06, -6.33457e-06]
z=[-0.63, -0.02, -1.05, -0.22, -0.51, -1.26, -0.53, -4.97, -7.32, -0.44, -0.30, -0.19, -0.55, -1.48, -0.58, -0.20, -0.57, -5.00, -3.84, -0.01, -1.17, -3.31, -8.13, -22.59, -30.52, -13.69, -0.27, -2.88, -5.48, -0.49, -0.08, -0.86, -2.94, -0.01, -4.32, -5.49, -0.69, -1.15, -4.70, -3.73, -0.44, -0.14, -1.49, -0.77, -0.07, -1.08, -1.04, -0.93, -1.42, -1.76, -0.51, -0.25, -0.07, -0.18, -0.09, -0.00, -1.08, -5.03, -2.64, -0.65, -1.65]
xi = np.linspace(min(x), max(x), 30)
yi = np.linspace(min(y), max(y), 30)
[X, Y] = np.meshgrid(xi, yi)
Z = np.griddata(x, y, z, X, Y)
plt.figure()
levels = [0.5, -2.3, -4.61, -9.21]
contour = plt.contour(X, Y, Z, levels)

But the step of griddata give me this error
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/usr/lib64/python2.7/site-packages/matplotlib/mlab.py", line 2775, in griddata
tri = delaunay.Triangulation(x,y)
File "/usr/lib64/python2.7/site-packages/matplotlib/delaunay/triangulate.py", line 123, in __init__
self.hull = self._compute_convex_hull()
File "/usr/lib64/python2.7/site-packages/matplotlib/delaunay/triangulate.py", line 158, in _compute_convex_hull
hull.append(edges.pop(hull[-1]))
KeyError: 0

What could be the problem?

[afvincent]

Just a few comments (while waiting for someone more expert than myself in contour and griddata):

• is it normal that you call np.griddata while you imported griddata from matplotlib.mlab? (I guess it is a typo)
• in plt.contour docstring, it is written that levels have to be increasing. When I test your example (with np.griddata <- griddata and also passing interp='linear' to griddata, because otherwise mpl complains about natgrid not being installed on my computer), on mpl 1.5.1, Python 2.7.12 & Linux, it raises indeed ValueError: Contour levels must be increasing, while it seems to be working with levels=sorted([0.5, -2.3, -4.61, -9.21]).

[Rehamm]

• Iam sorry I meant griddata not np.griddata.
• thanks a lot for your note about using levels but also with adding interp='linear' it still give me the same error from griddata. When you tested my code on your computer, Did it work without error from griddata?

[afvincent]

Running

import numpy as np
import matplotlib.pyplot as plt

from matplotlib.mlab import griddata

x = [4, 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8, 5.9, 6, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7, 6.8, 6.9, 7, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 8, 8.1, 8.2, 8.3, 8.4, 8.5, 8.6, 8.7, 8.8, 8.9, 9, 9.1, 9.2, 9.3, 9.4, 9.5, 9.6, 9.7, 9.8, 9.9, 10]
y = [-6.95138e-06, -9.09998e-07, 8.24384e-06, 3.4941e-06, 5.08276e-06, 7.82652e-06, -4.7378e-06, -1.40027e-05, -1.62638e-05, -3.97604e-06, 3.19294e-06, 2.50123e-06, -4.13063e-06, -6.60289e-06, -4.02982e-06, -2.32882e-06, -3.86464e-06, -1.09167e-05, -9.42387e-06, -3.48118e-07, 5.22e-06, 8.74445e-06, 1.35842e-05, 2.33632e-05, 2.71328e-05, 1.747e-05, 2.32177e-06, -7.25386e-06, -9.75881e-06, -2.99633e-06, 1.19281e-06, -4.24077e-06, -7.4252e-06, -4.54435e-07, 1.03078e-05, 1.14579e-05, 3.90613e-06, -4.77174e-06, -9.25321e-06, -8.36579e-06, -3.0257e-06, -1.69309e-06, -5.36534e-06, -4.01092e-06, 1.20577e-06, 5.13284e-06, 5.06792e-06, 4.81178e-06, 5.9607e-06, 6.70492e-06, 3.45118e-06, 2.51942e-06, 1.23012e-06, 2.09802e-06, 1.44658e-06, -8.93274e-08, -5.14753e-06, -9.93717e-06, -7.91692e-06, -4.12816e-06, -6.33457e-06]
z = [-0.63, -0.02, -1.05, -0.22, -0.51, -1.26, -0.53, -4.97, -7.32, -0.44, -0.30, -0.19, -0.55, -1.48, -0.58, -0.20, -0.57, -5.00, -3.84, -0.01, -1.17, -3.31, -8.13, -22.59, -30.52, -13.69, -0.27, -2.88, -5.48, -0.49, -0.08, -0.86, -2.94, -0.01, -4.32, -5.49, -0.69, -1.15, -4.70, -3.73, -0.44, -0.14, -1.49, -0.77, -0.07, -1.08, -1.04, -0.93, -1.42, -1.76, -0.51, -0.25, -0.07, -0.18, -0.09, -0.00, -1.08, -5.03, -2.64, -0.65, -1.65]

xi = np.linspace(min(x), max(x), 30)
yi = np.linspace(min(y), max(y), 30)
X, Y = np.meshgrid(xi, yi)
Z = griddata(x, y, z, X, Y, interp='linear')  # switched to linear interpolation

plt.style.use('classic')  # just to make it look a bit more like older matplotlib plots
plt.figure()
levels = [-9.21, -4.61, -2.3, 0.5]  # levels sorted by increasing values
contour = plt.contour(X, Y, Z, levels)
plt.clabel(contour, fontsize='medium', inline_spacing=2, fmt='%0.3g')  # adding labels on the contour lines

plt.show()

seems fine with matplotlib 1.5.1 and produces

So I guess something has changed in griddata between mpl 1.2 and mpl 1.5.1…

[ianthomas23]

Matplotlib version 1.2.0 is quite old. Since then we have changed the Delaunay triangulator which is used by both mlab.griddata and the triangulation module. The old Delaunay triangulator used to give errors just like the one you have experienced; the new Delaunay triangulator is much more robust. You should upgrade matplotlib and try again.

Also, you have an unstructured (triangular) grid, but you don't need to create a rectangular gridded version of the data using griddata to produce useful output. Instead you can use the various triangular grid functions, like tricontour, directly.

For example:

import numpy as np
import matplotlib.pyplot as plt

x = [4, 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8, 5.9, 6, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7, 6.8, 6.9, 7, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 8, 8.1, 8.2, 8.3, 8.4, 8.5, 8.6, 8.7, 8.8, 8.9, 9, 9.1, 9.2, 9.3, 9.4, 9.5, 9.6, 9.7, 9.8, 9.9, 10]
y = [-6.95138e-06, -9.09998e-07, 8.24384e-06, 3.4941e-06, 5.08276e-06, 7.82652e-06, -4.7378e-06, -1.40027e-05, -1.62638e-05, -3.97604e-06, 3.19294e-06, 2.50123e-06, -4.13063e-06, -6.60289e-06, -4.02982e-06, -2.32882e-06, -3.86464e-06, -1.09167e-05, -9.42387e-06, -3.48118e-07, 5.22e-06, 8.74445e-06, 1.35842e-05, 2.33632e-05, 2.71328e-05, 1.747e-05, 2.32177e-06, -7.25386e-06, -9.75881e-06, -2.99633e-06, 1.19281e-06, -4.24077e-06, -7.4252e-06, -4.54435e-07, 1.03078e-05, 1.14579e-05, 3.90613e-06, -4.77174e-06, -9.25321e-06, -8.36579e-06, -3.0257e-06, -1.69309e-06, -5.36534e-06, -4.01092e-06, 1.20577e-06, 5.13284e-06, 5.06792e-06, 4.81178e-06, 5.9607e-06, 6.70492e-06, 3.45118e-06, 2.51942e-06, 1.23012e-06, 2.09802e-06, 1.44658e-06, -8.93274e-08, -5.14753e-06, -9.93717e-06, -7.91692e-06, -4.12816e-06, -6.33457e-06]
z = [-0.63, -0.02, -1.05, -0.22, -0.51, -1.26, -0.53, -4.97, -7.32, -0.44, -0.30, -0.19, -0.55, -1.48, -0.58, -0.20, -0.57, -5.00, -3.84, -0.01, -1.17, -3.31, -8.13, -22.59, -30.52, -13.69, -0.27, -2.88, -5.48, -0.49, -0.08, -0.86, -2.94, -0.01, -4.32, -5.49, -0.69, -1.15, -4.70, -3.73, -0.44, -0.14, -1.49, -0.77, -0.07, -1.08, -1.04, -0.93, -1.42, -1.76, -0.51, -0.25, -0.07, -0.18, -0.09, -0.00, -1.08, -5.03, -2.64, -0.65, -1.65]

levels = [-9.21, -4.61, -2.3, 0.5]  # Levels must be increasing.
plt.tricontour(x, y, z, levels)
plt.show()

Using this approach you write less code and the results are more accurate, not depending on the resolution of the griddata grid.

If you want to see the grid that is used for the tricontour call, add something like the following line:

plt.triplot(x, y, 'bo-', alpha=0.1)

[Rehamm]

Thanks a lot after upgrade the matplotlib it works.