scikit-image

1
import inspect
2
import sys
3
from functools import wraps
4
from math import atan2
5
from math import pi as PI
6
from math import sqrt
7
from warnings import warn
8

9
import numpy as np
10
from scipy import ndimage as ndi
11
from scipy.spatial.distance import pdist
12

13
from . import _moments
14
from ._find_contours import find_contours
15
from ._marching_cubes_lewiner import marching_cubes
16
from ._regionprops_utils import (
17
    euler_number,
18
    perimeter,
19
    perimeter_crofton,
20
    _normalize_spacing,
21
)
22

23
__all__ = ['regionprops', 'euler_number', 'perimeter', 'perimeter_crofton']
24

25

26
# All values in this PROPS dict correspond to current scikit-image property
27
# names. The keys in this PROPS dict correspond to older names used in prior
28
# releases. For backwards compatibility, these older names will continue to
29
# work, but will not be documented.
30
PROPS = {
31
    'Area': 'area',
32
    'BoundingBox': 'bbox',
33
    'BoundingBoxArea': 'area_bbox',
34
    'bbox_area': 'area_bbox',
35
    'CentralMoments': 'moments_central',
36
    'Centroid': 'centroid',
37
    'ConvexArea': 'area_convex',
38
    'convex_area': 'area_convex',
39
    # 'ConvexHull',
40
    'ConvexImage': 'image_convex',
41
    'convex_image': 'image_convex',
42
    'Coordinates': 'coords',
43
    'Eccentricity': 'eccentricity',
44
    'EquivDiameter': 'equivalent_diameter_area',
45
    'equivalent_diameter': 'equivalent_diameter_area',
46
    'EulerNumber': 'euler_number',
47
    'Extent': 'extent',
48
    # 'Extrema',
49
    'FeretDiameter': 'feret_diameter_max',
50
    'FeretDiameterMax': 'feret_diameter_max',
51
    'FilledArea': 'area_filled',
52
    'filled_area': 'area_filled',
53
    'FilledImage': 'image_filled',
54
    'filled_image': 'image_filled',
55
    'HuMoments': 'moments_hu',
56
    'Image': 'image',
57
    'InertiaTensor': 'inertia_tensor',
58
    'InertiaTensorEigvals': 'inertia_tensor_eigvals',
59
    'IntensityImage': 'image_intensity',
60
    'intensity_image': 'image_intensity',
61
    'Label': 'label',
62
    'LocalCentroid': 'centroid_local',
63
    'local_centroid': 'centroid_local',
64
    'MajorAxisLength': 'axis_major_length',
65
    'major_axis_length': 'axis_major_length',
66
    'MaxIntensity': 'intensity_max',
67
    'max_intensity': 'intensity_max',
68
    'MeanIntensity': 'intensity_mean',
69
    'mean_intensity': 'intensity_mean',
70
    'MinIntensity': 'intensity_min',
71
    'min_intensity': 'intensity_min',
72
    'std_intensity': 'intensity_std',
73
    'MinorAxisLength': 'axis_minor_length',
74
    'minor_axis_length': 'axis_minor_length',
75
    'Moments': 'moments',
76
    'NormalizedMoments': 'moments_normalized',
77
    'Orientation': 'orientation',
78
    'Perimeter': 'perimeter',
79
    'CroftonPerimeter': 'perimeter_crofton',
80
    # 'PixelIdxList',
81
    # 'PixelList',
82
    'Slice': 'slice',
83
    'Solidity': 'solidity',
84
    # 'SubarrayIdx'
85
    'WeightedCentralMoments': 'moments_weighted_central',
86
    'weighted_moments_central': 'moments_weighted_central',
87
    'WeightedCentroid': 'centroid_weighted',
88
    'weighted_centroid': 'centroid_weighted',
89
    'WeightedHuMoments': 'moments_weighted_hu',
90
    'weighted_moments_hu': 'moments_weighted_hu',
91
    'WeightedLocalCentroid': 'centroid_weighted_local',
92
    'weighted_local_centroid': 'centroid_weighted_local',
93
    'WeightedMoments': 'moments_weighted',
94
    'weighted_moments': 'moments_weighted',
95
    'WeightedNormalizedMoments': 'moments_weighted_normalized',
96
    'weighted_moments_normalized': 'moments_weighted_normalized',
97
}
98

99
COL_DTYPES = {
100
    'area': float,
101
    'area_bbox': float,
102
    'area_convex': float,
103
    'area_filled': float,
104
    'axis_major_length': float,
105
    'axis_minor_length': float,
106
    'bbox': int,
107
    'centroid': float,
108
    'centroid_local': float,
109
    'centroid_weighted': float,
110
    'centroid_weighted_local': float,
111
    'coords': object,
112
    'coords_scaled': object,
113
    'eccentricity': float,
114
    'equivalent_diameter_area': float,
115
    'euler_number': int,
116
    'extent': float,
117
    'feret_diameter_max': float,
118
    'image': object,
119
    'image_convex': object,
120
    'image_filled': object,
121
    'image_intensity': object,
122
    'inertia_tensor': float,
123
    'inertia_tensor_eigvals': float,
124
    'intensity_max': float,
125
    'intensity_mean': float,
126
    'intensity_min': float,
127
    'intensity_std': float,
128
    'label': int,
129
    'moments': float,
130
    'moments_central': float,
131
    'moments_hu': float,
132
    'moments_normalized': float,
133
    'moments_weighted': float,
134
    'moments_weighted_central': float,
135
    'moments_weighted_hu': float,
136
    'moments_weighted_normalized': float,
137
    'num_pixels': int,
138
    'orientation': float,
139
    'perimeter': float,
140
    'perimeter_crofton': float,
141
    'slice': object,
142
    'solidity': float,
143
}
144

145
OBJECT_COLUMNS = [col for col, dtype in COL_DTYPES.items() if dtype == object]
146

147
PROP_VALS = set(PROPS.values())
148

149
_require_intensity_image = (
150
    'image_intensity',
151
    'intensity_max',
152
    'intensity_mean',
153
    'intensity_min',
154
    'intensity_std',
155
    'moments_weighted',
156
    'moments_weighted_central',
157
    'centroid_weighted',
158
    'centroid_weighted_local',
159
    'moments_weighted_hu',
160
    'moments_weighted_normalized',
161
)
162

163

164
def _infer_number_of_required_args(func):
165
    """Infer the number of required arguments for a function
166

167
    Parameters
168
    ----------
169
    func : callable
170
        The function that is being inspected.
171

172
    Returns
173
    -------
174
    n_args : int
175
        The number of required arguments of func.
176
    """
177
    argspec = inspect.getfullargspec(func)
178
    n_args = len(argspec.args)
179
    if argspec.defaults is not None:
180
        n_args -= len(argspec.defaults)
181
    return n_args
182

183

184
def _infer_regionprop_dtype(func, *, intensity, ndim):
185
    """Infer the dtype of a region property calculated by func.
186

187
    If a region property function always returns the same shape and type of
188
    output regardless of input size, then the dtype is the dtype of the
189
    returned array. Otherwise, the property has object dtype.
190

191
    Parameters
192
    ----------
193
    func : callable
194
        Function to be tested. The signature should be array[bool] -> Any if
195
        intensity is False, or *(array[bool], array[float]) -> Any otherwise.
196
    intensity : bool
197
        Whether the regionprop is calculated on an intensity image.
198
    ndim : int
199
        The number of dimensions for which to check func.
200

201
    Returns
202
    -------
203
    dtype : NumPy data type
204
        The data type of the returned property.
205
    """
206
    mask_1 = np.ones((1,) * ndim, dtype=bool)
207
    mask_1 = np.pad(mask_1, (0, 1), constant_values=False)
208
    mask_2 = np.ones((2,) * ndim, dtype=bool)
209
    mask_2 = np.pad(mask_2, (1, 0), constant_values=False)
210
    propmasks = [mask_1, mask_2]
211

212
    rng = np.random.default_rng()
213

214
    if intensity and _infer_number_of_required_args(func) == 2:
215

216
        def _func(mask):
217
            return func(mask, rng.random(mask.shape))
218

219
    else:
220
        _func = func
221
    props1, props2 = map(_func, propmasks)
222
    if (
223
        np.isscalar(props1)
224
        and np.isscalar(props2)
225
        or np.array(props1).shape == np.array(props2).shape
226
    ):
227
        dtype = np.array(props1).dtype.type
228
    else:
229
        dtype = np.object_
230
    return dtype
231

232

233
def _cached(f):
234
    @wraps(f)
235
    def wrapper(obj):
236
        cache = obj._cache
237
        prop = f.__name__
238

239
        if not obj._cache_active:
240
            return f(obj)
241

242
        if prop not in cache:
243
            cache[prop] = f(obj)
244

245
        return cache[prop]
246

247
    return wrapper
248

249

250
def only2d(method):
251
    @wraps(method)
252
    def func2d(self, *args, **kwargs):
253
        if self._ndim > 2:
254
            raise NotImplementedError(
255
                f"Property {method.__name__} is not implemented for 3D images"
256
            )
257
        return method(self, *args, **kwargs)
258

259
    return func2d
260

261

262
def _inertia_eigvals_to_axes_lengths_3D(inertia_tensor_eigvals):
263
    """Compute ellipsoid axis lengths from inertia tensor eigenvalues.
264

265
    Parameters
266
    ---------
267
    inertia_tensor_eigvals : sequence of float
268
        A sequence of 3 floating point eigenvalues, sorted in descending order.
269

270
    Returns
271
    -------
272
    axis_lengths : list of float
273
        The ellipsoid axis lengths sorted in descending order.
274

275
    Notes
276
    -----
277
    Let a >= b >= c be the ellipsoid semi-axes and s1 >= s2 >= s3 be the
278
    inertia tensor eigenvalues.
279

280
    The inertia tensor eigenvalues are given for a solid ellipsoid in [1]_.
281
    s1 = 1 / 5 * (a**2 + b**2)
282
    s2 = 1 / 5 * (a**2 + c**2)
283
    s3 = 1 / 5 * (b**2 + c**2)
284

285
    Rearranging to solve for a, b, c in terms of s1, s2, s3 gives
286
    a = math.sqrt(5 / 2 * ( s1 + s2 - s3))
287
    b = math.sqrt(5 / 2 * ( s1 - s2 + s3))
288
    c = math.sqrt(5 / 2 * (-s1 + s2 + s3))
289

290
    We can then simply replace sqrt(5/2) by sqrt(10) to get the full axes
291
    lengths rather than the semi-axes lengths.
292

293
    References
294
    ----------
295
    ..[1] https://en.wikipedia.org/wiki/List_of_moments_of_inertia#List_of_3D_inertia_tensors  # noqa
296
    """
297
    axis_lengths = []
298
    for ax in range(2, -1, -1):
299
        w = sum(v * -1 if i == ax else v for i, v in enumerate(inertia_tensor_eigvals))
300
        axis_lengths.append(sqrt(10 * w))
301
    return axis_lengths
302

303

304
class RegionProperties:
305
    """Please refer to `skimage.measure.regionprops` for more information
306
    on the available region properties.
307
    """
308

309
    def __init__(
310
        self,
311
        slice,
312
        label,
313
        label_image,
314
        intensity_image,
315
        cache_active,
316
        *,
317
        extra_properties=None,
318
        spacing=None,
319
        offset=None,
320
    ):
321
        if intensity_image is not None:
322
            ndim = label_image.ndim
323
            if not (
324
                intensity_image.shape[:ndim] == label_image.shape
325
                and intensity_image.ndim in [ndim, ndim + 1]
326
            ):
327
                raise ValueError(
328
                    'Label and intensity image shapes must match,'
329
                    ' except for channel (last) axis.'
330
                )
331
            multichannel = label_image.shape < intensity_image.shape
332
        else:
333
            multichannel = False
334

335
        self.label = label
336
        if offset is None:
337
            offset = np.zeros((label_image.ndim,), dtype=int)
338
        self._offset = np.array(offset)
339

340
        self._slice = slice
341
        self.slice = slice
342
        self._label_image = label_image
343
        self._intensity_image = intensity_image
344

345
        self._cache_active = cache_active
346
        self._cache = {}
347
        self._ndim = label_image.ndim
348
        self._multichannel = multichannel
349
        self._spatial_axes = tuple(range(self._ndim))
350
        if spacing is None:
351
            spacing = np.full(self._ndim, 1.0)
352
        self._spacing = _normalize_spacing(spacing, self._ndim)
353
        self._pixel_area = np.prod(self._spacing)
354

355
        self._extra_properties = {}
356
        if extra_properties is not None:
357
            for func in extra_properties:
358
                name = func.__name__
359
                if hasattr(self, name):
360
                    msg = (
361
                        f"Extra property '{name}' is shadowed by existing "
362
                        f"property and will be inaccessible. Consider "
363
                        f"renaming it."
364
                    )
365
                    warn(msg)
366
            self._extra_properties = {func.__name__: func for func in extra_properties}
367

368
    def __getattr__(self, attr):
369
        if self._intensity_image is None and attr in _require_intensity_image:
370
            raise AttributeError(
371
                f"Attribute '{attr}' unavailable when `intensity_image` "
372
                f"has not been specified."
373
            )
374
        if attr in self._extra_properties:
375
            func = self._extra_properties[attr]
376
            n_args = _infer_number_of_required_args(func)
377
            # determine whether func requires intensity image
378
            if n_args == 2:
379
                if self._intensity_image is not None:
380
                    if self._multichannel:
381
                        multichannel_list = [
382
                            func(self.image, self.image_intensity[..., i])
383
                            for i in range(self.image_intensity.shape[-1])
384
                        ]
385
                        return np.stack(multichannel_list, axis=-1)
386
                    else:
387
                        return func(self.image, self.image_intensity)
388
                else:
389
                    raise AttributeError(
390
                        f'intensity image required to calculate {attr}'
391
                    )
392
            elif n_args == 1:
393
                return func(self.image)
394
            else:
395
                raise AttributeError(
396
                    f'Custom regionprop function\'s number of arguments must '
397
                    f'be 1 or 2, but {attr} takes {n_args} arguments.'
398
                )
399
        elif attr in PROPS and attr.lower() == attr:
400
            if (
401
                self._intensity_image is None
402
                and PROPS[attr] in _require_intensity_image
403
            ):
404
                raise AttributeError(
405
                    f"Attribute '{attr}' unavailable when `intensity_image` "
406
                    f"has not been specified."
407
                )
408
            # retrieve deprecated property (excluding old CamelCase ones)
409
            return getattr(self, PROPS[attr])
410
        else:
411
            raise AttributeError(f"'{type(self)}' object has no attribute '{attr}'")
412

413
    def __setattr__(self, name, value):
414
        if name in PROPS:
415
            super().__setattr__(PROPS[name], value)
416
        else:
417
            super().__setattr__(name, value)
418

419
    @property
420
    @_cached
421
    def num_pixels(self):
422
        return np.sum(self.image)
423

424
    @property
425
    @_cached
426
    def area(self):
427
        return np.sum(self.image) * self._pixel_area
428

429
    @property
430
    def bbox(self):
431
        """
432
        Returns
433
        -------
434
        A tuple of the bounding box's start coordinates for each dimension,
435
        followed by the end coordinates for each dimension
436
        """
437
        return tuple(
438
            [self.slice[i].start for i in range(self._ndim)]
439
            + [self.slice[i].stop for i in range(self._ndim)]
440
        )
441

442
    @property
443
    def area_bbox(self):
444
        return self.image.size * self._pixel_area
445

446
    @property
447
    def centroid(self):
448
        return tuple(self.coords_scaled.mean(axis=0))
449

450
    @property
451
    @_cached
452
    def area_convex(self):
453
        return np.sum(self.image_convex) * self._pixel_area
454

455
    @property
456
    @_cached
457
    def image_convex(self):
458
        from ..morphology.convex_hull import convex_hull_image
459

460
        return convex_hull_image(self.image)
461

462
    @property
463
    def coords_scaled(self):
464
        indices = np.argwhere(self.image)
465
        object_offset = np.array([self.slice[i].start for i in range(self._ndim)])
466
        return (object_offset + indices) * self._spacing + self._offset
467

468
    @property
469
    def coords(self):
470
        indices = np.argwhere(self.image)
471
        object_offset = np.array([self.slice[i].start for i in range(self._ndim)])
472
        return object_offset + indices + self._offset
473

474
    @property
475
    @only2d
476
    def eccentricity(self):
477
        l1, l2 = self.inertia_tensor_eigvals
478
        if l1 == 0:
479
            return 0
480
        return sqrt(1 - l2 / l1)
481

482
    @property
483
    def equivalent_diameter_area(self):
484
        return (2 * self._ndim * self.area / PI) ** (1 / self._ndim)
485

486
    @property
487
    def euler_number(self):
488
        if self._ndim not in [2, 3]:
489
            raise NotImplementedError(
490
                'Euler number is implemented for ' '2D or 3D images only'
491
            )
492
        return euler_number(self.image, self._ndim)
493

494
    @property
495
    def extent(self):
496
        return self.area / self.area_bbox
497

498
    @property
499
    def feret_diameter_max(self):
500
        identity_convex_hull = np.pad(
501
            self.image_convex, 2, mode='constant', constant_values=0
502
        )
503
        if self._ndim == 2:
504
            coordinates = np.vstack(
505
                find_contours(identity_convex_hull, 0.5, fully_connected='high')
506
            )
507
        elif self._ndim == 3:
508
            coordinates, _, _, _ = marching_cubes(identity_convex_hull, level=0.5)
509
        distances = pdist(coordinates * self._spacing, 'sqeuclidean')
510
        return sqrt(np.max(distances))
511

512
    @property
513
    def area_filled(self):
514
        return np.sum(self.image_filled) * self._pixel_area
515

516
    @property
517
    @_cached
518
    def image_filled(self):
519
        structure = np.ones((3,) * self._ndim)
520
        return ndi.binary_fill_holes(self.image, structure)
521

522
    @property
523
    @_cached
524
    def image(self):
525
        return self._label_image[self.slice] == self.label
526

527
    @property
528
    @_cached
529
    def inertia_tensor(self):
530
        mu = self.moments_central
531
        return _moments.inertia_tensor(self.image, mu, spacing=self._spacing)
532

533
    @property
534
    @_cached
535
    def inertia_tensor_eigvals(self):
536
        return _moments.inertia_tensor_eigvals(self.image, T=self.inertia_tensor)
537

538
    @property
539
    @_cached
540
    def image_intensity(self):
541
        if self._intensity_image is None:
542
            raise AttributeError('No intensity image specified.')
543
        image = (
544
            self.image
545
            if not self._multichannel
546
            else np.expand_dims(self.image, self._ndim)
547
        )
548
        return self._intensity_image[self.slice] * image
549

550
    def _image_intensity_double(self):
551
        return self.image_intensity.astype(np.float64, copy=False)
552

553
    @property
554
    def centroid_local(self):
555
        M = self.moments
556
        M0 = M[(0,) * self._ndim]
557

558
        def _get_element(axis):
559
            return (0,) * axis + (1,) + (0,) * (self._ndim - 1 - axis)
560

561
        return np.asarray(
562
            tuple(M[_get_element(axis)] / M0 for axis in range(self._ndim))
563
        )
564

565
    @property
566
    def intensity_max(self):
567
        vals = self.image_intensity[self.image]
568
        return np.max(vals, axis=0).astype(np.float64, copy=False)
569

570
    @property
571
    def intensity_mean(self):
572
        return np.mean(self.image_intensity[self.image], axis=0)
573

574
    @property
575
    def intensity_min(self):
576
        vals = self.image_intensity[self.image]
577
        return np.min(vals, axis=0).astype(np.float64, copy=False)
578

579
    @property
580
    def intensity_std(self):
581
        vals = self.image_intensity[self.image]
582
        return np.std(vals, axis=0)
583

584
    @property
585
    def axis_major_length(self):
586
        if self._ndim == 2:
587
            l1 = self.inertia_tensor_eigvals[0]
588
            return 4 * sqrt(l1)
589
        elif self._ndim == 3:
590
            # equivalent to _inertia_eigvals_to_axes_lengths_3D(ev)[0]
591
            ev = self.inertia_tensor_eigvals
592
            return sqrt(10 * (ev[0] + ev[1] - ev[2]))
593
        else:
594
            raise ValueError("axis_major_length only available in 2D and 3D")
595

596
    @property
597
    def axis_minor_length(self):
598
        if self._ndim == 2:
599
            l2 = self.inertia_tensor_eigvals[-1]
600
            return 4 * sqrt(l2)
601
        elif self._ndim == 3:
602
            # equivalent to _inertia_eigvals_to_axes_lengths_3D(ev)[-1]
603
            ev = self.inertia_tensor_eigvals
604
            return sqrt(10 * (-ev[0] + ev[1] + ev[2]))
605
        else:
606
            raise ValueError("axis_minor_length only available in 2D and 3D")
607

608
    @property
609
    @_cached
610
    def moments(self):
611
        M = _moments.moments(self.image.astype(np.uint8), 3, spacing=self._spacing)
612
        return M
613

614
    @property
615
    @_cached
616
    def moments_central(self):
617
        mu = _moments.moments_central(
618
            self.image.astype(np.uint8),
619
            self.centroid_local,
620
            order=3,
621
            spacing=self._spacing,
622
        )
623
        return mu
624

625
    @property
626
    @only2d
627
    def moments_hu(self):
628
        if any(s != 1.0 for s in self._spacing):
629
            raise NotImplementedError('`moments_hu` supports spacing = (1, 1) only')
630
        return _moments.moments_hu(self.moments_normalized)
631

632
    @property
633
    @_cached
634
    def moments_normalized(self):
635
        return _moments.moments_normalized(
636
            self.moments_central, 3, spacing=self._spacing
637
        )
638

639
    @property
640
    @only2d
641
    def orientation(self):
642
        a, b, b, c = self.inertia_tensor.flat
643
        if a - c == 0:
644
            if b < 0:
645
                return PI / 4.0
646
            else:
647
                return -PI / 4.0
648
        else:
649
            return 0.5 * atan2(-2 * b, c - a)
650

651
    @property
652
    @only2d
653
    def perimeter(self):
654
        if len(np.unique(self._spacing)) != 1:
655
            raise NotImplementedError('`perimeter` supports isotropic spacings only')
656
        return perimeter(self.image, 4) * self._spacing[0]
657

658
    @property
659
    @only2d
660
    def perimeter_crofton(self):
661
        if len(np.unique(self._spacing)) != 1:
662
            raise NotImplementedError('`perimeter` supports isotropic spacings only')
663
        return perimeter_crofton(self.image, 4) * self._spacing[0]
664

665
    @property
666
    def solidity(self):
667
        return self.area / self.area_convex
668

669
    @property
670
    def centroid_weighted(self):
671
        ctr = self.centroid_weighted_local
672
        return tuple(
673
            idx + slc.start * spc
674
            for idx, slc, spc in zip(ctr, self.slice, self._spacing)
675
        )
676

677
    @property
678
    def centroid_weighted_local(self):
679
        M = self.moments_weighted
680
        M0 = M[(0,) * self._ndim]
681

682
        def _get_element(axis):
683
            return (0,) * axis + (1,) + (0,) * (self._ndim - 1 - axis)
684

685
        return np.asarray(
686
            tuple(M[_get_element(axis)] / M0 for axis in range(self._ndim))
687
        )
688

689
    @property
690
    @_cached
691
    def moments_weighted(self):
692
        image = self._image_intensity_double()
693
        if self._multichannel:
694
            moments = np.stack(
695
                [
696
                    _moments.moments(image[..., i], order=3, spacing=self._spacing)
697
                    for i in range(image.shape[-1])
698
                ],
699
                axis=-1,
700
            )
701
        else:
702
            moments = _moments.moments(image, order=3, spacing=self._spacing)
703
        return moments
704

705
    @property
706
    @_cached
707
    def moments_weighted_central(self):
708
        ctr = self.centroid_weighted_local
709
        image = self._image_intensity_double()
710
        if self._multichannel:
711
            moments_list = [
712
                _moments.moments_central(
713
                    image[..., i], center=ctr[..., i], order=3, spacing=self._spacing
714
                )
715
                for i in range(image.shape[-1])
716
            ]
717
            moments = np.stack(moments_list, axis=-1)
718
        else:
719
            moments = _moments.moments_central(
720
                image, ctr, order=3, spacing=self._spacing
721
            )
722
        return moments
723

724
    @property
725
    @only2d
726
    def moments_weighted_hu(self):
727
        if not (np.array(self._spacing) == np.array([1, 1])).all():
728
            raise NotImplementedError('`moments_hu` supports spacing = (1, 1) only')
729
        nu = self.moments_weighted_normalized
730
        if self._multichannel:
731
            nchannels = self._intensity_image.shape[-1]
732
            return np.stack(
733
                [_moments.moments_hu(nu[..., i]) for i in range(nchannels)],
734
                axis=-1,
735
            )
736
        else:
737
            return _moments.moments_hu(nu)
738

739
    @property
740
    @_cached
741
    def moments_weighted_normalized(self):
742
        mu = self.moments_weighted_central
743
        if self._multichannel:
744
            nchannels = self._intensity_image.shape[-1]
745
            return np.stack(
746
                [
747
                    _moments.moments_normalized(
748
                        mu[..., i], order=3, spacing=self._spacing
749
                    )
750
                    for i in range(nchannels)
751
                ],
752
                axis=-1,
753
            )
754
        else:
755
            return _moments.moments_normalized(mu, order=3, spacing=self._spacing)
756

757
    def __iter__(self):
758
        props = PROP_VALS
759

760
        if self._intensity_image is None:
761
            unavailable_props = _require_intensity_image
762
            props = props.difference(unavailable_props)
763

764
        return iter(sorted(props))
765

766
    def __getitem__(self, key):
767
        value = getattr(self, key, None)
768
        if value is not None:
769
            return value
770
        else:  # backwards compatibility
771
            return getattr(self, PROPS[key])
772

773
    def __eq__(self, other):
774
        if not isinstance(other, RegionProperties):
775
            return False
776

777
        for key in PROP_VALS:
778
            try:
779
                # so that NaNs are equal
780
                np.testing.assert_equal(
781
                    getattr(self, key, None), getattr(other, key, None)
782
                )
783
            except AssertionError:
784
                return False
785

786
        return True
787

788

789
# For compatibility with code written prior to 0.16
790
_RegionProperties = RegionProperties
791

792

793
def _props_to_dict(regions, properties=('label', 'bbox'), separator='-'):
794
    """Convert image region properties list into a column dictionary.
795

796
    Parameters
797
    ----------
798
    regions : (K,) list
799
        List of RegionProperties objects as returned by :func:`regionprops`.
800
    properties : tuple or list of str, optional
801
        Properties that will be included in the resulting dictionary
802
        For a list of available properties, please see :func:`regionprops`.
803
        Users should remember to add "label" to keep track of region
804
        identities.
805
    separator : str, optional
806
        For non-scalar properties not listed in OBJECT_COLUMNS, each element
807
        will appear in its own column, with the index of that element separated
808
        from the property name by this separator. For example, the inertia
809
        tensor of a 2D region will appear in four columns:
810
        ``inertia_tensor-0-0``, ``inertia_tensor-0-1``, ``inertia_tensor-1-0``,
811
        and ``inertia_tensor-1-1`` (where the separator is ``-``).
812

813
        Object columns are those that cannot be split in this way because the
814
        number of columns would change depending on the object. For example,
815
        ``image`` and ``coords``.
816

817
    Returns
818
    -------
819
    out_dict : dict
820
        Dictionary mapping property names to an array of values of that
821
        property, one value per region. This dictionary can be used as input to
822
        pandas ``DataFrame`` to map property names to columns in the frame and
823
        regions to rows.
824

825
    Notes
826
    -----
827
    Each column contains either a scalar property, an object property, or an
828
    element in a multidimensional array.
829

830
    Properties with scalar values for each region, such as "eccentricity", will
831
    appear as a float or int array with that property name as key.
832

833
    Multidimensional properties *of fixed size* for a given image dimension,
834
    such as "centroid" (every centroid will have three elements in a 3D image,
835
    no matter the region size), will be split into that many columns, with the
836
    name {property_name}{separator}{element_num} (for 1D properties),
837
    {property_name}{separator}{elem_num0}{separator}{elem_num1} (for 2D
838
    properties), and so on.
839

840
    For multidimensional properties that don't have a fixed size, such as
841
    "image" (the image of a region varies in size depending on the region
842
    size), an object array will be used, with the corresponding property name
843
    as the key.
844

845
    Examples
846
    --------
847
    >>> from skimage import data, util, measure
848
    >>> image = data.coins()
849
    >>> label_image = measure.label(image > 110, connectivity=image.ndim)
850
    >>> proplist = regionprops(label_image, image)
851
    >>> props = _props_to_dict(proplist, properties=['label', 'inertia_tensor',
852
    ...                                              'inertia_tensor_eigvals'])
853
    >>> props  # doctest: +ELLIPSIS +SKIP
854
    {'label': array([ 1,  2, ...]), ...
855
     'inertia_tensor-0-0': array([  4.012...e+03,   8.51..., ...]), ...
856
     ...,
857
     'inertia_tensor_eigvals-1': array([  2.67...e+02,   2.83..., ...])}
858

859
    The resulting dictionary can be directly passed to pandas, if installed, to
860
    obtain a clean DataFrame:
861

862
    >>> import pandas as pd  # doctest: +SKIP
863
    >>> data = pd.DataFrame(props)  # doctest: +SKIP
864
    >>> data.head()  # doctest: +SKIP
865
       label  inertia_tensor-0-0  ...  inertia_tensor_eigvals-1
866
    0      1         4012.909888  ...                267.065503
867
    1      2            8.514739  ...                  2.834806
868
    2      3            0.666667  ...                  0.000000
869
    3      4            0.000000  ...                  0.000000
870
    4      5            0.222222  ...                  0.111111
871

872
    """
873

874
    out = {}
875
    n = len(regions)
876
    for prop in properties:
877
        r = regions[0]
878
        # Copy the original property name so the output will have the
879
        # user-provided property name in the case of deprecated names.
880
        orig_prop = prop
881
        # determine the current property name for any deprecated property.
882
        prop = PROPS.get(prop, prop)
883
        rp = getattr(r, prop)
884
        if prop in COL_DTYPES:
885
            dtype = COL_DTYPES[prop]
886
        else:
887
            func = r._extra_properties[prop]
888
            dtype = _infer_regionprop_dtype(
889
                func,
890
                intensity=r._intensity_image is not None,
891
                ndim=r.image.ndim,
892
            )
893

894
        # scalars and objects are dedicated one column per prop
895
        # array properties are raveled into multiple columns
896
        # for more info, refer to notes 1
897
        if np.isscalar(rp) or prop in OBJECT_COLUMNS or dtype is np.object_:
898
            column_buffer = np.empty(n, dtype=dtype)
899
            for i in range(n):
900
                column_buffer[i] = regions[i][prop]
901
            out[orig_prop] = np.copy(column_buffer)
902
        else:
903
            # precompute property column names and locations
904
            modified_props = []
905
            locs = []
906
            for ind in np.ndindex(np.shape(rp)):
907
                modified_props.append(separator.join(map(str, (orig_prop,) + ind)))
908
                locs.append(ind if len(ind) > 1 else ind[0])
909

910
            # fill temporary column data_array
911
            n_columns = len(locs)
912
            column_data = np.empty((n, n_columns), dtype=dtype)
913
            for k in range(n):
914
                # we coerce to a numpy array to ensure structures like
915
                # tuple-of-arrays expand correctly into columns
916
                rp = np.asarray(regions[k][prop])
917
                for i, loc in enumerate(locs):
918
                    column_data[k, i] = rp[loc]
919

920
            # add the columns to the output dictionary
921
            for i, modified_prop in enumerate(modified_props):
922
                out[modified_prop] = column_data[:, i]
923
    return out
924

925

926
def regionprops_table(
927
    label_image,
928
    intensity_image=None,
929
    properties=('label', 'bbox'),
930
    *,
931
    cache=True,
932
    separator='-',
933
    extra_properties=None,
934
    spacing=None,
935
):
936
    """Compute image properties and return them as a pandas-compatible table.
937

938
    The table is a dictionary mapping column names to value arrays. See Notes
939
    section below for details.
940

941
    .. versionadded:: 0.16
942

943
    Parameters
944
    ----------
945
    label_image : (M, N[, P]) ndarray
946
        Labeled input image. Labels with value 0 are ignored.
947
    intensity_image : (M, N[, P][, C]) ndarray, optional
948
        Intensity (i.e., input) image with same size as labeled image, plus
949
        optionally an extra dimension for multichannel data. The channel dimension,
950
        if present, must be the last axis. Default is None.
951

952
        .. versionchanged:: 0.18.0
953
            The ability to provide an extra dimension for channels was added.
954
    properties : tuple or list of str, optional
955
        Properties that will be included in the resulting dictionary
956
        For a list of available properties, please see :func:`regionprops`.
957
        Users should remember to add "label" to keep track of region
958
        identities.
959
    cache : bool, optional
960
        Determine whether to cache calculated properties. The computation is
961
        much faster for cached properties, whereas the memory consumption
962
        increases.
963
    separator : str, optional
964
        For non-scalar properties not listed in OBJECT_COLUMNS, each element
965
        will appear in its own column, with the index of that element separated
966
        from the property name by this separator. For example, the inertia
967
        tensor of a 2D region will appear in four columns:
968
        ``inertia_tensor-0-0``, ``inertia_tensor-0-1``, ``inertia_tensor-1-0``,
969
        and ``inertia_tensor-1-1`` (where the separator is ``-``).
970

971
        Object columns are those that cannot be split in this way because the
972
        number of columns would change depending on the object. For example,
973
        ``image`` and ``coords``.
974
    extra_properties : Iterable of callables
975
        Add extra property computation functions that are not included with
976
        skimage. The name of the property is derived from the function name,
977
        the dtype is inferred by calling the function on a small sample.
978
        If the name of an extra property clashes with the name of an existing
979
        property the extra property will not be visible and a UserWarning is
980
        issued. A property computation function must take a region mask as its
981
        first argument. If the property requires an intensity image, it must
982
        accept the intensity image as the second argument.
983
    spacing: tuple of float, shape (ndim,)
984
        The pixel spacing along each axis of the image.
985

986
    Returns
987
    -------
988
    out_dict : dict
989
        Dictionary mapping property names to an array of values of that
990
        property, one value per region. This dictionary can be used as input to
991
        pandas ``DataFrame`` to map property names to columns in the frame and
992
        regions to rows. If the image has no regions,
993
        the arrays will have length 0, but the correct type.
994

995
    Notes
996
    -----
997
    Each column contains either a scalar property, an object property, or an
998
    element in a multidimensional array.
999

1000
    Properties with scalar values for each region, such as "eccentricity", will
1001
    appear as a float or int array with that property name as key.
1002

1003
    Multidimensional properties *of fixed size* for a given image dimension,
1004
    such as "centroid" (every centroid will have three elements in a 3D image,
1005
    no matter the region size), will be split into that many columns, with the
1006
    name {property_name}{separator}{element_num} (for 1D properties),
1007
    {property_name}{separator}{elem_num0}{separator}{elem_num1} (for 2D
1008
    properties), and so on.
1009

1010
    For multidimensional properties that don't have a fixed size, such as
1011
    "image" (the image of a region varies in size depending on the region
1012
    size), an object array will be used, with the corresponding property name
1013
    as the key.
1014

1015
    Examples
1016
    --------
1017
    >>> from skimage import data, util, measure
1018
    >>> image = data.coins()
1019
    >>> label_image = measure.label(image > 110, connectivity=image.ndim)
1020
    >>> props = measure.regionprops_table(label_image, image,
1021
    ...                           properties=['label', 'inertia_tensor',
1022
    ...                                       'inertia_tensor_eigvals'])
1023
    >>> props  # doctest: +ELLIPSIS +SKIP
1024
    {'label': array([ 1,  2, ...]), ...
1025
     'inertia_tensor-0-0': array([  4.012...e+03,   8.51..., ...]), ...
1026
     ...,
1027
     'inertia_tensor_eigvals-1': array([  2.67...e+02,   2.83..., ...])}
1028

1029
    The resulting dictionary can be directly passed to pandas, if installed, to
1030
    obtain a clean DataFrame:
1031

1032
    >>> import pandas as pd  # doctest: +SKIP
1033
    >>> data = pd.DataFrame(props)  # doctest: +SKIP
1034
    >>> data.head()  # doctest: +SKIP
1035
       label  inertia_tensor-0-0  ...  inertia_tensor_eigvals-1
1036
    0      1         4012.909888  ...                267.065503
1037
    1      2            8.514739  ...                  2.834806
1038
    2      3            0.666667  ...                  0.000000
1039
    3      4            0.000000  ...                  0.000000
1040
    4      5            0.222222  ...                  0.111111
1041

1042
    [5 rows x 7 columns]
1043

1044
    If we want to measure a feature that does not come as a built-in
1045
    property, we can define custom functions and pass them as
1046
    ``extra_properties``. For example, we can create a custom function
1047
    that measures the intensity quartiles in a region:
1048

1049
    >>> from skimage import data, util, measure
1050
    >>> import numpy as np
1051
    >>> def quartiles(regionmask, intensity):
1052
    ...     return np.percentile(intensity[regionmask], q=(25, 50, 75))
1053
    >>>
1054
    >>> image = data.coins()
1055
    >>> label_image = measure.label(image > 110, connectivity=image.ndim)
1056
    >>> props = measure.regionprops_table(label_image, intensity_image=image,
1057
    ...                                   properties=('label',),
1058
    ...                                   extra_properties=(quartiles,))
1059
    >>> import pandas as pd # doctest: +SKIP
1060
    >>> pd.DataFrame(props).head() # doctest: +SKIP
1061
           label  quartiles-0  quartiles-1  quartiles-2
1062
    0      1       117.00        123.0        130.0
1063
    1      2       111.25        112.0        114.0
1064
    2      3       111.00        111.0        111.0
1065
    3      4       111.00        111.5        112.5
1066
    4      5       112.50        113.0        114.0
1067

1068
    """
1069
    regions = regionprops(
1070
        label_image,
1071
        intensity_image=intensity_image,
1072
        cache=cache,
1073
        extra_properties=extra_properties,
1074
        spacing=spacing,
1075
    )
1076
    if extra_properties is not None:
1077
        properties = list(properties) + [prop.__name__ for prop in extra_properties]
1078
    if len(regions) == 0:
1079
        ndim = label_image.ndim
1080
        label_image = np.zeros((3,) * ndim, dtype=int)
1081
        label_image[(1,) * ndim] = 1
1082
        if intensity_image is not None:
1083
            intensity_image = np.zeros(
1084
                label_image.shape + intensity_image.shape[ndim:],
1085
                dtype=intensity_image.dtype,
1086
            )
1087
        regions = regionprops(
1088
            label_image,
1089
            intensity_image=intensity_image,
1090
            cache=cache,
1091
            extra_properties=extra_properties,
1092
            spacing=spacing,
1093
        )
1094

1095
        out_d = _props_to_dict(regions, properties=properties, separator=separator)
1096
        return {k: v[:0] for k, v in out_d.items()}
1097

1098
    return _props_to_dict(regions, properties=properties, separator=separator)
1099

1100

1101
def regionprops(
1102
    label_image,
1103
    intensity_image=None,
1104
    cache=True,
1105
    *,
1106
    extra_properties=None,
1107
    spacing=None,
1108
    offset=None,
1109
):
1110
    r"""Measure properties of labeled image regions.
1111

1112
    Parameters
1113
    ----------
1114
    label_image : (M, N[, P]) ndarray
1115
        Labeled input image. Labels with value 0 are ignored.
1116

1117
        .. versionchanged:: 0.14.1
1118
            Previously, ``label_image`` was processed by ``numpy.squeeze`` and
1119
            so any number of singleton dimensions was allowed. This resulted in
1120
            inconsistent handling of images with singleton dimensions. To
1121
            recover the old behaviour, use
1122
            ``regionprops(np.squeeze(label_image), ...)``.
1123
    intensity_image : (M, N[, P][, C]) ndarray, optional
1124
        Intensity (i.e., input) image with same size as labeled image, plus
1125
        optionally an extra dimension for multichannel data. Currently,
1126
        this extra channel dimension, if present, must be the last axis.
1127
        Default is None.
1128

1129
        .. versionchanged:: 0.18.0
1130
            The ability to provide an extra dimension for channels was added.
1131
    cache : bool, optional
1132
        Determine whether to cache calculated properties. The computation is
1133
        much faster for cached properties, whereas the memory consumption
1134
        increases.
1135
    extra_properties : Iterable of callables
1136
        Add extra property computation functions that are not included with
1137
        skimage. The name of the property is derived from the function name,
1138
        the dtype is inferred by calling the function on a small sample.
1139
        If the name of an extra property clashes with the name of an existing
1140
        property the extra property will not be visible and a UserWarning is
1141
        issued. A property computation function must take a region mask as its
1142
        first argument. If the property requires an intensity image, it must
1143
        accept the intensity image as the second argument.
1144
    spacing: tuple of float, shape (ndim,)
1145
        The pixel spacing along each axis of the image.
1146
    offset : array-like of int, shape `(label_image.ndim,)`, optional
1147
        Coordinates of the origin ("top-left" corner) of the label image.
1148
        Normally this is ([0, ]0, 0), but it might be different if one wants
1149
        to obtain regionprops of subvolumes within a larger volume.
1150

1151
    Returns
1152
    -------
1153
    properties : list of RegionProperties
1154
        Each item describes one labeled region, and can be accessed using the
1155
        attributes listed below.
1156

1157
    Notes
1158
    -----
1159
    The following properties can be accessed as attributes or keys:
1160

1161
    **area** : float
1162
        Area of the region i.e. number of pixels of the region scaled by pixel-area.
1163
    **area_bbox** : float
1164
        Area of the bounding box i.e. number of pixels of bounding box scaled by pixel-area.
1165
    **area_convex** : float
1166
        Area of the convex hull image, which is the smallest convex
1167
        polygon that encloses the region.
1168
    **area_filled** : float
1169
        Area of the region with all the holes filled in.
1170
    **axis_major_length** : float
1171
        The length of the major axis of the ellipse that has the same
1172
        normalized second central moments as the region.
1173
    **axis_minor_length** : float
1174
        The length of the minor axis of the ellipse that has the same
1175
        normalized second central moments as the region.
1176
    **bbox** : tuple
1177
        Bounding box ``(min_row, min_col, max_row, max_col)``.
1178
        Pixels belonging to the bounding box are in the half-open interval
1179
        ``[min_row; max_row)`` and ``[min_col; max_col)``.
1180
    **centroid** : array
1181
        Centroid coordinate tuple ``(row, col)``.
1182
    **centroid_local** : array
1183
        Centroid coordinate tuple ``(row, col)``, relative to region bounding
1184
        box.
1185
    **centroid_weighted** : array
1186
        Centroid coordinate tuple ``(row, col)`` weighted with intensity
1187
        image.
1188
    **centroid_weighted_local** : array
1189
        Centroid coordinate tuple ``(row, col)``, relative to region bounding
1190
        box, weighted with intensity image.
1191
    **coords_scaled** : (K, 2) ndarray
1192
        Coordinate list ``(row, col)`` of the region scaled by ``spacing``.
1193
    **coords** : (K, 2) ndarray
1194
        Coordinate list ``(row, col)`` of the region.
1195
    **eccentricity** : float
1196
        Eccentricity of the ellipse that has the same second-moments as the
1197
        region. The eccentricity is the ratio of the focal distance
1198
        (distance between focal points) over the major axis length.
1199
        The value is in the interval [0, 1).
1200
        When it is 0, the ellipse becomes a circle.
1201
    **equivalent_diameter_area** : float
1202
        The diameter of a circle with the same area as the region.
1203
    **euler_number** : int
1204
        Euler characteristic of the set of non-zero pixels.
1205
        Computed as number of connected components subtracted by number of
1206
        holes (input.ndim connectivity). In 3D, number of connected
1207
        components plus number of holes subtracted by number of tunnels.
1208
    **extent** : float
1209
        Ratio of pixels in the region to pixels in the total bounding box.
1210
        Computed as ``area / (rows * cols)``
1211
    **feret_diameter_max** : float
1212
        Maximum Feret's diameter computed as the longest distance between
1213
        points around a region's convex hull contour as determined by
1214
        ``find_contours``. [5]_
1215
    **image** : (H, J) ndarray
1216
        Sliced binary region image which has the same size as bounding box.
1217
    **image_convex** : (H, J) ndarray
1218
        Binary convex hull image which has the same size as bounding box.
1219
    **image_filled** : (H, J) ndarray
1220
        Binary region image with filled holes which has the same size as
1221
        bounding box.
1222
    **image_intensity** : ndarray
1223
        Image inside region bounding box.
1224
    **inertia_tensor** : ndarray
1225
        Inertia tensor of the region for the rotation around its mass.
1226
    **inertia_tensor_eigvals** : tuple
1227
        The eigenvalues of the inertia tensor in decreasing order.
1228
    **intensity_max** : float
1229
        Value with the greatest intensity in the region.
1230
    **intensity_mean** : float
1231
        Value with the mean intensity in the region.
1232
    **intensity_min** : float
1233
        Value with the least intensity in the region.
1234
    **intensity_std** : float
1235
        Standard deviation of the intensity in the region.
1236
    **label** : int
1237
        The label in the labeled input image.
1238
    **moments** : (3, 3) ndarray
1239
        Spatial moments up to 3rd order::
1240

1241
            m_ij = sum{ array(row, col) * row^i * col^j }
1242

1243
        where the sum is over the `row`, `col` coordinates of the region.
1244
    **moments_central** : (3, 3) ndarray
1245
        Central moments (translation invariant) up to 3rd order::
1246

1247
            mu_ij = sum{ array(row, col) * (row - row_c)^i * (col - col_c)^j }
1248

1249
        where the sum is over the `row`, `col` coordinates of the region,
1250
        and `row_c` and `col_c` are the coordinates of the region's centroid.
1251
    **moments_hu** : tuple
1252
        Hu moments (translation, scale and rotation invariant).
1253
    **moments_normalized** : (3, 3) ndarray
1254
        Normalized moments (translation and scale invariant) up to 3rd order::
1255

1256
            nu_ij = mu_ij / m_00^[(i+j)/2 + 1]
1257

1258
        where `m_00` is the zeroth spatial moment.
1259
    **moments_weighted** : (3, 3) ndarray
1260
        Spatial moments of intensity image up to 3rd order::
1261

1262
            wm_ij = sum{ array(row, col) * row^i * col^j }
1263

1264
        where the sum is over the `row`, `col` coordinates of the region.
1265
    **moments_weighted_central** : (3, 3) ndarray
1266
        Central moments (translation invariant) of intensity image up to
1267
        3rd order::
1268

1269
            wmu_ij = sum{ array(row, col) * (row - row_c)^i * (col - col_c)^j }
1270

1271
        where the sum is over the `row`, `col` coordinates of the region,
1272
        and `row_c` and `col_c` are the coordinates of the region's weighted
1273
        centroid.
1274
    **moments_weighted_hu** : tuple
1275
        Hu moments (translation, scale and rotation invariant) of intensity
1276
        image.
1277
    **moments_weighted_normalized** : (3, 3) ndarray
1278
        Normalized moments (translation and scale invariant) of intensity
1279
        image up to 3rd order::
1280

1281
            wnu_ij = wmu_ij / wm_00^[(i+j)/2 + 1]
1282

1283
        where ``wm_00`` is the zeroth spatial moment (intensity-weighted area).
1284
    **num_pixels** : int
1285
        Number of foreground pixels.
1286
    **orientation** : float
1287
        Angle between the 0th axis (rows) and the major
1288
        axis of the ellipse that has the same second moments as the region,
1289
        ranging from `-pi/2` to `pi/2` counter-clockwise.
1290
    **perimeter** : float
1291
        Perimeter of object which approximates the contour as a line
1292
        through the centers of border pixels using a 4-connectivity.
1293
    **perimeter_crofton** : float
1294
        Perimeter of object approximated by the Crofton formula in 4
1295
        directions.
1296
    **slice** : tuple of slices
1297
        A slice to extract the object from the source image.
1298
    **solidity** : float
1299
        Ratio of pixels in the region to pixels of the convex hull image.
1300

1301
    Each region also supports iteration, so that you can do::
1302

1303
      for prop in region:
1304
          print(prop, region[prop])
1305

1306
    See Also
1307
    --------
1308
    label
1309

1310
    References
1311
    ----------
1312
    .. [1] Wilhelm Burger, Mark Burge. Principles of Digital Image Processing:
1313
           Core Algorithms. Springer-Verlag, London, 2009.
1314
    .. [2] B. Jähne. Digital Image Processing. Springer-Verlag,
1315
           Berlin-Heidelberg, 6. edition, 2005.
1316
    .. [3] T. H. Reiss. Recognizing Planar Objects Using Invariant Image
1317
           Features, from Lecture notes in computer science, p. 676. Springer,
1318
           Berlin, 1993.
1319
    .. [4] https://en.wikipedia.org/wiki/Image_moment
1320
    .. [5] W. Pabst, E. Gregorová. Characterization of particles and particle
1321
           systems, pp. 27-28. ICT Prague, 2007.
1322
           https://old.vscht.cz/sil/keramika/Characterization_of_particles/CPPS%20_English%20version_.pdf
1323

1324
    Examples
1325
    --------
1326
    >>> from skimage import data, util
1327
    >>> from skimage.measure import label, regionprops
1328
    >>> img = util.img_as_ubyte(data.coins()) > 110
1329
    >>> label_img = label(img, connectivity=img.ndim)
1330
    >>> props = regionprops(label_img)
1331
    >>> # centroid of first labeled object
1332
    >>> props[0].centroid
1333
    (22.72987986048314, 81.91228523446583)
1334
    >>> # centroid of first labeled object
1335
    >>> props[0]['centroid']
1336
    (22.72987986048314, 81.91228523446583)
1337

1338
    Add custom measurements by passing functions as ``extra_properties``
1339

1340
    >>> from skimage import data, util
1341
    >>> from skimage.measure import label, regionprops
1342
    >>> import numpy as np
1343
    >>> img = util.img_as_ubyte(data.coins()) > 110
1344
    >>> label_img = label(img, connectivity=img.ndim)
1345
    >>> def pixelcount(regionmask):
1346
    ...     return np.sum(regionmask)
1347
    >>> props = regionprops(label_img, extra_properties=(pixelcount,))
1348
    >>> props[0].pixelcount
1349
    7741
1350
    >>> props[1]['pixelcount']
1351
    42
1352

1353
    """
1354

1355
    if label_image.ndim not in (2, 3):
1356
        raise TypeError('Only 2-D and 3-D images supported.')
1357

1358
    if not np.issubdtype(label_image.dtype, np.integer):
1359
        if np.issubdtype(label_image.dtype, bool):
1360
            raise TypeError(
1361
                'Non-integer image types are ambiguous: '
1362
                'use skimage.measure.label to label the connected '
1363
                'components of label_image, '
1364
                'or label_image.astype(np.uint8) to interpret '
1365
                'the True values as a single label.'
1366
            )
1367
        else:
1368
            raise TypeError('Non-integer label_image types are ambiguous')
1369

1370
    if offset is None:
1371
        offset_arr = np.zeros((label_image.ndim,), dtype=int)
1372
    else:
1373
        offset_arr = np.asarray(offset)
1374
        if offset_arr.ndim != 1 or offset_arr.size != label_image.ndim:
1375
            raise ValueError(
1376
                'Offset should be an array-like of integers '
1377
                'of shape (label_image.ndim,); '
1378
                f'{offset} was provided.'
1379
            )
1380

1381
    regions = []
1382

1383
    objects = ndi.find_objects(label_image)
1384
    for i, sl in enumerate(objects):
1385
        if sl is None:
1386
            continue
1387

1388
        label = i + 1
1389

1390
        props = RegionProperties(
1391
            sl,
1392
            label,
1393
            label_image,
1394
            intensity_image,
1395
            cache,
1396
            spacing=spacing,
1397
            extra_properties=extra_properties,
1398
            offset=offset_arr,
1399
        )
1400
        regions.append(props)
1401

1402
    return regions
1403

1404

1405
def _parse_docs():
1406
    import re
1407
    import textwrap
1408

1409
    doc = regionprops.__doc__ or ''
1410
    arg_regex = r'\*\*(\w+)\*\* \:.*?\n(.*?)(?=\n    [\*\S]+)'
1411
    if sys.version_info >= (3, 13):
1412
        arg_regex = r'\*\*(\w+)\*\* \:.*?\n(.*?)(?=\n[\*\S]+)'
1413

1414
    matches = re.finditer(arg_regex, doc, flags=re.DOTALL)
1415
    prop_doc = {m.group(1): textwrap.dedent(m.group(2)) for m in matches}
1416

1417
    return prop_doc
1418

1419

1420
def _install_properties_docs():
1421
    prop_doc = _parse_docs()
1422

1423
    for p in [member for member in dir(RegionProperties) if not member.startswith('_')]:
1424
        getattr(RegionProperties, p).__doc__ = prop_doc[p]
1425

1426

1427
if __debug__:
1428
    # don't install docstrings when in optimized/non-debug mode
1429
    _install_properties_docs()
1430