System and method for GPU-based 3D nonrigid registration

• US 7,483,034 B2
• Filed: 02/22/2005
• Issued: 01/27/2009
• Est. Priority Date: 02/25/2004
• Status: Expired due to Fees

• Current Assignee: Siemens Medical Solutions USA Incorporated (Siemens AG)
• Original Assignee: Siemens Medical Solutions USA Incorporated (Siemens AG)
• Inventors: Saddi, Kinda Anna, Chefd'hotel, Christophe
• Primary Examiner(s): Tung; Kee M
• Assistant Examiner(s): AMIN, JWALANT B
• Application Number:
  US11/062,962
• Publication Number:
  US 20050190189A1
• Time in Patent Office: 1,435 Days
• Field of Search: 345/582, 345/646, 382/300, 382294-295
• US Class Current: 345/582
• CPC Class Codes: G06T 15/005 General purpose rendering architectures
• CPC Class Codes: G06T 15/04 Texture mapping
• CPC Class Codes: G06T 7/30 Determination of transform parameters for the alignment of images, i.e. image registration

Abstract

A method of registering two images using a graphics processing unit includes providing a pair of images with a first and second image, calculating a gradient of the second image, initializing a displacement field on the grid point domain of the pair of images, generating textures for the first image, the second image, the gradient, and the displacement field, and loading said textures into the graphics processing unit. A pixel buffer is created and initialized with the texture containing the displacement field. The displacement field is updated from the first image, the second image, and the gradient for one or more iterations in one or more rendering passes performed by the graphics processing unit.

Claims

1. 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₁−

I₂∘

(Id+u_k))∇

I₂∘

(Id+u_k)]+α

Δ

u_k,

wherein u_k and u_k+1 are the displacement field and updated displacement field, respectively, I₁ is the first image, I₂ is the second image, ∇

I₂ 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).

2. The method of claim 1, wherein creating a pixel buffer further comprises creating a rendering context, initializing the pixel buffer attributes, and generating a new texture where the pixel buffer is bound.

3. The method of claim 1, further comprising the steps of setting texture parameters, and setting uniform variables.

4. The method of claim 1, further comprising creating a rendering context for the image registration, and compiling and loading fragment programs for updating the displacement field into the graphics processing unit.

5. The method of claim 1, further comprising expanding the resolution of the displacement field after reading the displacement field from the pixel buffer.

6. The method of claim 5, wherein expanding the resolution includes linearly interpolating displacement field values for points between grid points in the displacement field domain.

9. The method of claim 1, wherein the Laplacian of the displacement field U_k is defined in 3-dimensions as

Δ

u_k(x,y,z)=u_k(x+1,y,z)+u_k(x−

1,y,z)+u_k(x,y+1,z)+u_k(x,y−

1,z)+u_k(x,y,z+1)+u_k(x,y,z−

1)−

6.0×

u_k(x,y,x),

and in 2-dimensions as

Δ

u_k(p,q)=u_k(p+1,q)+u_k(p−

1,q)+u_k(p,q+1)+u_k(p,q−

1)+u_k(r,s)+u_k(t,u)−

60.0×

u_k(p,q),

wherein (x, y, z) is a 3-dimensional grid point, (p, q) is a 2-dimensional grid point, R_x is the number of columns, R_y is the number of rows, R_z is the number of slices, and N is an integer defined as N=√

{square root over (R_z)}, if the result is an integer, N={square root over (2×

R_z)}, otherwise,

10. The method of claim 1, further comprising smoothing the displacement field and reading the updated displacement field from the pixel buffer.

7. The method of claim I. wherein the domain of grid points is a 3-dimensional space. and wherein the displacement field is converted from 3 dimensions to 2 dimensions for the rendering pass according to the transformations

8. The method of claim 7, wherein the displacement field is converted from a 2-dimensional grid to a 3-dimensional grid according to the transformations

11. A method of registering two images, said method comprising the steps of:

providing a graphics processing unit;

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, including creating a rendering context, initializing the pixel buffer attributes, generating a new texture where the pixel buffer is bound and initializing it with the texture containing the displacement field;

updating the displacement field from a rendering pass performed by the graphics processing unit according to the equation

u_k+1=u_k+δ

[(I₁−

I₂∘

(Id+u_k))∇

I₂∘

(Id+u_k)]+α

Δ

u_k,

wherein u_k and u_k+1 are the displacement field and updated displacement field, respectively, I₁ is the first image, I₂ is the second image, ∇

I₂ is the gradient, Id is an identity map, Δ

U_k is the Laplacian of the displacement field, α

is a regularization parameter, and δ

is the time step, wherein the Laplacian of the displacement field u_k is defined in 3-dimensions as

Δ

u_k(X,Y,Z)=u_k(X+1,y,z)+u_k(X,y+1,Z)+U_k(X,Y=1,Z) +u_k(x, y, z+1)+u_k(x, y, z=1)=6.0×

u_k (x, y, x),

and in 2-dimensions as

Δ

u_k(p,q)=u_k(p+1,q)+u_k(t,u)=6.0×

u_k(p,q),

wherein (x, y, z) is a 3-dimensional grid point, (p,q) is a 2-dimensional grid point, R_x is the number of columns, R_y, is the number of rows, R_z is the number of slices, and N is an integer defined as N=√

{square root over (R_z)}, if the result is an integer, N=√

{square root over (2×

R_z)}, otherwise,

12. The method of claim 11, further comprising smoothing the displacement field and reading the updated displacement field from the pixel buffer.

13. A program storage device readable by a computer, tangibly embodying a program of instructions executable by the computer to perform the method steps for 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₁−

I₂∘

(Id+u_k))∇

I₂∘

(Id+u_k)]+α

Δ

u_k,

wherein u_k and u_k+1 are the displacement field and updated displacement field, respectively, I₁ is the first image, I₂ is the second image, ∇

I₂ 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).

14. The computer readable program storage device of claim 13, wherein creating a pixel buffer further comprises creating a rendering context, initializing the pixel buffer attributes, and generating a new texture where the pixel buffer is bound.

15. The computer readable program storage device of claim 14, wherein the method further comprising the steps of setting texture parameters, and setting uniform variables.

16. The computer readable program storage device of claim 13, wherein the method further comprises creating a rendering context for the image registration. and compiling and loading fragment programs for updating the displacement field into the graphics processing unit.

17. The computer readable program storage device of claim 13, wherein the method further comprises expanding the resolution of the displacement field after reading the displacement field from the pixel buffer.

18. The computer readable program storage device of claim 17, wherein expanding the resolution includes linearly interpolating displacement field values for points between grid points in the displacement field domain.

19. The computer readable program storage device of claim 14, wherein the domain of grid points is a 3-dimensional space, and wherein the displacement field is converted from 3 dimensions to 2 dimensions for the rendering pass according to the transformations

20. The computer readable program storage device of claim 19, wherein the displacement field is converted from a 2-dimensional grid to a 3-dimensional grid according to the transformations

21. The computer readable program storage device of claim 13, wherein the Laplacian of the displacement field u_k is defined in 3-dimensions as

Δ

u_k (X,Y,Z)=u_k(X+1,y,z)+u_k(X,y+1,Z)+U_k(X,Y=1,Z) +u_k(x, y, z+1)+u_k (x, y, z=1)=6.0×

u_k(x, y, x),

and in 2-dimensions as

Δ

u_k(p,q)=u_k(p+1,q)+u_k(p−

1,q)+u_k(p,q+1)+u_k(p,q−

1)+u_k(r,s)+u_k(t,u)−

60.0×

u_k(p,q),

wherein (x, y, z) is a 3-dimensional grid point, (p, q) is a 2-dimensional grid point, R_x is the number of columns, R_y is the number of rows, R_z is the number of slices, and N is an integer defined as N=√

{square root over (R_z)}, if the result is an integer, N=√

{square root over (2×

R_z)}, otherwise,

22. The computer readable program storage device of claim 13, further comprising smoothing the displacement field and reading the updated displacement field from the pixel buffer.