Crop non-overlapping regions with a 3D mask


I’m working on python 3.9.
I would like the crop the skull of fMRI scans based on a mask created based on anatomical data. The dimensions of the mask is 181 x 217 X 181. My fMRI scans are 4D time-series data with dimensions 61 X 73 X 61 X 180.
First, I resampled the mask using:

resampled_mask = nilearn.image.resample_to_img(fmri, mask, interpolation = ‘nearest’).

Is this a good method to compute a mask or should I recompute a mask directly on epi images?

Then, I’m looking for cropping regions of my fMRI images that are outside of the mask to recreate 3D sequence of scans but I did not find any functions on nilearn to do that (nilearn.masking.apply_mask outputs a 2D matrix of features as a function of time but I would like an output of dimension x, y, z, time)

Thanks per advance.

I found a solution using nilearn to compute mask epi for all the time series, apply them to all volumes iteratively and then concatenate all volumes into a time series. The last part is to crop 0 values and to resample all time series to the same dimensions.

import numpy as np
import nibabel as nib
import nilearn
from nilearn.image import resample_to_img as res
from nilearn.masking import compute_epi_mask
from nilearn.image import iter_img
from nilearn.image import math_img
from nilearn.image import concat_imgs
from nilearn.image import crop_img

for f in file_list:
    fmri = nib.load(f)
    mask = nilearn.masking.compute_epi_mask(fmri)
    filename = f[:-17] + "masked.nii"
    masked_list = []
    for img in iter_img(fmri):
        masked = math_img('img*mask', img = img, mask=mask)
    final = concat_imgs(masked_list)
    crop = crop_img(final)
    crop_resampled = res(crop, example, interpolation = 'nearest')