Building the ijk to RAS transform from a nrrd file

Hi Slicers,

I am working on a project with Slicer generated nrrd files and after looking at the nrrd header I notice there is no info about the ijk to RAS matrix yet when I load that nrrd file into Slicer I obtain the matrix. So I think I figured out how to build it in code and I am looking for validation from someone who is more experienced than me. This is what I am doing:

  1. I read the nrrd header and obtain the space direction vectors (say a,b,c)

  2. I build a 3x3 matrix (A) where the columns correspond to a, b, and c.

  3. Given the t translation column vector (as defined in the nrrd space origin header), I can obtain world coordinates x for the index column vector v by:

x = A*v + t

  1. since my nrrd volume appears to be in LPS frame (as per the space header), I have to multiply x by the LPS to RAS transform.

Say that’s T = [-1 0 0; 0 -1 0; 0 0 1], so:

x_ras = T * x

Is this correct?

Now, the IJK to RAS direction matrix reported in Slicer (Volumes module) seems to be:

T*A

is this correct?

Thanks!

Almost correct. IJKToRAS is a homogeneous transformation matrix, which contains axis directions in the top-left 3x3 sub-matrix, origin in the first 3 values of 4th column, and [0,0,0,1] in the 4th row (so you don’t “add” the translation component but concatenate). In Matlab syntax:

ijk_to_lps = [[space_directions, space_origin]; [0 0 0 1]]
lps_to_ras = diag([-1 -1 1 1])
ijk_to_ras = lps_to_ras * ijk_to_lps

Ok, the confusion was because in the volumes module the matrix appears as a 3x3 matrix. But it doesn’t matter, anyways the 4x4 form is equivalent as you show since you are operating in homogeneous coordinates.
Thank you Andras.

Hello Slicers

OS: Windows11
Slicer5.2.2

I’m beginner of 3DSlicer. I am working on a project with Slicer generated nrrd files too. This is what I want to do:

Step1. Forearm segmentation in Slicer and output as nrrd data.

Step2. Calculate the centre of gravity for a particular slice using python and obtain its coordinate ijk

import nrrd
tensor, header = nrrd.read("C:\\Users.....\\Segmentation.seg.nrrd", index_order ='C')
import cv2
def calc_center_gravity(tensor, gk):
    gravity = []
    img = tensor[gk]
    contours, _ = cv2.findContours(img, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    for cnt in contours:
        M = cv2.moments(cnt)
        if M["m00"] == 0:
            continue
        gx = int(M["m10"] / M["m00"])
        gy = int(M["m01"] / M["m00"])
        plt.plot(gi, gj, marker='.')
        gravity.append((gi, gj, gk))
    plt.imshow(img, cmap="gray")
    plt.title(f"centre of gravity index\n (i, j, k)=({gi}, {gj}, {gk})")
    plt.show()
    return gravity
gravity = calc_center_gravity(tensor, 24)

スクリーンショット 2023-06-25 145035

space_directions = header["space directions"]
space_origin = header["space origin"]
ijk2lps = np.vstack([np.hstack([space_directions, space_origin.reshape(-1, 1)]), np.array([0, 0, 0, 1])])
lps2ras = np.diag([-1, -1, 1, 1])
ijk2ras = np.dot(lps2ras, ijk2lps)
​
gravity_ijk = np.append(np.array(gravity[0]), 1).reshape(-1, 1)
gravity_ras = np.dot(ijk2ras, gravity_ijk)
gravity_ras

Step3. transform the acquired coordinate ijk to ras and markup it again in Slicer.

Problem
Cannot convert nrrd centre-of-gravity coordinates ijk to slicer coordinates ras between step2 and 3
Could you tell me what is wrong? Thank you!

You should be able to work from this example:
Script repository — 3D Slicer documentation

Thank you. I will have a try.