Coverage for /usr/local/lib/python3.8/site-packages/spam/measurements/globalDescriptors.py: 100%
156 statements
« prev ^ index » next coverage.py v7.2.3, created at 2023-11-22 13:26 +0000
1# Library of SPAM functions for measuring global descriptors
2# Copyright (C) 2020 SPAM Contributors
3#
4# This program is free software: you can redistribute it and/or modify it
5# under the terms of the GNU General Public License as published by the Free
6# Software Foundation, either version 3 of the License, or (at your option)
7# any later version.
8#
9# This program is distributed in the hope that it will be useful, but WITHOUT
10# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
12# more details.
13#
14# You should have received a copy of the GNU General Public License along with
15# this program. If not, see <http://www.gnu.org/licenses/>.
18import numpy
19import spam.helpers
22def volume(segmented, phase=1, aspectRatio=(1.0, 1.0, 1.0), specific=False, ROI=None):
23 """
24 Computes the (specific) volume of a given phase for a certain *Region Of
25 Interest*.
27 Parameters
28 -----------
29 segmented: array of integers
30 Array where integers represent the phases.
31 The value 0 is stricly reserved for **not** Region Of Interest.
33 phase: integer, default=1
34 The phase for which the fraction volume is computed.
35 The phase should be greater than 0 if ROI is used.
37 aspectRatio: tuple of floats, default=(1.0, 1.0, 1.0)
38 Length between two nodes in every direction e.g. size of a cell.
39 The length of the tuple is the same as the spatial dimension.
41 specific: boolean, default=False
42 Returns the specific volume if True.
44 ROI: array of boolean, default=None
45 If not None, boolean array which represents the Region Of Interest.
46 Segmented will be 0 where ROI is False.
48 Returns
49 --------
50 float
51 The (specific) volume of the phase.
52 """
54 if phase == 0:
55 print(
56 "spam.measurements.globalDescriptors.volume: phase={}. Should be > 0.".format(
57 phase
58 )
59 )
61 if ROI is not None:
62 segmented = segmented * ROI
64 voxelSize = 1.0
65 for a in aspectRatio:
66 voxelSize *= a
68 if specific:
69 return float((segmented == phase).sum()) / float(segmented.size)
70 else:
71 return voxelSize * float((segmented == phase).sum())
74def eulerCharacteristic(segmented, phase=1, ROI=None, verbose=False):
75 """
76 Compute the Euler Characteristic of a given phase by counting
77 n-dimensionnal topological features using the formula:
79 .. code-block:: text
81 1D: X = #{ Vertices } - #{ Edges }
82 2D: X = #{ Vertices } - #{ Edges } + #{ Faces }
83 3D: X = #{ Vertices } - #{ Edges } + #{ Faces } - #{ Polyedra }
85 This function as been implemented using a optimized algorithm using rolled
86 views with unsymetrical structuring elements to count, thus avoiding
87 n nested loops.
89 Parameters
90 -----------
91 segmented: array of integers
92 Array where integers represent the phases.
93 The value 0 is stricly reserved for **not** Region Of Interest.
95 phase: integer, default=1
96 The phase for which the fraction volume is computed.
97 The phase should be greater than 0 if ROI is used.
99 ROI: array of boolean, default=None
100 If not None, boolean array which represents the Region Of Interest.
101 Segmented will be 0 where ROI is False.
103 verbose: boolean, default=False,
104 Verbose mode.
106 Returns
107 --------
108 int:
109 The Euler Characteristic
110 """
112 if phase == 0:
113 print(
114 "spam.measurements.globalDescriptors.eulerCharacteristic: phase={}. Should be > 0.".format(
115 phase
116 )
117 )
119 # get dimension
120 dim = len(segmented.shape)
122 if ROI is not None:
123 segmented = segmented * ROI
125 if dim == 1: # DIMENSION 1
127 # Count the vertices
128 counter = (segmented == phase).astype("<u1")
129 V = counter.sum()
131 # Count the edges
132 s_x = numpy.roll(counter, 1)
133 s_x[0] = 0
134 counter_x = numpy.logical_and(counter, s_x)
135 E = counter_x.sum()
137 # Euler characteristic
138 X = int(V - E)
140 if verbose:
141 print("Compute Euler Characteristic in 1D")
142 print("* V = #(vertices) = (+) {}".format(V))
143 print("* E = #(edges) = (-) {}".format(E))
144 print("* X = V - E = {}".format(X))
146 elif dim == 2: # DIMENSION 2
148 # Count the vertices
149 counter = (segmented == phase).astype("<u1")
150 V = counter.sum()
152 # Count the edges
153 s_x = numpy.roll(counter, 1, axis=0)
154 s_x[0, :] = 0
155 counter_x = numpy.logical_and(counter, s_x)
156 s_y = numpy.roll(counter, 1, axis=1)
157 s_y[:, 0] = 0
158 counter_y = numpy.logical_and(counter, s_y)
159 E = counter_x.sum() + counter_y.sum()
161 # Count the faces
162 s_xy = numpy.roll(s_x, 1, axis=1)
163 s_xy[:, 0] = 0
164 f_xy = numpy.logical_and(numpy.logical_and(counter_x, s_xy), counter_y)
165 F = f_xy.sum()
167 # Euler characteristic
168 X = int(V - E + F)
170 if verbose:
171 print("Compute Euler Characteristic in 2D")
172 print("* V = #(vertices) = (+) {}".format(V))
173 print("* E = #(edges) = (-) {}".format(E))
174 print("* F = #(faces) = (+) {}".format(F))
175 print("* X = V - E + F = {}".format(X))
177 elif dim == 3: # DIMENSION 3
179 # Count the vertices
180 counter = (segmented == phase).astype("<u1")
181 V = counter.sum()
183 # Count the edges
184 s_x = spam.helpers.singleShift(counter, 1, 0, sub=0)
185 counter_x = numpy.logical_and(counter, s_x)
186 s_y = spam.helpers.singleShift(counter, 1, 1, sub=0)
187 counter_y = numpy.logical_and(counter, s_y)
188 s_z = spam.helpers.singleShift(counter, 1, 2, sub=0)
189 counter_z = numpy.logical_and(counter, s_z)
190 E = counter_x.sum() + counter_y.sum() + counter_z.sum()
192 # Count the faces
193 s_xy = spam.helpers.singleShift(s_x, 1, 1, sub=0)
194 f_xy = numpy.logical_and(numpy.logical_and(counter_x, s_xy), counter_y)
195 s_xz = spam.helpers.singleShift(s_x, 1, 2, sub=0)
196 f_xz = numpy.logical_and(numpy.logical_and(counter_x, s_xz), counter_z)
197 s_yz = spam.helpers.singleShift(s_y, 1, 2, sub=0)
198 f_yz = numpy.logical_and(numpy.logical_and(counter_y, s_yz), counter_z)
199 F = f_xy.sum() + f_xz.sum() + f_yz.sum()
201 # Count the polyhedra
202 s_xyz = spam.helpers.singleShift(s_xy, 1, 2, sub=0)
203 f_xyz = numpy.logical_and(numpy.logical_and(f_xy, f_xz), f_yz)
204 p_xyz = numpy.logical_and(f_xyz, s_xyz)
205 P = p_xyz.sum()
207 # Euler characteristic
208 X = int(V - E + F - P)
210 if verbose:
211 print("* V = #(vertices) = (+) {}".format(V))
212 print("* E = #(edges) = (-) {}".format(E))
213 print("* F = #(faces) = (+) {}".format(F))
214 print("* P = #(polyedra) = (-) {}".format(P))
215 print("* X = V - E + F - P = {}".format(X))
217 return X
220def surfaceArea(field, level=0.0, aspectRatio=(1.0, 1.0, 1.0)):
221 """
222 Computes the surface area of a continuous field over a given level set.
224 This function is based on a marching cubes algorithm of skimage.
226 Parameters
227 -----------
228 field: array of floats
229 Array where some level sets represent the interface between phases.
231 level: float, default=0.0
232 The level set.
233 See ``skimage.measure.marching_cubes``.
234 Contour value to search for isosurfaces in volume.
235 If not given or None, the average of the min and max of vol is used.
237 aspectRatio: length 3 tuple of floats, default=(1.0, 1.0, 1.0)
238 Length between two nodes in every direction e.g. size of a cell.
239 It corresponds to ``spacing`` in ``skimage.measure.marching_cubes``.
241 Returns
242 --------
243 float:
244 The surface area.
245 """
246 import skimage.measure
248 # add a margin of 3 voxels
249 margedField = numpy.zeros([n + 6 for n in field.shape])
250 margedField[:, :, :] = field.min()
251 margedField[3:-3, 3:-3, 3:-3] = field
253 # http://scikit-image.org/docs/dev/api/skimage.measure.html?highlight=mesh_surface_area#skimage.measure.marching_cubes
254 verts, faces, _, _ = skimage.measure.marching_cubes(
255 margedField, level=level, spacing=aspectRatio
256 )
258 return skimage.measure.mesh_surface_area(verts, faces)
261# def totalCurvature(field, level=0.0, aspectRatio=(1.0, 1.0, 1.0)):
262# """
263# Computes the total curvature of a continuous field over a given level set.
264#
265# This function is based on a marching cubes algorithm of skimage.
266#
267# Parameters
268# -----------
269# field: array of floats
270# Array where some level sets represent the interface between phases.
271# level: float, default=0.0
272# The level set.
273# See ``skimage.measure.marching_cubes``.
274# Contour value to search for isosurfaces in volume.
275# If not given or None, the average of the min and max of vol is used.
276# aspectRatio: length 3 tuple of floats, default=(1.0, 1.0, 1.0)
277# Length between two nodes in every direction e.g. size of a cell.
278# It corresponds to ``spacing`` in ``skimage.measure.marching_cubes``.
279#
280# Returns
281# --------
282# float:
283# The total curvature.
284# """
285#
286# # specific function
287# def faceNormal(tri3d):
288# f_normal=numpy.zeros((tri3d.shape[0],3))
289# p1=tri3d[:,0,:]-tri3d[:,1,:]
290# p2=tri3d[:,0,:]-tri3d[:,2,:]
291# n=numpy.cross(p1, p2)
292# norm=LA.norm(n,axis=1)
293# f_normal[:,0]=n[:,0]/norm[:]
294# f_normal[:,1]=n[:,1]/norm[:]
295# f_normal[:,2]=n[:,2]/norm[:]
296#
297# return f_normal
298#
299# # import
300# import skimage.measure
301# from numpy import linalg as LA
302# import math
303#
304# # marching cubes
305#
306# # add a margin of 3 voxels
307# margedField = numpy.zeros([n+6 for n in field.shape])
308# margedField[:, :, :] = field.min()
309# margedField[3:-3, 3:-3, 3:-3] = field
310#
311# # http://scikit-image.org/docs/dev/api/skimage.measure.html?highlight=mesh_surface_area#skimage.measure.marching_cubes
312# verts, faces, _, _ = skimage.measure.marching_cubes(margedField, level=level, spacing=aspectRatio)
313#
314# # compute curvature
315#
316# tri3d=verts[faces]
317# nb_tri = tri3d.shape[0]
318# nb_p = verts.shape[0]
319#
320# # freeBoundary ### not yet
321# f_normal = faceNormal(tri3d)
322#
323# ### vectors, areas, angles, edges for every triangle
324# l_edg=numpy.zeros((nb_tri,3))
325# ang_tri=numpy.zeros((nb_tri,3))
326# f_center=numpy.zeros((nb_tri,3))
327# p=numpy.zeros((nb_tri))
328# v0=numpy.zeros((nb_tri,3))
329# v1=numpy.zeros((nb_tri,3))
330# v2=numpy.zeros((nb_tri,3))
331# area_tri=numpy.zeros((nb_tri))
332# for i in range(nb_tri):
333# # print("Boucle 1: {}/{}".format(i, nb_tri))
334# # points
335# p0=faces[i,0]
336# p1=faces[i,1]
337# p2=faces[i,2]
338#
339# # edges (vectors)
340# v0[i,:]=tri3d[i,1,:]-tri3d[i,0,:]
341# v1[i,:]=tri3d[i,2,:]-tri3d[i,1,:]
342# v2[i,:]=tri3d[i,0,:]-tri3d[i,2,:]
343#
344# # length of edges
345# l_edg[i,0]= LA.norm(v0[i,:])
346# l_edg[i,1]= LA.norm(v1[i,:])
347# l_edg[i,2]= LA.norm(v2[i,:])
348#
349# # angle
350# ang_tri[i,0]=math.acos(numpy.dot(v0[i,:]/l_edg[i,0], -v2[i,:]/l_edg[i,2]))
351# ang_tri[i,1]=math.acos(numpy.dot(-v0[i,:]/l_edg[i,0], v1[i,:]/l_edg[i,1]))
352# ang_tri[i,2]=math.pi-ang_tri[i,0]-ang_tri[i,1]
353#
354# # incenter
355# p[i]=numpy.sum(l_edg[i,:])
356# f_center[i,0]=(l_edg[i,0]*tri3d[i,2,0]+l_edg[i,1]*tri3d[i,0,0]+l_edg[i,2]*tri3d[i,1,0])/p[i]
357# f_center[i,1]=(l_edg[i,0]*tri3d[i,2,1]+l_edg[i,1]*tri3d[i,0,1]+l_edg[i,2]*tri3d[i,1,1])/p[i]
358# f_center[i,2]=(l_edg[i,0]*tri3d[i,2,2]+l_edg[i,1]*tri3d[i,0,2]+l_edg[i,2]*tri3d[i,1,2])/p[i]
359#
360# # surface
361# area_tri[i]=((numpy.cross(v0[i,:],v1[i,:])**2).sum()**0.5).sum() /2.0
362#
363# a_mixed=numpy.zeros((nb_p))
364# alf=numpy.zeros((nb_p))
365# MC =numpy.zeros((nb_p))
366# GC =numpy.zeros((nb_p))
367#
368# for i in range(nb_p):
369# # print("Boucle 2: {}/{}".format(i, nb_p))
370# mc_vec=[0,0,0]
371# n_vec=[0,0,0]
372#
373# # if find(bndry_edge) ### not yet
374#
375# # vertex attachment ? find?
376# find=numpy.where(faces==i)
377# neib_tri=find[0]
378#
379# for j in range(len(neib_tri)):
380# neib=neib_tri[j]
381# # sum of angles around point i ==> GC
382# for k in range(3):
383# if faces[neib,k]==i:
384# alf[i]=alf[i]+ang_tri[neib,k]
385# break
386# # mean curvature operator
387# if k==0:
388# mc_vec=mc_vec+(v0[neib,:]/math.tan(ang_tri[neib,2])-v2[neib,:]/math.tan(ang_tri[neib,1]))
389# elif k==1:
390# mc_vec=mc_vec+(v1[neib,:]/math.tan(ang_tri[neib,0])-v0[neib,:]/math.tan(ang_tri[neib,2]))
391# elif k==2:
392# mc_vec=mc_vec+(v2[neib,:]/math.tan(ang_tri[neib,1])-v1[neib,:]/math.tan(ang_tri[neib,0]))
393#
394# # A_mixed calculation
395# if (ang_tri[neib,k]>=math.pi/2.0):
396# a_mixed[i]=a_mixed[i]+area_tri[neib]/2.0
397# else:
398# if any(ang >=math.pi/2.0 for ang in ang_tri[neib,:]):
399# a_mixed[i]=a_mixed[i]+area_tri[neib]/4.0
400# else:
401# sum=0
402# for m in range(3):
403# if m!=k:
404# ll=m+1
405# if ll==3:
406# ll=0
407# sum=sum+(l_edg[neib,ll]**2/math.tan(ang_tri[neib,m]))
408# a_mixed[i]=a_mixed[i]+sum/8.0
409#
410#
411# # normal vector at each vertex
412# # weighted average of normal vectors of neighbour triangles
413# wi=1.0/LA.norm( [f_center[neib,0]-verts[i,0],f_center[neib,1]-verts[i,1],f_center[neib,2]-verts[i,2]] )
414# n_vec=n_vec+wi*f_normal[neib,:]
415#
416# GC[i]=(2.0*math.pi-alf[i])/a_mixed[i]
417#
418# mc_vec=0.25*mc_vec/a_mixed[i]
419# n_vec=n_vec/LA.norm(n_vec)
420#
421# # sign of MC
422# if numpy.dot(mc_vec, n_vec)<0:
423# MC[i]=-LA.norm(mc_vec)
424# else:
425# MC[i]=LA.norm(mc_vec)
426#
427#
428# totalCurvature = numpy.multiply(MC, a_mixed).sum()
429#
430# return totalCurvature
433def totalCurvature(
434 field,
435 level=0.0,
436 aspectRatio=(1.0, 1.0, 1.0),
437 stepSize=None,
438 getMeshValues=False,
439 fileName=None,
440):
441 """
442 Computes the total curvature of a continuous field over a given level set.
444 This function is based on a marching cubes algorithm of skimage.
446 Parameters
447 -----------
448 field: array of floats
449 Array where some level sets represent the interface between phases
451 level: float
452 The level set
453 See ``skimage.measure.marching_cubes``.
454 Contour value to search for isosurfaces in volume.
455 If None is given, the average of the min and max of vol is used
456 Default = 0.0
458 aspectRatio: length 3 tuple of floats
459 Length between two nodes in every direction `e.g.`, size of a cell.
460 It corresponds to ``spacing`` in ``skimage.measure.marching_cubes``
461 Default = (1.0, 1.0, 1.0)
463 step_size: int
464 Step size in voxels
465 Larger steps yield faster but coarser results. The result will always be topologically correct though.
466 Default = 1
468 getMeshValues: boolean
469 Return the mean and gauss curvatures and the areas of the mesh vertices
470 Default = False
472 fileName: string
473 Name of the output vtk file
474 Only if `getMeshValues` = True
475 Default = None
477 Returns
478 --------
479 totalCurvature: float
480 If ``getMeshValues`` = True,
481 total curvature and a (nx3) array of n vertices
482 where 1st column is the mean curvature, 2nd column is the gaussian curvature and 3rd column is the area
483 """
484 # import
485 import skimage.measure
486 import spam.measurements.measurementsToolkit as mtk
488 # marching cubes
489 # add a margin of 3 voxels
490 margedField = numpy.zeros([n + 6 for n in field.shape])
491 margedField[:, :, :] = field.min()
492 margedField[3:-3, 3:-3, 3:-3] = field
494 # http://scikit-image.org/docs/dev/api/skimage.measure.html?highlight=mesh_surface_area#skimage.measure.marching_cubes
495 if stepSize is not None:
496 verts, faces, _, _ = skimage.measure.marching_cubes(
497 margedField, level=level, spacing=aspectRatio, step_size=stepSize
498 )
499 else:
500 verts, faces, _, _ = skimage.measure.marching_cubes(
501 margedField, level=level, spacing=aspectRatio
502 )
504 # call CPP function return gaussian curvature, mean curvature and areas
505 meanGaussArea = mtk.computeCurvatures(
506 faces.astype("<u8").tolist(), verts.astype("<f8").tolist()
507 )
508 meanGaussArea = numpy.array(meanGaussArea)
510 totalCurvature = numpy.multiply(meanGaussArea[:, 0], meanGaussArea[:, 2]).sum()
512 if fileName is not None:
513 pointData = {
514 "meanCurvature": meanGaussArea[:, 0],
515 "gaussCurvature": meanGaussArea[:, 1],
516 "area": meanGaussArea[:, 2],
517 }
519 import spam.helpers
521 spam.helpers.writeUnstructuredVTK(
522 verts, faces, pointData=pointData, elementType="triangle", fileName=fileName
523 )
525 if getMeshValues:
526 return (
527 totalCurvature,
528 meanGaussArea[:, 0],
529 meanGaussArea[:, 1],
530 meanGaussArea[:, 2],
531 )
532 else:
533 return totalCurvature
536def perimeter(segmented, phase=1, aspectRatio=(1.0, 1.0), ROI=None):
537 """
538 Computes the perimeter of a given phase.
540 Parameters
541 -----------
542 segmented: array of integers
543 Array where integers represent the phases.
544 The value 0 is stricly reserved for **not** Region Of Interest.
546 phase: integer, default=1
547 The phase for which the fraction volume is computed.
548 The phase should be greater than 0 if ROI is used.
550 aspectRatio: length 2 tuple of floats, default=(1.0, 1.0)
551 Length between two nodes in every direction e.g. size of a cell.
553 ROI: array of boolean, default=None
554 If not None, boolean array which represents the Region Of Interest.
555 Segmented will be 0 where ROI is False.
557 Returns
558 --------
559 float:
560 The perimeter.
562 """
563 import spam.filters
565 if phase == 0:
566 print(
567 "spam.measurements.globalDescriptors.perimeter: phase={}. Should be > 0.".format(
568 phase
569 )
570 )
572 if ROI is not None:
573 segmented = segmented * ROI
575 # add border
576 im_l = numpy.zeros((segmented.shape[0] + 2, segmented.shape[1] + 2), dtype=bool)
577 im_l[1 : segmented.shape[0] + 1, 1 : segmented.shape[1] + 1] = segmented == phase
579 # take perimeters
580 peri = spam.filters.binaryMorphologicalGradient(im_l)
582 # compute number of each voxel types
583 late = peri & spam.helpers.multipleShifts(peri, shifts=[+1, 0])
584 late = late + (peri & spam.helpers.multipleShifts(peri, shifts=[0, +1]))
585 late = late + (peri & spam.helpers.multipleShifts(peri, shifts=[-1, 0]))
586 late = late + (peri & spam.helpers.multipleShifts(peri, shifts=[0, -1]))
587 diag = peri & spam.helpers.multipleShifts(peri, shifts=[+1, +1])
588 diag = diag + (peri & spam.helpers.multipleShifts(peri, shifts=[+1, -1]))
589 diag = diag + (peri & spam.helpers.multipleShifts(peri, shifts=[-1, +1]))
590 diag = diag + (peri & spam.helpers.multipleShifts(peri, shifts=[-1, -1]))
591 inte = diag & late
592 jdiag = diag & numpy.logical_not(inte)
593 jlate = late & numpy.logical_not(inte)
595 voxelSize = aspectRatio[0]
597 # compute global coefficient
598 alpha = (
599 (1.0 + 0.5**2) ** 0.5 * numpy.sum(inte)
600 + numpy.sum(jlate)
601 + 2.0**0.5 * numpy.sum(jdiag)
602 )
604 return alpha * voxelSize
607def generic(
608 segmented,
609 n,
610 phase=1,
611 level=0.0,
612 aspectRatio=(1.0, 1.0, 1.0),
613 specific=False,
614 ROI=None,
615 verbose=False,
616):
617 """
618 Maps all measures in a homogeneous format
620 Parameters
621 -----------
622 segmented: array of integers
623 Array where integers represent the phases.
624 The value 0 is stricly reserved for **not** Region Of Interest.
626 n: int
627 The number of the measure. The meaning depends of the dimension.
629 .. code-block:: text
631 * In 1D:
632 * n=0: Euler characteristic
633 * n=1: (specific) length
634 * In 2D:
635 * n=0: Euler characteristic
636 * n=1: perimeter
637 * n=2: (specific) surface
638 * In 3D:
639 * n=0: Euler Characteristic
640 * n=1: ?
641 * n=2: surface
642 * n=3: (specific) volume
644 phase: integer, default=1
645 The phase for which the fraction volume is computed.
646 The phase should be greater than 0 if ROI is used.
648 level: float, default=0.0
649 The level set.
650 See ``skimage.measure.marching_cubes``.
651 Contour value to search for isosurfaces in volume.
652 If not given or None, the average of the min and max of vol is used.
654 specific: boolean, default=False
655 Returns the specific volume if True.
657 aspectRatio: length 3 tuple of floats, default=(1.0, 1.0, 1.0)
658 Length between two nodes in every direction e.g. size of a cell.
660 ROI: array of boolean, default=None
661 If not None, boolean array which represents the Region Of Interest.
662 Segmented will be 0 where ROI is False.
664 verbose: boolean, default=False,
665 Verbose mode.
666 Returns
667 --------
668 float:
669 The measure
670 """
672 # get dimension
673 dim = len(segmented.shape)
675 if dim == 1:
676 if n == 0:
677 return eulerCharacteristic(segmented, phase=phase, ROI=ROI, verbose=verbose)
678 elif n == 1:
679 return volume(
680 segmented,
681 phase=phase,
682 aspectRatio=aspectRatio,
683 specific=specific,
684 ROI=ROI,
685 )
686 elif dim == 2:
687 if n == 0:
688 return eulerCharacteristic(segmented, phase=phase, ROI=ROI, verbose=verbose)
689 if n == 1:
690 return perimeter(segmented, phase=phase, ROI=ROI)
691 if n == 2:
692 return volume(
693 segmented,
694 phase=phase,
695 aspectRatio=aspectRatio,
696 specific=specific,
697 ROI=ROI,
698 )
699 elif dim == 3:
700 if n == 0:
701 return eulerCharacteristic(segmented, phase=phase, ROI=ROI, verbose=verbose)
702 if n == 1:
703 return -1
704 if n == 2:
705 return surfaceArea(segmented, level=level, aspectRatio=aspectRatio)
706 if n == 3:
707 return volume(
708 segmented,
709 phase=phase,
710 aspectRatio=aspectRatio,
711 specific=specific,
712 ROI=ROI,
713 )
715 print(
716 "spam.measurements.globalDescriptors.generic: dim={} (should be between 1 and 3) and n={} (should be between 0 and dim). Return NaN".format(
717 dim, n
718 )
719 )
720 return numpy.nan