You could use the boolean array produced by in1d to index an arange. Reversing a so that the indices are different from the values:
>>> a[::-1]
array([10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0])
>>> a = a[::-1]
intersect1d still returns the same values...
>>> numpy.intersect1d(a, b)
array([ 2, 7, 10])
But in1d returns a boolean array:
>>> numpy.in1d(a, b)
array([ True, False, False, True, False, False, False, False, True,
False, False], dtype=bool)
Which can be used to index a range:
>>> numpy.arange(a.shape[0])[numpy.in1d(a, b)]
array([0, 3, 8])
>>> indices = numpy.arange(a.shape[0])[numpy.in1d(a, b)]
>>> a[indices]
array([10, 7, 2])
To simplify the above, though, you could use nonzero -- this is probably the most correct approach, because it returns a tuple of uniform lists of X, Y... coordinates:
>>> numpy.nonzero(numpy.in1d(a, b))
(array([0, 3, 8]),)
Or, equivalently:
>>> numpy.in1d(a, b).nonzero()
(array([0, 3, 8]),)
The result can be used as an index to arrays of the same shape as a with no problems.
>>> a[numpy.nonzero(numpy.in1d(a, b))]
array([10, 7, 2])
But note that under many circumstances, it makes sense just to use the boolean array itself, rather than converting it into a set of non-boolean indices.
Finally, you can also pass the boolean array to argwhere, which produces a slightly differently-shaped result that's not as suitable for indexing, but might be useful for other purposes.
>>> numpy.argwhere(numpy.in1d(a, b))
array([[0],
[3],
[8]])
If you need to get unique values as given by intersect1d:
import numpy as np
a = np.array([range(11,21), range(11,21)]).reshape(20)
b = np.array([12, 17, 20])
print(np.intersect1d(a,b))
#unique values
inter = np.in1d(a, b)
print(a[inter])
#you can see these values are not unique
indices=np.array(range(len(a)))[inter]
#These are the non-unique indices
_,unique=np.unique(a[inter], return_index=True)
uniqueIndices=indices[unique]
#this grabs the unique indices
print(uniqueIndices)
print(a[uniqueIndices])
#now they are unique as you would get from np.intersect1d()
Output:
[12 17 20]
[12 17 20 12 17 20]
[1 6 9]
[12 17 20]
Stolen directly from https://web.archive.org/web/20111108065352/https://www.cs.mun.ca/~rod/2500/notes/numpy-arrays/numpy-arrays.html
#
# line segment intersection using vectors
# see Computer Graphics by F.S. Hill
#
from numpy import *
def perp( a ) :
b = empty_like(a)
b[0] = -a[1]
b[1] = a[0]
return b
# line segment a given by endpoints a1, a2
# line segment b given by endpoints b1, b2
# return
def seg_intersect(a1,a2, b1,b2) :
da = a2-a1
db = b2-b1
dp = a1-b1
dap = perp(da)
denom = dot( dap, db)
num = dot( dap, dp )
return (num / denom.astype(float))*db + b1
p1 = array( [0.0, 0.0] )
p2 = array( [1.0, 0.0] )
p3 = array( [4.0, -5.0] )
p4 = array( [4.0, 2.0] )
print seg_intersect( p1,p2, p3,p4)
p1 = array( [2.0, 2.0] )
p2 = array( [4.0, 3.0] )
p3 = array( [6.0, 0.0] )
p4 = array( [6.0, 3.0] )
print seg_intersect( p1,p2, p3,p4)
import numpy as np
def get_intersect(a1, a2, b1, b2):
"""
Returns the point of intersection of the lines passing through a2,a1 and b2,b1.
a1: [x, y] a point on the first line
a2: [x, y] another point on the first line
b1: [x, y] a point on the second line
b2: [x, y] another point on the second line
"""
s = np.vstack([a1,a2,b1,b2]) # s for stacked
h = np.hstack((s, np.ones((4, 1)))) # h for homogeneous
l1 = np.cross(h[0], h[1]) # get first line
l2 = np.cross(h[2], h[3]) # get second line
x, y, z = np.cross(l1, l2) # point of intersection
if z == 0: # lines are parallel
return (float('inf'), float('inf'))
return (x/z, y/z)
if __name__ == "__main__":
print get_intersect((0, 1), (0, 2), (1, 10), (1, 9)) # parallel lines
print get_intersect((0, 1), (0, 2), (1, 10), (2, 10)) # vertical and horizontal lines
print get_intersect((0, 1), (1, 2), (0, 10), (1, 9)) # another line for fun
Explanation
Note that the equation of a line is ax+by+c=0. So if a point is on this line, then it is a solution to (a,b,c).(x,y,1)=0 (. is the dot product)
let l1=(a1,b1,c1), l2=(a2,b2,c2) be two lines and p1=(x1,y1,1), p2=(x2,y2,1) be two points.
Finding the line passing through two points:
let t=p1xp2 (the cross product of two points) be a vector representing a line.
We know that p1 is on the line t because t.p1 = (p1xp2).p1=0.
We also know that p2 is on t because t.p2 = (p1xp2).p2=0. So t must be the line passing through p1 and p2.
This means that we can get the vector representation of a line by taking the cross product of two points on that line.
Finding the point of intersection:
Now let r=l1xl2 (the cross product of two lines) be a vector representing a point
We know r lies on l1 because r.l1=(l1xl2).l1=0. We also know r lies on l2 because r.l2=(l1xl2).l2=0. So r must be the point of intersection of the lines l1 and l2.
Interestingly, we can find the point of intersection by taking the cross product of two lines.