Patent ID: 7483034

Claim:
A method of registering two images using a graphics processing unit, said method comprising the steps of: providing a pair of images with a first image and a second image, wherein said images comprise a plurality of intensities corresponding to a domain of grid points in a 3-dimensional space; calculating a gradient of the second image; initializing a displacement field on the grid point domain of the pair of images, wherein said displacement field transforms said second image into said first image; generating textures for the first image, the second image, the gradient, and the displacement field, and loading said textures into the graphics processing unit; creating a pixel buffer and initializing it with the texture containing the displacement field; and updating the displacement field in a rendering pass performed by the graphics processing unit according to the equation u k+1 =u k +δ[( I 1 −I 2 ∘(Id+ u k ))∇ I 2 ∘(Id+ u k )]+αΔ u k , wherein u k and u k+1 are the displacement field and updated displacement field, respectively, I 1 is the first image, I 2 is the second image, ∇ I 2 is the gradient, Id is the identity map, ∇u k is the Laplacian of the displacement field, α is a regularization parameter, and δ is the time step for one or more iterations, wherein the boundary conditions of the 3-dimensional displacement field are clamp-to-edge, and further comprising decomposing the 3-dimensional displacement field into a 2-dimensional displacement field and controlling said boundary conditions of the 2-D displacement field by imposing a 2D mask, wherein said mask is a precomputed vector comprising 7 coefficients (a,b,c,d,e,f,g) that multiply each term of the displacement field according to the equation, a×u k (x+1,y)+b×u k (x−1,y)+c×u k (x,y+1)+d×u k (x,y−1)+e×u k (r,s)+f×u k (t,u)+g×u k (x,y), wherein u k (x,y) is the 2D displacement field, (a,b,c,d,e,f) are either 0 or 1, wherein 0 indicates a displacement field value outside the boundary, and 1 a point inside the boundary, g can vary between −6 and 0 depending on the number of nearest neighbor interior points about (x,y).