docs: def get_missing_bins(original, trimmed):
      """Retrieve indices of a trimmed matrix with respect to the original matrix.
      Fairly fast but is only correct if diagonal values are different, which is
      always the case in practice.
      """
  
      original_diag = np.diag(original)
      trimmed_diag = np.diag(trimmed)
      index = []
      m = min(original.shape)
      for j in range(min(trimmed.shape)):
          k = 0
          while original_diag[j + k] != trimmed_diag[j] and k < 2 * m:
              k += 1
          index.append(k + j)
      return np.array(index)