Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

After smoothing, the image is divided by the smoothed mask to correct… #321

Merged
merged 2 commits into from
Sep 28, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion rabies/analysis_pkg/diagnosis_pkg/analysis_QC.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,10 +57,11 @@ def get_maps(prior, prior_list, corr_variable, mask_file, smoothing=False, non_p
maps.append(corr_map)

if smoothing:
mask_img = sitk.ReadImage(mask_file)
import nibabel as nb
affine = nb.load(mask_file).affine[:3,:3]
for i in range(len(maps)):
maps[i] = sitk.GetArrayFromImage(smooth_image(recover_3D(mask_file,maps[i]), affine, 0.3))[volume_indices]
maps[i] = sitk.GetArrayFromImage(smooth_image(recover_3D(mask_file,maps[i]), affine, 0.3, mask_img))[volume_indices]

return maps

Expand Down
3 changes: 2 additions & 1 deletion rabies/confound_correction_pkg/confound_correction.py
Original file line number Diff line number Diff line change
Expand Up @@ -490,7 +490,8 @@ def _run_interface(self, runtime):
'''
import nibabel as nb
affine = nb.load(bold_file).affine[:3,:3] # still not sure how to match nibabel's affine reliably
timeseries_img = smooth_image(timeseries_img, affine, cr_opts.smoothing_filter)
mask_img = sitk.ReadImage(brain_mask_file, sitk.sitkFloat32)
timeseries_img = smooth_image(timeseries_img, affine, cr_opts.smoothing_filter, mask_img)

cleaned_path = cr_out+'/'+filename_split[0]+'_cleaned.nii.gz'
sitk.WriteImage(timeseries_img, cleaned_path)
Expand Down
18 changes: 17 additions & 1 deletion rabies/confound_correction_pkg/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -447,7 +447,7 @@ def phase_randomized_regressors(confounds_array, frame_mask, TR):
return randomized_confounds_array


def smooth_image(img, affine, fwhm):
def smooth_image(img, affine, fwhm, mask_img):
# apply nilearn's Gaussian smoothing on a SITK image
from nilearn.image.image import _smooth_array
from rabies.utils import copyInfo_4DImage, copyInfo_3DImage
Expand All @@ -466,12 +466,28 @@ def smooth_image(img, affine, fwhm):
array_3d = sitk.GetArrayFromImage(img)
arr = array_3d.transpose(2,1,0) # re-orient the array to match the affine
smoothed_arr = _smooth_array(arr, affine, fwhm=fwhm, ensure_finite=True, copy=True)

# smoothing creates leakage around mask boundaries
# correct for edge effects by dividing by the smoothed mask like FSL https://johnmuschelli.com/fslr/reference/fslsmooth.html
mask_arr = sitk.GetArrayFromImage(mask_img).transpose(2,1,0) # re-orient the array to match the affine
smoothed_mask = _smooth_array(mask_arr, affine, fwhm=fwhm, ensure_finite=True, copy=True)
smoothed_mask = smoothed_mask.transpose(2,1,0) # recover SITK orientation

# recover SITK img
if dim==4:
smoothed_arr = smoothed_arr.transpose(3,2,1,0)
smoothed_arr /= smoothed_mask # correct for edge effect
smoothed_arr *= sitk.GetArrayFromImage(mask_img) # re-apply mask
smoothed_arr[np.isnan(smoothed_arr)] = 0

smoothed_img = copyInfo_4DImage(sitk.GetImageFromArray(
smoothed_arr, isVector=False), img, img)
elif dim==3:
smoothed_arr = smoothed_arr.transpose(2,1,0)
smoothed_arr /= smoothed_mask # correct for edge effect
smoothed_arr *= sitk.GetArrayFromImage(mask_img) # re-apply mask
smoothed_arr[np.isnan(smoothed_arr)] = 0

smoothed_img = copyInfo_3DImage(sitk.GetImageFromArray(
smoothed_arr, isVector=False), img)

Expand Down