Arbitrary image slice with python/numpy

Question

I would like to plot a 1D profile of a 2D image along an arbitrary line. The code below loads the image data hosted on github and plots it:

import urllib
import numpy as np
import matplotlib.pyplot as plt

url = "https://gist.github.com/andreiberceanu/7141843/raw/0b9d50d3d417b1cbe651560470c098700df5a1fc/image.dat"
f = urllib.urlopen(url)
data = np.loadtxt(f)

plt.imshow(data)

waves with line

The red line in the plot above was drawn by hand, as an example. I suppose one can parametrize it in the form a*x + b. I am also guessing some sort of interpolation is necessary, because the line passes though points which may not be part of the original 2D array of data.

See this related question: stackoverflow.com/questions/7878398/… — aganders3
– aganders3, Commented Oct 24, 2013 at 22:15

Ben · Accepted Answer · 2013-10-28 12:18:30Z

2

You want to use scipy.ndimage.map_coordinates. You need to build up a 2xn array that is the coordinates at which to sample and then do map_coordinates(im, samples).

I think this is it:

def sliceImage(I, a, b, *arg, **kws):
    from scipy import linspace, asarray
    from scipy.ndimage import map_coordinates
    from scipy.linalg import norm
    dst = norm(asarray(b) - a) + 1
    return map_coordinates(I, [linspace(strt, end, dst) 
                               for strt, end in zip(a, b)],
                           *arg, **kws)

Edit: On further consideration, I think this is more elegant:

def sliceImage(I, a, b, *arg, **kws):
    from scipy import linspace, asarray
    from scipy.ndimage import map_coordinates
    from scipy.linalg import norm
    a = asarray(a)
    b = asarray(b)
    dst = norm(b - a) + 1
    return map_coordinates(I, (a[:,newaxis] * linspace(1, 0, dst) +
                               b[:,newaxis] * linspace(0, 1, dst)),
                           *arg, **kws)

Edit: Thanks tcaswell: Added 1 to dst.

edited Oct 28, 2013 at 12:18

answered Oct 24, 2013 at 22:57

Ben

9,8933 gold badges50 silver badges62 bronze badges

Sign up to request clarification or add additional context in comments.

3 Comments

tacaswell Over a year ago

I feel there are clearer ways to generate the list of points to sample at.

Ben Over a year ago

I agree. There's this approach, which is more numpyish: a = asarray(a); b = asarray(b); dst = norm(b - a); return map_coordinates(I, (a[:,newaxis] * linspace(1, 0, dst) + b[:,newaxis] * linspace(0, 1, dst)), *arg, **kws)

tacaswell Over a year ago

I would use something like linspace(0, 1, int(np.ceil(dst)) + 1) to make sure you get enough points. (make sure a, b =(0, 0), (0, 1) gives back what you think it should)

Collectives™ on Stack Overflow

Arbitrary image slice with python/numpy

1 Answer 1

3 Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

3 Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related