diff --git a/HypVINN/README.md b/HypVINN/README.md index 2317c7047..0937e3e9e 100644 --- a/HypVINN/README.md +++ b/HypVINN/README.md @@ -7,7 +7,7 @@ Hypothalamic subfields segmentation pipeline ### Requirements * Same as FastSurfer. -* If the T1w and T2w images are available and not co-registered, FreeSurfer should be sourced to run the registration code, and the mri_coreg and mri_vol2vol binaries should also be available. +* If the T1w and T2w images are available and not co-registered, HypVINN can register the T2w image to the T1w reference internally. ### Model weights * EUDAT (FZ Jülich) data repository: https://b2share.fz-juelich.de/records/2af6da63d5c1414b832c1f606bbd068a @@ -30,7 +30,7 @@ Note: These weights (version 1.1) are retrained compared to paper ([version 1.0] * `--t2 ` : T2 image path (optional) * `--seg_log` : Path to file in which run logs will be saved. If not set logs will be stored in `/scripts/hypvinn_seg.log` ### Image processing options - * `--reg_mode` : Ignored, if no T2 image is passed. Specifies the registration method used to register T1 and T2 images. Options are 'coreg' (default) for mri_coreg, 'robust' for mri_robust_register, and 'none' to skip registration (this requires T1 and T2 are externally co-registered). + * `--reg_mode` : Ignored, if no T2 image is passed. Specifies the registration method used to register T1 and T2 images. Options are 'coreg' (default) for `neuroreg.coreg` and 'none' to skip registration (this requires T1 and T2 are externally co-registered). * `--qc_snap`: Activate the creation of QC snapshots of the predicted HypVINN segmentation in `/qc_snapshots`. The created QC snapshots are created to simplify the visual quality control process. ### FastSurfer Technical parameters (see FastSurfer documentation) * `--device` @@ -116,7 +116,7 @@ Note: These weights (version 1.1) are retrained compared to paper ([version 1.0] |--hypothalamus.HypVINN.nii.gz(Hypothalamus Segmentation) |-- hypothalamus_mask.HypVINN.nii.gz (Hypothalamus Segmentation Mask) |-- transforms - |-- t2tot1.lta (FreeSurfer registration file, only available if registration is performed) + |-- t2tot1.lta (T2-to-T1 registration transform, only available if registration is performed) |-- qc_snapshots : QC outputs (optional) |-- hypothalamus.HypVINN_qc_screenshoot.png (Coronal quality control image) |-- stats : Statistics outputs diff --git a/HypVINN/run_prediction.py b/HypVINN/run_prediction.py index 1e6fe4014..6086c8cf0 100644 --- a/HypVINN/run_prediction.py +++ b/HypVINN/run_prediction.py @@ -109,11 +109,11 @@ def option_parse() -> argparse.ArgumentParser: "--reg_mode", type=str, default="coreg", - choices=["none", "coreg", "robust"], - help="Freesurfer Registration type to run. coreg: mri_coreg, " - "robust : mri_robust_register, none: entirely deactivates " - "registration of T2 to T1, if both images are passed, " - "images need to be register properly externally.", + choices=["none", "coreg"], + help="Registration type to run for T2-to-T1 alignment when both images " + "are passed. coreg: neuroreg.coreg, none: entirely deactivates " + "registration of T2 to T1, so images need to be registered " + "properly externally.", ) parser.add_argument( @@ -170,7 +170,7 @@ def main( hypo_statsfile: str = HYPVINN_STATS_NAME, qc_snapshots: bool = False, threads: int | None = None, - reg_mode: Literal["coreg", "robust", "none"] = "coreg", + reg_mode: Literal["coreg", "none"] = "coreg", batch_size: int = 1, async_io: bool = False, device: str = "auto", @@ -211,7 +211,7 @@ def main( Whether to create QC snapshots. threads : int, optional If not None, updates the FastSurfer global setting in `FastSurfer.utils.parallel`. - reg_mode : "coreg", "robust", "none", default="coreg" + reg_mode : "coreg", "none", default="coreg" The registration mode to use. batch_size : int, default=1 The batch size to use. @@ -229,6 +229,10 @@ def main( """ if threads is not None and threads > 1: set_num_threads(threads) + # here on the main thread, not in Inference, because the T2 registration below is submitted to + # a worker before Inference exists and neuroreg does no thread handling of its own. Without + # this it would run on torch's own default, ignoring --threads. + torch.set_num_threads(get_num_threads()) from concurrent.futures import Future diff --git a/HypVINN/utils/__init__.py b/HypVINN/utils/__init__.py index 64cd0b1c4..0753985d1 100644 --- a/HypVINN/utils/__init__.py +++ b/HypVINN/utils/__init__.py @@ -14,4 +14,4 @@ class ViewOperationDefinition(TypedDict): ViewOperations = dict[Plane, ViewOperationDefinition | None] ModalityMode = Literal["t1", "t2", "t1t2"] ModalityDict = dict[Literal["t1", "t2"], ndarray] -RegistrationMode = Literal["robust", "coreg", "none"] +RegistrationMode = Literal["coreg", "none"] diff --git a/HypVINN/utils/preproc.py b/HypVINN/utils/preproc.py index f983ee636..31ab09e30 100644 --- a/HypVINN/utils/preproc.py +++ b/HypVINN/utils/preproc.py @@ -36,7 +36,7 @@ def t1_to_t2_registration( threads: int = -1, ) -> Path: """ - Register the T2 image to the T1 reference using either mri_coreg or mri_robust_register. + Register the T2 image to the T1 reference using neuroreg. Parameters ---------- @@ -49,7 +49,7 @@ def t1_to_t2_registration( lta_path : Path The path to the lta transform. registration_type : RegistrationMode, default="coreg" - The type of registration to be used. It can be either "coreg" or "robust". + The type of registration to be used. It can be either "coreg" or "none". threads : int, default=-1 The number of threads to be used. If it is less than or equal to 0, the number of threads will be automatically determined. @@ -61,50 +61,25 @@ def t1_to_t2_registration( Raises ------ - RuntimeError - If mri_coreg, mri_vol2vol, or mri_robust_register fails to run or if they cannot - be found. + ValueError + If `registration_type` names a mode this function does not implement. """ - import shutil - - from FastSurferCNN.utils.parallel import get_num_threads - from FastSurferCNN.utils.run_tools import Popen - - if threads <= 0: - threads = get_num_threads() - - def from_freesurfer_home(fs_binary: str) -> str: - if not os.environ.get("FREESURFER_HOME", ""): - raise RuntimeError( - f"Could not find {fs_binary}, source FreeSurfer or set the FREESURFER_HOME environment variable" - ) - return os.environ["FREESURFER_HOME"] + "/bin/" + fs_binary - - def run_fs_binary(fs_binary: str, args: list[str]) -> int: - fs_binary = shutil.which(fs_binary) or from_freesurfer_home(fs_binary) - args = [fs_binary] + list(map(str, args)) - LOGGER.info("Running " + " ".join(args)) - retval = Popen(args).finish() - if retval.retcode != 0: - LOGGER.error(f"{fs_binary} failed with error code {retval.retcode}.") - raise RuntimeError(f"{fs_binary} failed") - - LOGGER.info(f"{fs_binary} finished in {retval.runtime}!") + from neuroreg import coreg if registration_type == "coreg": - run_fs_binary( - "mri_coreg", - ["--mov", t2_path, "--targ", t1_path, "--reg", lta_path, "--threads", str(threads)], - ) - run_fs_binary( - "mri_vol2vol", - ["--mov", t2_path, "--targ", t1_path, "--reg", lta_path, "--o", output_path, "--cubic", "--keep-precision"], + LOGGER.info("Running neuroreg.coreg for T2-to-T1 registration.") + # on the CPU, coreg's default. Whether a GPU would be faster here has not been measured. + # The thread count is the process-wide torch one, set in run_prediction before this is + # submitted, since neuroreg does no thread handling of its own. + coreg( + str(t2_path), + str(t1_path), + lta_name=str(lta_path), + mapped_name=str(output_path), + keep_dtype=True, ) else: - run_fs_binary( - "mri_robust_register", - ["--mov", t2_path, "--dst", t1_path, "--lta", lta_path, "--mapmov", output_path, "--cost", "NMI"], - ) + raise ValueError(f"Unknown registration type: {registration_type}") return output_path diff --git a/doc/overview/LONG.md b/doc/overview/LONG.md index 374c404cd..d18495404 100644 --- a/doc/overview/LONG.md +++ b/doc/overview/LONG.md @@ -58,7 +58,7 @@ $FASTSURFER_HOME/long_fastsurfer.sh \ Note, with a few exceptions, you can add additional flags that can be understood by `run_fastsurfer.sh`, which will be passed through, e.g. the `--3T` when working with 3T images. The exceptions: `--sid`, `--t1` and `--t2` are set per time point from `--tpids`, `--t1s` and `--t2s`; `--seg_only` and `--surf_only` are replaced by stages (see [Running Individual Stages](#running-individual-stages) below); and `--remove_suffix`, `--keepgeom` and `--native_image` are not supported in longitudinal processing. -If you have T2-weighted images, pass them with `--t2s ...`, one per time point in the same order as `--tpids`. They are used by the hypothalamus module (HypVINN) in the segmentation of each time point, where each T2 is registered directly to that time point's T1 in template space. The template itself does not need them, as the hypothalamus module does not run on it. Either every time point gets a T2 or none does: a T2 changes what the hypothalamus module computes, so a series where only some time points had one would not be comparable over time. `--reg_mode` selects how each T2 reaches template space: with `coreg` (the default) or `robust`, it is registered to its time point's T1 there; with `none`, it must already be co-registered with the T1 you passed in `--t1s`, and it is then mapped into template space with the same transform as that T1. +If you have T2-weighted images, pass them with `--t2s ...`, one per time point in the same order as `--tpids`. They are used by the hypothalamus module (HypVINN) in the segmentation of each time point, where each T2 is registered directly to that time point's T1 in template space. The template itself does not need them, as the hypothalamus module does not run on it. Either every time point gets a T2 or none does: a T2 changes what the hypothalamus module computes, so a series where only some time points had one would not be comparable over time. `--reg_mode` selects how each T2 reaches template space: with `coreg` (the default), it is registered to its time point's T1 there; with `none`, it must already be co-registered with the T1 you passed in `--t1s`, and it is then mapped into template space with the same transform as that T1. The above command will, of course, be slightly different when using your preferred installation in Singularity or Docker. For example, for Singularity: @@ -161,7 +161,7 @@ Behind the Scenes --t1s ... \ --tpids ``` - This segments each time point once, only to obtain a brain mask, and then registers (aligns) all time point images into the unbiased mid-space using `mri_robust_template`. It will also create the template image, the voxel-wise median across time. For single time point cases, it will align the input into a standard upright position. Finally, it maps each time point's input into template space, where the time point stages below pick it up. + This segments each time point once, only to obtain a brain mask, and then registers (aligns) all time point images into the unbiased mid-space using neuroreg's robust multi-image registration. It will also create the template image, the voxel-wise median across time. For single time point cases, it will align the input into a standard upright position, obtained from a rigid registration of its segmentation to MNI centroids. Finally, it maps each time point's input into template space, where the time point stages below pick it up. 2. **Template Seg** (stage `template_seg`): Next, the template image will be segmented via a call to `run_fastsurfer.sh --sid --base --seg_only ...` where the `--base` flag indicates that the input image will be taken from the already existing template directory. The cerebellum and hypothalamus modules do not run on the template. 3. **Template Surf** (stage `template_surf`): This is followed by the surface processing of the template `run_fastsurfer.sh --sid --base --surf_only ...`. 4. **Long Seg** (stage `long_seg`): Next, the segmentation of each time point is performed `run_fastsurfer.sh --sid --long --seg_only ...`. It only needs step 1, so it can run at the same time as steps 2 and 3. diff --git a/pyproject.toml b/pyproject.toml index 20067df77..4d1f8a590 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -38,7 +38,7 @@ dependencies = [ 'matplotlib>=3.7.1', 'meshpy>=2025.1.1', # needed for FastSurfer-CC 'monai>=1.4.0', # needed for FastSurfer-CC - 'neuroreg>=0.6.2', # needed for registration / etiv + 'neuroreg>=0.12.0', # needed for registration / etiv 'nibabel>=5.4.0', # needed to fix a bug in nibabel 'numpy>=1.25', 'packaging', diff --git a/recon_surf/long_prepare_template.sh b/recon_surf/long_prepare_template.sh index 72af9a14c..98f518327 100755 --- a/recon_surf/long_prepare_template.sh +++ b/recon_surf/long_prepare_template.sh @@ -24,7 +24,7 @@ # to setup a base template for processing with a slightly modified version of run_fastsurfer. # # 1. It obtains a brainmask on all time points for a specific subject -# 2. It uses mri_robust_template to co-register all time points into a mid-space +# 2. It uses neuroreg's multireg to co-register all time points into a mid-space # and sets up the base (subject template) directory with the aligned images the median image # and the forward and backward rigid transformations # @@ -38,8 +38,7 @@ # - Maybe use intersection of tp masks as brainmask for base, as done in FreeSurfer. # - Add flag for adding a new time point to an existing base/template. # -# FreeSurfer requirements: -# mri_convert, mri_robust_template, mri_mask, mri_concatenate_lta, make_upright +# neuroreg requirements: multireg, segreg, vol2vol, lta, mri (mask, diff, info) # ################################################################################################### @@ -74,6 +73,7 @@ robust_template_avg_arg=1 # median for template creation (if more than 1 time p # default arguments python="python3 -s" # avoid user-directory package inclusion sd="$SUBJECTS_DIR" +device="" # torch device for neuroreg (empty/auto = let neuroreg pick GPU else CPU) # init variables that need to be passed tid="" @@ -175,7 +175,9 @@ case $key in --t1s) while [[ $# -gt 0 ]] && [[ $1 != -* ]] ; do t1s+=("$1") ; shift ; done ;; --sd) export SUBJECTS_DIR="$1" ; shift ;; # these flags are passed through to run_prediction.py - --vox_size|--device|--viewagg_device|--conform_to_1mm_threshold) run_pred_flags+=("$key" "$1") ; shift ;; + --vox_size|--viewagg_device|--conform_to_1mm_threshold) run_pred_flags+=("$key" "$1") ; shift ;; + # device is passed through to run_prediction.py and also forwarded to neuroreg + --device) run_pred_flags+=("$key" "$1") ; device="$1" ; shift ;; --threads|--threads_seg) threads="$1" ; shift ;; --batch_size|--batch) if [[ "$key" == "--batch" ]] @@ -243,10 +245,6 @@ fi # check that SUBJECTS_DIR exists check_create_subjects_dir_properties "$SUBJECTS_DIR" -# the FreeSurfer tools below (mri_convert, mri_mask, mri_diff, mri_robust_template, ...) refuse to -# run without a license, so fail here rather than after the segmentation of every time point -auto_detect_fs_license "the longitudinal template preparation" || exit $? - ################################## SETUP and LOGFILE ############################## @@ -274,7 +272,7 @@ if [[ "$code" != 0 ]] ; then echo "ERROR: Getting the version failed (code=$code echo "Version: $VERSION" | tee -a "$LF" echo "Log file for long_prepare_template" >> "$LF" -# the segmentation per time point and mri_robust_template below share this budget; exported +# the segmentation per time point and neuroreg's multireg below share this budget; exported # before the host block, so that it logs the limits in force if ! resolve_threads "$threads" auto ; then echo "$threads_note" | tee -a "$LF" ; exit 1 ; fi set_thread_env "$threads_budget" @@ -315,7 +313,7 @@ do # check if geometry differs across time if [[ "$s" != "${t1s[0]}" ]] then - cmda=(mri_diff --notallow-pix --notallow-geo "$s" "${t1s[0]}" --res-thresh "0.000001") + cmda=($python -m neuroreg.cli.mri diff --notallow-pix --notallow-geo "$s" "${t1s[0]}" --res-thresh "0.000001") difftext=$("${cmda[@]}") retcode=${PIPESTATUS[0]} if [[ "$retcode" != 0 ]] ; then geodiff+="Comparing $s and ${t1s[0]} (code $retcode):\n$difftext\n" ; fi @@ -328,7 +326,7 @@ then echo "*******************************************************************************" echo "WARNING: Image parameters differ across time, maybe due to acquisition changes?" echo " Consistent changes in, e.g., resolution can potentially bias a " - echo " longitudinal study! You can check image parameters by running mri_info" + echo " longitudinal study! You can check image parameters by running neuroreg mri info" echo " on each input image." echo "*******************************************************************************" echo "$geodiff" @@ -391,7 +389,7 @@ for ((i=0;i<${#tpids[@]};++i)); do run_it "$LF" "${cmda[@]}" # mask is binary, we need to use on conformed image: - cmda=(mri_mask "$conformed_name" "$mask_name" "$mdir/cross_brainmask${extension}") + cmda=($python -m neuroreg.cli.mri mask "$conformed_name" "$mask_name" "$mdir/cross_brainmask${extension}") run_it "$LF" "${cmda[@]}" done @@ -430,40 +428,142 @@ do done +reference_centroids="mni_icbm152_t1_tal_nlin_asym_09c" + +# The base grid, built from the conformed standard rather than copied from a time point: a 256 mm +# field of view, axis aligned, centred on the world origin. The voxel size is the largest across +# time points of each one's smallest side, which is the rule multireg's create_template_geometry +# applies when it derives a grid itself, and mri_robust_template before it. Passing the grid is +# only needed to keep it axis aligned, since the pre-posed volumes below are oblique, so it should +# not quietly decide the resolution differently. It also does not depend on the order the time +# points were given. Stating the orientation rather than inheriting it drops the float dust a +# conformed header can carry. +base_vox=$(for v in "${normInVols[@]}" ; do + $python -m neuroreg.cli.mri info --res "$v" + done | awk '{s=$1; for(j=2;j<=NF;j++) if($jg) g=s} END{printf "%.10g\n", g}') +echo "base grid: 256 mm field of view at ${base_vox}mm" | tee -a "$LF" +base_geom="${SUBJECTS_DIR}/$tid/mri/base_geom${extension}" +cmd="$python -m neuroreg.cli.mri geom --fov 256 --vox-size $base_vox" +cmd="$cmd --orientation LIA --cras 0,0,0 --o $base_geom" +RunIt "$cmd" "$LF" + if [ ${#tpids[@]} == 1 ] then - # if only a single time point, create fake 'base' by making the image upright - # this assures that also subjects with a single time point get processed as the other - # subjects in the longitudinal stream + # If only a single time point, we still create a 'base' so that single-tp subjects are + # processed like the others. Instead of FreeSurfer's make_upright (mirror-symmetry + # uprighting), we rigidly register this time point's segmentation to the reference centroid + # template - the same way the talairach initialization prealigns - which yields a + # reproducible, template-aligned (upright) base pose without scaling the anatomy. The + # segmentation always exists here (it was created above to obtain the brainmask). The base + # lives in the time point's own (conformed) geometry. + + seg0="$SUBJECTS_DIR/$tid/long-inputs/${tpids[0]}/cross_aparc+aseg.orig${extension}" + + # 1. rigid segmentation-to-template registration: tp -> base pose (RAS-to-RAS transform) + cmd="$python -m neuroreg.cli.segreg --seg $seg0 --centroids $reference_centroids" + cmd="$cmd --dof 6 --lta ${ltaXforms[0]}" + RunIt "$cmd" "$LF" - # 1. make the norm upright (base space) - cmd="make_upright ${normInVols[0]} \ - ${SUBJECTS_DIR}/$tid/mri/base_brainmask${extension} ${ltaXforms[0]}" + # 2. create the base brainmask by mapping the norm into the base pose. The base grid is the one + # built above, centred on the world origin, which is where the uprighting puts the brain; + # keeping the time point's own centre instead would rotate the head towards the edge of its + # field of view and could clip it. --keep-dtype because the inputs are conformed uchar and the + # default is float32, which the base run would not accept as conformed. + cmd="$python -m neuroreg.cli.vol2vol --in ${normInVols[0]} --transform ${ltaXforms[0]}" + cmd="$cmd --ref $base_geom --interp cubic --keep-dtype" + cmd="$cmd --out ${SUBJECTS_DIR}/$tid/mri/base_brainmask${extension}" RunIt "$cmd" "$LF" - # 2. create the upright orig volume - cmd="mri_convert -rt cubic \ - -at ${ltaXforms[0]} ${subjInVols[0]} ${SUBJECTS_DIR}/$tid/mri/orig.mgz" + # 3. create the base orig volume the same way, so that it lands on the same grid + cmd="$python -m neuroreg.cli.vol2vol --in ${subjInVols[0]} --transform ${ltaXforms[0]}" + cmd="$cmd --ref $base_geom --interp cubic --keep-dtype" + cmd="$cmd --out ${SUBJECTS_DIR}/$tid/mri/orig.mgz" RunIt "$cmd" "$LF" -else #more than 1 time point: + # 4. record the base geometry in the transform. segreg registers to a centroid set rather than + # to an image, so it leaves the destination geometry empty, while multireg fills it in for + # more than one time point. Consumers such as mri_fuse_segmentations read it, so both paths + # should describe the base the same way. + # Written to a side file and moved into place, so a failure here cannot leave a transform + # behind that is half rewritten. + tmplta="${SUBJECTS_DIR}/$tid/mri/transforms/${tpids[0]}_withgeom.lta" + cmd="$python -m neuroreg.cli.lta convert ${ltaXforms[0]} $tmplta" + cmd="$cmd --dst-img ${SUBJECTS_DIR}/$tid/mri/base_brainmask${extension}" + RunIt "$cmd" "$LF" + cmd="mv -f $tmplta ${ltaXforms[0]}" + RunIt "$cmd" "$LF" +else #more than 1 time point: - # create the 'mean/median' norm volume: - cmd="mri_robust_template --mov ${normInVols[*]}" - cmd="$cmd --lta ${ltaXforms[*]}" + # forward the device to neuroreg if explicitly requested (empty/auto = neuroreg auto-selects) + device_opt="" + if [[ -n "$device" ]] && [[ "$device" != "auto" ]] ; then device_opt="--device $device" ; fi + + # Put every time point into the reference centroid pose before co-registering them. Otherwise + # multireg opens by registering everything to one pseudo-randomly chosen time point, and the + # base inherits that time point's pose, so the same subject lands differently depending on + # which one was drawn. Pre-posing them removes that and makes the base upright, the same pose + # the single time point branch above produces. + # + # The pose is applied to the header only, so no interpolation happens here and multireg still + # does all of its own registering, including the opening round. The transforms are composed + # back onto the original time point geometry afterwards. + hdrInVols=() + preXforms=() + mrgXforms=() + for ((i=0;i<${#tpids[@]};++i)) + do + mdir="$SUBJECTS_DIR/$tid/long-inputs/${tpids[i]}" + pre="${SUBJECTS_DIR}/$tid/mri/transforms/${tpids[i]}_prereg.lta" + preXforms+=("$pre") + mrgXforms+=("${SUBJECTS_DIR}/$tid/mri/transforms/${tpids[i]}_prereg_to_${tid}.lta") + hdr="$mdir/cross_brainmask_prereg${extension}" + hdrInVols+=("$hdr") + + cmd="$python -m neuroreg.cli.segreg --seg $mdir/cross_aparc+aseg.orig${extension}" + cmd="$cmd --centroids $reference_centroids --dof 6 --lta $pre" + RunIt "$cmd" "$LF" + + cmd="$python -m neuroreg.cli.vol2vol --in ${normInVols[i]} --transform $pre" + cmd="$cmd --header-only --out $hdr" + RunIt "$cmd" "$LF" + done + + # robust co-registration of all time points into an unbiased mid-space, creating the + # 'mean/median' norm (brainmask) volume and the forward transforms (tp -> base). The base grid + # has to be passed explicitly because the pre-posed volumes are oblique, and the one multireg + # would derive by averaging their direction cosines would be oblique too, which the base run + # would reject as unconformed: + cmd="$python -m neuroreg.cli.multireg --mov ${hdrInVols[*]}" + cmd="$cmd --lta ${mrgXforms[*]}" cmd="$cmd --template ${SUBJECTS_DIR}/$tid/mri/base_brainmask${extension}" + cmd="$cmd --template-geom $base_geom" cmd="$cmd --average ${robust_template_avg_arg}" - cmd="$cmd --sat 4.685" + cmd="$cmd --sat 4.685 --keep-dtype $device_opt" RunIt "$cmd" "$LF" - # create the 'mean/median' input (orig) volume: - cmd="mri_robust_template --mov ${subjInVols[*]}" + # compose the pre-pose back in, so the transforms map the original time point and not the + # pre-posed copy. 'lta concat A B' applies A first, and carries the geometry of A's source and + # B's destination. The pre-posed volumes have served their purpose once this is done. + for ((i=0;i<${#tpids[@]};++i)) + do + cmd="$python -m neuroreg.cli.lta concat ${preXforms[i]} ${mrgXforms[i]} ${ltaXforms[i]}" + RunIt "$cmd" "$LF" + cmd="rm -f ${mrgXforms[i]} ${hdrInVols[i]}" + RunIt "$cmd" "$LF" + done + + # create the 'mean/median' input (orig) volume by reusing the transforms above (no + # registration, --noit); the per-tp LTAs carry the base geometry so orig.mgz lands on + # the same grid as base_brainmask. --keep-dtype because the time points are conformed + # uchar and the base run checks orig.mgz against exactly that: + cmd="$python -m neuroreg.cli.multireg --mov ${subjInVols[*]}" cmd="$cmd --average ${robust_template_avg_arg}" cmd="$cmd --ixforms ${ltaXforms[*]}" - cmd="$cmd --noit" + cmd="$cmd --noit --keep-dtype" t1=${SUBJECTS_DIR}/$tid/mri/orig.mgz - cmd="$cmd --template $t1" + cmd="$cmd --template $t1 $device_opt" RunIt "$cmd" "$LF" fi # more than one time point @@ -482,8 +582,13 @@ done for ((i=0;i<${#tpids[@]};++i)) do mdir="$SUBJECTS_DIR/$tid/long-inputs/${tpids[i]}" - # map orig to base space - cmd="mri_convert -at ${ltaXforms[$i]} -rt $interpol $mdir/mri/rawavg.mgz $mdir/long_conform${extension}" + # map each time point input into base space. The raw input, not the cross_conform the base + # template was built from: the --long run conforms this itself, so it should see the original + # intensities rather than a second rescale of an already uchar image. The transforms are + # RAS-to-RAS, so the base geometry comes from --ref and not from the transform. + cmd="$python -m neuroreg.cli.vol2vol --in $mdir/mri/rawavg.mgz" + cmd="$cmd --transform ${ltaXforms[$i]} --ref ${SUBJECTS_DIR}/$tid/mri/base_brainmask${extension}" + cmd="$cmd --interp $interpol --out $mdir/long_conform${extension}" RunIt "$cmd" "$LF" done diff --git a/recon_surf/recon-surf.sh b/recon_surf/recon-surf.sh index 474baf493..339e2ddc5 100755 --- a/recon_surf/recon-surf.sh +++ b/recon_surf/recon-surf.sh @@ -533,7 +533,7 @@ fi # create orig.mgz and aparc.DKTatlas+aseg.orig.mgz (copy of T1 and segmentation) # also ensures .mgz format (in case inputs are nifti) -cmd="mri_convert $t1 $mdir/orig.mgz" +cmd="$python -m neuroreg.cli.vol2vol --in $t1 --out $mdir/orig.mgz" RunIt "$cmd" "$LF" asegdkt_segfile_manedit=$(add_file_suffix "$asegdkt_segfile" "manedit") @@ -542,7 +542,7 @@ if [[ ! "$asegdkt_segfile_manedit" =~ (\.manedit){2,}\. ]] && [[ -f "$asegdkt_se then asegdkt_segfile="$asegdkt_segfile_manedit" # use the manedit file fi -cmd="mri_convert $asegdkt_segfile $mdir/aparc.DKTatlas+aseg.orig.mgz" +cmd="$python -m neuroreg.cli.vol2vol --in $asegdkt_segfile --out $mdir/aparc.DKTatlas+aseg.orig.mgz" RunIt "$cmd" "$LF" # pctsurfcon reads mri/rawavg.mgz and hardcodes that path. run_fastsurfer.sh writes it from the @@ -672,7 +672,7 @@ fi # the difference between nu and orig_nu is the fact that nu has the talairach-registration header # create norm by masking nu (supports manedit-ed mask) -cmda=(mri_mask "$mdir/nu.mgz" "$mask" "$mdir/norm.mgz") +cmda=($python -m neuroreg.cli.mri mask "$mdir/nu.mgz" "$mask" "$mdir/norm.mgz") run_it "$LF" "${cmda[@]}" if [[ "$get_t1" == "true" ]] then @@ -693,7 +693,7 @@ then cmda=(mri_normalize -g 1 -seed 1234 -mprage "$mdir/nu.mgz" "$mdir/T1.mgz" "${noconform_if_hires[@]}") run_it "$LF" "${cmda[@]}" # create brainmask by masking T1 (supports manedit-ed mask) - cmda=(mri_mask "$mdir/T1.mgz" "$mask" "$mdir/brainmask.mgz") + cmda=($python -m neuroreg.cli.mri mask "$mdir/T1.mgz" "$mask" "$mdir/brainmask.mgz") run_it "$LF" "${cmda[@]}" else # create brainmask by linkage to norm.mgz (masked nu.mgz) diff --git a/recon_surf/talairach-reg.sh b/recon_surf/talairach-reg.sh index c2c81e4f2..3dd2efe8f 100755 --- a/recon_surf/talairach-reg.sh +++ b/recon_surf/talairach-reg.sh @@ -214,7 +214,10 @@ else prealigned_name=$mdir/segreg_prealigned.mgz target_geom="$FREESURFER_HOME/average/mni305.cor.mgz" # only used for its conformed 1mm header; does not imply registration to this template - cmd=(mri_convert --apply_transform "$prealigned_lta" --reslice_like "$target_geom" "$norm_name" "$prealigned_name") + # --keep-dtype: norm.mgz is uchar and talairach_avi below converts this to Analyze, so the + # float32 default would change the intensities mpr2mni305 sees, and with them talairach.xfm + cmd=($python -m neuroreg.cli.vol2vol --in "$norm_name" --transform "$prealigned_lta" + --ref "$target_geom" --interp linear --keep-dtype --out "$prealigned_name") run_it "$LF" "${cmd[@]}" # talairach.xfm: compute talairach full head (25sec) diff --git a/run_fastsurfer.sh b/run_fastsurfer.sh index b20b9be77..83b8403c4 100755 --- a/run_fastsurfer.sh +++ b/run_fastsurfer.sh @@ -285,16 +285,15 @@ SEGMENTATION PIPELINE: --t2 *Optional* T2 full head input (must be externally biasfield corrected when called with --no_biasfield). Requires an ABSOLUTE Path! - --reg_mode + --reg_mode Ignored, if no T2 image is passed. Specifies the registration method used to register T1 and T2 images. Options are 'coreg' (default) for - mri_coreg, 'robust' for mri_robust_register, and 'none' - to skip registration (this requires T1 and T2 are - externally co-registered). With --long, 'none' means - the T2 is co-registered with the T1 this time point - was built from, and it is mapped into template space - with the same transform as that T1. + neuroreg.coreg and 'none' to skip registration (this + requires T1 and T2 are externally co-registered). With + --long, 'none' means the T2 is co-registered with the T1 + this time point was built from, and it is mapped into + template space with the same transform as that T1. --qc_snap Create QC snapshots in \$SUBJECTS_DIR/\$sid/qc_snapshots to simplify the QC process. @@ -587,8 +586,8 @@ case $key in --hypo_statsfile) hypo_statsfile="$1" ; shift ;; --reg_mode) mode=$(echo "$1" | tr "[:upper:]" "[:lower:]") - if [[ "$mode" =~ ^(none|coreg|robust)$ ]] ; then hypvinn_regmode="$mode" - else echo "Invalid --reg_mode option, must be 'none', 'coreg' or 'robust'." ; exit 1 + if [[ "$mode" =~ ^(none|coreg)$ ]] ; then hypvinn_regmode="$mode" + else echo "Invalid --reg_mode option, must be 'none' or 'coreg'." ; exit 1 fi shift # past value ;; @@ -911,9 +910,6 @@ if [[ "$run_seg_pipeline" == "true" ]] ; then if [[ "$run_talairach_registration" == "true" ]] ; then what_needs_license+=" and the talairach-registration in the segmentation pipeline" fi - if [[ -n "$t2" ]] && [[ "$hypvinn_regmode" != "none" ]] ; then - what_needs_license+=" and the T1-T2 registration in the segmentation pipeline" - fi fi if [[ -n "$what_needs_license" ]] then