Find the intervals which have a non-empty intersection with a given interval The problem is the following: I have a list intervals which consists of tuples of the form (start, end) [with start &lt;= end]. Each tuple represents an interval (of the real line). We assume that the intervals in intervals are not overlapping each other. Given a new interval (s,e), I would like to write a Python function which checks if (s, e) overlaps any of the intervals in intervals. If (s, e) has a non-empty intersection with at least one of the intervals in intervals, the function should return the indices of these intervals in the list intervals. 

Say that the function is called find_intersections. Then, given that intervals = [(1, 3.5), (5.5, 8.7), (10.2, 22.6), (22.7, 23.1)], expected outputs would be: 


find_intersection(intervals, (3.2, 5.)) returns array([0])
find_intersection(intervals, (6.1, 7.3)) returns array([1])
find_intersection(intervals, (9.1, 10.2)) returns No intersection.
find_intersection(intervals, (5.8, 22.9)) returns array([1, 2, 3]).


The code for find_intersection I have written is:

import itertools

def find_intersection(intervals, new_interval):
    _times = sorted(list(itertools.chain.from_iterable(intervals)))
    ind = np.searchsorted(_times, np.asarray(new_interval))
    parity = np.mod(ind, 2)
    if (not np.any(parity)) and ind[1] == ind[0]:
        print('No intersection.')
    elif parity[0] == 1:
        ub = ind[1] if parity[1] == 1 else ind[1] - 1
        return np.arange((ind[0] - 1) / 2, (ub - 1) / 2 + 1)
    elif parity[1] == 1:
        lb = ind[0] if parity[0] == 1 else ind[0] + 1
        return np.arange((lb - 1) / 2, (ind[1] - 1) / 2 + 1)
    else:
        lb = ind[0] if parity[0] == 1 else ind[0] + 1
        ub = ind[1] if parity[1] == 1 else ind[1] - 1
        return np.arange((lb - 1) / 2, (ub - 1) / 2 + 1)


I believe that the code does the job.

Is there an easier/smarter way to address this problem?

Edit: this question was cross-posted on stackoverflow.com.
1. Review


There is no docstring. What does the function find_intersections do? The text in the post is pretty clear, so you could use this as the basis of a docstring. I would make it explicit that the intervals are open (if that's what you really intend), since in Python intervals are normally half-open.
There seems to be some uncertainty about whether the input intervals is supposed to be already sorted. On the one hand, if it's already sorted then there's no point in calling sorted. But on the other, if it's not sorted then returning a range of indexes is not very useful, as these are indexes into _times, not into intervals.

I recommend picking one possibility or the other: if intervals is supposed to be sorted, then don't sort it (you might assert that it is already sorted); but if not, then translate the returned indexes so that they refer to intervals.
If there is no intersection, find_intersections prints a message. But if the caller is already printing some output, it would be inconvenient (at best) to have this kind of message appear in the middle of the output. Better to leave the printing of messages (or not) up to the caller.
If there is an intersection, find_intersections returns a range of indexes. But if there are no intersection, find_intersections returns None. It is not a good idea to return an exceptional value like this, because it would be easy for the caller to forget to check for the exceptional value, and to write something like this:

result = find_intersection(intervals, new_interval)
print(f"Number of intersections = {len(result)}.")


which would fail with:

TypeError: object of type 'NoneType' has no len()


What they should have written is:

result = find_intersection(intervals, new_interval)
if result is None:
    print("Number of intersections = 0.")
else:
    print(f"Number of intersections = {len(result)}.")


but this is very long-winded. In this case a simple improvement would be to return an empty range of indexes.
If you're going to use NumPy, then you might as well use it wherever possible. So instead of itertools.chain.from_iterable, use numpy.reshape; instead of sorted, use numpy.sort; and so on.


2. Revised code if intervals are unsorted

If intervals is unsorted (as suggested by the call to sorted), then there's no point in sorting it. Sorting costs \$O(n\log n)\$ whereas comparing the query interval with every element of intervals only costs \$O(n)\$.

So in this case we can take advantage of the fact that the open interval \$(a, b)\$ overlaps with the open interval \$(c, d)\$ if and only if \$a &lt; d\$ and \$c &lt; b\$, and write:

def find_intersection(intervals, query):
    """Find intersections between intervals.
    Intervals are open and are represented as pairs (lower bound,
    upper bound).

    Arguments:
    intervals: array_like, shape=(N, 2) -- Array of intervals.
    query: array_like, shape=(2,) -- Interval to query.

    Returns:
    Array of indexes of intervals that overlap with query.

    """
    intervals = np.asarray(intervals)
    lower, upper = query
    return np.argwhere((lower &lt; intervals[:, 1]) &amp; (intervals[:, 0] &lt; upper))


3. Revised code if intervals are sorted

If intervals is already sorted, as suggested by returning indexes into the sorted array, then take advantage of the side argument to numpy.searchsorted and the floor division operator, //:

def find_intersection(intervals, query):
    """Find intersections between intervals.
    Intervals are open and are represented as pairs (lower bound,
    upper bound).

    Arguments:
    intervals: array_like, shape=(N, 2) -- Array of intervals.
    query: array_like, shape=(2,) -- Interval to query.

    Returns:
    Array of indexes of intervals that overlap with query.

    """
    endpoints = np.reshape(intervals, -1)
    lower, upper = query
    i = np.searchsorted(endpoints, lower, side='right')
    j = np.searchsorted(endpoints, upper, side='left')
    return np.arange(i // 2, (j + 1) // 2)