3D Slicer:从入门到精通(四)——坐标系统
文章目录
坐标系统
引言
在处理医学图像和应用时,一个问题是坐标系统之间的差异。在成像应用中通常使用三种坐标系统:可以区分世界、解剖和图像坐标系统。
下图示意了这三个空间及其对应的轴:
每个坐标系统都有其特定的用途,并以特定的方式表示其数据。
解剖图像基于My MS组织分享的图像。
注意,斯坦福的Chand John创建了一个详细介绍Slicer中坐标处理方式的演示。
世界坐标系统
世界坐标系统通常是一个笛卡尔坐标系统,模型(例如MRI扫描仪或患者)被定位在其中。每个模型都有自己的坐标系统,但只有一个世界坐标系统来定义每个模型的位置和方向。
解剖坐标系统
对于医学成像技术来说,最重要的模型坐标系统是解剖空间(也称为患者坐标系统)。这个空间由三个平面组成,用以描述人体的标准解剖位置:
- 轴向平面平行于地面,将头部(上方)与脚部(下方)分开。
- 冠状平面垂直于地面,将前部(前方)与后部(后方)分开。
- 矢状平面垂直于地面,将左侧与右侧分开。
由这些平面可知,所有轴都在正方向有其标记(例如,负上方轴由下方轴表示)。
解剖坐标系统是一个连续的三维空间,在该空间中采样了图像。在神经成像中,通常根据被扫描的人的解剖轴定义这个空间。因此,3D基础沿着前后、上下和左右的解剖轴定义。
然而,不同的医学应用使用不同的这个3D基础的定义。最常见的基础有:
- LPS(左、后、上)用于DICOM图像
L P S = { 从右向左 从前向后 从下向上 } LPS = \begin{Bmatrix} \text{从右向左} \\ \text{从前向后} \\ \text{从下向上} \end{Bmatrix} LPS=⎩ ⎨ ⎧从右向左从前向后从下向上⎭ ⎬ ⎫
- RAS(右、前、上)与LPS相似,但前两个轴颠倒
R A S = { 从左向右 从后向前 从下向上 } RAS = \begin{Bmatrix} \text{从左向右} \\ \text{从后向前} \\ \text{从下向上} \end{Bmatrix} RAS=⎩ ⎨ ⎧从左向右从后向前从下向上⎭ ⎬ ⎫
因此,这两种约定之间的唯一区别是前两个坐标的符号被反转。
这两种基础都同样有用且合逻辑。只需知道图像是参照哪个基础的。
图像坐标系统
图像坐标系统描述了图像相对于解剖结构的获取方式。医疗扫描仪创建的是规则的、矩形的点和单元阵列,从左上角开始。 i i i轴向右增加, j j j轴向下增加, k k k轴向后增加。
除了每个体素 ( i j k ) (ijk) (ijk)的强度值外,解剖坐标的原点和间距也被存储。
- 原点表示第一个体素(0,0,0)在解剖坐标系统中的位置,例如(100, 50, -25)毫米。
- 间距表示体素之间的距离,例如,(1.5, 0.5, 0.5)毫米。
以下2D示例展示了原点和间距的含义:
利用原点和间距,可以计算出每个(图像坐标)体素在解剖坐标中对应的位置。
图像变换
从图像空间向量 ( i j k ) ′ (i j k)' (ijk)′到解剖空间向量 x ⃗ \vec{x} x的变换是一种仿射变换,由一个线性变换 A \mathbf{A} A和一个平移 t ⃗ \vec{t} t组成。
x ⃗ = A ( i j k ) ′ + t ⃗ \vec{x} = A \begin{pmatrix} i & j & k \end{pmatrix}' + \vec{t} x=A(ijk)′+t
变换矩阵 A \mathbf{A} A是一个 3 × 3 3 \times 3 3×3矩阵,包含所有关于空间方向和轴向缩放的信息。
t ⃗ \vec{t} t是一个 3 × 1 3 \times 1 3×1向量,包含第一个体素的几何位置信息。
( x 1 x 2 x 3 ) = ( A 11 A 12 A 13 A 21 A 22 A 23 A 31 A 32 A 33 ) ( i j k ) + ( t 1 t 2 t 3 ) \begin{pmatrix} x_1 \\ x_2 \\ x_3 \end{pmatrix} = \begin{pmatrix} A_{11} & A_{12} & A_{13} \\ A_{21} & A_{22} & A_{23} \\ A_{31} & A_{32} & A_{33} \end{pmatrix} \begin{pmatrix} i \\ j \\ k \end{pmatrix} + \begin{pmatrix} t_1 \\ t_2 \\ t_3 \end{pmatrix} x1x2x3 = A11A21A31A12A22A32A13A23A33 ijk + t1t2t3
最后一个等式表明,通过矩阵乘法进行线性变换,通过向量加法进行平移。为了通过矩阵乘法同时表示变换和平移,必须使用增广矩阵。这种技术要求在 A \mathbf{A} A矩阵的底部增加一行零,在右侧增加一个额外的列—平移向量,以及在右下角增加一个 1 1 1。此外,所有向量都必须表示为齐次坐标,这意味着在末尾增加一个 1 1 1。
( x 1 x 2 x 3 1 ) = ( A 11 A 12 A 13 t 1 A 21 A 22 A 23 t 2 A 31 A 32 A 33 t 3 0 0 0 1 ) ( i j k 1 ) \begin{pmatrix} x_1 \\ x_2 \\ x_3 \\ 1 \end{pmatrix} = \begin{pmatrix} A_{11} & A_{12} & A_{13} & t_1 \\ A_{21} & A_{22} & A_{23} & t_2 \\ A_{31} & A_{32} & A_{33} & t_3 \\ 0 & 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} i \\ j \\ k \\ 1 \end{pmatrix} x1x2x31 = A11A21A310A12A22A320A13A23A330t1t2t31 ijk1
根据使用的解剖空间(LPS或RAS), 4 × 4 4 \times 4 4×4矩阵被称为IJKtoLPS或IJKtoRAS矩阵,因为它表示从IJK到LPS或RAS的变换。
计算IJtoLS-矩阵的2D示例
以下图像显示了左侧的解剖空间与L§S基础,以及右侧相应的图像坐标。
原点(解剖空间中第一个像素的坐标)是(50, 300)mm,间距(两个像素之间的距离)是(50, 50)mm。
因为这是一个2D示例, A \mathbf{A} A是一个 2 × 2 2 \times 2 2×2矩阵, t ⃗ \vec{t} t是一个 2 × 1 2 \times 1 2×1向量。因此,仿射变换的方程为:
( L S 1 ) = ( A 11 A 12 t 1 A 21 A 22 t 2 0 0 1 ) ( i j 1 ) \begin{pmatrix} L \\ S \\ 1 \end{pmatrix} = \begin{pmatrix} A_{11} & A_{12} & t_1 \\ A_{21} & A_{22} & t_2 \\ 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} i \\ j \\ 1 \end{pmatrix} LS1 = A11A210A12A220t1t21 ij1
通过乘以IJtoLS矩阵,我们可以将图像坐标转换为解剖坐标。
这种变换允许我们在图像处理和医学成像分析中精确地定位和量化结构。
通过将IJtoLS矩阵和右侧的向量相乘,可以得到以下乘积:
最后一个等式和矩阵乘积显示,共有6个未知变量
(
A
11
,
A
12
,
A
21
,
A
22
,
t
1
,
t
2
)
(A_{11}, A_{12}, A_{21}, A_{22}, t_1, t_2)
(A11,A12,A21,A22,t1,t2)需要确定。然而,原点和间距的知识允许图像和解剖空间之间存在以下关系:
( L S ) ≡ ( i j ) ( 50 300 ) ≡ ( 0 0 ) ( 100 300 ) ≡ ( 1 0 ) ( 50 250 ) ≡ ( 0 1 ) … \begin{pmatrix} L \\ S \end{pmatrix} \equiv \begin{pmatrix} i \\ j \end{pmatrix} \qquad \begin{pmatrix} 50 \\ 300 \end{pmatrix} \equiv \begin{pmatrix} 0 \\ 0 \end{pmatrix} \qquad \begin{pmatrix} 100 \\ 300 \end{pmatrix} \equiv \begin{pmatrix} 1 \\ 0 \end{pmatrix} \qquad \begin{pmatrix} 50 \\ 250 \end{pmatrix} \equiv \begin{pmatrix} 0 \\ 1 \end{pmatrix} \qquad \dots (LS)≡(ij)(50300)≡(00)(100300)≡(10)(50250)≡(01)…
因此,至少可以得到六个方程式:
50 = A 11 ⋅ 0 + A 12 ⋅ 0 + t 1 ⋅ 1 300 = A 21 ⋅ 0 + A 22 ⋅ 0 + t 2 ⋅ 1 100 = A 11 ⋅ 1 + A 12 ⋅ 0 + t 1 ⋅ 1 300 = A 21 ⋅ 1 + A 22 ⋅ 0 + t 2 ⋅ 1 50 = A 11 ⋅ 0 + A 12 ⋅ 1 + t 1 ⋅ 1 250 = A 21 ⋅ 0 + A 22 ⋅ 1 + t 2 ⋅ 1 \begin{align*} 50 &= A_{11} \cdot 0 + A_{12} \cdot 0 + t_1 \cdot 1 \\ 300 &= A_{21} \cdot 0 + A_{22} \cdot 0 + t_2 \cdot 1 \\ 100 &= A_{11} \cdot 1 + A_{12} \cdot 0 + t_1 \cdot 1 \\ 300 &= A_{21} \cdot 1 + A_{22} \cdot 0 + t_2 \cdot 1 \\ 50 &= A_{11} \cdot 0 + A_{12} \cdot 1 + t_1 \cdot 1 \\ 250 &= A_{21} \cdot 0 + A_{22} \cdot 1 + t_2 \cdot 1 \end{align*} 5030010030050250=A11⋅0+A12⋅0+t1⋅1=A21⋅0+A22⋅0+t2⋅1=A11⋅1+A12⋅0+t1⋅1=A21⋅1+A22⋅0+t2⋅1=A11⋅0+A12⋅1+t1⋅1=A21⋅0+A22⋅1+t2⋅1
如上所述,平移 t ⃗ \vec{t} t包含了第一个像素的几何位置信息,因此等同于原点。这个结果也被第一个方程式证实。
其他方程式的解决方案导致以下IJtoLS矩阵:
I J t o L S = ( 50 0 50 0 − 50 300 0 0 1 ) IJtoLS = \begin{pmatrix} 50 & 0 & 50 \\ 0 & -50 & 300 \\ 0 & 0 & 1 \end{pmatrix} IJtoLS= 50000−500503001
在使用R(A)S基准的情况下,只需翻转解剖空间的左和前轴,图像坐标系统就会以与L(P)S情况相同的方式出现。
对于这个2D示例,IJtoRS矩阵将是:
I J t o R S = ( − 50 0 250 0 − 50 300 0 0 1 ) IJtoRS = \begin{pmatrix} -50 & 0 & 250 \\ 0 & -50 & 300 \\ 0 & 0 & 1 \end{pmatrix} IJtoRS= −50000−5002503001
这个矩阵与IJtoLS矩阵非常相似,但有2个不同之处:
- 由于原点不同,平移 t ⃗ \vec{t} t发生了变化。
- 右轴被翻转,所以IJtoRS矩阵的第一列只是符号相反。
Slicer中的坐标系统约定
DICOM和大多数医学影像软件使用LPS坐标系统来存储所有数据。因为只有相对差异有意义,所以原点的选择是任意的,没有通用标准,但通常设置为成像系统的某个几何中心,或者选择靠近感兴趣对象的中心。
早在2000年代初期,LPS和RAS两种坐标系统就已广泛使用,当时Slicer的开发者选择使用RAS坐标系统。历史上,GE设备扫描使用RAS,而西门子和其他公司使用LPS。由于几位GE研究员是Slicer的早期贡献者,因此采用了RAS作为内部表示。
Slicer仍然内部使用RAS坐标系统来存储坐标值对于所有数据类型,但为了与其他软件兼容,它假设所有文件中的数据都存储在LPS坐标系统中(除非文件中的坐标系统明确声明为RAS)。为了实现这一点,每当Slicer读取或写入一个文件时,它可能需要翻转前两个坐标轴的符号,以将数据转换为RAS坐标系统。
与其他软件/约定的关系
ITK
ITK 使用LPS约定。
使用MATLAB将Slicer RAS坐标(例如,标志点)映射到NIfTI图像的体素空间
要从NIFTI文件头部提取“体素到世界”变换矩阵(条目:qto_xyz:1-4
)在MATLAB中:
d = inv(M) * [ R A S 1 ]'
其中M
是矩阵,R A S
是Slicer中的坐标,然后d
给出一个体素坐标的向量。
(解决方案由Debrecen大学的András Jakab提供)
参考资料
- https://people.cs.uchicago.edu/~glk/unlinked/nrrd-iomf.pdf
- http://www.grahamwideman.com/gw/brain/orientation/orientterms.htm
- https://nifti.nimh.nih.gov/nifti-1/documentation/faq
- https://teem.sourceforge.net/nrrd/format.html
- DICOM 2013 PS3.3 图像位置和图像方向