def heaviest_interpolation(self, coords, n_jobs=-<NUM_LIT:1>):
li = self.linear_interpolation(coords, n_jobs=n_jobs)<EOL>mx = flattest(li.argmax(axis=<NUM_LIT:1>))<EOL>wh = np.where(~np.isclose(flattest(li.sum(axis=<NUM_LIT:1>)), <NUM_LIT:0>))[<NUM_LIT:0>]<EOL>mx = mx[wh]<EOL>return sps.csr_matrix((np.ones(len(mx)), (wh, mx)), shape=li.shape, dtype=np.int)<EOL>
mesh.heaviest_interpolation(x) yields an interpolation matrix for the given coordinate or 
  coordinate matrix x. The heaviest interpolation matrix is like a linear interpolation
  except that it assigns to each value the vertex with the maximum weight.