voxelmorph / voxelmorph

Unsupervised Learning for Image Registration
Apache License 2.0
2.25k stars 574 forks source link

How to properly crop the freesurfer images #515

Open dyhan316 opened 1 year ago

dyhan316 commented 1 year ago

Task (what are you trying to do/register?)

I am trying to prepare my data to be used as according to the README.md, but I am stuck in knowing how/where to crop exactly

What have you tried

By following the author's comments on this repository's issues, I have been able to register the freesurfer result to the Talarich space using talarich.xfm. Now all that's left is the crop the images, but I can't figure out how.

for reference, the README.md states,

In particular, we perform FreeSurfer recon-all steps up to skull stripping and affine normalization to Talairach space, and crop the images via ((48, 48), (31, 33), (3, 29)).

I am confused as to what (48,48) and so on mean in this context. The information does not seem to be enough to make a plane (only three 2d coordinates). Or am I missing something?

It would be great if anyone could write a code that they used or more description about how the cropping was exactly done!

Thank you in advance for your help 👍

dyhan316 commented 1 year ago

After digging around with the image, it seems that ((48, 48), (31, 33), (3, 29)) means to crop the 48,-48 pixels on the x-axis, and so on. For anyone who might come across this same issue, I have made a code below using nibabel that does this (I think it's right) (any corrections are welcome!)

import pandas as pd 
import numpy as np 
import shutil
import nibabel as nib
import os 
import subprocess
import glob
from pathlib import Path

#file_test
tgz_dir = "/storage/bigdata/ABCD/freesurfer/smri/0.ALL_FS_ABCD"
save_dir = "/scratch/connectome/dyhan316/test/save_here"
final_product_dir = "/scratch/connectome/dyhan316/test/save_here/final_product" 

os.makedirs(save_dir, exist_ok= True)
os.makedirs(final_product_dir, exist_ok= True)

os.chdir(save_dir)
import time 
time_0 = time.time()

for tar_file in glob.glob(tgz_dir + '/*.tar') : 
    try : 
        os.chdir(save_dir)
        tar_file = Path(tar_file)
        path, file = tar_file.parent.absolute(), tar_file.parts[-1]

        #step1 : extract all the brain.mgz and talarich.xfm
        subprocess.run(f'tar -xvf {tar_file} {file[:-4]}/mri/brain.mgz {file[:-4]}/mri/transforms/talairach.xfm',
                      capture_output= True, shell = True)

        #step2 : perform talarich transform
        os.chdir(save_dir + f"/{file[:-4]}/mri")
        subprocess.run(f"mri_convert brain.mgz --apply_transform transforms/talairach.xfm -oc 0 0 0 {final_product_dir}/sub-{file[:-4].split(sep = '_')[0]}.mgz",
                      capture_output=True, shell = True)

        #step3 : perform crop and saving it 
        os.chdir(final_product_dir)
        img = nib.load(f"./sub-{file[:-4].split(sep = '_')[0]}.mgz")
        img_slice = img.slicer[48:-48, 31:-33, 3:-29]
        #import pdb ; pdb.set_trace()
        nib.save(img_slice, f"./sub-{file[:-4].split(sep = '_')[0]}.nii")

    except : 
        print(f"{tar_file} has failed")
    #os.chdir(save_dir)
    #subprocess.run(f"rm -r file[:-4]", capture_output= True, shell = True)

print(time.time() - time_0)