vedo.utils
Utilities submodule.
1#!/usr/bin/env python3 2# -*- coding: utf-8 -*- 3import os 4import time 5from typing import Union, Tuple, MutableSequence, List 6import numpy as np 7 8from vtkmodules.util.numpy_support import numpy_to_vtk, vtk_to_numpy 9from vtkmodules.util.numpy_support import numpy_to_vtkIdTypeArray 10import vedo.vtkclasses as vtki 11 12import vedo 13 14 15__docformat__ = "google" 16 17__doc__ = "Utilities submodule." 18 19__all__ = [ 20 "OperationNode", 21 "ProgressBar", 22 "progressbar", 23 "Minimizer", 24 "geometry", 25 "is_sequence", 26 "lin_interpolate", 27 "vector", 28 "mag", 29 "mag2", 30 "versor", 31 "precision", 32 "round_to_digit", 33 "point_in_triangle", 34 "point_line_distance", 35 "closest", 36 "grep", 37 "make_bands", 38 "pack_spheres", 39 "humansort", 40 "print_histogram", 41 "print_inheritance_tree", 42 "camera_from_quaternion", 43 "camera_from_neuroglancer", 44 "camera_from_dict", 45 "camera_to_dict", 46 "oriented_camera", 47 "vedo2trimesh", 48 "trimesh2vedo", 49 "vedo2meshlab", 50 "meshlab2vedo", 51 "vedo2open3d", 52 "open3d2vedo", 53 "vtk2numpy", 54 "numpy2vtk", 55 "get_uv", 56 "andrews_curves", 57] 58 59########################################################################### 60class OperationNode: 61 """ 62 Keep track of the operations which led to a final state. 63 """ 64 # https://www.graphviz.org/doc/info/shapes.html#html 65 # Mesh #e9c46a 66 # Follower #d9ed92 67 # Volume, UnstructuredGrid #4cc9f0 68 # TetMesh #9e2a2b 69 # File #8a817c 70 # Image #f28482 71 # Assembly #f08080 72 73 def __init__( 74 self, operation, parents=(), comment="", shape="none", c="#e9c46a", style="filled" 75 ) -> None: 76 """ 77 Keep track of the operations which led to a final object. 78 This allows to show the `pipeline` tree for any `vedo` object with e.g.: 79 80 ```python 81 from vedo import * 82 sp = Sphere() 83 sp.clean().subdivide() 84 sp.pipeline.show() 85 ``` 86 87 Arguments: 88 operation : (str, class) 89 descriptor label, if a class is passed then grab its name 90 parents : (list) 91 list of the parent classes the object comes from 92 comment : (str) 93 a second-line text description 94 shape : (str) 95 shape of the frame, check out [this link.](https://graphviz.org/doc/info/shapes.html) 96 c : (hex) 97 hex color 98 style : (str) 99 comma-separated list of styles 100 101 Example: 102 ```python 103 from vedo.utils import OperationNode 104 105 op_node1 = OperationNode("Operation1", c="lightblue") 106 op_node2 = OperationNode("Operation2") 107 op_node3 = OperationNode("Operation3", shape='diamond') 108 op_node4 = OperationNode("Operation4") 109 op_node5 = OperationNode("Operation5") 110 op_node6 = OperationNode("Result", c="lightgreen") 111 112 op_node3.add_parent(op_node1) 113 op_node4.add_parent(op_node1) 114 op_node3.add_parent(op_node2) 115 op_node5.add_parent(op_node2) 116 op_node6.add_parent(op_node3) 117 op_node6.add_parent(op_node5) 118 op_node6.add_parent(op_node1) 119 120 op_node6.show(orientation="TB") 121 ``` 122  123 """ 124 if not vedo.settings.enable_pipeline: 125 return 126 127 if isinstance(operation, str): 128 self.operation = operation 129 else: 130 self.operation = operation.__class__.__name__ 131 self.operation_plain = str(self.operation) 132 133 pp = [] # filter out invalid stuff 134 for p in parents: 135 if hasattr(p, "pipeline"): 136 pp.append(p.pipeline) 137 self.parents = pp 138 139 if comment: 140 self.operation = f"<{self.operation}<BR/><SUB><I>{comment}</I></SUB>>" 141 142 self.dot = None 143 self.time = time.time() 144 self.shape = shape 145 self.style = style 146 self.color = c 147 self.counts = 0 148 149 def add_parent(self, parent) -> None: 150 self.parents.append(parent) 151 152 def _build_tree(self, dot): 153 dot.node( 154 str(id(self)), 155 label=self.operation, 156 shape=self.shape, 157 color=self.color, 158 style=self.style, 159 ) 160 for parent in self.parents: 161 if parent: 162 t = f"{self.time - parent.time: .1f}s" 163 dot.edge(str(id(parent)), str(id(self)), label=t) 164 parent._build_tree(dot) 165 166 def __repr__(self): 167 try: 168 from treelib import Tree 169 except ImportError: 170 vedo.logger.error( 171 "To use this functionality please install treelib:" 172 "\n pip install treelib" 173 ) 174 return "" 175 176 def _build_tree(parent): 177 for par in parent.parents: 178 if par: 179 op = par.operation_plain 180 tree.create_node( 181 op, op + str(par.time), parent=parent.operation_plain + str(parent.time) 182 ) 183 _build_tree(par) 184 try: 185 tree = Tree() 186 tree.create_node(self.operation_plain, self.operation_plain + str(self.time)) 187 _build_tree(self) 188 out = tree.show(stdout=False) 189 except: 190 out = f"Sorry treelib failed to build the tree for '{self.operation_plain}()'." 191 return out 192 193 def print(self) -> None: 194 """Print the tree of operations.""" 195 print(self.__str__()) 196 197 def show(self, orientation="LR", popup=True) -> None: 198 """Show the graphviz output for the pipeline of this object""" 199 if not vedo.settings.enable_pipeline: 200 return 201 202 try: 203 from graphviz import Digraph 204 except ImportError: 205 vedo.logger.error("please install graphviz with command\n pip install graphviz") 206 return 207 208 # visualize the entire tree 209 dot = Digraph( 210 node_attr={"fontcolor": "#201010", "fontname": "Helvetica", "fontsize": "12"}, 211 edge_attr={"fontname": "Helvetica", "fontsize": "6", "arrowsize": "0.4"}, 212 ) 213 dot.attr(rankdir=orientation) 214 215 self.counts = 0 216 self._build_tree(dot) 217 self.dot = dot 218 219 home_dir = os.path.expanduser("~") 220 gpath = os.path.join( 221 home_dir, vedo.settings.cache_directory, "vedo", "pipeline_graphviz") 222 223 dot.render(gpath, view=popup) 224 225 226########################################################################### 227class ProgressBar: 228 """ 229 Class to print a progress bar. 230 """ 231 232 def __init__( 233 self, 234 start, 235 stop, 236 step=1, 237 c=None, 238 bold=True, 239 italic=False, 240 title="", 241 eta=True, 242 delay=-1, 243 width=25, 244 char="\U00002501", 245 char_back="\U00002500", 246 ) -> None: 247 """ 248 Class to print a progress bar with optional text message. 249 250 Check out also function `progressbar()`. 251 252 Arguments: 253 start : (int) 254 starting value 255 stop : (int) 256 stopping value 257 step : (int) 258 step value 259 c : (str) 260 color in hex format 261 title : (str) 262 title text 263 eta : (bool) 264 estimate time of arrival 265 delay : (float) 266 minimum time before printing anything, 267 if negative use the default value 268 as set in `vedo.settings.progressbar_delay` 269 width : (int) 270 width of the progress bar 271 char : (str) 272 character to use for the progress bar 273 char_back : (str) 274 character to use for the background of the progress bar 275 276 Example: 277 ```python 278 import time 279 from vedo import ProgressBar 280 pb = ProgressBar(0,40, c='r') 281 for i in pb.range(): 282 time.sleep(0.1) 283 pb.print() 284 ``` 285  286 """ 287 self.char = char 288 self.char_back = char_back 289 290 self.title = title + " " 291 if title: 292 self.title = " " + self.title 293 294 if delay < 0: 295 delay = vedo.settings.progressbar_delay 296 297 self.start = start 298 self.stop = stop 299 self.step = step 300 301 self.color = c 302 self.bold = bold 303 self.italic = italic 304 self.width = width 305 self.pbar = "" 306 self.percent = 0.0 307 self.percent_int = 0 308 self.eta = eta 309 self.delay = delay 310 311 self.t0 = time.time() 312 self._remaining = 1e10 313 314 self._update(0) 315 316 self._counts = 0 317 self._oldbar = "" 318 self._lentxt = 0 319 self._range = np.arange(start, stop, step) 320 321 def print(self, txt="", c=None) -> None: 322 """Print the progress bar with an optional message.""" 323 if not c: 324 c = self.color 325 326 self._update(self._counts + self.step) 327 328 if self.delay: 329 if time.time() - self.t0 < self.delay: 330 return 331 332 if self.pbar != self._oldbar: 333 self._oldbar = self.pbar 334 335 if self.eta and self._counts > 1: 336 337 tdenom = time.time() - self.t0 338 if tdenom: 339 vel = self._counts / tdenom 340 self._remaining = (self.stop - self._counts) / vel 341 else: 342 vel = 1 343 self._remaining = 0.0 344 345 if self._remaining > 60: 346 _mins = int(self._remaining / 60) 347 _secs = self._remaining - 60 * _mins 348 mins = f"{_mins}m" 349 secs = f"{int(_secs + 0.5)}s " 350 else: 351 mins = "" 352 secs = f"{int(self._remaining + 0.5)}s " 353 354 vel = round(vel, 1) 355 eta = f"eta: {mins}{secs}({vel} it/s) " 356 if self._remaining < 0.5: 357 dt = time.time() - self.t0 358 if dt > 60: 359 _mins = int(dt / 60) 360 _secs = dt - 60 * _mins 361 mins = f"{_mins}m" 362 secs = f"{int(_secs + 0.5)}s " 363 else: 364 mins = "" 365 secs = f"{int(dt + 0.5)}s " 366 eta = f"elapsed: {mins}{secs}({vel} it/s) " 367 txt = "" 368 else: 369 eta = "" 370 371 eraser = " " * self._lentxt + "\b" * self._lentxt 372 373 s = f"{self.pbar} {eraser}{eta}{txt}\r" 374 vedo.printc(s, c=c, bold=self.bold, italic=self.italic, end="") 375 if self.percent > 99.999: 376 print("") 377 378 self._lentxt = len(txt) 379 380 def range(self) -> np.ndarray: 381 """Return the range iterator.""" 382 return self._range 383 384 def _update(self, counts): 385 if counts < self.start: 386 counts = self.start 387 elif counts > self.stop: 388 counts = self.stop 389 self._counts = counts 390 391 self.percent = (self._counts - self.start) * 100.0 392 393 delta = self.stop - self.start 394 if delta: 395 self.percent /= delta 396 else: 397 self.percent = 0.0 398 399 self.percent_int = int(round(self.percent)) 400 af = self.width - 2 401 nh = int(round(self.percent_int / 100 * af)) 402 pbar_background = "\x1b[2m" + self.char_back * (af - nh) 403 self.pbar = f"{self.title}{self.char * (nh-1)}{pbar_background}" 404 if self.percent < 100.0: 405 ps = f" {self.percent_int}%" 406 else: 407 ps = "" 408 self.pbar += ps 409 410 411##################################### 412def progressbar( 413 iterable, 414 c=None, bold=True, italic=False, title="", 415 eta=True, width=25, delay=-1, 416 ): 417 """ 418 Function to print a progress bar with optional text message. 419 420 Use delay to set a minimum time before printing anything. 421 If delay is negative, then use the default value 422 as set in `vedo.settings.progressbar_delay`. 423 424 Arguments: 425 start : (int) 426 starting value 427 stop : (int) 428 stopping value 429 step : (int) 430 step value 431 c : (str) 432 color in hex format 433 title : (str) 434 title text 435 eta : (bool) 436 estimate time of arrival 437 delay : (float) 438 minimum time before printing anything, 439 if negative use the default value 440 set in `vedo.settings.progressbar_delay` 441 width : (int) 442 width of the progress bar 443 char : (str) 444 character to use for the progress bar 445 char_back : (str) 446 character to use for the background of the progress bar 447 448 Example: 449 ```python 450 import time 451 for i in progressbar(range(100), c='r'): 452 time.sleep(0.1) 453 ``` 454  455 """ 456 try: 457 if is_number(iterable): 458 total = int(iterable) 459 iterable = range(total) 460 else: 461 total = len(iterable) 462 except TypeError: 463 iterable = list(iterable) 464 total = len(iterable) 465 466 pb = ProgressBar( 467 0, total, c=c, bold=bold, italic=italic, title=title, 468 eta=eta, delay=delay, width=width, 469 ) 470 for item in iterable: 471 pb.print() 472 yield item 473 474 475########################################################### 476class Minimizer: 477 """ 478 A function minimizer that uses the Nelder-Mead method. 479 480 The algorithm constructs an n-dimensional simplex in parameter 481 space (i.e. a tetrahedron if the number or parameters is 3) 482 and moves the vertices around parameter space until 483 a local minimum is found. The amoeba method is robust, 484 reasonably efficient, but is not guaranteed to find 485 the global minimum if several local minima exist. 486 487 Arguments: 488 function : (callable) 489 the function to minimize 490 max_iterations : (int) 491 the maximum number of iterations 492 contraction_ratio : (float) 493 The contraction ratio. 494 The default value of 0.5 gives fast convergence, 495 but larger values such as 0.6 or 0.7 provide greater stability. 496 expansion_ratio : (float) 497 The expansion ratio. 498 The default value is 2.0, which provides rapid expansion. 499 Values between 1.1 and 2.0 are valid. 500 tol : (float) 501 the tolerance for convergence 502 503 Example: 504 - [nelder-mead.py](https://github.com/marcomusy/vedo/blob/master/examples/others/nelder-mead.py) 505 """ 506 def __init__( 507 self, 508 function=None, 509 max_iterations=10000, 510 contraction_ratio=0.5, 511 expansion_ratio=2.0, 512 tol=1e-5, 513 ) -> None: 514 self.function = function 515 self.tolerance = tol 516 self.contraction_ratio = contraction_ratio 517 self.expansion_ratio = expansion_ratio 518 self.max_iterations = max_iterations 519 self.minimizer = vtki.new("AmoebaMinimizer") 520 self.minimizer.SetFunction(self._vtkfunc) 521 self.results = {} 522 self.parameters_path = [] 523 self.function_path = [] 524 525 def _vtkfunc(self): 526 n = self.minimizer.GetNumberOfParameters() 527 ain = [self.minimizer.GetParameterValue(i) for i in range(n)] 528 r = self.function(ain) 529 self.minimizer.SetFunctionValue(r) 530 self.parameters_path.append(ain) 531 self.function_path.append(r) 532 return r 533 534 def eval(self, parameters=()) -> float: 535 """ 536 Evaluate the function at the current or given parameters. 537 """ 538 if len(parameters) == 0: 539 return self.minimizer.EvaluateFunction() 540 self.set_parameters(parameters) 541 return self.function(parameters) 542 543 def set_parameter(self, name, value, scale=1.0) -> None: 544 """ 545 Set the parameter value. 546 The initial amount by which the parameter 547 will be modified during the search for the minimum. 548 """ 549 self.minimizer.SetParameterValue(name, value) 550 self.minimizer.SetParameterScale(name, scale) 551 552 def set_parameters(self, parameters) -> None: 553 """ 554 Set the parameters names and values from a dictionary. 555 """ 556 for name, value in parameters.items(): 557 if len(value) == 2: 558 self.set_parameter(name, value[0], value[1]) 559 else: 560 self.set_parameter(name, value) 561 562 def minimize(self) -> dict: 563 """ 564 Minimize the input function. 565 566 Returns: 567 dict : 568 the minimization results 569 init_parameters : (dict) 570 the initial parameters 571 parameters : (dict) 572 the final parameters 573 min_value : (float) 574 the minimum value 575 iterations : (int) 576 the number of iterations 577 max_iterations : (int) 578 the maximum number of iterations 579 tolerance : (float) 580 the tolerance for convergence 581 convergence_flag : (int) 582 zero if the tolerance stopping 583 criterion has been met. 584 parameters_path : (np.array) 585 the path of the minimization 586 algorithm in parameter space 587 function_path : (np.array) 588 the path of the minimization 589 algorithm in function space 590 hessian : (np.array) 591 the Hessian matrix of the 592 function at the minimum 593 parameter_errors : (np.array) 594 the errors on the parameters 595 """ 596 n = self.minimizer.GetNumberOfParameters() 597 out = [( 598 self.minimizer.GetParameterName(i), 599 (self.minimizer.GetParameterValue(i), 600 self.minimizer.GetParameterScale(i)) 601 ) for i in range(n)] 602 self.results["init_parameters"] = dict(out) 603 604 self.minimizer.SetTolerance(self.tolerance) 605 self.minimizer.SetContractionRatio(self.contraction_ratio) 606 self.minimizer.SetExpansionRatio(self.expansion_ratio) 607 self.minimizer.SetMaxIterations(self.max_iterations) 608 609 self.minimizer.Minimize() 610 self.results["convergence_flag"] = not bool(self.minimizer.Iterate()) 611 612 out = [( 613 self.minimizer.GetParameterName(i), 614 self.minimizer.GetParameterValue(i), 615 ) for i in range(n)] 616 617 self.results["parameters"] = dict(out) 618 self.results["min_value"] = self.minimizer.GetFunctionValue() 619 self.results["iterations"] = self.minimizer.GetIterations() 620 self.results["max_iterations"] = self.minimizer.GetMaxIterations() 621 self.results["tolerance"] = self.minimizer.GetTolerance() 622 self.results["expansion_ratio"] = self.expansion_ratio 623 self.results["contraction_ratio"] = self.contraction_ratio 624 self.results["parameters_path"] = np.array(self.parameters_path) 625 self.results["function_path"] = np.array(self.function_path) 626 self.results["hessian"] = np.zeros((n,n)) 627 self.results["parameter_errors"] = np.zeros(n) 628 return self.results 629 630 def compute_hessian(self, epsilon=0) -> np.array: 631 """ 632 Compute the Hessian matrix of `function` at the 633 minimum numerically. 634 635 Arguments: 636 epsilon : (float) 637 Step size used for numerical approximation. 638 639 Returns: 640 array: Hessian matrix of `function` at minimum. 641 """ 642 if not epsilon: 643 epsilon = self.tolerance * 10 644 n = self.minimizer.GetNumberOfParameters() 645 x0 = [self.minimizer.GetParameterValue(i) for i in range(n)] 646 hessian = np.zeros((n, n)) 647 for i in vedo.progressbar(n, title="Computing Hessian", delay=2): 648 for j in range(n): 649 xijp = np.copy(x0) 650 xijp[i] += epsilon 651 xijp[j] += epsilon 652 xijm = np.copy(x0) 653 xijm[i] += epsilon 654 xijm[j] -= epsilon 655 xjip = np.copy(x0) 656 xjip[i] -= epsilon 657 xjip[j] += epsilon 658 xjim = np.copy(x0) 659 xjim[i] -= epsilon 660 xjim[j] -= epsilon 661 # Second derivative approximation 662 fijp = self.function(xijp) 663 fijm = self.function(xijm) 664 fjip = self.function(xjip) 665 fjim = self.function(xjim) 666 hessian[i, j] = (fijp - fijm - fjip + fjim) / (2 * epsilon**2) 667 self.results["hessian"] = hessian 668 try: 669 ihess = np.linalg.inv(hessian) 670 self.results["parameter_errors"] = np.sqrt(np.diag(ihess)) 671 except: 672 vedo.logger.warning("Cannot compute hessian for parameter errors") 673 self.results["parameter_errors"] = np.zeros(n) 674 return hessian 675 676 def __str__(self) -> str: 677 out = vedo.printc( 678 f"vedo.utils.Minimizer at ({hex(id(self))})".ljust(75), 679 bold=True, invert=True, return_string=True, 680 ) 681 out += "Function name".ljust(20) + self.function.__name__ + "()\n" 682 out += "-------- parameters initial value -----------\n" 683 out += "Name".ljust(20) + "Value".ljust(20) + "Scale\n" 684 for name, value in self.results["init_parameters"].items(): 685 out += name.ljust(20) + str(value[0]).ljust(20) + str(value[1]) + "\n" 686 out += "-------- parameters final value --------------\n" 687 for name, value in self.results["parameters"].items(): 688 out += name.ljust(20) + f"{value:.6f}" 689 ierr = list(self.results["parameters"]).index(name) 690 err = self.results["parameter_errors"][ierr] 691 if err: 692 out += f" ± {err:.4f}" 693 out += "\n" 694 out += "Value at minimum".ljust(20)+ f'{self.results["min_value"]}\n' 695 out += "Iterations".ljust(20) + f'{self.results["iterations"]}\n' 696 out += "Max iterations".ljust(20) + f'{self.results["max_iterations"]}\n' 697 out += "Convergence flag".ljust(20)+ f'{self.results["convergence_flag"]}\n' 698 out += "Tolerance".ljust(20) + f'{self.results["tolerance"]}\n' 699 try: 700 arr = np.array2string( 701 self.compute_hessian(), 702 separator=', ', precision=6, suppress_small=True, 703 ) 704 out += "Hessian Matrix:\n" + arr 705 except: 706 out += "Hessian Matrix: (not available)" 707 return out 708 709 710########################################################### 711def andrews_curves(M, res=100) -> np.ndarray: 712 """ 713 Computes the [Andrews curves](https://en.wikipedia.org/wiki/Andrews_plot) 714 for the provided data. 715 716 The input array is an array of shape (n,m) where n is the number of 717 features and m is the number of observations. 718 719 Arguments: 720 M : (ndarray) 721 the data matrix (or data vector). 722 res : (int) 723 the resolution (n. of points) of the output curve. 724 725 Example: 726 - [andrews_cluster.py](https://github.com/marcomusy/vedo/blob/master/examples/pyplot/andrews_cluster.py) 727 728  729 """ 730 # Credits: 731 # https://gist.github.com/ryuzakyl/12c221ff0e54d8b1ac171c69ea552c0a 732 M = np.asarray(M) 733 m = int(res + 0.5) 734 735 # getting data vectors 736 X = np.reshape(M, (1, -1)) if len(M.shape) == 1 else M.copy() 737 _rows, n = X.shape 738 739 # andrews curve dimension (n. theta angles) 740 t = np.linspace(-np.pi, np.pi, m) 741 742 # m: range of values for angle theta 743 # n: amount of components of the Fourier expansion 744 A = np.empty((m, n)) 745 746 # setting first column of A 747 A[:, 0] = [1/np.sqrt(2)] * m 748 749 # filling columns of A 750 for i in range(1, n): 751 # computing the scaling coefficient for angle theta 752 c = np.ceil(i / 2) 753 # computing i-th column of matrix A 754 col = np.sin(c * t) if i % 2 == 1 else np.cos(c * t) 755 # setting column in matrix A 756 A[:, i] = col[:] 757 758 # computing Andrews curves for provided data 759 andrew_curves = np.dot(A, X.T).T 760 761 # returning the Andrews Curves (raveling if needed) 762 return np.ravel(andrew_curves) if andrew_curves.shape[0] == 1 else andrew_curves 763 764 765########################################################### 766def numpy2vtk(arr, dtype=None, deep=True, name=""): 767 """ 768 Convert a numpy array into a `vtkDataArray`. 769 Use `dtype='id'` for `vtkIdTypeArray` objects. 770 """ 771 # https://github.com/Kitware/VTK/blob/master/Wrapping/Python/vtkmodules/util/numpy_support.py 772 if arr is None: 773 return None 774 775 arr = np.ascontiguousarray(arr) 776 777 if dtype == "id": 778 if vtki.vtkIdTypeArray().GetDataTypeSize() != 4: 779 ast = np.int64 780 else: 781 ast = np.int32 782 varr = numpy_to_vtkIdTypeArray(arr.astype(ast), deep=deep) 783 elif dtype: 784 varr = numpy_to_vtk(arr.astype(dtype), deep=deep) 785 else: 786 # let numpy_to_vtk() decide what is best type based on arr type 787 if arr.dtype == np.bool_: 788 arr = arr.astype(np.uint8) 789 varr = numpy_to_vtk(arr, deep=deep) 790 791 if name: 792 varr.SetName(name) 793 return varr 794 795def vtk2numpy(varr): 796 """Convert a `vtkDataArray`, `vtkIdList` or `vtTransform` into a numpy array.""" 797 if varr is None: 798 return np.array([]) 799 if isinstance(varr, vtki.vtkIdList): 800 return np.array([varr.GetId(i) for i in range(varr.GetNumberOfIds())]) 801 elif isinstance(varr, vtki.vtkBitArray): 802 carr = vtki.vtkCharArray() 803 carr.DeepCopy(varr) 804 varr = carr 805 elif isinstance(varr, vtki.vtkHomogeneousTransform): 806 try: 807 varr = varr.GetMatrix() 808 except AttributeError: 809 pass 810 M = [[varr.GetElement(i, j) for j in range(4)] for i in range(4)] 811 return np.array(M) 812 return vtk_to_numpy(varr) 813 814 815def make3d(pts) -> np.ndarray: 816 """ 817 Make an array which might be 2D to 3D. 818 819 Array can also be in the form `[allx, ally, allz]`. 820 """ 821 if pts is None: 822 return np.array([]) 823 pts = np.asarray(pts) 824 825 if pts.dtype == "object": 826 raise ValueError("Cannot form a valid numpy array, input may be non-homogenous") 827 828 if pts.size == 0: # empty list 829 return pts 830 831 if pts.ndim == 1: 832 if pts.shape[0] == 2: 833 return np.hstack([pts, [0]]).astype(pts.dtype) 834 elif pts.shape[0] == 3: 835 return pts 836 else: 837 raise ValueError 838 839 if pts.shape[1] == 3: 840 return pts 841 842 # if 2 <= pts.shape[0] <= 3 and pts.shape[1] > 3: 843 # pts = pts.T 844 845 if pts.shape[1] == 2: 846 return np.c_[pts, np.zeros(pts.shape[0], dtype=pts.dtype)] 847 848 if pts.shape[1] != 3: 849 raise ValueError(f"input shape is not supported: {pts.shape}") 850 return pts 851 852 853def geometry(obj, extent=None) -> "vedo.Mesh": 854 """ 855 Apply the `vtkGeometryFilter` to the input object. 856 This is a general-purpose filter to extract geometry (and associated data) 857 from any type of dataset. 858 This filter also may be used to convert any type of data to polygonal type. 859 The conversion process may be less than satisfactory for some 3D datasets. 860 For example, this filter will extract the outer surface of a volume 861 or structured grid dataset. 862 863 Returns a `vedo.Mesh` object. 864 865 Set `extent` as the `[xmin,xmax, ymin,ymax, zmin,zmax]` bounding box to clip data. 866 """ 867 gf = vtki.new("GeometryFilter") 868 gf.SetInputData(obj) 869 if extent is not None: 870 gf.SetExtent(extent) 871 gf.Update() 872 return vedo.Mesh(gf.GetOutput()) 873 874 875def buildPolyData(vertices, faces=None, lines=None, strips=None, index_offset=0) -> vtki.vtkPolyData: 876 """ 877 Build a `vtkPolyData` object from a list of vertices 878 where faces represents the connectivity of the polygonal mesh. 879 Lines and triangle strips can also be specified. 880 881 E.g. : 882 - `vertices=[[x1,y1,z1],[x2,y2,z2], ...]` 883 - `faces=[[0,1,2], [1,2,3], ...]` 884 - `lines=[[0,1], [1,2,3,4], ...]` 885 - `strips=[[0,1,2,3,4,5], [2,3,9,7,4], ...]` 886 887 A flat list of faces can be passed as `faces=[3, 0,1,2, 4, 1,2,3,4, ...]`. 888 For lines use `lines=[2, 0,1, 4, 1,2,3,4, ...]`. 889 890 Use `index_offset=1` if face numbering starts from 1 instead of 0. 891 """ 892 if is_sequence(faces) and len(faces) == 0: 893 faces=None 894 if is_sequence(lines) and len(lines) == 0: 895 lines=None 896 if is_sequence(strips) and len(strips) == 0: 897 strips=None 898 899 poly = vtki.vtkPolyData() 900 901 if len(vertices) == 0: 902 return poly 903 904 vertices = make3d(vertices) 905 source_points = vtki.vtkPoints() 906 if vedo.settings.force_single_precision_points: 907 source_points.SetData(numpy2vtk(vertices, dtype=np.float32)) 908 else: 909 source_points.SetData(numpy2vtk(vertices)) 910 poly.SetPoints(source_points) 911 912 if lines is not None: 913 # Create a cell array to store the lines in and add the lines to it 914 linesarr = vtki.vtkCellArray() 915 if is_sequence(lines[0]): # assume format [(id0,id1),..] 916 for iline in lines: 917 for i in range(0, len(iline) - 1): 918 i1, i2 = iline[i], iline[i + 1] 919 if i1 != i2: 920 vline = vtki.vtkLine() 921 vline.GetPointIds().SetId(0, i1) 922 vline.GetPointIds().SetId(1, i2) 923 linesarr.InsertNextCell(vline) 924 else: # assume format [id0,id1,...] 925 # print("buildPolyData: assuming lines format [id0,id1,...]", lines) 926 # TODO CORRECT THIS CASE, MUST BE [2, id0,id1,...] 927 for i in range(0, len(lines) - 1): 928 vline = vtki.vtkLine() 929 vline.GetPointIds().SetId(0, lines[i]) 930 vline.GetPointIds().SetId(1, lines[i + 1]) 931 linesarr.InsertNextCell(vline) 932 poly.SetLines(linesarr) 933 934 if faces is not None: 935 936 source_polygons = vtki.vtkCellArray() 937 938 if isinstance(faces, np.ndarray) or not is_ragged(faces): 939 ##### all faces are composed of equal nr of vtxs, FAST 940 faces = np.asarray(faces) 941 if vtki.vtkIdTypeArray().GetDataTypeSize() != 4: 942 ast = np.int64 943 else: 944 ast = np.int32 945 946 if faces.ndim > 1: 947 nf, nc = faces.shape 948 hs = np.hstack((np.zeros(nf)[:, None] + nc, faces)) 949 else: 950 nf = faces.shape[0] 951 hs = faces 952 arr = numpy_to_vtkIdTypeArray(hs.astype(ast).ravel(), deep=True) 953 source_polygons.SetCells(nf, arr) 954 955 else: 956 ############################# manually add faces, SLOW 957 for f in faces: 958 n = len(f) 959 960 if n == 3: 961 tri = vtki.vtkTriangle() 962 pids = tri.GetPointIds() 963 for i in range(3): 964 pids.SetId(i, f[i] - index_offset) 965 source_polygons.InsertNextCell(tri) 966 967 else: 968 ele = vtki.vtkPolygon() 969 pids = ele.GetPointIds() 970 pids.SetNumberOfIds(n) 971 for i in range(n): 972 pids.SetId(i, f[i] - index_offset) 973 source_polygons.InsertNextCell(ele) 974 975 poly.SetPolys(source_polygons) 976 977 if strips is not None: 978 tscells = vtki.vtkCellArray() 979 for strip in strips: 980 # create a triangle strip 981 # https://vtk.org/doc/nightly/html/classvtkTriangleStrip.html 982 n = len(strip) 983 tstrip = vtki.vtkTriangleStrip() 984 tstrip_ids = tstrip.GetPointIds() 985 tstrip_ids.SetNumberOfIds(n) 986 for i in range(n): 987 tstrip_ids.SetId(i, strip[i] - index_offset) 988 tscells.InsertNextCell(tstrip) 989 poly.SetStrips(tscells) 990 991 if faces is None and lines is None and strips is None: 992 source_vertices = vtki.vtkCellArray() 993 for i in range(len(vertices)): 994 source_vertices.InsertNextCell(1) 995 source_vertices.InsertCellPoint(i) 996 poly.SetVerts(source_vertices) 997 998 # print("buildPolyData \n", 999 # poly.GetNumberOfPoints(), 1000 # poly.GetNumberOfCells(), # grand total 1001 # poly.GetNumberOfLines(), 1002 # poly.GetNumberOfPolys(), 1003 # poly.GetNumberOfStrips(), 1004 # poly.GetNumberOfVerts(), 1005 # ) 1006 return poly 1007 1008 1009############################################################################## 1010def get_font_path(font: str) -> str: 1011 """Internal use.""" 1012 if font in vedo.settings.font_parameters.keys(): 1013 if vedo.settings.font_parameters[font]["islocal"]: 1014 fl = os.path.join(vedo.fonts_path, f"{font}.ttf") 1015 else: 1016 try: 1017 fl = vedo.file_io.download(f"https://vedo.embl.es/fonts/{font}.ttf", verbose=False) 1018 except: 1019 vedo.logger.warning(f"Could not download https://vedo.embl.es/fonts/{font}.ttf") 1020 fl = os.path.join(vedo.fonts_path, "Normografo.ttf") 1021 else: 1022 if font.startswith("https://"): 1023 fl = vedo.file_io.download(font, verbose=False) 1024 elif os.path.isfile(font): 1025 fl = font # assume user is passing a valid file 1026 else: 1027 if font.endswith(".ttf"): 1028 vedo.logger.error( 1029 f"Could not set font file {font}" 1030 f"-> using default: {vedo.settings.default_font}" 1031 ) 1032 else: 1033 vedo.settings.default_font = "Normografo" 1034 vedo.logger.error( 1035 f"Could not set font name {font}" 1036 f" -> using default: Normografo\n" 1037 f"Check out https://vedo.embl.es/fonts for additional fonts\n" 1038 f"Type 'vedo -r fonts' to see available fonts" 1039 ) 1040 fl = get_font_path(vedo.settings.default_font) 1041 return fl 1042 1043 1044def is_sequence(arg) -> bool: 1045 """Check if the input is iterable.""" 1046 if hasattr(arg, "strip"): 1047 return False 1048 if hasattr(arg, "__getslice__"): 1049 return True 1050 if hasattr(arg, "__iter__"): 1051 return True 1052 return False 1053 1054 1055def is_ragged(arr, deep=False) -> bool: 1056 """ 1057 A ragged or inhomogeneous array in Python is an array 1058 with arrays of different lengths as its elements. 1059 To check if an array is ragged, we iterate through the elements 1060 and check if their lengths are the same. 1061 1062 Example: 1063 ```python 1064 arr = [[1, 2, 3], [[4, 5], [6], 1], [7, 8, 9]] 1065 print(is_ragged(arr, deep=True)) # output: True 1066 ``` 1067 """ 1068 n = len(arr) 1069 if n == 0: 1070 return False 1071 if is_sequence(arr[0]): 1072 length = len(arr[0]) 1073 for i in range(1, n): 1074 if len(arr[i]) != length or (deep and is_ragged(arr[i])): 1075 return True 1076 return False 1077 return False 1078 1079 1080def flatten(list_to_flatten) -> list: 1081 """Flatten out a list.""" 1082 1083 def _genflatten(lst): 1084 for elem in lst: 1085 if isinstance(elem, (list, tuple)): 1086 for x in flatten(elem): 1087 yield x 1088 else: 1089 yield elem 1090 1091 return list(_genflatten(list_to_flatten)) 1092 1093 1094def humansort(alist) -> list: 1095 """ 1096 Sort in place a given list the way humans expect. 1097 1098 E.g. `['file11', 'file1'] -> ['file1', 'file11']` 1099 1100 .. warning:: input list is modified in-place by this function. 1101 """ 1102 import re 1103 1104 def alphanum_key(s): 1105 # Turn a string into a list of string and number chunks. 1106 # e.g. "z23a" -> ["z", 23, "a"] 1107 def tryint(s): 1108 if s.isdigit(): 1109 return int(s) 1110 return s 1111 1112 return [tryint(c) for c in re.split("([0-9]+)", s)] 1113 1114 alist.sort(key=alphanum_key) 1115 return alist # NB: input list is modified 1116 1117 1118def sort_by_column(arr, nth, invert=False) -> np.ndarray: 1119 """Sort a numpy array by its `n-th` column.""" 1120 arr = np.asarray(arr) 1121 arr = arr[arr[:, nth].argsort()] 1122 if invert: 1123 return np.flip(arr, axis=0) 1124 return arr 1125 1126 1127def point_in_triangle(p, p1, p2, p3) -> Union[bool, None]: 1128 """ 1129 Return True if a point is inside (or above/below) 1130 a triangle defined by 3 points in space. 1131 """ 1132 p1 = np.array(p1) 1133 u = p2 - p1 1134 v = p3 - p1 1135 n = np.cross(u, v) 1136 w = p - p1 1137 ln = np.dot(n, n) 1138 if not ln: 1139 return None # degenerate triangle 1140 gamma = (np.dot(np.cross(u, w), n)) / ln 1141 if 0 < gamma < 1: 1142 beta = (np.dot(np.cross(w, v), n)) / ln 1143 if 0 < beta < 1: 1144 alpha = 1 - gamma - beta 1145 if 0 < alpha < 1: 1146 return True 1147 return False 1148 1149 1150def intersection_ray_triangle(P0, P1, V0, V1, V2) -> Union[bool, None, np.ndarray]: 1151 """ 1152 Fast intersection between a directional ray defined by `P0,P1` 1153 and triangle `V0, V1, V2`. 1154 1155 Returns the intersection point or 1156 - `None` if triangle is degenerate, or ray is parallel to triangle plane. 1157 - `False` if no intersection, or ray direction points away from triangle. 1158 """ 1159 # Credits: http://geomalgorithms.com/a06-_intersect-2.html 1160 # Get triangle edge vectors and plane normal 1161 # todo : this is slow should check 1162 # https://vtk.org/doc/nightly/html/classvtkCell.html 1163 V0 = np.asarray(V0, dtype=float) 1164 P0 = np.asarray(P0, dtype=float) 1165 u = V1 - V0 1166 v = V2 - V0 1167 n = np.cross(u, v) 1168 if not np.abs(v).sum(): # triangle is degenerate 1169 return None # do not deal with this case 1170 1171 rd = P1 - P0 # ray direction vector 1172 w0 = P0 - V0 1173 a = -np.dot(n, w0) 1174 b = np.dot(n, rd) 1175 if not b: # ray is parallel to triangle plane 1176 return None 1177 1178 # Get intersect point of ray with triangle plane 1179 r = a / b 1180 if r < 0.0: # ray goes away from triangle 1181 return False # => no intersect 1182 1183 # Gor a segment, also test if (r > 1.0) => no intersect 1184 I = P0 + r * rd # intersect point of ray and plane 1185 1186 # is I inside T? 1187 uu = np.dot(u, u) 1188 uv = np.dot(u, v) 1189 vv = np.dot(v, v) 1190 w = I - V0 1191 wu = np.dot(w, u) 1192 wv = np.dot(w, v) 1193 D = uv * uv - uu * vv 1194 1195 # Get and test parametric coords 1196 s = (uv * wv - vv * wu) / D 1197 if s < 0.0 or s > 1.0: # I is outside T 1198 return False 1199 t = (uv * wu - uu * wv) / D 1200 if t < 0.0 or (s + t) > 1.0: # I is outside T 1201 return False 1202 return I # I is in T 1203 1204 1205def triangle_solver(**input_dict): 1206 """ 1207 Solve a triangle from any 3 known elements. 1208 (Note that there might be more than one solution or none). 1209 Angles are in radians. 1210 1211 Example: 1212 ```python 1213 print(triangle_solver(a=3, b=4, c=5)) 1214 print(triangle_solver(a=3, ac=0.9273, ab=1.5716)) 1215 print(triangle_solver(a=3, b=4, ab=1.5716)) 1216 print(triangle_solver(b=4, bc=.64, ab=1.5716)) 1217 print(triangle_solver(c=5, ac=.9273, bc=0.6435)) 1218 print(triangle_solver(a=3, c=5, bc=0.6435)) 1219 print(triangle_solver(b=4, c=5, ac=0.927)) 1220 ``` 1221 """ 1222 a = input_dict.get("a") 1223 b = input_dict.get("b") 1224 c = input_dict.get("c") 1225 ab = input_dict.get("ab") 1226 bc = input_dict.get("bc") 1227 ac = input_dict.get("ac") 1228 1229 if ab and bc: 1230 ac = np.pi - bc - ab 1231 elif bc and ac: 1232 ab = np.pi - bc - ac 1233 elif ab and ac: 1234 bc = np.pi - ab - ac 1235 1236 if a is not None and b is not None and c is not None: 1237 ab = np.arccos((a ** 2 + b ** 2 - c ** 2) / (2 * a * b)) 1238 sinab = np.sin(ab) 1239 ac = np.arcsin(a / c * sinab) 1240 bc = np.arcsin(b / c * sinab) 1241 1242 elif a is not None and b is not None and ab is not None: 1243 c = np.sqrt(a ** 2 + b ** 2 - 2 * a * b * np.cos(ab)) 1244 sinab = np.sin(ab) 1245 ac = np.arcsin(a / c * sinab) 1246 bc = np.arcsin(b / c * sinab) 1247 1248 elif a is not None and ac is not None and ab is not None: 1249 h = a * np.sin(ac) 1250 b = h / np.sin(bc) 1251 c = b * np.cos(bc) + a * np.cos(ac) 1252 1253 elif b is not None and bc is not None and ab is not None: 1254 h = b * np.sin(bc) 1255 a = h / np.sin(ac) 1256 c = np.sqrt(a * a + b * b) 1257 1258 elif c is not None and ac is not None and bc is not None: 1259 h = c * np.sin(bc) 1260 b1 = c * np.cos(bc) 1261 b2 = h / np.tan(ab) 1262 b = b1 + b2 1263 a = np.sqrt(b2 * b2 + h * h) 1264 1265 elif a is not None and c is not None and bc is not None: 1266 # double solution 1267 h = c * np.sin(bc) 1268 k = np.sqrt(a * a - h * h) 1269 omega = np.arcsin(k / a) 1270 cosbc = np.cos(bc) 1271 b = c * cosbc - k 1272 phi = np.pi / 2 - bc - omega 1273 ac = phi 1274 ab = np.pi - ac - bc 1275 if k: 1276 b2 = c * cosbc + k 1277 ac2 = phi + 2 * omega 1278 ab2 = np.pi - ac2 - bc 1279 return [ 1280 {"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}, 1281 {"a": a, "b": b2, "c": c, "ab": ab2, "bc": bc, "ac": ac2}, 1282 ] 1283 1284 elif b is not None and c is not None and ac is not None: 1285 # double solution 1286 h = c * np.sin(ac) 1287 k = np.sqrt(b * b - h * h) 1288 omega = np.arcsin(k / b) 1289 cosac = np.cos(ac) 1290 a = c * cosac - k 1291 phi = np.pi / 2 - ac - omega 1292 bc = phi 1293 ab = np.pi - bc - ac 1294 if k: 1295 a2 = c * cosac + k 1296 bc2 = phi + 2 * omega 1297 ab2 = np.pi - ac - bc2 1298 return [ 1299 {"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}, 1300 {"a": a2, "b": b, "c": c, "ab": ab2, "bc": bc2, "ac": ac}, 1301 ] 1302 1303 else: 1304 vedo.logger.error(f"Case {input_dict} is not supported.") 1305 return [] 1306 1307 return [{"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}] 1308 1309 1310############################################################################# 1311def circle_from_3points(p1, p2, p3) -> np.ndarray: 1312 """ 1313 Find the center and radius of a circle given 3 points in 3D space. 1314 1315 Returns the center of the circle. 1316 1317 Example: 1318 ```python 1319 from vedo.utils import mag, circle_from_3points 1320 p1 = [0,1,1] 1321 p2 = [3,0,1] 1322 p3 = [1,2,0] 1323 c = circle_from_3points(p1, p2, p3) 1324 print(mag(c-p1), mag(c-p2), mag(c-p3)) 1325 ``` 1326 """ 1327 p1 = np.asarray(p1) 1328 p2 = np.asarray(p2) 1329 p3 = np.asarray(p3) 1330 v1 = p2 - p1 1331 v2 = p3 - p1 1332 v11 = np.dot(v1, v1) 1333 v22 = np.dot(v2, v2) 1334 v12 = np.dot(v1, v2) 1335 f = 1.0 / (2 * (v11 * v22 - v12 * v12)) 1336 k1 = f * v22 * (v11-v12) 1337 k2 = f * v11 * (v22-v12) 1338 return p1 + k1 * v1 + k2 * v2 1339 1340def point_line_distance(p, p1, p2) -> float: 1341 """ 1342 Compute the distance of a point to a line (not the segment) 1343 defined by `p1` and `p2`. 1344 """ 1345 return np.sqrt(vtki.vtkLine.DistanceToLine(p, p1, p2)) 1346 1347def line_line_distance(p1, p2, q1, q2) -> Tuple[float, np.ndarray, np.ndarray, float, float]: 1348 """ 1349 Compute the distance of a line to a line (not the segment) 1350 defined by `p1` and `p2` and `q1` and `q2`. 1351 1352 Returns the distance, 1353 the closest point on line 1, the closest point on line 2. 1354 Their parametric coords (-inf <= t0, t1 <= inf) are also returned. 1355 """ 1356 closest_pt1: MutableSequence[float] = [0,0,0] 1357 closest_pt2: MutableSequence[float] = [0,0,0] 1358 t1, t2 = 0.0, 0.0 1359 d = vtki.vtkLine.DistanceBetweenLines( 1360 p1, p2, q1, q2, closest_pt1, closest_pt2, t1, t2) 1361 return np.sqrt(d), closest_pt1, closest_pt2, t1, t2 1362 1363def segment_segment_distance(p1, p2, q1, q2): 1364 """ 1365 Compute the distance of a segment to a segment 1366 defined by `p1` and `p2` and `q1` and `q2`. 1367 1368 Returns the distance, 1369 the closest point on line 1, the closest point on line 2. 1370 Their parametric coords (-inf <= t0, t1 <= inf) are also returned. 1371 """ 1372 closest_pt1 = [0,0,0] 1373 closest_pt2 = [0,0,0] 1374 t1, t2 = 0.0, 0.0 1375 d = vtki.vtkLine.DistanceBetweenLineSegments( 1376 p1, p2, q1, q2, closest_pt1, closest_pt2, t1, t2) 1377 return np.sqrt(d), closest_pt1, closest_pt2, t1, t2 1378 1379 1380def closest(point, points, n=1, return_ids=False, use_tree=False): 1381 """ 1382 Returns the distances and the closest point(s) to the given set of points. 1383 Needs `scipy.spatial` library. 1384 1385 Arguments: 1386 n : (int) 1387 the nr of closest points to return 1388 return_ids : (bool) 1389 return the ids instead of the points coordinates 1390 use_tree : (bool) 1391 build a `scipy.spatial.KDTree`. 1392 An already existing one can be passed to avoid rebuilding. 1393 """ 1394 from scipy.spatial import distance, KDTree 1395 1396 points = np.asarray(points) 1397 if n == 1: 1398 dists = distance.cdist([point], points) 1399 closest_idx = np.argmin(dists) 1400 else: 1401 if use_tree: 1402 if isinstance(use_tree, KDTree): # reuse 1403 tree = use_tree 1404 else: 1405 tree = KDTree(points) 1406 dists, closest_idx = tree.query([point], k=n) 1407 closest_idx = closest_idx[0] 1408 else: 1409 dists = distance.cdist([point], points) 1410 closest_idx = np.argsort(dists)[0][:n] 1411 if return_ids: 1412 return dists, closest_idx 1413 else: 1414 return dists, points[closest_idx] 1415 1416 1417############################################################################# 1418def lin_interpolate(x, rangeX, rangeY): 1419 """ 1420 Interpolate linearly the variable `x` in `rangeX` onto the new `rangeY`. 1421 If `x` is a 3D vector the linear weight is the distance to the two 3D `rangeX` vectors. 1422 1423 E.g. if `x` runs in `rangeX=[x0,x1]` and I want it to run in `rangeY=[y0,y1]` then 1424 1425 `y = lin_interpolate(x, rangeX, rangeY)` will interpolate `x` onto `rangeY`. 1426 1427 Examples: 1428 - [lin_interpolate.py](https://github.com/marcomusy/vedo/tree/master/examples/basic/lin_interpolate.py) 1429 1430  1431 """ 1432 if is_sequence(x): 1433 x = np.asarray(x) 1434 x0, x1 = np.asarray(rangeX) 1435 y0, y1 = np.asarray(rangeY) 1436 dx = x1 - x0 1437 dxn = np.linalg.norm(dx) 1438 if not dxn: 1439 return y0 1440 s = np.linalg.norm(x - x0) / dxn 1441 t = np.linalg.norm(x - x1) / dxn 1442 st = s + t 1443 out = y0 * (t / st) + y1 * (s / st) 1444 1445 else: # faster 1446 1447 x0 = rangeX[0] 1448 dx = rangeX[1] - x0 1449 if not dx: 1450 return rangeY[0] 1451 s = (x - x0) / dx 1452 out = rangeY[0] * (1 - s) + rangeY[1] * s 1453 return out 1454 1455 1456def get_uv(p, x, v): 1457 """ 1458 Obtain the texture uv-coords of a point p belonging to a face that has point 1459 coordinates (x0, x1, x2) with the corresponding uv-coordinates v=(v0, v1, v2). 1460 All p and x0,x1,x2 are 3D-vectors, while v are their 2D uv-coordinates. 1461 1462 Example: 1463 ```python 1464 from vedo import * 1465 1466 pic = Image(dataurl+"coloured_cube_faces.jpg") 1467 cb = Mesh(dataurl+"coloured_cube.obj").lighting("off").texture(pic) 1468 1469 cbpts = cb.vertices 1470 faces = cb.cells 1471 uv = cb.pointdata["Material"] 1472 1473 pt = [-0.2, 0.75, 2] 1474 pr = cb.closest_point(pt) 1475 1476 idface = cb.closest_point(pt, return_cell_id=True) 1477 idpts = faces[idface] 1478 uv_face = uv[idpts] 1479 1480 uv_pr = utils.get_uv(pr, cbpts[idpts], uv_face) 1481 print("interpolated uv =", uv_pr) 1482 1483 sx, sy = pic.dimensions() 1484 i_interp_uv = uv_pr * [sy, sx] 1485 ix, iy = i_interp_uv.astype(int) 1486 mpic = pic.tomesh() 1487 rgba = mpic.pointdata["RGBA"].reshape(sy, sx, 3) 1488 print("color =", rgba[ix, iy]) 1489 1490 show( 1491 [[cb, Point(pr), cb.labels("Material")], 1492 [pic, Point(i_interp_uv)]], 1493 N=2, axes=1, sharecam=False, 1494 ).close() 1495 ``` 1496  1497 """ 1498 # Vector vp=p-x0 is representable as alpha*s + beta*t, 1499 # where s = x1-x0 and t = x2-x0, in matrix form 1500 # vp = [alpha, beta] . matrix(s,t) 1501 # M = matrix(s,t) is 2x3 matrix, so (alpha, beta) can be found by 1502 # inverting any of its minor A with non-zero determinant. 1503 # Once found, uv-coords of p are vt0 + alpha (vt1-v0) + beta (vt2-v0) 1504 1505 p = np.asarray(p) 1506 x0, x1, x2 = np.asarray(x)[:3] 1507 vt0, vt1, vt2 = np.asarray(v)[:3] 1508 1509 s = x1 - x0 1510 t = x2 - x0 1511 vs = vt1 - vt0 1512 vt = vt2 - vt0 1513 vp = p - x0 1514 1515 # finding a minor with independent rows 1516 M = np.matrix([s, t]) 1517 mnr = [0, 1] 1518 A = M[:, mnr] 1519 if np.abs(np.linalg.det(A)) < 0.000001: 1520 mnr = [0, 2] 1521 A = M[:, mnr] 1522 if np.abs(np.linalg.det(A)) < 0.000001: 1523 mnr = [1, 2] 1524 A = M[:, mnr] 1525 Ainv = np.linalg.inv(A) 1526 alpha_beta = vp[mnr].dot(Ainv) # [alpha, beta] 1527 return np.asarray(vt0 + alpha_beta.dot(np.matrix([vs, vt])))[0] 1528 1529 1530def vector(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1531 """ 1532 Return a 3D numpy array representing a vector. 1533 1534 If `y` is `None`, assume input is already in the form `[x,y,z]`. 1535 """ 1536 if y is None: # assume x is already [x,y,z] 1537 return np.asarray(x, dtype=dtype) 1538 return np.array([x, y, z], dtype=dtype) 1539 1540 1541def versor(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1542 """Return the unit vector. Input can be a list of vectors.""" 1543 v = vector(x, y, z, dtype) 1544 if isinstance(v[0], np.ndarray): 1545 return np.divide(v, mag(v)[:, None]) 1546 return v / mag(v) 1547 1548 1549def mag(v): 1550 """Get the magnitude of a vector or array of vectors.""" 1551 v = np.asarray(v) 1552 if v.ndim == 1: 1553 return np.linalg.norm(v) 1554 return np.linalg.norm(v, axis=1) 1555 1556 1557def mag2(v) -> np.ndarray: 1558 """Get the squared magnitude of a vector or array of vectors.""" 1559 v = np.asarray(v) 1560 if v.ndim == 1: 1561 return np.square(v).sum() 1562 return np.square(v).sum(axis=1) 1563 1564 1565def is_integer(n) -> bool: 1566 """Check if input is an integer.""" 1567 try: 1568 float(n) 1569 except (ValueError, TypeError): 1570 return False 1571 else: 1572 return float(n).is_integer() 1573 1574 1575def is_number(n) -> bool: 1576 """Check if input is a number""" 1577 try: 1578 float(n) 1579 return True 1580 except (ValueError, TypeError): 1581 return False 1582 1583 1584def round_to_digit(x, p) -> float: 1585 """Round a real number to the specified number of significant digits.""" 1586 if not x: 1587 return 0 1588 r = np.round(x, p - int(np.floor(np.log10(abs(x)))) - 1) 1589 if int(r) == r: 1590 return int(r) 1591 return r 1592 1593 1594def pack_spheres(bounds, radius) -> np.ndarray: 1595 """ 1596 Packing spheres into a bounding box. 1597 Returns a numpy array of sphere centers. 1598 """ 1599 h = 0.8164965 / 2 1600 d = 0.8660254 1601 a = 0.288675135 1602 1603 if is_sequence(bounds): 1604 x0, x1, y0, y1, z0, z1 = bounds 1605 else: 1606 x0, x1, y0, y1, z0, z1 = bounds.bounds() 1607 1608 x = np.arange(x0, x1, radius) 1609 nul = np.zeros_like(x) 1610 nz = int((z1 - z0) / radius / h / 2 + 1.5) 1611 ny = int((y1 - y0) / radius / d + 1.5) 1612 1613 pts = [] 1614 for iz in range(nz): 1615 z = z0 + nul + iz * h * radius 1616 dx, dy, dz = [radius * 0.5, radius * a, iz * h * radius] 1617 for iy in range(ny): 1618 y = y0 + nul + iy * d * radius 1619 if iy % 2: 1620 xs = x 1621 else: 1622 xs = x + radius * 0.5 1623 if iz % 2: 1624 p = np.c_[xs, y, z] + [dx, dy, dz] 1625 else: 1626 p = np.c_[xs, y, z] + [0, 0, dz] 1627 pts += p.tolist() 1628 return np.array(pts) 1629 1630 1631def precision(x, p: int, vrange=None, delimiter="e") -> str: 1632 """ 1633 Returns a string representation of `x` formatted to precision `p`. 1634 1635 Set `vrange` to the range in which x exists (to snap x to '0' if below precision). 1636 """ 1637 # Based on the webkit javascript implementation 1638 # `from here <https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp>`_, 1639 # and implemented by `randlet <https://github.com/randlet/to-precision>`_. 1640 # Modified for vedo by M.Musy 2020 1641 1642 if isinstance(x, str): # do nothing 1643 return x 1644 1645 if is_sequence(x): 1646 out = "(" 1647 nn = len(x) - 1 1648 for i, ix in enumerate(x): 1649 1650 try: 1651 if np.isnan(ix): 1652 return "NaN" 1653 except: 1654 # cannot handle list of list 1655 continue 1656 1657 out += precision(ix, p) 1658 if i < nn: 1659 out += ", " 1660 return out + ")" ############ <-- 1661 1662 try: 1663 if np.isnan(x): 1664 return "NaN" 1665 except TypeError: 1666 return "NaN" 1667 1668 x = float(x) 1669 1670 if x == 0.0 or (vrange is not None and abs(x) < vrange / pow(10, p)): 1671 return "0" 1672 1673 out = [] 1674 if x < 0: 1675 out.append("-") 1676 x = -x 1677 1678 e = int(np.log10(x)) 1679 # tens = np.power(10, e - p + 1) 1680 tens = 10 ** (e - p + 1) 1681 n = np.floor(x / tens) 1682 1683 # if n < np.power(10, p - 1): 1684 if n < 10 ** (p - 1): 1685 e = e - 1 1686 # tens = np.power(10, e - p + 1) 1687 tens = 10 ** (e - p + 1) 1688 n = np.floor(x / tens) 1689 1690 if abs((n + 1.0) * tens - x) <= abs(n * tens - x): 1691 n = n + 1 1692 1693 # if n >= np.power(10, p): 1694 if n >= 10 ** p: 1695 n = n / 10.0 1696 e = e + 1 1697 1698 m = "%.*g" % (p, n) 1699 if e < -2 or e >= p: 1700 out.append(m[0]) 1701 if p > 1: 1702 out.append(".") 1703 out.extend(m[1:p]) 1704 out.append(delimiter) 1705 if e > 0: 1706 out.append("+") 1707 out.append(str(e)) 1708 elif e == (p - 1): 1709 out.append(m) 1710 elif e >= 0: 1711 out.append(m[: e + 1]) 1712 if e + 1 < len(m): 1713 out.append(".") 1714 out.extend(m[e + 1 :]) 1715 else: 1716 out.append("0.") 1717 out.extend(["0"] * -(e + 1)) 1718 out.append(m) 1719 return "".join(out) 1720 1721 1722################################################################################## 1723def grep(filename: str, tag: str, column=None, first_occurrence_only=False) -> list: 1724 """Greps the line in a file that starts with a specific `tag` string inside the file.""" 1725 import re 1726 1727 with open(filename, "r", encoding="UTF-8") as afile: 1728 content = [] 1729 for line in afile: 1730 if re.search(tag, line): 1731 c = line.split() 1732 c[-1] = c[-1].replace("\n", "") 1733 if column is not None: 1734 c = c[column] 1735 content.append(c) 1736 if first_occurrence_only: 1737 break 1738 return content 1739 1740def parse_pattern(query, strings_to_parse) -> list: 1741 """ 1742 Parse a pattern query to a list of strings. 1743 The query string can contain wildcards like * and ?. 1744 1745 Arguments: 1746 query : (str) 1747 the query to parse 1748 strings_to_parse : (str/list) 1749 the string or list of strings to parse 1750 1751 Returns: 1752 a list of booleans, one for each string in strings_to_parse 1753 1754 Example: 1755 >>> query = r'*Sphere 1?3*' 1756 >>> strings = ["Sphere 143 red", "Sphere 13 red", "Sphere 123", "ASphere 173"] 1757 >>> parse_pattern(query, strings) 1758 [True, True, False, False] 1759 """ 1760 from re import findall as re_findall 1761 if not isinstance(query, str): 1762 return [False] 1763 1764 if not is_sequence(strings_to_parse): 1765 strings_to_parse = [strings_to_parse] 1766 1767 outs = [] 1768 for sp in strings_to_parse: 1769 if not isinstance(sp, str): 1770 outs.append(False) 1771 continue 1772 1773 s = query 1774 if s.startswith("*"): 1775 s = s[1:] 1776 else: 1777 s = "^" + s 1778 1779 t = "" 1780 if not s.endswith("*"): 1781 t = "$" 1782 else: 1783 s = s[:-1] 1784 1785 pattern = s.replace('?', r'\w').replace(' ', r'\s').replace("*", r"\w+") + t 1786 1787 # Search for the pattern in the input string 1788 match = re_findall(pattern, sp) 1789 out = bool(match) 1790 outs.append(out) 1791 # Print the matches for debugging 1792 print("pattern", pattern, "in:", strings_to_parse) 1793 print("matches", match, "result:", out) 1794 return outs 1795 1796def print_histogram( 1797 data, 1798 bins=10, 1799 height=10, 1800 logscale=False, 1801 minbin=0, 1802 horizontal=True, 1803 char="\U00002589", 1804 c=None, 1805 bold=True, 1806 title="histogram", 1807 spacer="", 1808) -> np.ndarray: 1809 """ 1810 Ascii histogram printing. 1811 1812 Input can be a `vedo.Volume` or `vedo.Mesh`. 1813 Returns the raw data before binning (useful when passing vtk objects). 1814 1815 Arguments: 1816 bins : (int) 1817 number of histogram bins 1818 height : (int) 1819 height of the histogram in character units 1820 logscale : (bool) 1821 use logscale for frequencies 1822 minbin : (int) 1823 ignore bins before minbin 1824 horizontal : (bool) 1825 show histogram horizontally 1826 char : (str) 1827 character to be used 1828 bold : (bool) 1829 use boldface 1830 title : (str) 1831 histogram title 1832 spacer : (str) 1833 horizontal spacer 1834 1835 Example: 1836 ```python 1837 from vedo import print_histogram 1838 import numpy as np 1839 d = np.random.normal(size=1000) 1840 data = print_histogram(d, c='b', logscale=True, title='my scalars') 1841 data = print_histogram(d, c='o') 1842 print(np.mean(data)) # data here is same as d 1843 ``` 1844  1845 """ 1846 # credits: http://pyinsci.blogspot.com/2009/10/ascii-histograms.html 1847 # adapted for vedo by M.Musy, 2019 1848 1849 if not horizontal: # better aspect ratio 1850 bins *= 2 1851 1852 try: 1853 data = vtk2numpy(data.dataset.GetPointData().GetScalars()) 1854 except AttributeError: 1855 # already an array 1856 data = np.asarray(data) 1857 1858 if isinstance(data, vtki.vtkImageData): 1859 dims = data.GetDimensions() 1860 nvx = min(100000, dims[0] * dims[1] * dims[2]) 1861 idxs = np.random.randint(0, min(dims), size=(nvx, 3)) 1862 data = [] 1863 for ix, iy, iz in idxs: 1864 d = data.GetScalarComponentAsFloat(ix, iy, iz, 0) 1865 data.append(d) 1866 data = np.array(data) 1867 1868 elif isinstance(data, vtki.vtkPolyData): 1869 arr = data.GetPointData().GetScalars() 1870 if not arr: 1871 arr = data.GetCellData().GetScalars() 1872 if not arr: 1873 return np.array([]) 1874 data = vtk2numpy(arr) 1875 1876 try: 1877 h = np.histogram(data, bins=bins) 1878 except TypeError as e: 1879 vedo.logger.error(f"cannot compute histogram: {e}") 1880 return np.array([]) 1881 1882 if minbin: 1883 hi = h[0][minbin:-1] 1884 else: 1885 hi = h[0] 1886 1887 if char == "\U00002589" and horizontal: 1888 char = "\U00002586" 1889 1890 title = title.ljust(14) + ":" 1891 entrs = " entries=" + str(len(data)) 1892 if logscale: 1893 h0 = np.log10(hi + 1) 1894 maxh0 = int(max(h0) * 100) / 100 1895 title = title + entrs + " (logscale)" 1896 else: 1897 h0 = hi 1898 maxh0 = max(h0) 1899 title = title + entrs 1900 1901 def _v(): 1902 his = "" 1903 if title: 1904 his += title + "\n" 1905 bars = h0 / maxh0 * height 1906 for l in reversed(range(1, height + 1)): 1907 line = "" 1908 if l == height: 1909 line = "%s " % maxh0 1910 else: 1911 line = " |" + " " * (len(str(maxh0)) - 3) 1912 for c in bars: 1913 if c >= np.ceil(l): 1914 line += char 1915 else: 1916 line += " " 1917 line += "\n" 1918 his += line 1919 his += "%.2f" % h[1][0] + "." * (bins) + "%.2f" % h[1][-1] + "\n" 1920 return his 1921 1922 def _h(): 1923 his = "" 1924 if title: 1925 his += title + "\n" 1926 xl = ["%.2f" % n for n in h[1]] 1927 lxl = [len(l) for l in xl] 1928 bars = h0 / maxh0 * height 1929 his += spacer + " " * int(max(bars) + 2 + max(lxl)) + "%s\n" % maxh0 1930 for i, c in enumerate(bars): 1931 line = xl[i] + " " * int(max(lxl) - lxl[i]) + "| " + char * int(c) + "\n" 1932 his += spacer + line 1933 return his 1934 1935 if horizontal: 1936 height *= 2 1937 vedo.printc(_h(), c=c, bold=bold) 1938 else: 1939 vedo.printc(_v(), c=c, bold=bold) 1940 return data 1941 1942 1943def print_table(*columns, headers=None, c="g") -> None: 1944 """ 1945 Print lists as tables. 1946 1947 Example: 1948 ```python 1949 from vedo.utils import print_table 1950 list1 = ["A", "B", "C"] 1951 list2 = [142, 220, 330] 1952 list3 = [True, False, True] 1953 headers = ["First Column", "Second Column", "Third Column"] 1954 print_table(list1, list2, list3, headers=headers) 1955 ``` 1956 1957  1958 """ 1959 # If headers is not provided, use default header names 1960 corner = "─" 1961 if headers is None: 1962 headers = [f"Column {i}" for i in range(1, len(columns) + 1)] 1963 assert len(headers) == len(columns) 1964 1965 # Find the maximum length of the elements in each column and header 1966 max_lens = [max(len(str(x)) for x in column) for column in columns] 1967 max_len_headers = [max(len(str(header)), max_len) for header, max_len in zip(headers, max_lens)] 1968 1969 # Construct the table header 1970 header = ( 1971 "│ " 1972 + " │ ".join(header.ljust(max_len) for header, max_len in zip(headers, max_len_headers)) 1973 + " │" 1974 ) 1975 1976 # Construct the line separator 1977 line1 = "┌" + corner.join("─" * (max_len + 2) for max_len in max_len_headers) + "┐" 1978 line2 = "└" + corner.join("─" * (max_len + 2) for max_len in max_len_headers) + "┘" 1979 1980 # Print the table header 1981 vedo.printc(line1, c=c) 1982 vedo.printc(header, c=c) 1983 vedo.printc(line2, c=c) 1984 1985 # Print the data rows 1986 for row in zip(*columns): 1987 row = ( 1988 "│ " 1989 + " │ ".join(str(col).ljust(max_len) for col, max_len in zip(row, max_len_headers)) 1990 + " │" 1991 ) 1992 vedo.printc(row, bold=False, c=c) 1993 1994 # Print the line separator again to close the table 1995 vedo.printc(line2, c=c) 1996 1997def print_inheritance_tree(C) -> None: 1998 """Prints the inheritance tree of class C.""" 1999 # Adapted from: https://stackoverflow.com/questions/26568976/ 2000 def class_tree(cls): 2001 subc = [class_tree(sub_class) for sub_class in cls.__subclasses__()] 2002 return {cls.__name__: subc} 2003 2004 def print_tree(tree, indent=8, current_ind=0): 2005 for k, v in tree.items(): 2006 if current_ind: 2007 before_dashes = current_ind - indent 2008 m = " " * before_dashes + "└" + "─" * (indent - 1) + " " + k 2009 vedo.printc(m) 2010 else: 2011 vedo.printc(k) 2012 for sub_tree in v: 2013 print_tree(sub_tree, indent=indent, current_ind=current_ind + indent) 2014 2015 if str(C.__class__) != "<class 'type'>": 2016 C = C.__class__ 2017 ct = class_tree(C) 2018 print_tree(ct) 2019 2020 2021def make_bands(inputlist, n): 2022 """ 2023 Group values of a list into bands of equal value, where 2024 `n` is the number of bands, a positive integer > 2. 2025 2026 Returns a binned list of the same length as the input. 2027 """ 2028 if n < 2: 2029 return inputlist 2030 vmin = np.min(inputlist) 2031 vmax = np.max(inputlist) 2032 bb = np.linspace(vmin, vmax, n, endpoint=0) 2033 dr = bb[1] - bb[0] 2034 bb += dr / 2 2035 tol = dr / 2 * 1.001 2036 newlist = [] 2037 for s in inputlist: 2038 for b in bb: 2039 if abs(s - b) < tol: 2040 newlist.append(b) 2041 break 2042 return np.array(newlist) 2043 2044 2045################################################################# 2046# Functions adapted from: 2047# https://github.com/sdorkenw/MeshParty/blob/master/meshparty/trimesh_vtk.py 2048def camera_from_quaternion(pos, quaternion, distance=10000, ngl_correct=True) -> vtki.vtkCamera: 2049 """ 2050 Define a `vtkCamera` with a particular orientation. 2051 2052 Arguments: 2053 pos: (np.array, list, tuple) 2054 an iterator of length 3 containing the focus point of the camera 2055 quaternion: (np.array, list, tuple) 2056 a len(4) quaternion `(x,y,z,w)` describing the rotation of the camera 2057 such as returned by neuroglancer `x,y,z,w` all in `[0,1]` range 2058 distance: (float) 2059 the desired distance from pos to the camera (default = 10000 nm) 2060 2061 Returns: 2062 `vtki.vtkCamera`, a vtk camera setup according to these rules. 2063 """ 2064 camera = vtki.vtkCamera() 2065 # define the quaternion in vtk, note the swapped order 2066 # w,x,y,z instead of x,y,z,w 2067 quat_vtk = vtki.get_class("Quaternion")( 2068 quaternion[3], quaternion[0], quaternion[1], quaternion[2]) 2069 # use this to define a rotation matrix in x,y,z 2070 # right handed units 2071 M = np.zeros((3, 3), dtype=np.float32) 2072 quat_vtk.ToMatrix3x3(M) 2073 # the default camera orientation is y up 2074 up = [0, 1, 0] 2075 # calculate default camera position is backed off in positive z 2076 pos = [0, 0, distance] 2077 2078 # set the camera rototation by applying the rotation matrix 2079 camera.SetViewUp(*np.dot(M, up)) 2080 # set the camera position by applying the rotation matrix 2081 camera.SetPosition(*np.dot(M, pos)) 2082 if ngl_correct: 2083 # neuroglancer has positive y going down 2084 # so apply these azimuth and roll corrections 2085 # to fix orientatins 2086 camera.Azimuth(-180) 2087 camera.Roll(180) 2088 2089 # shift the camera posiiton and focal position 2090 # to be centered on the desired location 2091 p = camera.GetPosition() 2092 p_new = np.array(p) + pos 2093 camera.SetPosition(*p_new) 2094 camera.SetFocalPoint(*pos) 2095 return camera 2096 2097 2098def camera_from_neuroglancer(state, zoom) -> vtki.vtkCamera: 2099 """ 2100 Define a `vtkCamera` from a neuroglancer state dictionary. 2101 2102 Arguments: 2103 state: (dict) 2104 an neuroglancer state dictionary. 2105 zoom: (float) 2106 how much to multiply zoom by to get camera backoff distance 2107 2108 Returns: 2109 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2110 """ 2111 orient = state.get("perspectiveOrientation", [0.0, 0.0, 0.0, 1.0]) 2112 pzoom = state.get("perspectiveZoom", 10.0) 2113 position = state["navigation"]["pose"]["position"] 2114 pos_nm = np.array(position["voxelCoordinates"]) * position["voxelSize"] 2115 return camera_from_quaternion(pos_nm, orient, pzoom * zoom, ngl_correct=True) 2116 2117 2118def oriented_camera(center, up_vector, backoff_vector, backoff) -> vtki.vtkCamera: 2119 """ 2120 Generate a `vtkCamera` pointed at a specific location, 2121 oriented with a given up direction, set to a backoff. 2122 """ 2123 vup = np.array(up_vector) 2124 vup = vup / np.linalg.norm(vup) 2125 pt_backoff = center - backoff * np.array(backoff_vector) 2126 camera = vtki.vtkCamera() 2127 camera.SetFocalPoint(center[0], center[1], center[2]) 2128 camera.SetViewUp(vup[0], vup[1], vup[2]) 2129 camera.SetPosition(pt_backoff[0], pt_backoff[1], pt_backoff[2]) 2130 return camera 2131 2132 2133def camera_from_dict(camera, modify_inplace=None) -> vtki.vtkCamera: 2134 """ 2135 Generate a `vtkCamera` object from a python dictionary. 2136 2137 Parameters of the camera are: 2138 - `position` or `pos` (3-tuple) 2139 - `focal_point` (3-tuple) 2140 - `viewup` (3-tuple) 2141 - `distance` (float) 2142 - `clipping_range` (2-tuple) 2143 - `parallel_scale` (float) 2144 - `thickness` (float) 2145 - `view_angle` (float) 2146 - `roll` (float) 2147 2148 Exaplanation of the parameters can be found in the 2149 [vtkCamera documentation](https://vtk.org/doc/nightly/html/classvtkCamera.html). 2150 2151 Arguments: 2152 camera: (dict) 2153 a python dictionary containing camera parameters. 2154 modify_inplace: (vtkCamera) 2155 an existing `vtkCamera` object to modify in place. 2156 2157 Returns: 2158 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2159 """ 2160 if modify_inplace: 2161 vcam = modify_inplace 2162 else: 2163 vcam = vtki.vtkCamera() 2164 2165 camera = dict(camera) # make a copy so input is not emptied by pop() 2166 2167 cm_pos = camera.pop("position", camera.pop("pos", None)) 2168 cm_focal_point = camera.pop("focal_point", camera.pop("focalPoint", None)) 2169 cm_viewup = camera.pop("viewup", None) 2170 cm_distance = camera.pop("distance", None) 2171 cm_clipping_range = camera.pop("clipping_range", camera.pop("clippingRange", None)) 2172 cm_parallel_scale = camera.pop("parallel_scale", camera.pop("parallelScale", None)) 2173 cm_thickness = camera.pop("thickness", None) 2174 cm_view_angle = camera.pop("view_angle", camera.pop("viewAngle", None)) 2175 cm_roll = camera.pop("roll", None) 2176 2177 if len(camera.keys()) > 0: 2178 vedo.logger.warning(f"in camera_from_dict, key(s) not recognized: {camera.keys()}") 2179 if cm_pos is not None: vcam.SetPosition(cm_pos) 2180 if cm_focal_point is not None: vcam.SetFocalPoint(cm_focal_point) 2181 if cm_viewup is not None: vcam.SetViewUp(cm_viewup) 2182 if cm_distance is not None: vcam.SetDistance(cm_distance) 2183 if cm_clipping_range is not None: vcam.SetClippingRange(cm_clipping_range) 2184 if cm_parallel_scale is not None: vcam.SetParallelScale(cm_parallel_scale) 2185 if cm_thickness is not None: vcam.SetThickness(cm_thickness) 2186 if cm_view_angle is not None: vcam.SetViewAngle(cm_view_angle) 2187 if cm_roll is not None: vcam.SetRoll(cm_roll) 2188 return vcam 2189 2190def camera_to_dict(vtkcam) -> dict: 2191 """ 2192 Convert a [vtkCamera](https://vtk.org/doc/nightly/html/classvtkCamera.html) 2193 object into a python dictionary. 2194 2195 Parameters of the camera are: 2196 - `position` (3-tuple) 2197 - `focal_point` (3-tuple) 2198 - `viewup` (3-tuple) 2199 - `distance` (float) 2200 - `clipping_range` (2-tuple) 2201 - `parallel_scale` (float) 2202 - `thickness` (float) 2203 - `view_angle` (float) 2204 - `roll` (float) 2205 2206 Arguments: 2207 vtkcam: (vtkCamera) 2208 a `vtkCamera` object to convert. 2209 """ 2210 cam = dict() 2211 cam["position"] = np.array(vtkcam.GetPosition()) 2212 cam["focal_point"] = np.array(vtkcam.GetFocalPoint()) 2213 cam["viewup"] = np.array(vtkcam.GetViewUp()) 2214 cam["distance"] = vtkcam.GetDistance() 2215 cam["clipping_range"] = np.array(vtkcam.GetClippingRange()) 2216 cam["parallel_scale"] = vtkcam.GetParallelScale() 2217 cam["thickness"] = vtkcam.GetThickness() 2218 cam["view_angle"] = vtkcam.GetViewAngle() 2219 cam["roll"] = vtkcam.GetRoll() 2220 return cam 2221 2222 2223def vtkCameraToK3D(vtkcam) -> np.ndarray: 2224 """ 2225 Convert a `vtkCamera` object into a 9-element list to be used by the K3D backend. 2226 2227 Output format is: 2228 `[posx,posy,posz, targetx,targety,targetz, upx,upy,upz]`. 2229 """ 2230 cpos = np.array(vtkcam.GetPosition()) 2231 kam = [cpos.tolist()] 2232 kam.append(vtkcam.GetFocalPoint()) 2233 kam.append(vtkcam.GetViewUp()) 2234 return np.array(kam).ravel() 2235 2236 2237def make_ticks( 2238 x0: float, 2239 x1: float, 2240 n=None, 2241 labels=None, 2242 digits=None, 2243 logscale=False, 2244 useformat="", 2245 ) -> Tuple[np.ndarray, List[str]]: 2246 """ 2247 Generate numeric labels for the `[x0, x1]` range. 2248 2249 The format specifier could be expressed in the format: 2250 `:[[fill]align][sign][#][0][width][,][.precision][type]` 2251 2252 where, the options are: 2253 ``` 2254 fill = any character 2255 align = < | > | = | ^ 2256 sign = + | - | " " 2257 width = integer 2258 precision = integer 2259 type = b | c | d | e | E | f | F | g | G | n | o | s | x | X | % 2260 ``` 2261 2262 E.g.: useformat=":.2f" 2263 """ 2264 # Copyright M. Musy, 2021, license: MIT. 2265 # 2266 # useformat eg: ":.2f", check out: 2267 # https://mkaz.blog/code/python-string-format-cookbook/ 2268 # https://www.programiz.com/python-programming/methods/built-in/format 2269 2270 if x1 <= x0: 2271 # vedo.printc("Error in make_ticks(): x0 >= x1", x0,x1, c='r') 2272 return np.array([0.0, 1.0]), ["", ""] 2273 2274 ticks_str, ticks_float = [], [] 2275 baseline = (1, 2, 5, 10, 20, 50) 2276 2277 if logscale: 2278 if x0 <= 0 or x1 <= 0: 2279 vedo.logger.error("make_ticks: zero or negative range with log scale.") 2280 raise RuntimeError 2281 if n is None: 2282 n = int(abs(np.log10(x1) - np.log10(x0))) + 1 2283 x0, x1 = np.log10([x0, x1]) 2284 2285 if not n: 2286 n = 5 2287 2288 if labels is not None: 2289 # user is passing custom labels 2290 2291 ticks_float.append(0) 2292 ticks_str.append("") 2293 for tp, ts in labels: 2294 if tp == x1: 2295 continue 2296 ticks_str.append(str(ts)) 2297 tickn = lin_interpolate(tp, [x0, x1], [0, 1]) 2298 ticks_float.append(tickn) 2299 2300 else: 2301 # ..here comes one of the shortest and most painful pieces of code: 2302 # automatically choose the best natural axis subdivision based on multiples of 1,2,5 2303 dstep = (x1 - x0) / n # desired step size, begin of the nightmare 2304 2305 basestep = pow(10, np.floor(np.log10(dstep))) 2306 steps = np.array([basestep * i for i in baseline]) 2307 idx = (np.abs(steps - dstep)).argmin() 2308 s = steps[idx] # chosen step size 2309 2310 low_bound, up_bound = 0, 0 2311 if x0 < 0: 2312 low_bound = -pow(10, np.ceil(np.log10(-x0))) 2313 if x1 > 0: 2314 up_bound = pow(10, np.ceil(np.log10(x1))) 2315 2316 if low_bound < 0: 2317 if up_bound < 0: 2318 negaxis = np.arange(low_bound, int(up_bound / s) * s) 2319 else: 2320 if -low_bound / s > 1.0e06: 2321 return np.array([0.0, 1.0]), ["", ""] 2322 negaxis = np.arange(low_bound, 0, s) 2323 else: 2324 negaxis = np.array([]) 2325 2326 if up_bound > 0: 2327 if low_bound > 0: 2328 posaxis = np.arange(int(low_bound / s) * s, up_bound, s) 2329 else: 2330 if up_bound / s > 1.0e06: 2331 return np.array([0.0, 1.0]), ["", ""] 2332 posaxis = np.arange(0, up_bound, s) 2333 else: 2334 posaxis = np.array([]) 2335 2336 fulaxis = np.unique(np.clip(np.concatenate([negaxis, posaxis]), x0, x1)) 2337 # end of the nightmare 2338 2339 if useformat: 2340 sf = "{" + f"{useformat}" + "}" 2341 sas = "" 2342 for x in fulaxis: 2343 sas += sf.format(x) + " " 2344 elif digits is None: 2345 np.set_printoptions(suppress=True) # avoid zero precision 2346 sas = str(fulaxis).replace("[", "").replace("]", "") 2347 sas = sas.replace(".e", "e").replace("e+0", "e+").replace("e-0", "e-") 2348 np.set_printoptions(suppress=None) # set back to default 2349 else: 2350 sas = precision(fulaxis, digits, vrange=(x0, x1)) 2351 sas = sas.replace("[", "").replace("]", "").replace(")", "").replace(",", "") 2352 2353 sas2 = [] 2354 for s in sas.split(): 2355 if s.endswith("."): 2356 s = s[:-1] 2357 if s == "-0": 2358 s = "0" 2359 if digits is not None and "e" in s: 2360 s += " " # add space to terminate modifiers 2361 sas2.append(s) 2362 2363 for ts, tp in zip(sas2, fulaxis): 2364 if tp == x1: 2365 continue 2366 tickn = lin_interpolate(tp, [x0, x1], [0, 1]) 2367 ticks_float.append(tickn) 2368 if logscale: 2369 val = np.power(10, tp) 2370 if useformat: 2371 sf = "{" + f"{useformat}" + "}" 2372 ticks_str.append(sf.format(val)) 2373 else: 2374 if val >= 10: 2375 val = int(val + 0.5) 2376 else: 2377 val = round_to_digit(val, 2) 2378 ticks_str.append(str(val)) 2379 else: 2380 ticks_str.append(ts) 2381 2382 ticks_str.append("") 2383 ticks_float.append(1) 2384 return np.array(ticks_float), ticks_str 2385 2386 2387def grid_corners(i: int, nm: list, size: list, margin=0, yflip=True) -> Tuple[np.ndarray, np.ndarray]: 2388 """ 2389 Compute the 2 corners coordinates of the i-th box in a grid of shape n*m. 2390 The top-left square is square number 1. 2391 2392 Arguments: 2393 i : (int) 2394 input index of the desired grid square (to be used in `show(..., at=...)`). 2395 nm : (list) 2396 grid shape as (n,m). 2397 size : (list) 2398 total size of the grid along x and y. 2399 margin : (float) 2400 keep a small margin between boxes. 2401 yflip : (bool) 2402 y-coordinate points downwards 2403 2404 Returns: 2405 Two 2D points representing the bottom-left corner and the top-right corner 2406 of the `i`-nth box in the grid. 2407 2408 Example: 2409 ```python 2410 from vedo import * 2411 acts=[] 2412 n,m = 5,7 2413 for i in range(1, n*m + 1): 2414 c1,c2 = utils.grid_corners(i, [n,m], [1,1], 0.01) 2415 t = Text3D(i, (c1+c2)/2, c='k', s=0.02, justify='center').z(0.01) 2416 r = Rectangle(c1, c2, c=i) 2417 acts += [t,r] 2418 show(acts, axes=1).close() 2419 ``` 2420  2421 """ 2422 i -= 1 2423 n, m = nm 2424 sx, sy = size 2425 dx, dy = sx / n, sy / m 2426 nx = i % n 2427 ny = int((i - nx) / n) 2428 if yflip: 2429 ny = n - ny 2430 c1 = (dx * nx + margin, dy * ny + margin) 2431 c2 = (dx * (nx + 1) - margin, dy * (ny + 1) - margin) 2432 return np.array(c1), np.array(c2) 2433 2434 2435############################################################################ 2436# Trimesh support 2437# 2438# Install trimesh with: 2439# 2440# sudo apt install python3-rtree 2441# pip install rtree shapely 2442# conda install trimesh 2443# 2444# Check the example gallery in: examples/other/trimesh> 2445########################################################################### 2446def vedo2trimesh(mesh): 2447 """ 2448 Convert `vedo.mesh.Mesh` to `Trimesh.Mesh` object. 2449 """ 2450 if is_sequence(mesh): 2451 tms = [] 2452 for a in mesh: 2453 tms.append(vedo2trimesh(a)) 2454 return tms 2455 2456 try: 2457 from trimesh import Trimesh 2458 except ModuleNotFoundError: 2459 vedo.logger.error("Need trimesh to run:\npip install trimesh") 2460 return None 2461 2462 tris = mesh.cells 2463 carr = mesh.celldata["CellIndividualColors"] 2464 ccols = carr 2465 2466 points = mesh.vertices 2467 varr = mesh.pointdata["VertexColors"] 2468 vcols = varr 2469 2470 if len(tris) == 0: 2471 tris = None 2472 2473 return Trimesh(vertices=points, faces=tris, face_colors=ccols, vertex_colors=vcols, process=False) 2474 2475 2476def trimesh2vedo(inputobj): 2477 """ 2478 Convert a `Trimesh` object to `vedo.Mesh` or `vedo.Assembly` object. 2479 """ 2480 if is_sequence(inputobj): 2481 vms = [] 2482 for ob in inputobj: 2483 vms.append(trimesh2vedo(ob)) 2484 return vms 2485 2486 inputobj_type = str(type(inputobj)) 2487 2488 if "Trimesh" in inputobj_type or "primitives" in inputobj_type: 2489 faces = inputobj.faces 2490 poly = buildPolyData(inputobj.vertices, faces) 2491 tact = vedo.Mesh(poly) 2492 if inputobj.visual.kind == "face": 2493 trim_c = inputobj.visual.face_colors 2494 elif inputobj.visual.kind == "texture": 2495 trim_c = inputobj.visual.to_color().vertex_colors 2496 else: 2497 trim_c = inputobj.visual.vertex_colors 2498 2499 if is_sequence(trim_c): 2500 if is_sequence(trim_c[0]): 2501 same_color = len(np.unique(trim_c, axis=0)) < 2 # all vtxs have same color 2502 2503 if same_color: 2504 tact.c(trim_c[0, [0, 1, 2]]).alpha(trim_c[0, 3]) 2505 else: 2506 if inputobj.visual.kind == "face": 2507 tact.cellcolors = trim_c 2508 return tact 2509 2510 if "PointCloud" in inputobj_type: 2511 2512 vdpts = vedo.shapes.Points(inputobj.vertices, r=8, c='k') 2513 if hasattr(inputobj, "vertices_color"): 2514 vcols = (inputobj.vertices_color * 1).astype(np.uint8) 2515 vdpts.pointcolors = vcols 2516 return vdpts 2517 2518 if "path" in inputobj_type: 2519 2520 lines = [] 2521 for e in inputobj.entities: 2522 # print('trimesh entity', e.to_dict()) 2523 l = vedo.shapes.Line(inputobj.vertices[e.points], c="k", lw=2) 2524 lines.append(l) 2525 return vedo.Assembly(lines) 2526 2527 return None 2528 2529 2530def vedo2meshlab(vmesh): 2531 """Convert a `vedo.Mesh` to a Meshlab object.""" 2532 try: 2533 import pymeshlab as mlab 2534 except ModuleNotFoundError: 2535 vedo.logger.error("Need pymeshlab to run:\npip install pymeshlab") 2536 2537 vertex_matrix = vmesh.vertices.astype(np.float64) 2538 2539 try: 2540 face_matrix = np.asarray(vmesh.cells, dtype=np.float64) 2541 except: 2542 print("WARNING: in vedo2meshlab(), need to triangulate mesh first!") 2543 face_matrix = np.array(vmesh.clone().triangulate().cells, dtype=np.float64) 2544 2545 # v_normals_matrix = vmesh.normals(cells=False, recompute=False) 2546 v_normals_matrix = vmesh.vertex_normals 2547 if not v_normals_matrix.shape[0]: 2548 v_normals_matrix = np.empty((0, 3), dtype=np.float64) 2549 2550 # f_normals_matrix = vmesh.normals(cells=True, recompute=False) 2551 f_normals_matrix = vmesh.cell_normals 2552 if not f_normals_matrix.shape[0]: 2553 f_normals_matrix = np.empty((0, 3), dtype=np.float64) 2554 2555 v_color_matrix = vmesh.pointdata["RGBA"] 2556 if v_color_matrix is None: 2557 v_color_matrix = np.empty((0, 4), dtype=np.float64) 2558 else: 2559 v_color_matrix = v_color_matrix.astype(np.float64) / 255 2560 if v_color_matrix.shape[1] == 3: 2561 v_color_matrix = np.c_[ 2562 v_color_matrix, np.ones(v_color_matrix.shape[0], dtype=np.float64) 2563 ] 2564 2565 f_color_matrix = vmesh.celldata["RGBA"] 2566 if f_color_matrix is None: 2567 f_color_matrix = np.empty((0, 4), dtype=np.float64) 2568 else: 2569 f_color_matrix = f_color_matrix.astype(np.float64) / 255 2570 if f_color_matrix.shape[1] == 3: 2571 f_color_matrix = np.c_[ 2572 f_color_matrix, np.ones(f_color_matrix.shape[0], dtype=np.float64) 2573 ] 2574 2575 m = mlab.Mesh( 2576 vertex_matrix=vertex_matrix, 2577 face_matrix=face_matrix, 2578 v_normals_matrix=v_normals_matrix, 2579 f_normals_matrix=f_normals_matrix, 2580 v_color_matrix=v_color_matrix, 2581 f_color_matrix=f_color_matrix, 2582 ) 2583 2584 for k in vmesh.pointdata.keys(): 2585 data = vmesh.pointdata[k] 2586 if data is not None: 2587 if data.ndim == 1: # scalar 2588 m.add_vertex_custom_scalar_attribute(data.astype(np.float64), k) 2589 elif data.ndim == 2: # vectorial data 2590 if "tcoord" not in k.lower() and k not in ["Normals", "TextureCoordinates"]: 2591 m.add_vertex_custom_point_attribute(data.astype(np.float64), k) 2592 2593 for k in vmesh.celldata.keys(): 2594 data = vmesh.celldata[k] 2595 if data is not None: 2596 if data.ndim == 1: # scalar 2597 m.add_face_custom_scalar_attribute(data.astype(np.float64), k) 2598 elif data.ndim == 2 and k != "Normals": # vectorial data 2599 m.add_face_custom_point_attribute(data.astype(np.float64), k) 2600 2601 m.update_bounding_box() 2602 return m 2603 2604 2605def meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()): 2606 """Convert a Meshlab object to `vedo.Mesh`.""" 2607 inputtype = str(type(mmesh)) 2608 2609 if "MeshSet" in inputtype: 2610 mmesh = mmesh.current_mesh() 2611 2612 mpoints, mcells = mmesh.vertex_matrix(), mmesh.face_matrix() 2613 if len(mcells) > 0: 2614 polydata = buildPolyData(mpoints, mcells) 2615 else: 2616 polydata = buildPolyData(mpoints, None) 2617 2618 if mmesh.has_vertex_scalar(): 2619 parr = mmesh.vertex_scalar_array() 2620 parr_vtk = numpy_to_vtk(parr) 2621 parr_vtk.SetName("MeshLabScalars") 2622 polydata.GetPointData().AddArray(parr_vtk) 2623 polydata.GetPointData().SetActiveScalars("MeshLabScalars") 2624 2625 if mmesh.has_face_scalar(): 2626 carr = mmesh.face_scalar_array() 2627 carr_vtk = numpy_to_vtk(carr) 2628 carr_vtk.SetName("MeshLabScalars") 2629 x0, x1 = carr_vtk.GetRange() 2630 polydata.GetCellData().AddArray(carr_vtk) 2631 polydata.GetCellData().SetActiveScalars("MeshLabScalars") 2632 2633 for k in pointdata_keys: 2634 parr = mmesh.vertex_custom_scalar_attribute_array(k) 2635 parr_vtk = numpy_to_vtk(parr) 2636 parr_vtk.SetName(k) 2637 polydata.GetPointData().AddArray(parr_vtk) 2638 polydata.GetPointData().SetActiveScalars(k) 2639 2640 for k in celldata_keys: 2641 carr = mmesh.face_custom_scalar_attribute_array(k) 2642 carr_vtk = numpy_to_vtk(carr) 2643 carr_vtk.SetName(k) 2644 polydata.GetCellData().AddArray(carr_vtk) 2645 polydata.GetCellData().SetActiveScalars(k) 2646 2647 pnorms = mmesh.vertex_normal_matrix() 2648 if len(pnorms) > 0: 2649 polydata.GetPointData().SetNormals(numpy2vtk(pnorms, name="Normals")) 2650 2651 cnorms = mmesh.face_normal_matrix() 2652 if len(cnorms) > 0: 2653 polydata.GetCellData().SetNormals(numpy2vtk(cnorms, name="Normals")) 2654 return vedo.Mesh(polydata) 2655 2656 2657def open3d2vedo(o3d_mesh): 2658 """Convert `open3d.geometry.TriangleMesh` to a `vedo.Mesh`.""" 2659 m = vedo.Mesh([np.array(o3d_mesh.vertices), np.array(o3d_mesh.triangles)]) 2660 # TODO: could also check whether normals and color are present in 2661 # order to port with the above vertices/faces 2662 return m 2663 2664 2665def vedo2open3d(vedo_mesh): 2666 """ 2667 Return an `open3d.geometry.TriangleMesh` version of the current mesh. 2668 """ 2669 try: 2670 import open3d as o3d 2671 except RuntimeError: 2672 vedo.logger.error("Need open3d to run:\npip install open3d") 2673 2674 # create from numpy arrays 2675 o3d_mesh = o3d.geometry.TriangleMesh( 2676 vertices=o3d.utility.Vector3dVector(vedo_mesh.vertices), 2677 triangles=o3d.utility.Vector3iVector(vedo_mesh.cells), 2678 ) 2679 # TODO: need to add some if check here in case color and normals 2680 # info are not existing 2681 # o3d_mesh.vertex_colors = o3d.utility.Vector3dVector(vedo_mesh.pointdata["RGB"]/255) 2682 # o3d_mesh.vertex_normals= o3d.utility.Vector3dVector(vedo_mesh.pointdata["Normals"]) 2683 return o3d_mesh 2684 2685def vedo2madcad(vedo_mesh): 2686 """ 2687 Convert a `vedo.Mesh` to a `madcad.Mesh`. 2688 """ 2689 try: 2690 import madcad 2691 import numbers 2692 except ModuleNotFoundError: 2693 vedo.logger.error("Need madcad to run:\npip install pymadcad") 2694 2695 points = [madcad.vec3(*pt) for pt in vedo_mesh.vertices] 2696 faces = [madcad.vec3(*fc) for fc in vedo_mesh.cells] 2697 2698 options = {} 2699 for key, val in vedo_mesh.pointdata.items(): 2700 vec_type = f"vec{val.shape[-1]}" 2701 is_float = np.issubdtype(val.dtype, np.floating) 2702 madcad_dtype = getattr(madcad, f"f{vec_type}" if is_float else vec_type) 2703 options[key] = [madcad_dtype(v) for v in val] 2704 2705 madcad_mesh = madcad.Mesh(points=points, faces=faces, options=options) 2706 2707 return madcad_mesh 2708 2709 2710def madcad2vedo(madcad_mesh): 2711 """ 2712 Convert a `madcad.Mesh` to a `vedo.Mesh`. 2713 2714 A pointdata or celldata array named "tracks" is added to the output mesh, indicating 2715 the mesh region each point belongs to. 2716 2717 A metadata array named "madcad_groups" is added to the output mesh, indicating 2718 the mesh groups. 2719 2720 See [pymadcad website](https://pymadcad.readthedocs.io/en/latest/index.html) 2721 for more info. 2722 """ 2723 try: 2724 madcad_mesh = madcad_mesh["part"] 2725 except: 2726 pass 2727 2728 madp = [] 2729 for p in madcad_mesh.points: 2730 madp.append([float(p[0]), float(p[1]), float(p[2])]) 2731 madp = np.array(madp) 2732 2733 madf = [] 2734 try: 2735 for f in madcad_mesh.faces: 2736 madf.append([int(f[0]), int(f[1]), int(f[2])]) 2737 madf = np.array(madf).astype(np.uint16) 2738 except AttributeError: 2739 # print("no faces") 2740 pass 2741 2742 made = [] 2743 try: 2744 edges = madcad_mesh.edges 2745 for e in edges: 2746 made.append([int(e[0]), int(e[1])]) 2747 made = np.array(made).astype(np.uint16) 2748 except (AttributeError, TypeError): 2749 # print("no edges") 2750 pass 2751 2752 try: 2753 line = np.array(madcad_mesh.indices).astype(np.uint16) 2754 made.append(line) 2755 except AttributeError: 2756 # print("no indices") 2757 pass 2758 2759 madt = [] 2760 try: 2761 for t in madcad_mesh.tracks: 2762 madt.append(int(t)) 2763 madt = np.array(madt).astype(np.uint16) 2764 except AttributeError: 2765 # print("no tracks") 2766 pass 2767 2768 ############################### 2769 poly = vedo.utils.buildPolyData(madp, madf, made) 2770 if len(madf) == 0 and len(made) == 0: 2771 m = vedo.Points(poly) 2772 else: 2773 m = vedo.Mesh(poly) 2774 2775 if len(madt) == len(madf): 2776 m.celldata["tracks"] = madt 2777 maxt = np.max(madt) 2778 m.mapper.SetScalarRange(0, np.max(madt)) 2779 if maxt==0: m.mapper.SetScalarVisibility(0) 2780 elif len(madt) == len(madp): 2781 m.pointdata["tracks"] = madt 2782 maxt = np.max(madt) 2783 m.mapper.SetScalarRange(0, maxt) 2784 if maxt==0: m.mapper.SetScalarVisibility(0) 2785 2786 try: 2787 m.info["madcad_groups"] = madcad_mesh.groups 2788 except AttributeError: 2789 # print("no groups") 2790 pass 2791 2792 try: 2793 options = dict(madcad_mesh.options) 2794 if "display_wire" in options and options["display_wire"]: 2795 m.lw(1).lc(madcad_mesh.c()) 2796 if "display_faces" in options and not options["display_faces"]: 2797 m.alpha(0.2) 2798 if "color" in options: 2799 m.c(options["color"]) 2800 2801 for key, val in options.items(): 2802 m.pointdata[key] = val 2803 2804 except AttributeError: 2805 # print("no options") 2806 pass 2807 2808 return m 2809 2810 2811def vtk_version_at_least(major, minor=0, build=0) -> bool: 2812 """ 2813 Check the installed VTK version. 2814 2815 Return `True` if the requested VTK version is greater or equal to the actual VTK version. 2816 2817 Arguments: 2818 major : (int) 2819 Major version. 2820 minor : (int) 2821 Minor version. 2822 build : (int) 2823 Build version. 2824 """ 2825 needed_version = 10000000000 * int(major) + 100000000 * int(minor) + int(build) 2826 try: 2827 vtk_version_number = vtki.VTK_VERSION_NUMBER 2828 except AttributeError: # as error: 2829 ver = vtki.vtkVersion() 2830 vtk_version_number = ( 2831 10000000000 * ver.GetVTKMajorVersion() 2832 + 100000000 * ver.GetVTKMinorVersion() 2833 + ver.GetVTKBuildVersion() 2834 ) 2835 return vtk_version_number >= needed_version 2836 2837 2838def ctf2lut(vol, logscale=False): 2839 """Internal use.""" 2840 # build LUT from a color transfer function for tmesh or volume 2841 2842 ctf = vol.properties.GetRGBTransferFunction() 2843 otf = vol.properties.GetScalarOpacity() 2844 x0, x1 = vol.dataset.GetScalarRange() 2845 cols, alphas = [], [] 2846 for x in np.linspace(x0, x1, 256): 2847 cols.append(ctf.GetColor(x)) 2848 alphas.append(otf.GetValue(x)) 2849 2850 if logscale: 2851 lut = vtki.vtkLogLookupTable() 2852 else: 2853 lut = vtki.vtkLookupTable() 2854 2855 lut.SetRange(x0, x1) 2856 lut.SetNumberOfTableValues(len(cols)) 2857 for i, col in enumerate(cols): 2858 r, g, b = col 2859 lut.SetTableValue(i, r, g, b, alphas[i]) 2860 lut.Build() 2861 return lut 2862 2863 2864def get_vtk_name_event(name: str) -> str: 2865 """ 2866 Return the name of a VTK event. 2867 2868 Frequently used events are: 2869 - KeyPress, KeyRelease: listen to keyboard events 2870 - LeftButtonPress, LeftButtonRelease: listen to mouse clicks 2871 - MiddleButtonPress, MiddleButtonRelease 2872 - RightButtonPress, RightButtonRelease 2873 - MouseMove: listen to mouse pointer changing position 2874 - MouseWheelForward, MouseWheelBackward 2875 - Enter, Leave: listen to mouse entering or leaving the window 2876 - Pick, StartPick, EndPick: listen to object picking 2877 - ResetCamera, ResetCameraClippingRange 2878 - Error, Warning, Char, Timer 2879 2880 Check the complete list of events here: 2881 https://vtk.org/doc/nightly/html/classvtkCommand.html 2882 """ 2883 # as vtk names are ugly and difficult to remember: 2884 ln = name.lower() 2885 if "click" in ln or "button" in ln: 2886 event_name = "LeftButtonPress" 2887 if "right" in ln: 2888 event_name = "RightButtonPress" 2889 elif "mid" in ln: 2890 event_name = "MiddleButtonPress" 2891 if "release" in ln: 2892 event_name = event_name.replace("Press", "Release") 2893 else: 2894 event_name = name 2895 if "key" in ln: 2896 if "release" in ln: 2897 event_name = "KeyRelease" 2898 else: 2899 event_name = "KeyPress" 2900 2901 if ("mouse" in ln and "mov" in ln) or "over" in ln: 2902 event_name = "MouseMove" 2903 2904 words = [ 2905 "pick", "timer", "reset", "enter", "leave", "char", 2906 "error", "warning", "start", "end", "wheel", "clipping", 2907 "range", "camera", "render", "interaction", "modified", 2908 ] 2909 for w in words: 2910 if w in ln: 2911 event_name = event_name.replace(w, w.capitalize()) 2912 2913 event_name = event_name.replace("REnd ", "Rend") 2914 event_name = event_name.replace("the ", "") 2915 event_name = event_name.replace(" of ", "").replace(" ", "") 2916 2917 if not event_name.endswith("Event"): 2918 event_name += "Event" 2919 2920 if vtki.vtkCommand.GetEventIdFromString(event_name) == 0: 2921 vedo.printc( 2922 f"Error: '{name}' is not a valid event name.", c='r') 2923 vedo.printc("Check the list of events here:", c='r') 2924 vedo.printc("\thttps://vtk.org/doc/nightly/html/classvtkCommand.html", c='r') 2925 # raise RuntimeError 2926 2927 return event_name
class
OperationNode:
61class OperationNode: 62 """ 63 Keep track of the operations which led to a final state. 64 """ 65 # https://www.graphviz.org/doc/info/shapes.html#html 66 # Mesh #e9c46a 67 # Follower #d9ed92 68 # Volume, UnstructuredGrid #4cc9f0 69 # TetMesh #9e2a2b 70 # File #8a817c 71 # Image #f28482 72 # Assembly #f08080 73 74 def __init__( 75 self, operation, parents=(), comment="", shape="none", c="#e9c46a", style="filled" 76 ) -> None: 77 """ 78 Keep track of the operations which led to a final object. 79 This allows to show the `pipeline` tree for any `vedo` object with e.g.: 80 81 ```python 82 from vedo import * 83 sp = Sphere() 84 sp.clean().subdivide() 85 sp.pipeline.show() 86 ``` 87 88 Arguments: 89 operation : (str, class) 90 descriptor label, if a class is passed then grab its name 91 parents : (list) 92 list of the parent classes the object comes from 93 comment : (str) 94 a second-line text description 95 shape : (str) 96 shape of the frame, check out [this link.](https://graphviz.org/doc/info/shapes.html) 97 c : (hex) 98 hex color 99 style : (str) 100 comma-separated list of styles 101 102 Example: 103 ```python 104 from vedo.utils import OperationNode 105 106 op_node1 = OperationNode("Operation1", c="lightblue") 107 op_node2 = OperationNode("Operation2") 108 op_node3 = OperationNode("Operation3", shape='diamond') 109 op_node4 = OperationNode("Operation4") 110 op_node5 = OperationNode("Operation5") 111 op_node6 = OperationNode("Result", c="lightgreen") 112 113 op_node3.add_parent(op_node1) 114 op_node4.add_parent(op_node1) 115 op_node3.add_parent(op_node2) 116 op_node5.add_parent(op_node2) 117 op_node6.add_parent(op_node3) 118 op_node6.add_parent(op_node5) 119 op_node6.add_parent(op_node1) 120 121 op_node6.show(orientation="TB") 122 ``` 123  124 """ 125 if not vedo.settings.enable_pipeline: 126 return 127 128 if isinstance(operation, str): 129 self.operation = operation 130 else: 131 self.operation = operation.__class__.__name__ 132 self.operation_plain = str(self.operation) 133 134 pp = [] # filter out invalid stuff 135 for p in parents: 136 if hasattr(p, "pipeline"): 137 pp.append(p.pipeline) 138 self.parents = pp 139 140 if comment: 141 self.operation = f"<{self.operation}<BR/><SUB><I>{comment}</I></SUB>>" 142 143 self.dot = None 144 self.time = time.time() 145 self.shape = shape 146 self.style = style 147 self.color = c 148 self.counts = 0 149 150 def add_parent(self, parent) -> None: 151 self.parents.append(parent) 152 153 def _build_tree(self, dot): 154 dot.node( 155 str(id(self)), 156 label=self.operation, 157 shape=self.shape, 158 color=self.color, 159 style=self.style, 160 ) 161 for parent in self.parents: 162 if parent: 163 t = f"{self.time - parent.time: .1f}s" 164 dot.edge(str(id(parent)), str(id(self)), label=t) 165 parent._build_tree(dot) 166 167 def __repr__(self): 168 try: 169 from treelib import Tree 170 except ImportError: 171 vedo.logger.error( 172 "To use this functionality please install treelib:" 173 "\n pip install treelib" 174 ) 175 return "" 176 177 def _build_tree(parent): 178 for par in parent.parents: 179 if par: 180 op = par.operation_plain 181 tree.create_node( 182 op, op + str(par.time), parent=parent.operation_plain + str(parent.time) 183 ) 184 _build_tree(par) 185 try: 186 tree = Tree() 187 tree.create_node(self.operation_plain, self.operation_plain + str(self.time)) 188 _build_tree(self) 189 out = tree.show(stdout=False) 190 except: 191 out = f"Sorry treelib failed to build the tree for '{self.operation_plain}()'." 192 return out 193 194 def print(self) -> None: 195 """Print the tree of operations.""" 196 print(self.__str__()) 197 198 def show(self, orientation="LR", popup=True) -> None: 199 """Show the graphviz output for the pipeline of this object""" 200 if not vedo.settings.enable_pipeline: 201 return 202 203 try: 204 from graphviz import Digraph 205 except ImportError: 206 vedo.logger.error("please install graphviz with command\n pip install graphviz") 207 return 208 209 # visualize the entire tree 210 dot = Digraph( 211 node_attr={"fontcolor": "#201010", "fontname": "Helvetica", "fontsize": "12"}, 212 edge_attr={"fontname": "Helvetica", "fontsize": "6", "arrowsize": "0.4"}, 213 ) 214 dot.attr(rankdir=orientation) 215 216 self.counts = 0 217 self._build_tree(dot) 218 self.dot = dot 219 220 home_dir = os.path.expanduser("~") 221 gpath = os.path.join( 222 home_dir, vedo.settings.cache_directory, "vedo", "pipeline_graphviz") 223 224 dot.render(gpath, view=popup)
Keep track of the operations which led to a final state.
OperationNode( operation, parents=(), comment='', shape='none', c='#e9c46a', style='filled')
74 def __init__( 75 self, operation, parents=(), comment="", shape="none", c="#e9c46a", style="filled" 76 ) -> None: 77 """ 78 Keep track of the operations which led to a final object. 79 This allows to show the `pipeline` tree for any `vedo` object with e.g.: 80 81 ```python 82 from vedo import * 83 sp = Sphere() 84 sp.clean().subdivide() 85 sp.pipeline.show() 86 ``` 87 88 Arguments: 89 operation : (str, class) 90 descriptor label, if a class is passed then grab its name 91 parents : (list) 92 list of the parent classes the object comes from 93 comment : (str) 94 a second-line text description 95 shape : (str) 96 shape of the frame, check out [this link.](https://graphviz.org/doc/info/shapes.html) 97 c : (hex) 98 hex color 99 style : (str) 100 comma-separated list of styles 101 102 Example: 103 ```python 104 from vedo.utils import OperationNode 105 106 op_node1 = OperationNode("Operation1", c="lightblue") 107 op_node2 = OperationNode("Operation2") 108 op_node3 = OperationNode("Operation3", shape='diamond') 109 op_node4 = OperationNode("Operation4") 110 op_node5 = OperationNode("Operation5") 111 op_node6 = OperationNode("Result", c="lightgreen") 112 113 op_node3.add_parent(op_node1) 114 op_node4.add_parent(op_node1) 115 op_node3.add_parent(op_node2) 116 op_node5.add_parent(op_node2) 117 op_node6.add_parent(op_node3) 118 op_node6.add_parent(op_node5) 119 op_node6.add_parent(op_node1) 120 121 op_node6.show(orientation="TB") 122 ``` 123  124 """ 125 if not vedo.settings.enable_pipeline: 126 return 127 128 if isinstance(operation, str): 129 self.operation = operation 130 else: 131 self.operation = operation.__class__.__name__ 132 self.operation_plain = str(self.operation) 133 134 pp = [] # filter out invalid stuff 135 for p in parents: 136 if hasattr(p, "pipeline"): 137 pp.append(p.pipeline) 138 self.parents = pp 139 140 if comment: 141 self.operation = f"<{self.operation}<BR/><SUB><I>{comment}</I></SUB>>" 142 143 self.dot = None 144 self.time = time.time() 145 self.shape = shape 146 self.style = style 147 self.color = c 148 self.counts = 0
Keep track of the operations which led to a final object.
This allows to show the pipeline tree for any vedo object with e.g.:
from vedo import *
sp = Sphere()
sp.clean().subdivide()
sp.pipeline.show()
Arguments:
- operation : (str, class) descriptor label, if a class is passed then grab its name
- parents : (list) list of the parent classes the object comes from
- comment : (str) a second-line text description
- shape : (str) shape of the frame, check out this link.
- c : (hex) hex color
- style : (str) comma-separated list of styles
Example:
from vedo.utils import OperationNode op_node1 = OperationNode("Operation1", c="lightblue") op_node2 = OperationNode("Operation2") op_node3 = OperationNode("Operation3", shape='diamond') op_node4 = OperationNode("Operation4") op_node5 = OperationNode("Operation5") op_node6 = OperationNode("Result", c="lightgreen") op_node3.add_parent(op_node1) op_node4.add_parent(op_node1) op_node3.add_parent(op_node2) op_node5.add_parent(op_node2) op_node6.add_parent(op_node3) op_node6.add_parent(op_node5) op_node6.add_parent(op_node1) op_node6.show(orientation="TB")
def
show(self, orientation='LR', popup=True) -> None:
198 def show(self, orientation="LR", popup=True) -> None: 199 """Show the graphviz output for the pipeline of this object""" 200 if not vedo.settings.enable_pipeline: 201 return 202 203 try: 204 from graphviz import Digraph 205 except ImportError: 206 vedo.logger.error("please install graphviz with command\n pip install graphviz") 207 return 208 209 # visualize the entire tree 210 dot = Digraph( 211 node_attr={"fontcolor": "#201010", "fontname": "Helvetica", "fontsize": "12"}, 212 edge_attr={"fontname": "Helvetica", "fontsize": "6", "arrowsize": "0.4"}, 213 ) 214 dot.attr(rankdir=orientation) 215 216 self.counts = 0 217 self._build_tree(dot) 218 self.dot = dot 219 220 home_dir = os.path.expanduser("~") 221 gpath = os.path.join( 222 home_dir, vedo.settings.cache_directory, "vedo", "pipeline_graphviz") 223 224 dot.render(gpath, view=popup)
Show the graphviz output for the pipeline of this object
class
ProgressBar:
228class ProgressBar: 229 """ 230 Class to print a progress bar. 231 """ 232 233 def __init__( 234 self, 235 start, 236 stop, 237 step=1, 238 c=None, 239 bold=True, 240 italic=False, 241 title="", 242 eta=True, 243 delay=-1, 244 width=25, 245 char="\U00002501", 246 char_back="\U00002500", 247 ) -> None: 248 """ 249 Class to print a progress bar with optional text message. 250 251 Check out also function `progressbar()`. 252 253 Arguments: 254 start : (int) 255 starting value 256 stop : (int) 257 stopping value 258 step : (int) 259 step value 260 c : (str) 261 color in hex format 262 title : (str) 263 title text 264 eta : (bool) 265 estimate time of arrival 266 delay : (float) 267 minimum time before printing anything, 268 if negative use the default value 269 as set in `vedo.settings.progressbar_delay` 270 width : (int) 271 width of the progress bar 272 char : (str) 273 character to use for the progress bar 274 char_back : (str) 275 character to use for the background of the progress bar 276 277 Example: 278 ```python 279 import time 280 from vedo import ProgressBar 281 pb = ProgressBar(0,40, c='r') 282 for i in pb.range(): 283 time.sleep(0.1) 284 pb.print() 285 ``` 286  287 """ 288 self.char = char 289 self.char_back = char_back 290 291 self.title = title + " " 292 if title: 293 self.title = " " + self.title 294 295 if delay < 0: 296 delay = vedo.settings.progressbar_delay 297 298 self.start = start 299 self.stop = stop 300 self.step = step 301 302 self.color = c 303 self.bold = bold 304 self.italic = italic 305 self.width = width 306 self.pbar = "" 307 self.percent = 0.0 308 self.percent_int = 0 309 self.eta = eta 310 self.delay = delay 311 312 self.t0 = time.time() 313 self._remaining = 1e10 314 315 self._update(0) 316 317 self._counts = 0 318 self._oldbar = "" 319 self._lentxt = 0 320 self._range = np.arange(start, stop, step) 321 322 def print(self, txt="", c=None) -> None: 323 """Print the progress bar with an optional message.""" 324 if not c: 325 c = self.color 326 327 self._update(self._counts + self.step) 328 329 if self.delay: 330 if time.time() - self.t0 < self.delay: 331 return 332 333 if self.pbar != self._oldbar: 334 self._oldbar = self.pbar 335 336 if self.eta and self._counts > 1: 337 338 tdenom = time.time() - self.t0 339 if tdenom: 340 vel = self._counts / tdenom 341 self._remaining = (self.stop - self._counts) / vel 342 else: 343 vel = 1 344 self._remaining = 0.0 345 346 if self._remaining > 60: 347 _mins = int(self._remaining / 60) 348 _secs = self._remaining - 60 * _mins 349 mins = f"{_mins}m" 350 secs = f"{int(_secs + 0.5)}s " 351 else: 352 mins = "" 353 secs = f"{int(self._remaining + 0.5)}s " 354 355 vel = round(vel, 1) 356 eta = f"eta: {mins}{secs}({vel} it/s) " 357 if self._remaining < 0.5: 358 dt = time.time() - self.t0 359 if dt > 60: 360 _mins = int(dt / 60) 361 _secs = dt - 60 * _mins 362 mins = f"{_mins}m" 363 secs = f"{int(_secs + 0.5)}s " 364 else: 365 mins = "" 366 secs = f"{int(dt + 0.5)}s " 367 eta = f"elapsed: {mins}{secs}({vel} it/s) " 368 txt = "" 369 else: 370 eta = "" 371 372 eraser = " " * self._lentxt + "\b" * self._lentxt 373 374 s = f"{self.pbar} {eraser}{eta}{txt}\r" 375 vedo.printc(s, c=c, bold=self.bold, italic=self.italic, end="") 376 if self.percent > 99.999: 377 print("") 378 379 self._lentxt = len(txt) 380 381 def range(self) -> np.ndarray: 382 """Return the range iterator.""" 383 return self._range 384 385 def _update(self, counts): 386 if counts < self.start: 387 counts = self.start 388 elif counts > self.stop: 389 counts = self.stop 390 self._counts = counts 391 392 self.percent = (self._counts - self.start) * 100.0 393 394 delta = self.stop - self.start 395 if delta: 396 self.percent /= delta 397 else: 398 self.percent = 0.0 399 400 self.percent_int = int(round(self.percent)) 401 af = self.width - 2 402 nh = int(round(self.percent_int / 100 * af)) 403 pbar_background = "\x1b[2m" + self.char_back * (af - nh) 404 self.pbar = f"{self.title}{self.char * (nh-1)}{pbar_background}" 405 if self.percent < 100.0: 406 ps = f" {self.percent_int}%" 407 else: 408 ps = "" 409 self.pbar += ps
Class to print a progress bar.
ProgressBar( start, stop, step=1, c=None, bold=True, italic=False, title='', eta=True, delay=-1, width=25, char='━', char_back='─')
233 def __init__( 234 self, 235 start, 236 stop, 237 step=1, 238 c=None, 239 bold=True, 240 italic=False, 241 title="", 242 eta=True, 243 delay=-1, 244 width=25, 245 char="\U00002501", 246 char_back="\U00002500", 247 ) -> None: 248 """ 249 Class to print a progress bar with optional text message. 250 251 Check out also function `progressbar()`. 252 253 Arguments: 254 start : (int) 255 starting value 256 stop : (int) 257 stopping value 258 step : (int) 259 step value 260 c : (str) 261 color in hex format 262 title : (str) 263 title text 264 eta : (bool) 265 estimate time of arrival 266 delay : (float) 267 minimum time before printing anything, 268 if negative use the default value 269 as set in `vedo.settings.progressbar_delay` 270 width : (int) 271 width of the progress bar 272 char : (str) 273 character to use for the progress bar 274 char_back : (str) 275 character to use for the background of the progress bar 276 277 Example: 278 ```python 279 import time 280 from vedo import ProgressBar 281 pb = ProgressBar(0,40, c='r') 282 for i in pb.range(): 283 time.sleep(0.1) 284 pb.print() 285 ``` 286  287 """ 288 self.char = char 289 self.char_back = char_back 290 291 self.title = title + " " 292 if title: 293 self.title = " " + self.title 294 295 if delay < 0: 296 delay = vedo.settings.progressbar_delay 297 298 self.start = start 299 self.stop = stop 300 self.step = step 301 302 self.color = c 303 self.bold = bold 304 self.italic = italic 305 self.width = width 306 self.pbar = "" 307 self.percent = 0.0 308 self.percent_int = 0 309 self.eta = eta 310 self.delay = delay 311 312 self.t0 = time.time() 313 self._remaining = 1e10 314 315 self._update(0) 316 317 self._counts = 0 318 self._oldbar = "" 319 self._lentxt = 0 320 self._range = np.arange(start, stop, step)
Class to print a progress bar with optional text message.
Check out also function progressbar().
Arguments:
- start : (int) starting value
- stop : (int) stopping value
- step : (int) step value
- c : (str) color in hex format
- title : (str) title text
- eta : (bool) estimate time of arrival
- delay : (float)
minimum time before printing anything,
if negative use the default value
as set in
vedo.settings.progressbar_delay - width : (int) width of the progress bar
- char : (str) character to use for the progress bar
- char_back : (str) character to use for the background of the progress bar
Example:
import time from vedo import ProgressBar pb = ProgressBar(0,40, c='r') for i in pb.range(): time.sleep(0.1) pb.print()
def
print(self, txt='', c=None) -> None:
322 def print(self, txt="", c=None) -> None: 323 """Print the progress bar with an optional message.""" 324 if not c: 325 c = self.color 326 327 self._update(self._counts + self.step) 328 329 if self.delay: 330 if time.time() - self.t0 < self.delay: 331 return 332 333 if self.pbar != self._oldbar: 334 self._oldbar = self.pbar 335 336 if self.eta and self._counts > 1: 337 338 tdenom = time.time() - self.t0 339 if tdenom: 340 vel = self._counts / tdenom 341 self._remaining = (self.stop - self._counts) / vel 342 else: 343 vel = 1 344 self._remaining = 0.0 345 346 if self._remaining > 60: 347 _mins = int(self._remaining / 60) 348 _secs = self._remaining - 60 * _mins 349 mins = f"{_mins}m" 350 secs = f"{int(_secs + 0.5)}s " 351 else: 352 mins = "" 353 secs = f"{int(self._remaining + 0.5)}s " 354 355 vel = round(vel, 1) 356 eta = f"eta: {mins}{secs}({vel} it/s) " 357 if self._remaining < 0.5: 358 dt = time.time() - self.t0 359 if dt > 60: 360 _mins = int(dt / 60) 361 _secs = dt - 60 * _mins 362 mins = f"{_mins}m" 363 secs = f"{int(_secs + 0.5)}s " 364 else: 365 mins = "" 366 secs = f"{int(dt + 0.5)}s " 367 eta = f"elapsed: {mins}{secs}({vel} it/s) " 368 txt = "" 369 else: 370 eta = "" 371 372 eraser = " " * self._lentxt + "\b" * self._lentxt 373 374 s = f"{self.pbar} {eraser}{eta}{txt}\r" 375 vedo.printc(s, c=c, bold=self.bold, italic=self.italic, end="") 376 if self.percent > 99.999: 377 print("") 378 379 self._lentxt = len(txt)
Print the progress bar with an optional message.
def
progressbar( iterable, c=None, bold=True, italic=False, title='', eta=True, width=25, delay=-1):
413def progressbar( 414 iterable, 415 c=None, bold=True, italic=False, title="", 416 eta=True, width=25, delay=-1, 417 ): 418 """ 419 Function to print a progress bar with optional text message. 420 421 Use delay to set a minimum time before printing anything. 422 If delay is negative, then use the default value 423 as set in `vedo.settings.progressbar_delay`. 424 425 Arguments: 426 start : (int) 427 starting value 428 stop : (int) 429 stopping value 430 step : (int) 431 step value 432 c : (str) 433 color in hex format 434 title : (str) 435 title text 436 eta : (bool) 437 estimate time of arrival 438 delay : (float) 439 minimum time before printing anything, 440 if negative use the default value 441 set in `vedo.settings.progressbar_delay` 442 width : (int) 443 width of the progress bar 444 char : (str) 445 character to use for the progress bar 446 char_back : (str) 447 character to use for the background of the progress bar 448 449 Example: 450 ```python 451 import time 452 for i in progressbar(range(100), c='r'): 453 time.sleep(0.1) 454 ``` 455  456 """ 457 try: 458 if is_number(iterable): 459 total = int(iterable) 460 iterable = range(total) 461 else: 462 total = len(iterable) 463 except TypeError: 464 iterable = list(iterable) 465 total = len(iterable) 466 467 pb = ProgressBar( 468 0, total, c=c, bold=bold, italic=italic, title=title, 469 eta=eta, delay=delay, width=width, 470 ) 471 for item in iterable: 472 pb.print() 473 yield item
Function to print a progress bar with optional text message.
Use delay to set a minimum time before printing anything.
If delay is negative, then use the default value
as set in vedo.settings.progressbar_delay.
Arguments:
- start : (int) starting value
- stop : (int) stopping value
- step : (int) step value
- c : (str) color in hex format
- title : (str) title text
- eta : (bool) estimate time of arrival
- delay : (float)
minimum time before printing anything,
if negative use the default value
set in
vedo.settings.progressbar_delay - width : (int) width of the progress bar
- char : (str) character to use for the progress bar
- char_back : (str) character to use for the background of the progress bar
Example:
import time for i in progressbar(range(100), c='r'): time.sleep(0.1)
class
Minimizer:
477class Minimizer: 478 """ 479 A function minimizer that uses the Nelder-Mead method. 480 481 The algorithm constructs an n-dimensional simplex in parameter 482 space (i.e. a tetrahedron if the number or parameters is 3) 483 and moves the vertices around parameter space until 484 a local minimum is found. The amoeba method is robust, 485 reasonably efficient, but is not guaranteed to find 486 the global minimum if several local minima exist. 487 488 Arguments: 489 function : (callable) 490 the function to minimize 491 max_iterations : (int) 492 the maximum number of iterations 493 contraction_ratio : (float) 494 The contraction ratio. 495 The default value of 0.5 gives fast convergence, 496 but larger values such as 0.6 or 0.7 provide greater stability. 497 expansion_ratio : (float) 498 The expansion ratio. 499 The default value is 2.0, which provides rapid expansion. 500 Values between 1.1 and 2.0 are valid. 501 tol : (float) 502 the tolerance for convergence 503 504 Example: 505 - [nelder-mead.py](https://github.com/marcomusy/vedo/blob/master/examples/others/nelder-mead.py) 506 """ 507 def __init__( 508 self, 509 function=None, 510 max_iterations=10000, 511 contraction_ratio=0.5, 512 expansion_ratio=2.0, 513 tol=1e-5, 514 ) -> None: 515 self.function = function 516 self.tolerance = tol 517 self.contraction_ratio = contraction_ratio 518 self.expansion_ratio = expansion_ratio 519 self.max_iterations = max_iterations 520 self.minimizer = vtki.new("AmoebaMinimizer") 521 self.minimizer.SetFunction(self._vtkfunc) 522 self.results = {} 523 self.parameters_path = [] 524 self.function_path = [] 525 526 def _vtkfunc(self): 527 n = self.minimizer.GetNumberOfParameters() 528 ain = [self.minimizer.GetParameterValue(i) for i in range(n)] 529 r = self.function(ain) 530 self.minimizer.SetFunctionValue(r) 531 self.parameters_path.append(ain) 532 self.function_path.append(r) 533 return r 534 535 def eval(self, parameters=()) -> float: 536 """ 537 Evaluate the function at the current or given parameters. 538 """ 539 if len(parameters) == 0: 540 return self.minimizer.EvaluateFunction() 541 self.set_parameters(parameters) 542 return self.function(parameters) 543 544 def set_parameter(self, name, value, scale=1.0) -> None: 545 """ 546 Set the parameter value. 547 The initial amount by which the parameter 548 will be modified during the search for the minimum. 549 """ 550 self.minimizer.SetParameterValue(name, value) 551 self.minimizer.SetParameterScale(name, scale) 552 553 def set_parameters(self, parameters) -> None: 554 """ 555 Set the parameters names and values from a dictionary. 556 """ 557 for name, value in parameters.items(): 558 if len(value) == 2: 559 self.set_parameter(name, value[0], value[1]) 560 else: 561 self.set_parameter(name, value) 562 563 def minimize(self) -> dict: 564 """ 565 Minimize the input function. 566 567 Returns: 568 dict : 569 the minimization results 570 init_parameters : (dict) 571 the initial parameters 572 parameters : (dict) 573 the final parameters 574 min_value : (float) 575 the minimum value 576 iterations : (int) 577 the number of iterations 578 max_iterations : (int) 579 the maximum number of iterations 580 tolerance : (float) 581 the tolerance for convergence 582 convergence_flag : (int) 583 zero if the tolerance stopping 584 criterion has been met. 585 parameters_path : (np.array) 586 the path of the minimization 587 algorithm in parameter space 588 function_path : (np.array) 589 the path of the minimization 590 algorithm in function space 591 hessian : (np.array) 592 the Hessian matrix of the 593 function at the minimum 594 parameter_errors : (np.array) 595 the errors on the parameters 596 """ 597 n = self.minimizer.GetNumberOfParameters() 598 out = [( 599 self.minimizer.GetParameterName(i), 600 (self.minimizer.GetParameterValue(i), 601 self.minimizer.GetParameterScale(i)) 602 ) for i in range(n)] 603 self.results["init_parameters"] = dict(out) 604 605 self.minimizer.SetTolerance(self.tolerance) 606 self.minimizer.SetContractionRatio(self.contraction_ratio) 607 self.minimizer.SetExpansionRatio(self.expansion_ratio) 608 self.minimizer.SetMaxIterations(self.max_iterations) 609 610 self.minimizer.Minimize() 611 self.results["convergence_flag"] = not bool(self.minimizer.Iterate()) 612 613 out = [( 614 self.minimizer.GetParameterName(i), 615 self.minimizer.GetParameterValue(i), 616 ) for i in range(n)] 617 618 self.results["parameters"] = dict(out) 619 self.results["min_value"] = self.minimizer.GetFunctionValue() 620 self.results["iterations"] = self.minimizer.GetIterations() 621 self.results["max_iterations"] = self.minimizer.GetMaxIterations() 622 self.results["tolerance"] = self.minimizer.GetTolerance() 623 self.results["expansion_ratio"] = self.expansion_ratio 624 self.results["contraction_ratio"] = self.contraction_ratio 625 self.results["parameters_path"] = np.array(self.parameters_path) 626 self.results["function_path"] = np.array(self.function_path) 627 self.results["hessian"] = np.zeros((n,n)) 628 self.results["parameter_errors"] = np.zeros(n) 629 return self.results 630 631 def compute_hessian(self, epsilon=0) -> np.array: 632 """ 633 Compute the Hessian matrix of `function` at the 634 minimum numerically. 635 636 Arguments: 637 epsilon : (float) 638 Step size used for numerical approximation. 639 640 Returns: 641 array: Hessian matrix of `function` at minimum. 642 """ 643 if not epsilon: 644 epsilon = self.tolerance * 10 645 n = self.minimizer.GetNumberOfParameters() 646 x0 = [self.minimizer.GetParameterValue(i) for i in range(n)] 647 hessian = np.zeros((n, n)) 648 for i in vedo.progressbar(n, title="Computing Hessian", delay=2): 649 for j in range(n): 650 xijp = np.copy(x0) 651 xijp[i] += epsilon 652 xijp[j] += epsilon 653 xijm = np.copy(x0) 654 xijm[i] += epsilon 655 xijm[j] -= epsilon 656 xjip = np.copy(x0) 657 xjip[i] -= epsilon 658 xjip[j] += epsilon 659 xjim = np.copy(x0) 660 xjim[i] -= epsilon 661 xjim[j] -= epsilon 662 # Second derivative approximation 663 fijp = self.function(xijp) 664 fijm = self.function(xijm) 665 fjip = self.function(xjip) 666 fjim = self.function(xjim) 667 hessian[i, j] = (fijp - fijm - fjip + fjim) / (2 * epsilon**2) 668 self.results["hessian"] = hessian 669 try: 670 ihess = np.linalg.inv(hessian) 671 self.results["parameter_errors"] = np.sqrt(np.diag(ihess)) 672 except: 673 vedo.logger.warning("Cannot compute hessian for parameter errors") 674 self.results["parameter_errors"] = np.zeros(n) 675 return hessian 676 677 def __str__(self) -> str: 678 out = vedo.printc( 679 f"vedo.utils.Minimizer at ({hex(id(self))})".ljust(75), 680 bold=True, invert=True, return_string=True, 681 ) 682 out += "Function name".ljust(20) + self.function.__name__ + "()\n" 683 out += "-------- parameters initial value -----------\n" 684 out += "Name".ljust(20) + "Value".ljust(20) + "Scale\n" 685 for name, value in self.results["init_parameters"].items(): 686 out += name.ljust(20) + str(value[0]).ljust(20) + str(value[1]) + "\n" 687 out += "-------- parameters final value --------------\n" 688 for name, value in self.results["parameters"].items(): 689 out += name.ljust(20) + f"{value:.6f}" 690 ierr = list(self.results["parameters"]).index(name) 691 err = self.results["parameter_errors"][ierr] 692 if err: 693 out += f" ± {err:.4f}" 694 out += "\n" 695 out += "Value at minimum".ljust(20)+ f'{self.results["min_value"]}\n' 696 out += "Iterations".ljust(20) + f'{self.results["iterations"]}\n' 697 out += "Max iterations".ljust(20) + f'{self.results["max_iterations"]}\n' 698 out += "Convergence flag".ljust(20)+ f'{self.results["convergence_flag"]}\n' 699 out += "Tolerance".ljust(20) + f'{self.results["tolerance"]}\n' 700 try: 701 arr = np.array2string( 702 self.compute_hessian(), 703 separator=', ', precision=6, suppress_small=True, 704 ) 705 out += "Hessian Matrix:\n" + arr 706 except: 707 out += "Hessian Matrix: (not available)" 708 return out
A function minimizer that uses the Nelder-Mead method.
The algorithm constructs an n-dimensional simplex in parameter space (i.e. a tetrahedron if the number or parameters is 3) and moves the vertices around parameter space until a local minimum is found. The amoeba method is robust, reasonably efficient, but is not guaranteed to find the global minimum if several local minima exist.
Arguments:
- function : (callable) the function to minimize
- max_iterations : (int) the maximum number of iterations
- contraction_ratio : (float) The contraction ratio. The default value of 0.5 gives fast convergence, but larger values such as 0.6 or 0.7 provide greater stability.
- expansion_ratio : (float) The expansion ratio. The default value is 2.0, which provides rapid expansion. Values between 1.1 and 2.0 are valid.
- tol : (float) the tolerance for convergence
Example:
Minimizer( function=None, max_iterations=10000, contraction_ratio=0.5, expansion_ratio=2.0, tol=1e-05)
507 def __init__( 508 self, 509 function=None, 510 max_iterations=10000, 511 contraction_ratio=0.5, 512 expansion_ratio=2.0, 513 tol=1e-5, 514 ) -> None: 515 self.function = function 516 self.tolerance = tol 517 self.contraction_ratio = contraction_ratio 518 self.expansion_ratio = expansion_ratio 519 self.max_iterations = max_iterations 520 self.minimizer = vtki.new("AmoebaMinimizer") 521 self.minimizer.SetFunction(self._vtkfunc) 522 self.results = {} 523 self.parameters_path = [] 524 self.function_path = []
def
eval(self, parameters=()) -> float:
535 def eval(self, parameters=()) -> float: 536 """ 537 Evaluate the function at the current or given parameters. 538 """ 539 if len(parameters) == 0: 540 return self.minimizer.EvaluateFunction() 541 self.set_parameters(parameters) 542 return self.function(parameters)
Evaluate the function at the current or given parameters.
def
set_parameter(self, name, value, scale=1.0) -> None:
544 def set_parameter(self, name, value, scale=1.0) -> None: 545 """ 546 Set the parameter value. 547 The initial amount by which the parameter 548 will be modified during the search for the minimum. 549 """ 550 self.minimizer.SetParameterValue(name, value) 551 self.minimizer.SetParameterScale(name, scale)
Set the parameter value. The initial amount by which the parameter will be modified during the search for the minimum.
def
set_parameters(self, parameters) -> None:
553 def set_parameters(self, parameters) -> None: 554 """ 555 Set the parameters names and values from a dictionary. 556 """ 557 for name, value in parameters.items(): 558 if len(value) == 2: 559 self.set_parameter(name, value[0], value[1]) 560 else: 561 self.set_parameter(name, value)
Set the parameters names and values from a dictionary.
def
minimize(self) -> dict:
563 def minimize(self) -> dict: 564 """ 565 Minimize the input function. 566 567 Returns: 568 dict : 569 the minimization results 570 init_parameters : (dict) 571 the initial parameters 572 parameters : (dict) 573 the final parameters 574 min_value : (float) 575 the minimum value 576 iterations : (int) 577 the number of iterations 578 max_iterations : (int) 579 the maximum number of iterations 580 tolerance : (float) 581 the tolerance for convergence 582 convergence_flag : (int) 583 zero if the tolerance stopping 584 criterion has been met. 585 parameters_path : (np.array) 586 the path of the minimization 587 algorithm in parameter space 588 function_path : (np.array) 589 the path of the minimization 590 algorithm in function space 591 hessian : (np.array) 592 the Hessian matrix of the 593 function at the minimum 594 parameter_errors : (np.array) 595 the errors on the parameters 596 """ 597 n = self.minimizer.GetNumberOfParameters() 598 out = [( 599 self.minimizer.GetParameterName(i), 600 (self.minimizer.GetParameterValue(i), 601 self.minimizer.GetParameterScale(i)) 602 ) for i in range(n)] 603 self.results["init_parameters"] = dict(out) 604 605 self.minimizer.SetTolerance(self.tolerance) 606 self.minimizer.SetContractionRatio(self.contraction_ratio) 607 self.minimizer.SetExpansionRatio(self.expansion_ratio) 608 self.minimizer.SetMaxIterations(self.max_iterations) 609 610 self.minimizer.Minimize() 611 self.results["convergence_flag"] = not bool(self.minimizer.Iterate()) 612 613 out = [( 614 self.minimizer.GetParameterName(i), 615 self.minimizer.GetParameterValue(i), 616 ) for i in range(n)] 617 618 self.results["parameters"] = dict(out) 619 self.results["min_value"] = self.minimizer.GetFunctionValue() 620 self.results["iterations"] = self.minimizer.GetIterations() 621 self.results["max_iterations"] = self.minimizer.GetMaxIterations() 622 self.results["tolerance"] = self.minimizer.GetTolerance() 623 self.results["expansion_ratio"] = self.expansion_ratio 624 self.results["contraction_ratio"] = self.contraction_ratio 625 self.results["parameters_path"] = np.array(self.parameters_path) 626 self.results["function_path"] = np.array(self.function_path) 627 self.results["hessian"] = np.zeros((n,n)) 628 self.results["parameter_errors"] = np.zeros(n) 629 return self.results
Minimize the input function.
Returns:
dict : the minimization results init_parameters : (dict) the initial parameters parameters : (dict) the final parameters min_value : (float) the minimum value iterations : (int) the number of iterations max_iterations : (int) the maximum number of iterations tolerance : (float) the tolerance for convergence convergence_flag : (int) zero if the tolerance stopping criterion has been met. parameters_path : (np.array) the path of the minimization algorithm in parameter space function_path : (np.array) the path of the minimization algorithm in function space hessian : (np.array) the Hessian matrix of the function at the minimum parameter_errors : (np.array) the errors on the parameters
def
compute_hessian(self, epsilon=0) -> <built-in function array>:
631 def compute_hessian(self, epsilon=0) -> np.array: 632 """ 633 Compute the Hessian matrix of `function` at the 634 minimum numerically. 635 636 Arguments: 637 epsilon : (float) 638 Step size used for numerical approximation. 639 640 Returns: 641 array: Hessian matrix of `function` at minimum. 642 """ 643 if not epsilon: 644 epsilon = self.tolerance * 10 645 n = self.minimizer.GetNumberOfParameters() 646 x0 = [self.minimizer.GetParameterValue(i) for i in range(n)] 647 hessian = np.zeros((n, n)) 648 for i in vedo.progressbar(n, title="Computing Hessian", delay=2): 649 for j in range(n): 650 xijp = np.copy(x0) 651 xijp[i] += epsilon 652 xijp[j] += epsilon 653 xijm = np.copy(x0) 654 xijm[i] += epsilon 655 xijm[j] -= epsilon 656 xjip = np.copy(x0) 657 xjip[i] -= epsilon 658 xjip[j] += epsilon 659 xjim = np.copy(x0) 660 xjim[i] -= epsilon 661 xjim[j] -= epsilon 662 # Second derivative approximation 663 fijp = self.function(xijp) 664 fijm = self.function(xijm) 665 fjip = self.function(xjip) 666 fjim = self.function(xjim) 667 hessian[i, j] = (fijp - fijm - fjip + fjim) / (2 * epsilon**2) 668 self.results["hessian"] = hessian 669 try: 670 ihess = np.linalg.inv(hessian) 671 self.results["parameter_errors"] = np.sqrt(np.diag(ihess)) 672 except: 673 vedo.logger.warning("Cannot compute hessian for parameter errors") 674 self.results["parameter_errors"] = np.zeros(n) 675 return hessian
Compute the Hessian matrix of function at the
minimum numerically.
Arguments:
- epsilon : (float) Step size used for numerical approximation.
Returns:
array: Hessian matrix of
functionat minimum.
854def geometry(obj, extent=None) -> "vedo.Mesh": 855 """ 856 Apply the `vtkGeometryFilter` to the input object. 857 This is a general-purpose filter to extract geometry (and associated data) 858 from any type of dataset. 859 This filter also may be used to convert any type of data to polygonal type. 860 The conversion process may be less than satisfactory for some 3D datasets. 861 For example, this filter will extract the outer surface of a volume 862 or structured grid dataset. 863 864 Returns a `vedo.Mesh` object. 865 866 Set `extent` as the `[xmin,xmax, ymin,ymax, zmin,zmax]` bounding box to clip data. 867 """ 868 gf = vtki.new("GeometryFilter") 869 gf.SetInputData(obj) 870 if extent is not None: 871 gf.SetExtent(extent) 872 gf.Update() 873 return vedo.Mesh(gf.GetOutput())
Apply the vtkGeometryFilter to the input object.
This is a general-purpose filter to extract geometry (and associated data)
from any type of dataset.
This filter also may be used to convert any type of data to polygonal type.
The conversion process may be less than satisfactory for some 3D datasets.
For example, this filter will extract the outer surface of a volume
or structured grid dataset.
Returns a vedo.Mesh object.
Set extent as the [xmin,xmax, ymin,ymax, zmin,zmax] bounding box to clip data.
def
is_sequence(arg) -> bool:
1045def is_sequence(arg) -> bool: 1046 """Check if the input is iterable.""" 1047 if hasattr(arg, "strip"): 1048 return False 1049 if hasattr(arg, "__getslice__"): 1050 return True 1051 if hasattr(arg, "__iter__"): 1052 return True 1053 return False
Check if the input is iterable.
def
lin_interpolate(x, rangeX, rangeY):
1419def lin_interpolate(x, rangeX, rangeY): 1420 """ 1421 Interpolate linearly the variable `x` in `rangeX` onto the new `rangeY`. 1422 If `x` is a 3D vector the linear weight is the distance to the two 3D `rangeX` vectors. 1423 1424 E.g. if `x` runs in `rangeX=[x0,x1]` and I want it to run in `rangeY=[y0,y1]` then 1425 1426 `y = lin_interpolate(x, rangeX, rangeY)` will interpolate `x` onto `rangeY`. 1427 1428 Examples: 1429 - [lin_interpolate.py](https://github.com/marcomusy/vedo/tree/master/examples/basic/lin_interpolate.py) 1430 1431  1432 """ 1433 if is_sequence(x): 1434 x = np.asarray(x) 1435 x0, x1 = np.asarray(rangeX) 1436 y0, y1 = np.asarray(rangeY) 1437 dx = x1 - x0 1438 dxn = np.linalg.norm(dx) 1439 if not dxn: 1440 return y0 1441 s = np.linalg.norm(x - x0) / dxn 1442 t = np.linalg.norm(x - x1) / dxn 1443 st = s + t 1444 out = y0 * (t / st) + y1 * (s / st) 1445 1446 else: # faster 1447 1448 x0 = rangeX[0] 1449 dx = rangeX[1] - x0 1450 if not dx: 1451 return rangeY[0] 1452 s = (x - x0) / dx 1453 out = rangeY[0] * (1 - s) + rangeY[1] * s 1454 return out
Interpolate linearly the variable x in rangeX onto the new rangeY.
If x is a 3D vector the linear weight is the distance to the two 3D rangeX vectors.
E.g. if x runs in rangeX=[x0,x1] and I want it to run in rangeY=[y0,y1] then
y = lin_interpolate(x, rangeX, rangeY) will interpolate x onto rangeY.
Examples:
def
vector(x, y=None, z=0.0, dtype=<class 'numpy.float64'>) -> numpy.ndarray:
1531def vector(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1532 """ 1533 Return a 3D numpy array representing a vector. 1534 1535 If `y` is `None`, assume input is already in the form `[x,y,z]`. 1536 """ 1537 if y is None: # assume x is already [x,y,z] 1538 return np.asarray(x, dtype=dtype) 1539 return np.array([x, y, z], dtype=dtype)
Return a 3D numpy array representing a vector.
If y is None, assume input is already in the form [x,y,z].
def
mag(v):
1550def mag(v): 1551 """Get the magnitude of a vector or array of vectors.""" 1552 v = np.asarray(v) 1553 if v.ndim == 1: 1554 return np.linalg.norm(v) 1555 return np.linalg.norm(v, axis=1)
Get the magnitude of a vector or array of vectors.
def
mag2(v) -> numpy.ndarray:
1558def mag2(v) -> np.ndarray: 1559 """Get the squared magnitude of a vector or array of vectors.""" 1560 v = np.asarray(v) 1561 if v.ndim == 1: 1562 return np.square(v).sum() 1563 return np.square(v).sum(axis=1)
Get the squared magnitude of a vector or array of vectors.
def
versor(x, y=None, z=0.0, dtype=<class 'numpy.float64'>) -> numpy.ndarray:
1542def versor(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1543 """Return the unit vector. Input can be a list of vectors.""" 1544 v = vector(x, y, z, dtype) 1545 if isinstance(v[0], np.ndarray): 1546 return np.divide(v, mag(v)[:, None]) 1547 return v / mag(v)
Return the unit vector. Input can be a list of vectors.
def
precision(x, p: int, vrange=None, delimiter='e') -> str:
1632def precision(x, p: int, vrange=None, delimiter="e") -> str: 1633 """ 1634 Returns a string representation of `x` formatted to precision `p`. 1635 1636 Set `vrange` to the range in which x exists (to snap x to '0' if below precision). 1637 """ 1638 # Based on the webkit javascript implementation 1639 # `from here <https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp>`_, 1640 # and implemented by `randlet <https://github.com/randlet/to-precision>`_. 1641 # Modified for vedo by M.Musy 2020 1642 1643 if isinstance(x, str): # do nothing 1644 return x 1645 1646 if is_sequence(x): 1647 out = "(" 1648 nn = len(x) - 1 1649 for i, ix in enumerate(x): 1650 1651 try: 1652 if np.isnan(ix): 1653 return "NaN" 1654 except: 1655 # cannot handle list of list 1656 continue 1657 1658 out += precision(ix, p) 1659 if i < nn: 1660 out += ", " 1661 return out + ")" ############ <-- 1662 1663 try: 1664 if np.isnan(x): 1665 return "NaN" 1666 except TypeError: 1667 return "NaN" 1668 1669 x = float(x) 1670 1671 if x == 0.0 or (vrange is not None and abs(x) < vrange / pow(10, p)): 1672 return "0" 1673 1674 out = [] 1675 if x < 0: 1676 out.append("-") 1677 x = -x 1678 1679 e = int(np.log10(x)) 1680 # tens = np.power(10, e - p + 1) 1681 tens = 10 ** (e - p + 1) 1682 n = np.floor(x / tens) 1683 1684 # if n < np.power(10, p - 1): 1685 if n < 10 ** (p - 1): 1686 e = e - 1 1687 # tens = np.power(10, e - p + 1) 1688 tens = 10 ** (e - p + 1) 1689 n = np.floor(x / tens) 1690 1691 if abs((n + 1.0) * tens - x) <= abs(n * tens - x): 1692 n = n + 1 1693 1694 # if n >= np.power(10, p): 1695 if n >= 10 ** p: 1696 n = n / 10.0 1697 e = e + 1 1698 1699 m = "%.*g" % (p, n) 1700 if e < -2 or e >= p: 1701 out.append(m[0]) 1702 if p > 1: 1703 out.append(".") 1704 out.extend(m[1:p]) 1705 out.append(delimiter) 1706 if e > 0: 1707 out.append("+") 1708 out.append(str(e)) 1709 elif e == (p - 1): 1710 out.append(m) 1711 elif e >= 0: 1712 out.append(m[: e + 1]) 1713 if e + 1 < len(m): 1714 out.append(".") 1715 out.extend(m[e + 1 :]) 1716 else: 1717 out.append("0.") 1718 out.extend(["0"] * -(e + 1)) 1719 out.append(m) 1720 return "".join(out)
Returns a string representation of x formatted to precision p.
Set vrange to the range in which x exists (to snap x to '0' if below precision).
def
round_to_digit(x, p) -> float:
1585def round_to_digit(x, p) -> float: 1586 """Round a real number to the specified number of significant digits.""" 1587 if not x: 1588 return 0 1589 r = np.round(x, p - int(np.floor(np.log10(abs(x)))) - 1) 1590 if int(r) == r: 1591 return int(r) 1592 return r
Round a real number to the specified number of significant digits.
def
point_in_triangle(p, p1, p2, p3) -> Optional[bool]:
1128def point_in_triangle(p, p1, p2, p3) -> Union[bool, None]: 1129 """ 1130 Return True if a point is inside (or above/below) 1131 a triangle defined by 3 points in space. 1132 """ 1133 p1 = np.array(p1) 1134 u = p2 - p1 1135 v = p3 - p1 1136 n = np.cross(u, v) 1137 w = p - p1 1138 ln = np.dot(n, n) 1139 if not ln: 1140 return None # degenerate triangle 1141 gamma = (np.dot(np.cross(u, w), n)) / ln 1142 if 0 < gamma < 1: 1143 beta = (np.dot(np.cross(w, v), n)) / ln 1144 if 0 < beta < 1: 1145 alpha = 1 - gamma - beta 1146 if 0 < alpha < 1: 1147 return True 1148 return False
Return True if a point is inside (or above/below) a triangle defined by 3 points in space.
def
point_line_distance(p, p1, p2) -> float:
1341def point_line_distance(p, p1, p2) -> float: 1342 """ 1343 Compute the distance of a point to a line (not the segment) 1344 defined by `p1` and `p2`. 1345 """ 1346 return np.sqrt(vtki.vtkLine.DistanceToLine(p, p1, p2))
Compute the distance of a point to a line (not the segment)
defined by p1 and p2.
def
closest(point, points, n=1, return_ids=False, use_tree=False):
1381def closest(point, points, n=1, return_ids=False, use_tree=False): 1382 """ 1383 Returns the distances and the closest point(s) to the given set of points. 1384 Needs `scipy.spatial` library. 1385 1386 Arguments: 1387 n : (int) 1388 the nr of closest points to return 1389 return_ids : (bool) 1390 return the ids instead of the points coordinates 1391 use_tree : (bool) 1392 build a `scipy.spatial.KDTree`. 1393 An already existing one can be passed to avoid rebuilding. 1394 """ 1395 from scipy.spatial import distance, KDTree 1396 1397 points = np.asarray(points) 1398 if n == 1: 1399 dists = distance.cdist([point], points) 1400 closest_idx = np.argmin(dists) 1401 else: 1402 if use_tree: 1403 if isinstance(use_tree, KDTree): # reuse 1404 tree = use_tree 1405 else: 1406 tree = KDTree(points) 1407 dists, closest_idx = tree.query([point], k=n) 1408 closest_idx = closest_idx[0] 1409 else: 1410 dists = distance.cdist([point], points) 1411 closest_idx = np.argsort(dists)[0][:n] 1412 if return_ids: 1413 return dists, closest_idx 1414 else: 1415 return dists, points[closest_idx]
Returns the distances and the closest point(s) to the given set of points.
Needs scipy.spatial library.
Arguments:
- n : (int) the nr of closest points to return
- return_ids : (bool) return the ids instead of the points coordinates
- use_tree : (bool)
build a
scipy.spatial.KDTree. An already existing one can be passed to avoid rebuilding.
def
grep( filename: str, tag: str, column=None, first_occurrence_only=False) -> list:
1724def grep(filename: str, tag: str, column=None, first_occurrence_only=False) -> list: 1725 """Greps the line in a file that starts with a specific `tag` string inside the file.""" 1726 import re 1727 1728 with open(filename, "r", encoding="UTF-8") as afile: 1729 content = [] 1730 for line in afile: 1731 if re.search(tag, line): 1732 c = line.split() 1733 c[-1] = c[-1].replace("\n", "") 1734 if column is not None: 1735 c = c[column] 1736 content.append(c) 1737 if first_occurrence_only: 1738 break 1739 return content
Greps the line in a file that starts with a specific tag string inside the file.
def
make_bands(inputlist, n):
2022def make_bands(inputlist, n): 2023 """ 2024 Group values of a list into bands of equal value, where 2025 `n` is the number of bands, a positive integer > 2. 2026 2027 Returns a binned list of the same length as the input. 2028 """ 2029 if n < 2: 2030 return inputlist 2031 vmin = np.min(inputlist) 2032 vmax = np.max(inputlist) 2033 bb = np.linspace(vmin, vmax, n, endpoint=0) 2034 dr = bb[1] - bb[0] 2035 bb += dr / 2 2036 tol = dr / 2 * 1.001 2037 newlist = [] 2038 for s in inputlist: 2039 for b in bb: 2040 if abs(s - b) < tol: 2041 newlist.append(b) 2042 break 2043 return np.array(newlist)
Group values of a list into bands of equal value, where
n is the number of bands, a positive integer > 2.
Returns a binned list of the same length as the input.
def
pack_spheres(bounds, radius) -> numpy.ndarray:
1595def pack_spheres(bounds, radius) -> np.ndarray: 1596 """ 1597 Packing spheres into a bounding box. 1598 Returns a numpy array of sphere centers. 1599 """ 1600 h = 0.8164965 / 2 1601 d = 0.8660254 1602 a = 0.288675135 1603 1604 if is_sequence(bounds): 1605 x0, x1, y0, y1, z0, z1 = bounds 1606 else: 1607 x0, x1, y0, y1, z0, z1 = bounds.bounds() 1608 1609 x = np.arange(x0, x1, radius) 1610 nul = np.zeros_like(x) 1611 nz = int((z1 - z0) / radius / h / 2 + 1.5) 1612 ny = int((y1 - y0) / radius / d + 1.5) 1613 1614 pts = [] 1615 for iz in range(nz): 1616 z = z0 + nul + iz * h * radius 1617 dx, dy, dz = [radius * 0.5, radius * a, iz * h * radius] 1618 for iy in range(ny): 1619 y = y0 + nul + iy * d * radius 1620 if iy % 2: 1621 xs = x 1622 else: 1623 xs = x + radius * 0.5 1624 if iz % 2: 1625 p = np.c_[xs, y, z] + [dx, dy, dz] 1626 else: 1627 p = np.c_[xs, y, z] + [0, 0, dz] 1628 pts += p.tolist() 1629 return np.array(pts)
Packing spheres into a bounding box. Returns a numpy array of sphere centers.
def
humansort(alist) -> list:
1095def humansort(alist) -> list: 1096 """ 1097 Sort in place a given list the way humans expect. 1098 1099 E.g. `['file11', 'file1'] -> ['file1', 'file11']` 1100 1101 .. warning:: input list is modified in-place by this function. 1102 """ 1103 import re 1104 1105 def alphanum_key(s): 1106 # Turn a string into a list of string and number chunks. 1107 # e.g. "z23a" -> ["z", 23, "a"] 1108 def tryint(s): 1109 if s.isdigit(): 1110 return int(s) 1111 return s 1112 1113 return [tryint(c) for c in re.split("([0-9]+)", s)] 1114 1115 alist.sort(key=alphanum_key) 1116 return alist # NB: input list is modified
Sort in place a given list the way humans expect.
E.g. ['file11', 'file1'] -> ['file1', 'file11']
input list is modified in-place by this function.
def
print_histogram( data, bins=10, height=10, logscale=False, minbin=0, horizontal=True, char='▉', c=None, bold=True, title='histogram', spacer='') -> numpy.ndarray:
1797def print_histogram( 1798 data, 1799 bins=10, 1800 height=10, 1801 logscale=False, 1802 minbin=0, 1803 horizontal=True, 1804 char="\U00002589", 1805 c=None, 1806 bold=True, 1807 title="histogram", 1808 spacer="", 1809) -> np.ndarray: 1810 """ 1811 Ascii histogram printing. 1812 1813 Input can be a `vedo.Volume` or `vedo.Mesh`. 1814 Returns the raw data before binning (useful when passing vtk objects). 1815 1816 Arguments: 1817 bins : (int) 1818 number of histogram bins 1819 height : (int) 1820 height of the histogram in character units 1821 logscale : (bool) 1822 use logscale for frequencies 1823 minbin : (int) 1824 ignore bins before minbin 1825 horizontal : (bool) 1826 show histogram horizontally 1827 char : (str) 1828 character to be used 1829 bold : (bool) 1830 use boldface 1831 title : (str) 1832 histogram title 1833 spacer : (str) 1834 horizontal spacer 1835 1836 Example: 1837 ```python 1838 from vedo import print_histogram 1839 import numpy as np 1840 d = np.random.normal(size=1000) 1841 data = print_histogram(d, c='b', logscale=True, title='my scalars') 1842 data = print_histogram(d, c='o') 1843 print(np.mean(data)) # data here is same as d 1844 ``` 1845  1846 """ 1847 # credits: http://pyinsci.blogspot.com/2009/10/ascii-histograms.html 1848 # adapted for vedo by M.Musy, 2019 1849 1850 if not horizontal: # better aspect ratio 1851 bins *= 2 1852 1853 try: 1854 data = vtk2numpy(data.dataset.GetPointData().GetScalars()) 1855 except AttributeError: 1856 # already an array 1857 data = np.asarray(data) 1858 1859 if isinstance(data, vtki.vtkImageData): 1860 dims = data.GetDimensions() 1861 nvx = min(100000, dims[0] * dims[1] * dims[2]) 1862 idxs = np.random.randint(0, min(dims), size=(nvx, 3)) 1863 data = [] 1864 for ix, iy, iz in idxs: 1865 d = data.GetScalarComponentAsFloat(ix, iy, iz, 0) 1866 data.append(d) 1867 data = np.array(data) 1868 1869 elif isinstance(data, vtki.vtkPolyData): 1870 arr = data.GetPointData().GetScalars() 1871 if not arr: 1872 arr = data.GetCellData().GetScalars() 1873 if not arr: 1874 return np.array([]) 1875 data = vtk2numpy(arr) 1876 1877 try: 1878 h = np.histogram(data, bins=bins) 1879 except TypeError as e: 1880 vedo.logger.error(f"cannot compute histogram: {e}") 1881 return np.array([]) 1882 1883 if minbin: 1884 hi = h[0][minbin:-1] 1885 else: 1886 hi = h[0] 1887 1888 if char == "\U00002589" and horizontal: 1889 char = "\U00002586" 1890 1891 title = title.ljust(14) + ":" 1892 entrs = " entries=" + str(len(data)) 1893 if logscale: 1894 h0 = np.log10(hi + 1) 1895 maxh0 = int(max(h0) * 100) / 100 1896 title = title + entrs + " (logscale)" 1897 else: 1898 h0 = hi 1899 maxh0 = max(h0) 1900 title = title + entrs 1901 1902 def _v(): 1903 his = "" 1904 if title: 1905 his += title + "\n" 1906 bars = h0 / maxh0 * height 1907 for l in reversed(range(1, height + 1)): 1908 line = "" 1909 if l == height: 1910 line = "%s " % maxh0 1911 else: 1912 line = " |" + " " * (len(str(maxh0)) - 3) 1913 for c in bars: 1914 if c >= np.ceil(l): 1915 line += char 1916 else: 1917 line += " " 1918 line += "\n" 1919 his += line 1920 his += "%.2f" % h[1][0] + "." * (bins) + "%.2f" % h[1][-1] + "\n" 1921 return his 1922 1923 def _h(): 1924 his = "" 1925 if title: 1926 his += title + "\n" 1927 xl = ["%.2f" % n for n in h[1]] 1928 lxl = [len(l) for l in xl] 1929 bars = h0 / maxh0 * height 1930 his += spacer + " " * int(max(bars) + 2 + max(lxl)) + "%s\n" % maxh0 1931 for i, c in enumerate(bars): 1932 line = xl[i] + " " * int(max(lxl) - lxl[i]) + "| " + char * int(c) + "\n" 1933 his += spacer + line 1934 return his 1935 1936 if horizontal: 1937 height *= 2 1938 vedo.printc(_h(), c=c, bold=bold) 1939 else: 1940 vedo.printc(_v(), c=c, bold=bold) 1941 return data
Ascii histogram printing.
Input can be a vedo.Volume or vedo.Mesh.
Returns the raw data before binning (useful when passing vtk objects).
Arguments:
- bins : (int) number of histogram bins
- height : (int) height of the histogram in character units
- logscale : (bool) use logscale for frequencies
- minbin : (int) ignore bins before minbin
- horizontal : (bool) show histogram horizontally
- char : (str) character to be used
- bold : (bool) use boldface
- title : (str) histogram title
- spacer : (str) horizontal spacer
Example:
from vedo import print_histogram import numpy as np d = np.random.normal(size=1000) data = print_histogram(d, c='b', logscale=True, title='my scalars') data = print_histogram(d, c='o') print(np.mean(data)) # data here is same as d
def
print_inheritance_tree(C) -> None:
1998def print_inheritance_tree(C) -> None: 1999 """Prints the inheritance tree of class C.""" 2000 # Adapted from: https://stackoverflow.com/questions/26568976/ 2001 def class_tree(cls): 2002 subc = [class_tree(sub_class) for sub_class in cls.__subclasses__()] 2003 return {cls.__name__: subc} 2004 2005 def print_tree(tree, indent=8, current_ind=0): 2006 for k, v in tree.items(): 2007 if current_ind: 2008 before_dashes = current_ind - indent 2009 m = " " * before_dashes + "└" + "─" * (indent - 1) + " " + k 2010 vedo.printc(m) 2011 else: 2012 vedo.printc(k) 2013 for sub_tree in v: 2014 print_tree(sub_tree, indent=indent, current_ind=current_ind + indent) 2015 2016 if str(C.__class__) != "<class 'type'>": 2017 C = C.__class__ 2018 ct = class_tree(C) 2019 print_tree(ct)
Prints the inheritance tree of class C.
def
camera_from_quaternion( pos, quaternion, distance=10000, ngl_correct=True) -> vtkmodules.vtkRenderingCore.vtkCamera:
2049def camera_from_quaternion(pos, quaternion, distance=10000, ngl_correct=True) -> vtki.vtkCamera: 2050 """ 2051 Define a `vtkCamera` with a particular orientation. 2052 2053 Arguments: 2054 pos: (np.array, list, tuple) 2055 an iterator of length 3 containing the focus point of the camera 2056 quaternion: (np.array, list, tuple) 2057 a len(4) quaternion `(x,y,z,w)` describing the rotation of the camera 2058 such as returned by neuroglancer `x,y,z,w` all in `[0,1]` range 2059 distance: (float) 2060 the desired distance from pos to the camera (default = 10000 nm) 2061 2062 Returns: 2063 `vtki.vtkCamera`, a vtk camera setup according to these rules. 2064 """ 2065 camera = vtki.vtkCamera() 2066 # define the quaternion in vtk, note the swapped order 2067 # w,x,y,z instead of x,y,z,w 2068 quat_vtk = vtki.get_class("Quaternion")( 2069 quaternion[3], quaternion[0], quaternion[1], quaternion[2]) 2070 # use this to define a rotation matrix in x,y,z 2071 # right handed units 2072 M = np.zeros((3, 3), dtype=np.float32) 2073 quat_vtk.ToMatrix3x3(M) 2074 # the default camera orientation is y up 2075 up = [0, 1, 0] 2076 # calculate default camera position is backed off in positive z 2077 pos = [0, 0, distance] 2078 2079 # set the camera rototation by applying the rotation matrix 2080 camera.SetViewUp(*np.dot(M, up)) 2081 # set the camera position by applying the rotation matrix 2082 camera.SetPosition(*np.dot(M, pos)) 2083 if ngl_correct: 2084 # neuroglancer has positive y going down 2085 # so apply these azimuth and roll corrections 2086 # to fix orientatins 2087 camera.Azimuth(-180) 2088 camera.Roll(180) 2089 2090 # shift the camera posiiton and focal position 2091 # to be centered on the desired location 2092 p = camera.GetPosition() 2093 p_new = np.array(p) + pos 2094 camera.SetPosition(*p_new) 2095 camera.SetFocalPoint(*pos) 2096 return camera
Define a vtkCamera with a particular orientation.
Arguments:
- pos: (np.array, list, tuple) an iterator of length 3 containing the focus point of the camera
- quaternion: (np.array, list, tuple)
a len(4) quaternion
(x,y,z,w)describing the rotation of the camera such as returned by neuroglancerx,y,z,wall in[0,1]range - distance: (float) the desired distance from pos to the camera (default = 10000 nm)
Returns:
vtki.vtkCamera, a vtk camera setup according to these rules.
def
camera_from_neuroglancer(state, zoom) -> vtkmodules.vtkRenderingCore.vtkCamera:
2099def camera_from_neuroglancer(state, zoom) -> vtki.vtkCamera: 2100 """ 2101 Define a `vtkCamera` from a neuroglancer state dictionary. 2102 2103 Arguments: 2104 state: (dict) 2105 an neuroglancer state dictionary. 2106 zoom: (float) 2107 how much to multiply zoom by to get camera backoff distance 2108 2109 Returns: 2110 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2111 """ 2112 orient = state.get("perspectiveOrientation", [0.0, 0.0, 0.0, 1.0]) 2113 pzoom = state.get("perspectiveZoom", 10.0) 2114 position = state["navigation"]["pose"]["position"] 2115 pos_nm = np.array(position["voxelCoordinates"]) * position["voxelSize"] 2116 return camera_from_quaternion(pos_nm, orient, pzoom * zoom, ngl_correct=True)
Define a vtkCamera from a neuroglancer state dictionary.
Arguments:
- state: (dict) an neuroglancer state dictionary.
- zoom: (float) how much to multiply zoom by to get camera backoff distance
Returns:
vtki.vtkCamera, a vtk camera setup that matches this state.
def
camera_from_dict(camera, modify_inplace=None) -> vtkmodules.vtkRenderingCore.vtkCamera:
2134def camera_from_dict(camera, modify_inplace=None) -> vtki.vtkCamera: 2135 """ 2136 Generate a `vtkCamera` object from a python dictionary. 2137 2138 Parameters of the camera are: 2139 - `position` or `pos` (3-tuple) 2140 - `focal_point` (3-tuple) 2141 - `viewup` (3-tuple) 2142 - `distance` (float) 2143 - `clipping_range` (2-tuple) 2144 - `parallel_scale` (float) 2145 - `thickness` (float) 2146 - `view_angle` (float) 2147 - `roll` (float) 2148 2149 Exaplanation of the parameters can be found in the 2150 [vtkCamera documentation](https://vtk.org/doc/nightly/html/classvtkCamera.html). 2151 2152 Arguments: 2153 camera: (dict) 2154 a python dictionary containing camera parameters. 2155 modify_inplace: (vtkCamera) 2156 an existing `vtkCamera` object to modify in place. 2157 2158 Returns: 2159 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2160 """ 2161 if modify_inplace: 2162 vcam = modify_inplace 2163 else: 2164 vcam = vtki.vtkCamera() 2165 2166 camera = dict(camera) # make a copy so input is not emptied by pop() 2167 2168 cm_pos = camera.pop("position", camera.pop("pos", None)) 2169 cm_focal_point = camera.pop("focal_point", camera.pop("focalPoint", None)) 2170 cm_viewup = camera.pop("viewup", None) 2171 cm_distance = camera.pop("distance", None) 2172 cm_clipping_range = camera.pop("clipping_range", camera.pop("clippingRange", None)) 2173 cm_parallel_scale = camera.pop("parallel_scale", camera.pop("parallelScale", None)) 2174 cm_thickness = camera.pop("thickness", None) 2175 cm_view_angle = camera.pop("view_angle", camera.pop("viewAngle", None)) 2176 cm_roll = camera.pop("roll", None) 2177 2178 if len(camera.keys()) > 0: 2179 vedo.logger.warning(f"in camera_from_dict, key(s) not recognized: {camera.keys()}") 2180 if cm_pos is not None: vcam.SetPosition(cm_pos) 2181 if cm_focal_point is not None: vcam.SetFocalPoint(cm_focal_point) 2182 if cm_viewup is not None: vcam.SetViewUp(cm_viewup) 2183 if cm_distance is not None: vcam.SetDistance(cm_distance) 2184 if cm_clipping_range is not None: vcam.SetClippingRange(cm_clipping_range) 2185 if cm_parallel_scale is not None: vcam.SetParallelScale(cm_parallel_scale) 2186 if cm_thickness is not None: vcam.SetThickness(cm_thickness) 2187 if cm_view_angle is not None: vcam.SetViewAngle(cm_view_angle) 2188 if cm_roll is not None: vcam.SetRoll(cm_roll) 2189 return vcam
Generate a vtkCamera object from a python dictionary.
Parameters of the camera are:
positionorpos(3-tuple)focal_point(3-tuple)viewup(3-tuple)distance(float)clipping_range(2-tuple)parallel_scale(float)thickness(float)view_angle(float)roll(float)
Exaplanation of the parameters can be found in the vtkCamera documentation.
Arguments:
- camera: (dict) a python dictionary containing camera parameters.
- modify_inplace: (vtkCamera)
an existing
vtkCameraobject to modify in place.
Returns:
vtki.vtkCamera, a vtk camera setup that matches this state.
def
camera_to_dict(vtkcam) -> dict:
2191def camera_to_dict(vtkcam) -> dict: 2192 """ 2193 Convert a [vtkCamera](https://vtk.org/doc/nightly/html/classvtkCamera.html) 2194 object into a python dictionary. 2195 2196 Parameters of the camera are: 2197 - `position` (3-tuple) 2198 - `focal_point` (3-tuple) 2199 - `viewup` (3-tuple) 2200 - `distance` (float) 2201 - `clipping_range` (2-tuple) 2202 - `parallel_scale` (float) 2203 - `thickness` (float) 2204 - `view_angle` (float) 2205 - `roll` (float) 2206 2207 Arguments: 2208 vtkcam: (vtkCamera) 2209 a `vtkCamera` object to convert. 2210 """ 2211 cam = dict() 2212 cam["position"] = np.array(vtkcam.GetPosition()) 2213 cam["focal_point"] = np.array(vtkcam.GetFocalPoint()) 2214 cam["viewup"] = np.array(vtkcam.GetViewUp()) 2215 cam["distance"] = vtkcam.GetDistance() 2216 cam["clipping_range"] = np.array(vtkcam.GetClippingRange()) 2217 cam["parallel_scale"] = vtkcam.GetParallelScale() 2218 cam["thickness"] = vtkcam.GetThickness() 2219 cam["view_angle"] = vtkcam.GetViewAngle() 2220 cam["roll"] = vtkcam.GetRoll() 2221 return cam
Convert a vtkCamera object into a python dictionary.
Parameters of the camera are:
position(3-tuple)focal_point(3-tuple)viewup(3-tuple)distance(float)clipping_range(2-tuple)parallel_scale(float)thickness(float)view_angle(float)roll(float)
Arguments:
- vtkcam: (vtkCamera)
a
vtkCameraobject to convert.
def
oriented_camera( center, up_vector, backoff_vector, backoff) -> vtkmodules.vtkRenderingCore.vtkCamera:
2119def oriented_camera(center, up_vector, backoff_vector, backoff) -> vtki.vtkCamera: 2120 """ 2121 Generate a `vtkCamera` pointed at a specific location, 2122 oriented with a given up direction, set to a backoff. 2123 """ 2124 vup = np.array(up_vector) 2125 vup = vup / np.linalg.norm(vup) 2126 pt_backoff = center - backoff * np.array(backoff_vector) 2127 camera = vtki.vtkCamera() 2128 camera.SetFocalPoint(center[0], center[1], center[2]) 2129 camera.SetViewUp(vup[0], vup[1], vup[2]) 2130 camera.SetPosition(pt_backoff[0], pt_backoff[1], pt_backoff[2]) 2131 return camera
Generate a vtkCamera pointed at a specific location,
oriented with a given up direction, set to a backoff.
def
vedo2trimesh(mesh):
2447def vedo2trimesh(mesh): 2448 """ 2449 Convert `vedo.mesh.Mesh` to `Trimesh.Mesh` object. 2450 """ 2451 if is_sequence(mesh): 2452 tms = [] 2453 for a in mesh: 2454 tms.append(vedo2trimesh(a)) 2455 return tms 2456 2457 try: 2458 from trimesh import Trimesh 2459 except ModuleNotFoundError: 2460 vedo.logger.error("Need trimesh to run:\npip install trimesh") 2461 return None 2462 2463 tris = mesh.cells 2464 carr = mesh.celldata["CellIndividualColors"] 2465 ccols = carr 2466 2467 points = mesh.vertices 2468 varr = mesh.pointdata["VertexColors"] 2469 vcols = varr 2470 2471 if len(tris) == 0: 2472 tris = None 2473 2474 return Trimesh(vertices=points, faces=tris, face_colors=ccols, vertex_colors=vcols, process=False)
Convert vedo.mesh.Mesh to Trimesh.Mesh object.
def
trimesh2vedo(inputobj):
2477def trimesh2vedo(inputobj): 2478 """ 2479 Convert a `Trimesh` object to `vedo.Mesh` or `vedo.Assembly` object. 2480 """ 2481 if is_sequence(inputobj): 2482 vms = [] 2483 for ob in inputobj: 2484 vms.append(trimesh2vedo(ob)) 2485 return vms 2486 2487 inputobj_type = str(type(inputobj)) 2488 2489 if "Trimesh" in inputobj_type or "primitives" in inputobj_type: 2490 faces = inputobj.faces 2491 poly = buildPolyData(inputobj.vertices, faces) 2492 tact = vedo.Mesh(poly) 2493 if inputobj.visual.kind == "face": 2494 trim_c = inputobj.visual.face_colors 2495 elif inputobj.visual.kind == "texture": 2496 trim_c = inputobj.visual.to_color().vertex_colors 2497 else: 2498 trim_c = inputobj.visual.vertex_colors 2499 2500 if is_sequence(trim_c): 2501 if is_sequence(trim_c[0]): 2502 same_color = len(np.unique(trim_c, axis=0)) < 2 # all vtxs have same color 2503 2504 if same_color: 2505 tact.c(trim_c[0, [0, 1, 2]]).alpha(trim_c[0, 3]) 2506 else: 2507 if inputobj.visual.kind == "face": 2508 tact.cellcolors = trim_c 2509 return tact 2510 2511 if "PointCloud" in inputobj_type: 2512 2513 vdpts = vedo.shapes.Points(inputobj.vertices, r=8, c='k') 2514 if hasattr(inputobj, "vertices_color"): 2515 vcols = (inputobj.vertices_color * 1).astype(np.uint8) 2516 vdpts.pointcolors = vcols 2517 return vdpts 2518 2519 if "path" in inputobj_type: 2520 2521 lines = [] 2522 for e in inputobj.entities: 2523 # print('trimesh entity', e.to_dict()) 2524 l = vedo.shapes.Line(inputobj.vertices[e.points], c="k", lw=2) 2525 lines.append(l) 2526 return vedo.Assembly(lines) 2527 2528 return None
Convert a Trimesh object to vedo.Mesh or vedo.Assembly object.
def
vedo2meshlab(vmesh):
2531def vedo2meshlab(vmesh): 2532 """Convert a `vedo.Mesh` to a Meshlab object.""" 2533 try: 2534 import pymeshlab as mlab 2535 except ModuleNotFoundError: 2536 vedo.logger.error("Need pymeshlab to run:\npip install pymeshlab") 2537 2538 vertex_matrix = vmesh.vertices.astype(np.float64) 2539 2540 try: 2541 face_matrix = np.asarray(vmesh.cells, dtype=np.float64) 2542 except: 2543 print("WARNING: in vedo2meshlab(), need to triangulate mesh first!") 2544 face_matrix = np.array(vmesh.clone().triangulate().cells, dtype=np.float64) 2545 2546 # v_normals_matrix = vmesh.normals(cells=False, recompute=False) 2547 v_normals_matrix = vmesh.vertex_normals 2548 if not v_normals_matrix.shape[0]: 2549 v_normals_matrix = np.empty((0, 3), dtype=np.float64) 2550 2551 # f_normals_matrix = vmesh.normals(cells=True, recompute=False) 2552 f_normals_matrix = vmesh.cell_normals 2553 if not f_normals_matrix.shape[0]: 2554 f_normals_matrix = np.empty((0, 3), dtype=np.float64) 2555 2556 v_color_matrix = vmesh.pointdata["RGBA"] 2557 if v_color_matrix is None: 2558 v_color_matrix = np.empty((0, 4), dtype=np.float64) 2559 else: 2560 v_color_matrix = v_color_matrix.astype(np.float64) / 255 2561 if v_color_matrix.shape[1] == 3: 2562 v_color_matrix = np.c_[ 2563 v_color_matrix, np.ones(v_color_matrix.shape[0], dtype=np.float64) 2564 ] 2565 2566 f_color_matrix = vmesh.celldata["RGBA"] 2567 if f_color_matrix is None: 2568 f_color_matrix = np.empty((0, 4), dtype=np.float64) 2569 else: 2570 f_color_matrix = f_color_matrix.astype(np.float64) / 255 2571 if f_color_matrix.shape[1] == 3: 2572 f_color_matrix = np.c_[ 2573 f_color_matrix, np.ones(f_color_matrix.shape[0], dtype=np.float64) 2574 ] 2575 2576 m = mlab.Mesh( 2577 vertex_matrix=vertex_matrix, 2578 face_matrix=face_matrix, 2579 v_normals_matrix=v_normals_matrix, 2580 f_normals_matrix=f_normals_matrix, 2581 v_color_matrix=v_color_matrix, 2582 f_color_matrix=f_color_matrix, 2583 ) 2584 2585 for k in vmesh.pointdata.keys(): 2586 data = vmesh.pointdata[k] 2587 if data is not None: 2588 if data.ndim == 1: # scalar 2589 m.add_vertex_custom_scalar_attribute(data.astype(np.float64), k) 2590 elif data.ndim == 2: # vectorial data 2591 if "tcoord" not in k.lower() and k not in ["Normals", "TextureCoordinates"]: 2592 m.add_vertex_custom_point_attribute(data.astype(np.float64), k) 2593 2594 for k in vmesh.celldata.keys(): 2595 data = vmesh.celldata[k] 2596 if data is not None: 2597 if data.ndim == 1: # scalar 2598 m.add_face_custom_scalar_attribute(data.astype(np.float64), k) 2599 elif data.ndim == 2 and k != "Normals": # vectorial data 2600 m.add_face_custom_point_attribute(data.astype(np.float64), k) 2601 2602 m.update_bounding_box() 2603 return m
Convert a vedo.Mesh to a Meshlab object.
def
meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()):
2606def meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()): 2607 """Convert a Meshlab object to `vedo.Mesh`.""" 2608 inputtype = str(type(mmesh)) 2609 2610 if "MeshSet" in inputtype: 2611 mmesh = mmesh.current_mesh() 2612 2613 mpoints, mcells = mmesh.vertex_matrix(), mmesh.face_matrix() 2614 if len(mcells) > 0: 2615 polydata = buildPolyData(mpoints, mcells) 2616 else: 2617 polydata = buildPolyData(mpoints, None) 2618 2619 if mmesh.has_vertex_scalar(): 2620 parr = mmesh.vertex_scalar_array() 2621 parr_vtk = numpy_to_vtk(parr) 2622 parr_vtk.SetName("MeshLabScalars") 2623 polydata.GetPointData().AddArray(parr_vtk) 2624 polydata.GetPointData().SetActiveScalars("MeshLabScalars") 2625 2626 if mmesh.has_face_scalar(): 2627 carr = mmesh.face_scalar_array() 2628 carr_vtk = numpy_to_vtk(carr) 2629 carr_vtk.SetName("MeshLabScalars") 2630 x0, x1 = carr_vtk.GetRange() 2631 polydata.GetCellData().AddArray(carr_vtk) 2632 polydata.GetCellData().SetActiveScalars("MeshLabScalars") 2633 2634 for k in pointdata_keys: 2635 parr = mmesh.vertex_custom_scalar_attribute_array(k) 2636 parr_vtk = numpy_to_vtk(parr) 2637 parr_vtk.SetName(k) 2638 polydata.GetPointData().AddArray(parr_vtk) 2639 polydata.GetPointData().SetActiveScalars(k) 2640 2641 for k in celldata_keys: 2642 carr = mmesh.face_custom_scalar_attribute_array(k) 2643 carr_vtk = numpy_to_vtk(carr) 2644 carr_vtk.SetName(k) 2645 polydata.GetCellData().AddArray(carr_vtk) 2646 polydata.GetCellData().SetActiveScalars(k) 2647 2648 pnorms = mmesh.vertex_normal_matrix() 2649 if len(pnorms) > 0: 2650 polydata.GetPointData().SetNormals(numpy2vtk(pnorms, name="Normals")) 2651 2652 cnorms = mmesh.face_normal_matrix() 2653 if len(cnorms) > 0: 2654 polydata.GetCellData().SetNormals(numpy2vtk(cnorms, name="Normals")) 2655 return vedo.Mesh(polydata)
Convert a Meshlab object to vedo.Mesh.
def
vedo2open3d(vedo_mesh):
2666def vedo2open3d(vedo_mesh): 2667 """ 2668 Return an `open3d.geometry.TriangleMesh` version of the current mesh. 2669 """ 2670 try: 2671 import open3d as o3d 2672 except RuntimeError: 2673 vedo.logger.error("Need open3d to run:\npip install open3d") 2674 2675 # create from numpy arrays 2676 o3d_mesh = o3d.geometry.TriangleMesh( 2677 vertices=o3d.utility.Vector3dVector(vedo_mesh.vertices), 2678 triangles=o3d.utility.Vector3iVector(vedo_mesh.cells), 2679 ) 2680 # TODO: need to add some if check here in case color and normals 2681 # info are not existing 2682 # o3d_mesh.vertex_colors = o3d.utility.Vector3dVector(vedo_mesh.pointdata["RGB"]/255) 2683 # o3d_mesh.vertex_normals= o3d.utility.Vector3dVector(vedo_mesh.pointdata["Normals"]) 2684 return o3d_mesh
Return an open3d.geometry.TriangleMesh version of the current mesh.
def
open3d2vedo(o3d_mesh):
2658def open3d2vedo(o3d_mesh): 2659 """Convert `open3d.geometry.TriangleMesh` to a `vedo.Mesh`.""" 2660 m = vedo.Mesh([np.array(o3d_mesh.vertices), np.array(o3d_mesh.triangles)]) 2661 # TODO: could also check whether normals and color are present in 2662 # order to port with the above vertices/faces 2663 return m
Convert open3d.geometry.TriangleMesh to a vedo.Mesh.
def
vtk2numpy(varr):
796def vtk2numpy(varr): 797 """Convert a `vtkDataArray`, `vtkIdList` or `vtTransform` into a numpy array.""" 798 if varr is None: 799 return np.array([]) 800 if isinstance(varr, vtki.vtkIdList): 801 return np.array([varr.GetId(i) for i in range(varr.GetNumberOfIds())]) 802 elif isinstance(varr, vtki.vtkBitArray): 803 carr = vtki.vtkCharArray() 804 carr.DeepCopy(varr) 805 varr = carr 806 elif isinstance(varr, vtki.vtkHomogeneousTransform): 807 try: 808 varr = varr.GetMatrix() 809 except AttributeError: 810 pass 811 M = [[varr.GetElement(i, j) for j in range(4)] for i in range(4)] 812 return np.array(M) 813 return vtk_to_numpy(varr)
Convert a vtkDataArray, vtkIdList or vtTransform into a numpy array.
def
numpy2vtk(arr, dtype=None, deep=True, name=''):
767def numpy2vtk(arr, dtype=None, deep=True, name=""): 768 """ 769 Convert a numpy array into a `vtkDataArray`. 770 Use `dtype='id'` for `vtkIdTypeArray` objects. 771 """ 772 # https://github.com/Kitware/VTK/blob/master/Wrapping/Python/vtkmodules/util/numpy_support.py 773 if arr is None: 774 return None 775 776 arr = np.ascontiguousarray(arr) 777 778 if dtype == "id": 779 if vtki.vtkIdTypeArray().GetDataTypeSize() != 4: 780 ast = np.int64 781 else: 782 ast = np.int32 783 varr = numpy_to_vtkIdTypeArray(arr.astype(ast), deep=deep) 784 elif dtype: 785 varr = numpy_to_vtk(arr.astype(dtype), deep=deep) 786 else: 787 # let numpy_to_vtk() decide what is best type based on arr type 788 if arr.dtype == np.bool_: 789 arr = arr.astype(np.uint8) 790 varr = numpy_to_vtk(arr, deep=deep) 791 792 if name: 793 varr.SetName(name) 794 return varr
Convert a numpy array into a vtkDataArray.
Use dtype='id' for vtkIdTypeArray objects.
def
get_uv(p, x, v):
1457def get_uv(p, x, v): 1458 """ 1459 Obtain the texture uv-coords of a point p belonging to a face that has point 1460 coordinates (x0, x1, x2) with the corresponding uv-coordinates v=(v0, v1, v2). 1461 All p and x0,x1,x2 are 3D-vectors, while v are their 2D uv-coordinates. 1462 1463 Example: 1464 ```python 1465 from vedo import * 1466 1467 pic = Image(dataurl+"coloured_cube_faces.jpg") 1468 cb = Mesh(dataurl+"coloured_cube.obj").lighting("off").texture(pic) 1469 1470 cbpts = cb.vertices 1471 faces = cb.cells 1472 uv = cb.pointdata["Material"] 1473 1474 pt = [-0.2, 0.75, 2] 1475 pr = cb.closest_point(pt) 1476 1477 idface = cb.closest_point(pt, return_cell_id=True) 1478 idpts = faces[idface] 1479 uv_face = uv[idpts] 1480 1481 uv_pr = utils.get_uv(pr, cbpts[idpts], uv_face) 1482 print("interpolated uv =", uv_pr) 1483 1484 sx, sy = pic.dimensions() 1485 i_interp_uv = uv_pr * [sy, sx] 1486 ix, iy = i_interp_uv.astype(int) 1487 mpic = pic.tomesh() 1488 rgba = mpic.pointdata["RGBA"].reshape(sy, sx, 3) 1489 print("color =", rgba[ix, iy]) 1490 1491 show( 1492 [[cb, Point(pr), cb.labels("Material")], 1493 [pic, Point(i_interp_uv)]], 1494 N=2, axes=1, sharecam=False, 1495 ).close() 1496 ``` 1497  1498 """ 1499 # Vector vp=p-x0 is representable as alpha*s + beta*t, 1500 # where s = x1-x0 and t = x2-x0, in matrix form 1501 # vp = [alpha, beta] . matrix(s,t) 1502 # M = matrix(s,t) is 2x3 matrix, so (alpha, beta) can be found by 1503 # inverting any of its minor A with non-zero determinant. 1504 # Once found, uv-coords of p are vt0 + alpha (vt1-v0) + beta (vt2-v0) 1505 1506 p = np.asarray(p) 1507 x0, x1, x2 = np.asarray(x)[:3] 1508 vt0, vt1, vt2 = np.asarray(v)[:3] 1509 1510 s = x1 - x0 1511 t = x2 - x0 1512 vs = vt1 - vt0 1513 vt = vt2 - vt0 1514 vp = p - x0 1515 1516 # finding a minor with independent rows 1517 M = np.matrix([s, t]) 1518 mnr = [0, 1] 1519 A = M[:, mnr] 1520 if np.abs(np.linalg.det(A)) < 0.000001: 1521 mnr = [0, 2] 1522 A = M[:, mnr] 1523 if np.abs(np.linalg.det(A)) < 0.000001: 1524 mnr = [1, 2] 1525 A = M[:, mnr] 1526 Ainv = np.linalg.inv(A) 1527 alpha_beta = vp[mnr].dot(Ainv) # [alpha, beta] 1528 return np.asarray(vt0 + alpha_beta.dot(np.matrix([vs, vt])))[0]
Obtain the texture uv-coords of a point p belonging to a face that has point coordinates (x0, x1, x2) with the corresponding uv-coordinates v=(v0, v1, v2). All p and x0,x1,x2 are 3D-vectors, while v are their 2D uv-coordinates.
Example:
from vedo import * pic = Image(dataurl+"coloured_cube_faces.jpg") cb = Mesh(dataurl+"coloured_cube.obj").lighting("off").texture(pic) cbpts = cb.vertices faces = cb.cells uv = cb.pointdata["Material"] pt = [-0.2, 0.75, 2] pr = cb.closest_point(pt) idface = cb.closest_point(pt, return_cell_id=True) idpts = faces[idface] uv_face = uv[idpts] uv_pr = utils.get_uv(pr, cbpts[idpts], uv_face) print("interpolated uv =", uv_pr) sx, sy = pic.dimensions() i_interp_uv = uv_pr * [sy, sx] ix, iy = i_interp_uv.astype(int) mpic = pic.tomesh() rgba = mpic.pointdata["RGBA"].reshape(sy, sx, 3) print("color =", rgba[ix, iy]) show( [[cb, Point(pr), cb.labels("Material")], [pic, Point(i_interp_uv)]], N=2, axes=1, sharecam=False, ).close()
def
andrews_curves(M, res=100) -> numpy.ndarray:
712def andrews_curves(M, res=100) -> np.ndarray: 713 """ 714 Computes the [Andrews curves](https://en.wikipedia.org/wiki/Andrews_plot) 715 for the provided data. 716 717 The input array is an array of shape (n,m) where n is the number of 718 features and m is the number of observations. 719 720 Arguments: 721 M : (ndarray) 722 the data matrix (or data vector). 723 res : (int) 724 the resolution (n. of points) of the output curve. 725 726 Example: 727 - [andrews_cluster.py](https://github.com/marcomusy/vedo/blob/master/examples/pyplot/andrews_cluster.py) 728 729  730 """ 731 # Credits: 732 # https://gist.github.com/ryuzakyl/12c221ff0e54d8b1ac171c69ea552c0a 733 M = np.asarray(M) 734 m = int(res + 0.5) 735 736 # getting data vectors 737 X = np.reshape(M, (1, -1)) if len(M.shape) == 1 else M.copy() 738 _rows, n = X.shape 739 740 # andrews curve dimension (n. theta angles) 741 t = np.linspace(-np.pi, np.pi, m) 742 743 # m: range of values for angle theta 744 # n: amount of components of the Fourier expansion 745 A = np.empty((m, n)) 746 747 # setting first column of A 748 A[:, 0] = [1/np.sqrt(2)] * m 749 750 # filling columns of A 751 for i in range(1, n): 752 # computing the scaling coefficient for angle theta 753 c = np.ceil(i / 2) 754 # computing i-th column of matrix A 755 col = np.sin(c * t) if i % 2 == 1 else np.cos(c * t) 756 # setting column in matrix A 757 A[:, i] = col[:] 758 759 # computing Andrews curves for provided data 760 andrew_curves = np.dot(A, X.T).T 761 762 # returning the Andrews Curves (raveling if needed) 763 return np.ravel(andrew_curves) if andrew_curves.shape[0] == 1 else andrew_curves
Computes the Andrews curves for the provided data.
The input array is an array of shape (n,m) where n is the number of features and m is the number of observations.
Arguments:
- M : (ndarray) the data matrix (or data vector).
- res : (int) the resolution (n. of points) of the output curve.
