Python 2D z-order curve

zorder.py

def zorder(N):
    if N == 0:
        return np.zeros((0, 0), dtype=np.int)
    elif N == 1:
        return np.zeros((1, 2), dtype=np.int)
    else:
        pnts = zorder(N - 1)
        l = pnts.shape[0]
        offs = 2 ** (N - 2)
        res = np.zeros((4 * l, 2), dtype=np.int)
        for i in range(2):
            for j in range(2):
                i0, i1 = (2 * i + j) * l, (2 * i + j + 1) * l
                res[i0:i1, 0] = pnts[:, 0] + j * offs
                res[i0:i1, 1] = pnts[:, 1] + i * offs
        return res