torchio.visualization.make_video() - Code Metrics - TorchIO-project/torchio - Measure and Improve Code Quality continuously with Scrutinizer

torchio.visualization.make_video() F
last analyzed 2025-08-12 20:07 UTC

↳ Parent: torchio.visualization

Complexity

Conditions

Size

Total Lines	117
Code Lines	88

Duplication

Lines	0
Ratio	0 %

Importance

Changes

Metric	Value
cc	21
eloc	88
nop	6
dl	0
loc	117
rs	0
c	0
b	0
f	0

How to fix Long Method Complexity

from __future__ import annotations

import warnings
from itertools import cycle
from pathlib import Path
from typing import TYPE_CHECKING
from typing import Any

import numpy as np
import torch
from einops import rearrange

from .data.image import Image
from .data.image import LabelMap
from .data.image import ScalarImage
from .data.subject import Subject
from .external.imports import get_ffmpeg
from .transforms.preprocessing.intensity.rescale import RescaleIntensity
from .transforms.preprocessing.intensity.to import To
from .transforms.preprocessing.spatial.ensure_shape_multiple import EnsureShapeMultiple
from .transforms.preprocessing.spatial.resample import Resample
from .transforms.preprocessing.spatial.to_canonical import ToCanonical
from .transforms.preprocessing.spatial.to_orientation import ToOrientation
from .types import TypePath

if TYPE_CHECKING:
    from matplotlib.colors import BoundaryNorm
    from matplotlib.colors import ListedColormap
    from matplotlib.figure import Figure


def import_mpl_plt():
    try:
        import matplotlib as mpl
        import matplotlib.pyplot as plt
    except ImportError as e:
        raise ImportError('Install matplotlib for plotting support') from e
    return mpl, plt


def rotate(image: np.ndarray, *, radiological: bool = True, n: int = -1) -> np.ndarray:
    # Rotate for visualization purposes
    image = np.rot90(image, n, axes=(0, 1))
    if radiological:
        image = np.fliplr(image)
    return image


def _create_categorical_colormap(
    data: torch.Tensor,
    cmap_name: str = 'glasbey_category10',
) -> tuple[ListedColormap, BoundaryNorm]:
    num_classes = int(data.max())
    mpl, _ = import_mpl_plt()

    colors = [
        (0, 0, 0),  # black for background
        (1, 1, 1),  # white for class 1
    ]
    if num_classes > 1:
        from .external.imports import get_colorcet

        colorcet = get_colorcet()
        cmap = getattr(colorcet.cm, cmap_name)
        color_cycle = cycle(cmap.colors)
        distinct_colors = [next(color_cycle) for _ in range(num_classes - 1)]
        colors.extend(distinct_colors)
    boundaries = np.arange(-0.5, num_classes + 1.5, 1)
    colormap = mpl.colors.ListedColormap(colors)
    boundary_norm = mpl.colors.BoundaryNorm(boundaries, ncolors=colormap.N)
    return colormap, boundary_norm


def plot_volume(
    image: Image,
    radiological=True,
    channel=None,
    axes=None,
    cmap=None,
    output_path=None,
    show=True,
    xlabels=True,
    percentiles: tuple[float, float] = (0, 100),
    figsize=None,
    title=None,
    reorient=True,
    indices=None,
    rgb=True,
    savefig_kwargs: dict[str, Any] | None = None,
    **imshow_kwargs,
) -> Figure | None:
    _, plt = import_mpl_plt()
    fig: Figure | None = None

    if axes is None:
        fig, axes = plt.subplots(1, 3, figsize=figsize)

    if reorient:
        image = ToCanonical()(image)  # type: ignore[assignment]

    is_label = isinstance(image, LabelMap)
    if is_label:  # probabilistic label map
        data = image.data[np.newaxis, -1]
    elif rgb and image.num_channels == 3:
        data = image.data  # keep image as it is
    elif channel is None:
        data = image.data[0:1]  # just use the first channel
    else:
        data = image.data[np.newaxis, channel]
    data = rearrange(data, 'c x y z -> x y z c')
    data_numpy: np.ndarray = data.cpu().numpy()

    if indices is None:
        indices = np.array(data_numpy.shape[:3]) // 2
    i, j, k = indices
    slice_x = rotate(data_numpy[i, :, :], radiological=radiological)
    slice_y = rotate(data_numpy[:, j, :], radiological=radiological)
    slice_z = rotate(data_numpy[:, :, k], radiological=radiological)

    if isinstance(cmap, dict):
        slices = slice_x, slice_y, slice_z
        slice_x, slice_y, slice_z = color_labels(slices, cmap)
    else:
        boundary_norm = None
        if cmap is None:
            if is_label:
                cmap, boundary_norm = _create_categorical_colormap(data)
            else:
                cmap = 'gray'
        imshow_kwargs['cmap'] = cmap
        imshow_kwargs['norm'] = boundary_norm

    if is_label:
        imshow_kwargs['interpolation'] = 'none'
    else:
        if 'interpolation' not in imshow_kwargs:
            imshow_kwargs['interpolation'] = 'bicubic'

    imshow_kwargs['origin'] = 'lower'

    if not is_label:
        displayed_data = np.concatenate(
            [
                slice_x.flatten(),
                slice_y.flatten(),
                slice_z.flatten(),
            ]
        )
        p1, p2 = np.percentile(displayed_data, percentiles)
        imshow_kwargs['vmin'] = p1
        imshow_kwargs['vmax'] = p2

    spacing_r, spacing_a, spacing_s = image.spacing
    sag_axis, cor_axis, axi_axis = axes
    slices_dict = {
        'Sagittal': {
            'aspect': spacing_s / spacing_a,
            'slice': slice_x,
            'xlabel': 'A',
            'ylabel': 'S',
            'axis': sag_axis,
        },
        'Coronal': {
            'aspect': spacing_s / spacing_r,
            'slice': slice_y,
            'xlabel': 'R',
            'ylabel': 'S',
            'axis': cor_axis,
        },
        'Axial': {
            'aspect': spacing_a / spacing_r,
            'slice': slice_z,
            'xlabel': 'R',
            'ylabel': 'A',
            'axis': axi_axis,
        },
    }

    for axis_title, info in slices_dict.items():
        axis = info['axis']
        axis.imshow(info['slice'], aspect=info['aspect'], **imshow_kwargs)
        if xlabels:
            axis.set_xlabel(info['xlabel'])
        axis.set_ylabel(info['ylabel'])
        axis.invert_xaxis()
        axis.set_title(axis_title)

    plt.tight_layout()
    if title is not None:
        plt.suptitle(title)

    if output_path is not None and fig is not None:
        if savefig_kwargs is None:
            savefig_kwargs = {}
        fig.savefig(output_path, **savefig_kwargs)
    if show:
        plt.show()
    return fig


def plot_subject(
    subject: Subject,
    cmap_dict=None,
    show=True,
    output_path=None,
    figsize=None,
    clear_axes=True,
    **plot_volume_kwargs,
):
    _, plt = import_mpl_plt()
    num_images = len(subject)
    many_images = num_images > 2
    subplots_kwargs = {'figsize': figsize}
    try:
        if clear_axes:
            subject.check_consistent_spatial_shape()
            subplots_kwargs['sharex'] = 'row' if many_images else 'col'
            subplots_kwargs['sharey'] = 'row' if many_images else 'col'
    except RuntimeError:  # different shapes in subject
        pass
    args = (3, num_images) if many_images else (num_images, 3)
    fig, axes = plt.subplots(*args, **subplots_kwargs)
    # The array of axes must be 2D so that it can be indexed correctly within
    # the plot_volume() function
    axes = axes.T if many_images else axes.reshape(-1, 3)
    iterable = enumerate(subject.get_images_dict(intensity_only=False).items())
    axes_names = 'sagittal', 'coronal', 'axial'
    for image_index, (name, image) in iterable:
        image_axes = axes[image_index]
        cmap = None
        if cmap_dict is not None and name in cmap_dict:
            cmap = cmap_dict[name]
        last_row = image_index == len(axes) - 1
        plot_volume(
            image,
            axes=image_axes,
            show=False,
            cmap=cmap,
            xlabels=last_row,
            **plot_volume_kwargs,
        )
        for axis, axis_name in zip(image_axes, axes_names):
            axis.set_title(f'{name} ({axis_name})')
    plt.tight_layout()
    if output_path is not None:
        fig.savefig(output_path)
    if show:
        plt.show()


def get_num_bins(x: np.ndarray) -> int:
    """Get the optimal number of bins for a histogram.

    This method uses the Freedman–Diaconis rule to compute the histogram that
    minimizes "the integral of the squared difference between the histogram
    (i.e., relative frequency density) and the density of the theoretical
    probability distribution" (`Wikipedia <https://en.wikipedia.org/wiki/Freedman%E2%80%93Diaconis_rule>`_).

    Args:
        x: Input values.
    """
    # Freedman–Diaconis number of bins
    q25, q75 = np.percentile(x, [25, 75])
    bin_width = 2 * (q75 - q25) * len(x) ** (-1 / 3)
    bins = round((x.max() - x.min()) / bin_width)
    return bins


def plot_histogram(x: np.ndarray, show=True, **kwargs) -> None:
    _, plt = import_mpl_plt()
    plt.hist(x, bins=get_num_bins(x), **kwargs)
    plt.xlabel('Intensity')
    density = kwargs.pop('density', False)
    ylabel = 'Density' if density else 'Frequency'
    plt.ylabel(ylabel)
    if show:
        plt.show()


def color_labels(arrays, cmap_dict):
    results = []
    for slice_array in arrays:
        si, sj, _ = slice_array.shape
        rgb = np.zeros((si, sj, 3), dtype=np.uint8)
        for label, color in cmap_dict.items():
            if isinstance(color, str):
                mpl, _ = import_mpl_plt()
                color = mpl.colors.to_rgb(color)
                color = [255 * n for n in color]
            rgb[slice_array[..., 0] == label] = color
        results.append(rgb)
    return results


def make_gif(
    tensor: torch.Tensor,
    axis: int,
    duration: float,  # of full gif
    output_path: TypePath,
    loop: int = 0,
    optimize: bool = True,
    rescale: bool = True,
    reverse: bool = False,
) -> None:
    try:
        from PIL import Image as ImagePIL
    except ModuleNotFoundError as e:
        message = 'Please install Pillow to use Image.to_gif(): pip install Pillow'
        raise RuntimeError(message) from e
    transform = RescaleIntensity((0, 255))
    tensor = transform(tensor) if rescale else tensor  # type: ignore[assignment]
    single_channel = len(tensor) == 1

    # Move channels dimension to the end and bring selected axis to 0
    axes = np.roll(range(1, 4), -axis)
    tensor = tensor.permute(*axes, 0)

    if single_channel:
        mode = 'P'
        tensor = tensor[..., 0]
    else:
        mode = 'RGB'
    array = tensor.byte().numpy()
    n = 2 if axis == 1 else 1
    images = [ImagePIL.fromarray(rotate(i, n=n)).convert(mode) for i in array]
    num_images = len(images)
    images = list(reversed(images)) if reverse else images
    frame_duration_ms = duration / num_images * 1000
    if frame_duration_ms < 10:
        fps = round(1000 / frame_duration_ms)
        frame_duration_ms = 10
        new_duration = frame_duration_ms * num_images / 1000
        message = (
            'The computed frame rate from the given duration is too high'
            f' ({fps} fps). The highest possible frame rate in the GIF'
            ' file format specification is 100 fps. The duration has been set'
            f' to {new_duration:.1f} seconds, instead of {duration:.1f}'
        )
        warnings.warn(message, RuntimeWarning, stacklevel=2)
    images[0].save(
        output_path,
        save_all=True,
        append_images=images[1:],
        optimize=optimize,
        duration=frame_duration_ms,
        loop=loop,
    )


def make_video(
    image: ScalarImage,
    output_path: TypePath,
    seconds: float | None = None,
    frame_rate: float | None = None,
    direction: str = 'I',
    verbosity: str = 'error',
) -> None:
    ffmpeg = get_ffmpeg()

    if seconds is None and frame_rate is None:
        message = 'Either seconds or frame_rate must be provided.'
        raise ValueError(message)
    if seconds is not None and frame_rate is not None:
        message = 'Provide either seconds or frame_rate, not both.'
        raise ValueError(message)
    if image.num_channels > 1:
        message = 'Only single-channel tensors are supported for video output for now.'
        raise ValueError(message)
    tmin, tmax = image.data.min(), image.data.max()
    if tmin < 0 or tmax > 255:
        message = (
            'The tensor must be in the range [0, 256) for video output.'
            ' The image data will be rescaled to this range.'
        )
        warnings.warn(message, RuntimeWarning, stacklevel=2)
        image = RescaleIntensity((0, 255))(image)
    if image.data.dtype != torch.uint8:
        message = (
            'Only uint8 tensors are supported for video output. The image data'
            ' will be cast to uint8.'
        )
        warnings.warn(message, RuntimeWarning, stacklevel=2)
        image = To(torch.uint8)(image)

    # Reorient so the output looks like in typical visualization software
    direction = direction.upper()
    if direction == 'I':  # axial top to bottom
        target = 'IPL'
    elif direction == 'S':  # axial bottom to top
        target = 'SPL'
    elif direction == 'A':  # coronal back to front
        target = 'AIL'
    elif direction == 'P':  # coronal front to back
        target = 'PIL'
    elif direction == 'R':  # sagittal left to right
        target = 'RIP'
    elif direction == 'L':  # sagittal right to left
        target = 'LIP'
    else:
        message = (
            'Direction must be one of "I", "S", "P", "A", "R" or "L".'
            f' Got {direction!r}.'
        )
        raise ValueError(message)
    image = ToOrientation(target)(image)

    # Check isotropy
    spacing_f, spacing_h, spacing_w = image.spacing
    if spacing_h != spacing_w:
        message = (
            'The height and width spacings should be the same video output.'
            f' Got {spacing_h:.2f} and {spacing_w:.2f}.'
            f' Resampling both to {spacing_f:.2f}.'
        )
        warnings.warn(message, RuntimeWarning, stacklevel=2)
        spacing_iso = min(spacing_h, spacing_w)
        target_spacing = spacing_f, spacing_iso, spacing_iso
        image = Resample(target_spacing)(image)  # type: ignore[assignment]

    # Check that height and width are multiples of 2 for H.265 encoding
    num_frames, height, width = image.spatial_shape
    if height % 2 != 0 or width % 2 != 0:
        message = (
            f'The height ({height}) and width ({width}) must be even.'
            ' The image will be cropped to the nearest even number.'
        )
        warnings.warn(message, RuntimeWarning, stacklevel=2)
        image = EnsureShapeMultiple((1, 2, 2), method='crop')(image)

    if seconds is not None:
        frame_rate = num_frames / seconds

    output_path = Path(output_path)
    if output_path.suffix.lower() != '.mp4':
        message = 'Only .mp4 files are supported for video output.'
        raise NotImplementedError(message)

    frames = image.numpy()[0]
    first = frames[0]
    height, width = first.shape

    process = (
        ffmpeg.input(
            'pipe:',
            format='rawvideo',
            pix_fmt='gray',
            s=f'{width}x{height}',
            framerate=frame_rate,
        )
        .output(
            str(output_path),
            vcodec='libx265',
            pix_fmt='yuv420p',
            loglevel=verbosity,
            **{'x265-params': f'log-level={verbosity}'},
        )
        .overwrite_output()
        .run_async(pipe_stdin=True)
    )

    for array in frames:
        buffer = array.tobytes()
        process.stdin.write(buffer)

    process.stdin.close()
    process.wait()


1			from __future__ import annotations
2
3			import warnings
4			from itertools import cycle
5			from pathlib import Path
6			from typing import TYPE_CHECKING
7			from typing import Any
8
9			import numpy as np
10			import torch
11			from einops import rearrange
12
13			from .data.image import Image
14			from .data.image import LabelMap
15			from .data.image import ScalarImage
16			from .data.subject import Subject
17			from .external.imports import get_ffmpeg
18			from .transforms.preprocessing.intensity.rescale import RescaleIntensity
19			from .transforms.preprocessing.intensity.to import To
20			from .transforms.preprocessing.spatial.ensure_shape_multiple import EnsureShapeMultiple
21			from .transforms.preprocessing.spatial.resample import Resample
22			from .transforms.preprocessing.spatial.to_canonical import ToCanonical
23			from .transforms.preprocessing.spatial.to_orientation import ToOrientation
24			from .types import TypePath
25
26			if TYPE_CHECKING:
27			from matplotlib.colors import BoundaryNorm
28			from matplotlib.colors import ListedColormap
29			from matplotlib.figure import Figure
30
31
32			def import_mpl_plt():
33			try:
34			import matplotlib as mpl
35			import matplotlib.pyplot as plt
36			except ImportError as e:
37			raise ImportError('Install matplotlib for plotting support') from e
38			return mpl, plt
39
40
41			def rotate(image: np.ndarray, *, radiological: bool = True, n: int = -1) -> np.ndarray:
42			# Rotate for visualization purposes
43			image = np.rot90(image, n, axes=(0, 1))
44			if radiological:
45			image = np.fliplr(image)
46			return image
47
48
49			def _create_categorical_colormap(
50			data: torch.Tensor,
51			cmap_name: str = 'glasbey_category10',
52			) -> tuple[ListedColormap, BoundaryNorm]:
53			num_classes = int(data.max())
54			mpl, _ = import_mpl_plt()
55
56			colors = [
57			(0, 0, 0), # black for background
58			(1, 1, 1), # white for class 1
59			]
60			if num_classes > 1:
61			from .external.imports import get_colorcet
62
63			colorcet = get_colorcet()
64			cmap = getattr(colorcet.cm, cmap_name)
65			color_cycle = cycle(cmap.colors)
66			distinct_colors = [next(color_cycle) for _ in range(num_classes - 1)]
67			colors.extend(distinct_colors)
68			boundaries = np.arange(-0.5, num_classes + 1.5, 1)
69			colormap = mpl.colors.ListedColormap(colors)
70			boundary_norm = mpl.colors.BoundaryNorm(boundaries, ncolors=colormap.N)
71			return colormap, boundary_norm
72
73
74			def plot_volume(
75			image: Image,
76			radiological=True,
77			channel=None,
78			axes=None,
79			cmap=None,
80			output_path=None,
81			show=True,
82			xlabels=True,
83			percentiles: tuple[float, float] = (0, 100),
84			figsize=None,
85			title=None,
86			reorient=True,
87			indices=None,
88			rgb=True,
89			savefig_kwargs: dict[str, Any] \| None = None,
90			**imshow_kwargs,
91			) -> Figure \| None:
92			_, plt = import_mpl_plt()
93			fig: Figure \| None = None
			0 ignored issues – show introduced 2025-08-08 10:48 UTC by Report Bug Copy Issue Report The variable `Figure` does not seem to be defined in case `TYPE_CHECKING` on line `26` is `False`. Are you sure this can never be the case? Loading history...
94			if axes is None:
95			fig, axes = plt.subplots(1, 3, figsize=figsize)
96
97			if reorient:
98			image = ToCanonical()(image) # type: ignore[assignment]
99
100			is_label = isinstance(image, LabelMap)
101			if is_label: # probabilistic label map
102			data = image.data[np.newaxis, -1]
103			elif rgb and image.num_channels == 3:
104			data = image.data # keep image as it is
105			elif channel is None:
106			data = image.data[0:1] # just use the first channel
107			else:
108			data = image.data[np.newaxis, channel]
109			data = rearrange(data, 'c x y z -> x y z c')
110			data_numpy: np.ndarray = data.cpu().numpy()
111
112			if indices is None:
113			indices = np.array(data_numpy.shape[:3]) // 2
114			i, j, k = indices
115			slice_x = rotate(data_numpy[i, :, :], radiological=radiological)
116			slice_y = rotate(data_numpy[:, j, :], radiological=radiological)
117			slice_z = rotate(data_numpy[:, :, k], radiological=radiological)
118
119			if isinstance(cmap, dict):
120			slices = slice_x, slice_y, slice_z
121			slice_x, slice_y, slice_z = color_labels(slices, cmap)
122			else:
123			boundary_norm = None
124			if cmap is None:
125			if is_label:
126			cmap, boundary_norm = _create_categorical_colormap(data)
127			else:
128			cmap = 'gray'
129			imshow_kwargs['cmap'] = cmap
130			imshow_kwargs['norm'] = boundary_norm
131
132			if is_label:
133			imshow_kwargs['interpolation'] = 'none'
134			else:
135			if 'interpolation' not in imshow_kwargs:
136			imshow_kwargs['interpolation'] = 'bicubic'
137
138			imshow_kwargs['origin'] = 'lower'
139
140			if not is_label:
141			displayed_data = np.concatenate(
142			[
143			slice_x.flatten(),
144			slice_y.flatten(),
145			slice_z.flatten(),
146			]
147			)
148			p1, p2 = np.percentile(displayed_data, percentiles)
149			imshow_kwargs['vmin'] = p1
150			imshow_kwargs['vmax'] = p2
151
152			spacing_r, spacing_a, spacing_s = image.spacing
153			sag_axis, cor_axis, axi_axis = axes
154			slices_dict = {
155			'Sagittal': {
156			'aspect': spacing_s / spacing_a,
157			'slice': slice_x,
158			'xlabel': 'A',
159			'ylabel': 'S',
160			'axis': sag_axis,
161			},
162			'Coronal': {
163			'aspect': spacing_s / spacing_r,
164			'slice': slice_y,
165			'xlabel': 'R',
166			'ylabel': 'S',
167			'axis': cor_axis,
168			},
169			'Axial': {
170			'aspect': spacing_a / spacing_r,
171			'slice': slice_z,
172			'xlabel': 'R',
173			'ylabel': 'A',
174			'axis': axi_axis,
175			},
176			}
177
178			for axis_title, info in slices_dict.items():
179			axis = info['axis']
180			axis.imshow(info['slice'], aspect=info['aspect'], **imshow_kwargs)
181			if xlabels:
182			axis.set_xlabel(info['xlabel'])
183			axis.set_ylabel(info['ylabel'])
184			axis.invert_xaxis()
185			axis.set_title(axis_title)
186
187			plt.tight_layout()
188			if title is not None:
189			plt.suptitle(title)
190
191			if output_path is not None and fig is not None:
192			if savefig_kwargs is None:
193			savefig_kwargs = {}
194			fig.savefig(output_path, **savefig_kwargs)
195			if show:
196			plt.show()
197			return fig
198
199
200			def plot_subject(
201			subject: Subject,
202			cmap_dict=None,
203			show=True,
204			output_path=None,
205			figsize=None,
206			clear_axes=True,
207			**plot_volume_kwargs,
208			):
209			_, plt = import_mpl_plt()
210			num_images = len(subject)
211			many_images = num_images > 2
212			subplots_kwargs = {'figsize': figsize}
213			try:
214			if clear_axes:
215			subject.check_consistent_spatial_shape()
216			subplots_kwargs['sharex'] = 'row' if many_images else 'col'
217			subplots_kwargs['sharey'] = 'row' if many_images else 'col'
218			except RuntimeError: # different shapes in subject
219			pass
220			args = (3, num_images) if many_images else (num_images, 3)
221			fig, axes = plt.subplots(args, *subplots_kwargs)
222			# The array of axes must be 2D so that it can be indexed correctly within
223			# the plot_volume() function
224			axes = axes.T if many_images else axes.reshape(-1, 3)
225			iterable = enumerate(subject.get_images_dict(intensity_only=False).items())
226			axes_names = 'sagittal', 'coronal', 'axial'
227			for image_index, (name, image) in iterable:
228			image_axes = axes[image_index]
229			cmap = None
230			if cmap_dict is not None and name in cmap_dict:
231			cmap = cmap_dict[name]
232			last_row = image_index == len(axes) - 1
233			plot_volume(
234			image,
235			axes=image_axes,
236			show=False,
237			cmap=cmap,
238			xlabels=last_row,
239			**plot_volume_kwargs,
240			)
241			for axis, axis_name in zip(image_axes, axes_names):
242			axis.set_title(f'{name} ({axis_name})')
243			plt.tight_layout()
244			if output_path is not None:
245			fig.savefig(output_path)
246			if show:
247			plt.show()
248
249
250			def get_num_bins(x: np.ndarray) -> int:
251			"""Get the optimal number of bins for a histogram.
252
253			This method uses the Freedman–Diaconis rule to compute the histogram that
254			minimizes "the integral of the squared difference between the histogram
255			(i.e., relative frequency density) and the density of the theoretical
256			probability distribution" (`Wikipedia <https://en.wikipedia.org/wiki/Freedman%E2%80%93Diaconis_rule>`_).
257
258			Args:
259			x: Input values.
260			"""
261			# Freedman–Diaconis number of bins
262			q25, q75 = np.percentile(x, [25, 75])
263			bin_width = 2 * (q75 - q25) * len(x) ** (-1 / 3)
264			bins = round((x.max() - x.min()) / bin_width)
265			return bins
266
267
268			def plot_histogram(x: np.ndarray, show=True, **kwargs) -> None:
269			_, plt = import_mpl_plt()
270			plt.hist(x, bins=get_num_bins(x), **kwargs)
271			plt.xlabel('Intensity')
272			density = kwargs.pop('density', False)
273			ylabel = 'Density' if density else 'Frequency'
274			plt.ylabel(ylabel)
275			if show:
276			plt.show()
277
278
279			def color_labels(arrays, cmap_dict):
280			results = []
281			for slice_array in arrays:
282			si, sj, _ = slice_array.shape
283			rgb = np.zeros((si, sj, 3), dtype=np.uint8)
284			for label, color in cmap_dict.items():
285			if isinstance(color, str):
286			mpl, _ = import_mpl_plt()
287			color = mpl.colors.to_rgb(color)
288			color = [255 * n for n in color]
289			rgb[slice_array[..., 0] == label] = color
290			results.append(rgb)
291			return results
292
293
294			def make_gif(
295			tensor: torch.Tensor,
296			axis: int,
297			duration: float, # of full gif
298			output_path: TypePath,
299			loop: int = 0,
300			optimize: bool = True,
301			rescale: bool = True,
302			reverse: bool = False,
303			) -> None:
304			try:
305			from PIL import Image as ImagePIL
306			except ModuleNotFoundError as e:
307			message = 'Please install Pillow to use Image.to_gif(): pip install Pillow'
308			raise RuntimeError(message) from e
309			transform = RescaleIntensity((0, 255))
310			tensor = transform(tensor) if rescale else tensor # type: ignore[assignment]
311			single_channel = len(tensor) == 1
312
313			# Move channels dimension to the end and bring selected axis to 0
314			axes = np.roll(range(1, 4), -axis)
315			tensor = tensor.permute(*axes, 0)
316
317			if single_channel:
318			mode = 'P'
319			tensor = tensor[..., 0]
320			else:
321			mode = 'RGB'
322			array = tensor.byte().numpy()
323			n = 2 if axis == 1 else 1
324			images = [ImagePIL.fromarray(rotate(i, n=n)).convert(mode) for i in array]
325			num_images = len(images)
326			images = list(reversed(images)) if reverse else images
327			frame_duration_ms = duration / num_images * 1000
328			if frame_duration_ms < 10:
329			fps = round(1000 / frame_duration_ms)
330			frame_duration_ms = 10
331			new_duration = frame_duration_ms * num_images / 1000
332			message = (
333			'The computed frame rate from the given duration is too high'
334			f' ({fps} fps). The highest possible frame rate in the GIF'
335			' file format specification is 100 fps. The duration has been set'
336			f' to {new_duration:.1f} seconds, instead of {duration:.1f}'
337			)
338			warnings.warn(message, RuntimeWarning, stacklevel=2)
339			images[0].save(
340			output_path,
341			save_all=True,
342			append_images=images[1:],
343			optimize=optimize,
344			duration=frame_duration_ms,
345			loop=loop,
346			)
347
348
349			def make_video(
350			image: ScalarImage,
351			output_path: TypePath,
352			seconds: float \| None = None,
353			frame_rate: float \| None = None,
354			direction: str = 'I',
355			verbosity: str = 'error',
356			) -> None:
357			ffmpeg = get_ffmpeg()
358
359			if seconds is None and frame_rate is None:
360			message = 'Either seconds or frame_rate must be provided.'
361			raise ValueError(message)
362			if seconds is not None and frame_rate is not None:
363			message = 'Provide either seconds or frame_rate, not both.'
364			raise ValueError(message)
365			if image.num_channels > 1:
366			message = 'Only single-channel tensors are supported for video output for now.'
367			raise ValueError(message)
368			tmin, tmax = image.data.min(), image.data.max()
369			if tmin < 0 or tmax > 255:
370			message = (
371			'The tensor must be in the range [0, 256) for video output.'
372			' The image data will be rescaled to this range.'
373			)
374			warnings.warn(message, RuntimeWarning, stacklevel=2)
375			image = RescaleIntensity((0, 255))(image)
376			if image.data.dtype != torch.uint8:
377			message = (
378			'Only uint8 tensors are supported for video output. The image data'
379			' will be cast to uint8.'
380			)
381			warnings.warn(message, RuntimeWarning, stacklevel=2)
382			image = To(torch.uint8)(image)
383
384			# Reorient so the output looks like in typical visualization software
385			direction = direction.upper()
386			if direction == 'I': # axial top to bottom
387			target = 'IPL'
388			elif direction == 'S': # axial bottom to top
389			target = 'SPL'
390			elif direction == 'A': # coronal back to front
391			target = 'AIL'
392			elif direction == 'P': # coronal front to back
393			target = 'PIL'
394			elif direction == 'R': # sagittal left to right
395			target = 'RIP'
396			elif direction == 'L': # sagittal right to left
397			target = 'LIP'
398			else:
399			message = (
400			'Direction must be one of "I", "S", "P", "A", "R" or "L".'
401			f' Got {direction!r}.'
402			)
403			raise ValueError(message)
404			image = ToOrientation(target)(image)
405
406			# Check isotropy
407			spacing_f, spacing_h, spacing_w = image.spacing
408			if spacing_h != spacing_w:
409			message = (
410			'The height and width spacings should be the same video output.'
411			f' Got {spacing_h:.2f} and {spacing_w:.2f}.'
412			f' Resampling both to {spacing_f:.2f}.'
413			)
414			warnings.warn(message, RuntimeWarning, stacklevel=2)
415			spacing_iso = min(spacing_h, spacing_w)
416			target_spacing = spacing_f, spacing_iso, spacing_iso
417			image = Resample(target_spacing)(image) # type: ignore[assignment]
418
419			# Check that height and width are multiples of 2 for H.265 encoding
420			num_frames, height, width = image.spatial_shape
421			if height % 2 != 0 or width % 2 != 0:
422			message = (
423			f'The height ({height}) and width ({width}) must be even.'
424			' The image will be cropped to the nearest even number.'
425			)
426			warnings.warn(message, RuntimeWarning, stacklevel=2)
427			image = EnsureShapeMultiple((1, 2, 2), method='crop')(image)
428
429			if seconds is not None:
430			frame_rate = num_frames / seconds
431
432			output_path = Path(output_path)
433			if output_path.suffix.lower() != '.mp4':
434			message = 'Only .mp4 files are supported for video output.'
435			raise NotImplementedError(message)
436
437			frames = image.numpy()[0]
438			first = frames[0]
439			height, width = first.shape
440
441			process = (
442			ffmpeg.input(
443			'pipe:',
444			format='rawvideo',
445			pix_fmt='gray',
446			s=f'{width}x{height}',
447			framerate=frame_rate,
448			)
449			.output(
450			str(output_path),
451			vcodec='libx265',
452			pix_fmt='yuv420p',
453			loglevel=verbosity,
454			**{'x265-params': f'log-level={verbosity}'},
455			)
456			.overwrite_output()
457			.run_async(pipe_stdin=True)
458			)
459
460			for array in frames:
461			buffer = array.tobytes()
462			process.stdin.write(buffer)
463
464			process.stdin.close()
465			process.wait()
466

TorchIO-project / torchio

torchio.visualization.make_video() F last analyzed 2025-08-12 20:07 UTC

Complexity

Size

Duplication

Importance

How to fix Long Method Complexity

Long Method

Complexity

Duplication Side-by-Side

Filter issues like

torchio.visualization.make_video() F
last analyzed 2025-08-12 20:07 UTC