Question

I've been trying to get the perimeters and areas of a bunch of contours (I'm using a level-set method to model burnback in solid/hybrid rocket motors). I've taken the approach recomended here:

matplotlib - extracting values from contour lines

Here's a minimal example that reproduces the error:

import numpy as np
import pylab as pl
import skfmm

n_pts = 1000
r_f = 98 #mm (final radius)
fin_l = 20 #mm (fin length)
fin_w = 10 #mm (fin length)

#set up grid
X, Y = np.meshgrid(np.linspace(-r_f,r_f,n_pts), np.linspace(-r_f,r_f,n_pts))
initial_grid = -1 * np.ones_like(X)

#set up initial shape
initial_grid[np.logical_and(np.abs(Y) < fin_l, np.abs(X) < fin_w/2.)] = 1
initial_grid[np.logical_and(np.abs(Y) < fin_w/2., np.abs(X) < fin_l)] = 1

#obtain 3D Array
contour_grid = skfmm.distance(initial_grid, dx=1e-2)

perims = np.array([])
webs = np.linspace(0, r_f, n_pts/10)
for web in webs:
    web_ctr = pl.contour(X, Y, contour_grid, [web])
    no_error_here = p = web_ctr.collections[0].get_paths()
    out_of_range_error = web_ctr.collections[0].get_paths()[0]

I get a "List index out of range" error, and it turns out that no_error_here is an empty list. I looked around a little, and I can't find any reason that matplotlib does this. Tips?

Was it helpful?

Solution

Basically, contour returns an empty collection if you contour a value outside the range of your data.

Either check that the value you're contouring is inside the bounds of your data or check to see if you get an empty collection.

Also, you're currently assuming that there will be only a single contour at the given value. That's probably a bad assumption. Because you're contouring a single values, you do want to use collections[0], but that contour is likely to have more than one path.


If you don't want to plot the results, why not just extract the contour values directly? skimage.measure.find_contours sounds like a better fit to what you're doing.

For example:

import numpy as np
from skimage.measure import find_contours

data = np.random.random((10,10))
contours = find_contours(data, 0.5)
for xy in contours:
    x, y = xy.T
    perim = np.hypot(np.diff(x), np.diff(y)).sum()
Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top