diff --git a/.github/workflows/run_script_and_create_release.yml b/.github/workflows/run_script_and_create_release.yml
index 805c58e..1d5e66c 100644
--- a/.github/workflows/run_script_and_create_release.yml
+++ b/.github/workflows/run_script_and_create_release.yml
@@ -26,6 +26,12 @@ on:
         description: 'Release title (e.g. rYYYYMMDD)'
         required: true
 
+env:
+  # Hardcoding this variable just to make accessing the output files a bit simpler (for testing)
+  # See: https://github.com/spinalcordtoolbox/sct_tutorial_data/actions/runs/28119017640/job/83265838588?pr=35#step:3:43
+  # Also: https://github.com/orgs/community/discussions/25576#discussioncomment-3248364
+  SCT_DIR: "${{ github.workspace }}/spinalcordtoolbox"
+
 jobs:
   run-course-script:
     runs-on: ubuntu-latest
@@ -34,6 +40,7 @@ jobs:
     - name: Checkout spinalcordtoolbox
       uses: actions/checkout@v4
       with:
+        ref: jn/3508-merge-example-datasets-and-batch-scripts
         repository: spinalcordtoolbox/spinalcordtoolbox
         path: spinalcordtoolbox
 
@@ -56,6 +63,9 @@ jobs:
         path: ${{ github.event.repository.name }}
 
     - name: Run batch_single_subject.sh to generate intermediate files
+      # NB: SCT's copy of batch_single_subject.sh is standalone and cds into SCT_DIR
+      # So, even though we run the script here, the output files will be generated elsewhere
+      # This is just temporary for testing purposes, so don't commit these changes.
       run: |
         cd "${{ github.event.repository.name }}/single_subject"
         ./batch_single_subject.sh
@@ -70,11 +80,11 @@ jobs:
       uses: actions/upload-artifact@v4
       with:
         name: Single Subject CSV Files
-        path: ${{ github.event.repository.name }}/single_subject/data/**/*.csv
+        path: ${{ env.SCT_DIR }}/data/sct_course_data/single_subject/data/**/*.csv
 
     - name: "Package data into tutorial-specific .zip files"
       run: |
-        cd ${{ github.event.repository.name }}
+        cd "${{ env.SCT_DIR }}/data/sct_course_data/"
         awk -F, '{ print $1,$2 }' tutorial-datasets.csv | xargs -l zip -ur
 
     - uses: ncipollo/release-action@v1
@@ -92,7 +102,7 @@ jobs:
       env:
         GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
       with:
-        file: "${{ github.event.repository.name }}/*.zip"
+        file: "${{ env.SCT_DIR }}/data/sct_course_data/*.zip"
         draft: false
         tag_name: ${{ github.event.inputs.release_title }}
         release_id: ${{ steps.create_release.outputs.id }}
diff --git a/.github/workflows/scripts/process_csvs.py b/.github/workflows/scripts/process_csvs.py
index e0be780..bc29c11 100644
--- a/.github/workflows/scripts/process_csvs.py
+++ b/.github/workflows/scripts/process_csvs.py
@@ -124,8 +124,8 @@ def process_ap_ratio_norm_pam50(ci_root: Path, tutorials_root: Path) -> None:
 def process_csa_c3c4(ci_root: Path, tutorials_root: Path) -> None:
     out = tutorials_root / "shape-analysis/compute-csa-and-other-shape-metrics"
     process(
-        ci_root / "t2/csa_c3c4.csv",
-        out / "csa_c3c4.csv",
+        ci_root / "t2/csa_c2c3.csv",
+        out / "csa_c2c3.csv",
         keep_cols=["Filename", "Slice (I->S)", "VertLevel", "MEAN(area)", "STD(area)"],
         conditional_vertlevel=True,
     )
diff --git a/single_subject/batch_single_subject.sh b/single_subject/batch_single_subject.sh
index f752a2f..dab946e 100755
--- a/single_subject/batch_single_subject.sh
+++ b/single_subject/batch_single_subject.sh
@@ -18,7 +18,26 @@
 # ======================================================================================================================
 
 # If a command fails, set -e will make the whole script exit, instead of just resuming on the next line
-set -e
+set -ve
+
+# For full verbose, uncomment the next line
+# set -x
+
+# get starting time:
+start=$(date +%s)
+
+# Fetch OS type
+if uname -a | grep -i  darwin > /dev/null 2>&1; then
+  # OSX
+  open_command="open"
+elif uname -a | grep -i  linux > /dev/null 2>&1; then
+  # Linux
+  open_command="xdg-open"
+fi
+
+# download example data and enter our data directory
+sct_download_data -d sct_course_data
+cd "$SCT_DIR/data/sct_course_data/single_subject/data/"
 
 # Exit if user presses CTRL+C (Linux) or CMD+C (OSX)
 trap "echo Caught Keyboard Interrupt within script. Exiting now.; exit" INT
@@ -44,7 +63,7 @@ fi
 # ======================================================================================================================
 
 # Go to T2 folder
-cd data/t2
+cd t2
 # Spinal cord segmentation (using the new 2024 contrast-agnostic method)
 sct_deepseg spinalcord -i t2.nii.gz -qc ~/qc_singleSubj
 # The default output is t2_seg.nii.gz
@@ -74,9 +93,21 @@ fsleyes t2.nii.gz -cm greyscale t2_step1_canal.nii.gz -cm YlOrRd -a 70.0 t2_step
 # Check QC report: Go to your browser and do "refresh".
 
 # Optionally, you can use the generated disc labels to create a labeled segmentation
+# sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -discfile t2_totalspineseg_discs.nii.gz
+
+# If you wish to bypass the disc labeling model entirely and instead wish to use the legacy
+# vertebral labeling method (`sct_label_vertebrae`), you can do so.
 # Note: This approach is no longer recommended. Instead, use the disc labels directly in subsequent commands (e.g. `sct_process_segmentation`).
-sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -discfile t2_totalspineseg_discs.nii.gz
-# FIXME: Remove this command once the web tutorials are updated to no longer use labeled segmentations
+# sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -qc ~/qc_singleSubj
+
+# Optionally, if the legacy `sct_label_vertebrae` fails as well, you can initialize the command
+# by manually labeling the c2-c3 disc.
+# sct_label_utils -i t2.nii.gz -create-viewer 3 -o label_c2c3.nii.gz -msg "Click at the posterior tip of C2/C3 inter-vertebral disc"
+# sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -initlabel label_c2c3.nii.gz -qc ~/qc_singleSubj
+
+# If generating a labeled segmentation using `sct_label_vertebrae`, you can then extract vertebral
+# body labels from the segmentation using `-vert-body`.
+# sct_label_utils -i t2_seg_labeled.nii.gz -vert-body 3,9 -o t2_labels_vert.nii.gz
 
 
 
@@ -84,9 +115,9 @@ sct_label_vertebrae -i t2.nii.gz -s t2_seg.nii.gz -c t2 -discfile t2_totalspines
 # ======================================================================================================================
 
 # Compute cross-sectional area (CSA) of spinal cord and average it across levels C3 and C4
-sct_process_segmentation -i t2_seg.nii.gz -vert 3:4 -discfile t2_totalspineseg_discs.nii.gz -o csa_c3c4.csv
+sct_process_segmentation -i t2_seg.nii.gz -vert 2:3 -discfile t2_totalspineseg_discs.nii.gz -o csa_c2c3.csv
 # Aggregate CSA value per level (including new anat-based symmetry metrics)
-sct_process_segmentation -i t2_seg.nii.gz -anat t2.nii.gz -vert 3:4 -discfile t2_totalspineseg_discs.nii.gz -perlevel 1 -o csa_perlevel.csv
+sct_process_segmentation -i t2_seg.nii.gz -anat t2.nii.gz -vert 2:3 -discfile t2_totalspineseg_discs.nii.gz -perlevel 1 -o csa_perlevel.csv
 # Aggregate CSA value per slices
 sct_process_segmentation -i t2_seg.nii.gz -z 30:35 -discfile t2_totalspineseg_discs.nii.gz -perslice 1 -o csa_perslice.csv
 
@@ -96,11 +127,11 @@ sct_process_segmentation -i t2_seg.nii.gz -z 30:35 -discfile t2_totalspineseg_di
 # of the spinal cord will vary depending on the position of the neck.
 sct_detect_pmj -i t2.nii.gz -c t2 -qc ~/qc_singleSubj
 # Check the QC to make sure PMJ was properly detected, then compute CSA using the distance from the PMJ:
-sct_process_segmentation -i t2_seg.nii.gz -pmj t2_pmj.nii.gz -pmj-distance 64 -pmj-extent 30 -o csa_pmj.csv -qc ~/qc_singleSubj -qc-image t2.nii.gz
+sct_process_segmentation -i t2_seg.nii.gz -pmj t2_pmj.nii.gz -pmj-distance 60 -pmj-extent 30 -o csa_pmj.csv -qc ~/qc_singleSubj -qc-image t2.nii.gz
 
 # The above commands will output the metrics in the subject space (with the original image's slice numbers)
 # However, you can get the corresponding slice number in the PAM50 space by using the flag `-normalize-PAM50 1`
-sct_process_segmentation -i t2_seg.nii.gz -discfile t2_totalspineseg_discs.nii.gz -perslice 1 -normalize-PAM50 1 -o csa_PAM50.csv
+sct_process_segmentation -i t2_seg.nii.gz -discfile t2_totalspineseg_discs.nii.gz -perslice 1 -normalize-PAM50 1 -o csa_pam50.csv
 
 
 
@@ -124,6 +155,14 @@ sct_compute_compression -i t2_compressed_seg.nii.gz -vertfile t2_compressed_seg_
 # Compute ratio of AP diameter, normalized with healthy controls using `-normalize-hc 1`.
 sct_compute_compression -i t2_compressed_seg.nii.gz -vertfile t2_compressed_seg_labeled.nii.gz -l t2_compressed_labels-compression.nii.gz -metric diameter_AP -normalize-hc 1 -o ap_ratio_norm_PAM50.csv
 
+# Canal segmentation
+sct_deepseg canal -i t2_compressed.nii.gz -o t2_compressed_canal_seg.nii.gz -qc ~/qc_singleSubj
+# Check results using FSLeyes
+fsleyes t2.nii.gz -cm greyscale t2_canal_seg_seg.nii.gz -cm red -a 70.0 &
+# Compute aSCOR (Adapted Spinal Cord Occupation Ratio)
+# i.e. Spinal cord to canal ratio using the canal seg
+sct_compute_ascor -i-SC t2_compressed_seg.nii.gz -i-canal t2_compressed_canal_seg.nii.gz -perlevel 1 -o ascor.csv
+
 
 
 # Registration to template
@@ -181,7 +220,7 @@ sct_create_mask -i mt1.nii.gz -p centerline,mt1_seg.nii.gz -size 35mm -f cylinde
 
 # Register template->mt1. The flag -initwarp ../t2/warp_template2anat.nii.gz initializes the registration using the
 # template->t2 transformation which was previously estimated
-sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t2.nii.gz -iseg "${SCT_DIR}"/data/PAM50/template/PAM50_cord.nii.gz -d mt1.nii.gz -dseg mt1_seg.nii.gz -m mask_mt1.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -owarp warp_template2mt.nii.gz -qc ~/qc_singleSubj
+sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t2.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_cord.nii.gz" -d mt1.nii.gz -dseg mt1_seg.nii.gz -m mask_mt1.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -owarp warp_template2mt.nii.gz -qc ~/qc_singleSubj
 # Tips: Here we only use the segmentations (type=seg) to minimize the sensitivity of the registration procedure to
 #       image artifacts.
 # Tips: Step 1: algo=centermass to align source and destination segmentations, then Step 2: algo=bpslinesyn to adapt the
@@ -194,7 +233,7 @@ sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t2.nii.gz -ise
 # sct_register_to_template -i mt1.nii.gz -s mt1_seg.nii.gz -ldisc label_c3c4.nii.gz -ref subject -param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=bsplinesyn,slicewise=1
 
 # Warp template
-sct_warp_template -d mt1.nii.gz -w warp_template2mt.nii.gz -a 1 -qc ~/qc_singleSubj
+sct_warp_template -d mt1.nii.gz -w warp_template2mt.nii.gz -a 0 -qc ~/qc_singleSubj
 # Check results using FSLeyes
 fsleyes mt1.nii.gz -cm greyscale -a 100.0 label/template/PAM50_t2.nii.gz -cm greyscale -dr 0 4000 -a 100.0 label/template/PAM50_gm.nii.gz -cm red-yellow -dr 0.4 1 -a 50.0 label/template/PAM50_wm.nii.gz -cm blue-lightblue -dr 0.4 1 -a 50.0 &
 
@@ -204,7 +243,9 @@ fsleyes mt1.nii.gz -cm greyscale -a 100.0 label/template/PAM50_t2.nii.gz -cm gre
 # ======================================================================================================================
 
 # Register mt0->mt1 using z-regularized slicewise translations (algo=slicereg)
-# Note: Segmentation and mask can be re-used from "MT registration" section
+# Note: Segmentation and mask can be re-used from "Register MT to PAM50" section (or reproduced using commands below)
+# sct_deepseg spinalcord -i mt1.nii.gz -qc ~/qc_singleSubj
+# sct_create_mask -i mt1.nii.gz -p centerline,mt1_seg.nii.gz -size 35mm -f cylinder -o mask_mt1.nii.gz
 sct_register_multimodal -i mt0.nii.gz -d mt1.nii.gz -dseg mt1_seg.nii.gz -m mask_mt1.nii.gz -param step=1,type=im,algo=slicereg,metric=CC -x spline -qc ~/qc_singleSubj
 # Check results using FSLeyes
 fsleyes mt1.nii.gz mt0_reg.nii.gz &
@@ -214,11 +255,42 @@ sct_compute_mtr -mt0 mt0_reg.nii.gz -mt1 mt1.nii.gz
 
 
 
+# Registering additional contrasts (contrast-agnostic registration with deep learning between T2 and T1)
+# ======================================================================================================================
+
+# Go to T2 folder
+cd ../t2
+
+# Segment the spinal cord on T2-weighted data
+sct_deepseg spinalcord -i t2.nii.gz -qc ~/qc_singleSubj
+# Create a mask around the spinal cord
+sct_create_mask -i t2.nii.gz -p centerline,t2_seg.nii.gz -size 35mm -o mask_t2.nii.gz
+# Crop around the spinal cord to speed up and improve the subsequent registration
+sct_crop_image -i t2.nii.gz -m mask_t2.nii.gz
+
+# Go to T1 folder
+cd ../t1
+
+# Segment the spinal cord on T1-weighted data
+sct_deepseg spinalcord -i t1.nii.gz -qc ~/qc_singleSubj
+# Create a mask around the spinal cord
+sct_create_mask -i t1.nii.gz -p centerline,t1_seg.nii.gz -size 35mm -o mask_t1.nii.gz
+# Crop around the spinal cord to speed up and improve the subsequent registration
+sct_crop_image -i t1.nii.gz -m mask_t1.nii.gz
+
+# Register T1->T2 using the contrast-agnostic deep learning method
+# Note: We use the cropped images to focus the registration on the spinal cord region
+# Note: The destination segmentation is provided for QC reporting only
+sct_register_multimodal -i t1_crop.nii.gz -d ../t2/t2_crop.nii.gz -param step=1,type=im,algo=dl -qc ~/qc_singleSubj -dseg ../t2/t2_seg.nii.gz
+
+
+
 # Registering additional contrasts (T2 lumbar data)
 # ======================================================================================================================
 cd ../t2_lumbar
 
 # Use lumbar-specific `sct_deepseg` model to segment the spinal cord
+sct_deepseg sc_lumbar_t2 -install
 sct_deepseg sc_lumbar_t2 -i t2_lumbar.nii.gz -qc ~/qc_singleSubj
 
 # Generate labels for the 2 spinal cord landmarks: cauda equinea ('99') and T9-T10 disc ('17')
@@ -250,7 +322,7 @@ sct_deepseg spinalcord -i t2s.nii.gz -qc ~/qc_singleSubj
 # Subtract GM segmentation from cord segmentation to obtain WM segmentation
 # Note that we use the flag -thr 0 in case some voxels in the GM segmentation are *not* included in the cord
 # segmentation. That would results in voxels in the WM segmentation having the value “-1”, which would cause issues
-# with the registration. 
+# with the registration.
 sct_maths -i t2s_seg.nii.gz -sub t2s_gmseg.nii.gz -thr 0 -o t2s_wmseg.nii.gz
 
 
@@ -261,8 +333,8 @@ sct_maths -i t2s_seg.nii.gz -sub t2s_gmseg.nii.gz -thr 0 -o t2s_wmseg.nii.gz
 # Compute cross-sectional area (CSA) of the gray and white matter for all slices in the volume.
 # Note: Here we use the flag -angle-corr 0, because we do not want to correct the computed CSA by the cosine of the
 # angle between the cord centerline and the S-I axis: we assume that slices were acquired orthogonally to the cord.
-sct_process_segmentation -i t2s_wmseg.nii.gz -o csa_wm.csv -perslice 1 -angle-corr 0
-sct_process_segmentation -i t2s_gmseg.nii.gz -o csa_gm.csv -perslice 1 -angle-corr 0
+sct_process_segmentation -i t2s_wmseg.nii.gz -o csa_wm_perslice.csv -perslice 1 -angle-corr 0
+sct_process_segmentation -i t2s_gmseg.nii.gz -o csa_gm_perslice.csv -perslice 1 -angle-corr 0
 
 # You can also use the binary masks to extract signal intensity from MRI data.
 # The example below will show how to use the GM and WM segmentations to quantify T2* signal intensity, as done in
@@ -279,16 +351,19 @@ sct_extract_metric -i t2s.nii.gz -f t2s_gmseg.nii.gz -method bin -z 2:12 -o t2s_
 
 # Register template->t2s (using warping field generated from template<->t2 registration)
 # Tips: Here we use the WM seg for the iseg/dseg fields in order to account for both the cord and the GM shape.
-sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t2s.nii.gz -iseg "${SCT_DIR}"/data/PAM50/template/PAM50_wm.nii.gz -d t2s.nii.gz -dseg t2s_wmseg.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2t2s.nii.gz -owarpinv warp_t2s2template.nii.gz -param step=1,type=seg,algo=rigid:step=2,type=seg,metric=CC,algo=bsplinesyn,slicewise=1,iter=3:step=3,type=im,metric=CC,algo=syn,slicewise=1,iter=2 -qc ~/qc_singleSubj
+sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t2s.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_wm.nii.gz" -d t2s.nii.gz -dseg t2s_wmseg.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2t2s.nii.gz -owarpinv warp_t2s2template.nii.gz -param step=1,type=seg,algo=rigid:step=2,type=seg,metric=CC,algo=bsplinesyn,slicewise=1,iter=3:step=3,type=im,metric=CC,algo=syn,slicewise=1,iter=2 -qc ~/qc_singleSubj
 # Warp template
 sct_warp_template -d t2s.nii.gz -w warp_template2t2s.nii.gz -qc ~/qc_singleSubj
+# Compute vertebral level-based metrics using warped template (needed for the template's vertlevel file)
+sct_process_segmentation -i t2s_gmseg.nii.gz -vert 2:5 -perlevel 1 -o csa_gm.csv -centerline t2s_seg.nii.gz -centerline-exclude-missing 1
+sct_process_segmentation -i t2s_wmseg.nii.gz -vert 2:5 -perlevel 1 -o csa_wm.csv -centerline t2s_seg.nii.gz -centerline-exclude-missing 1
 
 # Register another metric while reusing newly-created GM-informed warping fields
 cd ../mt
 # Register template->mt using `-initwarp` with t2s to account for GM segmentation
-sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t2.nii.gz -iseg "${SCT_DIR}"/data/PAM50/template/PAM50_cord.nii.gz -d mt1.nii.gz -dseg mt1_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -m mask_mt1.nii.gz -initwarp ../t2s/warp_template2t2s.nii.gz -owarp warp_template2mt.nii.gz -qc ~/qc_singleSubj
+sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t2.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_cord.nii.gz" -d mt1.nii.gz -dseg mt1_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -m mask_mt1.nii.gz -initwarp ../t2s/warp_template2t2s.nii.gz -owarp warp_template2mt.nii.gz -qc ~/qc_singleSubj
 # Warp template
-sct_warp_template -d mt1.nii.gz -w warp_template2mt.nii.gz -qc ~/qc_singleSubj
+sct_warp_template -d mt1.nii.gz -w warp_template2mt.nii.gz -a 1 -qc ~/qc_singleSubj
 # Check results
 fsleyes mt1.nii.gz -cm greyscale -a 100.0 label/template/PAM50_t2.nii.gz -cm greyscale -dr 0 4000 -a 100.0 label/template/PAM50_gm.nii.gz -cm red-yellow -dr 0.4 1 -a 100.0 label/template/PAM50_wm.nii.gz -cm blue-lightblue -dr 0.4 1 -a 100.0 &
 
@@ -299,7 +374,7 @@ fsleyes mt1.nii.gz -cm greyscale -a 100.0 label/template/PAM50_t2.nii.gz -cm gre
 
 # Extract MTR for each slice within the white matter (combined label: #51)
 # Tips: To list all available labels, type: "sct_extract_metric"
-sct_extract_metric -i mtr.nii.gz -f label/atlas -method map -l 51 -o mtr_in_wm.csv
+sct_extract_metric -i mtr.nii.gz -f label/atlas -method map -l 51 -vert 2:5 -o mtr_in_wm.csv
 
 # Extract MTR within the right and left corticospinal tract and aggregate across specific slices
 sct_extract_metric -i mtr.nii.gz -f label/atlas -method map -l 4,5 -z 5:15 -o mtr_in_cst.csv
@@ -315,7 +390,7 @@ sct_extract_metric -i mtr.nii.gz -f label/atlas -method map -l 53 -vert 2:4 -ver
 cd ../dmri
 # Preprocessing steps
 # Compute mean dMRI from dMRI data
-sct_dmri_separate_b0_and_dwi -i dmri.nii.gz -bvec bvecs.txt 
+sct_dmri_separate_b0_and_dwi -i dmri.nii.gz -bvec bvecs.txt
 # Segment SC on mean dMRI data
 # Note: This segmentation does not need to be accurate-- it is only used to create a mask around the cord
 sct_deepseg spinalcord -i dmri_dwi_mean.nii.gz -qc ~/qc_singleSubj
@@ -334,7 +409,7 @@ sct_deepseg spinalcord -i dmri_moco_dwi_mean.nii.gz -qc ~/qc_singleSubj
 #       -param, so it will not make a difference here)
 # Note: the flag “-initwarpinv" provides a transformation dmri->template, in case you would like to bring all your DTI
 #       metrics in the PAM50 space (e.g. group averaging of FA maps)
-sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t1.nii.gz -iseg "${SCT_DIR}"/data/PAM50/template/PAM50_cord.nii.gz -d dmri_moco_dwi_mean.nii.gz -dseg dmri_moco_dwi_mean_seg.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2dmri.nii.gz -owarpinv warp_dmri2template.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -qc ~/qc_singleSubj
+sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t1.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_cord.nii.gz" -d dmri_moco_dwi_mean.nii.gz -dseg dmri_moco_dwi_mean_seg.nii.gz -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2dmri.nii.gz -owarpinv warp_dmri2template.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -qc ~/qc_singleSubj
 # Warp template (so 'label/atlas' can be used to extract metrics)
 sct_warp_template -d dmri_moco_dwi_mean.nii.gz -w warp_template2dmri.nii.gz -qc ~/qc_singleSubj
 # Check results in the QC report
@@ -345,12 +420,18 @@ sct_dmri_compute_dti -i dmri_moco.nii.gz -bval bvals.txt -bvec bvecs.txt
 
 # Compute FA within the white matter from individual level 2 to 5
 sct_extract_metric -i dti_FA.nii.gz -f label/atlas -l 51 -method map -vert 2:5 -vertfile label/template/PAM50_levels.nii.gz -perlevel 1 -o fa_in_wm.csv
+# Compute FA within the CST, aggregated across z slices, using the weighted average method
+sct_extract_metric -i dti_FA.nii.gz -f label/atlas -l 4,5 -method wa -z 2:14 -o fa_in_cst.csv
 
 
 
 # Functional MRI
 # ======================================================================================================================
 
+# the T2 segmentation will be reused, but it can also be generated using the commands below:
+# cd ../t2
+# sct_deepseg spinalcord -i t2.nii.gz
+
 cd ../fmri
 # Preprocessing steps
 # Average all fMRI time series to make it a 3D volume (needed by the next command)
@@ -364,18 +445,18 @@ sct_create_mask -i fmri.nii.gz -p centerline,t2_seg_reg.nii.gz -size 35mm -f cyl
 sct_fmri_moco -i fmri.nii.gz -m mask_fmri.nii.gz -qc ~/qc_singleSubj -qc-seg t2_seg_reg.nii.gz
 
 # Cord segmentation on motion-corrected averaged time series
-sct_deepseg spinalcord -i fmri_moco_mean.nii.gz -qc ~/qc_singleSubj/
+sct_deepseg spinalcord -i fmri_moco_mean.nii.gz -qc ~/qc_singleSubj
 # TSNR before/after motion correction with QC report
 sct_fmri_compute_tsnr -i fmri.nii.gz
 sct_fmri_compute_tsnr -i fmri_moco.nii.gz
-sct_qc -i fmri_tsnr.nii.gz -d fmri_moco_tsnr.nii.gz -s fmri_moco_mean_seg.nii.gz -p sct_fmri_compute_tsnr -qc ~/qc_singleSubj/
+sct_qc -i fmri_tsnr.nii.gz -d fmri_moco_tsnr.nii.gz -s fmri_moco_mean_seg.nii.gz -p sct_fmri_compute_tsnr -qc ~/qc_singleSubj
 
 # Register the template to the fMRI scan.
 # Note: here we don't rely on the segmentation because it is difficult to obtain one automatically. Instead, we rely on
 #       ANTs_SyN superpower to find a suitable transformation between the PAM50_t2s and the fMRI scan. We don't want to
 #       put too many iterations because this registration is very sensitive to the artifacts (drop out) in the image.
 #       Also, we want a 3D transformation (not 2D) because we need the through-z regularization.
-sct_register_multimodal -i "${SCT_DIR}"/data/PAM50/template/PAM50_t2s.nii.gz -iseg "${SCT_DIR}"/data/PAM50/template/PAM50_cord.nii.gz -d fmri_moco_mean.nii.gz -dseg fmri_moco_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3:step=3,type=im,algo=syn,metric=CC,iter=3,slicewise=1 -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2fmri.nii.gz -owarpinv warp_fmri2template.nii.gz -qc ~/qc_singleSubj
+sct_register_multimodal -i "${SCT_DIR}/data/PAM50/template/PAM50_t2s.nii.gz" -iseg "${SCT_DIR}/data/PAM50/template/PAM50_cord.nii.gz" -d fmri_moco_mean.nii.gz -dseg fmri_moco_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3:step=3,type=im,algo=syn,metric=CC,iter=3,slicewise=1 -initwarp ../t2/warp_template2anat.nii.gz -initwarpinv ../t2/warp_anat2template.nii.gz -owarp warp_template2fmri.nii.gz -owarpinv warp_fmri2template.nii.gz -qc ~/qc_singleSubj
 # Check results in the QC report
 
 # Warp template with the spinal levels (can be found at $SCT_DIR/data/PAM50/template/)
@@ -386,8 +467,8 @@ sct_warp_template -d fmri_moco_mean.nii.gz -w warp_template2fmri.nii.gz -a 0 -qc
 # ======================================================================================================================
 
 cd ../t1
-# Segment T1-weighted image (to be used in later steps)
-sct_deepseg spinalcord -i t1.nii.gz -qc ~/qc_singleSubj/
+# Segment T1-weighted image (but, reuse the previous segmentation to save processing time)
+## sct_deepseg spinalcord -i t1.nii.gz -qc ~/qc_singleSubj
 
 # Smooth spinal cord along centerline (extracted from the segmentation)
 sct_smooth_spinalcord -i t1.nii.gz -s t1_seg.nii.gz
@@ -424,6 +505,9 @@ sct_analyze_lesion -m t2_lesion_seg.nii.gz -s t2_sc_seg.nii.gz -qc ~/qc_singleSu
 sct_warp_template -d t2.nii.gz -w ../t2/warp_template2anat.nii.gz
 sct_analyze_lesion -m t2_lesion_seg.nii.gz -s t2_sc_seg.nii.gz -f label -qc ~/qc_singleSubj
 
+# You can also use the legacy method if the new methods fail for your data (`-c t2s` is also supported)
+sct_deepseg_lesion -i t2.nii.gz -c t2
+
 # Segment the spinal cord on gradient echo EPI data
 cd ../fmri/
 # Crop extraneous tissue using the t2-based mask generated earlier
@@ -431,20 +515,13 @@ sct_crop_image -i fmri_moco_mean.nii.gz -m mask_fmri.nii.gz -b 0
 # Segment the cord using the cropped image
 sct_deepseg sc_epi -i fmri_moco_mean_crop.nii.gz -qc ~/qc_singleSubj
 
-# Canal segmentation
-cd ../t2
-sct_deepseg sc_canal_t2 -i t2.nii.gz -qc ~/qc_singleSubj
-# Check results using FSLeyes
-fsleyes t2.nii.gz -cm greyscale t2_canal_seg_seg.nii.gz -cm red -a 70.0 &
-
-# Compute aSCOR (Adapted Spinal Cord Occupation Ratio)
-# i.e. Spinal cord to canal ratio using the canal seg
-sct_compute_ascor -i-SC t2_seg.nii.gz -i-canal t2_canal_seg.nii.gz -perlevel 1 -o ascor.csv
-
 # Segment the spinal nerve rootlets
-sct_deepseg rootlets -i t2.nii.gz -qc ~/qc_singleSubj
+cd ../t2/
+sct_deepseg rootlets -i t2.nii.gz -o t2_rootlets.nii.gz -qc ~/qc_singleSubj
 # Check results using FSLeyes
 fsleyes t2.nii.gz -cm greyscale t2_rootlets.nii.gz -cm subcortical -a 70.0 &
+# Rootlets-based registration
+sct_register_to_template -i t2.nii.gz -s t2_seg.nii.gz -lrootlet t2_rootlets.nii.gz -c t2 -ofolder rootlets-reg -qc ~/qc_singleSubj
 
 # Multiple sclerosis lesion segmentation on T2-weighted images
 cd ../t2_ms/
@@ -457,3 +534,23 @@ fsleyes t2.nii.gz -cm greyscale t2_lesion_seg.nii.gz -cm red -a 70.0 &
 
 # Return to parent directory
 cd ..
+
+
+# Display results (to easily compare integrity across SCT versions)
+# ===========================================================================================
+set +v
+end=$(date +%s)
+runtime=$((end-start))
+echo "~~~"  # these are used to format as code when copy/pasting in github's markdown
+echo "Version:         $(sct_version)"
+echo "Ran on:          $(uname -nsr)"
+echo "Duration:        $((runtime / 3600))hrs $(((runtime / 60) % 60))min $((runtime % 60))sec"
+echo "---"
+# The file `test_batch_processing.py` will output tested values when run as a script
+"$SCT_DIR"/python/envs/venv_sct/bin/python "$SCT_DIR"/testing/batch_single_subject/test_batch_single_subject.py ||
+"$SCT_DIR"/python/envs/venv_sct/python.exe "$SCT_DIR"/testing/batch_single_subject/test_batch_single_subject.py
+echo "~~~"
+
+# Display syntax to open QC report on web browser
+echo "To open Quality Control (QC) report on a web-browser, run the following:"
+echo "$open_command $SCT_BP_QC_FOLDER/index.html"