docs: def convex_hull_ijv(pixel_labels, indexes, fast=True):
      '''Return the convex hull for each label using an ijv labeling
      
      pixel_labels: the labeling of the pixels in i,j,v form where
                    i & j are the coordinates of a pixel and v is
                    the pixel's label number
      indexes: the indexes at which to measure the convex hull
  
      Returns a matrix and a vector. The matrix consists of one row per
      point in the convex hull. Each row has three columns, the label #,
      the i coordinate of the point and the j coordinate of the point. The
      result is organized first by label, then the points are arranged
      counter-clockwise around the perimeter.
      The vector is a vector of #s of points in the convex hull per label
      '''
      if fast:
          return _convex_hull.convex_hull_ijv(pixel_labels, indexes)
  
      # We keep this code for testing.
      if len(indexes) == 0:
          return np.zeros((0,3),int),np.zeros((0,),int)
      #
      # An array that converts from label # to index in "indexes"
      anti_indexes = np.zeros((np.max(indexes)+1,),int)
      anti_indexes[indexes] = list(range(len(indexes)))
  
      coords = pixel_labels[:,:2]
      i = coords[:, 0]
      j = coords[:, 1]
      # This disgusting copy spooge appears to be needed for scipy 0.7.0
      labels_per_point = np.zeros(len(pixel_labels), int)
      labels_per_point[:] = pixel_labels[:,2]
      #
      # Calculate the centers for each label
      #
      center_i = fixup_scipy_ndimage_result(
          scind.mean(i.astype(float), labels_per_point, indexes))
      center_j = fixup_scipy_ndimage_result(
          scind.mean(j.astype(float), labels_per_point, indexes))
      centers = np.column_stack((center_i, center_j))
      #
      # Now make an array with one outline point per row and the following
      # columns:
      #
      # index of label # in indexes array
      # angle of the point relative to the center
      # i coordinate of the point
      # j coordinate of the point
      #
      anti_indexes_per_point = anti_indexes[labels_per_point]
      centers_per_point = centers[anti_indexes_per_point]
      angle = np.arctan2(i-centers_per_point[:,0],j-centers_per_point[:,1])
      a = np.zeros((len(i),3), np.int32)
      a[:,0] = anti_indexes_per_point
      a[:,1:] = coords
      #
      # Sort the array first by label # (sort of), then by angle
      #
      order = np.lexsort((angle,anti_indexes_per_point))
      #
      # Make unique
      #
      same_as_next = np.hstack([np.all(a[order[:-1],:] == a[order[1:],:], 1), [False]])
      order = order[~same_as_next]
      a=a[order]
      anti_indexes_per_point = anti_indexes_per_point[order]
      angle = angle[order]
      centers_per_point = centers_per_point[order]
      #
      # Make the result matrix, leaving enough space so that all points might
      # be on the convex hull.
      #
      result = np.zeros((len(order),3), np.int32)
      result[:,0] = labels_per_point[order]
      #
      # Create an initial count vector
      #
      v = np.ones((a.shape[0],),dtype = np.int32)
      result_counts = scipy.sparse.coo_matrix((v,(a[:,0],v*0)),
                                              shape=(len(indexes),1))
      result_counts = result_counts.toarray().flatten()
      r_anti_indexes_per_point = anti_indexes_per_point # save this
      #
      # Create a vector that indexes into the results for each label
      #
      result_index = np.zeros(result_counts.shape, np.int32)
      result_index[1:]=np.cumsum(result_counts[:-1])
      #
      # Initialize the counts of convex hull points to a ridiculous number
      #
      counts = np.iinfo(np.int32).max
      first_pass = True
      while True:
          #
          # Figure out how many putative convex hull points there are for
          # each label.
          #
          # If the count for a label is 3 or less, it's a convex hull or
          # degenerate case.
          #
          # If the count hasn't changed in an iteration, then we've done
          # as well as we can hope to do.
          #
          v = np.ones((a.shape[0],),dtype = np.int32)
          new_counts = scipy.sparse.coo_matrix((v,(a[:,0],v*0)),
                                               shape=(len(indexes),1))
          new_counts = new_counts.toarray().flatten()
          done_count = (2 if first_pass else 3)
          finish_me = ((new_counts > 0) & 
                       ((new_counts <= done_count) | 
                        (new_counts == counts)))
          indexes_to_finish = np.argwhere(finish_me).astype(np.int32).squeeze()
          keep_me = (new_counts > done_count) & (new_counts < counts)
          indexes_to_keep = np.argwhere(keep_me).astype(np.int32).squeeze()
          if len(indexes_to_finish):
              result_counts[finish_me] = new_counts[finish_me]
              #
              # Store the coordinates of each of the points to finish
              #
              finish_this_row = finish_me[a[:,0]]
              rows_to_finish = np.argwhere(finish_this_row).flatten()
              a_to_finish = a[rows_to_finish]
              atf_indexes = a_to_finish[:,0]
              #
              # Map label #s to the index into indexes_to_finish of that label #
              #
              anti_indexes_to_finish = np.zeros((len(indexes),), np.int32)
              anti_indexes_to_finish[indexes_to_finish] = list(range(len(indexes_to_finish)))
              #
              # Figure out the indices of each point in a label to be finished.
              # We figure out how much to subtract for each label, then
              # subtract that much from 0:N to get successive indexes at
              # each label.
              # Then we add the result_index to figure out where to store it
              # in the result table.
              #
              finish_idx_base = np.zeros((len(indexes_to_finish),), np.int32)
              finish_idx_base[1:]=np.cumsum(new_counts[indexes_to_finish])[:-1]
              finish_idx_bases = finish_idx_base[anti_indexes_to_finish[atf_indexes]]
              finish_idx = (np.array(list(range(a_to_finish.shape[0])))-
                            finish_idx_bases)
              finish_idx = finish_idx + result_index[atf_indexes]
              result[finish_idx,1:] = a_to_finish[:,1:]
          if len(indexes_to_keep) == 0:
              break
          #
          # Figure out which points are still available
          #
          rows_to_keep = np.argwhere(keep_me[a[:,0].astype(np.int32)]).flatten()
          rows_to_keep = rows_to_keep.astype(np.int32)
          a = a[rows_to_keep]
          centers_per_point = centers_per_point[rows_to_keep]
          counts = new_counts
          #
          # The rule is that the area of the triangle from the center to
          # point N-1 to point N plus the area of the triangle from the center
          # to point N to point N+1 must be greater than the area of the
          # triangle from the center to point N-1 to point N+1 for a point
          # to be on the convex hull.
          # N-1 and N+1 have to be modulo "counts", so we make special arrays
          # to address those situations.
          #
          anti_indexes_to_keep = np.zeros((len(indexes),), np.int32)
          anti_indexes_to_keep[indexes_to_keep] = list(range(len(indexes_to_keep)))
          idx_base = np.zeros((len(indexes_to_keep),), np.int32)
          idx_base[1:]=np.cumsum(counts[keep_me])[0:-1]
          idx_bases = idx_base[anti_indexes_to_keep[a[:,0]]]
          counts_per_pt = counts[a[:,0]]
          idx = np.array(list(range(a.shape[0])), np.int32)-idx_bases
          n_minus_one = np.mod(idx+counts_per_pt-1,counts_per_pt)+idx_bases
          n_plus_one  = np.mod(idx+1,counts_per_pt)+idx_bases
          #
          # Compute the triangle areas
          #
          t_left_plus_right = (triangle_areas(centers_per_point,
                                              a[n_minus_one, 1:],
                                              a[:, 1:]) +
                               triangle_areas(centers_per_point,
                                              a[:, 1:],
                                              a[n_plus_one, 1:]))
          t_lr = triangle_areas(centers_per_point,
                                a[n_minus_one, 1:], a[n_plus_one, 1:])
          #
          # Keep the points where the area of the left triangle plus the
          # area of the right triangle is bigger than the area of the triangle
          # composed of the points to the left and right. This means that
          # there's a little triangle sitting on top of t_lr with our point
          # on top and convex in relation to its neighbors.
          #
          keep_me = t_left_plus_right > t_lr
          #
          # If all points on a line are co-linear with the center, then the
          # whole line goes away. Special handling for this to find the points
          # most distant from the center and on the same side
          #
          consider_me = t_left_plus_right == 0
          del t_lr, t_left_plus_right
          if np.any(consider_me):
              diff_i = a[:,1]-centers_per_point[:,0]
              diff_j = a[:,2]-centers_per_point[:,1]
              #
              # The manhattan distance is good enough
              #
              dist = np.abs(diff_i)+np.abs(diff_j)
              # The sign is different on different sides of a line including
              # the center. Multiply j by 2 to keep from colliding with i
              #
              # If both signs are zero, then the point is in the center
              #
              # Compute np.sign(diff_i, out=diff_i) + np.sign(diff_j) * 2,
              # but without any temporaries.
              sign = np.sign(diff_j, diff_j)
              sign *= 2
              sign += np.sign(diff_i)
              del diff_i, diff_j
  
              n_minus_one_consider = n_minus_one[consider_me]
              n_plus_one_consider = n_plus_one[consider_me]
              left_is_worse = ((dist[consider_me] > dist[n_minus_one_consider]) |
                               (sign[consider_me] != sign[n_minus_one_consider]))
              right_is_worse = ((dist[consider_me] > dist[n_plus_one_consider]) |
                                (sign[consider_me] != sign[n_plus_one_consider]))
              to_keep = left_is_worse & right_is_worse & (sign[consider_me] != 0)
              keep_me[consider_me] = to_keep
              del dist, sign
          a = a[keep_me,:]
          centers_per_point = centers_per_point[keep_me]
          del keep_me
          first_pass = False
      #
      # Finally, we have to shrink the results. We number each of the
      # points for a label, then only keep those whose indexes are
      # less than the count for their label.
      #
      within_label_index = np.array(list(range(result.shape[0])), np.int32)
      counts_per_point = result_counts[r_anti_indexes_per_point]
      result_indexes_per_point = result_index[r_anti_indexes_per_point] 
      within_label_index = (within_label_index - result_indexes_per_point)
      result = result[within_label_index < counts_per_point]
      return result, result_counts