csttool.metrics

CST metric computation and PDF/HTML reporting. The CLI counterpart is metrics.

from csttool.metrics import compute_metrics

compute_metrics(
    cst_left="extract/cst_left.trk",
    cst_right="extract/cst_right.trk",
    fa="tracking/dti_FA.nii.gz",
    out="./metrics",
    generate_pdf=True,
)

csttool.metrics

CST Metrics Module

Modular metrics computation for bilateral CST analysis.

Structure: - modules/unilateral_analysis: Analyze individual hemispheres - modules/bilateral_analysis: Compare left vs right - modules/visualizations: Generate all plots

Functions

analyze_cst_hemisphere(streamlines, fa_map=None, md_map=None, rd_map=None, ad_map=None, affine=None, hemisphere='unknown')

Complete analysis of a single CST hemisphere.

Parameters:

Name	Type	Description	Default
streamlines	Streamlines	Streamlines for one hemisphere (left OR right)	required
fa_map	ndarray	3D fractional anisotropy map	None
md_map	ndarray	3D mean diffusivity map	None
rd_map	ndarray	3D radial diffusivity map	None
ad_map	ndarray	3D axial diffusivity map	None
affine	ndarray	4x4 affine transformation matrix	None
hemisphere	str	'left' or 'right' for identification	'unknown'

Returns:

Name	Type	Description
metrics	dict	Comprehensive metrics dictionary containing: - morphology: streamline count, length stats, volume - fa: mean, std, median, profile, all sampled values - md: mean, std, median, profile, all sampled values - rd: mean, std, median, profile, all sampled values - ad: mean, std, median, profile, all sampled values - hemisphere: identification string

Source code in src/csttool/metrics/modules/unilateral_analysis.py

def analyze_cst_hemisphere(
    streamlines,
    fa_map=None,
    md_map=None,
    rd_map=None,
    ad_map=None,
    affine=None,
    hemisphere='unknown'
):
    """
    Complete analysis of a single CST hemisphere.

    Parameters
    ----------
    streamlines : Streamlines
        Streamlines for one hemisphere (left OR right)
    fa_map : ndarray, optional
        3D fractional anisotropy map
    md_map : ndarray, optional
        3D mean diffusivity map
    rd_map : ndarray, optional
        3D radial diffusivity map
    ad_map : ndarray, optional
        3D axial diffusivity map
    affine : ndarray, optional
        4x4 affine transformation matrix
    hemisphere : str
        'left' or 'right' for identification

    Returns
    -------
    metrics : dict
        Comprehensive metrics dictionary containing:
        - morphology: streamline count, length stats, volume
        - fa: mean, std, median, profile, all sampled values
        - md: mean, std, median, profile, all sampled values
        - rd: mean, std, median, profile, all sampled values
        - ad: mean, std, median, profile, all sampled values
        - hemisphere: identification string
    """

    print(f"\nAnalyzing {hemisphere.upper()} CST...")

    metrics = {
        'hemisphere': hemisphere,
        'morphology': compute_morphology(streamlines, affine)
    }

    # Microstructural analysis requires affine
    if fa_map is not None and affine is not None:
        fa_values = sample_scalar_along_tract(streamlines, fa_map, affine)
        if len(fa_values) > 0:
            fa_profile = compute_tract_profile(streamlines, fa_map, affine, n_points=20)
            fa_localized = compute_localized_metrics(fa_profile)
            metrics['fa'] = {
                'mean': float(np.mean(fa_values)),
                'std': float(np.std(fa_values)),
                'median': float(np.median(fa_values)),
                'min': float(np.min(fa_values)),
                'max': float(np.max(fa_values)),
                'profile': fa_profile,
                'n_samples': len(fa_values),
                'pontine': fa_localized['pontine'],
                'plic': fa_localized['plic'],
                'precentral': fa_localized['precentral']
            }
            print(f"  FA: {metrics['fa']['mean']:.3f} ± {metrics['fa']['std']:.3f}")
        else:
            # Handle empty case
            metrics['fa'] = {
                'mean': 0.0, 'std': 0.0, 'median': 0.0, 'min': 0.0, 'max': 0.0,
                'profile': [], 'n_samples': 0,
                'pontine': 0.0, 'plic': 0.0, 'precentral': 0.0
            }

    if md_map is not None and affine is not None:
        md_values = sample_scalar_along_tract(streamlines, md_map, affine)
        if len(md_values) > 0:
            md_profile = compute_tract_profile(streamlines, md_map, affine, n_points=20)
            md_localized = compute_localized_metrics(md_profile)
            metrics['md'] = {
                'mean': float(np.mean(md_values)),
                'std': float(np.std(md_values)),
                'median': float(np.median(md_values)),
                'min': float(np.min(md_values)),
                'max': float(np.max(md_values)),
                'profile': md_profile,
                'n_samples': len(md_values),
                'pontine': md_localized['pontine'],
                'plic': md_localized['plic'],
                'precentral': md_localized['precentral']
            }
            print(f"  MD: {metrics['md']['mean']:.3e} ± {metrics['md']['std']:.3e}")
        else:
            # Handle empty case
            metrics['md'] = {
                'mean': 0.0, 'std': 0.0, 'median': 0.0, 'min': 0.0, 'max': 0.0,
                'profile': [], 'n_samples': 0,
                'pontine': 0.0, 'plic': 0.0, 'precentral': 0.0
            }

    if rd_map is not None and affine is not None:
        rd_values = sample_scalar_along_tract(streamlines, rd_map, affine)
        if len(rd_values) > 0:
            rd_profile = compute_tract_profile(streamlines, rd_map, affine, n_points=20)
            rd_localized = compute_localized_metrics(rd_profile)
            metrics['rd'] = {
                'mean': float(np.mean(rd_values)),
                'std': float(np.std(rd_values)),
                'median': float(np.median(rd_values)),
                'min': float(np.min(rd_values)),
                'max': float(np.max(rd_values)),
                'profile': rd_profile,
                'n_samples': len(rd_values),
                'pontine': rd_localized['pontine'],
                'plic': rd_localized['plic'],
                'precentral': rd_localized['precentral']
            }
            print(f"  RD: {metrics['rd']['mean']:.3e} ± {metrics['rd']['std']:.3e}")
        else:
            metrics['rd'] = {
                'mean': 0.0, 'std': 0.0, 'median': 0.0, 'min': 0.0, 'max': 0.0,
                'profile': [], 'n_samples': 0,
                'pontine': 0.0, 'plic': 0.0, 'precentral': 0.0
            }

    if ad_map is not None and affine is not None:
        ad_values = sample_scalar_along_tract(streamlines, ad_map, affine)
        if len(ad_values) > 0:
            ad_profile = compute_tract_profile(streamlines, ad_map, affine, n_points=20)
            ad_localized = compute_localized_metrics(ad_profile)
            metrics['ad'] = {
                'mean': float(np.mean(ad_values)),
                'std': float(np.std(ad_values)),
                'median': float(np.median(ad_values)),
                'min': float(np.min(ad_values)),
                'max': float(np.max(ad_values)),
                'profile': ad_profile,
                'n_samples': len(ad_values),
                'pontine': ad_localized['pontine'],
                'plic': ad_localized['plic'],
                'precentral': ad_localized['precentral']
            }
            print(f"  AD: {metrics['ad']['mean']:.3e} ± {metrics['ad']['std']:.3e}")
        else:
            metrics['ad'] = {
                'mean': 0.0, 'std': 0.0, 'median': 0.0, 'min': 0.0, 'max': 0.0,
                'profile': [], 'n_samples': 0,
                'pontine': 0.0, 'plic': 0.0, 'precentral': 0.0
            }

    return metrics

compute_morphology(streamlines, affine)

Compute morphological properties of a streamline bundle.

Parameters:

Name	Type	Description	Default
streamlines	Streamlines	Input streamlines	required
affine	ndarray	4x4 affine transformation matrix	required

Returns:

Name	Type	Description
morphology	dict	Dictionary containing: - n_streamlines: number of streamlines - mean_length: average streamline length in mm - std_length: standard deviation of lengths - min_length: minimum streamline length - max_length: maximum streamline length - tract_volume: volume in mm³

Source code in src/csttool/metrics/modules/unilateral_analysis.py

def compute_morphology(streamlines, affine):
    """
    Compute morphological properties of a streamline bundle.

    Parameters
    ----------
    streamlines : Streamlines
        Input streamlines
    affine : ndarray
        4x4 affine transformation matrix

    Returns
    -------
    morphology : dict
        Dictionary containing:
        - n_streamlines: number of streamlines
        - mean_length: average streamline length in mm
        - std_length: standard deviation of lengths
        - min_length: minimum streamline length
        - max_length: maximum streamline length
        - tract_volume: volume in mm³
    """

    if len(streamlines) == 0:
        return {
            'n_streamlines': 0,
            'mean_length': 0.0,
            'std_length': 0.0,
            'min_length': 0.0,
            'max_length': 0.0,
            'tract_volume': 0.0
        }

    # Compute streamline lengths
    lengths = np.array([length(s) for s in streamlines])

    # Compute tract volume by counting unique voxels
    # Get all points from all streamlines
    all_points = np.vstack(streamlines)

    # Transform to voxel coordinates using inverse affine
    inv_affine = np.linalg.inv(affine)
    all_points_hom = np.c_[all_points, np.ones(len(all_points))]
    voxel_coords = (inv_affine @ all_points_hom.T).T[:, :3]

    # Round to integer voxel indices
    voxel_indices = np.round(voxel_coords).astype(int)

    # Count unique voxels
    unique_voxels = np.unique(voxel_indices, axis=0)
    n_voxels = len(unique_voxels)

    # Compute voxel volume
    voxel_size = np.sqrt(np.sum(affine[:3, :3]**2, axis=0))
    voxel_volume = np.prod(voxel_size)

    # Total tract volume
    tract_volume = n_voxels * voxel_volume

    morphology = {
        'n_streamlines': len(streamlines),
        'mean_length': float(np.mean(lengths)),
        'std_length': float(np.std(lengths)),
        'min_length': float(np.min(lengths)),
        'max_length': float(np.max(lengths)),
        'tract_volume': float(tract_volume)
    }

    print(f"  Morphology: {morphology['n_streamlines']} streamlines, "
          f"{morphology['mean_length']:.1f} mm mean length, "
          f"{morphology['tract_volume']:.0f} mm³ volume")

    return morphology

sample_scalar_along_tract(streamlines, scalar_map, affine)

Sample scalar values at every point along all streamlines.

Parameters:

Name	Type	Description	Default
streamlines	Streamlines	Input streamlines in world coordinates (mm)	required
scalar_map	ndarray	3D scalar map (e.g., FA or MD)	required
affine	ndarray	4x4 affine transformation matrix	required

Returns:

Name	Type	Description
scalar_values	ndarray	Array of all sampled scalar values (flattened across all streamlines)

Source code in src/csttool/metrics/modules/unilateral_analysis.py

def sample_scalar_along_tract(streamlines, scalar_map, affine):
    """
    Sample scalar values at every point along all streamlines.

    Parameters
    ----------
    streamlines : Streamlines
        Input streamlines in world coordinates (mm)
    scalar_map : ndarray
        3D scalar map (e.g., FA or MD)
    affine : ndarray
        4x4 affine transformation matrix

    Returns
    -------
    scalar_values : ndarray
        Array of all sampled scalar values (flattened across all streamlines)
    """

    if len(streamlines) == 0:
        return np.array([])

    scalar_values = []

    for streamline in streamlines:
        for point in streamline:
            # Convert world coordinates to voxel coordinates
            voxel_coord = world_to_voxel(point, affine)

            # Check bounds
            if (0 <= voxel_coord[0] < scalar_map.shape[0] and
                0 <= voxel_coord[1] < scalar_map.shape[1] and
                0 <= voxel_coord[2] < scalar_map.shape[2]):

                scalar_value = scalar_map[voxel_coord[0], voxel_coord[1], voxel_coord[2]]
                scalar_values.append(scalar_value)

    return np.array(scalar_values)

compute_tract_profile(streamlines, scalar_map, affine, n_points=20)

Compute normalized tract profile (average scalar along tract).

This function samples scalar values along each streamline, normalizes each streamline to the same number of points, and averages across all streamlines to create a representative profile.

Parameters:

Name	Type	Description	Default
streamlines	Streamlines	Input streamlines	required
scalar_map	ndarray	3D scalar map (FA or MD)	required
affine	ndarray	4x4 affine transformation matrix	required
n_points	int	Number of points in output profile (default: 20)	20

Returns:

Name	Type	Description
profile	ndarray	Average scalar values at n_points positions along the tract

Source code in src/csttool/metrics/modules/unilateral_analysis.py

def compute_tract_profile(streamlines, scalar_map, affine, n_points=20):
    """
    Compute normalized tract profile (average scalar along tract).

    This function samples scalar values along each streamline, normalizes
    each streamline to the same number of points, and averages across all
    streamlines to create a representative profile.

    Parameters
    ----------
    streamlines : Streamlines
        Input streamlines
    scalar_map : ndarray
        3D scalar map (FA or MD)
    affine : ndarray
        4x4 affine transformation matrix
    n_points : int
        Number of points in output profile (default: 20)

    Returns
    -------
    profile : ndarray
        Average scalar values at n_points positions along the tract
    """

    if len(streamlines) == 0:
        return np.zeros(n_points)

    all_profiles = []

    for streamline in streamlines:
        if len(streamline) < 2:
            continue

        # Sample scalar values at each point
        streamline_scalars = []
        for point in streamline:
            voxel_coord = world_to_voxel(point, affine)

            if (0 <= voxel_coord[0] < scalar_map.shape[0] and
                0 <= voxel_coord[1] < scalar_map.shape[1] and
                0 <= voxel_coord[2] < scalar_map.shape[2]):

                scalar_value = scalar_map[voxel_coord[0], voxel_coord[1], voxel_coord[2]]
                streamline_scalars.append(scalar_value)

        if len(streamline_scalars) < 5:  # Need minimum points
            continue

        # Normalize to n_points
        if len(streamline_scalars) >= n_points:
            # Downsample
            indices = np.linspace(0, len(streamline_scalars)-1, n_points).astype(int)
            normalized_profile = np.array(streamline_scalars)[indices]
        else:
            # Upsample with interpolation
            x_original = np.linspace(0, 1, len(streamline_scalars))
            x_target = np.linspace(0, 1, n_points)
            normalized_profile = np.interp(x_target, x_original, streamline_scalars)

        all_profiles.append(normalized_profile)

    if len(all_profiles) == 0:
        return np.zeros(n_points)

    # Average across all streamlines
    final_profile = np.mean(all_profiles, axis=0)

    return final_profile.tolist()  # Convert to list for JSON serialization

print_hemisphere_summary(metrics)

Print human-readable summary of hemisphere metrics.

Parameters:

Name	Type	Description	Default
metrics	dict	Metrics dictionary from analyze_cst_hemisphere()	required

Source code in src/csttool/metrics/modules/unilateral_analysis.py

def print_hemisphere_summary(metrics):
    """
    Print human-readable summary of hemisphere metrics.

    Parameters
    ----------
    metrics : dict
        Metrics dictionary from analyze_cst_hemisphere()
    """

    hemisphere = metrics['hemisphere'].upper()
    print(f"\n{'='*60}")
    print(f"{hemisphere} CST ANALYSIS SUMMARY")
    print(f"{'='*60}")

    morph = metrics['morphology']
    print(f"Streamline Count: {morph['n_streamlines']}")
    print(f"Mean Length: {morph['mean_length']:.1f} ± {morph['std_length']:.1f} mm")
    print(f"Length Range: [{morph['min_length']:.1f}, {morph['max_length']:.1f}] mm")
    print(f"Tract Volume: {morph['tract_volume']:.0f} mm³")

    if 'fa' in metrics:
        fa = metrics['fa']
        print(f"\nFractional Anisotropy:")
        print(f"  Mean: {fa['mean']:.3f} ± {fa['std']:.3f}")
        print(f"  Range: [{fa['min']:.3f}, {fa['max']:.3f}]")
        print(f"  Samples: {fa['n_samples']}")

    if 'md' in metrics:
        md = metrics['md']
        print(f"\nMean Diffusivity:")
        print(f"  Mean: {md['mean']:.3e} ± {md['std']:.3e}")
        print(f"  Range: [{md['min']:.3e}, {md['max']:.3e}]")
        print(f"  Samples: {md['n_samples']}")

    print(f"{'='*60}\n")

compare_bilateral_cst(left_metrics, right_metrics)

Compare left and right CST metrics and compute asymmetry measures.

Parameters:

Name	Type	Description	Default
left_metrics	dict	Metrics from analyze_cst_hemisphere() for left CST	required
right_metrics	dict	Metrics from analyze_cst_hemisphere() for right CST	required

Returns:

Name	Type	Description
comparison	dict	Comprehensive bilateral comparison containing: - left: complete left hemisphere metrics - right: complete right hemisphere metrics - asymmetry: laterality indices and differences

Source code in src/csttool/metrics/modules/bilateral_analysis.py

def compare_bilateral_cst(left_metrics, right_metrics):
    """
    Compare left and right CST metrics and compute asymmetry measures.

    Parameters
    ----------
    left_metrics : dict
        Metrics from analyze_cst_hemisphere() for left CST
    right_metrics : dict
        Metrics from analyze_cst_hemisphere() for right CST

    Returns
    -------
    comparison : dict
        Comprehensive bilateral comparison containing:
        - left: complete left hemisphere metrics
        - right: complete right hemisphere metrics
        - asymmetry: laterality indices and differences
    """

    print("\nComputing bilateral comparison...")

    comparison = {
        'left': left_metrics,
        'right': right_metrics,
        'asymmetry': compute_laterality_indices(left_metrics, right_metrics)
    }

    print_bilateral_summary(comparison)

    return comparison

compute_laterality_indices(left_metrics, right_metrics)

Compute laterality indices for all available metrics.

Laterality Index (LI) = (L - R) / (L + R)

LI interpretation: - LI > 0: Left hemisphere larger/higher - LI < 0: Right hemisphere larger/higher
- LI ≈ 0: Symmetric - |LI| > 0.1: Potentially significant asymmetry

Parameters:

Name	Type	Description	Default
left_metrics	dict	Left hemisphere metrics	required
right_metrics	dict	Right hemisphere metrics	required

Returns:

Name	Type	Description
asymmetry	dict	Laterality indices and absolute differences for: - volume - streamline count - mean length - mean FA (if available) - mean MD (if available) - mean RD (if available) - mean AD (if available)

Source code in src/csttool/metrics/modules/bilateral_analysis.py

def compute_laterality_indices(left_metrics, right_metrics):
    """
    Compute laterality indices for all available metrics.

    Laterality Index (LI) = (L - R) / (L + R)

    LI interpretation:
    - LI > 0: Left hemisphere larger/higher
    - LI < 0: Right hemisphere larger/higher  
    - LI ≈ 0: Symmetric
    - |LI| > 0.1: Potentially significant asymmetry

    Parameters
    ----------
    left_metrics : dict
        Left hemisphere metrics
    right_metrics : dict
        Right hemisphere metrics

    Returns
    -------
    asymmetry : dict
        Laterality indices and absolute differences for:
        - volume
        - streamline count
        - mean length
        - mean FA (if available)
        - mean MD (if available)
        - mean RD (if available)
        - mean AD (if available)
    """

    asymmetry = {}

    # Morphological asymmetry
    left_morph = left_metrics['morphology']
    right_morph = right_metrics['morphology']

    asymmetry['volume'] = compute_li(
        left_morph['tract_volume'],
        right_morph['tract_volume']
    )

    asymmetry['streamline_count'] = compute_li(
        left_morph['n_streamlines'],
        right_morph['n_streamlines']
    )

    asymmetry['mean_length'] = compute_li(
        left_morph['mean_length'],
        right_morph['mean_length']
    )

    # Microstructural asymmetry (global)
    if 'fa' in left_metrics and 'fa' in right_metrics:
        asymmetry['fa'] = compute_li(
            left_metrics['fa']['mean'],
            right_metrics['fa']['mean']
        )

    if 'md' in left_metrics and 'md' in right_metrics:
        asymmetry['md'] = compute_li(
            left_metrics['md']['mean'],
            right_metrics['md']['mean']
        )

    if 'rd' in left_metrics and 'rd' in right_metrics:
        asymmetry['rd'] = compute_li(
            left_metrics['rd']['mean'],
            right_metrics['rd']['mean']
        )

    if 'ad' in left_metrics and 'ad' in right_metrics:
        asymmetry['ad'] = compute_li(
            left_metrics['ad']['mean'],
            right_metrics['ad']['mean']
        )

    # Localized microstructural asymmetry (per region)
    regions = ['pontine', 'plic', 'precentral']
    scalars = ['fa', 'md', 'rd', 'ad']

    for scalar in scalars:
        if scalar in left_metrics and scalar in right_metrics:
            for region in regions:
                if region in left_metrics[scalar] and region in right_metrics[scalar]:
                    key = f'{scalar}_{region}'
                    asymmetry[key] = compute_li(
                        left_metrics[scalar][region],
                        right_metrics[scalar][region]
                    )

    return asymmetry

compute_li(left_value, right_value)

Compute laterality index and absolute difference.

Parameters:

Name	Type	Description	Default
left_value	float	Metric value from left hemisphere	required
right_value	float	Metric value from right hemisphere	required

Returns:

Name	Type	Description
li_info	dict	Dictionary containing: - laterality_index: (L-R)/(L+R) - absolute_difference: |L-R| - percent_difference: 100 * |L-R| / mean(L,R) - interpretation: text description

Source code in src/csttool/metrics/modules/bilateral_analysis.py

def compute_li(left_value, right_value):
    """
    Compute laterality index and absolute difference.

    Parameters
    ----------
    left_value : float
        Metric value from left hemisphere
    right_value : float
        Metric value from right hemisphere

    Returns
    -------
    li_info : dict
        Dictionary containing:
        - laterality_index: (L-R)/(L+R)
        - absolute_difference: |L-R|
        - percent_difference: 100 * |L-R| / mean(L,R)
        - interpretation: text description
    """

    total = left_value + right_value

    if total == 0:
        return {
            'laterality_index': 0.0,
            'absolute_difference': 0.0,
            'percent_difference': 0.0,
            'interpretation': 'no data'
        }

    li = (left_value - right_value) / total
    abs_diff = abs(left_value - right_value)
    mean_value = (left_value + right_value) / 2
    pct_diff = 100 * abs_diff / mean_value if mean_value > 0 else 0.0

    # Interpret laterality
    if abs(li) < 0.05:
        interpretation = 'symmetric'
    elif abs(li) < 0.10:
        interpretation = 'mild asymmetry'
    elif abs(li) < 0.20:
        interpretation = 'moderate asymmetry'
    else:
        interpretation = 'strong asymmetry'

    if li > 0:
        interpretation += ' (left > right)'
    elif li < 0:
        interpretation += ' (right > left)'

    return {
        'laterality_index': float(li),
        'absolute_difference': float(abs_diff),
        'percent_difference': float(pct_diff),
        'interpretation': interpretation,
        'left_value': float(left_value),
        'right_value': float(right_value)
    }

print_bilateral_summary(comparison)

Print human-readable bilateral comparison summary.

Parameters:

Name	Type	Description	Default
comparison	dict	Output from compare_bilateral_cst()	required

Source code in src/csttool/metrics/modules/bilateral_analysis.py

def print_bilateral_summary(comparison):
    """
    Print human-readable bilateral comparison summary.

    Parameters
    ----------
    comparison : dict
        Output from compare_bilateral_cst()
    """

    print("\n" + "=" * 60)
    print("BILATERAL CST COMPARISON")
    print("=" * 60)

    left = comparison['left']
    right = comparison['right']
    asym = comparison['asymmetry']

    # Morphology comparison
    print("\nMORPHOLOGY:")
    print(f"  Streamline Count:")
    print(f"    Left:  {left['morphology']['n_streamlines']}")
    print(f"    Right: {right['morphology']['n_streamlines']}")
    print(f"    LI:    {asym['streamline_count']['laterality_index']:+.3f} "
          f"({asym['streamline_count']['interpretation']})")

    print(f"\n  Tract Volume:")
    print(f"    Left:  {left['morphology']['tract_volume']:.0f} mm³")
    print(f"    Right: {right['morphology']['tract_volume']:.0f} mm³")
    print(f"    LI:    {asym['volume']['laterality_index']:+.3f} "
          f"({asym['volume']['interpretation']})")
    print(f"    Diff:  {asym['volume']['absolute_difference']:.0f} mm³ "
          f"({asym['volume']['percent_difference']:.1f}%)")

    print(f"\n  Mean Length:")
    print(f"    Left:  {left['morphology']['mean_length']:.1f} mm")
    print(f"    Right: {right['morphology']['mean_length']:.1f} mm")
    print(f"    LI:    {asym['mean_length']['laterality_index']:+.3f} "
          f"({asym['mean_length']['interpretation']})")

    # Microstructural comparison
    if 'fa' in asym:
        print(f"\nFRACTIONAL ANISOTROPY:")
        print(f"    Left:  {left['fa']['mean']:.3f} ± {left['fa']['std']:.3f}")
        print(f"    Right: {right['fa']['mean']:.3f} ± {right['fa']['std']:.3f}")
        print(f"    LI:    {asym['fa']['laterality_index']:+.3f} "
              f"({asym['fa']['interpretation']})")
        print(f"    Diff:  {asym['fa']['absolute_difference']:.3f} "
              f"({asym['fa']['percent_difference']:.1f}%)")

    if 'md' in asym:
        print(f"\nMEAN DIFFUSIVITY:")
        print(f"    Left:  {left['md']['mean']:.3e} ± {left['md']['std']:.3e}")
        print(f"    Right: {right['md']['mean']:.3e} ± {right['md']['std']:.3e}")
        print(f"    LI:    {asym['md']['laterality_index']:+.3f} "
              f"({asym['md']['interpretation']})")
        print(f"    Diff:  {asym['md']['absolute_difference']:.3e} "
              f"({asym['md']['percent_difference']:.1f}%)")

    print("=" * 60)

assess_clinical_significance(asymmetry, thresholds=None)

Assess clinical significance of asymmetry based on thresholds.

Parameters:

Name	Type	Description	Default
asymmetry	dict	Asymmetry metrics from compute_laterality_indices()	required
thresholds	dict	Custom thresholds for each metric Default thresholds are based on literature (if available)	None

Returns:

Name	Type	Description
assessment	dict	Clinical assessment with flags for significant asymmetries

Source code in src/csttool/metrics/modules/bilateral_analysis.py

def assess_clinical_significance(asymmetry, thresholds=None):
    """
    Assess clinical significance of asymmetry based on thresholds.

    Parameters
    ----------
    asymmetry : dict
        Asymmetry metrics from compute_laterality_indices()
    thresholds : dict, optional
        Custom thresholds for each metric
        Default thresholds are based on literature (if available)

    Returns
    -------
    assessment : dict
        Clinical assessment with flags for significant asymmetries
    """

    if thresholds is None:
        # Default thresholds (can be refined based on normative data)
        thresholds = {
            'volume': 0.15,        # 15% asymmetry is potentially significant
            'fa': 0.10,            # 10% FA asymmetry
            'md': 0.10,            # 10% MD asymmetry
            'streamline_count': 0.20  # 20% streamline asymmetry
        }

    assessment = {
        'flags': [],
        'severity': 'normal',
        'recommendations': []
    }

    # Check each metric
    for metric, threshold in thresholds.items():
        if metric in asymmetry:
            li = abs(asymmetry[metric]['laterality_index'])

            if li > threshold:
                assessment['flags'].append({
                    'metric': metric,
                    'laterality_index': asymmetry[metric]['laterality_index'],
                    'threshold': threshold,
                    'interpretation': asymmetry[metric]['interpretation']
                })

    # Determine overall severity
    if len(assessment['flags']) == 0:
        assessment['severity'] = 'normal'
        assessment['recommendations'].append('No significant asymmetries detected')
    elif len(assessment['flags']) == 1:
        assessment['severity'] = 'mild'
        assessment['recommendations'].append('Single metric shows asymmetry - consider follow-up')
    elif len(assessment['flags']) == 2:
        assessment['severity'] = 'moderate'
        assessment['recommendations'].append('Multiple metrics show asymmetry - recommend clinical correlation')
    else:
        assessment['severity'] = 'significant'
        assessment['recommendations'].append('Significant asymmetries detected - recommend further clinical evaluation')

    return assessment

plot_tract_profiles(left_metrics, right_metrics, output_dir, subject_id, scalar='fa', anatomical_labels=True)

Plot along-tract profiles for bilateral comparison.

Creates a figure showing FA or MD profiles along normalized tract length for both left and right CST.

Parameters:

Name	Type	Description	Default
left_metrics	dict	Left hemisphere metrics with 'fa' or 'md' profile	required
right_metrics	dict	Right hemisphere metrics with 'fa' or 'md' profile	required
output_dir	str or Path	Output directory for saving figure	required
subject_id	str	Subject identifier for filename	required
scalar	str	'fa' or 'md' - which scalar to plot	'fa'
anatomical_labels	bool	If True, add anatomical labels to x-axis (Pontine Level, PLIC, Precentral Gyrus)	True

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_tract_profiles(
    left_metrics,
    right_metrics,
    output_dir,
    subject_id,
    scalar='fa',
    anatomical_labels=True
):
    """
    Plot along-tract profiles for bilateral comparison.

    Creates a figure showing FA or MD profiles along normalized tract length
    for both left and right CST.

    Parameters
    ----------
    left_metrics : dict
        Left hemisphere metrics with 'fa' or 'md' profile
    right_metrics : dict
        Right hemisphere metrics with 'fa' or 'md' profile
    output_dir : str or Path
        Output directory for saving figure
    subject_id : str
        Subject identifier for filename
    scalar : str
        'fa' or 'md' - which scalar to plot
    anatomical_labels : bool
        If True, add anatomical labels to x-axis (Pontine Level, PLIC, Precentral Gyrus)

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    # Get profiles
    if scalar not in left_metrics or scalar not in right_metrics:
        print(f"Warning: {scalar.upper()} not available in metrics")
        return None

    left_profile = np.array(left_metrics[scalar]['profile'])
    right_profile = np.array(right_metrics[scalar]['profile'])

    n_points = len(left_profile)
    x = np.linspace(0, 100, n_points)  # Normalized position (0-100%)

    # Create figure
    fig, ax = plt.subplots(figsize=(10, 6))

    # Plot profiles
    ax.plot(x, left_profile, 'b-', linewidth=2, label='Left CST', marker='o', markersize=4)
    ax.plot(x, right_profile, 'r-', linewidth=2, label='Right CST', marker='s', markersize=4)

    # Add mean lines
    left_mean = left_metrics[scalar]['mean']
    right_mean = right_metrics[scalar]['mean']
    ax.axhline(left_mean, color='b', linestyle='--', alpha=0.5, label=f'Left mean: {left_mean:.3f}')
    ax.axhline(right_mean, color='r', linestyle='--', alpha=0.5, label=f'Right mean: {right_mean:.3f}')

    # Labels and formatting
    scalar_label = 'Fractional Anisotropy' if scalar == 'fa' else 'Mean Diffusivity (×10⁻³ mm²/s)'

    # Set x-axis with anatomical labels
    if anatomical_labels:
        ax.set_xticks([0, 50, 100])
        ax.set_xticklabels(['0%', '50%', '100%'])
        # Add anatomical labels as secondary text
        ax.text(0, -0.12, 'Pontine Level', transform=ax.get_xaxis_transform(), 
                ha='center', fontsize=9, style='italic')
        ax.text(50, -0.12, 'PLIC', transform=ax.get_xaxis_transform(), 
                ha='center', fontsize=9, style='italic')
        ax.text(100, -0.12, 'Precentral Gyrus', transform=ax.get_xaxis_transform(), 
                ha='center', fontsize=9, style='italic')
        ax.set_xlabel('Normalized Tract Position', fontsize=12)
    else:
        ax.set_xlabel('Normalized Tract Position (%)', fontsize=12)

    ax.set_ylabel(scalar_label, fontsize=12)
    ax.set_title(f'{scalar_label.split(" (")[0]} Profile - {subject_id}', fontsize=14, fontweight='bold')
    ax.legend(loc='best', fontsize=10)
    ax.grid(True, alpha=0.3)

    plt.tight_layout()
    plt.subplots_adjust(bottom=0.15)  # Extra space for anatomical labels

    # Save figure
    fig_path = output_dir / f"{subject_id}_tract_profile_{scalar}.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight')
    plt.close()

    print(f"✓ Tract profile saved: {fig_path}")
    return fig_path

plot_bilateral_comparison(comparison, output_dir, subject_id)

Create bar charts comparing left vs right CST metrics.

Parameters:

Name	Type	Description	Default
comparison	dict	Output from compare_bilateral_cst()	required
output_dir	str or Path	Output directory for saving figure	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_bilateral_comparison(
    comparison,
    output_dir,
    subject_id
):
    """
    Create bar charts comparing left vs right CST metrics.

    Parameters
    ----------
    comparison : dict
        Output from compare_bilateral_cst()
    output_dir : str or Path
        Output directory for saving figure
    subject_id : str
        Subject identifier

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    left = comparison['left']
    right = comparison['right']
    asym = comparison['asymmetry']

    # Prepare data for plotting
    metrics_to_plot = []

    # Morphology
    metrics_to_plot.append({
        'name': 'Streamline\nCount',
        'left': left['morphology']['n_streamlines'],
        'right': right['morphology']['n_streamlines'],
        'unit': '',
        'li': asym['streamline_count']['laterality_index']
    })

    metrics_to_plot.append({
        'name': 'Tract Volume\n(mm³)',
        'left': left['morphology']['tract_volume'],
        'right': right['morphology']['tract_volume'],
        'unit': 'mm³',
        'li': asym['volume']['laterality_index']
    })

    metrics_to_plot.append({
        'name': 'Mean Length\n(mm)',
        'left': left['morphology']['mean_length'],
        'right': right['morphology']['mean_length'],
        'unit': 'mm',
        'li': asym['mean_length']['laterality_index']
    })

    # Microstructure
    if 'fa' in left:
        metrics_to_plot.append({
            'name': 'Mean FA',
            'left': left['fa']['mean'],
            'right': right['fa']['mean'],
            'unit': '',
            'li': asym['fa']['laterality_index']
        })

    if 'md' in left:
        metrics_to_plot.append({
            'name': 'Mean MD\n(×10⁻³)',
            'left': left['md']['mean'] * 1000,  # Convert to 10^-3
            'right': right['md']['mean'] * 1000,
            'unit': '×10⁻³',
            'li': asym['md']['laterality_index']
        })

    # Create figure with subplots
    n_metrics = len(metrics_to_plot)
    fig, axes = plt.subplots(1, n_metrics, figsize=(4*n_metrics, 6))

    if n_metrics == 1:
        axes = [axes]

    # Plot each metric
    for ax, metric in zip(axes, metrics_to_plot):
        x_pos = [0, 1]
        values = [metric['left'], metric['right']]
        colors = ['#2196F3', '#F44336']  # Blue for left, red for right

        bars = ax.bar(x_pos, values, color=colors, alpha=0.7, edgecolor='black', linewidth=1.5)

        # Add value labels on bars
        for bar, val in zip(bars, values):
            height = bar.get_height()
            if 'FA' in metric['name']:
                label_text = f'{val:.3f}'
            elif 'MD' in metric['name']:
                label_text = f'{val:.2f}'
            else:
                label_text = f'{val:.0f}'

            ax.text(bar.get_x() + bar.get_width()/2., height,
                   label_text,
                   ha='center', va='bottom', fontsize=10, fontweight='bold')

        # Add laterality index
        li_text = f'LI = {metric["li"]:+.3f}'
        ax.text(0.5, ax.get_ylim()[1]*0.9, li_text,
               ha='center', va='top', fontsize=10,
               bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))

        ax.set_xticks(x_pos)
        ax.set_xticklabels(['Left', 'Right'], fontsize=11)
        ax.set_ylabel(metric['unit'], fontsize=10)
        ax.set_title(metric['name'], fontsize=12, fontweight='bold')
        ax.grid(True, axis='y', alpha=0.3)

    plt.suptitle(f'Bilateral CST Comparison - {subject_id}', fontsize=14, fontweight='bold', y=1.02)
    plt.tight_layout()

    # Save figure
    fig_path = output_dir / f"{subject_id}_bilateral_comparison.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight')
    plt.close()

    print(f"✓ Bilateral comparison saved: {fig_path}")
    return fig_path

plot_3d_streamlines(streamlines_left, streamlines_right, fa_map, affine, output_dir, subject_id)

Create 3D visualization of CST streamlines colored by FA.

Parameters:

Name	Type	Description	Default
streamlines_left	Streamlines	Left CST streamlines	required
streamlines_right	Streamlines	Right CST streamlines	required
fa_map	ndarray	3D FA map for coloring	required
affine	ndarray	4x4 affine transformation matrix	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_3d_streamlines(
    streamlines_left,
    streamlines_right,
    fa_map,
    affine,
    output_dir,
    subject_id
):
    """
    Create 3D visualization of CST streamlines colored by FA.

    Parameters
    ----------
    streamlines_left : Streamlines
        Left CST streamlines
    streamlines_right : Streamlines
        Right CST streamlines
    fa_map : ndarray
        3D FA map for coloring
    affine : ndarray
        4x4 affine transformation matrix
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    try:
        from dipy.viz import window, actor

        # Create renderer
        renderer = window.Renderer()
        renderer.SetBackground(1, 1, 1)  # White background

        # Add left CST (blue tones)
        if len(streamlines_left) > 0:
            renderer.add(actor.line(
                streamlines_left,
                colors=(0.2, 0.4, 0.8),
                linewidth=2,
                opacity=0.8
            ))

        # Add right CST (red tones)
        if len(streamlines_right) > 0:
            renderer.add(actor.line(
                streamlines_right,
                colors=(0.8, 0.2, 0.2),
                linewidth=2,
                opacity=0.8
            ))

        # Set camera
        renderer.set_camera(position=(200, 200, 200), focal_point=(0, 0, 0))

        # Save figure
        fig_path = output_dir / f"{subject_id}_3d_streamlines.png"
        window.record(renderer, out_path=str(fig_path), size=(800, 800))

        print(f"  ✓ 3D streamlines saved: {fig_path}")
        return fig_path

    except Exception as e:
        print(f"  ⚠️ 3D visualization failed: {e}")
        return None

create_summary_figure(comparison, streamlines_left, streamlines_right, fa_map, affine, output_dir, subject_id)

Create multi-panel summary figure with all key visualizations.

Parameters:

Name	Type	Description	Default
comparison	dict	Bilateral comparison metrics	required
streamlines_left	Streamlines	Left CST streamlines	required
streamlines_right	Streamlines	Right CST streamlines	required
fa_map	ndarray	3D FA map	required
affine	ndarray	4x4 affine transformation matrix	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
fig_path	Path	Path to saved summary figure

Source code in src/csttool/metrics/modules/visualizations.py

def create_summary_figure(
    comparison,
    streamlines_left,
    streamlines_right,
    fa_map,
    affine,
    output_dir,
    subject_id
):
    """
    Create multi-panel summary figure with all key visualizations.

    Parameters
    ----------
    comparison : dict
        Bilateral comparison metrics
    streamlines_left : Streamlines
        Left CST streamlines
    streamlines_right : Streamlines
        Right CST streamlines
    fa_map : ndarray
        3D FA map
    affine : ndarray
        4x4 affine transformation matrix
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier

    Returns
    -------
    fig_path : Path
        Path to saved summary figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    fig = plt.figure(figsize=(16, 10))

    # Create grid layout
    gs = fig.add_gridspec(3, 3, hspace=0.3, wspace=0.3)

    left = comparison['left']
    right = comparison['right']
    asym = comparison['asymmetry']

    # Panel 1: FA Tract Profile
    ax1 = fig.add_subplot(gs[0, :2])
    if 'fa' in left:
        left_profile = np.array(left['fa']['profile'])
        right_profile = np.array(right['fa']['profile'])
        x = np.linspace(0, 100, len(left_profile))

        ax1.plot(x, left_profile, 'b-', linewidth=2, label='Left', marker='o', markersize=3)
        ax1.plot(x, right_profile, 'r-', linewidth=2, label='Right', marker='s', markersize=3)
        ax1.set_xlabel('Normalized Position (%)')
        ax1.set_ylabel('Fractional Anisotropy')
        ax1.set_title('FA Tract Profile', fontweight='bold')
        ax1.legend()
        ax1.grid(True, alpha=0.3)

    # Panel 2: Key Metrics Table
    ax2 = fig.add_subplot(gs[0, 2])
    ax2.axis('off')

    table_data = [
        ['Metric', 'Left', 'Right', 'LI'],
        ['Streamlines', f"{left['morphology']['n_streamlines']}", 
         f"{right['morphology']['n_streamlines']}", 
         f"{asym['streamline_count']['laterality_index']:+.2f}"],
        ['Volume (mm³)', f"{left['morphology']['tract_volume']:.0f}", 
         f"{right['morphology']['tract_volume']:.0f}", 
         f"{asym['volume']['laterality_index']:+.2f}"]
    ]

    if 'fa' in left:
        table_data.append(['Mean FA', f"{left['fa']['mean']:.3f}", 
                          f"{right['fa']['mean']:.3f}", 
                          f"{asym['fa']['laterality_index']:+.2f}"])

    table = ax2.table(cellText=table_data, cellLoc='center', loc='center',
                     colWidths=[0.3, 0.2, 0.2, 0.2])
    table.auto_set_font_size(False)
    table.set_fontsize(9)
    table.scale(1, 2)

    # Style header row
    for i in range(4):
        table[(0, i)].set_facecolor('#4CAF50')
        table[(0, i)].set_text_props(weight='bold', color='white')

    ax2.set_title('Summary Metrics', fontweight='bold', pad=20)

    # Panel 3: Streamline Count Comparison
    ax3 = fig.add_subplot(gs[1, 0])
    ax3.bar(['Left', 'Right'], 
           [left['morphology']['n_streamlines'], right['morphology']['n_streamlines']],
           color=['#2196F3', '#F44336'], alpha=0.7, edgecolor='black')
    ax3.set_ylabel('Count')
    ax3.set_title('Streamline Count', fontweight='bold')
    ax3.grid(True, axis='y', alpha=0.3)

    # Panel 4: Volume Comparison
    ax4 = fig.add_subplot(gs[1, 1])
    ax4.bar(['Left', 'Right'],
           [left['morphology']['tract_volume'], right['morphology']['tract_volume']],
           color=['#2196F3', '#F44336'], alpha=0.7, edgecolor='black')
    ax4.set_ylabel('Volume (mm³)')
    ax4.set_title('Tract Volume', fontweight='bold')
    ax4.grid(True, axis='y', alpha=0.3)

    # Panel 5: FA Comparison
    ax5 = fig.add_subplot(gs[1, 2])
    if 'fa' in left:
        ax5.bar(['Left', 'Right'],
               [left['fa']['mean'], right['fa']['mean']],
               color=['#2196F3', '#F44336'], alpha=0.7, edgecolor='black')
        ax5.set_ylabel('Mean FA')
        ax5.set_title('Fractional Anisotropy', fontweight='bold')
        ax5.grid(True, axis='y', alpha=0.3)

    # Panel 6: Laterality Indices
    ax6 = fig.add_subplot(gs[2, :])

    metrics_names = ['Volume', 'Streamlines', 'Length']
    li_values = [
        asym['volume']['laterality_index'],
        asym['streamline_count']['laterality_index'],
        asym['mean_length']['laterality_index']
    ]

    if 'fa' in asym:
        metrics_names.append('FA')
        li_values.append(asym['fa']['laterality_index'])

    if 'md' in asym:
        metrics_names.append('MD')
        li_values.append(asym['md']['laterality_index'])

    colors = ['#2196F3' if li > 0 else '#F44336' for li in li_values]
    ax6.barh(metrics_names, li_values, color=colors, alpha=0.7, edgecolor='black')
    ax6.axvline(0, color='black', linewidth=1)
    ax6.axvline(-0.1, color='gray', linestyle='--', alpha=0.5)
    ax6.axvline(0.1, color='gray', linestyle='--', alpha=0.5)
    ax6.set_xlabel('Laterality Index (LI)')
    ax6.set_title('Asymmetry Analysis', fontweight='bold')
    ax6.grid(True, axis='x', alpha=0.3)
    ax6.text(0.12, 0.5, 'Left > Right', transform=ax6.transAxes, fontsize=9, color='blue')
    ax6.text(-0.12, 0.5, 'Right > Left', transform=ax6.transAxes, fontsize=9, color='red', ha='right')

    # Overall title
    plt.suptitle(f'CST Analysis Summary - {subject_id}', fontsize=16, fontweight='bold', y=0.98)

    # Save figure
    fig_path = output_dir / f"{subject_id}_summary.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight')
    plt.close()

    print(f"✓ Summary figure saved: {fig_path}")
    return fig_path

plot_asymmetry_radar(asymmetry, output_dir, subject_id)

Create radar plot showing asymmetry across multiple metrics.

Parameters:

Name	Type	Description	Default
asymmetry	dict	Asymmetry metrics from bilateral comparison	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_asymmetry_radar(asymmetry, output_dir, subject_id):
    """
    Create radar plot showing asymmetry across multiple metrics.

    Parameters
    ----------
    asymmetry : dict
        Asymmetry metrics from bilateral comparison
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    # Prepare data
    metrics = []
    values = []

    if 'volume' in asymmetry:
        metrics.append('Volume')
        values.append(abs(asymmetry['volume']['laterality_index']))

    if 'streamline_count' in asymmetry:
        metrics.append('Streamline\nCount')
        values.append(abs(asymmetry['streamline_count']['laterality_index']))

    if 'fa' in asymmetry:
        metrics.append('FA')
        values.append(abs(asymmetry['fa']['laterality_index']))

    if 'md' in asymmetry:
        metrics.append('MD')
        values.append(abs(asymmetry['md']['laterality_index']))

    if len(metrics) < 3:
        print("Not enough metrics for radar plot")
        return None

    # Create radar plot
    angles = np.linspace(0, 2 * np.pi, len(metrics), endpoint=False).tolist()
    values += values[:1]  # Close the plot
    angles += angles[:1]

    fig, ax = plt.subplots(figsize=(8, 8), subplot_kw=dict(projection='polar'))
    ax.plot(angles, values, 'o-', linewidth=2, color='#2196F3')
    ax.fill(angles, values, alpha=0.25, color='#2196F3')

    # Add threshold circle
    threshold = [0.1] * (len(metrics) + 1)
    ax.plot(angles, threshold, '--', linewidth=1, color='red', alpha=0.5, label='Asymmetry threshold (0.1)')

    ax.set_xticks(angles[:-1])
    ax.set_xticklabels(metrics, fontsize=11)
    ax.set_ylim(0, max(0.3, max(values)))
    ax.set_title(f'Asymmetry Profile - {subject_id}', fontsize=14, fontweight='bold', pad=20)
    ax.legend(loc='upper right', bbox_to_anchor=(1.3, 1.1))
    ax.grid(True)

    # Save
    fig_path = output_dir / f"{subject_id}_asymmetry_radar.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight')
    plt.close()

    print(f"✓ Asymmetry radar plot saved: {fig_path}")
    return fig_path

plot_stacked_profiles(left_metrics, right_metrics, output_dir, subject_id)

Create stacked FA, MD, RD, and AD profile plots for PDF report.

Creates a vertically stacked figure with 4 subplots (if data available): - FA profile - MD profile - RD profile - AD profile

All have shared x-axis (anatomical labels only on bottom).

Parameters:

Name	Type	Description	Default
left_metrics	dict	Left hemisphere metrics	required
right_metrics	dict	Right hemisphere metrics	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_stacked_profiles(
    left_metrics,
    right_metrics,
    output_dir,
    subject_id
):
    """
    Create stacked FA, MD, RD, and AD profile plots for PDF report.

    Creates a vertically stacked figure with 4 subplots (if data available):
    - FA profile
    - MD profile
    - RD profile
    - AD profile

    All have shared x-axis (anatomical labels only on bottom).

    Parameters
    ----------
    left_metrics : dict
        Left hemisphere metrics
    right_metrics : dict
        Right hemisphere metrics
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    # Define metrics to plot in order
    metrics_config = [
        {'key': 'fa', 'title': 'Fractional Anisotropy', 'ylabel': 'FA', 'ylim': (0, 0.8), 'scale': 1},
        {'key': 'md', 'title': 'Mean Diffusivity', 'ylabel': 'MD (×10⁻³)', 'ylim': (0.5, 1.2), 'scale': 1000},
        {'key': 'rd', 'title': 'Radial Diffusivity', 'ylabel': 'RD (×10⁻³)', 'ylim': (0.3, 1.0), 'scale': 1000},
        {'key': 'ad', 'title': 'Axial Diffusivity', 'ylabel': 'AD (×10⁻³)', 'ylim': (0.8, 1.8), 'scale': 1000}
    ]

    # Filter available metrics
    available_metrics = []
    for m in metrics_config:
        if m['key'] in left_metrics and m['key'] in right_metrics:
            available_metrics.append(m)

    if not available_metrics:
        print("Warning: No profiles available for stacking")
        return None

    n_plots = len(available_metrics)
    # Fixed height per plot
    fig, axes = plt.subplots(n_plots, 1, figsize=(6, 2.2 * n_plots), sharex=True)

    if n_plots == 1:
        axes = [axes]

    for i, (ax, m) in enumerate(zip(axes, available_metrics)):
        key = m['key']
        scale = m['scale']

        left_profile = np.array(left_metrics[key]['profile']) * scale
        right_profile = np.array(right_metrics[key]['profile']) * scale
        n_points = len(left_profile)
        x = np.linspace(0, 100, n_points)

        ax.plot(x, left_profile, 'b-', linewidth=2, label='Left CST', marker='o', markersize=3)
        ax.plot(x, right_profile, 'r-', linewidth=2, label='Right CST', marker='s', markersize=3)

        ax.set_ylabel(m['ylabel'], fontsize=10)
        # Auto stats for ylim might be better, but keeping fixed range as starting point logic
        # If 'ylim' is provided, use it, else auto
        if 'ylim' in m:
             # Basic check to see if data fits in default range, if not, auto-scale
             all_data = np.concatenate([left_profile, right_profile])
             if np.min(all_data) < m['ylim'][0] or np.max(all_data) > m['ylim'][1]:
                 # Auto scale with margin
                 margin = (np.max(all_data) - np.min(all_data)) * 0.1
                 ax.set_ylim(max(0, np.min(all_data) - margin), np.max(all_data) + margin)
             else:
                 ax.set_ylim(m['ylim'])

        ax.text(0.5, 0.9, m['title'], transform=ax.transAxes, fontsize=10, fontweight='bold', ha='center')

        ax.grid(True, alpha=0.3)

        # Only add legend to first plot
        if i == 0:
            ax.legend(loc='upper right', fontsize=9, framealpha=0.9)

    # Shared X-axis formatting (only on bottom plot)
    axes[-1].set_xticks([0, 50, 100])
    axes[-1].set_xticklabels(['0%', '50%', '100%'])

    # Add anatomical annotations
    # Create transform for the bottom axis
    trans = axes[-1].get_xaxis_transform()

    axes[-1].text(0, -0.3, 'Pontine\nLevel', transform=trans, 
            ha='center', fontsize=9, style='italic')
    axes[-1].text(50, -0.3, 'PLIC', transform=trans, 
            ha='center', fontsize=9, style='italic')
    axes[-1].text(100, -0.3, 'Precentral\nGyrus', transform=trans, 
            ha='center', fontsize=9, style='italic')

    plt.tight_layout()
    # Increase bottom margin to prevent x-axis overlap
    plt.subplots_adjust(hspace=0.2, bottom=0.15)  # Minimize vertical space between plots

    fig_path = output_dir / f"{subject_id}_stacked_profiles.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight')
    plt.close()

    print(f"✓ Stacked profiles saved: {fig_path}")
    return fig_path

plot_tractogram_qc_preview(streamlines_left, streamlines_right, background_image, affine, output_dir, subject_id, slice_type='axial', max_streamlines=500, set_title=True)

Create compact 3D tractogram QC preview for PDF report.

Renders left (blue) and right (red) CST streamlines overlaid on a brain slice at the internal capsule level.

Parameters:

Name	Type	Description	Default
streamlines_left	Streamlines	Left CST streamlines	required
streamlines_right	Streamlines	Right CST streamlines	required
background_image	ndarray	3D T1 or FA image for background	required
affine	ndarray	4x4 affine transformation matrix	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required
slice_type	str	'axial', 'sagittal', or 'coronal'	'axial'
max_streamlines	int	Maximum streamlines to render per hemisphere (for performance)	500

Returns:

Name	Type	Description
fig_path	Path	Path to saved figure

Source code in src/csttool/metrics/modules/visualizations.py

def plot_tractogram_qc_preview(
    streamlines_left,
    streamlines_right,
    background_image,
    affine,
    output_dir,
    subject_id,
    slice_type='axial',
    max_streamlines=500,
    set_title=True
):
    """
    Create compact 3D tractogram QC preview for PDF report.

    Renders left (blue) and right (red) CST streamlines overlaid
    on a brain slice at the internal capsule level.

    Parameters
    ----------
    streamlines_left : Streamlines
        Left CST streamlines
    streamlines_right : Streamlines
        Right CST streamlines
    background_image : ndarray
        3D T1 or FA image for background
    affine : ndarray
        4x4 affine transformation matrix
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier
    slice_type : str
        'axial', 'sagittal', or 'coronal'
    max_streamlines : int
        Maximum streamlines to render per hemisphere (for performance)

    Returns
    -------
    fig_path : Path
        Path to saved figure
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    fig, ax = plt.subplots(figsize=(4, 4))

    # Get slice at center (internal capsule level)
    shape = background_image.shape

    if slice_type == 'axial':
        slice_idx = shape[2] // 2 + 5  # Slightly above center for IC
        bg_slice = background_image[:, :, slice_idx].T
    elif slice_type == 'sagittal':
        slice_idx = shape[0] // 2
        bg_slice = background_image[slice_idx, :, :].T
    else:  # coronal
        slice_idx = shape[1] // 2
        bg_slice = background_image[:, slice_idx, :].T

    # Display background
    ax.imshow(bg_slice, cmap='gray', origin='lower', aspect='equal')

    # Adjust figure size to match data aspect ratio
    # logical_height / logical_width
    data_ratio = bg_slice.shape[0] / bg_slice.shape[1]

    # Set fixed width of 4 inches and adjust height
    fig.set_size_inches(4, 4 * data_ratio)

    # Project streamlines onto slice
    def project_streamlines(streamlines, color, alpha=0.6):
        """Project streamlines onto 2D and plot."""
        count = 0
        for sl in streamlines:
            if count >= max_streamlines:
                break
            # Convert to voxel coordinates
            inv_affine = np.linalg.inv(affine)
            voxel_coords = np.dot(sl, inv_affine[:3, :3].T) + inv_affine[:3, 3]

            if slice_type == 'axial':
                # Filter points near the slice
                near_slice = np.abs(voxel_coords[:, 2] - slice_idx) < 5
                if np.any(near_slice):
                    ax.plot(voxel_coords[near_slice, 0], 
                           voxel_coords[near_slice, 1], 
                           color=color, linewidth=0.5, alpha=alpha)
                    count += 1
            elif slice_type == 'sagittal':
                near_slice = np.abs(voxel_coords[:, 0] - slice_idx) < 5
                if np.any(near_slice):
                    ax.plot(voxel_coords[near_slice, 1], 
                           voxel_coords[near_slice, 2], 
                           color=color, linewidth=0.5, alpha=alpha)
                    count += 1
            else:  # coronal
                near_slice = np.abs(voxel_coords[:, 1] - slice_idx) < 5
                if np.any(near_slice):
                    ax.plot(voxel_coords[near_slice, 0], 
                           voxel_coords[near_slice, 2], 
                           color=color, linewidth=0.5, alpha=alpha)
                    count += 1

    # Plot streamlines
    if len(streamlines_left) > 0:
        project_streamlines(streamlines_left, '#2196F3')  # Blue
    if len(streamlines_right) > 0:
        project_streamlines(streamlines_right, '#F44336')  # Red

    # Add legend
    from matplotlib.lines import Line2D
    legend_elements = [
        Line2D([0], [0], color='#2196F3', linewidth=2, label='Left CST'),
        Line2D([0], [0], color='#F44336', linewidth=2, label='Right CST')
    ]
    ax.legend(handles=legend_elements, loc='upper right', fontsize=8)

    if set_title:
        ax.set_title(f'CST Tractogram ({slice_type.title()})', fontsize=10, fontweight='bold')

    ax.axis('off')

    plt.tight_layout()

    fig_path = output_dir / f"{subject_id}_tractogram_qc_{slice_type}.png"
    plt.savefig(fig_path, dpi=150, bbox_inches='tight', facecolor='white')
    plt.close()

    print(f"✓ Tractogram QC preview saved: {fig_path}")
    return fig_path

save_json_report(comparison, output_dir, subject_id, metadata=None)

Save comprehensive metrics report as JSON.

Parameters:

Name	Type	Description	Default
comparison	dict	Output from compare_bilateral_cst()	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required
metadata	dict	Additional metadata including: - acquisition: dict with protocol, b_values, n_directions, resolution - processing: dict with denoising_method, tracking_method, etc. - qc_thresholds: dict with fa_threshold, min_length, max_length	None

Returns:

Name	Type	Description
json_path	Path	Path to saved JSON file

Source code in src/csttool/metrics/modules/reports.py

def save_json_report(comparison, output_dir, subject_id, metadata=None):
    """
    Save comprehensive metrics report as JSON.

    Parameters
    ----------
    comparison : dict
        Output from compare_bilateral_cst()
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier
    metadata : dict, optional
        Additional metadata including:
        - acquisition: dict with protocol, b_values, n_directions, resolution
        - processing: dict with denoising_method, tracking_method, etc.
        - qc_thresholds: dict with fa_threshold, min_length, max_length

    Returns
    -------
    json_path : Path
        Path to saved JSON file
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    # Initialize metadata if not provided
    if metadata is None:
        metadata = {}

    # Build report with extended schema
    report = {
        'subject_id': subject_id,
        'processing_date': datetime.now().strftime("%Y-%m-%d %H:%M:%S"),
        'csttool_version': __version__,
        'acquisition': metadata.get('acquisition', {}),
        'processing': metadata.get('processing', {}),
        'qc_thresholds': metadata.get('qc_thresholds', {}),
        'metrics': comparison
    }

    # Save JSON
    json_path = output_dir / f"{subject_id}_bilateral_metrics.json"
    with open(json_path, 'w') as f:
        json.dump(report, f, indent=2)

    print(f"  ✓ JSON report saved: {json_path}")
    return json_path

save_csv_summary(comparison, output_dir, subject_id)

Save metrics summary as CSV table.

Creates a flat CSV file with key metrics suitable for group analysis.

Parameters:

Name	Type	Description	Default
comparison	dict	Output from compare_bilateral_cst()	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required

Returns:

Name	Type	Description
csv_path	Path	Path to saved CSV file

Source code in src/csttool/metrics/modules/reports.py

def save_csv_summary(comparison, output_dir, subject_id):
    """
    Save metrics summary as CSV table.

    Creates a flat CSV file with key metrics suitable for group analysis.

    Parameters
    ----------
    comparison : dict
        Output from compare_bilateral_cst()
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier

    Returns
    -------
    csv_path : Path
        Path to saved CSV file
    """

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    left = comparison['left']
    right = comparison['right']
    asym = comparison['asymmetry']

    # Prepare data row
    data = {
        'subject_id': subject_id,
        'processing_date': datetime.now().strftime("%Y-%m-%d"),

        # Left morphology
        'left_n_streamlines': left['morphology']['n_streamlines'],
        'left_mean_length_mm': left['morphology']['mean_length'],
        'left_tract_volume_mm3': left['morphology']['tract_volume'],

        # Right morphology
        'right_n_streamlines': right['morphology']['n_streamlines'],
        'right_mean_length_mm': right['morphology']['mean_length'],
        'right_tract_volume_mm3': right['morphology']['tract_volume'],

        # Asymmetry
        'volume_laterality_index': asym['volume']['laterality_index'],
        'streamline_count_laterality_index': asym['streamline_count']['laterality_index'],
    }

    # Add FA if available
    if 'fa' in left:
        data.update({
            'left_fa_mean': left['fa']['mean'],
            'left_fa_std': left['fa']['std'],
            'right_fa_mean': right['fa']['mean'],
            'right_fa_std': right['fa']['std'],
            'fa_laterality_index': asym['fa']['laterality_index'],
        })

    # Add MD if available
    if 'md' in left:
        data.update({
            'left_md_mean': left['md']['mean'],
            'left_md_std': left['md']['std'],
            'right_md_mean': right['md']['mean'],
            'right_md_std': right['md']['std'],
            'md_laterality_index': asym['md']['laterality_index'],
        })

    # Add RD if available
    if 'rd' in left:
        data.update({
            'left_rd_mean': left['rd']['mean'],
            'left_rd_std': left['rd']['std'],
            'right_rd_mean': right['rd']['mean'],
            'right_rd_std': right['rd']['std'],
            'rd_laterality_index': asym['rd']['laterality_index'],
        })

    # Add AD if available
    if 'ad' in left:
        data.update({
            'left_ad_mean': left['ad']['mean'],
            'left_ad_std': left['ad']['std'],
            'right_ad_mean': right['ad']['mean'],
            'right_ad_std': right['ad']['std'],
            'ad_laterality_index': asym['ad']['laterality_index'],
        })

    # Add localized metrics (pontine, plic, precentral) for each scalar
    regions = ['pontine', 'plic', 'precentral']
    scalars = ['fa', 'md', 'rd', 'ad']

    for scalar in scalars:
        if scalar in left and 'pontine' in left[scalar]:
            for region in regions:
                data[f'left_{scalar}_{region}'] = left[scalar].get(region, 0.0)
                data[f'right_{scalar}_{region}'] = right[scalar].get(region, 0.0)
                # Add LI for localized metric
                li_key = f'{scalar}_{region}'
                if li_key in asym:
                    data[f'{scalar}_{region}_laterality_index'] = asym[li_key]['laterality_index']

    # Save CSV
    csv_path = output_dir / f"{subject_id}_metrics_summary.csv"
    with open(csv_path, 'w', newline='') as f:
        writer = csv.DictWriter(f, fieldnames=data.keys())
        writer.writeheader()
        writer.writerow(data)

    print(f"  ✓ CSV summary saved: {csv_path}")
    return csv_path

save_html_report(comparison, visualization_paths, output_dir, subject_id, version=None, space='Native Space', metadata=None)

Generate HTML report using Jinja2 template.

Parameters:

Name	Type	Description	Default
comparison	dict	Output from compare_bilateral_cst()	required
visualization_paths	dict	Paths to generated visualizations	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required
version	str	csttool version (defaults to package version)	None
space	str	Space declaration (e.g., "Native Space")	'Native Space'
metadata	dict	Acquisition and processing metadata	None

Returns:

Name	Type	Description
html_path	Path	Path to saved HTML file

Source code in src/csttool/metrics/modules/reports.py

def save_html_report(
    comparison,
    visualization_paths,
    output_dir,
    subject_id,
    version=None,
    space="Native Space",
    metadata=None
):
    """
    Generate HTML report using Jinja2 template.

    Parameters
    ----------
    comparison : dict
        Output from compare_bilateral_cst()
    visualization_paths : dict
        Paths to generated visualizations
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier
    version : str, optional
        csttool version (defaults to package version)
    space : str
        Space declaration (e.g., "Native Space")
    metadata : dict, optional
        Acquisition and processing metadata

    Returns
    -------
    html_path : Path
        Path to saved HTML file
    """
    if version is None:
        version = __version__

    if metadata is None:
        metadata = {}

    output_dir = Path(output_dir)
    output_dir.mkdir(parents=True, exist_ok=True)

    # Extract comparison data
    left = comparison['left']
    right = comparison['right']
    asym = comparison['asymmetry']

    # Helper functions for formatting
    def fmt_mean_sd(mean, std, is_diffusivity=False):
        """Format mean ± SD string."""
        if is_diffusivity:
            return f"{mean*1000:.2f} ± {std*1000:.2f}"
        return f"{mean:.3f} ± {std:.3f}"

    def fmt_med_range(median, min_val, max_val, is_diffusivity=False):
        """Format median (min-max) string."""
        if is_diffusivity:
            return f"{median*1000:.2f} ({min_val*1000:.2f}-{max_val*1000:.2f})"
        return f"{median:.3f} ({min_val:.3f}-{max_val:.3f})"

    # Build metrics list for template (6-column format)
    metrics = []

    # Streamlines (no SD/range available, use simple format)
    metrics.append({
        "label": "Streamlines",
        "left_mean_sd": str(left['morphology']['n_streamlines']),
        "left_med_range": "-",
        "right_mean_sd": str(right['morphology']['n_streamlines']),
        "right_med_range": "-",
        "li": asym['streamline_count']['laterality_index']
    })

    # Volume (convert mm³ to cm³, no SD/range available)
    metrics.append({
        "label": "Volume (cm³)",
        "left_mean_sd": f"{left['morphology']['tract_volume'] / 1000.0:.2f}",
        "left_med_range": "-",
        "right_mean_sd": f"{right['morphology']['tract_volume'] / 1000.0:.2f}",
        "right_med_range": "-",
        "li": asym['volume']['laterality_index']
    })

    # Length
    lm = left['morphology']
    rm = right['morphology']
    metrics.append({
        "label": "Length (mm)",
        "left_mean_sd": f"{lm['mean_length']:.1f} ± {lm['std_length']:.1f}",
        "left_med_range": f"({lm['min_length']:.1f}-{lm['max_length']:.1f})",
        "right_mean_sd": f"{rm['mean_length']:.1f} ± {rm['std_length']:.1f}",
        "right_med_range": f"({rm['min_length']:.1f}-{rm['max_length']:.1f})",
        "li": asym['mean_length']['laterality_index']
    })

    # FA
    if 'fa' in left:
        # Backward compatibility: use mean as fallback for median if not present
        left_fa_median = left['fa'].get('median', left['fa']['mean'])
        left_fa_min = left['fa'].get('min', max(0.0, left['fa']['mean'] - 3*left['fa']['std']))
        left_fa_max = left['fa'].get('max', min(1.0, left['fa']['mean'] + 3*left['fa']['std']))
        right_fa_median = right['fa'].get('median', right['fa']['mean'])
        right_fa_min = right['fa'].get('min', max(0.0, right['fa']['mean'] - 3*right['fa']['std']))
        right_fa_max = right['fa'].get('max', min(1.0, right['fa']['mean'] + 3*right['fa']['std']))

        metrics.append({
            "label": "FA",
            "left_mean_sd": fmt_mean_sd(left['fa']['mean'], left['fa']['std']),
            "left_med_range": fmt_med_range(left_fa_median, left_fa_min, left_fa_max),
            "right_mean_sd": fmt_mean_sd(right['fa']['mean'], right['fa']['std']),
            "right_med_range": fmt_med_range(right_fa_median, right_fa_min, right_fa_max),
            "li": asym['fa']['laterality_index']
        })

    # MD (×10⁻³ mm²/s)
    if 'md' in left:
        # Backward compatibility: use mean as fallback for median if not present
        left_md_median = left['md'].get('median', left['md']['mean'])
        left_md_min = left['md'].get('min', max(0.0, left['md']['mean'] - 3*left['md']['std']))
        left_md_max = left['md'].get('max', left['md']['mean'] + 3*left['md']['std'])
        right_md_median = right['md'].get('median', right['md']['mean'])
        right_md_min = right['md'].get('min', max(0.0, right['md']['mean'] - 3*right['md']['std']))
        right_md_max = right['md'].get('max', right['md']['mean'] + 3*right['md']['std'])

        metrics.append({
            "label": "MD (×10⁻³)",
            "left_mean_sd": fmt_mean_sd(left['md']['mean'], left['md']['std'], is_diffusivity=True),
            "left_med_range": fmt_med_range(left_md_median, left_md_min, left_md_max, is_diffusivity=True),
            "right_mean_sd": fmt_mean_sd(right['md']['mean'], right['md']['std'], is_diffusivity=True),
            "right_med_range": fmt_med_range(right_md_median, right_md_min, right_md_max, is_diffusivity=True),
            "li": asym['md']['laterality_index']
        })

    # RD (×10⁻³ mm²/s)
    if 'rd' in left:
        # Backward compatibility: use mean as fallback for median if not present
        left_rd_median = left['rd'].get('median', left['rd']['mean'])
        left_rd_min = left['rd'].get('min', max(0.0, left['rd']['mean'] - 3*left['rd']['std']))
        left_rd_max = left['rd'].get('max', left['rd']['mean'] + 3*left['rd']['std'])
        right_rd_median = right['rd'].get('median', right['rd']['mean'])
        right_rd_min = right['rd'].get('min', max(0.0, right['rd']['mean'] - 3*right['rd']['std']))
        right_rd_max = right['rd'].get('max', right['rd']['mean'] + 3*right['rd']['std'])

        metrics.append({
            "label": "RD (×10⁻³)",
            "left_mean_sd": fmt_mean_sd(left['rd']['mean'], left['rd']['std'], is_diffusivity=True),
            "left_med_range": fmt_med_range(left_rd_median, left_rd_min, left_rd_max, is_diffusivity=True),
            "right_mean_sd": fmt_mean_sd(right['rd']['mean'], right['rd']['std'], is_diffusivity=True),
            "right_med_range": fmt_med_range(right_rd_median, right_rd_min, right_rd_max, is_diffusivity=True),
            "li": asym['rd']['laterality_index']
        })

    # AD (×10⁻³ mm²/s)
    if 'ad' in left:
        # Backward compatibility: use mean as fallback for median if not present
        left_ad_median = left['ad'].get('median', left['ad']['mean'])
        left_ad_min = left['ad'].get('min', max(0.0, left['ad']['mean'] - 3*left['ad']['std']))
        left_ad_max = left['ad'].get('max', left['ad']['mean'] + 3*left['ad']['std'])
        right_ad_median = right['ad'].get('median', right['ad']['mean'])
        right_ad_min = right['ad'].get('min', max(0.0, right['ad']['mean'] - 3*right['ad']['std']))
        right_ad_max = right['ad'].get('max', right['ad']['mean'] + 3*right['ad']['std'])

        metrics.append({
            "label": "AD (×10⁻³)",
            "left_mean_sd": fmt_mean_sd(left['ad']['mean'], left['ad']['std'], is_diffusivity=True),
            "left_med_range": fmt_med_range(left_ad_median, left_ad_min, left_ad_max, is_diffusivity=True),
            "right_mean_sd": fmt_mean_sd(right['ad']['mean'], right['ad']['std'], is_diffusivity=True),
            "right_med_range": fmt_med_range(right_ad_median, right_ad_min, right_ad_max, is_diffusivity=True),
            "li": asym['ad']['laterality_index']
        })

    # Build localized metrics for template
    localized_metrics = []
    regions = [('Pontine', 'pontine'), ('PLIC', 'plic'), ('Precentral', 'precentral')]

    def fmt_localized(scalar, region):
        """Format L / R / LI string for localized metric."""
        if scalar not in left or region not in left[scalar]:
            return "-"
        l_val = left[scalar][region]
        r_val = right[scalar][region]
        li_key = f'{scalar}_{region}'
        li_val = asym.get(li_key, {}).get('laterality_index', 0.0)
        if scalar == 'fa':
            return f"{l_val:.3f} / {r_val:.3f} / {li_val:+.3f}"
        else:
            # Diffusivity values (×10⁻³)
            return f"{l_val*1000:.2f} / {r_val*1000:.2f} / {li_val:+.3f}"

    for region_name, region_key in regions:
        localized_metrics.append({
            'name': region_name,
            'fa': fmt_localized('fa', region_key),
            'md': fmt_localized('md', region_key),
            'rd': fmt_localized('rd', region_key),
            'ad': fmt_localized('ad', region_key)
        })

    # Build visualization data for template (coronal only)
    viz = {
        "stacked_profiles": _embed_image(visualization_paths.get('stacked_profiles')),
        "profiles": [
            {"title": "Fractional Anisotropy", "data": _embed_image(visualization_paths.get('profile_fa'))},
            {"title": "Mean Diffusivity", "data": _embed_image(visualization_paths.get('profile_md'))},
            {"title": "Radial Diffusivity", "data": _embed_image(visualization_paths.get('profile_rd'))},
            {"title": "Axial Diffusivity", "data": _embed_image(visualization_paths.get('profile_ad'))},
        ],
        "tractogram_coronal": _embed_image(visualization_paths.get('tractogram_qc_coronal'))
    }

    # Get acquisition/processing metadata with defaults
    acquisition = metadata.get('acquisition', {})
    processing = metadata.get('processing', {})

    # Format whole brain streamline count
    whole_brain = processing.get('whole_brain_streamlines', 'N/A')
    if whole_brain != 'N/A' and isinstance(whole_brain, (int, float)):
        processing = {**processing, 'whole_brain': f'{int(whole_brain):,} streamlines'}
    else:
        processing = {**processing, 'whole_brain': whole_brain}

    # Load template and CSS
    template = _jinja_env.get_template("report.html.j2")
    css = (_TEMPLATE_DIR / "report.css").read_text()

    # Render template with context
    html_content = template.render(
        subject_id=subject_id,
        date=datetime.now().strftime("%Y-%m-%d"),
        version=version,
        space=space,
        css=css,
        metrics=metrics,
        localized_metrics=localized_metrics,
        viz=viz,
        acquisition=acquisition,
        processing=processing,
    )

    html_path = output_dir / f"{subject_id}_report.html"
    html_path.write_text(html_content, encoding="utf-8")

    print(f"  ✓ HTML report saved: {html_path}")
    return html_path

save_pdf_report(comparison, visualization_paths, output_dir, subject_id, version=None, space='Native Space', html_path=None)

Generate PDF report by converting HTML report.

Source code in src/csttool/metrics/modules/reports.py

def save_pdf_report(
    comparison,
    visualization_paths,
    output_dir,
    subject_id,
    version=None,
    space="Native Space",
    html_path=None
):
    """
    Generate PDF report by converting HTML report.
    """
    output_dir = Path(output_dir)
    pdf_path = output_dir / f"{subject_id}_report.pdf"

    if html_path is None:
        # Fallback: try to find the standard HTML report
        possible_html = output_dir / f"{subject_id}_report.html"
        if possible_html.exists():
            html_path = possible_html
        else:
            print("⚠️ HTML report not found for PDF conversion. Generating temporary HTML...")
            # We would need to call save_html_report here, but let's rely on caller
            return None

    return html_to_pdf(html_path, pdf_path)

generate_complete_report(comparison, streamlines_left, streamlines_right, fa_map, affine, output_dir, subject_id, background_image=None, version=None, space='Native Space', metadata=None)

Generate all report formats: JSON, CSV, and PDF with visualizations.

Parameters:

Name	Type	Description	Default
comparison	dict	Bilateral comparison metrics	required
streamlines_left	Streamlines	Left CST streamlines	required
streamlines_right	Streamlines	Right CST streamlines	required
fa_map	ndarray	3D FA map	required
affine	ndarray	4x4 affine transformation	required
output_dir	str or Path	Output directory	required
subject_id	str	Subject identifier	required
background_image	ndarray	3D T1 or FA image for tractogram QC background (defaults to fa_map)	None
version	str	csttool version string	None

Returns:

Name	Type	Description
report_paths	dict	Dictionary of all generated report file paths

Source code in src/csttool/metrics/modules/reports.py

def generate_complete_report(
    comparison,
    streamlines_left,
    streamlines_right,
    fa_map,
    affine,
    output_dir,
    subject_id,
    background_image=None,
    version=None,
    space="Native Space",
    metadata=None
):
    """
    Generate all report formats: JSON, CSV, and PDF with visualizations.

    Parameters
    ----------
    comparison : dict
        Bilateral comparison metrics
    streamlines_left : Streamlines
        Left CST streamlines
    streamlines_right : Streamlines
        Right CST streamlines
    fa_map : ndarray
        3D FA map
    affine : ndarray
        4x4 affine transformation
    output_dir : str or Path
        Output directory
    subject_id : str
        Subject identifier
    background_image : ndarray, optional
        3D T1 or FA image for tractogram QC background (defaults to fa_map)
    version : str
        csttool version string

    Returns
    -------
    report_paths : dict
        Dictionary of all generated report file paths
    """

    from .visualizations import (
        plot_tract_profiles,
        plot_bilateral_comparison,
        create_summary_figure,
        plot_stacked_profiles,
        plot_tractogram_qc_preview
    )

    # Use package version if not specified
    if version is None:
        version = __version__

    output_dir = Path(output_dir)
    viz_dir = output_dir / "visualizations"
    viz_dir.mkdir(parents=True, exist_ok=True)

    print("\nGenerating complete report package...")

    # Use FA map as background if no background image provided
    if background_image is None:
        background_image = fa_map

    # Generate visualizations for PDF (new single-page layout)
    pdf_viz_paths = {}

    # Stacked FA/MD profiles for PDF
    pdf_viz_paths['stacked_profiles'] = plot_stacked_profiles(
        comparison['left'],
        comparison['right'],
        viz_dir,
        subject_id
    )

    # Tractogram QC preview for PDF (Coronal only)
    pdf_viz_paths['tractogram_qc_coronal'] = plot_tractogram_qc_preview(
        streamlines_left,
        streamlines_right,
        background_image,
        affine,
        viz_dir,
        subject_id,
        slice_type='coronal',
        set_title=False
    )

    # Also generate individual plots for detailed analysis
    viz_paths = {}

    if 'fa' in comparison['left']:
        viz_paths['tract_profiles_fa'] = plot_tract_profiles(
            comparison['left'],
            comparison['right'],
            viz_dir,
            subject_id,
            scalar='fa'
        )

    if 'md' in comparison['left']:
        viz_paths['tract_profiles_md'] = plot_tract_profiles(
            comparison['left'],
            comparison['right'],
            viz_dir,
            subject_id,
            scalar='md'
        )

    viz_paths['bilateral_comparison'] = plot_bilateral_comparison(
        comparison,
        viz_dir,
        subject_id
    )

    viz_paths['summary'] = create_summary_figure(
        comparison,
        streamlines_left,
        streamlines_right,
        fa_map,
        affine,
        viz_dir,
        subject_id
    )

    # Merge all visualization paths
    viz_paths.update(pdf_viz_paths)

    # Generate reports
    # 1. HTML Report (now critical as PDF is derived from it)
    html_path = save_html_report(comparison, pdf_viz_paths, output_dir, subject_id, version, space, metadata)

    # 2. PDF Report (from HTML)
    pdf_path = save_pdf_report(comparison, pdf_viz_paths, output_dir, subject_id, version, space, html_path)

    report_paths = {
        'json': save_json_report(comparison, output_dir, subject_id, metadata=metadata),
        'csv': save_csv_summary(comparison, output_dir, subject_id),
        'html': html_path,
        'pdf': pdf_path,
        'visualizations': viz_paths
    }

    print("\n  ✓ Complete report package generated")
    return report_paths