docs: def cartesian_to_barycentric_3D(tri, xy):
      '''
      cartesian_to_barycentric_3D(tri,xy) is identical to cartesian_to_barycentric_2D(tri,xy) except
      it works on 3D data. Note that if tri is a 3 x 3 x n, a 3 x n x 3 or an n x 3 x 3 matrix, the
      first dimension must always be the triangle vertices and the second 3-sized dimension must be
      the (x,y,z) coordinates.
      '''
      xy = np.asarray(xy)
      tri = np.asarray(tri)
      if len(xy.shape) == 1:
          return cartesian_to_barycentric_3D(np.transpose(np.asarray([tri]), (1,2,0)),
                                             np.asarray([xy]).T)[:,0]
      xy = xy if xy.shape[0] == 3 else xy.T
      if tri.shape[0] == 3:
          tri = tri if tri.shape[1] == 3 else np.transpose(tri, (0,2,1))
      elif tri.shape[1] == 3:
          tri = tri.T if tri.shape[0] == 3 else np.transpose(tri, (1,2,0))
      elif tri.shape[2] == 3:
          tri = np.transpose(tri, (2,1,0) if tri.shape[1] == 3 else (2,0,1))
      if tri.shape[0] != 3 or tri.shape[1] != 3:
          raise ValueError('Triangle array did not have dimensions of sizes 3 and 3')
      if xy.shape[0] != 3:
          raise ValueError('coordinate matrix did not have a dimension of size 3')
      if tri.shape[2] != xy.shape[1]:
          raise ValueError('number of triangles and coordinates must match')
      # The algorithm here is borrowed from this stack-exchange post:
      # http://gamedev.stackexchange.com/questions/23743
      # in which it is attributed to Christer Ericson's book Real-Time Collision Detection.
      v0 = tri[1] - tri[0]
      v1 = tri[2] - tri[0]
      v2 = xy - tri[0]
      d00 = np.sum(v0 * v0, axis=0)
      d01 = np.sum(v0 * v1, axis=0)
      d11 = np.sum(v1 * v1, axis=0)
      d20 = np.sum(v2 * v0, axis=0)
      d21 = np.sum(v2 * v1, axis=0)
      den = d00*d11 - d01*d01
      zero = np.isclose(den, 0)
      unit = 1 - zero
      den += zero
      l2 = unit * (d11 * d20 - d01 * d21) / den
      l3 = unit * (d00 * d21 - d01 * d20) / den
      return np.asarray([1.0 - l2 - l3, l2])