ailist.pyx

# ctypedef struct ailist_t:
#     int64_t nr, mr  					# Number of regions
#     interval_t *interval_list			# Regions data
#     uint32_t first, last				# Record range of intervals
#     int nc, lenC[10], idxC[10]
#     uint32_t *maxE
# ...
#     uint32_t binary_search(interval_t* As, uint32_t idxS, uint32_t idxE, uint32_t qe) nogil


@cython.boundscheck(False)
@cython.wraparound(False)
@cython.initializedcheck(False)
def ifa(self, const long[::1] starts, const long[::1] ends, const long[::1] index):

    # Check if object is still open
    if self.is_closed:
        raise NameError("AIList object has been closed.")

    if self.is_constructed == False:
        self.construct()

    cdef int length = len(starts)
    cdef int nfound = 0

    output_arr = np.ones(length, dtype=np.int64) * -1
    output_arr_other = np.ones(length, dtype=np.int64) * -1

    cdef int64_t [::1] output
    cdef int64_t [::1] output_other

    output = output_arr
    output_other = output_arr_other

    cdef int i = 0
    cdef int k = 0
    cdef int t = 0
    cdef int32_t qs, qe, qx, cs, ce
    cdef ailist_t *ail = self.interval_list

    for i in range(len(starts)):

        qs, qe, qx = starts[i], ends[i], index[i]

        # print("-----" * 5)
        # print("qs, qe, qx", qs, qe, qx)

        for k in range(ail.nc):
            # print("  k", k)

            cs = ail.idxC[k]
            ce = cs + ail.lenC[k]

            if ail.lenC[k] > 15:
                # print("    if")

                t = binary_search(ail.interval_list, cs, ce, qe)

                while (t >= cs and ail.maxE[t] > qs):
                    # print("      t", t)

                    if (ail.interval_list[t].end > qs):
                        # print("        ail.interval_list[t].end", ail.interval_list[t].end)

                        if nfound >= length:

                            length = (length + 1) * 2
                            output_arr = np.resize(output_arr, length)
                            output_arr_other = np.resize(output_arr_other, length)
                            output = output_arr
                            output_other = output_arr_other

                        output[nfound] = i
                        output_other[nfound] = ail.interval_list[t].index
                        nfound += 1

                        # print("        nfound", nfound)
                        # print("        output_arr", output_arr)
                        # print("        output_arr_other", output_arr_other)

                    t -= 1
            else:
                # print("    else")
                for t in range(cs, ce):
                    # print("      t", t)
                    # print("      ail.interval_list[t].start", ail.interval_list[t].start)
                    # print("      ail.interval_list[t].end", ail.interval_list[t].end)

                    if (ail.interval_list[t].start < qe and ail.interval_list[t].end > qs):

                        if nfound >= length:

                            length = (length + 1) * 2
                            output_arr = np.resize(output_arr, length)
                            output_arr_other = np.resize(output_arr_other, length)
                            output = output_arr
                            output_other = output_arr_other

                        output[nfound] = i
                        output_other[nfound] = ail.interval_list[t].index
                        # print("        nfound", nfound)
                        # print("        output_arr", output_arr)
                        # print("        output_arr_other", output_arr_other)

                        nfound += 1

    return output_arr[:nfound], output_arr_other[:nfound]