title: Increasing efficiency of barycentric coordinate calculation in python
body: <p>Background: I'm attempting to warp one face to another of a different shape.</p>
  
  <p>In order to warp one image to another, I'm using a delaunay triangulation of facial landmarks and warping the triangles of one portrait to the corresponding triangles of the second portrait. I'm using a barycentric coordinate system to map a point within a triangle to its corresponding warped location on the other triangle.</p>
  
  <p>My first approach was to solve the system Ax = b with the inverse multiplication method, where A consists of the three corners of the triangle, b represents the current point, and x represents the barycentric coordinates of this point (alpha, beta, and gamma). I found the inverse of matrix A once per triangle, and then for every point within that triangle calculated the barycentric coordinates by finding the dot product of A^-1 and the point b. I found this to be very slow (the function takes 36 seconds to complete).</p>
  
  <p>Following the recommendation of other posts, I attempted to use a least squares solution to improve the efficiency of this process. However, the time increased to 154 seconds when I used numpy's lsq method. I believe this is due to the fact that the A matrix is factored every single time the inside loop runs, while before I was able to find the inverse only one time, before the two loops begin.</p>
  
  <p>My question is, how can I improve the efficiency of this function? Is there a way to store the factorization of A so that each time the least squares solution is calculated for a new point, it isn't repeating the same work?</p>
  
  <p>Pseudocode for this function:</p>
  
  <pre><code># Iterate through each triangle (and get corresponding warp triangle)
  for triangle in triangulation:
  
      # Extract corners of the unwarped triangle
      a = firstCornerUW
      b = secondCornerUW
      c = thirdCornerUW
  
      # Extract corners of the warp triangle
      a_prime = firstCornerW
      b_prime = secondCornerW
      c_prime = thirdCornerW
  
      # This matrix will be the same for all points within the triangle
      triMatrix = matrix of a, b, and c
  
      # Bounding box of the triangle
      xleft = min(ax, bx, cx)
      xright = max(ax, bx, cx)
      ytop = min(ay, by, cy)
      ybottom = max(ay, by, cy)
  
      for x in range(xleft, xright):
  
          for y in range(ytop, ybottom):
  
              # Store the current point as a matrix
              p = np.array([[x], [y], [1]])
  
              # Solve for least squares solution to get barycentric coordinates
              barycoor = np.linalg.lstsq(triMatrix, p)
  
              # Pull individual coordinates from the array
              alpha = barycoor[0]
              beta = barycoor[1]
              gamma = barycoor[2]
  
              # If any of these conditions are not met, the point is not inside the triangle
              if alpha, beta, gamma &gt; 0 and alpha + beta + gamma &lt;= 1:
  
                  # Now calculate the warped point by multiplying by alpha, beta, and gamma
                  # Warp the point from image to warped image
  </code></pre>