Skip to content

Reconstruction Module

ASL Mapping Reconstruction

This is the main module in asltk which aims to offer quantitative reconstruction algorithms for the ASLData.

In general, the ASL processing requires data fitting and manipulation to result in quantitative mappings. A classical example are the CBF and ATT mapping.

One of the major objective in asltk is to provide the state-of-the-art reconstruction methods available in the ASL research community, allowing an easy to use application to many studies.

Note

It is intended to the reconstruction module to be evolve by the contribution of scientific community. If you want to apply a new reconstruction method in the asltk project, please see more details at the How to Contribute section and into the reconstruction classes standards

Reconstruction class standard

In order to comply a standarized way to build the reconstruction classes, we assumed the following structure:

  1. All reconstruction classes must hierent from MRIParameters data class, which informs the base MRI parameters values for a generic data analysis and fitting
  2. The base reconstruction class must have a constructor (__init__ method) with an input data as the ASLData or other types (if needed).
  3. Even though there is not obligation (regulated by abstract methods) for determined class methods, it is expected that the reconstruction class can implement a create_map method, which is the actual reconstruction calculation that exposes (as an method output) the reconstructed data from ASLData.

Info

Depending on the specifity of the reconstruction method, it can vary the way that the inner reconstruction methods can be implemented. However, we try to maintain a basic API usage to get a more uniform implementation thoughout the community contribution.

Reconstruction Module API

CBFMapping

Bases: MRIParameters

Source code in asltk/reconstruction/cbf_mapping.py
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
class CBFMapping(MRIParameters):
    def __init__(self, asl_data: ASLData) -> None:
        """Basic CBFMapping constructor.

        Notes:
            The ASLData is the base data used in the object constructor.
            In order to create the CBF map correctly, a proper ASLData must be
            provided. Check whether the ASLData given as input is defined
            correctly

        Examples:
            The default MRIParameters are used as default in the object
            constructor
            >>> asl_data = ASLData(pcasl='./tests/files/pcasl_mte.nii.gz',m0='./tests/files/m0.nii.gz')
            >>> cbf = CBFMapping(asl_data)
            >>> cbf.get_constant('T1csf')
            1400.0

            If the user want to change the MRIParameter value, for a specific
            object, one can change it directly:
            >>> cbf.set_constant(1600.0, 'T1csf')
            >>> cbf.get_constant('T1csf')
            1600.0
            >>> default_param = MRIParameters()
            >>> default_param.get_constant('T1csf')
            1400.0

        Args:
            asl_data (ASLData): The ASL data object (ASLData)
        """
        super().__init__()
        self._asl_data = asl_data
        if self._asl_data('m0') is None:
            raise ValueError(
                'ASLData is incomplete. CBFMapping need pcasl and m0 images.'
            )

        self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
        self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

    def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
        """Defines a brain mask to limit CBF mapping calculations to specific regions.

        A brain mask significantly improves processing speed by limiting calculations
        to brain tissue voxels and excluding background regions. It also improves
        the quality of results by focusing the fitting algorithm on relevant tissue.

        A image mask is simply an image that defines the voxels where the ASL
        calculation should be made. The mask should have the same spatial dimensions
        as the M0 reference image.

        A most common approach is to use a binary image (zeros for background
        and 1 for brain tissues). However, the method can also handle multi-label
        masks by specifying which label value represents brain tissue.

        Args:
            brain_mask (np.ndarray): The image representing the brain mask.
                Must match the spatial dimensions of the M0 image.
            label (int, optional): The label value used to define the foreground
                tissue (brain). Defaults to 1. Voxels with this value will be
                included in processing.

        Examples:
            Use a binary brain mask:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import CBFMapping
            >>> import numpy as np
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     ld_values=[1.8], pld_values=[1.8]
            ... )
            >>> cbf_mapper = CBFMapping(asl_data)
            >>> # Create a simple brain mask (center region only)
            >>> mask_shape = asl_data('m0').get_as_numpy().shape  # Get M0 dimensions
            >>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
            >>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
            >>> adjusted_brain_mask[2:6, 10:25, 10:25] = 1  # Define brain region
            >>> brain_mask.update_image_data(adjusted_brain_mask)
            >>> cbf_mapper.set_brain_mask(brain_mask)

            Load and use an existing brain mask:
            >>> # Load pre-computed brain mask
            >>> from asltk.utils.io import ImageIO
            >>> brain_mask = ImageIO('./tests/files/m0_brain_mask.nii.gz')
            >>> cbf_mapper.set_brain_mask(brain_mask)

            Use multi-label mask (select specific region):
            >>> # Assuming a segmentation mask with different tissue labels
            >>> segmentation_mask = ImageIO(image_array=np.random.randint(0, 4, mask_shape))  # Example
            >>> # Use only label 2 (e.g., grey matter)
            >>> cbf_mapper.set_brain_mask(segmentation_mask, label=2)

            Automatic thresholding of M0 image as mask:
            >>> # Use M0 intensity to create brain mask
            >>> m0_data = asl_data('m0').get_as_numpy()
            >>> threshold = np.percentile(m0_data, 20)  # Bottom 20% as background
            >>> auto_mask = ImageIO(image_array=(m0_data > threshold).astype(np.uint8))
            >>> cbf_mapper.set_brain_mask(auto_mask)

        Raises:
            ValueError: If brain_mask dimensions don't match M0 image dimensions.
        """
        logger = get_logger('cbf_mapping')
        logger.info(f'Setting brain mask with label {label}')

        if not isinstance(brain_mask, ImageIO):
            raise ValueError(
                f'mask is not an ImageIO object. Type {type(brain_mask)}'
            )

        brain_mask_array = brain_mask.get_as_numpy()

        _check_mask_values(
            brain_mask, label, self._asl_data('m0').get_as_numpy().shape
        )

        binary_mask = (brain_mask_array == label).astype(np.uint8) * label
        self._brain_mask = binary_mask

        mask_volume = np.sum(binary_mask > 0)
        logger.info(f'Brain mask set successfully: {mask_volume} voxels')

    def get_brain_mask(self):
        """Get the current brain mask image being used for CBF calculations.

        Returns:
            np.ndarray: The brain mask image as a binary array where 1 indicates
                brain tissue voxels that will be processed, and 0 indicates
                background voxels that will be skipped.

        Examples:
            Check if a brain mask has been set:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import CBFMapping
            >>> import numpy as np
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     ld_values=[1.8], pld_values=[1.8]
            ... )
            >>> cbf_mapper = CBFMapping(asl_data)
            >>> # Initially, mask covers entire volume
            >>> current_mask = cbf_mapper.get_brain_mask()

            Verify brain mask after setting:
            >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
            >>> new_brain_mask = brain_mask.get_as_numpy().copy()
            >>> new_brain_mask[0:4, :, :] = 0  # Remove some slices
            >>> brain_mask.update_image_data(new_brain_mask)
            >>> cbf_mapper.set_brain_mask(brain_mask)
            >>> updated_mask = cbf_mapper.get_brain_mask()
        """
        return self._brain_mask

    def create_map(
        self,
        ub=[1.0, 5000.0],
        lb=[0.0, 0.0],
        par0=[1e-5, 1000],
        cores='auto',
        smoothing=None,
        smoothing_params=None,
    ):
        """Create the CBF and also ATT maps using the Buxton ASL model.

        This method performs voxel-wise non-linear fitting of the Buxton ASL model
        to generate Cerebral Blood Flow (CBF) and Arterial Transit Time (ATT) maps.
        The fitting is performed in parallel using multiple CPU cores for efficiency.

        Note:
            By default the ATT map is already calculated using the same Buxton
            formalism. Once the CBFMapping.create_map() method is called, both
            CBF and ATT maps are given in the output.

        Note:
            The CBF maps is given in two formats: the original pixel scale,
            resulted from the non-linear Buxton model fitting, and also
            a normalized version with the correct units of mL/100 g/min. In the
            output dictionary the user can select the 'cbf' and 'cbf_norm'
            options

        Args:
            ub (list, optional): The upper bounds for [CBF, ATT] fitting parameters.
                Defaults to [1.0, 5000.0]. CBF in relative units, ATT in ms.
            lb (list, optional): The lower bounds for [CBF, ATT] fitting parameters.
                Defaults to [0.0, 0.0]. Both parameters must be non-negative.
            par0 (list, optional): The initial guess for [CBF, ATT] parameters.
                Defaults to [1e-5, 1000]. Good starting values help convergence.
            cores (int, optional): Number of CPU threads to use for parallel processing.
                Defaults to using all available threads. Use fewer cores to preserve
                system resources.
            smoothing (str, optional): Type of spatial smoothing filter to apply.
                Options: None (default, no smoothing), 'gaussian', 'median'.
                Smoothing is applied to all output maps after reconstruction.
            smoothing_params (dict, optional): Parameters for the smoothing filter.
                For 'gaussian': {'sigma': float} (default: 1.0)
                For 'median': {'size': int} (default: 3)

        Returns:
            dict: A dictionary containing:
                - 'cbf': Raw CBF map in model units (numpy.ndarray)
                - 'cbf_norm': Normalized CBF map in mL/100g/min (numpy.ndarray)
                - 'att': ATT map in milliseconds (numpy.ndarray)
                All maps are smoothed if smoothing is enabled.

        Examples:  # doctest: +SKIP
            Basic CBF mapping with default parameters:
            >>> from asltk.asldata import ASLData
            >>> from asltk.utils.io import ImageIO
            >>> from asltk.reconstruction import CBFMapping
            >>> import numpy as np
            >>> # Load ASL data with LD/PLD values
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     ld_values=[1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8]
            ... )
            >>> cbf_mapper = CBFMapping(asl_data)
            >>> # Set brain mask (recommended for faster processing)
            >>> brain_mask = ImageIO(image_array=np.ones((5, 35, 35)))  # Example mask
            >>> cbf_mapper.set_brain_mask(brain_mask)
            >>> # Generate maps
            >>> results = cbf_mapper.create_map() # doctest: +SKIP

            Custom parameter bounds for specific tissue properties:
            >>> # For grey matter regions (higher CBF expected)
            >>> results_gm = cbf_mapper.create_map(
            ...     ub=[2.0, 3000.0],      # Higher CBF upper bound
            ...     lb=[0.1, 500.0],       # Reasonable lower bounds
            ...     par0=[0.5, 1200.0]     # Good initial guess for GM
            ... ) # doctest: +SKIP

            Apply spatial smoothing to reduce noise:
            >>> # Gaussian smoothing with default sigma=1.0
            >>> results_smooth = cbf_mapper.create_map(
            ...     smoothing='gaussian'
            ... ) # doctest: +SKIP
            >>> # Custom smoothing parameters
            >>> results_custom = cbf_mapper.create_map(
            ...     smoothing='gaussian',
            ...     smoothing_params={'sigma': 1.5}
            ... ) # doctest: +SKIP
            >>> # Median filtering for edge-preserving smoothing
            >>> results_median = cbf_mapper.create_map(
            ...     smoothing='median',
            ...     smoothing_params={'size': 5}
            ... ) # doctest: +SKIP

            Memory-efficient processing with limited cores:
            >>> # Use only 4 cores to preserve system resources
            >>> results = cbf_mapper.create_map(cores=4) # doctest: +SKIP

        Raises:
            ValueError: If cores parameter is invalid, or if LD/PLD values are missing.
        """
        logger = get_logger('cbf_mapping')
        logger.info('Starting CBF map creation')

        if not isinstance(cores, str):
            if (
                (cores < 0)
                or (cores > cpu_count())
                or not isinstance(cores, int)
            ):
                error_msg = 'Number of CPU cores must be at least 1 and less than maximum cores available.'
                logger.error(
                    f'{error_msg} Requested: {cores}, Available: {cpu_count()}'
                )
                raise ValueError(error_msg)
        elif isinstance(cores, str):
            if cores not in ['auto']:
                error_msg = (
                    'Cores parameter must be either "auto" or a integer.'
                )
                logger.error(error_msg)
                raise ValueError(error_msg)
        else:
            raise ValueError(
                'Cores parameter must be either "auto" or a integer.'
            )

        if (
            len(self._asl_data.get_ld()) == 0
            or len(self._asl_data.get_pld()) == 0
        ):
            error_msg = 'LD or PLD list of values must be provided.'
            logger.error(error_msg)
            raise ValueError(error_msg)

        logger.info(f'Using {cores} CPU cores for parallel processing')
        log_processing_step('Initializing CBF mapping computation')

        global asl_data, brain_mask
        asl_data = self._asl_data
        brain_mask = self._brain_mask

        BuxtonX = [self._asl_data.get_ld(), self._asl_data.get_pld()]

        x_axis, y_axis, z_axis = (
            self._asl_data('m0').get_as_numpy().shape[2],
            self._asl_data('m0').get_as_numpy().shape[1],
            self._asl_data('m0').get_as_numpy().shape[0],
        )

        logger.info(
            f'Processing volume dimensions: {z_axis}x{y_axis}x{x_axis}'
        )

        cbf_map_shared = Array('f', z_axis * y_axis * x_axis, lock=False)
        att_map_shared = Array('f', z_axis * y_axis * x_axis, lock=False)

        # Estimate all the memory usage needed for each core processing
        asldata_memory = estimate_memory_usage(
            self._asl_data('pcasl').get_as_numpy()
        )
        brain_mask_memory = estimate_memory_usage(self._brain_mask)
        cbf_memory = estimate_memory_usage(self._cbf_map)
        att_memory = estimate_memory_usage(self._att_map)

        actual_cores = get_optimal_core_count(
            cores,
            sum([asldata_memory, brain_mask_memory, cbf_memory, att_memory]),
        )

        # Make a copy of base information
        m0_array = asl_data('m0').get_as_numpy()
        pcasl_array = asl_data('pcasl').get_as_numpy()

        log_processing_step(
            'Running voxel-wise CBF fitting', 'this may take several minutes'
        )
        with Pool(
            processes=actual_cores,
            initializer=_cbf_init_globals,
            initargs=(cbf_map_shared, att_map_shared, brain_mask, asl_data),
        ) as pool:
            with Progress() as progress:
                task = progress.add_task('CBF/ATT processing...', total=x_axis)
                results = [
                    pool.apply_async(
                        _cbf_process_slice,
                        args=(
                            i,
                            x_axis,
                            y_axis,
                            z_axis,
                            BuxtonX,
                            par0,
                            lb,
                            ub,
                            m0_array,
                            pcasl_array,
                        ),
                        callback=lambda _: progress.update(task, advance=1),
                    )
                    for i in range(x_axis)
                ]
                for result in results:
                    result.wait()

        self._cbf_map = np.frombuffer(
            cbf_map_shared, dtype=np.float32
        ).reshape(z_axis, y_axis, x_axis)
        self._att_map = np.frombuffer(
            att_map_shared, dtype=np.float32
        ).reshape(z_axis, y_axis, x_axis)

        # Log completion statistics
        cbf_values = self._cbf_map[brain_mask > 0]
        att_values = self._att_map[brain_mask > 0]

        logger.info(f'CBF mapping completed successfully')
        logger.info(
            f'CBF statistics - Mean: {np.mean(cbf_values):.4f}, Std: {np.std(cbf_values):.4f}'
        )
        logger.info(
            f'ATT statistics - Mean: {np.mean(att_values):.4f}, Std: {np.std(att_values):.4f}'
        )

        # Prepare output maps
        base_volume = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_image = clone_image(base_volume)
        cbf_map_image.update_image_data(
            self._cbf_map, self._asl_data._asl_image._average_m0
        )

        cbf_map_norm_image = clone_image(base_volume)
        cbf_map_norm_image.update_image_data(
            self._cbf_map * (60 * 60 * 1000),
            self._asl_data._asl_image._average_m0,
        )

        att_map_image = clone_image(base_volume)
        att_map_image.update_image_data(
            self._att_map, self._asl_data._asl_image._average_m0
        )

        output_maps = {
            'cbf': cbf_map_image,
            'cbf_norm': cbf_map_norm_image,
            'att': att_map_image,
        }

        # Apply smoothing if requested
        return _apply_smoothing_to_maps(
            output_maps, smoothing, smoothing_params
        )

__init__(asl_data)

Basic CBFMapping constructor.

Notes

The ASLData is the base data used in the object constructor. In order to create the CBF map correctly, a proper ASLData must be provided. Check whether the ASLData given as input is defined correctly

Examples:

The default MRIParameters are used as default in the object constructor

>>> asl_data = ASLData(pcasl='./tests/files/pcasl_mte.nii.gz',m0='./tests/files/m0.nii.gz')
>>> cbf = CBFMapping(asl_data)
>>> cbf.get_constant('T1csf')
1400.0

If the user want to change the MRIParameter value, for a specific object, one can change it directly:

>>> cbf.set_constant(1600.0, 'T1csf')
>>> cbf.get_constant('T1csf')
1600.0
>>> default_param = MRIParameters()
>>> default_param.get_constant('T1csf')
1400.0

Parameters:

Name Type Description Default
asl_data ASLData

The ASL data object (ASLData)

required
Source code in asltk/reconstruction/cbf_mapping.py
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
def __init__(self, asl_data: ASLData) -> None:
    """Basic CBFMapping constructor.

    Notes:
        The ASLData is the base data used in the object constructor.
        In order to create the CBF map correctly, a proper ASLData must be
        provided. Check whether the ASLData given as input is defined
        correctly

    Examples:
        The default MRIParameters are used as default in the object
        constructor
        >>> asl_data = ASLData(pcasl='./tests/files/pcasl_mte.nii.gz',m0='./tests/files/m0.nii.gz')
        >>> cbf = CBFMapping(asl_data)
        >>> cbf.get_constant('T1csf')
        1400.0

        If the user want to change the MRIParameter value, for a specific
        object, one can change it directly:
        >>> cbf.set_constant(1600.0, 'T1csf')
        >>> cbf.get_constant('T1csf')
        1600.0
        >>> default_param = MRIParameters()
        >>> default_param.get_constant('T1csf')
        1400.0

    Args:
        asl_data (ASLData): The ASL data object (ASLData)
    """
    super().__init__()
    self._asl_data = asl_data
    if self._asl_data('m0') is None:
        raise ValueError(
            'ASLData is incomplete. CBFMapping need pcasl and m0 images.'
        )

    self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
    self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

create_map(ub=[1.0, 5000.0], lb=[0.0, 0.0], par0=[1e-05, 1000], cores='auto', smoothing=None, smoothing_params=None)

Create the CBF and also ATT maps using the Buxton ASL model.

This method performs voxel-wise non-linear fitting of the Buxton ASL model to generate Cerebral Blood Flow (CBF) and Arterial Transit Time (ATT) maps. The fitting is performed in parallel using multiple CPU cores for efficiency.

Note

By default the ATT map is already calculated using the same Buxton formalism. Once the CBFMapping.create_map() method is called, both CBF and ATT maps are given in the output.

Note

The CBF maps is given in two formats: the original pixel scale, resulted from the non-linear Buxton model fitting, and also a normalized version with the correct units of mL/100 g/min. In the output dictionary the user can select the 'cbf' and 'cbf_norm' options

Parameters:

Name Type Description Default
ub list

The upper bounds for [CBF, ATT] fitting parameters. Defaults to [1.0, 5000.0]. CBF in relative units, ATT in ms.

[1.0, 5000.0]
lb list

The lower bounds for [CBF, ATT] fitting parameters. Defaults to [0.0, 0.0]. Both parameters must be non-negative.

[0.0, 0.0]
par0 list

The initial guess for [CBF, ATT] parameters. Defaults to [1e-5, 1000]. Good starting values help convergence.

[1e-05, 1000]
cores int

Number of CPU threads to use for parallel processing. Defaults to using all available threads. Use fewer cores to preserve system resources.

'auto'
smoothing str

Type of spatial smoothing filter to apply. Options: None (default, no smoothing), 'gaussian', 'median'. Smoothing is applied to all output maps after reconstruction.

None
smoothing_params dict

Parameters for the smoothing filter. For 'gaussian': {'sigma': float} (default: 1.0) For 'median': {'size': int} (default: 3)

None

Returns:

Name Type Description
dict

A dictionary containing: - 'cbf': Raw CBF map in model units (numpy.ndarray) - 'cbf_norm': Normalized CBF map in mL/100g/min (numpy.ndarray) - 'att': ATT map in milliseconds (numpy.ndarray) All maps are smoothed if smoothing is enabled.

# doctest: +SKIP

Basic CBF mapping with default parameters:

>>> from asltk.asldata import ASLData
>>> from asltk.utils.io import ImageIO
>>> from asltk.reconstruction import CBFMapping
>>> import numpy as np
>>> # Load ASL data with LD/PLD values
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     ld_values=[1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8]
... )
>>> cbf_mapper = CBFMapping(asl_data)
>>> # Set brain mask (recommended for faster processing)
>>> brain_mask = ImageIO(image_array=np.ones((5, 35, 35)))  # Example mask
>>> cbf_mapper.set_brain_mask(brain_mask)
>>> # Generate maps
>>> results = cbf_mapper.create_map()

Custom parameter bounds for specific tissue properties:

>>> # For grey matter regions (higher CBF expected)
>>> results_gm = cbf_mapper.create_map(
...     ub=[2.0, 3000.0],      # Higher CBF upper bound
...     lb=[0.1, 500.0],       # Reasonable lower bounds
...     par0=[0.5, 1200.0]     # Good initial guess for GM
... )

Apply spatial smoothing to reduce noise:

>>> # Gaussian smoothing with default sigma=1.0
>>> results_smooth = cbf_mapper.create_map(
...     smoothing='gaussian'
... )
>>> # Custom smoothing parameters
>>> results_custom = cbf_mapper.create_map(
...     smoothing='gaussian',
...     smoothing_params={'sigma': 1.5}
... )
>>> # Median filtering for edge-preserving smoothing
>>> results_median = cbf_mapper.create_map(
...     smoothing='median',
...     smoothing_params={'size': 5}
... )

Memory-efficient processing with limited cores:

>>> # Use only 4 cores to preserve system resources
>>> results = cbf_mapper.create_map(cores=4)

Raises:

Type Description
ValueError

If cores parameter is invalid, or if LD/PLD values are missing.

Source code in asltk/reconstruction/cbf_mapping.py
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
def create_map(
    self,
    ub=[1.0, 5000.0],
    lb=[0.0, 0.0],
    par0=[1e-5, 1000],
    cores='auto',
    smoothing=None,
    smoothing_params=None,
):
    """Create the CBF and also ATT maps using the Buxton ASL model.

    This method performs voxel-wise non-linear fitting of the Buxton ASL model
    to generate Cerebral Blood Flow (CBF) and Arterial Transit Time (ATT) maps.
    The fitting is performed in parallel using multiple CPU cores for efficiency.

    Note:
        By default the ATT map is already calculated using the same Buxton
        formalism. Once the CBFMapping.create_map() method is called, both
        CBF and ATT maps are given in the output.

    Note:
        The CBF maps is given in two formats: the original pixel scale,
        resulted from the non-linear Buxton model fitting, and also
        a normalized version with the correct units of mL/100 g/min. In the
        output dictionary the user can select the 'cbf' and 'cbf_norm'
        options

    Args:
        ub (list, optional): The upper bounds for [CBF, ATT] fitting parameters.
            Defaults to [1.0, 5000.0]. CBF in relative units, ATT in ms.
        lb (list, optional): The lower bounds for [CBF, ATT] fitting parameters.
            Defaults to [0.0, 0.0]. Both parameters must be non-negative.
        par0 (list, optional): The initial guess for [CBF, ATT] parameters.
            Defaults to [1e-5, 1000]. Good starting values help convergence.
        cores (int, optional): Number of CPU threads to use for parallel processing.
            Defaults to using all available threads. Use fewer cores to preserve
            system resources.
        smoothing (str, optional): Type of spatial smoothing filter to apply.
            Options: None (default, no smoothing), 'gaussian', 'median'.
            Smoothing is applied to all output maps after reconstruction.
        smoothing_params (dict, optional): Parameters for the smoothing filter.
            For 'gaussian': {'sigma': float} (default: 1.0)
            For 'median': {'size': int} (default: 3)

    Returns:
        dict: A dictionary containing:
            - 'cbf': Raw CBF map in model units (numpy.ndarray)
            - 'cbf_norm': Normalized CBF map in mL/100g/min (numpy.ndarray)
            - 'att': ATT map in milliseconds (numpy.ndarray)
            All maps are smoothed if smoothing is enabled.

    Examples:  # doctest: +SKIP
        Basic CBF mapping with default parameters:
        >>> from asltk.asldata import ASLData
        >>> from asltk.utils.io import ImageIO
        >>> from asltk.reconstruction import CBFMapping
        >>> import numpy as np
        >>> # Load ASL data with LD/PLD values
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     ld_values=[1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8]
        ... )
        >>> cbf_mapper = CBFMapping(asl_data)
        >>> # Set brain mask (recommended for faster processing)
        >>> brain_mask = ImageIO(image_array=np.ones((5, 35, 35)))  # Example mask
        >>> cbf_mapper.set_brain_mask(brain_mask)
        >>> # Generate maps
        >>> results = cbf_mapper.create_map() # doctest: +SKIP

        Custom parameter bounds for specific tissue properties:
        >>> # For grey matter regions (higher CBF expected)
        >>> results_gm = cbf_mapper.create_map(
        ...     ub=[2.0, 3000.0],      # Higher CBF upper bound
        ...     lb=[0.1, 500.0],       # Reasonable lower bounds
        ...     par0=[0.5, 1200.0]     # Good initial guess for GM
        ... ) # doctest: +SKIP

        Apply spatial smoothing to reduce noise:
        >>> # Gaussian smoothing with default sigma=1.0
        >>> results_smooth = cbf_mapper.create_map(
        ...     smoothing='gaussian'
        ... ) # doctest: +SKIP
        >>> # Custom smoothing parameters
        >>> results_custom = cbf_mapper.create_map(
        ...     smoothing='gaussian',
        ...     smoothing_params={'sigma': 1.5}
        ... ) # doctest: +SKIP
        >>> # Median filtering for edge-preserving smoothing
        >>> results_median = cbf_mapper.create_map(
        ...     smoothing='median',
        ...     smoothing_params={'size': 5}
        ... ) # doctest: +SKIP

        Memory-efficient processing with limited cores:
        >>> # Use only 4 cores to preserve system resources
        >>> results = cbf_mapper.create_map(cores=4) # doctest: +SKIP

    Raises:
        ValueError: If cores parameter is invalid, or if LD/PLD values are missing.
    """
    logger = get_logger('cbf_mapping')
    logger.info('Starting CBF map creation')

    if not isinstance(cores, str):
        if (
            (cores < 0)
            or (cores > cpu_count())
            or not isinstance(cores, int)
        ):
            error_msg = 'Number of CPU cores must be at least 1 and less than maximum cores available.'
            logger.error(
                f'{error_msg} Requested: {cores}, Available: {cpu_count()}'
            )
            raise ValueError(error_msg)
    elif isinstance(cores, str):
        if cores not in ['auto']:
            error_msg = (
                'Cores parameter must be either "auto" or a integer.'
            )
            logger.error(error_msg)
            raise ValueError(error_msg)
    else:
        raise ValueError(
            'Cores parameter must be either "auto" or a integer.'
        )

    if (
        len(self._asl_data.get_ld()) == 0
        or len(self._asl_data.get_pld()) == 0
    ):
        error_msg = 'LD or PLD list of values must be provided.'
        logger.error(error_msg)
        raise ValueError(error_msg)

    logger.info(f'Using {cores} CPU cores for parallel processing')
    log_processing_step('Initializing CBF mapping computation')

    global asl_data, brain_mask
    asl_data = self._asl_data
    brain_mask = self._brain_mask

    BuxtonX = [self._asl_data.get_ld(), self._asl_data.get_pld()]

    x_axis, y_axis, z_axis = (
        self._asl_data('m0').get_as_numpy().shape[2],
        self._asl_data('m0').get_as_numpy().shape[1],
        self._asl_data('m0').get_as_numpy().shape[0],
    )

    logger.info(
        f'Processing volume dimensions: {z_axis}x{y_axis}x{x_axis}'
    )

    cbf_map_shared = Array('f', z_axis * y_axis * x_axis, lock=False)
    att_map_shared = Array('f', z_axis * y_axis * x_axis, lock=False)

    # Estimate all the memory usage needed for each core processing
    asldata_memory = estimate_memory_usage(
        self._asl_data('pcasl').get_as_numpy()
    )
    brain_mask_memory = estimate_memory_usage(self._brain_mask)
    cbf_memory = estimate_memory_usage(self._cbf_map)
    att_memory = estimate_memory_usage(self._att_map)

    actual_cores = get_optimal_core_count(
        cores,
        sum([asldata_memory, brain_mask_memory, cbf_memory, att_memory]),
    )

    # Make a copy of base information
    m0_array = asl_data('m0').get_as_numpy()
    pcasl_array = asl_data('pcasl').get_as_numpy()

    log_processing_step(
        'Running voxel-wise CBF fitting', 'this may take several minutes'
    )
    with Pool(
        processes=actual_cores,
        initializer=_cbf_init_globals,
        initargs=(cbf_map_shared, att_map_shared, brain_mask, asl_data),
    ) as pool:
        with Progress() as progress:
            task = progress.add_task('CBF/ATT processing...', total=x_axis)
            results = [
                pool.apply_async(
                    _cbf_process_slice,
                    args=(
                        i,
                        x_axis,
                        y_axis,
                        z_axis,
                        BuxtonX,
                        par0,
                        lb,
                        ub,
                        m0_array,
                        pcasl_array,
                    ),
                    callback=lambda _: progress.update(task, advance=1),
                )
                for i in range(x_axis)
            ]
            for result in results:
                result.wait()

    self._cbf_map = np.frombuffer(
        cbf_map_shared, dtype=np.float32
    ).reshape(z_axis, y_axis, x_axis)
    self._att_map = np.frombuffer(
        att_map_shared, dtype=np.float32
    ).reshape(z_axis, y_axis, x_axis)

    # Log completion statistics
    cbf_values = self._cbf_map[brain_mask > 0]
    att_values = self._att_map[brain_mask > 0]

    logger.info(f'CBF mapping completed successfully')
    logger.info(
        f'CBF statistics - Mean: {np.mean(cbf_values):.4f}, Std: {np.std(cbf_values):.4f}'
    )
    logger.info(
        f'ATT statistics - Mean: {np.mean(att_values):.4f}, Std: {np.std(att_values):.4f}'
    )

    # Prepare output maps
    base_volume = ImageIO(self._asl_data('m0').get_image_path())
    cbf_map_image = clone_image(base_volume)
    cbf_map_image.update_image_data(
        self._cbf_map, self._asl_data._asl_image._average_m0
    )

    cbf_map_norm_image = clone_image(base_volume)
    cbf_map_norm_image.update_image_data(
        self._cbf_map * (60 * 60 * 1000),
        self._asl_data._asl_image._average_m0,
    )

    att_map_image = clone_image(base_volume)
    att_map_image.update_image_data(
        self._att_map, self._asl_data._asl_image._average_m0
    )

    output_maps = {
        'cbf': cbf_map_image,
        'cbf_norm': cbf_map_norm_image,
        'att': att_map_image,
    }

    # Apply smoothing if requested
    return _apply_smoothing_to_maps(
        output_maps, smoothing, smoothing_params
    )

get_brain_mask()

Get the current brain mask image being used for CBF calculations.

Returns:

Type Description

np.ndarray: The brain mask image as a binary array where 1 indicates brain tissue voxels that will be processed, and 0 indicates background voxels that will be skipped.

Examples:

Check if a brain mask has been set:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import CBFMapping
>>> import numpy as np
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     ld_values=[1.8], pld_values=[1.8]
... )
>>> cbf_mapper = CBFMapping(asl_data)
>>> # Initially, mask covers entire volume
>>> current_mask = cbf_mapper.get_brain_mask()

Verify brain mask after setting:

>>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
>>> new_brain_mask = brain_mask.get_as_numpy().copy()
>>> new_brain_mask[0:4, :, :] = 0  # Remove some slices
>>> brain_mask.update_image_data(new_brain_mask)
>>> cbf_mapper.set_brain_mask(brain_mask)
>>> updated_mask = cbf_mapper.get_brain_mask()
Source code in asltk/reconstruction/cbf_mapping.py
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
def get_brain_mask(self):
    """Get the current brain mask image being used for CBF calculations.

    Returns:
        np.ndarray: The brain mask image as a binary array where 1 indicates
            brain tissue voxels that will be processed, and 0 indicates
            background voxels that will be skipped.

    Examples:
        Check if a brain mask has been set:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import CBFMapping
        >>> import numpy as np
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     ld_values=[1.8], pld_values=[1.8]
        ... )
        >>> cbf_mapper = CBFMapping(asl_data)
        >>> # Initially, mask covers entire volume
        >>> current_mask = cbf_mapper.get_brain_mask()

        Verify brain mask after setting:
        >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
        >>> new_brain_mask = brain_mask.get_as_numpy().copy()
        >>> new_brain_mask[0:4, :, :] = 0  # Remove some slices
        >>> brain_mask.update_image_data(new_brain_mask)
        >>> cbf_mapper.set_brain_mask(brain_mask)
        >>> updated_mask = cbf_mapper.get_brain_mask()
    """
    return self._brain_mask

set_brain_mask(brain_mask, label=1)

Defines a brain mask to limit CBF mapping calculations to specific regions.

A brain mask significantly improves processing speed by limiting calculations to brain tissue voxels and excluding background regions. It also improves the quality of results by focusing the fitting algorithm on relevant tissue.

A image mask is simply an image that defines the voxels where the ASL calculation should be made. The mask should have the same spatial dimensions as the M0 reference image.

A most common approach is to use a binary image (zeros for background and 1 for brain tissues). However, the method can also handle multi-label masks by specifying which label value represents brain tissue.

Parameters:

Name Type Description Default
brain_mask ndarray

The image representing the brain mask. Must match the spatial dimensions of the M0 image.

required
label int

The label value used to define the foreground tissue (brain). Defaults to 1. Voxels with this value will be included in processing.

1

Examples:

Use a binary brain mask:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import CBFMapping
>>> import numpy as np
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     ld_values=[1.8], pld_values=[1.8]
... )
>>> cbf_mapper = CBFMapping(asl_data)
>>> # Create a simple brain mask (center region only)
>>> mask_shape = asl_data('m0').get_as_numpy().shape  # Get M0 dimensions
>>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
>>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
>>> adjusted_brain_mask[2:6, 10:25, 10:25] = 1  # Define brain region
>>> brain_mask.update_image_data(adjusted_brain_mask)
>>> cbf_mapper.set_brain_mask(brain_mask)

Load and use an existing brain mask:

>>> # Load pre-computed brain mask
>>> from asltk.utils.io import ImageIO
>>> brain_mask = ImageIO('./tests/files/m0_brain_mask.nii.gz')
>>> cbf_mapper.set_brain_mask(brain_mask)

Use multi-label mask (select specific region):

>>> # Assuming a segmentation mask with different tissue labels
>>> segmentation_mask = ImageIO(image_array=np.random.randint(0, 4, mask_shape))  # Example
>>> # Use only label 2 (e.g., grey matter)
>>> cbf_mapper.set_brain_mask(segmentation_mask, label=2)

Automatic thresholding of M0 image as mask:

>>> # Use M0 intensity to create brain mask
>>> m0_data = asl_data('m0').get_as_numpy()
>>> threshold = np.percentile(m0_data, 20)  # Bottom 20% as background
>>> auto_mask = ImageIO(image_array=(m0_data > threshold).astype(np.uint8))
>>> cbf_mapper.set_brain_mask(auto_mask)

Raises:

Type Description
ValueError

If brain_mask dimensions don't match M0 image dimensions.

Source code in asltk/reconstruction/cbf_mapping.py
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
    """Defines a brain mask to limit CBF mapping calculations to specific regions.

    A brain mask significantly improves processing speed by limiting calculations
    to brain tissue voxels and excluding background regions. It also improves
    the quality of results by focusing the fitting algorithm on relevant tissue.

    A image mask is simply an image that defines the voxels where the ASL
    calculation should be made. The mask should have the same spatial dimensions
    as the M0 reference image.

    A most common approach is to use a binary image (zeros for background
    and 1 for brain tissues). However, the method can also handle multi-label
    masks by specifying which label value represents brain tissue.

    Args:
        brain_mask (np.ndarray): The image representing the brain mask.
            Must match the spatial dimensions of the M0 image.
        label (int, optional): The label value used to define the foreground
            tissue (brain). Defaults to 1. Voxels with this value will be
            included in processing.

    Examples:
        Use a binary brain mask:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import CBFMapping
        >>> import numpy as np
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     ld_values=[1.8], pld_values=[1.8]
        ... )
        >>> cbf_mapper = CBFMapping(asl_data)
        >>> # Create a simple brain mask (center region only)
        >>> mask_shape = asl_data('m0').get_as_numpy().shape  # Get M0 dimensions
        >>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
        >>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
        >>> adjusted_brain_mask[2:6, 10:25, 10:25] = 1  # Define brain region
        >>> brain_mask.update_image_data(adjusted_brain_mask)
        >>> cbf_mapper.set_brain_mask(brain_mask)

        Load and use an existing brain mask:
        >>> # Load pre-computed brain mask
        >>> from asltk.utils.io import ImageIO
        >>> brain_mask = ImageIO('./tests/files/m0_brain_mask.nii.gz')
        >>> cbf_mapper.set_brain_mask(brain_mask)

        Use multi-label mask (select specific region):
        >>> # Assuming a segmentation mask with different tissue labels
        >>> segmentation_mask = ImageIO(image_array=np.random.randint(0, 4, mask_shape))  # Example
        >>> # Use only label 2 (e.g., grey matter)
        >>> cbf_mapper.set_brain_mask(segmentation_mask, label=2)

        Automatic thresholding of M0 image as mask:
        >>> # Use M0 intensity to create brain mask
        >>> m0_data = asl_data('m0').get_as_numpy()
        >>> threshold = np.percentile(m0_data, 20)  # Bottom 20% as background
        >>> auto_mask = ImageIO(image_array=(m0_data > threshold).astype(np.uint8))
        >>> cbf_mapper.set_brain_mask(auto_mask)

    Raises:
        ValueError: If brain_mask dimensions don't match M0 image dimensions.
    """
    logger = get_logger('cbf_mapping')
    logger.info(f'Setting brain mask with label {label}')

    if not isinstance(brain_mask, ImageIO):
        raise ValueError(
            f'mask is not an ImageIO object. Type {type(brain_mask)}'
        )

    brain_mask_array = brain_mask.get_as_numpy()

    _check_mask_values(
        brain_mask, label, self._asl_data('m0').get_as_numpy().shape
    )

    binary_mask = (brain_mask_array == label).astype(np.uint8) * label
    self._brain_mask = binary_mask

    mask_volume = np.sum(binary_mask > 0)
    logger.info(f'Brain mask set successfully: {mask_volume} voxels')

MultiDW_ASLMapping

Bases: MRIParameters

Source code in asltk/reconstruction/multi_dw_mapping.py
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
class MultiDW_ASLMapping(MRIParameters):
    def __init__(self, asl_data: ASLData):
        """Multi-Diffusion-Weighted ASL mapping constructor for advanced perfusion analysis.

        MultiDW_ASLMapping enables sophisticated ASL analysis by incorporating multiple
        diffusion weightings (b-values) to separate intravascular and tissue
        compartments. This approach provides enhanced characterization of perfusion
        and can help differentiate between different vascular compartments.

        The class implements diffusion-weighted ASL analysis that can distinguish:
        - Fast-flowing blood (intravascular component)
        - Slow-flowing blood and tissue perfusion
        - Apparent diffusion coefficients for each compartment
        - Water exchange parameters between compartments

        Notes:
            The ASLData object must contain `dw_values` - a list of diffusion
            b-values used during ASL acquisition. These b-values are essential
            for the multi-compartment diffusion model fitting.

        Examples:
            Basic multi-DW ASL mapping setup:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import MultiDW_ASLMapping
            >>> # Create ASL data with diffusion weighting
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     dw_values=[0, 50, 100, 200],    # b-values in s/mm²
            ...     ld_values=[1.8, 1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8, 3.8]
            ... )
            >>> mdw_mapper = MultiDW_ASLMapping(asl_data)

            Access diffusion-related maps (after processing):
            >>> # These maps will be populated after create_map() is called
            >>> # A1: Signal amplitude for compartment 1
            >>> # D1: Apparent diffusion coefficient for compartment 1
            >>> # A2: Signal amplitude for compartment 2
            >>> # D2: Apparent diffusion coefficient for compartment 2
            >>> # kw: Water exchange parameter

        Args:
            asl_data (ASLData): The ASL data object containing multi-DW acquisition.
                Must include dw_values (b-values), ld_values, and pld_values.

        Raises:
            ValueError: If ASLData object lacks required DW values for
                diffusion-weighted analysis.

        See Also:
            CBFMapping: For basic CBF/ATT mapping without diffusion weighting
            MultiTE_ASLMapping: For multi-echo ASL analysis
        """
        super().__init__()
        self._asl_data = asl_data
        self._basic_maps = CBFMapping(asl_data)
        if self._asl_data.get_dw() is None:
            raise ValueError(
                'ASLData is incomplete. MultiDW_ASLMapping need a list of DW values.'
            )

        self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
        self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

        self._b_values = self._asl_data.get_dw()
        # self._A1 = np.zeros(tuple([len(self._b_values)]) + self._asl_data('m0').shape)
        self._A1 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        # self._D1 = np.zeros(tuple([1]) +self._asl_data('m0').shape)
        self._D1 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._A2 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        # self._A2 = np.zeros(tuple([len(self._b_values)])  + self._asl_data('m0').shape)
        # self._D2 = np.zeros(tuple([1]) +self._asl_data('m0').shape)
        self._D2 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._kw = np.zeros(self._asl_data('m0').get_as_numpy().shape)

    def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
        """Set brain mask for MultiDW-ASL processing (strongly recommended).

        A brain mask is especially important for multi-diffusion-weighted ASL
        processing as it significantly reduces computation time by limiting
        the intensive voxel-wise fitting to brain tissue regions only.

        Without a brain mask, processing time can be prohibitively long (hours)
        for whole-volume analysis. A proper brain mask can reduce processing
        time by 5-10x while maintaining analysis quality.

        Args:
            brain_mask (np.ndarray): The image representing the brain mask.
                Must match the spatial dimensions of the M0 image.
            label (int, optional): The label value used to define brain tissue.
                Defaults to 1. Voxels with this value will be processed.

        Examples:
            Set a brain mask for efficient processing:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import MultiDW_ASLMapping
            >>> import numpy as np
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     dw_values=[0, 50, 100], ld_values=[1.8]*3, pld_values=[1.8]*3
            ... )
            >>> mdw_mapper = MultiDW_ASLMapping(asl_data)
            >>> # Create conservative brain mask (center region only)
            >>> mask_shape = asl_data('m0').get_as_numpy().shape
            >>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
            >>> adjusted_brain_mask = brain_mask.get_as_numpy()
            >>> adjusted_brain_mask[1:4, 5:30, 5:30] = 1  # Conservative brain region
            >>> brain_mask.update_image_data(adjusted_brain_mask)
            >>> mdw_mapper.set_brain_mask(brain_mask)

        Note:
            For multi-DW ASL, consider using a more conservative (smaller) brain
            mask initially to test parameters and processing time, then expand
            to full brain analysis once satisfied with results.
        """
        if not isinstance(brain_mask, ImageIO):
            raise TypeError(
                'Brain mask must be an instance of ImageIO. '
                'Use ImageIO to load or create the mask.'
            )

        _check_mask_values(
            brain_mask, label, self._asl_data('m0').get_as_numpy().shape
        )

        binary_mask = (brain_mask.get_as_numpy() == label).astype(
            np.uint8
        ) * label
        self._brain_mask = binary_mask

    def get_brain_mask(self):
        """Get the brain mask image

        Returns:
            (np.ndarray): The brain mask image
        """
        return self._brain_mask

    def set_cbf_map(self, cbf_map: ImageIO):
        """Set the CBF map to the MultiDW_ASLMapping object.

        Note:
            The CBF maps must have the original scale in order to calculate the
            T1blGM map correclty. Hence, if the CBF map was made using
            CBFMapping class, one can use the 'cbf' output.

        Args:
            cbf_map (np.ndarray): The CBF map that is set in the MultiDW_ASLMapping object
        """
        self._cbf_map = cbf_map.get_as_numpy()

    def get_cbf_map(self) -> ImageIO:
        """Get the CBF map storaged at the MultiDW_ASLMapping object

        Returns:
            (np.ndarray): The CBF map that is storaged in the
            MultiDW_ASLMapping object
        """
        return self._cbf_map

    def set_att_map(self, att_map: ImageIO):
        """Set the ATT map to the MultiDW_ASLMapping object.

        Args:
            att_map (np.ndarray): The ATT map that is set in the MultiDW_ASLMapping object
        """
        self._att_map = att_map.get_as_numpy()

    def get_att_map(self):
        """Get the ATT map storaged at the MultiDW_ASLMapping object

        Returns:
            (np.ndarray): _description_
        """
        return self._att_map

    def create_map(
        self,
        lb: list = [0.0, 0.0, 0.0, 0.0],
        ub: list = [np.inf, np.inf, np.inf, np.inf],
        par0: list = [0.5, 0.000005, 0.5, 0.000005],
        smoothing=None,
        smoothing_params=None,
    ):
        """Create multi-diffusion-weighted ASL maps for compartment analysis.

        This method performs advanced diffusion-weighted ASL analysis to generate
        multi-compartment perfusion maps. The analysis uses multiple b-values to
        separate fast-flowing intravascular and slower tissue perfusion components.

        The method fits a bi-exponential diffusion model to estimate:
        - Signal amplitudes and apparent diffusion coefficients for two compartments
        - Water exchange parameters between vascular and tissue compartments
        - Enhanced CBF characterization with compartment specificity

        Note:
            The CBF and ATT maps can be provided before calling this method using
            set_cbf_map() and set_att_map() methods. If not provided, basic maps
            are automatically calculated using the CBFMapping class.

        Warning:
            This method is computationally intensive as it performs voxel-wise
            non-linear fitting without parallel processing. Consider using a brain
            mask to limit processing to relevant tissue areas.

        Args:
            lb (list, optional): Lower bounds for [A1, D1, A2, D2] parameters.
                Defaults to [0.0, 0.0, 0.0, 0.0]. All parameters should be non-negative.
                - A1, A2: Signal amplitudes (relative units)
                - D1, D2: Apparent diffusion coefficients (mm²/s)
            ub (list, optional): Upper bounds for [A1, D1, A2, D2] parameters.
                Defaults to [np.inf, np.inf, np.inf, np.inf].
            par0 (list, optional): Initial guess for [A1, D1, A2, D2] parameters.
                Defaults to [0.5, 0.000005, 0.5, 0.000005].
                - A1, A2: Typical values 0.1-1.0 (relative amplitudes)
                - D1, D2: Typical values 1e-6 to 1e-3 mm²/s
            smoothing (str, optional): Type of spatial smoothing filter to apply.
                Options: None (default, no smoothing), 'gaussian', 'median'.
                Smoothing is applied to all output maps after reconstruction.
            smoothing_params (dict, optional): Parameters for the smoothing filter.
                For 'gaussian': {'sigma': float} (default: 1.0)
                For 'median': {'size': int} (default: 3)

        Returns:
            dict: Dictionary containing diffusion-weighted ASL maps:
                - 'cbf': Basic CBF map in original units (numpy.ndarray)
                - 'cbf_norm': Normalized CBF in mL/100g/min (numpy.ndarray)
                - 'att': Arterial transit time in ms (numpy.ndarray)
                - 'A1': Signal amplitude for compartment 1 (numpy.ndarray)
                - 'D1': Apparent diffusion coefficient for compartment 1 in mm²/s (numpy.ndarray)
                - 'A2': Signal amplitude for compartment 2 (numpy.ndarray)
                - 'D2': Apparent diffusion coefficient for compartment 2 in mm²/s (numpy.ndarray)
                - 'kw': Water exchange parameter (numpy.ndarray)
                All maps are smoothed if smoothing is enabled.

        Examples:
            Basic multi-DW ASL analysis:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import MultiDW_ASLMapping
            >>> import numpy as np
            >>> # Load multi-DW ASL data
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     dw_values=[0, 50, 100, 200],  # b-values in s/mm²
            ...     ld_values=[1.8, 1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8, 3.8]
            ... )
            >>> mdw_mapper = MultiDW_ASLMapping(asl_data)
            >>> # Set brain mask for faster processing (recommended)
            >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
            >>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
            >>> adjusted_brain_mask[0:2, :, :] = 0  # Remove some background slices
            >>> brain_mask.update_image_data(adjusted_brain_mask)
            >>> mdw_mapper.set_brain_mask(brain_mask)
            >>> # Generate all maps (may take several minutes)
            >>> results = mdw_mapper.create_map() # doctest: +SKIP

            Custom parameters for specific tissue analysis:
            >>> # For analyzing fast vs slow perfusion components
            >>> results = mdw_mapper.create_map(
            ...     lb=[0.1, 1e-6, 0.1, 1e-7],      # Minimum realistic values
            ...     ub=[2.0, 1e-3, 2.0, 1e-4],      # Maximum realistic values
            ...     par0=[0.8, 5e-5, 0.3, 1e-5]     # Initial guesses
            ... ) # doctest: +SKIP

        Note:
            Processing time scales with brain mask size. For a full brain analysis,
            expect processing times of 30+ minutes depending on data size and
            hardware capabilities.

        See Also:
            set_cbf_map(): Provide pre-computed CBF map
            set_att_map(): Provide pre-computed ATT map
            CBFMapping: For basic CBF/ATT mapping
        """
        self._basic_maps.set_brain_mask(ImageIO(image_array=self._brain_mask))

        basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
        if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
            # If the CBF/ATT maps are zero (empty), then a new one is created
            print(
                '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
            )   # pragma: no cover
            basic_maps = self._basic_maps.create_map()   # pragma: no cover
            self._cbf_map = basic_maps[
                'cbf'
            ].get_as_numpy()   # pragma: no cover
            self._att_map = basic_maps[
                'att'
            ].get_as_numpy()   # pragma: no cover

        x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
        y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
        z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

        # TODO Fix the reconstruction method when ASL-DWI acquisition works properly
        print('multiDW-ASL processing...')
        for i in range(x_axis):
            for j in range(y_axis):
                for k in range(z_axis):
                    if self._brain_mask[k, j, i] != 0:
                        # Calculates the diffusion components for (A1, D1), (A2, D2)
                        def mod_diff(Xdata, par1, par2, par3, par4):
                            return asl_model_multi_dw(
                                b_values=Xdata,
                                A1=par1,
                                D1=par2,
                                A2=par3,
                                D2=par4,
                            )

                        # M(t,b)/M(t,0)
                        Ydata = (
                            self._asl_data('pcasl')
                            .get_as_numpy()[:, :, k, j, i]
                            .reshape(
                                (
                                    len(self._asl_data.get_ld())
                                    * len(self._asl_data.get_dw()),
                                    1,
                                )
                            )
                            .flatten()
                            / self._asl_data('m0').get_as_numpy()[k, j, i]
                        )

                        try:
                            # Xdata = self._b_values
                            Xdata = self._create_x_data(
                                self._asl_data.get_ld(),
                                self._asl_data.get_pld(),
                                self._asl_data.get_dw(),
                            )

                            par_fit, _ = curve_fit(
                                mod_diff,
                                Xdata[:, 2],
                                Ydata,
                                p0=par0,
                                bounds=(lb, ub),
                            )
                            self._A1[k, j, i] = par_fit[0]
                            self._D1[k, j, i] = par_fit[1]
                            self._A2[k, j, i] = par_fit[2]
                            self._D2[k, j, i] = par_fit[3]
                        except RuntimeError:
                            self._A1[k, j, i] = 0
                            self._D1[k, j, i] = 0
                            self._A2[k, j, i] = 0
                            self._D2[k, j, i] = 0

                        # Calculates the Mc fitting to alpha = kw + T1blood
                        m0_px = self._asl_data('m0').get_as_numpy()[k, j, i]

                        # def mod_2comp(Xdata, par1):
                        #     ...
                        #     # return asl_model_multi_te(
                        #     #     Xdata[:, 0],
                        #     #     Xdata[:, 1],
                        #     #     Xdata[:, 2],
                        #     #     m0_px,
                        #     #     basic_maps['cbf'][k, j, i],
                        #     #     basic_maps['att'][k, j, i],
                        #     #     par1,
                        #     #     self.T2bl,
                        #     #     self.T2gm,
                        #     # )

                        # Ydata = (
                        #     self._asl_data('pcasl')[:, :, k, j, i]
                        #     .reshape(
                        #         (
                        #             len(self._asl_data.get_ld())
                        #             * len(self._asl_data.get_te()),
                        #             1,
                        #         )
                        #     )
                        #     .flatten()
                        # )

                        # try:
                        #     Xdata = self._create_x_data(
                        #         self._asl_data.get_ld(),
                        #         self._asl_data.get_pld(),
                        #         self._asl_data.get_dw(),
                        #     )
                        #     par_fit, _ = curve_fit(
                        #         mod_2comp,
                        #         Xdata,
                        #         Ydata,
                        #         p0=par0,
                        #         bounds=(lb, ub),
                        #     )
                        #     self._kw[k, j, i] = par_fit[0]
                        # except RuntimeError:
                        #     self._kw[k, j, i] = 0.0

        # # Adjusting output image boundaries
        # self._kw = self._adjust_image_limits(self._kw, par0[0])

        # Prepare output maps
        cbf_map_image = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_image.update_image_data(self._cbf_map)

        cbf_map_norm_image = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_norm_image.update_image_data(
            self._cbf_map * (60 * 60 * 1000)
        )  # Convert to mL/100g/min

        att_map_image = ImageIO(self._asl_data('m0').get_image_path())
        att_map_image.update_image_data(self._att_map)

        a1_map_image = ImageIO(self._asl_data('m0').get_image_path())
        a1_map_image.update_image_data(self._A1)

        d1_map_image = ImageIO(self._asl_data('m0').get_image_path())
        d1_map_image.update_image_data(self._D1)

        a2_map_image = ImageIO(self._asl_data('m0').get_image_path())
        a2_map_image.update_image_data(self._A2)

        d2_map_image = ImageIO(self._asl_data('m0').get_image_path())
        d2_map_image.update_image_data(self._D2)

        kw_map_image = ImageIO(self._asl_data('m0').get_image_path())
        kw_map_image.update_image_data(self._kw)

        # Create output maps dictionary
        output_maps = {
            'cbf': cbf_map_image,
            'cbf_norm': cbf_map_norm_image,
            'att': att_map_image,
            'a1': a1_map_image,
            'd1': d1_map_image,
            'a2': a2_map_image,
            'd2': d2_map_image,
            'kw': kw_map_image,
        }

        # Apply smoothing if requested
        return _apply_smoothing_to_maps(
            output_maps, smoothing, smoothing_params
        )

    def _create_x_data(self, ld, pld, dw):
        # array for the x values, assuming an arbitrary size based on the PLD
        # and TE vector size
        Xdata = np.zeros((len(pld) * len(dw), 3))

        count = 0
        for i in range(len(pld)):
            for j in range(len(dw)):
                Xdata[count] = [ld[i], pld[i], dw[j]]
                count += 1

        return Xdata

__init__(asl_data)

Multi-Diffusion-Weighted ASL mapping constructor for advanced perfusion analysis.

MultiDW_ASLMapping enables sophisticated ASL analysis by incorporating multiple diffusion weightings (b-values) to separate intravascular and tissue compartments. This approach provides enhanced characterization of perfusion and can help differentiate between different vascular compartments.

The class implements diffusion-weighted ASL analysis that can distinguish: - Fast-flowing blood (intravascular component) - Slow-flowing blood and tissue perfusion - Apparent diffusion coefficients for each compartment - Water exchange parameters between compartments

Notes

The ASLData object must contain dw_values - a list of diffusion b-values used during ASL acquisition. These b-values are essential for the multi-compartment diffusion model fitting.

Examples:

Basic multi-DW ASL mapping setup:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import MultiDW_ASLMapping
>>> # Create ASL data with diffusion weighting
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mdw.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     dw_values=[0, 50, 100, 200],    # b-values in s/mm²
...     ld_values=[1.8, 1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8, 3.8]
... )
>>> mdw_mapper = MultiDW_ASLMapping(asl_data)

Access diffusion-related maps (after processing):

>>> # These maps will be populated after create_map() is called
>>> # A1: Signal amplitude for compartment 1
>>> # D1: Apparent diffusion coefficient for compartment 1
>>> # A2: Signal amplitude for compartment 2
>>> # D2: Apparent diffusion coefficient for compartment 2
>>> # kw: Water exchange parameter

Parameters:

Name Type Description Default
asl_data ASLData

The ASL data object containing multi-DW acquisition. Must include dw_values (b-values), ld_values, and pld_values.

required

Raises:

Type Description
ValueError

If ASLData object lacks required DW values for diffusion-weighted analysis.

See Also

CBFMapping: For basic CBF/ATT mapping without diffusion weighting MultiTE_ASLMapping: For multi-echo ASL analysis

Source code in asltk/reconstruction/multi_dw_mapping.py
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
def __init__(self, asl_data: ASLData):
    """Multi-Diffusion-Weighted ASL mapping constructor for advanced perfusion analysis.

    MultiDW_ASLMapping enables sophisticated ASL analysis by incorporating multiple
    diffusion weightings (b-values) to separate intravascular and tissue
    compartments. This approach provides enhanced characterization of perfusion
    and can help differentiate between different vascular compartments.

    The class implements diffusion-weighted ASL analysis that can distinguish:
    - Fast-flowing blood (intravascular component)
    - Slow-flowing blood and tissue perfusion
    - Apparent diffusion coefficients for each compartment
    - Water exchange parameters between compartments

    Notes:
        The ASLData object must contain `dw_values` - a list of diffusion
        b-values used during ASL acquisition. These b-values are essential
        for the multi-compartment diffusion model fitting.

    Examples:
        Basic multi-DW ASL mapping setup:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import MultiDW_ASLMapping
        >>> # Create ASL data with diffusion weighting
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     dw_values=[0, 50, 100, 200],    # b-values in s/mm²
        ...     ld_values=[1.8, 1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8, 3.8]
        ... )
        >>> mdw_mapper = MultiDW_ASLMapping(asl_data)

        Access diffusion-related maps (after processing):
        >>> # These maps will be populated after create_map() is called
        >>> # A1: Signal amplitude for compartment 1
        >>> # D1: Apparent diffusion coefficient for compartment 1
        >>> # A2: Signal amplitude for compartment 2
        >>> # D2: Apparent diffusion coefficient for compartment 2
        >>> # kw: Water exchange parameter

    Args:
        asl_data (ASLData): The ASL data object containing multi-DW acquisition.
            Must include dw_values (b-values), ld_values, and pld_values.

    Raises:
        ValueError: If ASLData object lacks required DW values for
            diffusion-weighted analysis.

    See Also:
        CBFMapping: For basic CBF/ATT mapping without diffusion weighting
        MultiTE_ASLMapping: For multi-echo ASL analysis
    """
    super().__init__()
    self._asl_data = asl_data
    self._basic_maps = CBFMapping(asl_data)
    if self._asl_data.get_dw() is None:
        raise ValueError(
            'ASLData is incomplete. MultiDW_ASLMapping need a list of DW values.'
        )

    self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
    self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

    self._b_values = self._asl_data.get_dw()
    # self._A1 = np.zeros(tuple([len(self._b_values)]) + self._asl_data('m0').shape)
    self._A1 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    # self._D1 = np.zeros(tuple([1]) +self._asl_data('m0').shape)
    self._D1 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._A2 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    # self._A2 = np.zeros(tuple([len(self._b_values)])  + self._asl_data('m0').shape)
    # self._D2 = np.zeros(tuple([1]) +self._asl_data('m0').shape)
    self._D2 = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._kw = np.zeros(self._asl_data('m0').get_as_numpy().shape)

create_map(lb=[0.0, 0.0, 0.0, 0.0], ub=[np.inf, np.inf, np.inf, np.inf], par0=[0.5, 5e-06, 0.5, 5e-06], smoothing=None, smoothing_params=None)

Create multi-diffusion-weighted ASL maps for compartment analysis.

This method performs advanced diffusion-weighted ASL analysis to generate multi-compartment perfusion maps. The analysis uses multiple b-values to separate fast-flowing intravascular and slower tissue perfusion components.

The method fits a bi-exponential diffusion model to estimate: - Signal amplitudes and apparent diffusion coefficients for two compartments - Water exchange parameters between vascular and tissue compartments - Enhanced CBF characterization with compartment specificity

Note

The CBF and ATT maps can be provided before calling this method using set_cbf_map() and set_att_map() methods. If not provided, basic maps are automatically calculated using the CBFMapping class.

Warning

This method is computationally intensive as it performs voxel-wise non-linear fitting without parallel processing. Consider using a brain mask to limit processing to relevant tissue areas.

Parameters:

Name Type Description Default
lb list

Lower bounds for [A1, D1, A2, D2] parameters. Defaults to [0.0, 0.0, 0.0, 0.0]. All parameters should be non-negative. - A1, A2: Signal amplitudes (relative units) - D1, D2: Apparent diffusion coefficients (mm²/s)

[0.0, 0.0, 0.0, 0.0]
ub list

Upper bounds for [A1, D1, A2, D2] parameters. Defaults to [np.inf, np.inf, np.inf, np.inf].

[inf, inf, inf, inf]
par0 list

Initial guess for [A1, D1, A2, D2] parameters. Defaults to [0.5, 0.000005, 0.5, 0.000005]. - A1, A2: Typical values 0.1-1.0 (relative amplitudes) - D1, D2: Typical values 1e-6 to 1e-3 mm²/s

[0.5, 5e-06, 0.5, 5e-06]
smoothing str

Type of spatial smoothing filter to apply. Options: None (default, no smoothing), 'gaussian', 'median'. Smoothing is applied to all output maps after reconstruction.

None
smoothing_params dict

Parameters for the smoothing filter. For 'gaussian': {'sigma': float} (default: 1.0) For 'median': {'size': int} (default: 3)

None

Returns:

Name Type Description
dict

Dictionary containing diffusion-weighted ASL maps: - 'cbf': Basic CBF map in original units (numpy.ndarray) - 'cbf_norm': Normalized CBF in mL/100g/min (numpy.ndarray) - 'att': Arterial transit time in ms (numpy.ndarray) - 'A1': Signal amplitude for compartment 1 (numpy.ndarray) - 'D1': Apparent diffusion coefficient for compartment 1 in mm²/s (numpy.ndarray) - 'A2': Signal amplitude for compartment 2 (numpy.ndarray) - 'D2': Apparent diffusion coefficient for compartment 2 in mm²/s (numpy.ndarray) - 'kw': Water exchange parameter (numpy.ndarray) All maps are smoothed if smoothing is enabled.

Examples:

Basic multi-DW ASL analysis:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import MultiDW_ASLMapping
>>> import numpy as np
>>> # Load multi-DW ASL data
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mdw.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     dw_values=[0, 50, 100, 200],  # b-values in s/mm²
...     ld_values=[1.8, 1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8, 3.8]
... )
>>> mdw_mapper = MultiDW_ASLMapping(asl_data)
>>> # Set brain mask for faster processing (recommended)
>>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
>>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
>>> adjusted_brain_mask[0:2, :, :] = 0  # Remove some background slices
>>> brain_mask.update_image_data(adjusted_brain_mask)
>>> mdw_mapper.set_brain_mask(brain_mask)
>>> # Generate all maps (may take several minutes)
>>> results = mdw_mapper.create_map()

Custom parameters for specific tissue analysis:

>>> # For analyzing fast vs slow perfusion components
>>> results = mdw_mapper.create_map(
...     lb=[0.1, 1e-6, 0.1, 1e-7],      # Minimum realistic values
...     ub=[2.0, 1e-3, 2.0, 1e-4],      # Maximum realistic values
...     par0=[0.8, 5e-5, 0.3, 1e-5]     # Initial guesses
... )
Note

Processing time scales with brain mask size. For a full brain analysis, expect processing times of 30+ minutes depending on data size and hardware capabilities.

See Also

set_cbf_map(): Provide pre-computed CBF map set_att_map(): Provide pre-computed ATT map CBFMapping: For basic CBF/ATT mapping

Source code in asltk/reconstruction/multi_dw_mapping.py
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
def create_map(
    self,
    lb: list = [0.0, 0.0, 0.0, 0.0],
    ub: list = [np.inf, np.inf, np.inf, np.inf],
    par0: list = [0.5, 0.000005, 0.5, 0.000005],
    smoothing=None,
    smoothing_params=None,
):
    """Create multi-diffusion-weighted ASL maps for compartment analysis.

    This method performs advanced diffusion-weighted ASL analysis to generate
    multi-compartment perfusion maps. The analysis uses multiple b-values to
    separate fast-flowing intravascular and slower tissue perfusion components.

    The method fits a bi-exponential diffusion model to estimate:
    - Signal amplitudes and apparent diffusion coefficients for two compartments
    - Water exchange parameters between vascular and tissue compartments
    - Enhanced CBF characterization with compartment specificity

    Note:
        The CBF and ATT maps can be provided before calling this method using
        set_cbf_map() and set_att_map() methods. If not provided, basic maps
        are automatically calculated using the CBFMapping class.

    Warning:
        This method is computationally intensive as it performs voxel-wise
        non-linear fitting without parallel processing. Consider using a brain
        mask to limit processing to relevant tissue areas.

    Args:
        lb (list, optional): Lower bounds for [A1, D1, A2, D2] parameters.
            Defaults to [0.0, 0.0, 0.0, 0.0]. All parameters should be non-negative.
            - A1, A2: Signal amplitudes (relative units)
            - D1, D2: Apparent diffusion coefficients (mm²/s)
        ub (list, optional): Upper bounds for [A1, D1, A2, D2] parameters.
            Defaults to [np.inf, np.inf, np.inf, np.inf].
        par0 (list, optional): Initial guess for [A1, D1, A2, D2] parameters.
            Defaults to [0.5, 0.000005, 0.5, 0.000005].
            - A1, A2: Typical values 0.1-1.0 (relative amplitudes)
            - D1, D2: Typical values 1e-6 to 1e-3 mm²/s
        smoothing (str, optional): Type of spatial smoothing filter to apply.
            Options: None (default, no smoothing), 'gaussian', 'median'.
            Smoothing is applied to all output maps after reconstruction.
        smoothing_params (dict, optional): Parameters for the smoothing filter.
            For 'gaussian': {'sigma': float} (default: 1.0)
            For 'median': {'size': int} (default: 3)

    Returns:
        dict: Dictionary containing diffusion-weighted ASL maps:
            - 'cbf': Basic CBF map in original units (numpy.ndarray)
            - 'cbf_norm': Normalized CBF in mL/100g/min (numpy.ndarray)
            - 'att': Arterial transit time in ms (numpy.ndarray)
            - 'A1': Signal amplitude for compartment 1 (numpy.ndarray)
            - 'D1': Apparent diffusion coefficient for compartment 1 in mm²/s (numpy.ndarray)
            - 'A2': Signal amplitude for compartment 2 (numpy.ndarray)
            - 'D2': Apparent diffusion coefficient for compartment 2 in mm²/s (numpy.ndarray)
            - 'kw': Water exchange parameter (numpy.ndarray)
            All maps are smoothed if smoothing is enabled.

    Examples:
        Basic multi-DW ASL analysis:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import MultiDW_ASLMapping
        >>> import numpy as np
        >>> # Load multi-DW ASL data
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     dw_values=[0, 50, 100, 200],  # b-values in s/mm²
        ...     ld_values=[1.8, 1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8, 3.8]
        ... )
        >>> mdw_mapper = MultiDW_ASLMapping(asl_data)
        >>> # Set brain mask for faster processing (recommended)
        >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
        >>> adjusted_brain_mask = brain_mask.get_as_numpy().copy()
        >>> adjusted_brain_mask[0:2, :, :] = 0  # Remove some background slices
        >>> brain_mask.update_image_data(adjusted_brain_mask)
        >>> mdw_mapper.set_brain_mask(brain_mask)
        >>> # Generate all maps (may take several minutes)
        >>> results = mdw_mapper.create_map() # doctest: +SKIP

        Custom parameters for specific tissue analysis:
        >>> # For analyzing fast vs slow perfusion components
        >>> results = mdw_mapper.create_map(
        ...     lb=[0.1, 1e-6, 0.1, 1e-7],      # Minimum realistic values
        ...     ub=[2.0, 1e-3, 2.0, 1e-4],      # Maximum realistic values
        ...     par0=[0.8, 5e-5, 0.3, 1e-5]     # Initial guesses
        ... ) # doctest: +SKIP

    Note:
        Processing time scales with brain mask size. For a full brain analysis,
        expect processing times of 30+ minutes depending on data size and
        hardware capabilities.

    See Also:
        set_cbf_map(): Provide pre-computed CBF map
        set_att_map(): Provide pre-computed ATT map
        CBFMapping: For basic CBF/ATT mapping
    """
    self._basic_maps.set_brain_mask(ImageIO(image_array=self._brain_mask))

    basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
    if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
        # If the CBF/ATT maps are zero (empty), then a new one is created
        print(
            '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
        )   # pragma: no cover
        basic_maps = self._basic_maps.create_map()   # pragma: no cover
        self._cbf_map = basic_maps[
            'cbf'
        ].get_as_numpy()   # pragma: no cover
        self._att_map = basic_maps[
            'att'
        ].get_as_numpy()   # pragma: no cover

    x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
    y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
    z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

    # TODO Fix the reconstruction method when ASL-DWI acquisition works properly
    print('multiDW-ASL processing...')
    for i in range(x_axis):
        for j in range(y_axis):
            for k in range(z_axis):
                if self._brain_mask[k, j, i] != 0:
                    # Calculates the diffusion components for (A1, D1), (A2, D2)
                    def mod_diff(Xdata, par1, par2, par3, par4):
                        return asl_model_multi_dw(
                            b_values=Xdata,
                            A1=par1,
                            D1=par2,
                            A2=par3,
                            D2=par4,
                        )

                    # M(t,b)/M(t,0)
                    Ydata = (
                        self._asl_data('pcasl')
                        .get_as_numpy()[:, :, k, j, i]
                        .reshape(
                            (
                                len(self._asl_data.get_ld())
                                * len(self._asl_data.get_dw()),
                                1,
                            )
                        )
                        .flatten()
                        / self._asl_data('m0').get_as_numpy()[k, j, i]
                    )

                    try:
                        # Xdata = self._b_values
                        Xdata = self._create_x_data(
                            self._asl_data.get_ld(),
                            self._asl_data.get_pld(),
                            self._asl_data.get_dw(),
                        )

                        par_fit, _ = curve_fit(
                            mod_diff,
                            Xdata[:, 2],
                            Ydata,
                            p0=par0,
                            bounds=(lb, ub),
                        )
                        self._A1[k, j, i] = par_fit[0]
                        self._D1[k, j, i] = par_fit[1]
                        self._A2[k, j, i] = par_fit[2]
                        self._D2[k, j, i] = par_fit[3]
                    except RuntimeError:
                        self._A1[k, j, i] = 0
                        self._D1[k, j, i] = 0
                        self._A2[k, j, i] = 0
                        self._D2[k, j, i] = 0

                    # Calculates the Mc fitting to alpha = kw + T1blood
                    m0_px = self._asl_data('m0').get_as_numpy()[k, j, i]

                    # def mod_2comp(Xdata, par1):
                    #     ...
                    #     # return asl_model_multi_te(
                    #     #     Xdata[:, 0],
                    #     #     Xdata[:, 1],
                    #     #     Xdata[:, 2],
                    #     #     m0_px,
                    #     #     basic_maps['cbf'][k, j, i],
                    #     #     basic_maps['att'][k, j, i],
                    #     #     par1,
                    #     #     self.T2bl,
                    #     #     self.T2gm,
                    #     # )

                    # Ydata = (
                    #     self._asl_data('pcasl')[:, :, k, j, i]
                    #     .reshape(
                    #         (
                    #             len(self._asl_data.get_ld())
                    #             * len(self._asl_data.get_te()),
                    #             1,
                    #         )
                    #     )
                    #     .flatten()
                    # )

                    # try:
                    #     Xdata = self._create_x_data(
                    #         self._asl_data.get_ld(),
                    #         self._asl_data.get_pld(),
                    #         self._asl_data.get_dw(),
                    #     )
                    #     par_fit, _ = curve_fit(
                    #         mod_2comp,
                    #         Xdata,
                    #         Ydata,
                    #         p0=par0,
                    #         bounds=(lb, ub),
                    #     )
                    #     self._kw[k, j, i] = par_fit[0]
                    # except RuntimeError:
                    #     self._kw[k, j, i] = 0.0

    # # Adjusting output image boundaries
    # self._kw = self._adjust_image_limits(self._kw, par0[0])

    # Prepare output maps
    cbf_map_image = ImageIO(self._asl_data('m0').get_image_path())
    cbf_map_image.update_image_data(self._cbf_map)

    cbf_map_norm_image = ImageIO(self._asl_data('m0').get_image_path())
    cbf_map_norm_image.update_image_data(
        self._cbf_map * (60 * 60 * 1000)
    )  # Convert to mL/100g/min

    att_map_image = ImageIO(self._asl_data('m0').get_image_path())
    att_map_image.update_image_data(self._att_map)

    a1_map_image = ImageIO(self._asl_data('m0').get_image_path())
    a1_map_image.update_image_data(self._A1)

    d1_map_image = ImageIO(self._asl_data('m0').get_image_path())
    d1_map_image.update_image_data(self._D1)

    a2_map_image = ImageIO(self._asl_data('m0').get_image_path())
    a2_map_image.update_image_data(self._A2)

    d2_map_image = ImageIO(self._asl_data('m0').get_image_path())
    d2_map_image.update_image_data(self._D2)

    kw_map_image = ImageIO(self._asl_data('m0').get_image_path())
    kw_map_image.update_image_data(self._kw)

    # Create output maps dictionary
    output_maps = {
        'cbf': cbf_map_image,
        'cbf_norm': cbf_map_norm_image,
        'att': att_map_image,
        'a1': a1_map_image,
        'd1': d1_map_image,
        'a2': a2_map_image,
        'd2': d2_map_image,
        'kw': kw_map_image,
    }

    # Apply smoothing if requested
    return _apply_smoothing_to_maps(
        output_maps, smoothing, smoothing_params
    )

get_att_map()

Get the ATT map storaged at the MultiDW_ASLMapping object

Returns:

Type Description
ndarray

description

Source code in asltk/reconstruction/multi_dw_mapping.py
201
202
203
204
205
206
207
def get_att_map(self):
    """Get the ATT map storaged at the MultiDW_ASLMapping object

    Returns:
        (np.ndarray): _description_
    """
    return self._att_map

get_brain_mask()

Get the brain mask image

Returns:

Type Description
ndarray

The brain mask image

Source code in asltk/reconstruction/multi_dw_mapping.py
163
164
165
166
167
168
169
def get_brain_mask(self):
    """Get the brain mask image

    Returns:
        (np.ndarray): The brain mask image
    """
    return self._brain_mask

get_cbf_map()

Get the CBF map storaged at the MultiDW_ASLMapping object

Returns:

Type Description
ndarray

The CBF map that is storaged in the

ImageIO

MultiDW_ASLMapping object

Source code in asltk/reconstruction/multi_dw_mapping.py
184
185
186
187
188
189
190
191
def get_cbf_map(self) -> ImageIO:
    """Get the CBF map storaged at the MultiDW_ASLMapping object

    Returns:
        (np.ndarray): The CBF map that is storaged in the
        MultiDW_ASLMapping object
    """
    return self._cbf_map

set_att_map(att_map)

Set the ATT map to the MultiDW_ASLMapping object.

Parameters:

Name Type Description Default
att_map ndarray

The ATT map that is set in the MultiDW_ASLMapping object

required
Source code in asltk/reconstruction/multi_dw_mapping.py
193
194
195
196
197
198
199
def set_att_map(self, att_map: ImageIO):
    """Set the ATT map to the MultiDW_ASLMapping object.

    Args:
        att_map (np.ndarray): The ATT map that is set in the MultiDW_ASLMapping object
    """
    self._att_map = att_map.get_as_numpy()

set_brain_mask(brain_mask, label=1)

Set brain mask for MultiDW-ASL processing (strongly recommended).

A brain mask is especially important for multi-diffusion-weighted ASL processing as it significantly reduces computation time by limiting the intensive voxel-wise fitting to brain tissue regions only.

Without a brain mask, processing time can be prohibitively long (hours) for whole-volume analysis. A proper brain mask can reduce processing time by 5-10x while maintaining analysis quality.

Parameters:

Name Type Description Default
brain_mask ndarray

The image representing the brain mask. Must match the spatial dimensions of the M0 image.

required
label int

The label value used to define brain tissue. Defaults to 1. Voxels with this value will be processed.

1

Examples:

Set a brain mask for efficient processing:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import MultiDW_ASLMapping
>>> import numpy as np
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mdw.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     dw_values=[0, 50, 100], ld_values=[1.8]*3, pld_values=[1.8]*3
... )
>>> mdw_mapper = MultiDW_ASLMapping(asl_data)
>>> # Create conservative brain mask (center region only)
>>> mask_shape = asl_data('m0').get_as_numpy().shape
>>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
>>> adjusted_brain_mask = brain_mask.get_as_numpy()
>>> adjusted_brain_mask[1:4, 5:30, 5:30] = 1  # Conservative brain region
>>> brain_mask.update_image_data(adjusted_brain_mask)
>>> mdw_mapper.set_brain_mask(brain_mask)
Note

For multi-DW ASL, consider using a more conservative (smaller) brain mask initially to test parameters and processing time, then expand to full brain analysis once satisfied with results.

Source code in asltk/reconstruction/multi_dw_mapping.py
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
    """Set brain mask for MultiDW-ASL processing (strongly recommended).

    A brain mask is especially important for multi-diffusion-weighted ASL
    processing as it significantly reduces computation time by limiting
    the intensive voxel-wise fitting to brain tissue regions only.

    Without a brain mask, processing time can be prohibitively long (hours)
    for whole-volume analysis. A proper brain mask can reduce processing
    time by 5-10x while maintaining analysis quality.

    Args:
        brain_mask (np.ndarray): The image representing the brain mask.
            Must match the spatial dimensions of the M0 image.
        label (int, optional): The label value used to define brain tissue.
            Defaults to 1. Voxels with this value will be processed.

    Examples:
        Set a brain mask for efficient processing:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import MultiDW_ASLMapping
        >>> import numpy as np
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mdw.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     dw_values=[0, 50, 100], ld_values=[1.8]*3, pld_values=[1.8]*3
        ... )
        >>> mdw_mapper = MultiDW_ASLMapping(asl_data)
        >>> # Create conservative brain mask (center region only)
        >>> mask_shape = asl_data('m0').get_as_numpy().shape
        >>> brain_mask = ImageIO(image_array=np.zeros(mask_shape))
        >>> adjusted_brain_mask = brain_mask.get_as_numpy()
        >>> adjusted_brain_mask[1:4, 5:30, 5:30] = 1  # Conservative brain region
        >>> brain_mask.update_image_data(adjusted_brain_mask)
        >>> mdw_mapper.set_brain_mask(brain_mask)

    Note:
        For multi-DW ASL, consider using a more conservative (smaller) brain
        mask initially to test parameters and processing time, then expand
        to full brain analysis once satisfied with results.
    """
    if not isinstance(brain_mask, ImageIO):
        raise TypeError(
            'Brain mask must be an instance of ImageIO. '
            'Use ImageIO to load or create the mask.'
        )

    _check_mask_values(
        brain_mask, label, self._asl_data('m0').get_as_numpy().shape
    )

    binary_mask = (brain_mask.get_as_numpy() == label).astype(
        np.uint8
    ) * label
    self._brain_mask = binary_mask

set_cbf_map(cbf_map)

Set the CBF map to the MultiDW_ASLMapping object.

Note

The CBF maps must have the original scale in order to calculate the T1blGM map correclty. Hence, if the CBF map was made using CBFMapping class, one can use the 'cbf' output.

Parameters:

Name Type Description Default
cbf_map ndarray

The CBF map that is set in the MultiDW_ASLMapping object

required
Source code in asltk/reconstruction/multi_dw_mapping.py
171
172
173
174
175
176
177
178
179
180
181
182
def set_cbf_map(self, cbf_map: ImageIO):
    """Set the CBF map to the MultiDW_ASLMapping object.

    Note:
        The CBF maps must have the original scale in order to calculate the
        T1blGM map correclty. Hence, if the CBF map was made using
        CBFMapping class, one can use the 'cbf' output.

    Args:
        cbf_map (np.ndarray): The CBF map that is set in the MultiDW_ASLMapping object
    """
    self._cbf_map = cbf_map.get_as_numpy()

MultiTE_ASLMapping

Bases: MRIParameters

Source code in asltk/reconstruction/multi_te_mapping.py
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
class MultiTE_ASLMapping(MRIParameters):
    def __init__(self, asl_data: ASLData) -> None:
        """Multi-Echo ASL mapping constructor for T1 tissue relaxometry.

        MultiTE_ASLMapping enables advanced ASL analysis by incorporating multiple
        echo times (TE) to estimate tissue-specific T1 relaxation times. This
        provides better characterization of blood vs. tissue compartments and
        improved CBF quantification.

        The class requires ASL data acquired with multiple echo times and performs:
        - Basic CBF and ATT mapping (via CBFMapping)
        - T1 relaxometry for blood-grey matter differentiation
        - Multi-TE model fitting for enhanced tissue characterization

        Notes:
            The ASLData object must contain `te_values` - a list of echo times
            used during ASL acquisition. These TE values are critical for the
            multi-echo model fitting and T1 estimation.

        Examples:
            Basic multi-TE ASL mapping setup:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import MultiTE_ASLMapping
            >>> # Create ASL data with multi-TE parameters
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
            ...     ld_values=[1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8]
            ... )
            >>> mte_mapper = MultiTE_ASLMapping(asl_data)
            >>> # Access default MRI parameters
            >>> mte_mapper.get_constant('T1csf')
            1400.0

            Custom MRI parameters for specific field strength:
            >>> # Adjust T1 values for 3T scanner
            >>> mte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
            >>> mte_mapper.get_constant('T1csf')
            1600.0
            >>> # Verify default parameters unchanged for other objects
            >>> from asltk.mri_parameters import MRIParameters
            >>> default_param = MRIParameters()
            >>> default_param.get_constant('T1csf')
            1400.0

        Args:
            asl_data (ASLData): The ASL data object containing multi-TE acquisition.
                Must include te_values, ld_values, and pld_values.

        Raises:
            ValueError: If ASLData object lacks required TE values for multi-echo analysis.

        See Also:
            CBFMapping: For basic CBF/ATT mapping without multi-TE analysis
            MultiDW_ASLMapping: For diffusion-weighted ASL analysis
        """
        super().__init__()
        self._asl_data = asl_data
        self._basic_maps = CBFMapping(asl_data)
        if self._asl_data.get_te() is None:
            raise ValueError(
                'ASLData is incomplete. MultiTE_ASLMapping need a list of TE values.'
            )

        self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
        self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._t1blgm_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

    def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
        """Defines whether a brain a mask is applied to the CBFMapping
        calculation

        A image mask is simply an image that defines the voxels where the ASL
        calculation should be made. Basically any integer value can be used as
        proper label mask.

        A most common approach is to use a binary image (zeros for background
        and 1 for the brain tissues). Anyway, the default behavior of the
        method can transform a integer-pixel values image to a binary mask with
        the `label` parameter provided by the user

        Args:
            brain_mask (np.ndarray): The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.
        """
        if not isinstance(brain_mask, ImageIO):
            raise TypeError(
                'The brain_mask parameter must be an instance of ImageIO.'
            )

        _check_mask_values(
            brain_mask, label, self._asl_data('m0').get_as_numpy().shape
        )

        binary_mask = (brain_mask.get_as_numpy() == label).astype(
            np.uint8
        ) * label
        self._brain_mask = binary_mask

    def get_brain_mask(self):
        """Get the brain mask image

        Returns:
            (np.ndarray): The brain mask image
        """
        return self._brain_mask

    def set_cbf_map(self, cbf_map: ImageIO):
        """Set the CBF map to the MultiTE_ASLMapping object.

        Note:
            The CBF maps must have the original scale in order to calculate the
            T1blGM map correclty. Hence, if the CBF map was made using
            CBFMapping class, one can use the 'cbf' output.

        Args:
            cbf_map (ImageIO): The CBF map that is set in the MultiTE_ASLMapping object
        """
        self._cbf_map = cbf_map.get_as_numpy()

    def get_cbf_map(self) -> np.ndarray:
        """Get the CBF map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The CBF map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._cbf_map

    def set_att_map(self, att_map: ImageIO):
        """Set the ATT map to the MultiTE_ASLMapping object.

        Args:
            att_map (ImageIO): The ATT map that is set in the MultiTE_ASLMapping object
        """
        self._att_map = att_map.get_as_numpy()

    def get_att_map(self):
        """Get the ATT map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The ATT map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._att_map

    def get_t1blgm_map(self):
        """Get the T1blGM map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The T1blGM map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._t1blgm_map

    def create_map(
        self,
        ub: list = [np.inf],
        lb: list = [0.0],
        par0: list = [400],
        cores=cpu_count(),
        smoothing=None,
        smoothing_params=None,
        suppress_warnings=True,
    ):
        """Create multi-TE ASL maps including T1 blood-gray matter exchange (T1blGM).

        This method performs advanced multi-echo ASL analysis to generate tissue-specific
        T1 relaxation maps that characterize blood-to-gray matter water exchange. The
        analysis uses multiple echo times to separate blood and tissue signal contributions.

        The method implements the multi-compartment TE ASL model described in:
        "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
        blood-to-CSF water transport in humans", NeuroImage, 2022.
        doi: 10.1016/j.neuroimage.2021.118755

        Note:
            The CBF and ATT maps can be provided before calling this method,
            using the set_cbf_map() and set_att_map() methods. If not provided,
            basic CBF/ATT maps are automatically calculated using the CBFMapping class.

        Note:
            The CBF map must be in original scale (not normalized) to perform the
            correct multiTE-ASL model fitting. Use the 'cbf' output from CBFMapping,
            not the 'cbf_norm' version.

        The method assumes the T1blGM values are well-characterized by the initial
        guess parameter. Results are filtered to include only positive values and
        values below 4 times the initial guess to remove unrealistic outliers.

        Note:
            Consider applying spatial smoothing to the output T1blGM map to improve
            SNR. The create_map() method does not apply filtering by default to
            preserve spatial resolution.

        Args:
            ub (list, optional): Upper bounds for T1blGM fitting. Defaults to [np.inf].
                Typically 800-1200 ms for healthy gray matter at 3T.
            lb (list, optional): Lower bounds for T1blGM fitting. Defaults to [0.0].
                Should be positive for realistic T1 values.
            par0 (list, optional): Initial guess for T1blGM in milliseconds.
                Defaults to [400]. Good starting values: 300-500 ms.
            cores (int, optional): Number of CPU threads for parallel processing.
                Defaults to all available cores.
            smoothing (str, optional): Type of spatial smoothing filter to apply.
                Options: None (default, no smoothing), 'gaussian', 'median'.
                Smoothing is applied to all output maps after reconstruction.
            smoothing_params (dict, optional): Parameters for the smoothing filter.
                For 'gaussian': {'sigma': float} (default: 1.0)
                For 'median': {'size': int} (default: 3)
            suppress_warnings (bool, optional): Whether to suppress warnings during
                processing. Defaults to True.

        Returns:
            dict: Dictionary containing:
                - 'cbf': Basic CBF map in original units (ImageIO)
                - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO)
                - 'att': Arterial transit time in ms (ImageIO)
                - 't1blgm': T1 blood-gray matter exchange time in ms (ImageIO)
                All maps are smoothed if smoothing is enabled.

        Examples:
            Basic multi-TE ASL analysis:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import MultiTE_ASLMapping
            >>> from asltk.utils.io import ImageIO
            >>> import numpy as np
            >>> # Load multi-TE ASL data
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
            ...     ld_values=[1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8]
            ... )
            >>> mte_mapper = MultiTE_ASLMapping(asl_data)
            >>> # Set brain mask for faster processing
            >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
            >>> mte_mapper.set_brain_mask(brain_mask)
            >>> # Generate all maps
            >>> results = mte_mapper.create_map() # doctest: +SKIP


            Custom parameters for specific analysis:
            >>> # For expected shorter T1blGM values (faster exchange)
            >>> results = mte_mapper.create_map(
            ...     ub=[600.0],        # Lower upper bound
            ...     lb=[50.0],         # Minimum realistic T1
            ...     par0=[300.0]       # Lower initial guess
            ... ) # doctest: +SKIP

            Apply spatial smoothing to improve SNR:
            >>> # Gaussian smoothing with default sigma=1.0
            >>> results_smooth = mte_mapper.create_map(
            ...     smoothing='gaussian'
            ... ) # doctest: +SKIP
            >>> # Custom smoothing parameters
            >>> results_custom = mte_mapper.create_map(
            ...     smoothing='gaussian',
            ...     smoothing_params={'sigma': 1.5}
            ... ) # doctest: +SKIP
            >>> # Median filtering for edge preservation
            >>> results_median = mte_mapper.create_map(
            ...     smoothing='median',
            ...     smoothing_params={'size': 5}
            ... ) # doctest: +SKIP

        Raises:
            ValueError: If cores parameter is invalid or required data is missing.

        See Also:
            set_cbf_map(): Provide pre-computed CBF map
            set_att_map(): Provide pre-computed ATT map
            CBFMapping: For basic CBF/ATT mapping
        """
        # Use context manager to suppress warnings if requested
        with warnings.catch_warnings():
            if suppress_warnings:
                # Filter common warnings that might appear during fitting and processing
                warnings.filterwarnings('ignore', category=RuntimeWarning)
                warnings.filterwarnings('ignore', category=UserWarning)
                warnings.filterwarnings(
                    'ignore', category=np.VisibleDeprecationWarning
                )

            self._basic_maps.set_brain_mask(
                ImageIO(image_array=self._brain_mask)
            )

            basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
            if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
                # If the CBF/ATT maps are zero (empty), then a new one is created
                print(
                    '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
                )
                basic_maps = self._basic_maps.create_map()
                self._cbf_map = basic_maps['cbf'].get_as_numpy()
                self._att_map = basic_maps['att'].get_as_numpy()

            global asl_data, brain_mask, cbf_map, att_map, t2bl, t2gm
            asl_data = self._asl_data
            brain_mask = self._brain_mask
            cbf_map = self._cbf_map
            att_map = self._att_map
            ld_arr = self._asl_data.get_ld()
            pld_arr = self._asl_data.get_pld()
            te_arr = self._asl_data.get_te()
            t2bl = self.T2bl
            t2gm = self.T2gm

            x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
            y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
            z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

            tblgm_map_shared = Array('d', z_axis * y_axis * x_axis, lock=False)

            with Pool(
                processes=cores,
                initializer=_multite_init_globals,
                initargs=(
                    cbf_map,
                    att_map,
                    brain_mask,
                    asl_data,
                    ld_arr,
                    pld_arr,
                    te_arr,
                    tblgm_map_shared,
                    t2bl,
                    t2gm,
                ),
            ) as pool:
                with Progress() as progress:
                    task = progress.add_task(
                        'multiTE-ASL processing...', total=x_axis
                    )
                    results = [
                        pool.apply_async(
                            _tblgm_multite_process_slice,
                            args=(i, x_axis, y_axis, z_axis, par0, lb, ub),
                            callback=lambda _: progress.update(
                                task, advance=1
                            ),
                        )
                        for i in range(x_axis)
                    ]
                    for result in results:
                        result.wait()

            self._t1blgm_map = np.frombuffer(tblgm_map_shared).reshape(
                z_axis, y_axis, x_axis
            )

            # Adjusting output image boundaries
            self._t1blgm_map = self._adjust_image_limits(
                self._t1blgm_map, par0[0]
            )

            # Prepare output maps
            base_volume = ImageIO(self._asl_data('m0').get_image_path())
            cbf_map_image = clone_image(base_volume)
            cbf_map_image.update_image_data(
                self._cbf_map, self._asl_data._asl_image._average_m0
            )

            cbf_map_norm_image = clone_image(base_volume)
            cbf_map_norm_image.update_image_data(
                self._cbf_map * (60 * 60 * 1000),
                self._asl_data._asl_image._average_m0,
            )

            att_map_image = clone_image(base_volume)
            att_map_image.update_image_data(
                self._att_map, self._asl_data._asl_image._average_m0
            )

            t1blgm_map_image = clone_image(base_volume)
            t1blgm_map_image.update_image_data(
                self._t1blgm_map, self._asl_data._asl_image._average_m0
            )

            # Create output maps dictionary
            output_maps = {
                'cbf': cbf_map_image,
                'cbf_norm': cbf_map_norm_image,
                'att': att_map_image,
                't1blgm': t1blgm_map_image,
            }

            # Apply smoothing if requested
            return _apply_smoothing_to_maps(
                output_maps, smoothing, smoothing_params
            )

    def _adjust_image_limits(self, map, init_guess):
        """Adjust image limits by rescaling values within realistic bounds.

        This method removes outliers and rescales T1csfGM values to a realistic
        physiological range based on the initial guess parameter.

        Args:
            map (np.ndarray): The T1csfGM map to adjust
            init_guess (float): Initial guess value used for determining bounds

        Returns:
            np.ndarray: Adjusted map with values rescaled to [0, 2*init_guess]
        """
        # Remove voxels that failed fitting (still at initial guess)
        img = sitk.GetImageFromArray(map * (map != init_guess))

        upper_limit = 2 * init_guess
        lower_limit = 0.0

        rescaler = sitk.RescaleIntensityImageFilter()
        rescaler.SetOutputMaximum(upper_limit)
        rescaler.SetOutputMinimum(lower_limit)

        img_rescaled = rescaler.Execute(img)

        return sitk.GetArrayFromImage(img_rescaled)

__init__(asl_data)

Multi-Echo ASL mapping constructor for T1 tissue relaxometry.

MultiTE_ASLMapping enables advanced ASL analysis by incorporating multiple echo times (TE) to estimate tissue-specific T1 relaxation times. This provides better characterization of blood vs. tissue compartments and improved CBF quantification.

The class requires ASL data acquired with multiple echo times and performs: - Basic CBF and ATT mapping (via CBFMapping) - T1 relaxometry for blood-grey matter differentiation - Multi-TE model fitting for enhanced tissue characterization

Notes

The ASLData object must contain te_values - a list of echo times used during ASL acquisition. These TE values are critical for the multi-echo model fitting and T1 estimation.

Examples:

Basic multi-TE ASL mapping setup:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import MultiTE_ASLMapping
>>> # Create ASL data with multi-TE parameters
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
...     ld_values=[1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8]
... )
>>> mte_mapper = MultiTE_ASLMapping(asl_data)
>>> # Access default MRI parameters
>>> mte_mapper.get_constant('T1csf')
1400.0

Custom MRI parameters for specific field strength:

>>> # Adjust T1 values for 3T scanner
>>> mte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
>>> mte_mapper.get_constant('T1csf')
1600.0
>>> # Verify default parameters unchanged for other objects
>>> from asltk.mri_parameters import MRIParameters
>>> default_param = MRIParameters()
>>> default_param.get_constant('T1csf')
1400.0

Parameters:

Name Type Description Default
asl_data ASLData

The ASL data object containing multi-TE acquisition. Must include te_values, ld_values, and pld_values.

required

Raises:

Type Description
ValueError

If ASLData object lacks required TE values for multi-echo analysis.

See Also

CBFMapping: For basic CBF/ATT mapping without multi-TE analysis MultiDW_ASLMapping: For diffusion-weighted ASL analysis

Source code in asltk/reconstruction/multi_te_mapping.py
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
def __init__(self, asl_data: ASLData) -> None:
    """Multi-Echo ASL mapping constructor for T1 tissue relaxometry.

    MultiTE_ASLMapping enables advanced ASL analysis by incorporating multiple
    echo times (TE) to estimate tissue-specific T1 relaxation times. This
    provides better characterization of blood vs. tissue compartments and
    improved CBF quantification.

    The class requires ASL data acquired with multiple echo times and performs:
    - Basic CBF and ATT mapping (via CBFMapping)
    - T1 relaxometry for blood-grey matter differentiation
    - Multi-TE model fitting for enhanced tissue characterization

    Notes:
        The ASLData object must contain `te_values` - a list of echo times
        used during ASL acquisition. These TE values are critical for the
        multi-echo model fitting and T1 estimation.

    Examples:
        Basic multi-TE ASL mapping setup:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import MultiTE_ASLMapping
        >>> # Create ASL data with multi-TE parameters
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
        ...     ld_values=[1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8]
        ... )
        >>> mte_mapper = MultiTE_ASLMapping(asl_data)
        >>> # Access default MRI parameters
        >>> mte_mapper.get_constant('T1csf')
        1400.0

        Custom MRI parameters for specific field strength:
        >>> # Adjust T1 values for 3T scanner
        >>> mte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
        >>> mte_mapper.get_constant('T1csf')
        1600.0
        >>> # Verify default parameters unchanged for other objects
        >>> from asltk.mri_parameters import MRIParameters
        >>> default_param = MRIParameters()
        >>> default_param.get_constant('T1csf')
        1400.0

    Args:
        asl_data (ASLData): The ASL data object containing multi-TE acquisition.
            Must include te_values, ld_values, and pld_values.

    Raises:
        ValueError: If ASLData object lacks required TE values for multi-echo analysis.

    See Also:
        CBFMapping: For basic CBF/ATT mapping without multi-TE analysis
        MultiDW_ASLMapping: For diffusion-weighted ASL analysis
    """
    super().__init__()
    self._asl_data = asl_data
    self._basic_maps = CBFMapping(asl_data)
    if self._asl_data.get_te() is None:
        raise ValueError(
            'ASLData is incomplete. MultiTE_ASLMapping need a list of TE values.'
        )

    self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
    self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._t1blgm_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

create_map(ub=[np.inf], lb=[0.0], par0=[400], cores=cpu_count(), smoothing=None, smoothing_params=None, suppress_warnings=True)

Create multi-TE ASL maps including T1 blood-gray matter exchange (T1blGM).

This method performs advanced multi-echo ASL analysis to generate tissue-specific T1 relaxation maps that characterize blood-to-gray matter water exchange. The analysis uses multiple echo times to separate blood and tissue signal contributions.

The method implements the multi-compartment TE ASL model described in: "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide blood-to-CSF water transport in humans", NeuroImage, 2022. doi: 10.1016/j.neuroimage.2021.118755

Note

The CBF and ATT maps can be provided before calling this method, using the set_cbf_map() and set_att_map() methods. If not provided, basic CBF/ATT maps are automatically calculated using the CBFMapping class.

Note

The CBF map must be in original scale (not normalized) to perform the correct multiTE-ASL model fitting. Use the 'cbf' output from CBFMapping, not the 'cbf_norm' version.

The method assumes the T1blGM values are well-characterized by the initial guess parameter. Results are filtered to include only positive values and values below 4 times the initial guess to remove unrealistic outliers.

Note

Consider applying spatial smoothing to the output T1blGM map to improve SNR. The create_map() method does not apply filtering by default to preserve spatial resolution.

Parameters:

Name Type Description Default
ub list

Upper bounds for T1blGM fitting. Defaults to [np.inf]. Typically 800-1200 ms for healthy gray matter at 3T.

[inf]
lb list

Lower bounds for T1blGM fitting. Defaults to [0.0]. Should be positive for realistic T1 values.

[0.0]
par0 list

Initial guess for T1blGM in milliseconds. Defaults to [400]. Good starting values: 300-500 ms.

[400]
cores int

Number of CPU threads for parallel processing. Defaults to all available cores.

cpu_count()
smoothing str

Type of spatial smoothing filter to apply. Options: None (default, no smoothing), 'gaussian', 'median'. Smoothing is applied to all output maps after reconstruction.

None
smoothing_params dict

Parameters for the smoothing filter. For 'gaussian': {'sigma': float} (default: 1.0) For 'median': {'size': int} (default: 3)

None
suppress_warnings bool

Whether to suppress warnings during processing. Defaults to True.

True

Returns:

Name Type Description
dict

Dictionary containing: - 'cbf': Basic CBF map in original units (ImageIO) - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO) - 'att': Arterial transit time in ms (ImageIO) - 't1blgm': T1 blood-gray matter exchange time in ms (ImageIO) All maps are smoothed if smoothing is enabled.

Examples:

Basic multi-TE ASL analysis:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import MultiTE_ASLMapping
>>> from asltk.utils.io import ImageIO
>>> import numpy as np
>>> # Load multi-TE ASL data
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
...     ld_values=[1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8]
... )
>>> mte_mapper = MultiTE_ASLMapping(asl_data)
>>> # Set brain mask for faster processing
>>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
>>> mte_mapper.set_brain_mask(brain_mask)
>>> # Generate all maps
>>> results = mte_mapper.create_map()

Custom parameters for specific analysis:

>>> # For expected shorter T1blGM values (faster exchange)
>>> results = mte_mapper.create_map(
...     ub=[600.0],        # Lower upper bound
...     lb=[50.0],         # Minimum realistic T1
...     par0=[300.0]       # Lower initial guess
... )

Apply spatial smoothing to improve SNR:

>>> # Gaussian smoothing with default sigma=1.0
>>> results_smooth = mte_mapper.create_map(
...     smoothing='gaussian'
... )
>>> # Custom smoothing parameters
>>> results_custom = mte_mapper.create_map(
...     smoothing='gaussian',
...     smoothing_params={'sigma': 1.5}
... )
>>> # Median filtering for edge preservation
>>> results_median = mte_mapper.create_map(
...     smoothing='median',
...     smoothing_params={'size': 5}
... )

Raises:

Type Description
ValueError

If cores parameter is invalid or required data is missing.

See Also

set_cbf_map(): Provide pre-computed CBF map set_att_map(): Provide pre-computed ATT map CBFMapping: For basic CBF/ATT mapping

Source code in asltk/reconstruction/multi_te_mapping.py
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
def create_map(
    self,
    ub: list = [np.inf],
    lb: list = [0.0],
    par0: list = [400],
    cores=cpu_count(),
    smoothing=None,
    smoothing_params=None,
    suppress_warnings=True,
):
    """Create multi-TE ASL maps including T1 blood-gray matter exchange (T1blGM).

    This method performs advanced multi-echo ASL analysis to generate tissue-specific
    T1 relaxation maps that characterize blood-to-gray matter water exchange. The
    analysis uses multiple echo times to separate blood and tissue signal contributions.

    The method implements the multi-compartment TE ASL model described in:
    "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
    blood-to-CSF water transport in humans", NeuroImage, 2022.
    doi: 10.1016/j.neuroimage.2021.118755

    Note:
        The CBF and ATT maps can be provided before calling this method,
        using the set_cbf_map() and set_att_map() methods. If not provided,
        basic CBF/ATT maps are automatically calculated using the CBFMapping class.

    Note:
        The CBF map must be in original scale (not normalized) to perform the
        correct multiTE-ASL model fitting. Use the 'cbf' output from CBFMapping,
        not the 'cbf_norm' version.

    The method assumes the T1blGM values are well-characterized by the initial
    guess parameter. Results are filtered to include only positive values and
    values below 4 times the initial guess to remove unrealistic outliers.

    Note:
        Consider applying spatial smoothing to the output T1blGM map to improve
        SNR. The create_map() method does not apply filtering by default to
        preserve spatial resolution.

    Args:
        ub (list, optional): Upper bounds for T1blGM fitting. Defaults to [np.inf].
            Typically 800-1200 ms for healthy gray matter at 3T.
        lb (list, optional): Lower bounds for T1blGM fitting. Defaults to [0.0].
            Should be positive for realistic T1 values.
        par0 (list, optional): Initial guess for T1blGM in milliseconds.
            Defaults to [400]. Good starting values: 300-500 ms.
        cores (int, optional): Number of CPU threads for parallel processing.
            Defaults to all available cores.
        smoothing (str, optional): Type of spatial smoothing filter to apply.
            Options: None (default, no smoothing), 'gaussian', 'median'.
            Smoothing is applied to all output maps after reconstruction.
        smoothing_params (dict, optional): Parameters for the smoothing filter.
            For 'gaussian': {'sigma': float} (default: 1.0)
            For 'median': {'size': int} (default: 3)
        suppress_warnings (bool, optional): Whether to suppress warnings during
            processing. Defaults to True.

    Returns:
        dict: Dictionary containing:
            - 'cbf': Basic CBF map in original units (ImageIO)
            - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO)
            - 'att': Arterial transit time in ms (ImageIO)
            - 't1blgm': T1 blood-gray matter exchange time in ms (ImageIO)
            All maps are smoothed if smoothing is enabled.

    Examples:
        Basic multi-TE ASL analysis:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import MultiTE_ASLMapping
        >>> from asltk.utils.io import ImageIO
        >>> import numpy as np
        >>> # Load multi-TE ASL data
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
        ...     ld_values=[1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8]
        ... )
        >>> mte_mapper = MultiTE_ASLMapping(asl_data)
        >>> # Set brain mask for faster processing
        >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
        >>> mte_mapper.set_brain_mask(brain_mask)
        >>> # Generate all maps
        >>> results = mte_mapper.create_map() # doctest: +SKIP


        Custom parameters for specific analysis:
        >>> # For expected shorter T1blGM values (faster exchange)
        >>> results = mte_mapper.create_map(
        ...     ub=[600.0],        # Lower upper bound
        ...     lb=[50.0],         # Minimum realistic T1
        ...     par0=[300.0]       # Lower initial guess
        ... ) # doctest: +SKIP

        Apply spatial smoothing to improve SNR:
        >>> # Gaussian smoothing with default sigma=1.0
        >>> results_smooth = mte_mapper.create_map(
        ...     smoothing='gaussian'
        ... ) # doctest: +SKIP
        >>> # Custom smoothing parameters
        >>> results_custom = mte_mapper.create_map(
        ...     smoothing='gaussian',
        ...     smoothing_params={'sigma': 1.5}
        ... ) # doctest: +SKIP
        >>> # Median filtering for edge preservation
        >>> results_median = mte_mapper.create_map(
        ...     smoothing='median',
        ...     smoothing_params={'size': 5}
        ... ) # doctest: +SKIP

    Raises:
        ValueError: If cores parameter is invalid or required data is missing.

    See Also:
        set_cbf_map(): Provide pre-computed CBF map
        set_att_map(): Provide pre-computed ATT map
        CBFMapping: For basic CBF/ATT mapping
    """
    # Use context manager to suppress warnings if requested
    with warnings.catch_warnings():
        if suppress_warnings:
            # Filter common warnings that might appear during fitting and processing
            warnings.filterwarnings('ignore', category=RuntimeWarning)
            warnings.filterwarnings('ignore', category=UserWarning)
            warnings.filterwarnings(
                'ignore', category=np.VisibleDeprecationWarning
            )

        self._basic_maps.set_brain_mask(
            ImageIO(image_array=self._brain_mask)
        )

        basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
        if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
            # If the CBF/ATT maps are zero (empty), then a new one is created
            print(
                '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
            )
            basic_maps = self._basic_maps.create_map()
            self._cbf_map = basic_maps['cbf'].get_as_numpy()
            self._att_map = basic_maps['att'].get_as_numpy()

        global asl_data, brain_mask, cbf_map, att_map, t2bl, t2gm
        asl_data = self._asl_data
        brain_mask = self._brain_mask
        cbf_map = self._cbf_map
        att_map = self._att_map
        ld_arr = self._asl_data.get_ld()
        pld_arr = self._asl_data.get_pld()
        te_arr = self._asl_data.get_te()
        t2bl = self.T2bl
        t2gm = self.T2gm

        x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
        y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
        z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

        tblgm_map_shared = Array('d', z_axis * y_axis * x_axis, lock=False)

        with Pool(
            processes=cores,
            initializer=_multite_init_globals,
            initargs=(
                cbf_map,
                att_map,
                brain_mask,
                asl_data,
                ld_arr,
                pld_arr,
                te_arr,
                tblgm_map_shared,
                t2bl,
                t2gm,
            ),
        ) as pool:
            with Progress() as progress:
                task = progress.add_task(
                    'multiTE-ASL processing...', total=x_axis
                )
                results = [
                    pool.apply_async(
                        _tblgm_multite_process_slice,
                        args=(i, x_axis, y_axis, z_axis, par0, lb, ub),
                        callback=lambda _: progress.update(
                            task, advance=1
                        ),
                    )
                    for i in range(x_axis)
                ]
                for result in results:
                    result.wait()

        self._t1blgm_map = np.frombuffer(tblgm_map_shared).reshape(
            z_axis, y_axis, x_axis
        )

        # Adjusting output image boundaries
        self._t1blgm_map = self._adjust_image_limits(
            self._t1blgm_map, par0[0]
        )

        # Prepare output maps
        base_volume = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_image = clone_image(base_volume)
        cbf_map_image.update_image_data(
            self._cbf_map, self._asl_data._asl_image._average_m0
        )

        cbf_map_norm_image = clone_image(base_volume)
        cbf_map_norm_image.update_image_data(
            self._cbf_map * (60 * 60 * 1000),
            self._asl_data._asl_image._average_m0,
        )

        att_map_image = clone_image(base_volume)
        att_map_image.update_image_data(
            self._att_map, self._asl_data._asl_image._average_m0
        )

        t1blgm_map_image = clone_image(base_volume)
        t1blgm_map_image.update_image_data(
            self._t1blgm_map, self._asl_data._asl_image._average_m0
        )

        # Create output maps dictionary
        output_maps = {
            'cbf': cbf_map_image,
            'cbf_norm': cbf_map_norm_image,
            'att': att_map_image,
            't1blgm': t1blgm_map_image,
        }

        # Apply smoothing if requested
        return _apply_smoothing_to_maps(
            output_maps, smoothing, smoothing_params
        )

get_att_map()

Get the ATT map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The ATT map that is storaged in the

MultiTE_ASLMapping object

Source code in asltk/reconstruction/multi_te_mapping.py
169
170
171
172
173
174
175
176
def get_att_map(self):
    """Get the ATT map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The ATT map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._att_map

get_brain_mask()

Get the brain mask image

Returns:

Type Description
ndarray

The brain mask image

Source code in asltk/reconstruction/multi_te_mapping.py
131
132
133
134
135
136
137
def get_brain_mask(self):
    """Get the brain mask image

    Returns:
        (np.ndarray): The brain mask image
    """
    return self._brain_mask

get_cbf_map()

Get the CBF map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The CBF map that is storaged in the

ndarray

MultiTE_ASLMapping object

Source code in asltk/reconstruction/multi_te_mapping.py
152
153
154
155
156
157
158
159
def get_cbf_map(self) -> np.ndarray:
    """Get the CBF map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The CBF map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._cbf_map

get_t1blgm_map()

Get the T1blGM map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The T1blGM map that is storaged in the

MultiTE_ASLMapping object

Source code in asltk/reconstruction/multi_te_mapping.py
178
179
180
181
182
183
184
185
def get_t1blgm_map(self):
    """Get the T1blGM map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The T1blGM map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._t1blgm_map

set_att_map(att_map)

Set the ATT map to the MultiTE_ASLMapping object.

Parameters:

Name Type Description Default
att_map ImageIO

The ATT map that is set in the MultiTE_ASLMapping object

required
Source code in asltk/reconstruction/multi_te_mapping.py
161
162
163
164
165
166
167
def set_att_map(self, att_map: ImageIO):
    """Set the ATT map to the MultiTE_ASLMapping object.

    Args:
        att_map (ImageIO): The ATT map that is set in the MultiTE_ASLMapping object
    """
    self._att_map = att_map.get_as_numpy()

set_brain_mask(brain_mask, label=1)

Defines whether a brain a mask is applied to the CBFMapping calculation

A image mask is simply an image that defines the voxels where the ASL calculation should be made. Basically any integer value can be used as proper label mask.

A most common approach is to use a binary image (zeros for background and 1 for the brain tissues). Anyway, the default behavior of the method can transform a integer-pixel values image to a binary mask with the label parameter provided by the user

Parameters:

Name Type Description Default
brain_mask ndarray

The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.

required
Source code in asltk/reconstruction/multi_te_mapping.py
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
    """Defines whether a brain a mask is applied to the CBFMapping
    calculation

    A image mask is simply an image that defines the voxels where the ASL
    calculation should be made. Basically any integer value can be used as
    proper label mask.

    A most common approach is to use a binary image (zeros for background
    and 1 for the brain tissues). Anyway, the default behavior of the
    method can transform a integer-pixel values image to a binary mask with
    the `label` parameter provided by the user

    Args:
        brain_mask (np.ndarray): The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.
    """
    if not isinstance(brain_mask, ImageIO):
        raise TypeError(
            'The brain_mask parameter must be an instance of ImageIO.'
        )

    _check_mask_values(
        brain_mask, label, self._asl_data('m0').get_as_numpy().shape
    )

    binary_mask = (brain_mask.get_as_numpy() == label).astype(
        np.uint8
    ) * label
    self._brain_mask = binary_mask

set_cbf_map(cbf_map)

Set the CBF map to the MultiTE_ASLMapping object.

Note

The CBF maps must have the original scale in order to calculate the T1blGM map correclty. Hence, if the CBF map was made using CBFMapping class, one can use the 'cbf' output.

Parameters:

Name Type Description Default
cbf_map ImageIO

The CBF map that is set in the MultiTE_ASLMapping object

required
Source code in asltk/reconstruction/multi_te_mapping.py
139
140
141
142
143
144
145
146
147
148
149
150
def set_cbf_map(self, cbf_map: ImageIO):
    """Set the CBF map to the MultiTE_ASLMapping object.

    Note:
        The CBF maps must have the original scale in order to calculate the
        T1blGM map correclty. Hence, if the CBF map was made using
        CBFMapping class, one can use the 'cbf' output.

    Args:
        cbf_map (ImageIO): The CBF map that is set in the MultiTE_ASLMapping object
    """
    self._cbf_map = cbf_map.get_as_numpy()

T2Scalar_ASLMapping

Bases: MRIParameters

Class for voxel-wise T2 mapping from multi-echo ASL data.

This class provides methods to calculate T2 relaxation maps from multi-echo ASL MRI data. It supports brain masking, multiprocessing for fast computation, and optional smoothing.

Main methods
  • set_brain_mask: Set a binary mask to restrict T2 fitting to brain voxels.
  • create_map: Compute T2 maps using multiprocessing (output shape: (N_PLDS, Z, Y, X)).
  • get_t2_maps: Retrieve the computed T2 maps.
  • get_mean_t2s: Retrieve mean T2 values per PLD.
Source code in asltk/reconstruction/t2_mapping.py
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
class T2Scalar_ASLMapping(MRIParameters):
    """
    Class for voxel-wise T2 mapping from multi-echo ASL data.

    This class provides methods to calculate T2 relaxation maps from multi-echo ASL MRI data.
    It supports brain masking, multiprocessing for fast computation, and optional smoothing.

    Main methods:
        - set_brain_mask: Set a binary mask to restrict T2 fitting to brain voxels.
        - create_map: Compute T2 maps using multiprocessing (output shape: (N_PLDS, Z, Y, X)).
        - get_t2_maps: Retrieve the computed T2 maps.
        - get_mean_t2s: Retrieve mean T2 values per PLD.
    """

    def __init__(self, asl_data: ASLData) -> None:
        super().__init__()
        self._asl_data = asl_data
        self._te_values = self._asl_data.get_te()
        self._pld_values = self._asl_data.get_pld()

        # Check if the ASLData has TE and PLD values
        if self._te_values is None or not self._pld_values:
            raise ValueError('ASLData must provide TE and PLD values.')

        # Check if the ASLData has DW values (not allowed for T2 mapping)
        if self._asl_data.get_dw() is not None:
            raise ValueError('ASLData must not include DW values.')

        self._brain_mask = ImageIO(
            image_array=np.ones(self._asl_data('m0').get_as_numpy().shape)
        )
        self._t2_maps = None  # Will be 4D: (N_PLDS, Z, Y, X)
        self._mean_t2s = None

    def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
        """
        Set a brain mask to restrict T2 fitting to specific voxels.

        Args:
            brain_mask (np.ndarray): Binary or integer mask with the same shape as the M0 image. Nonzero values indicate voxels to include.
            label (int, optional): The label value to use as foreground (default: 1).

        The mask should be a 3D numpy array matching the spatial dimensions of the ASL data.
        """
        _check_mask_values(
            brain_mask, label, self._asl_data('m0').get_as_numpy().shape
        )

        binary_mask = ImageIO(
            image_array=(brain_mask.get_as_numpy() == label).astype(np.uint8)
            * label
        )
        self._brain_mask = binary_mask

    def get_t2_maps(self):
        """Get the T2 maps storaged at the T2Scalar_ASLMapping object

        Returns:
            (np.ndarray): The T2 maps that is storaged in the
            T2Scalar_ASLMapping object
        """
        return self._t2_maps

    def get_mean_t2s(self):
        """Get the mean T2 values calculated from the T2 maps

        Returns:
            (list): The mean T2 values for each PLD
        """
        return self._mean_t2s

    def create_map(
        self,
        cores=cpu_count(),
        smoothing=None,
        smoothing_params=None,
        suppress_warnings=False,
    ):
        """
        Compute T2 maps using multi-echo ASL data and a brain mask, with multiprocessing.

        This method uses multiprocessing to accelerate voxel-wise T2 fitting. The output is a 4D array with shape (N_PLDS, Z, Y, X).

        Warning:
            For large datasets, memory usage can be significant due to parallel processing and storage of intermediate arrays.

        Args:
            cores (int, optional): Number of CPU cores for processing. Defaults to all available.
            smoothing (str, optional): Smoothing type ('gaussian', 'median', or None).
            smoothing_params (dict, optional): Smoothing parameters.

        Returns:
            dict: Dictionary with T2 maps ('t2', shape (N_PLDS, Z, Y, X)) and mean T2 values ('mean_t2').
        """
        logger = get_logger('t2_mapping')
        logger.info('Starting T2 map creation')

        # Optionally suppress warnings
        if suppress_warnings:
            warnings_context = warnings.catch_warnings()
            warnings_context.__enter__()
            warnings.filterwarnings('ignore', category=RuntimeWarning)
            warnings.filterwarnings('ignore', category=UserWarning)
            logger.info('Warnings suppressed during T2 mapping')

        try:
            data = self._asl_data('pcasl').get_as_numpy()
            mask = self._brain_mask.get_as_numpy()
            TEs = np.array(self._te_values)
            PLDs = np.array(self._pld_values)
            n_tes, n_plds, z_axis, y_axis, x_axis = data.shape

            t2_maps_all = []
            mean_t2s = []

            for pld_idx in range(n_plds):
                logger.info(
                    f'Processing PLD index {pld_idx} ({PLDs[pld_idx]} ms)'
                )
                t2_map_shared = Array(
                    'd', z_axis * y_axis * x_axis, lock=False
                )
                log_processing_step(
                    'Running voxel-wise T2 fitting',
                    'this may take several minutes',
                )
                with Pool(
                    processes=cores,
                    initializer=_t2_init_globals,
                    initargs=(t2_map_shared, mask, data, TEs),
                ) as pool:
                    with Progress() as progress:
                        task = progress.add_task(
                            f'T2 fitting (PLD {PLDs[pld_idx]} ms)...',
                            total=x_axis,
                        )
                        results = [
                            pool.apply_async(
                                _t2_process_slice,
                                args=(i, x_axis, y_axis, z_axis, pld_idx),
                                callback=lambda _: progress.update(
                                    task, advance=1
                                ),
                            )
                            for i in range(x_axis)
                        ]
                        for result in results:
                            result.wait()

                t2_map = np.frombuffer(t2_map_shared).reshape(
                    z_axis, y_axis, x_axis
                )
                t2_maps_all.append(t2_map)
                mean_t2s.append(np.nanmean(t2_map))

            t2_maps_stacked = np.array(t2_maps_all)  # shape: (N_PLDS, Z, Y, X)
            self._t2_maps = t2_maps_stacked
            self._mean_t2s = mean_t2s

            logger.info('T2 mapping completed successfully')
            logger.info(
                f'T2 statistics - Mean: {np.mean(self._t2_maps):.4f}, Std: {np.std(self._t2_maps):.4f}'
            )

            # Prepare output maps
            # TODO At the moment, the T2 maps and mean T2 maps are as ImageIO object, however, the Spacing, Dimension are not given as a 4D array. Check if can be imported from the m0 image is 3D.
            t2_maps_image = ImageIO(
                image_array=np.array(
                    [
                        self._asl_data('m0').get_as_numpy()
                        for _ in range(len(t2_maps_all))
                    ]
                )
            )
            t2_maps_image.update_image_data(self._t2_maps)

            # Update the _t2_maps attribute to be an ImageIO object
            self._t2_maps = t2_maps_image

            output_maps = {
                't2': t2_maps_image,
                'mean_t2': self._mean_t2s,
            }

            return _apply_smoothing_to_maps(
                output_maps, smoothing, smoothing_params
            )
        finally:
            # Ensure warnings are restored if suppressed
            if suppress_warnings:
                warnings_context.__exit__(None, None, None)

create_map(cores=cpu_count(), smoothing=None, smoothing_params=None, suppress_warnings=False)

Compute T2 maps using multi-echo ASL data and a brain mask, with multiprocessing.

This method uses multiprocessing to accelerate voxel-wise T2 fitting. The output is a 4D array with shape (N_PLDS, Z, Y, X).

Warning

For large datasets, memory usage can be significant due to parallel processing and storage of intermediate arrays.

Parameters:

Name Type Description Default
cores int

Number of CPU cores for processing. Defaults to all available.

cpu_count()
smoothing str

Smoothing type ('gaussian', 'median', or None).

None
smoothing_params dict

Smoothing parameters.

None

Returns:

Name Type Description
dict

Dictionary with T2 maps ('t2', shape (N_PLDS, Z, Y, X)) and mean T2 values ('mean_t2').

Source code in asltk/reconstruction/t2_mapping.py
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
def create_map(
    self,
    cores=cpu_count(),
    smoothing=None,
    smoothing_params=None,
    suppress_warnings=False,
):
    """
    Compute T2 maps using multi-echo ASL data and a brain mask, with multiprocessing.

    This method uses multiprocessing to accelerate voxel-wise T2 fitting. The output is a 4D array with shape (N_PLDS, Z, Y, X).

    Warning:
        For large datasets, memory usage can be significant due to parallel processing and storage of intermediate arrays.

    Args:
        cores (int, optional): Number of CPU cores for processing. Defaults to all available.
        smoothing (str, optional): Smoothing type ('gaussian', 'median', or None).
        smoothing_params (dict, optional): Smoothing parameters.

    Returns:
        dict: Dictionary with T2 maps ('t2', shape (N_PLDS, Z, Y, X)) and mean T2 values ('mean_t2').
    """
    logger = get_logger('t2_mapping')
    logger.info('Starting T2 map creation')

    # Optionally suppress warnings
    if suppress_warnings:
        warnings_context = warnings.catch_warnings()
        warnings_context.__enter__()
        warnings.filterwarnings('ignore', category=RuntimeWarning)
        warnings.filterwarnings('ignore', category=UserWarning)
        logger.info('Warnings suppressed during T2 mapping')

    try:
        data = self._asl_data('pcasl').get_as_numpy()
        mask = self._brain_mask.get_as_numpy()
        TEs = np.array(self._te_values)
        PLDs = np.array(self._pld_values)
        n_tes, n_plds, z_axis, y_axis, x_axis = data.shape

        t2_maps_all = []
        mean_t2s = []

        for pld_idx in range(n_plds):
            logger.info(
                f'Processing PLD index {pld_idx} ({PLDs[pld_idx]} ms)'
            )
            t2_map_shared = Array(
                'd', z_axis * y_axis * x_axis, lock=False
            )
            log_processing_step(
                'Running voxel-wise T2 fitting',
                'this may take several minutes',
            )
            with Pool(
                processes=cores,
                initializer=_t2_init_globals,
                initargs=(t2_map_shared, mask, data, TEs),
            ) as pool:
                with Progress() as progress:
                    task = progress.add_task(
                        f'T2 fitting (PLD {PLDs[pld_idx]} ms)...',
                        total=x_axis,
                    )
                    results = [
                        pool.apply_async(
                            _t2_process_slice,
                            args=(i, x_axis, y_axis, z_axis, pld_idx),
                            callback=lambda _: progress.update(
                                task, advance=1
                            ),
                        )
                        for i in range(x_axis)
                    ]
                    for result in results:
                        result.wait()

            t2_map = np.frombuffer(t2_map_shared).reshape(
                z_axis, y_axis, x_axis
            )
            t2_maps_all.append(t2_map)
            mean_t2s.append(np.nanmean(t2_map))

        t2_maps_stacked = np.array(t2_maps_all)  # shape: (N_PLDS, Z, Y, X)
        self._t2_maps = t2_maps_stacked
        self._mean_t2s = mean_t2s

        logger.info('T2 mapping completed successfully')
        logger.info(
            f'T2 statistics - Mean: {np.mean(self._t2_maps):.4f}, Std: {np.std(self._t2_maps):.4f}'
        )

        # Prepare output maps
        # TODO At the moment, the T2 maps and mean T2 maps are as ImageIO object, however, the Spacing, Dimension are not given as a 4D array. Check if can be imported from the m0 image is 3D.
        t2_maps_image = ImageIO(
            image_array=np.array(
                [
                    self._asl_data('m0').get_as_numpy()
                    for _ in range(len(t2_maps_all))
                ]
            )
        )
        t2_maps_image.update_image_data(self._t2_maps)

        # Update the _t2_maps attribute to be an ImageIO object
        self._t2_maps = t2_maps_image

        output_maps = {
            't2': t2_maps_image,
            'mean_t2': self._mean_t2s,
        }

        return _apply_smoothing_to_maps(
            output_maps, smoothing, smoothing_params
        )
    finally:
        # Ensure warnings are restored if suppressed
        if suppress_warnings:
            warnings_context.__exit__(None, None, None)

get_mean_t2s()

Get the mean T2 values calculated from the T2 maps

Returns:

Type Description
list

The mean T2 values for each PLD

Source code in asltk/reconstruction/t2_mapping.py
85
86
87
88
89
90
91
def get_mean_t2s(self):
    """Get the mean T2 values calculated from the T2 maps

    Returns:
        (list): The mean T2 values for each PLD
    """
    return self._mean_t2s

get_t2_maps()

Get the T2 maps storaged at the T2Scalar_ASLMapping object

Returns:

Type Description
ndarray

The T2 maps that is storaged in the

T2Scalar_ASLMapping object

Source code in asltk/reconstruction/t2_mapping.py
76
77
78
79
80
81
82
83
def get_t2_maps(self):
    """Get the T2 maps storaged at the T2Scalar_ASLMapping object

    Returns:
        (np.ndarray): The T2 maps that is storaged in the
        T2Scalar_ASLMapping object
    """
    return self._t2_maps

set_brain_mask(brain_mask, label=1)

Set a brain mask to restrict T2 fitting to specific voxels.

Parameters:

Name Type Description Default
brain_mask ndarray

Binary or integer mask with the same shape as the M0 image. Nonzero values indicate voxels to include.

required
label int

The label value to use as foreground (default: 1).

1

The mask should be a 3D numpy array matching the spatial dimensions of the ASL data.

Source code in asltk/reconstruction/t2_mapping.py
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
    """
    Set a brain mask to restrict T2 fitting to specific voxels.

    Args:
        brain_mask (np.ndarray): Binary or integer mask with the same shape as the M0 image. Nonzero values indicate voxels to include.
        label (int, optional): The label value to use as foreground (default: 1).

    The mask should be a 3D numpy array matching the spatial dimensions of the ASL data.
    """
    _check_mask_values(
        brain_mask, label, self._asl_data('m0').get_as_numpy().shape
    )

    binary_mask = ImageIO(
        image_array=(brain_mask.get_as_numpy() == label).astype(np.uint8)
        * label
    )
    self._brain_mask = binary_mask

UltraLongTE_ASLMapping

Bases: MRIParameters

Source code in asltk/reconstruction/ultralong_te_mapping.py
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
class UltraLongTE_ASLMapping(MRIParameters):
    def __init__(self, asl_data: ASLData) -> None:
        """UltraLongTE ASL mapping constructor for T1 time exchange tissue relaxometry.

        UltraLongTE_ASLMapping enables advanced ASL analysis by incorporating multiple
        echo times (TE) to estimate tissue-specific T1 relaxation times. This
        provides better characterization of blood vs. tissue compartments and
        improved CBF quantification.

        The class requires ASL data acquired with multiple echo times and performs:
        - Basic CBF and ATT mapping (via CBFMapping)
        - T1 relaxometry for blood-grey matter differentiation
        - Ultralong-TE model fitting for enhanced tissue characterization

        Notes:
            The ASLData object must contain `te_values` - a list of echo times
            used during ASL acquisition. These TE values are critical for the
            multi-echo model fitting and T1 estimation.

        Notes:
            This method is based from the original paper of:
            Leonie Petitclerc, Lydiane Hirschler, Jack A. Wells, David L. Thomas,
            Marianne A.A. van Walderveen, Mark A. van Buchem, Matthias J.P. van Osch,
            "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
            blood-to-CSF water transport in humans", NeuroImage, ISSN 1053-8119,
            doi: [10.1016/j.neuroimage.2021.118755](http://doi.org/10.1016/j.neuroimage.2021.118755).

        Important:
            Although this method applies a parallel CPU processing, it still a
            highly time and memory consuming procedure. Take this into account
            when working with large datasets or high-resolution images.

        Examples:
            Basic Ultralong-TE ASL mapping setup:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import UltraLongTE_ASLMapping
            >>> # Create ASL data with multi-TE parameters
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
            ...     ld_values=[1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8]
            ... )
            >>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
            >>> # Access default MRI parameters
            >>> ulte_mapper.get_constant('T1csf')
            1400.0

            Custom MRI parameters for specific field strength:
            >>> # Adjust T1 values for 3T scanner
            >>> ulte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
            >>> ulte_mapper.get_constant('T1csf')
            1600.0
            >>> # Verify default parameters unchanged for other objects
            >>> from asltk.mri_parameters import MRIParameters
            >>> default_param = MRIParameters()
            >>> default_param.get_constant('T1csf')
            1400.0

        Args:
            asl_data (ASLData): The ASL data object containing ultralong-TE acquisition.
                Must include te_values, ld_values, and pld_values.

        Raises:
            ValueError: If ASLData object lacks required TE values for ultralong-TE analysis.

        See Also:
            CBFMapping: For basic CBF/ATT mapping without multi-echo analysis
            MultiDW_ASLMapping: For diffusion-weighted ASL analysis
            MultiTE_ASLMapping: For the multi-echo TE ASL analysis (the predecessor of this method)
        """
        super().__init__()
        self._asl_data = asl_data
        self._basic_maps = CBFMapping(asl_data)
        if self._asl_data.get_te() is None:
            raise ValueError(
                'ASLData is incomplete. UltraLongTE_ASLMapping need a list of TE values.'
            )

        self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
        self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
        self._t1csfgm_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

        # Changing the T2csf and T2blood as requested in the original paper
        self.set_constant(1500.0, 'T2csf')  # T2 relaxation time for CSF in ms
        self.set_constant(100.0, 'T2bl')  # T2 relaxation time for blood in ms

    def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
        """Defines whether a brain a mask is applied to the CBFMapping
        calculation

        A image mask is simply an image that defines the voxels where the ASL
        calculation should be made. Basically any integer value can be used as
        proper label mask.

        A most common approach is to use a binary image (zeros for background
        and 1 for the brain tissues). Anyway, the default behavior of the
        method can transform a integer-pixel values image to a binary mask with
        the `label` parameter provided by the user

        Args:
            brain_mask (np.ndarray): The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.
        """
        if not isinstance(brain_mask, ImageIO):
            raise TypeError(
                'The brain_mask parameter must be an instance of ImageIO.'
            )

        _check_mask_values(
            brain_mask, label, self._asl_data('m0').get_as_numpy().shape
        )

        binary_mask = (brain_mask.get_as_numpy() == label).astype(
            np.uint8
        ) * label
        self._brain_mask = binary_mask

    def get_brain_mask(self):
        """Get the brain mask image

        Returns:
            (ImageIO): The brain mask image
        """
        return self._brain_mask

    def set_cbf_map(self, cbf_map: ImageIO):
        """Set the CBF map to the MultiTE_ASLMapping object.

        Note:
            The CBF maps must have the original scale in order to calculate the
            T1 CSF-GM map correclty. Hence, if the CBF map was made using
            CBFMapping class, one can use the 'cbf' output.

        Args:
            cbf_map (ImageIO): The CBF map that is set in the MultiTE_ASLMapping object
        """
        self._cbf_map = cbf_map.get_as_numpy()

    def get_cbf_map(self) -> np.ndarray:
        """Get the CBF map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The CBF map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._cbf_map

    def set_att_map(self, att_map: ImageIO):
        """Set the ATT map to the MultiTE_ASLMapping object.

        Args:
            att_map (ImageIO): The ATT map that is set in the MultiTE_ASLMapping object
        """
        self._att_map = att_map.get_as_numpy()

    def get_att_map(self):
        """Get the ATT map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The ATT map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._att_map

    def get_t1csfgm_map(self):
        """Get the T1csfGM map storaged at the MultiTE_ASLMapping object

        Returns:
            (ImageIO): The T1csfGM map that is storaged in the
            MultiTE_ASLMapping object
        """
        return self._t1csfgm_map

    def create_map(
        self,
        ub: list = [np.inf],
        lb: list = [0.0],
        par0: list = [80000],
        cores: Union[int, str] = 'auto',
        smoothing=None,
        smoothing_params=None,
        suppress_warnings=True,
    ):
        """Create ultra-long-TE ASL maps including T1 csf-gray matter exchange (T1csfGM).

        This method performs advanced multi-echo ASL analysis to generate tissue-specific
        T1 relaxation maps that characterize blood-to-gray matter water exchange. The
        analysis uses multiple echo times to separate blood and tissue signal contributions.

        Note:
            The method implements the multi-compartment TE ASL model described in:
            "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
            blood-to-CSF water transport in humans", NeuroImage, 2022.
            doi: [10.1016/j.neuroimage.2021.118755](http://doi.org/10.1016/j.neuroimage.2021.118755)

        Note:
            The CBF and ATT maps can be provided before calling this method,
            using the set_cbf_map() and set_att_map() methods. If not provided,
            basic CBF/ATT maps are automatically calculated using the CBFMapping class.

        Note:
            The CBF map must be in original scale (not normalized) to perform the
            correct ultralong-TE-ASL model fitting. Use the 'cbf' output from CBFMapping,
            not the 'cbf_norm' version.

        The method assumes the T1csfGM values are well-characterized by the initial
        guess parameter. Results are filtered to include only positive values and
        values below 4 times the initial guess to remove unrealistic outliers.

        Note:
            Consider applying spatial smoothing to the output T1csfGM map to improve
            SNR. The create_map() method does not apply filtering by default to
            preserve spatial resolution.

        Args:
            ub (list, optional): Upper bounds for T1csfGM fitting. Defaults to [np.inf].
                Typically 800-1200 ms for healthy gray matter at 3T.
            lb (list, optional): Lower bounds for T1csfGM fitting. Defaults to [0.0].
                Should be positive for realistic T1 values.
            par0 (list, optional): Initial guess for T1csfGM in milliseconds.
                Defaults to [400]. Good starting values: 300-500 ms.
            cores (int or str, optional): Number of CPU threads for parallel processing.
                If "auto" (default), automatically determines the optimal number based on
                available system memory. If an integer is provided, uses that specific number.
            smoothing (str, optional): Type of spatial smoothing filter to apply.
                Options: None (default, no smoothing), 'gaussian', 'median'.
                Smoothing is applied to all output maps after reconstruction.
            smoothing_params (dict, optional): Parameters for the smoothing filter.
                For 'gaussian': {'sigma': float} (default: 1.0)
                For 'median': {'size': int} (default: 3)
            suppress_warnings (bool, optional): Whether to suppress warnings during
                processing. Defaults to True.

        Returns:
            dict: Dictionary containing:
                - 'cbf': Basic CBF map in original units (ImageIO)
                - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO)
                - 'att': Arterial transit time in ms (ImageIO)
                - 't1csfgm': T1 csf-gray matter exchange time in ms (ImageIO)
                All maps are smoothed if smoothing is enabled.

        Examples:
            Basic multi-TE ASL analysis:
            >>> from asltk.asldata import ASLData
            >>> from asltk.reconstruction import UltraLongTE_ASLMapping
            >>> from asltk.utils.io import ImageIO
            >>> import numpy as np
            >>> # Load multi-TE ASL data
            >>> asl_data = ASLData(
            ...     pcasl='./tests/files/pcasl_mte.nii.gz',
            ...     m0='./tests/files/m0.nii.gz',
            ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
            ...     ld_values=[1.8, 1.8, 1.8],
            ...     pld_values=[0.8, 1.8, 2.8]
            ... )
            >>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
            >>> # Set brain mask for faster processing
            >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
            >>> ulte_mapper.set_brain_mask(brain_mask)
            >>> # Generate all maps
            >>> results = ulte_mapper.create_map() # doctest: +SKIP


            Custom parameters for specific analysis:
            >>> # For expected shorter T1csfGM values (faster exchange)
            >>> results = ulte_mapper.create_map(
            ...     ub=[600.0],        # Lower upper bound
            ...     lb=[50.0],         # Minimum realistic T1
            ...     par0=[300.0]       # Lower initial guess
            ... ) # doctest: +SKIP

            Apply spatial smoothing to improve SNR:
            >>> # Gaussian smoothing with default sigma=1.0
            >>> results_smooth = ulte_mapper.create_map(
            ...     smoothing='gaussian'
            ... ) # doctest: +SKIP
            >>> # Custom smoothing parameters
            >>> results_custom = ulte_mapper.create_map(
            ...     smoothing='gaussian',
            ...     smoothing_params={'sigma': 1.5}
            ... ) # doctest: +SKIP
            >>> # Median filtering for edge preservation
            >>> results_median = ulte_mapper.create_map(
            ...     smoothing='median',
            ...     smoothing_params={'size': 5}
            ... ) # doctest: +SKIP

        Raises:
            ValueError: If cores parameter is invalid or required data is missing.

        See Also:
            set_cbf_map(): Provide pre-computed CBF map
            set_att_map(): Provide pre-computed ATT map
            CBFMapping: For basic CBF/ATT mapping
        """
        # Determine optimal number of cores based on available memory
        actual_cores = get_optimal_core_count(cores)

        # Use context manager to suppress warnings if requested
        with warnings.catch_warnings():
            if suppress_warnings:
                # Filter common warnings that might appear during fitting and processing
                warnings.filterwarnings('ignore', category=RuntimeWarning)
                warnings.filterwarnings('ignore', category=UserWarning)
                warnings.filterwarnings(
                    'ignore', category=np.VisibleDeprecationWarning
                )

            self._basic_maps.set_brain_mask(
                ImageIO(image_array=self._brain_mask)
            )

            basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
            if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
                # If the CBF/ATT maps are zero (empty), then a new one is created
                print(
                    '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
                )
                basic_maps = self._basic_maps.create_map()
                self._cbf_map = basic_maps['cbf'].get_as_numpy()
                self._att_map = basic_maps['att'].get_as_numpy()

            global asl_data, brain_mask, cbf_map, att_map, t2bl, t2gm
            asl_data = self._asl_data
            brain_mask = self._brain_mask
            cbf_map = self._cbf_map
            att_map = self._att_map
            ld_arr = self._asl_data.get_ld()
            pld_arr = self._asl_data.get_pld()
            te_arr = self._asl_data.get_te()
            t2bl = self.T2bl
            t2gm = self.T2gm

            x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
            y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
            z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

            tcsfgm_map_shared = Array(
                'd', z_axis * y_axis * x_axis, lock=False
            )

            # Make a copy of base information
            m0_array = asl_data('m0').get_as_numpy()
            pcasl_array = asl_data('pcasl').get_as_numpy()

            with Pool(
                processes=actual_cores,
                initializer=_multite_init_globals,
                initargs=(
                    cbf_map,
                    att_map,
                    brain_mask,
                    asl_data,
                    ld_arr,
                    pld_arr,
                    te_arr,
                    tcsfgm_map_shared,
                    t2bl,
                    t2gm,
                ),
            ) as pool:
                with Progress() as progress:
                    task = progress.add_task(
                        'ultralongTE-ASL processing...', total=x_axis
                    )
                    results = [
                        pool.apply_async(
                            _tcsfgm_multite_process_slice,
                            args=(
                                i,
                                x_axis,
                                y_axis,
                                z_axis,
                                par0,
                                lb,
                                ub,
                                m0_array,
                                pcasl_array,
                            ),
                            callback=lambda _: progress.update(
                                task, advance=1
                            ),
                        )
                        for i in range(x_axis)
                    ]
                    for result in results:
                        result.wait()

            self._t1csfgm_map = np.frombuffer(tcsfgm_map_shared).reshape(
                z_axis, y_axis, x_axis
            )

            # Adjusting output image boundaries
            # TODO O ADJUST_LIMITS TEM QUE SER TESTADO E AJUSTADO PARA O CASO DO CSF-GM --> OK para a nova forma de adjust_image. Fazer a mesma coisa para multiTE e commitar o codigo
            self._t1csfgm_map = self._adjust_image_limits(
                self._t1csfgm_map, par0[0]
            )

            # Prepare output maps
            cbf_map_image = ImageIO(self._asl_data('m0').get_image_path())
            cbf_map_image.update_image_data(self._cbf_map)

            cbf_map_norm_image = ImageIO(self._asl_data('m0').get_image_path())
            cbf_map_norm_image.update_image_data(
                self._cbf_map * (60 * 60 * 1000)
            )

            att_map_image = ImageIO(self._asl_data('m0').get_image_path())
            att_map_image.update_image_data(self._att_map)

            t1csfgm_map_image = ImageIO(self._asl_data('m0').get_image_path())
            t1csfgm_map_image.update_image_data(self._t1csfgm_map)

            # Create output maps dictionary
            output_maps = {
                'cbf': cbf_map_image,
                'cbf_norm': cbf_map_norm_image,
                'att': att_map_image,
                't1csfgm': t1csfgm_map_image,
            }

            # Apply smoothing if requested
            return _apply_smoothing_to_maps(
                output_maps, smoothing, smoothing_params
            )

    def _adjust_image_limits(self, map, init_guess):
        """Adjust image limits by rescaling values within realistic bounds.

        This method removes outliers and rescales T1csfGM values to a realistic
        physiological range based on the initial guess parameter.

        Args:
            map (np.ndarray): The T1csfGM map to adjust
            init_guess (float): Initial guess value used for determining bounds

        Returns:
            np.ndarray: Adjusted map with values rescaled to [0, 2*init_guess]
        """
        # Remove voxels that failed fitting (still at initial guess)
        img = sitk.GetImageFromArray(map * (map != init_guess))

        upper_limit = 2 * init_guess
        lower_limit = 0.0

        rescaler = sitk.RescaleIntensityImageFilter()
        rescaler.SetOutputMaximum(upper_limit)
        rescaler.SetOutputMinimum(lower_limit)

        img_rescaled = rescaler.Execute(img)

        return sitk.GetArrayFromImage(img_rescaled)

__init__(asl_data)

UltraLongTE ASL mapping constructor for T1 time exchange tissue relaxometry.

UltraLongTE_ASLMapping enables advanced ASL analysis by incorporating multiple echo times (TE) to estimate tissue-specific T1 relaxation times. This provides better characterization of blood vs. tissue compartments and improved CBF quantification.

The class requires ASL data acquired with multiple echo times and performs: - Basic CBF and ATT mapping (via CBFMapping) - T1 relaxometry for blood-grey matter differentiation - Ultralong-TE model fitting for enhanced tissue characterization

Notes

The ASLData object must contain te_values - a list of echo times used during ASL acquisition. These TE values are critical for the multi-echo model fitting and T1 estimation.

Notes

This method is based from the original paper of: Leonie Petitclerc, Lydiane Hirschler, Jack A. Wells, David L. Thomas, Marianne A.A. van Walderveen, Mark A. van Buchem, Matthias J.P. van Osch, "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide blood-to-CSF water transport in humans", NeuroImage, ISSN 1053-8119, doi: 10.1016/j.neuroimage.2021.118755.

Important

Although this method applies a parallel CPU processing, it still a highly time and memory consuming procedure. Take this into account when working with large datasets or high-resolution images.

Examples:

Basic Ultralong-TE ASL mapping setup:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import UltraLongTE_ASLMapping
>>> # Create ASL data with multi-TE parameters
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
...     ld_values=[1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8]
... )
>>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
>>> # Access default MRI parameters
>>> ulte_mapper.get_constant('T1csf')
1400.0

Custom MRI parameters for specific field strength:

>>> # Adjust T1 values for 3T scanner
>>> ulte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
>>> ulte_mapper.get_constant('T1csf')
1600.0
>>> # Verify default parameters unchanged for other objects
>>> from asltk.mri_parameters import MRIParameters
>>> default_param = MRIParameters()
>>> default_param.get_constant('T1csf')
1400.0

Parameters:

Name Type Description Default
asl_data ASLData

The ASL data object containing ultralong-TE acquisition. Must include te_values, ld_values, and pld_values.

required

Raises:

Type Description
ValueError

If ASLData object lacks required TE values for ultralong-TE analysis.

See Also

CBFMapping: For basic CBF/ATT mapping without multi-echo analysis MultiDW_ASLMapping: For diffusion-weighted ASL analysis MultiTE_ASLMapping: For the multi-echo TE ASL analysis (the predecessor of this method)

Source code in asltk/reconstruction/ultralong_te_mapping.py
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
def __init__(self, asl_data: ASLData) -> None:
    """UltraLongTE ASL mapping constructor for T1 time exchange tissue relaxometry.

    UltraLongTE_ASLMapping enables advanced ASL analysis by incorporating multiple
    echo times (TE) to estimate tissue-specific T1 relaxation times. This
    provides better characterization of blood vs. tissue compartments and
    improved CBF quantification.

    The class requires ASL data acquired with multiple echo times and performs:
    - Basic CBF and ATT mapping (via CBFMapping)
    - T1 relaxometry for blood-grey matter differentiation
    - Ultralong-TE model fitting for enhanced tissue characterization

    Notes:
        The ASLData object must contain `te_values` - a list of echo times
        used during ASL acquisition. These TE values are critical for the
        multi-echo model fitting and T1 estimation.

    Notes:
        This method is based from the original paper of:
        Leonie Petitclerc, Lydiane Hirschler, Jack A. Wells, David L. Thomas,
        Marianne A.A. van Walderveen, Mark A. van Buchem, Matthias J.P. van Osch,
        "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
        blood-to-CSF water transport in humans", NeuroImage, ISSN 1053-8119,
        doi: [10.1016/j.neuroimage.2021.118755](http://doi.org/10.1016/j.neuroimage.2021.118755).

    Important:
        Although this method applies a parallel CPU processing, it still a
        highly time and memory consuming procedure. Take this into account
        when working with large datasets or high-resolution images.

    Examples:
        Basic Ultralong-TE ASL mapping setup:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import UltraLongTE_ASLMapping
        >>> # Create ASL data with multi-TE parameters
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
        ...     ld_values=[1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8]
        ... )
        >>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
        >>> # Access default MRI parameters
        >>> ulte_mapper.get_constant('T1csf')
        1400.0

        Custom MRI parameters for specific field strength:
        >>> # Adjust T1 values for 3T scanner
        >>> ulte_mapper.set_constant(1600.0, 'T1csf')  # CSF T1 at 3T
        >>> ulte_mapper.get_constant('T1csf')
        1600.0
        >>> # Verify default parameters unchanged for other objects
        >>> from asltk.mri_parameters import MRIParameters
        >>> default_param = MRIParameters()
        >>> default_param.get_constant('T1csf')
        1400.0

    Args:
        asl_data (ASLData): The ASL data object containing ultralong-TE acquisition.
            Must include te_values, ld_values, and pld_values.

    Raises:
        ValueError: If ASLData object lacks required TE values for ultralong-TE analysis.

    See Also:
        CBFMapping: For basic CBF/ATT mapping without multi-echo analysis
        MultiDW_ASLMapping: For diffusion-weighted ASL analysis
        MultiTE_ASLMapping: For the multi-echo TE ASL analysis (the predecessor of this method)
    """
    super().__init__()
    self._asl_data = asl_data
    self._basic_maps = CBFMapping(asl_data)
    if self._asl_data.get_te() is None:
        raise ValueError(
            'ASLData is incomplete. UltraLongTE_ASLMapping need a list of TE values.'
        )

    self._brain_mask = np.ones(self._asl_data('m0').get_as_numpy().shape)
    self._cbf_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._att_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)
    self._t1csfgm_map = np.zeros(self._asl_data('m0').get_as_numpy().shape)

    # Changing the T2csf and T2blood as requested in the original paper
    self.set_constant(1500.0, 'T2csf')  # T2 relaxation time for CSF in ms
    self.set_constant(100.0, 'T2bl')  # T2 relaxation time for blood in ms

create_map(ub=[np.inf], lb=[0.0], par0=[80000], cores='auto', smoothing=None, smoothing_params=None, suppress_warnings=True)

Create ultra-long-TE ASL maps including T1 csf-gray matter exchange (T1csfGM).

This method performs advanced multi-echo ASL analysis to generate tissue-specific T1 relaxation maps that characterize blood-to-gray matter water exchange. The analysis uses multiple echo times to separate blood and tissue signal contributions.

Note

The method implements the multi-compartment TE ASL model described in: "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide blood-to-CSF water transport in humans", NeuroImage, 2022. doi: 10.1016/j.neuroimage.2021.118755

Note

The CBF and ATT maps can be provided before calling this method, using the set_cbf_map() and set_att_map() methods. If not provided, basic CBF/ATT maps are automatically calculated using the CBFMapping class.

Note

The CBF map must be in original scale (not normalized) to perform the correct ultralong-TE-ASL model fitting. Use the 'cbf' output from CBFMapping, not the 'cbf_norm' version.

The method assumes the T1csfGM values are well-characterized by the initial guess parameter. Results are filtered to include only positive values and values below 4 times the initial guess to remove unrealistic outliers.

Note

Consider applying spatial smoothing to the output T1csfGM map to improve SNR. The create_map() method does not apply filtering by default to preserve spatial resolution.

Parameters:

Name Type Description Default
ub list

Upper bounds for T1csfGM fitting. Defaults to [np.inf]. Typically 800-1200 ms for healthy gray matter at 3T.

[inf]
lb list

Lower bounds for T1csfGM fitting. Defaults to [0.0]. Should be positive for realistic T1 values.

[0.0]
par0 list

Initial guess for T1csfGM in milliseconds. Defaults to [400]. Good starting values: 300-500 ms.

[80000]
cores int or str

Number of CPU threads for parallel processing. If "auto" (default), automatically determines the optimal number based on available system memory. If an integer is provided, uses that specific number.

'auto'
smoothing str

Type of spatial smoothing filter to apply. Options: None (default, no smoothing), 'gaussian', 'median'. Smoothing is applied to all output maps after reconstruction.

None
smoothing_params dict

Parameters for the smoothing filter. For 'gaussian': {'sigma': float} (default: 1.0) For 'median': {'size': int} (default: 3)

None
suppress_warnings bool

Whether to suppress warnings during processing. Defaults to True.

True

Returns:

Name Type Description
dict

Dictionary containing: - 'cbf': Basic CBF map in original units (ImageIO) - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO) - 'att': Arterial transit time in ms (ImageIO) - 't1csfgm': T1 csf-gray matter exchange time in ms (ImageIO) All maps are smoothed if smoothing is enabled.

Examples:

Basic multi-TE ASL analysis:

>>> from asltk.asldata import ASLData
>>> from asltk.reconstruction import UltraLongTE_ASLMapping
>>> from asltk.utils.io import ImageIO
>>> import numpy as np
>>> # Load multi-TE ASL data
>>> asl_data = ASLData(
...     pcasl='./tests/files/pcasl_mte.nii.gz',
...     m0='./tests/files/m0.nii.gz',
...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
...     ld_values=[1.8, 1.8, 1.8],
...     pld_values=[0.8, 1.8, 2.8]
... )
>>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
>>> # Set brain mask for faster processing
>>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
>>> ulte_mapper.set_brain_mask(brain_mask)
>>> # Generate all maps
>>> results = ulte_mapper.create_map()

Custom parameters for specific analysis:

>>> # For expected shorter T1csfGM values (faster exchange)
>>> results = ulte_mapper.create_map(
...     ub=[600.0],        # Lower upper bound
...     lb=[50.0],         # Minimum realistic T1
...     par0=[300.0]       # Lower initial guess
... )

Apply spatial smoothing to improve SNR:

>>> # Gaussian smoothing with default sigma=1.0
>>> results_smooth = ulte_mapper.create_map(
...     smoothing='gaussian'
... )
>>> # Custom smoothing parameters
>>> results_custom = ulte_mapper.create_map(
...     smoothing='gaussian',
...     smoothing_params={'sigma': 1.5}
... )
>>> # Median filtering for edge preservation
>>> results_median = ulte_mapper.create_map(
...     smoothing='median',
...     smoothing_params={'size': 5}
... )

Raises:

Type Description
ValueError

If cores parameter is invalid or required data is missing.

See Also

set_cbf_map(): Provide pre-computed CBF map set_att_map(): Provide pre-computed ATT map CBFMapping: For basic CBF/ATT mapping

Source code in asltk/reconstruction/ultralong_te_mapping.py
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
def create_map(
    self,
    ub: list = [np.inf],
    lb: list = [0.0],
    par0: list = [80000],
    cores: Union[int, str] = 'auto',
    smoothing=None,
    smoothing_params=None,
    suppress_warnings=True,
):
    """Create ultra-long-TE ASL maps including T1 csf-gray matter exchange (T1csfGM).

    This method performs advanced multi-echo ASL analysis to generate tissue-specific
    T1 relaxation maps that characterize blood-to-gray matter water exchange. The
    analysis uses multiple echo times to separate blood and tissue signal contributions.

    Note:
        The method implements the multi-compartment TE ASL model described in:
        "Ultra-long-TE arterial spin labeling reveals rapid and brain-wide
        blood-to-CSF water transport in humans", NeuroImage, 2022.
        doi: [10.1016/j.neuroimage.2021.118755](http://doi.org/10.1016/j.neuroimage.2021.118755)

    Note:
        The CBF and ATT maps can be provided before calling this method,
        using the set_cbf_map() and set_att_map() methods. If not provided,
        basic CBF/ATT maps are automatically calculated using the CBFMapping class.

    Note:
        The CBF map must be in original scale (not normalized) to perform the
        correct ultralong-TE-ASL model fitting. Use the 'cbf' output from CBFMapping,
        not the 'cbf_norm' version.

    The method assumes the T1csfGM values are well-characterized by the initial
    guess parameter. Results are filtered to include only positive values and
    values below 4 times the initial guess to remove unrealistic outliers.

    Note:
        Consider applying spatial smoothing to the output T1csfGM map to improve
        SNR. The create_map() method does not apply filtering by default to
        preserve spatial resolution.

    Args:
        ub (list, optional): Upper bounds for T1csfGM fitting. Defaults to [np.inf].
            Typically 800-1200 ms for healthy gray matter at 3T.
        lb (list, optional): Lower bounds for T1csfGM fitting. Defaults to [0.0].
            Should be positive for realistic T1 values.
        par0 (list, optional): Initial guess for T1csfGM in milliseconds.
            Defaults to [400]. Good starting values: 300-500 ms.
        cores (int or str, optional): Number of CPU threads for parallel processing.
            If "auto" (default), automatically determines the optimal number based on
            available system memory. If an integer is provided, uses that specific number.
        smoothing (str, optional): Type of spatial smoothing filter to apply.
            Options: None (default, no smoothing), 'gaussian', 'median'.
            Smoothing is applied to all output maps after reconstruction.
        smoothing_params (dict, optional): Parameters for the smoothing filter.
            For 'gaussian': {'sigma': float} (default: 1.0)
            For 'median': {'size': int} (default: 3)
        suppress_warnings (bool, optional): Whether to suppress warnings during
            processing. Defaults to True.

    Returns:
        dict: Dictionary containing:
            - 'cbf': Basic CBF map in original units (ImageIO)
            - 'cbf_norm': Normalized CBF in mL/100g/min (ImageIO)
            - 'att': Arterial transit time in ms (ImageIO)
            - 't1csfgm': T1 csf-gray matter exchange time in ms (ImageIO)
            All maps are smoothed if smoothing is enabled.

    Examples:
        Basic multi-TE ASL analysis:
        >>> from asltk.asldata import ASLData
        >>> from asltk.reconstruction import UltraLongTE_ASLMapping
        >>> from asltk.utils.io import ImageIO
        >>> import numpy as np
        >>> # Load multi-TE ASL data
        >>> asl_data = ASLData(
        ...     pcasl='./tests/files/pcasl_mte.nii.gz',
        ...     m0='./tests/files/m0.nii.gz',
        ...     te_values=[13.2, 25.7, 50.4],  # Multiple echo times
        ...     ld_values=[1.8, 1.8, 1.8],
        ...     pld_values=[0.8, 1.8, 2.8]
        ... )
        >>> ulte_mapper = UltraLongTE_ASLMapping(asl_data)
        >>> # Set brain mask for faster processing
        >>> brain_mask = ImageIO(image_array=np.ones(asl_data('m0').get_as_numpy().shape))
        >>> ulte_mapper.set_brain_mask(brain_mask)
        >>> # Generate all maps
        >>> results = ulte_mapper.create_map() # doctest: +SKIP


        Custom parameters for specific analysis:
        >>> # For expected shorter T1csfGM values (faster exchange)
        >>> results = ulte_mapper.create_map(
        ...     ub=[600.0],        # Lower upper bound
        ...     lb=[50.0],         # Minimum realistic T1
        ...     par0=[300.0]       # Lower initial guess
        ... ) # doctest: +SKIP

        Apply spatial smoothing to improve SNR:
        >>> # Gaussian smoothing with default sigma=1.0
        >>> results_smooth = ulte_mapper.create_map(
        ...     smoothing='gaussian'
        ... ) # doctest: +SKIP
        >>> # Custom smoothing parameters
        >>> results_custom = ulte_mapper.create_map(
        ...     smoothing='gaussian',
        ...     smoothing_params={'sigma': 1.5}
        ... ) # doctest: +SKIP
        >>> # Median filtering for edge preservation
        >>> results_median = ulte_mapper.create_map(
        ...     smoothing='median',
        ...     smoothing_params={'size': 5}
        ... ) # doctest: +SKIP

    Raises:
        ValueError: If cores parameter is invalid or required data is missing.

    See Also:
        set_cbf_map(): Provide pre-computed CBF map
        set_att_map(): Provide pre-computed ATT map
        CBFMapping: For basic CBF/ATT mapping
    """
    # Determine optimal number of cores based on available memory
    actual_cores = get_optimal_core_count(cores)

    # Use context manager to suppress warnings if requested
    with warnings.catch_warnings():
        if suppress_warnings:
            # Filter common warnings that might appear during fitting and processing
            warnings.filterwarnings('ignore', category=RuntimeWarning)
            warnings.filterwarnings('ignore', category=UserWarning)
            warnings.filterwarnings(
                'ignore', category=np.VisibleDeprecationWarning
            )

        self._basic_maps.set_brain_mask(
            ImageIO(image_array=self._brain_mask)
        )

        basic_maps = {'cbf': self._cbf_map, 'att': self._att_map}
        if np.mean(self._cbf_map) == 0 or np.mean(self._att_map) == 0:
            # If the CBF/ATT maps are zero (empty), then a new one is created
            print(
                '[blue][INFO] The CBF/ATT map were not provided. Creating these maps before next step...'
            )
            basic_maps = self._basic_maps.create_map()
            self._cbf_map = basic_maps['cbf'].get_as_numpy()
            self._att_map = basic_maps['att'].get_as_numpy()

        global asl_data, brain_mask, cbf_map, att_map, t2bl, t2gm
        asl_data = self._asl_data
        brain_mask = self._brain_mask
        cbf_map = self._cbf_map
        att_map = self._att_map
        ld_arr = self._asl_data.get_ld()
        pld_arr = self._asl_data.get_pld()
        te_arr = self._asl_data.get_te()
        t2bl = self.T2bl
        t2gm = self.T2gm

        x_axis = self._asl_data('m0').get_as_numpy().shape[2]   # height
        y_axis = self._asl_data('m0').get_as_numpy().shape[1]   # width
        z_axis = self._asl_data('m0').get_as_numpy().shape[0]   # depth

        tcsfgm_map_shared = Array(
            'd', z_axis * y_axis * x_axis, lock=False
        )

        # Make a copy of base information
        m0_array = asl_data('m0').get_as_numpy()
        pcasl_array = asl_data('pcasl').get_as_numpy()

        with Pool(
            processes=actual_cores,
            initializer=_multite_init_globals,
            initargs=(
                cbf_map,
                att_map,
                brain_mask,
                asl_data,
                ld_arr,
                pld_arr,
                te_arr,
                tcsfgm_map_shared,
                t2bl,
                t2gm,
            ),
        ) as pool:
            with Progress() as progress:
                task = progress.add_task(
                    'ultralongTE-ASL processing...', total=x_axis
                )
                results = [
                    pool.apply_async(
                        _tcsfgm_multite_process_slice,
                        args=(
                            i,
                            x_axis,
                            y_axis,
                            z_axis,
                            par0,
                            lb,
                            ub,
                            m0_array,
                            pcasl_array,
                        ),
                        callback=lambda _: progress.update(
                            task, advance=1
                        ),
                    )
                    for i in range(x_axis)
                ]
                for result in results:
                    result.wait()

        self._t1csfgm_map = np.frombuffer(tcsfgm_map_shared).reshape(
            z_axis, y_axis, x_axis
        )

        # Adjusting output image boundaries
        # TODO O ADJUST_LIMITS TEM QUE SER TESTADO E AJUSTADO PARA O CASO DO CSF-GM --> OK para a nova forma de adjust_image. Fazer a mesma coisa para multiTE e commitar o codigo
        self._t1csfgm_map = self._adjust_image_limits(
            self._t1csfgm_map, par0[0]
        )

        # Prepare output maps
        cbf_map_image = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_image.update_image_data(self._cbf_map)

        cbf_map_norm_image = ImageIO(self._asl_data('m0').get_image_path())
        cbf_map_norm_image.update_image_data(
            self._cbf_map * (60 * 60 * 1000)
        )

        att_map_image = ImageIO(self._asl_data('m0').get_image_path())
        att_map_image.update_image_data(self._att_map)

        t1csfgm_map_image = ImageIO(self._asl_data('m0').get_image_path())
        t1csfgm_map_image.update_image_data(self._t1csfgm_map)

        # Create output maps dictionary
        output_maps = {
            'cbf': cbf_map_image,
            'cbf_norm': cbf_map_norm_image,
            'att': att_map_image,
            't1csfgm': t1csfgm_map_image,
        }

        # Apply smoothing if requested
        return _apply_smoothing_to_maps(
            output_maps, smoothing, smoothing_params
        )

get_att_map()

Get the ATT map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The ATT map that is storaged in the

MultiTE_ASLMapping object

Source code in asltk/reconstruction/ultralong_te_mapping.py
192
193
194
195
196
197
198
199
def get_att_map(self):
    """Get the ATT map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The ATT map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._att_map

get_brain_mask()

Get the brain mask image

Returns:

Type Description
ImageIO

The brain mask image

Source code in asltk/reconstruction/ultralong_te_mapping.py
154
155
156
157
158
159
160
def get_brain_mask(self):
    """Get the brain mask image

    Returns:
        (ImageIO): The brain mask image
    """
    return self._brain_mask

get_cbf_map()

Get the CBF map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The CBF map that is storaged in the

ndarray

MultiTE_ASLMapping object

Source code in asltk/reconstruction/ultralong_te_mapping.py
175
176
177
178
179
180
181
182
def get_cbf_map(self) -> np.ndarray:
    """Get the CBF map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The CBF map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._cbf_map

get_t1csfgm_map()

Get the T1csfGM map storaged at the MultiTE_ASLMapping object

Returns:

Type Description
ImageIO

The T1csfGM map that is storaged in the

MultiTE_ASLMapping object

Source code in asltk/reconstruction/ultralong_te_mapping.py
201
202
203
204
205
206
207
208
def get_t1csfgm_map(self):
    """Get the T1csfGM map storaged at the MultiTE_ASLMapping object

    Returns:
        (ImageIO): The T1csfGM map that is storaged in the
        MultiTE_ASLMapping object
    """
    return self._t1csfgm_map

set_att_map(att_map)

Set the ATT map to the MultiTE_ASLMapping object.

Parameters:

Name Type Description Default
att_map ImageIO

The ATT map that is set in the MultiTE_ASLMapping object

required
Source code in asltk/reconstruction/ultralong_te_mapping.py
184
185
186
187
188
189
190
def set_att_map(self, att_map: ImageIO):
    """Set the ATT map to the MultiTE_ASLMapping object.

    Args:
        att_map (ImageIO): The ATT map that is set in the MultiTE_ASLMapping object
    """
    self._att_map = att_map.get_as_numpy()

set_brain_mask(brain_mask, label=1)

Defines whether a brain a mask is applied to the CBFMapping calculation

A image mask is simply an image that defines the voxels where the ASL calculation should be made. Basically any integer value can be used as proper label mask.

A most common approach is to use a binary image (zeros for background and 1 for the brain tissues). Anyway, the default behavior of the method can transform a integer-pixel values image to a binary mask with the label parameter provided by the user

Parameters:

Name Type Description Default
brain_mask ndarray

The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.

required
Source code in asltk/reconstruction/ultralong_te_mapping.py
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
def set_brain_mask(self, brain_mask: ImageIO, label: int = 1):
    """Defines whether a brain a mask is applied to the CBFMapping
    calculation

    A image mask is simply an image that defines the voxels where the ASL
    calculation should be made. Basically any integer value can be used as
    proper label mask.

    A most common approach is to use a binary image (zeros for background
    and 1 for the brain tissues). Anyway, the default behavior of the
    method can transform a integer-pixel values image to a binary mask with
    the `label` parameter provided by the user

    Args:
        brain_mask (np.ndarray): The image representing the brain mask label (int, optional): The label value used to define the foreground tissue (brain). Defaults to 1.
    """
    if not isinstance(brain_mask, ImageIO):
        raise TypeError(
            'The brain_mask parameter must be an instance of ImageIO.'
        )

    _check_mask_values(
        brain_mask, label, self._asl_data('m0').get_as_numpy().shape
    )

    binary_mask = (brain_mask.get_as_numpy() == label).astype(
        np.uint8
    ) * label
    self._brain_mask = binary_mask

set_cbf_map(cbf_map)

Set the CBF map to the MultiTE_ASLMapping object.

Note

The CBF maps must have the original scale in order to calculate the T1 CSF-GM map correclty. Hence, if the CBF map was made using CBFMapping class, one can use the 'cbf' output.

Parameters:

Name Type Description Default
cbf_map ImageIO

The CBF map that is set in the MultiTE_ASLMapping object

required
Source code in asltk/reconstruction/ultralong_te_mapping.py
162
163
164
165
166
167
168
169
170
171
172
173
def set_cbf_map(self, cbf_map: ImageIO):
    """Set the CBF map to the MultiTE_ASLMapping object.

    Note:
        The CBF maps must have the original scale in order to calculate the
        T1 CSF-GM map correclty. Hence, if the CBF map was made using
        CBFMapping class, one can use the 'cbf' output.

    Args:
        cbf_map (ImageIO): The CBF map that is set in the MultiTE_ASLMapping object
    """
    self._cbf_map = cbf_map.get_as_numpy()