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 varr = numpy_to_vtkIdTypeArray(arr.astype(np.int64), deep=deep) 779 elif dtype: 780 varr = numpy_to_vtk(arr.astype(dtype), deep=deep) 781 else: 782 # let numpy_to_vtk() decide what is best type based on arr type 783 if arr.dtype == np.bool_: 784 arr = arr.astype(np.uint8) 785 varr = numpy_to_vtk(arr, deep=deep) 786 787 if name: 788 varr.SetName(name) 789 return varr 790 791def vtk2numpy(varr): 792 """Convert a `vtkDataArray`, `vtkIdList` or `vtTransform` into a numpy array.""" 793 if varr is None: 794 return np.array([]) 795 if isinstance(varr, vtki.vtkIdList): 796 return np.array([varr.GetId(i) for i in range(varr.GetNumberOfIds())]) 797 elif isinstance(varr, vtki.vtkBitArray): 798 carr = vtki.vtkCharArray() 799 carr.DeepCopy(varr) 800 varr = carr 801 elif isinstance(varr, vtki.vtkHomogeneousTransform): 802 try: 803 varr = varr.GetMatrix() 804 except AttributeError: 805 pass 806 M = [[varr.GetElement(i, j) for j in range(4)] for i in range(4)] 807 return np.array(M) 808 return vtk_to_numpy(varr) 809 810 811def make3d(pts) -> np.ndarray: 812 """ 813 Make an array which might be 2D to 3D. 814 815 Array can also be in the form `[allx, ally, allz]`. 816 """ 817 if pts is None: 818 return np.array([]) 819 pts = np.asarray(pts) 820 821 if pts.dtype == "object": 822 raise ValueError("Cannot form a valid numpy array, input may be non-homogenous") 823 824 if pts.size == 0: # empty list 825 return pts 826 827 if pts.ndim == 1: 828 if pts.shape[0] == 2: 829 return np.hstack([pts, [0]]).astype(pts.dtype) 830 elif pts.shape[0] == 3: 831 return pts 832 else: 833 raise ValueError 834 835 if pts.shape[1] == 3: 836 return pts 837 838 # if 2 <= pts.shape[0] <= 3 and pts.shape[1] > 3: 839 # pts = pts.T 840 841 if pts.shape[1] == 2: 842 return np.c_[pts, np.zeros(pts.shape[0], dtype=pts.dtype)] 843 844 if pts.shape[1] != 3: 845 raise ValueError(f"input shape is not supported: {pts.shape}") 846 return pts 847 848 849def geometry(obj, extent=None) -> "vedo.Mesh": 850 """ 851 Apply the `vtkGeometryFilter` to the input object. 852 This is a general-purpose filter to extract geometry (and associated data) 853 from any type of dataset. 854 This filter also may be used to convert any type of data to polygonal type. 855 The conversion process may be less than satisfactory for some 3D datasets. 856 For example, this filter will extract the outer surface of a volume 857 or structured grid dataset. 858 859 Returns a `vedo.Mesh` object. 860 861 Set `extent` as the `[xmin,xmax, ymin,ymax, zmin,zmax]` bounding box to clip data. 862 """ 863 gf = vtki.new("GeometryFilter") 864 gf.SetInputData(obj) 865 if extent is not None: 866 gf.SetExtent(extent) 867 gf.Update() 868 return vedo.Mesh(gf.GetOutput()) 869 870 871def buildPolyData(vertices, faces=None, lines=None, strips=None, index_offset=0) -> vtki.vtkPolyData: 872 """ 873 Build a `vtkPolyData` object from a list of vertices 874 where faces represents the connectivity of the polygonal mesh. 875 Lines and triangle strips can also be specified. 876 877 E.g. : 878 - `vertices=[[x1,y1,z1],[x2,y2,z2], ...]` 879 - `faces=[[0,1,2], [1,2,3], ...]` 880 - `lines=[[0,1], [1,2,3,4], ...]` 881 - `strips=[[0,1,2,3,4,5], [2,3,9,7,4], ...]` 882 883 A flat list of faces can be passed as `faces=[3, 0,1,2, 4, 1,2,3,4, ...]`. 884 For lines use `lines=[2, 0,1, 4, 1,2,3,4, ...]`. 885 886 Use `index_offset=1` if face numbering starts from 1 instead of 0. 887 """ 888 if is_sequence(faces) and len(faces) == 0: 889 faces=None 890 if is_sequence(lines) and len(lines) == 0: 891 lines=None 892 if is_sequence(strips) and len(strips) == 0: 893 strips=None 894 895 poly = vtki.vtkPolyData() 896 897 if len(vertices) == 0: 898 return poly 899 900 vertices = make3d(vertices) 901 source_points = vtki.vtkPoints() 902 source_points.SetData(numpy2vtk(vertices, dtype=np.float32)) 903 poly.SetPoints(source_points) 904 905 if lines is not None: 906 # Create a cell array to store the lines in and add the lines to it 907 linesarr = vtki.vtkCellArray() 908 if is_sequence(lines[0]): # assume format [(id0,id1),..] 909 for iline in lines: 910 for i in range(0, len(iline) - 1): 911 i1, i2 = iline[i], iline[i + 1] 912 if i1 != i2: 913 vline = vtki.vtkLine() 914 vline.GetPointIds().SetId(0, i1) 915 vline.GetPointIds().SetId(1, i2) 916 linesarr.InsertNextCell(vline) 917 else: # assume format [id0,id1,...] 918 # print("buildPolyData: assuming lines format [id0,id1,...]", lines) 919 # TODO CORRECT THIS CASE, MUST BE [2, id0,id1,...] 920 for i in range(0, len(lines) - 1): 921 vline = vtki.vtkLine() 922 vline.GetPointIds().SetId(0, lines[i]) 923 vline.GetPointIds().SetId(1, lines[i + 1]) 924 linesarr.InsertNextCell(vline) 925 poly.SetLines(linesarr) 926 927 if faces is not None: 928 929 source_polygons = vtki.vtkCellArray() 930 931 if isinstance(faces, np.ndarray) or not is_ragged(faces): 932 ##### all faces are composed of equal nr of vtxs, FAST 933 faces = np.asarray(faces) 934 if vtki.vtkIdTypeArray().GetDataTypeSize() != 4: 935 ast = np.int64 936 else: 937 ast = np.int32 938 939 if faces.ndim > 1: 940 nf, nc = faces.shape 941 hs = np.hstack((np.zeros(nf)[:, None] + nc, faces)) 942 else: 943 nf = faces.shape[0] 944 hs = faces 945 arr = numpy_to_vtkIdTypeArray(hs.astype(ast).ravel(), deep=True) 946 source_polygons.SetCells(nf, arr) 947 948 else: 949 ############################# manually add faces, SLOW 950 for f in faces: 951 n = len(f) 952 953 if n == 3: 954 tri = vtki.vtkTriangle() 955 pids = tri.GetPointIds() 956 for i in range(3): 957 pids.SetId(i, f[i] - index_offset) 958 source_polygons.InsertNextCell(tri) 959 960 else: 961 ele = vtki.vtkPolygon() 962 pids = ele.GetPointIds() 963 pids.SetNumberOfIds(n) 964 for i in range(n): 965 pids.SetId(i, f[i] - index_offset) 966 source_polygons.InsertNextCell(ele) 967 968 poly.SetPolys(source_polygons) 969 970 if strips is not None: 971 tscells = vtki.vtkCellArray() 972 for strip in strips: 973 # create a triangle strip 974 # https://vtk.org/doc/nightly/html/classvtkTriangleStrip.html 975 n = len(strip) 976 tstrip = vtki.vtkTriangleStrip() 977 tstrip_ids = tstrip.GetPointIds() 978 tstrip_ids.SetNumberOfIds(n) 979 for i in range(n): 980 tstrip_ids.SetId(i, strip[i] - index_offset) 981 tscells.InsertNextCell(tstrip) 982 poly.SetStrips(tscells) 983 984 if faces is None and lines is None and strips is None: 985 source_vertices = vtki.vtkCellArray() 986 for i in range(len(vertices)): 987 source_vertices.InsertNextCell(1) 988 source_vertices.InsertCellPoint(i) 989 poly.SetVerts(source_vertices) 990 991 # print("buildPolyData \n", 992 # poly.GetNumberOfPoints(), 993 # poly.GetNumberOfCells(), # grand total 994 # poly.GetNumberOfLines(), 995 # poly.GetNumberOfPolys(), 996 # poly.GetNumberOfStrips(), 997 # poly.GetNumberOfVerts(), 998 # ) 999 return poly 1000 1001 1002############################################################################## 1003def get_font_path(font: str) -> str: 1004 """Internal use.""" 1005 if font in vedo.settings.font_parameters.keys(): 1006 if vedo.settings.font_parameters[font]["islocal"]: 1007 fl = os.path.join(vedo.fonts_path, f"{font}.ttf") 1008 else: 1009 try: 1010 fl = vedo.file_io.download(f"https://vedo.embl.es/fonts/{font}.ttf", verbose=False) 1011 except: 1012 vedo.logger.warning(f"Could not download https://vedo.embl.es/fonts/{font}.ttf") 1013 fl = os.path.join(vedo.fonts_path, "Normografo.ttf") 1014 else: 1015 if font.startswith("https://"): 1016 fl = vedo.file_io.download(font, verbose=False) 1017 elif os.path.isfile(font): 1018 fl = font # assume user is passing a valid file 1019 else: 1020 if font.endswith(".ttf"): 1021 vedo.logger.error( 1022 f"Could not set font file {font}" 1023 f"-> using default: {vedo.settings.default_font}" 1024 ) 1025 else: 1026 vedo.settings.default_font = "Normografo" 1027 vedo.logger.error( 1028 f"Could not set font name {font}" 1029 f" -> using default: Normografo\n" 1030 f"Check out https://vedo.embl.es/fonts for additional fonts\n" 1031 f"Type 'vedo -r fonts' to see available fonts" 1032 ) 1033 fl = get_font_path(vedo.settings.default_font) 1034 return fl 1035 1036 1037def is_sequence(arg) -> bool: 1038 """Check if the input is iterable.""" 1039 if hasattr(arg, "strip"): 1040 return False 1041 if hasattr(arg, "__getslice__"): 1042 return True 1043 if hasattr(arg, "__iter__"): 1044 return True 1045 return False 1046 1047 1048def is_ragged(arr, deep=False) -> bool: 1049 """ 1050 A ragged or inhomogeneous array in Python is an array 1051 with arrays of different lengths as its elements. 1052 To check if an array is ragged,we iterate through the elements 1053 and check if their lengths are the same. 1054 1055 Example: 1056 ```python 1057 arr = [[1, 2, 3], [[4, 5], [6], 1], [7, 8, 9]] 1058 print(is_ragged(arr, deep=True)) # output: True 1059 ``` 1060 """ 1061 n = len(arr) 1062 if n == 0: 1063 return False 1064 if is_sequence(arr[0]): 1065 length = len(arr[0]) 1066 for i in range(1, n): 1067 if len(arr[i]) != length or (deep and is_ragged(arr[i])): 1068 return True 1069 return False 1070 return False 1071 1072 1073def flatten(list_to_flatten) -> list: 1074 """Flatten out a list.""" 1075 1076 def _genflatten(lst): 1077 for elem in lst: 1078 if isinstance(elem, (list, tuple)): 1079 for x in flatten(elem): 1080 yield x 1081 else: 1082 yield elem 1083 1084 return list(_genflatten(list_to_flatten)) 1085 1086 1087def humansort(alist) -> list: 1088 """ 1089 Sort in place a given list the way humans expect. 1090 1091 E.g. `['file11', 'file1'] -> ['file1', 'file11']` 1092 1093 .. warning:: input list is modified in-place by this function. 1094 """ 1095 import re 1096 1097 def alphanum_key(s): 1098 # Turn a string into a list of string and number chunks. 1099 # e.g. "z23a" -> ["z", 23, "a"] 1100 def tryint(s): 1101 if s.isdigit(): 1102 return int(s) 1103 return s 1104 1105 return [tryint(c) for c in re.split("([0-9]+)", s)] 1106 1107 alist.sort(key=alphanum_key) 1108 return alist # NB: input list is modified 1109 1110 1111def sort_by_column(arr, nth, invert=False) -> np.ndarray: 1112 """Sort a numpy array by its `n-th` column.""" 1113 arr = np.asarray(arr) 1114 arr = arr[arr[:, nth].argsort()] 1115 if invert: 1116 return np.flip(arr, axis=0) 1117 return arr 1118 1119 1120def point_in_triangle(p, p1, p2, p3) -> Union[bool, None]: 1121 """ 1122 Return True if a point is inside (or above/below) 1123 a triangle defined by 3 points in space. 1124 """ 1125 p1 = np.array(p1) 1126 u = p2 - p1 1127 v = p3 - p1 1128 n = np.cross(u, v) 1129 w = p - p1 1130 ln = np.dot(n, n) 1131 if not ln: 1132 return None # degenerate triangle 1133 gamma = (np.dot(np.cross(u, w), n)) / ln 1134 if 0 < gamma < 1: 1135 beta = (np.dot(np.cross(w, v), n)) / ln 1136 if 0 < beta < 1: 1137 alpha = 1 - gamma - beta 1138 if 0 < alpha < 1: 1139 return True 1140 return False 1141 1142 1143def intersection_ray_triangle(P0, P1, V0, V1, V2) -> Union[bool, None, np.ndarray]: 1144 """ 1145 Fast intersection between a directional ray defined by `P0,P1` 1146 and triangle `V0, V1, V2`. 1147 1148 Returns the intersection point or 1149 - `None` if triangle is degenerate, or ray is parallel to triangle plane. 1150 - `False` if no intersection, or ray direction points away from triangle. 1151 """ 1152 # Credits: http://geomalgorithms.com/a06-_intersect-2.html 1153 # Get triangle edge vectors and plane normal 1154 # todo : this is slow should check 1155 # https://vtk.org/doc/nightly/html/classvtkCell.html 1156 V0 = np.asarray(V0, dtype=float) 1157 P0 = np.asarray(P0, dtype=float) 1158 u = V1 - V0 1159 v = V2 - V0 1160 n = np.cross(u, v) 1161 if not np.abs(v).sum(): # triangle is degenerate 1162 return None # do not deal with this case 1163 1164 rd = P1 - P0 # ray direction vector 1165 w0 = P0 - V0 1166 a = -np.dot(n, w0) 1167 b = np.dot(n, rd) 1168 if not b: # ray is parallel to triangle plane 1169 return None 1170 1171 # Get intersect point of ray with triangle plane 1172 r = a / b 1173 if r < 0.0: # ray goes away from triangle 1174 return False # => no intersect 1175 1176 # Gor a segment, also test if (r > 1.0) => no intersect 1177 I = P0 + r * rd # intersect point of ray and plane 1178 1179 # is I inside T? 1180 uu = np.dot(u, u) 1181 uv = np.dot(u, v) 1182 vv = np.dot(v, v) 1183 w = I - V0 1184 wu = np.dot(w, u) 1185 wv = np.dot(w, v) 1186 D = uv * uv - uu * vv 1187 1188 # Get and test parametric coords 1189 s = (uv * wv - vv * wu) / D 1190 if s < 0.0 or s > 1.0: # I is outside T 1191 return False 1192 t = (uv * wu - uu * wv) / D 1193 if t < 0.0 or (s + t) > 1.0: # I is outside T 1194 return False 1195 return I # I is in T 1196 1197 1198def triangle_solver(**input_dict): 1199 """ 1200 Solve a triangle from any 3 known elements. 1201 (Note that there might be more than one solution or none). 1202 Angles are in radians. 1203 1204 Example: 1205 ```python 1206 print(triangle_solver(a=3, b=4, c=5)) 1207 print(triangle_solver(a=3, ac=0.9273, ab=1.5716)) 1208 print(triangle_solver(a=3, b=4, ab=1.5716)) 1209 print(triangle_solver(b=4, bc=.64, ab=1.5716)) 1210 print(triangle_solver(c=5, ac=.9273, bc=0.6435)) 1211 print(triangle_solver(a=3, c=5, bc=0.6435)) 1212 print(triangle_solver(b=4, c=5, ac=0.927)) 1213 ``` 1214 """ 1215 a = input_dict.get("a") 1216 b = input_dict.get("b") 1217 c = input_dict.get("c") 1218 ab = input_dict.get("ab") 1219 bc = input_dict.get("bc") 1220 ac = input_dict.get("ac") 1221 1222 if ab and bc: 1223 ac = np.pi - bc - ab 1224 elif bc and ac: 1225 ab = np.pi - bc - ac 1226 elif ab and ac: 1227 bc = np.pi - ab - ac 1228 1229 if a is not None and b is not None and c is not None: 1230 ab = np.arccos((a ** 2 + b ** 2 - c ** 2) / (2 * a * b)) 1231 sinab = np.sin(ab) 1232 ac = np.arcsin(a / c * sinab) 1233 bc = np.arcsin(b / c * sinab) 1234 1235 elif a is not None and b is not None and ab is not None: 1236 c = np.sqrt(a ** 2 + b ** 2 - 2 * a * b * np.cos(ab)) 1237 sinab = np.sin(ab) 1238 ac = np.arcsin(a / c * sinab) 1239 bc = np.arcsin(b / c * sinab) 1240 1241 elif a is not None and ac is not None and ab is not None: 1242 h = a * np.sin(ac) 1243 b = h / np.sin(bc) 1244 c = b * np.cos(bc) + a * np.cos(ac) 1245 1246 elif b is not None and bc is not None and ab is not None: 1247 h = b * np.sin(bc) 1248 a = h / np.sin(ac) 1249 c = np.sqrt(a * a + b * b) 1250 1251 elif c is not None and ac is not None and bc is not None: 1252 h = c * np.sin(bc) 1253 b1 = c * np.cos(bc) 1254 b2 = h / np.tan(ab) 1255 b = b1 + b2 1256 a = np.sqrt(b2 * b2 + h * h) 1257 1258 elif a is not None and c is not None and bc is not None: 1259 # double solution 1260 h = c * np.sin(bc) 1261 k = np.sqrt(a * a - h * h) 1262 omega = np.arcsin(k / a) 1263 cosbc = np.cos(bc) 1264 b = c * cosbc - k 1265 phi = np.pi / 2 - bc - omega 1266 ac = phi 1267 ab = np.pi - ac - bc 1268 if k: 1269 b2 = c * cosbc + k 1270 ac2 = phi + 2 * omega 1271 ab2 = np.pi - ac2 - bc 1272 return [ 1273 {"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}, 1274 {"a": a, "b": b2, "c": c, "ab": ab2, "bc": bc, "ac": ac2}, 1275 ] 1276 1277 elif b is not None and c is not None and ac is not None: 1278 # double solution 1279 h = c * np.sin(ac) 1280 k = np.sqrt(b * b - h * h) 1281 omega = np.arcsin(k / b) 1282 cosac = np.cos(ac) 1283 a = c * cosac - k 1284 phi = np.pi / 2 - ac - omega 1285 bc = phi 1286 ab = np.pi - bc - ac 1287 if k: 1288 a2 = c * cosac + k 1289 bc2 = phi + 2 * omega 1290 ab2 = np.pi - ac - bc2 1291 return [ 1292 {"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}, 1293 {"a": a2, "b": b, "c": c, "ab": ab2, "bc": bc2, "ac": ac}, 1294 ] 1295 1296 else: 1297 vedo.logger.error(f"Case {input_dict} is not supported.") 1298 return [] 1299 1300 return [{"a": a, "b": b, "c": c, "ab": ab, "bc": bc, "ac": ac}] 1301 1302 1303############################################################################# 1304def circle_from_3points(p1, p2, p3) -> np.ndarray: 1305 """ 1306 Find the center and radius of a circle given 3 points in 3D space. 1307 1308 Returns the center of the circle. 1309 1310 Example: 1311 ```python 1312 from vedo.utils import mag, circle_from_3points 1313 p1 = [0,1,1] 1314 p2 = [3,0,1] 1315 p3 = [1,2,0] 1316 c = circle_from_3points(p1, p2, p3) 1317 print(mag(c-p1), mag(c-p2), mag(c-p3)) 1318 ``` 1319 """ 1320 p1 = np.asarray(p1) 1321 p2 = np.asarray(p2) 1322 p3 = np.asarray(p3) 1323 v1 = p2 - p1 1324 v2 = p3 - p1 1325 v11 = np.dot(v1, v1) 1326 v22 = np.dot(v2, v2) 1327 v12 = np.dot(v1, v2) 1328 f = 1.0 / (2 * (v11 * v22 - v12 * v12)) 1329 k1 = f * v22 * (v11-v12) 1330 k2 = f * v11 * (v22-v12) 1331 return p1 + k1 * v1 + k2 * v2 1332 1333def point_line_distance(p, p1, p2) -> float: 1334 """ 1335 Compute the distance of a point to a line (not the segment) 1336 defined by `p1` and `p2`. 1337 """ 1338 return np.sqrt(vtki.vtkLine.DistanceToLine(p, p1, p2)) 1339 1340def line_line_distance(p1, p2, q1, q2) -> Tuple[float, np.ndarray, np.ndarray, float, float]: 1341 """ 1342 Compute the distance of a line to a line (not the segment) 1343 defined by `p1` and `p2` and `q1` and `q2`. 1344 1345 Returns the distance, 1346 the closest point on line 1, the closest point on line 2. 1347 Their parametric coords (-inf <= t0, t1 <= inf) are also returned. 1348 """ 1349 closest_pt1: MutableSequence[float] = [0,0,0] 1350 closest_pt2: MutableSequence[float] = [0,0,0] 1351 t1, t2 = 0.0, 0.0 1352 d = vtki.vtkLine.DistanceBetweenLines( 1353 p1, p2, q1, q2, closest_pt1, closest_pt2, t1, t2) 1354 return np.sqrt(d), closest_pt1, closest_pt2, t1, t2 1355 1356def segment_segment_distance(p1, p2, q1, q2): 1357 """ 1358 Compute the distance of a segment to a segment 1359 defined by `p1` and `p2` and `q1` and `q2`. 1360 1361 Returns the distance, 1362 the closest point on line 1, the closest point on line 2. 1363 Their parametric coords (-inf <= t0, t1 <= inf) are also returned. 1364 """ 1365 closest_pt1 = [0,0,0] 1366 closest_pt2 = [0,0,0] 1367 t1, t2 = 0.0, 0.0 1368 d = vtki.vtkLine.DistanceBetweenLineSegments( 1369 p1, p2, q1, q2, closest_pt1, closest_pt2, t1, t2) 1370 return np.sqrt(d), closest_pt1, closest_pt2, t1, t2 1371 1372 1373def closest(point, points, n=1, return_ids=False, use_tree=False): 1374 """ 1375 Returns the distances and the closest point(s) to the given set of points. 1376 Needs `scipy.spatial` library. 1377 1378 Arguments: 1379 n : (int) 1380 the nr of closest points to return 1381 return_ids : (bool) 1382 return the ids instead of the points coordinates 1383 use_tree : (bool) 1384 build a `scipy.spatial.KDTree`. 1385 An already existing one can be passed to avoid rebuilding. 1386 """ 1387 from scipy.spatial import distance, KDTree 1388 1389 points = np.asarray(points) 1390 if n == 1: 1391 dists = distance.cdist([point], points) 1392 closest_idx = np.argmin(dists) 1393 else: 1394 if use_tree: 1395 if isinstance(use_tree, KDTree): # reuse 1396 tree = use_tree 1397 else: 1398 tree = KDTree(points) 1399 dists, closest_idx = tree.query([point], k=n) 1400 closest_idx = closest_idx[0] 1401 else: 1402 dists = distance.cdist([point], points) 1403 closest_idx = np.argsort(dists)[0][:n] 1404 if return_ids: 1405 return dists, closest_idx 1406 else: 1407 return dists, points[closest_idx] 1408 1409 1410############################################################################# 1411def lin_interpolate(x, rangeX, rangeY): 1412 """ 1413 Interpolate linearly the variable `x` in `rangeX` onto the new `rangeY`. 1414 If `x` is a 3D vector the linear weight is the distance to the two 3D `rangeX` vectors. 1415 1416 E.g. if `x` runs in `rangeX=[x0,x1]` and I want it to run in `rangeY=[y0,y1]` then 1417 1418 `y = lin_interpolate(x, rangeX, rangeY)` will interpolate `x` onto `rangeY`. 1419 1420 Examples: 1421 - [lin_interpolate.py](https://github.com/marcomusy/vedo/tree/master/examples/basic/lin_interpolate.py) 1422 1423  1424 """ 1425 if is_sequence(x): 1426 x = np.asarray(x) 1427 x0, x1 = np.asarray(rangeX) 1428 y0, y1 = np.asarray(rangeY) 1429 dx = x1 - x0 1430 dxn = np.linalg.norm(dx) 1431 if not dxn: 1432 return y0 1433 s = np.linalg.norm(x - x0) / dxn 1434 t = np.linalg.norm(x - x1) / dxn 1435 st = s + t 1436 out = y0 * (t / st) + y1 * (s / st) 1437 1438 else: # faster 1439 1440 x0 = rangeX[0] 1441 dx = rangeX[1] - x0 1442 if not dx: 1443 return rangeY[0] 1444 s = (x - x0) / dx 1445 out = rangeY[0] * (1 - s) + rangeY[1] * s 1446 return out 1447 1448 1449def get_uv(p, x, v): 1450 """ 1451 Obtain the texture uv-coords of a point p belonging to a face that has point 1452 coordinates (x0, x1, x2) with the corresponding uv-coordinates v=(v0, v1, v2). 1453 All p and x0,x1,x2 are 3D-vectors, while v are their 2D uv-coordinates. 1454 1455 Example: 1456 ```python 1457 from vedo import * 1458 1459 pic = Image(dataurl+"coloured_cube_faces.jpg") 1460 cb = Mesh(dataurl+"coloured_cube.obj").lighting("off").texture(pic) 1461 1462 cbpts = cb.vertices 1463 faces = cb.cells 1464 uv = cb.pointdata["Material"] 1465 1466 pt = [-0.2, 0.75, 2] 1467 pr = cb.closest_point(pt) 1468 1469 idface = cb.closest_point(pt, return_cell_id=True) 1470 idpts = faces[idface] 1471 uv_face = uv[idpts] 1472 1473 uv_pr = utils.get_uv(pr, cbpts[idpts], uv_face) 1474 print("interpolated uv =", uv_pr) 1475 1476 sx, sy = pic.dimensions() 1477 i_interp_uv = uv_pr * [sy, sx] 1478 ix, iy = i_interp_uv.astype(int) 1479 mpic = pic.tomesh() 1480 rgba = mpic.pointdata["RGBA"].reshape(sy, sx, 3) 1481 print("color =", rgba[ix, iy]) 1482 1483 show( 1484 [[cb, Point(pr), cb.labels("Material")], 1485 [pic, Point(i_interp_uv)]], 1486 N=2, axes=1, sharecam=False, 1487 ).close() 1488 ``` 1489  1490 """ 1491 # Vector vp=p-x0 is representable as alpha*s + beta*t, 1492 # where s = x1-x0 and t = x2-x0, in matrix form 1493 # vp = [alpha, beta] . matrix(s,t) 1494 # M = matrix(s,t) is 2x3 matrix, so (alpha, beta) can be found by 1495 # inverting any of its minor A with non-zero determinant. 1496 # Once found, uv-coords of p are vt0 + alpha (vt1-v0) + beta (vt2-v0) 1497 1498 p = np.asarray(p) 1499 x0, x1, x2 = np.asarray(x)[:3] 1500 vt0, vt1, vt2 = np.asarray(v)[:3] 1501 1502 s = x1 - x0 1503 t = x2 - x0 1504 vs = vt1 - vt0 1505 vt = vt2 - vt0 1506 vp = p - x0 1507 1508 # finding a minor with independent rows 1509 M = np.matrix([s, t]) 1510 mnr = [0, 1] 1511 A = M[:, mnr] 1512 if np.abs(np.linalg.det(A)) < 0.000001: 1513 mnr = [0, 2] 1514 A = M[:, mnr] 1515 if np.abs(np.linalg.det(A)) < 0.000001: 1516 mnr = [1, 2] 1517 A = M[:, mnr] 1518 Ainv = np.linalg.inv(A) 1519 alpha_beta = vp[mnr].dot(Ainv) # [alpha, beta] 1520 return np.asarray(vt0 + alpha_beta.dot(np.matrix([vs, vt])))[0] 1521 1522 1523def vector(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1524 """ 1525 Return a 3D numpy array representing a vector. 1526 1527 If `y` is `None`, assume input is already in the form `[x,y,z]`. 1528 """ 1529 if y is None: # assume x is already [x,y,z] 1530 return np.asarray(x, dtype=dtype) 1531 return np.array([x, y, z], dtype=dtype) 1532 1533 1534def versor(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1535 """Return the unit vector. Input can be a list of vectors.""" 1536 v = vector(x, y, z, dtype) 1537 if isinstance(v[0], np.ndarray): 1538 return np.divide(v, mag(v)[:, None]) 1539 return v / mag(v) 1540 1541 1542def mag(v): 1543 """Get the magnitude of a vector or array of vectors.""" 1544 v = np.asarray(v) 1545 if v.ndim == 1: 1546 return np.linalg.norm(v) 1547 return np.linalg.norm(v, axis=1) 1548 1549 1550def mag2(v) -> np.ndarray: 1551 """Get the squared magnitude of a vector or array of vectors.""" 1552 v = np.asarray(v) 1553 if v.ndim == 1: 1554 return np.square(v).sum() 1555 return np.square(v).sum(axis=1) 1556 1557 1558def is_integer(n) -> bool: 1559 """Check if input is an integer.""" 1560 try: 1561 float(n) 1562 except (ValueError, TypeError): 1563 return False 1564 else: 1565 return float(n).is_integer() 1566 1567 1568def is_number(n) -> bool: 1569 """Check if input is a number""" 1570 try: 1571 float(n) 1572 return True 1573 except (ValueError, TypeError): 1574 return False 1575 1576 1577def round_to_digit(x, p) -> float: 1578 """Round a real number to the specified number of significant digits.""" 1579 if not x: 1580 return 0 1581 r = np.round(x, p - int(np.floor(np.log10(abs(x)))) - 1) 1582 if int(r) == r: 1583 return int(r) 1584 return r 1585 1586 1587def pack_spheres(bounds, radius) -> np.ndarray: 1588 """ 1589 Packing spheres into a bounding box. 1590 Returns a numpy array of sphere centers. 1591 """ 1592 h = 0.8164965 / 2 1593 d = 0.8660254 1594 a = 0.288675135 1595 1596 if is_sequence(bounds): 1597 x0, x1, y0, y1, z0, z1 = bounds 1598 else: 1599 x0, x1, y0, y1, z0, z1 = bounds.bounds() 1600 1601 x = np.arange(x0, x1, radius) 1602 nul = np.zeros_like(x) 1603 nz = int((z1 - z0) / radius / h / 2 + 1.5) 1604 ny = int((y1 - y0) / radius / d + 1.5) 1605 1606 pts = [] 1607 for iz in range(nz): 1608 z = z0 + nul + iz * h * radius 1609 dx, dy, dz = [radius * 0.5, radius * a, iz * h * radius] 1610 for iy in range(ny): 1611 y = y0 + nul + iy * d * radius 1612 if iy % 2: 1613 xs = x 1614 else: 1615 xs = x + radius * 0.5 1616 if iz % 2: 1617 p = np.c_[xs, y, z] + [dx, dy, dz] 1618 else: 1619 p = np.c_[xs, y, z] + [0, 0, dz] 1620 pts += p.tolist() 1621 return np.array(pts) 1622 1623 1624def precision(x, p: int, vrange=None, delimiter="e") -> str: 1625 """ 1626 Returns a string representation of `x` formatted to precision `p`. 1627 1628 Set `vrange` to the range in which x exists (to snap x to '0' if below precision). 1629 """ 1630 # Based on the webkit javascript implementation 1631 # `from here <https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp>`_, 1632 # and implemented by `randlet <https://github.com/randlet/to-precision>`_. 1633 # Modified for vedo by M.Musy 2020 1634 1635 if isinstance(x, str): # do nothing 1636 return x 1637 1638 if is_sequence(x): 1639 out = "(" 1640 nn = len(x) - 1 1641 for i, ix in enumerate(x): 1642 1643 try: 1644 if np.isnan(ix): 1645 return "NaN" 1646 except: 1647 # cannot handle list of list 1648 continue 1649 1650 out += precision(ix, p) 1651 if i < nn: 1652 out += ", " 1653 return out + ")" ############ <-- 1654 1655 try: 1656 if np.isnan(x): 1657 return "NaN" 1658 except TypeError: 1659 return "NaN" 1660 1661 x = float(x) 1662 1663 if x == 0.0 or (vrange is not None and abs(x) < vrange / pow(10, p)): 1664 return "0" 1665 1666 out = [] 1667 if x < 0: 1668 out.append("-") 1669 x = -x 1670 1671 e = int(np.log10(x)) 1672 # tens = np.power(10, e - p + 1) 1673 tens = 10 ** (e - p + 1) 1674 n = np.floor(x / tens) 1675 1676 # if n < np.power(10, p - 1): 1677 if n < 10 ** (p - 1): 1678 e = e - 1 1679 # tens = np.power(10, e - p + 1) 1680 tens = 10 ** (e - p + 1) 1681 n = np.floor(x / tens) 1682 1683 if abs((n + 1.0) * tens - x) <= abs(n * tens - x): 1684 n = n + 1 1685 1686 # if n >= np.power(10, p): 1687 if n >= 10 ** p: 1688 n = n / 10.0 1689 e = e + 1 1690 1691 m = "%.*g" % (p, n) 1692 if e < -2 or e >= p: 1693 out.append(m[0]) 1694 if p > 1: 1695 out.append(".") 1696 out.extend(m[1:p]) 1697 out.append(delimiter) 1698 if e > 0: 1699 out.append("+") 1700 out.append(str(e)) 1701 elif e == (p - 1): 1702 out.append(m) 1703 elif e >= 0: 1704 out.append(m[: e + 1]) 1705 if e + 1 < len(m): 1706 out.append(".") 1707 out.extend(m[e + 1 :]) 1708 else: 1709 out.append("0.") 1710 out.extend(["0"] * -(e + 1)) 1711 out.append(m) 1712 return "".join(out) 1713 1714 1715################################################################################## 1716def grep(filename: str, tag: str, column=None, first_occurrence_only=False) -> list: 1717 """Greps the line in a file that starts with a specific `tag` string inside the file.""" 1718 import re 1719 1720 with open(filename, "r", encoding="UTF-8") as afile: 1721 content = [] 1722 for line in afile: 1723 if re.search(tag, line): 1724 c = line.split() 1725 c[-1] = c[-1].replace("\n", "") 1726 if column is not None: 1727 c = c[column] 1728 content.append(c) 1729 if first_occurrence_only: 1730 break 1731 return content 1732 1733def parse_pattern(query, strings_to_parse) -> list: 1734 """ 1735 Parse a pattern query to a list of strings. 1736 The query string can contain wildcards like * and ?. 1737 1738 Arguments: 1739 query : (str) 1740 the query to parse 1741 strings_to_parse : (str/list) 1742 the string or list of strings to parse 1743 1744 Returns: 1745 a list of booleans, one for each string in strings_to_parse 1746 1747 Example: 1748 >>> query = r'*Sphere 1?3*' 1749 >>> strings = ["Sphere 143 red", "Sphere 13 red", "Sphere 123", "ASphere 173"] 1750 >>> parse_pattern(query, strings) 1751 [True, True, False, False] 1752 """ 1753 from re import findall as re_findall 1754 if not isinstance(query, str): 1755 return [False] 1756 1757 if not is_sequence(strings_to_parse): 1758 strings_to_parse = [strings_to_parse] 1759 1760 outs = [] 1761 for sp in strings_to_parse: 1762 if not isinstance(sp, str): 1763 outs.append(False) 1764 continue 1765 1766 s = query 1767 if s.startswith("*"): 1768 s = s[1:] 1769 else: 1770 s = "^" + s 1771 1772 t = "" 1773 if not s.endswith("*"): 1774 t = "$" 1775 else: 1776 s = s[:-1] 1777 1778 pattern = s.replace('?', r'\w').replace(' ', r'\s').replace("*", r"\w+") + t 1779 1780 # Search for the pattern in the input string 1781 match = re_findall(pattern, sp) 1782 out = bool(match) 1783 outs.append(out) 1784 # Print the matches for debugging 1785 print("pattern", pattern, "in:", strings_to_parse) 1786 print("matches", match, "result:", out) 1787 return outs 1788 1789def print_histogram( 1790 data, 1791 bins=10, 1792 height=10, 1793 logscale=False, 1794 minbin=0, 1795 horizontal=True, 1796 char="\U00002589", 1797 c=None, 1798 bold=True, 1799 title="histogram", 1800 spacer="", 1801) -> np.ndarray: 1802 """ 1803 Ascii histogram printing. 1804 1805 Input can be a `vedo.Volume` or `vedo.Mesh`. 1806 Returns the raw data before binning (useful when passing vtk objects). 1807 1808 Arguments: 1809 bins : (int) 1810 number of histogram bins 1811 height : (int) 1812 height of the histogram in character units 1813 logscale : (bool) 1814 use logscale for frequencies 1815 minbin : (int) 1816 ignore bins before minbin 1817 horizontal : (bool) 1818 show histogram horizontally 1819 char : (str) 1820 character to be used 1821 bold : (bool) 1822 use boldface 1823 title : (str) 1824 histogram title 1825 spacer : (str) 1826 horizontal spacer 1827 1828 Example: 1829 ```python 1830 from vedo import print_histogram 1831 import numpy as np 1832 d = np.random.normal(size=1000) 1833 data = print_histogram(d, c='b', logscale=True, title='my scalars') 1834 data = print_histogram(d, c='o') 1835 print(np.mean(data)) # data here is same as d 1836 ``` 1837  1838 """ 1839 # credits: http://pyinsci.blogspot.com/2009/10/ascii-histograms.html 1840 # adapted for vedo by M.Musy, 2019 1841 1842 if not horizontal: # better aspect ratio 1843 bins *= 2 1844 1845 try: 1846 data = vtk2numpy(data.dataset.GetPointData().GetScalars()) 1847 except AttributeError: 1848 # already an array 1849 data = np.asarray(data) 1850 1851 if isinstance(data, vtki.vtkImageData): 1852 dims = data.GetDimensions() 1853 nvx = min(100000, dims[0] * dims[1] * dims[2]) 1854 idxs = np.random.randint(0, min(dims), size=(nvx, 3)) 1855 data = [] 1856 for ix, iy, iz in idxs: 1857 d = data.GetScalarComponentAsFloat(ix, iy, iz, 0) 1858 data.append(d) 1859 data = np.array(data) 1860 1861 elif isinstance(data, vtki.vtkPolyData): 1862 arr = data.GetPointData().GetScalars() 1863 if not arr: 1864 arr = data.GetCellData().GetScalars() 1865 if not arr: 1866 return np.array([]) 1867 data = vtk2numpy(arr) 1868 1869 try: 1870 h = np.histogram(data, bins=bins) 1871 except TypeError as e: 1872 vedo.logger.error(f"cannot compute histogram: {e}") 1873 return np.array([]) 1874 1875 if minbin: 1876 hi = h[0][minbin:-1] 1877 else: 1878 hi = h[0] 1879 1880 if char == "\U00002589" and horizontal: 1881 char = "\U00002586" 1882 1883 title = title.ljust(14) + ":" 1884 entrs = " entries=" + str(len(data)) 1885 if logscale: 1886 h0 = np.log10(hi + 1) 1887 maxh0 = int(max(h0) * 100) / 100 1888 title = title + entrs + " (logscale)" 1889 else: 1890 h0 = hi 1891 maxh0 = max(h0) 1892 title = title + entrs 1893 1894 def _v(): 1895 his = "" 1896 if title: 1897 his += title + "\n" 1898 bars = h0 / maxh0 * height 1899 for l in reversed(range(1, height + 1)): 1900 line = "" 1901 if l == height: 1902 line = "%s " % maxh0 1903 else: 1904 line = " |" + " " * (len(str(maxh0)) - 3) 1905 for c in bars: 1906 if c >= np.ceil(l): 1907 line += char 1908 else: 1909 line += " " 1910 line += "\n" 1911 his += line 1912 his += "%.2f" % h[1][0] + "." * (bins) + "%.2f" % h[1][-1] + "\n" 1913 return his 1914 1915 def _h(): 1916 his = "" 1917 if title: 1918 his += title + "\n" 1919 xl = ["%.2f" % n for n in h[1]] 1920 lxl = [len(l) for l in xl] 1921 bars = h0 / maxh0 * height 1922 his += spacer + " " * int(max(bars) + 2 + max(lxl)) + "%s\n" % maxh0 1923 for i, c in enumerate(bars): 1924 line = xl[i] + " " * int(max(lxl) - lxl[i]) + "| " + char * int(c) + "\n" 1925 his += spacer + line 1926 return his 1927 1928 if horizontal: 1929 height *= 2 1930 vedo.printc(_h(), c=c, bold=bold) 1931 else: 1932 vedo.printc(_v(), c=c, bold=bold) 1933 return data 1934 1935 1936def print_table(*columns, headers=None, c="g") -> None: 1937 """ 1938 Print lists as tables. 1939 1940 Example: 1941 ```python 1942 from vedo.utils import print_table 1943 list1 = ["A", "B", "C"] 1944 list2 = [142, 220, 330] 1945 list3 = [True, False, True] 1946 headers = ["First Column", "Second Column", "Third Column"] 1947 print_table(list1, list2, list3, headers=headers) 1948 ``` 1949 1950  1951 """ 1952 # If headers is not provided, use default header names 1953 corner = "─" 1954 if headers is None: 1955 headers = [f"Column {i}" for i in range(1, len(columns) + 1)] 1956 assert len(headers) == len(columns) 1957 1958 # Find the maximum length of the elements in each column and header 1959 max_lens = [max(len(str(x)) for x in column) for column in columns] 1960 max_len_headers = [max(len(str(header)), max_len) for header, max_len in zip(headers, max_lens)] 1961 1962 # Construct the table header 1963 header = ( 1964 "│ " 1965 + " │ ".join(header.ljust(max_len) for header, max_len in zip(headers, max_len_headers)) 1966 + " │" 1967 ) 1968 1969 # Construct the line separator 1970 line1 = "┌" + corner.join("─" * (max_len + 2) for max_len in max_len_headers) + "┐" 1971 line2 = "└" + corner.join("─" * (max_len + 2) for max_len in max_len_headers) + "┘" 1972 1973 # Print the table header 1974 vedo.printc(line1, c=c) 1975 vedo.printc(header, c=c) 1976 vedo.printc(line2, c=c) 1977 1978 # Print the data rows 1979 for row in zip(*columns): 1980 row = ( 1981 "│ " 1982 + " │ ".join(str(col).ljust(max_len) for col, max_len in zip(row, max_len_headers)) 1983 + " │" 1984 ) 1985 vedo.printc(row, bold=False, c=c) 1986 1987 # Print the line separator again to close the table 1988 vedo.printc(line2, c=c) 1989 1990def print_inheritance_tree(C) -> None: 1991 """Prints the inheritance tree of class C.""" 1992 # Adapted from: https://stackoverflow.com/questions/26568976/ 1993 def class_tree(cls): 1994 subc = [class_tree(sub_class) for sub_class in cls.__subclasses__()] 1995 return {cls.__name__: subc} 1996 1997 def print_tree(tree, indent=8, current_ind=0): 1998 for k, v in tree.items(): 1999 if current_ind: 2000 before_dashes = current_ind - indent 2001 m = " " * before_dashes + "└" + "─" * (indent - 1) + " " + k 2002 vedo.printc(m) 2003 else: 2004 vedo.printc(k) 2005 for sub_tree in v: 2006 print_tree(sub_tree, indent=indent, current_ind=current_ind + indent) 2007 2008 if str(C.__class__) != "<class 'type'>": 2009 C = C.__class__ 2010 ct = class_tree(C) 2011 print_tree(ct) 2012 2013 2014def make_bands(inputlist, n): 2015 """ 2016 Group values of a list into bands of equal value, where 2017 `n` is the number of bands, a positive integer > 2. 2018 2019 Returns a binned list of the same length as the input. 2020 """ 2021 if n < 2: 2022 return inputlist 2023 vmin = np.min(inputlist) 2024 vmax = np.max(inputlist) 2025 bb = np.linspace(vmin, vmax, n, endpoint=0) 2026 dr = bb[1] - bb[0] 2027 bb += dr / 2 2028 tol = dr / 2 * 1.001 2029 newlist = [] 2030 for s in inputlist: 2031 for b in bb: 2032 if abs(s - b) < tol: 2033 newlist.append(b) 2034 break 2035 return np.array(newlist) 2036 2037 2038################################################################# 2039# Functions adapted from: 2040# https://github.com/sdorkenw/MeshParty/blob/master/meshparty/trimesh_vtk.py 2041def camera_from_quaternion(pos, quaternion, distance=10000, ngl_correct=True) -> vtki.vtkCamera: 2042 """ 2043 Define a `vtkCamera` with a particular orientation. 2044 2045 Arguments: 2046 pos: (np.array, list, tuple) 2047 an iterator of length 3 containing the focus point of the camera 2048 quaternion: (np.array, list, tuple) 2049 a len(4) quaternion `(x,y,z,w)` describing the rotation of the camera 2050 such as returned by neuroglancer `x,y,z,w` all in `[0,1]` range 2051 distance: (float) 2052 the desired distance from pos to the camera (default = 10000 nm) 2053 2054 Returns: 2055 `vtki.vtkCamera`, a vtk camera setup according to these rules. 2056 """ 2057 camera = vtki.vtkCamera() 2058 # define the quaternion in vtk, note the swapped order 2059 # w,x,y,z instead of x,y,z,w 2060 quat_vtk = vtki.get_class("Quaternion")( 2061 quaternion[3], quaternion[0], quaternion[1], quaternion[2]) 2062 # use this to define a rotation matrix in x,y,z 2063 # right handed units 2064 M = np.zeros((3, 3), dtype=np.float32) 2065 quat_vtk.ToMatrix3x3(M) 2066 # the default camera orientation is y up 2067 up = [0, 1, 0] 2068 # calculate default camera position is backed off in positive z 2069 pos = [0, 0, distance] 2070 2071 # set the camera rototation by applying the rotation matrix 2072 camera.SetViewUp(*np.dot(M, up)) 2073 # set the camera position by applying the rotation matrix 2074 camera.SetPosition(*np.dot(M, pos)) 2075 if ngl_correct: 2076 # neuroglancer has positive y going down 2077 # so apply these azimuth and roll corrections 2078 # to fix orientatins 2079 camera.Azimuth(-180) 2080 camera.Roll(180) 2081 2082 # shift the camera posiiton and focal position 2083 # to be centered on the desired location 2084 p = camera.GetPosition() 2085 p_new = np.array(p) + pos 2086 camera.SetPosition(*p_new) 2087 camera.SetFocalPoint(*pos) 2088 return camera 2089 2090 2091def camera_from_neuroglancer(state, zoom=300) -> vtki.vtkCamera: 2092 """ 2093 Define a `vtkCamera` from a neuroglancer state dictionary. 2094 2095 Arguments: 2096 state: (dict) 2097 an neuroglancer state dictionary. 2098 zoom: (float) 2099 how much to multiply zoom by to get camera backoff distance 2100 default = 300 > ngl_zoom = 1 > 300 nm backoff distance. 2101 2102 Returns: 2103 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2104 """ 2105 orient = state.get("perspectiveOrientation", [0.0, 0.0, 0.0, 1.0]) 2106 pzoom = state.get("perspectiveZoom", 10.0) 2107 position = state["navigation"]["pose"]["position"] 2108 pos_nm = np.array(position["voxelCoordinates"]) * position["voxelSize"] 2109 return camera_from_quaternion(pos_nm, orient, pzoom * zoom, ngl_correct=True) 2110 2111 2112def oriented_camera(center=(0, 0, 0), up_vector=(0, 1, 0), backoff_vector=(0, 0, 1), backoff=1.0) -> vtki.vtkCamera: 2113 """ 2114 Generate a `vtkCamera` pointed at a specific location, 2115 oriented with a given up direction, set to a backoff. 2116 """ 2117 vup = np.array(up_vector) 2118 vup = vup / np.linalg.norm(vup) 2119 pt_backoff = center - backoff * np.array(backoff_vector) 2120 camera = vtki.vtkCamera() 2121 camera.SetFocalPoint(center[0], center[1], center[2]) 2122 camera.SetViewUp(vup[0], vup[1], vup[2]) 2123 camera.SetPosition(pt_backoff[0], pt_backoff[1], pt_backoff[2]) 2124 return camera 2125 2126 2127def camera_from_dict(camera, modify_inplace=None) -> vtki.vtkCamera: 2128 """ 2129 Generate a `vtkCamera` object from a python dictionary. 2130 2131 Parameters of the camera are: 2132 - `position` or `pos` (3-tuple) 2133 - `focal_point` (3-tuple) 2134 - `viewup` (3-tuple) 2135 - `distance` (float) 2136 - `clipping_range` (2-tuple) 2137 - `parallel_scale` (float) 2138 - `thickness` (float) 2139 - `view_angle` (float) 2140 - `roll` (float) 2141 2142 Exaplanation of the parameters can be found in the 2143 [vtkCamera documentation](https://vtk.org/doc/nightly/html/classvtkCamera.html). 2144 2145 Arguments: 2146 camera: (dict) 2147 a python dictionary containing camera parameters. 2148 modify_inplace: (vtkCamera) 2149 an existing `vtkCamera` object to modify in place. 2150 2151 Returns: 2152 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2153 """ 2154 if modify_inplace: 2155 vcam = modify_inplace 2156 else: 2157 vcam = vtki.vtkCamera() 2158 2159 camera = dict(camera) # make a copy so input is not emptied by pop() 2160 2161 cm_pos = camera.pop("position", camera.pop("pos", None)) 2162 cm_focal_point = camera.pop("focal_point", camera.pop("focalPoint", None)) 2163 cm_viewup = camera.pop("viewup", None) 2164 cm_distance = camera.pop("distance", None) 2165 cm_clipping_range = camera.pop("clipping_range", camera.pop("clippingRange", None)) 2166 cm_parallel_scale = camera.pop("parallel_scale", camera.pop("parallelScale", None)) 2167 cm_thickness = camera.pop("thickness", None) 2168 cm_view_angle = camera.pop("view_angle", camera.pop("viewAngle", None)) 2169 cm_roll = camera.pop("roll", None) 2170 2171 if len(camera.keys()) > 0: 2172 vedo.logger.warning(f"in camera_from_dict, key(s) not recognized: {camera.keys()}") 2173 if cm_pos is not None: vcam.SetPosition(cm_pos) 2174 if cm_focal_point is not None: vcam.SetFocalPoint(cm_focal_point) 2175 if cm_viewup is not None: vcam.SetViewUp(cm_viewup) 2176 if cm_distance is not None: vcam.SetDistance(cm_distance) 2177 if cm_clipping_range is not None: vcam.SetClippingRange(cm_clipping_range) 2178 if cm_parallel_scale is not None: vcam.SetParallelScale(cm_parallel_scale) 2179 if cm_thickness is not None: vcam.SetThickness(cm_thickness) 2180 if cm_view_angle is not None: vcam.SetViewAngle(cm_view_angle) 2181 if cm_roll is not None: vcam.SetRoll(cm_roll) 2182 return vcam 2183 2184def camera_to_dict(vtkcam) -> dict: 2185 """ 2186 Convert a [vtkCamera](https://vtk.org/doc/nightly/html/classvtkCamera.html) 2187 object into a python dictionary. 2188 2189 Parameters of the camera are: 2190 - `position` (3-tuple) 2191 - `focal_point` (3-tuple) 2192 - `viewup` (3-tuple) 2193 - `distance` (float) 2194 - `clipping_range` (2-tuple) 2195 - `parallel_scale` (float) 2196 - `thickness` (float) 2197 - `view_angle` (float) 2198 - `roll` (float) 2199 2200 Arguments: 2201 vtkcam: (vtkCamera) 2202 a `vtkCamera` object to convert. 2203 """ 2204 cam = dict() 2205 cam["position"] = np.array(vtkcam.GetPosition()) 2206 cam["focal_point"] = np.array(vtkcam.GetFocalPoint()) 2207 cam["viewup"] = np.array(vtkcam.GetViewUp()) 2208 cam["distance"] = vtkcam.GetDistance() 2209 cam["clipping_range"] = np.array(vtkcam.GetClippingRange()) 2210 cam["parallel_scale"] = vtkcam.GetParallelScale() 2211 cam["thickness"] = vtkcam.GetThickness() 2212 cam["view_angle"] = vtkcam.GetViewAngle() 2213 cam["roll"] = vtkcam.GetRoll() 2214 return cam 2215 2216 2217def vtkCameraToK3D(vtkcam) -> np.ndarray: 2218 """ 2219 Convert a `vtkCamera` object into a 9-element list to be used by the K3D backend. 2220 2221 Output format is: 2222 `[posx,posy,posz, targetx,targety,targetz, upx,upy,upz]`. 2223 """ 2224 cpos = np.array(vtkcam.GetPosition()) 2225 kam = [cpos.tolist()] 2226 kam.append(vtkcam.GetFocalPoint()) 2227 kam.append(vtkcam.GetViewUp()) 2228 return np.array(kam).ravel() 2229 2230 2231def make_ticks( 2232 x0: float, 2233 x1: float, 2234 n=None, 2235 labels=None, 2236 digits=None, 2237 logscale=False, 2238 useformat="", 2239 ) -> Tuple[np.ndarray, List[str]]: 2240 """ 2241 Generate numeric labels for the `[x0, x1]` range. 2242 2243 The format specifier could be expressed in the format: 2244 `:[[fill]align][sign][#][0][width][,][.precision][type]` 2245 2246 where, the options are: 2247 ``` 2248 fill = any character 2249 align = < | > | = | ^ 2250 sign = + | - | " " 2251 width = integer 2252 precision = integer 2253 type = b | c | d | e | E | f | F | g | G | n | o | s | x | X | % 2254 ``` 2255 2256 E.g.: useformat=":.2f" 2257 """ 2258 # Copyright M. Musy, 2021, license: MIT. 2259 # 2260 # useformat eg: ":.2f", check out: 2261 # https://mkaz.blog/code/python-string-format-cookbook/ 2262 # https://www.programiz.com/python-programming/methods/built-in/format 2263 2264 if x1 <= x0: 2265 # vedo.printc("Error in make_ticks(): x0 >= x1", x0,x1, c='r') 2266 return np.array([0.0, 1.0]), ["", ""] 2267 2268 ticks_str, ticks_float = [], [] 2269 baseline = (1, 2, 5, 10, 20, 50) 2270 2271 if logscale: 2272 if x0 <= 0 or x1 <= 0: 2273 vedo.logger.error("make_ticks: zero or negative range with log scale.") 2274 raise RuntimeError 2275 if n is None: 2276 n = int(abs(np.log10(x1) - np.log10(x0))) + 1 2277 x0, x1 = np.log10([x0, x1]) 2278 2279 if not n: 2280 n = 5 2281 2282 if labels is not None: 2283 # user is passing custom labels 2284 2285 ticks_float.append(0) 2286 ticks_str.append("") 2287 for tp, ts in labels: 2288 if tp == x1: 2289 continue 2290 ticks_str.append(str(ts)) 2291 tickn = lin_interpolate(tp, [x0, x1], [0, 1]) 2292 ticks_float.append(tickn) 2293 2294 else: 2295 # ..here comes one of the shortest and most painful pieces of code: 2296 # automatically choose the best natural axis subdivision based on multiples of 1,2,5 2297 dstep = (x1 - x0) / n # desired step size, begin of the nightmare 2298 2299 basestep = pow(10, np.floor(np.log10(dstep))) 2300 steps = np.array([basestep * i for i in baseline]) 2301 idx = (np.abs(steps - dstep)).argmin() 2302 s = steps[idx] # chosen step size 2303 2304 low_bound, up_bound = 0, 0 2305 if x0 < 0: 2306 low_bound = -pow(10, np.ceil(np.log10(-x0))) 2307 if x1 > 0: 2308 up_bound = pow(10, np.ceil(np.log10(x1))) 2309 2310 if low_bound < 0: 2311 if up_bound < 0: 2312 negaxis = np.arange(low_bound, int(up_bound / s) * s) 2313 else: 2314 if -low_bound / s > 1.0e06: 2315 return np.array([0.0, 1.0]), ["", ""] 2316 negaxis = np.arange(low_bound, 0, s) 2317 else: 2318 negaxis = np.array([]) 2319 2320 if up_bound > 0: 2321 if low_bound > 0: 2322 posaxis = np.arange(int(low_bound / s) * s, up_bound, s) 2323 else: 2324 if up_bound / s > 1.0e06: 2325 return np.array([0.0, 1.0]), ["", ""] 2326 posaxis = np.arange(0, up_bound, s) 2327 else: 2328 posaxis = np.array([]) 2329 2330 fulaxis = np.unique(np.clip(np.concatenate([negaxis, posaxis]), x0, x1)) 2331 # end of the nightmare 2332 2333 if useformat: 2334 sf = "{" + f"{useformat}" + "}" 2335 sas = "" 2336 for x in fulaxis: 2337 sas += sf.format(x) + " " 2338 elif digits is None: 2339 np.set_printoptions(suppress=True) # avoid zero precision 2340 sas = str(fulaxis).replace("[", "").replace("]", "") 2341 sas = sas.replace(".e", "e").replace("e+0", "e+").replace("e-0", "e-") 2342 np.set_printoptions(suppress=None) # set back to default 2343 else: 2344 sas = precision(fulaxis, digits, vrange=(x0, x1)) 2345 sas = sas.replace("[", "").replace("]", "").replace(")", "").replace(",", "") 2346 2347 sas2 = [] 2348 for s in sas.split(): 2349 if s.endswith("."): 2350 s = s[:-1] 2351 if s == "-0": 2352 s = "0" 2353 if digits is not None and "e" in s: 2354 s += " " # add space to terminate modifiers 2355 sas2.append(s) 2356 2357 for ts, tp in zip(sas2, fulaxis): 2358 if tp == x1: 2359 continue 2360 tickn = lin_interpolate(tp, [x0, x1], [0, 1]) 2361 ticks_float.append(tickn) 2362 if logscale: 2363 val = np.power(10, tp) 2364 if useformat: 2365 sf = "{" + f"{useformat}" + "}" 2366 ticks_str.append(sf.format(val)) 2367 else: 2368 if val >= 10: 2369 val = int(val + 0.5) 2370 else: 2371 val = round_to_digit(val, 2) 2372 ticks_str.append(str(val)) 2373 else: 2374 ticks_str.append(ts) 2375 2376 ticks_str.append("") 2377 ticks_float.append(1) 2378 return np.array(ticks_float), ticks_str 2379 2380 2381def grid_corners(i: int, nm: list, size: list, margin=0, yflip=True) -> Tuple[np.ndarray, np.ndarray]: 2382 """ 2383 Compute the 2 corners coordinates of the i-th box in a grid of shape n*m. 2384 The top-left square is square number 1. 2385 2386 Arguments: 2387 i : (int) 2388 input index of the desired grid square (to be used in `show(..., at=...)`). 2389 nm : (list) 2390 grid shape as (n,m). 2391 size : (list) 2392 total size of the grid along x and y. 2393 margin : (float) 2394 keep a small margin between boxes. 2395 yflip : (bool) 2396 y-coordinate points downwards 2397 2398 Returns: 2399 Two 2D points representing the bottom-left corner and the top-right corner 2400 of the `i`-nth box in the grid. 2401 2402 Example: 2403 ```python 2404 from vedo import * 2405 acts=[] 2406 n,m = 5,7 2407 for i in range(1, n*m + 1): 2408 c1,c2 = utils.grid_corners(i, [n,m], [1,1], 0.01) 2409 t = Text3D(i, (c1+c2)/2, c='k', s=0.02, justify='center').z(0.01) 2410 r = Rectangle(c1, c2, c=i) 2411 acts += [t,r] 2412 show(acts, axes=1).close() 2413 ``` 2414  2415 """ 2416 i -= 1 2417 n, m = nm 2418 sx, sy = size 2419 dx, dy = sx / n, sy / m 2420 nx = i % n 2421 ny = int((i - nx) / n) 2422 if yflip: 2423 ny = n - ny 2424 c1 = (dx * nx + margin, dy * ny + margin) 2425 c2 = (dx * (nx + 1) - margin, dy * (ny + 1) - margin) 2426 return np.array(c1), np.array(c2) 2427 2428 2429############################################################################ 2430# Trimesh support 2431# 2432# Install trimesh with: 2433# 2434# sudo apt install python3-rtree 2435# pip install rtree shapely 2436# conda install trimesh 2437# 2438# Check the example gallery in: examples/other/trimesh> 2439########################################################################### 2440def vedo2trimesh(mesh): 2441 """ 2442 Convert `vedo.mesh.Mesh` to `Trimesh.Mesh` object. 2443 """ 2444 if is_sequence(mesh): 2445 tms = [] 2446 for a in mesh: 2447 tms.append(vedo2trimesh(a)) 2448 return tms 2449 2450 try: 2451 from trimesh import Trimesh 2452 except ModuleNotFoundError: 2453 vedo.logger.error("Need trimesh to run:\npip install trimesh") 2454 return None 2455 2456 tris = mesh.cells 2457 carr = mesh.celldata["CellIndividualColors"] 2458 ccols = carr 2459 2460 points = mesh.vertices 2461 varr = mesh.pointdata["VertexColors"] 2462 vcols = varr 2463 2464 if len(tris) == 0: 2465 tris = None 2466 2467 return Trimesh(vertices=points, faces=tris, face_colors=ccols, vertex_colors=vcols) 2468 2469 2470def trimesh2vedo(inputobj): 2471 """ 2472 Convert a `Trimesh` object to `vedo.Mesh` or `vedo.Assembly` object. 2473 """ 2474 if is_sequence(inputobj): 2475 vms = [] 2476 for ob in inputobj: 2477 vms.append(trimesh2vedo(ob)) 2478 return vms 2479 2480 inputobj_type = str(type(inputobj)) 2481 2482 if "Trimesh" in inputobj_type or "primitives" in inputobj_type: 2483 faces = inputobj.faces 2484 poly = buildPolyData(inputobj.vertices, faces) 2485 tact = vedo.Mesh(poly) 2486 if inputobj.visual.kind == "face": 2487 trim_c = inputobj.visual.face_colors 2488 elif inputobj.visual.kind == "texture": 2489 trim_c = inputobj.visual.to_color().vertex_colors 2490 else: 2491 trim_c = inputobj.visual.vertex_colors 2492 2493 if is_sequence(trim_c): 2494 if is_sequence(trim_c[0]): 2495 same_color = len(np.unique(trim_c, axis=0)) < 2 # all vtxs have same color 2496 2497 if same_color: 2498 tact.c(trim_c[0, [0, 1, 2]]).alpha(trim_c[0, 3]) 2499 else: 2500 if inputobj.visual.kind == "face": 2501 tact.cellcolors = trim_c 2502 return tact 2503 2504 if "PointCloud" in inputobj_type: 2505 2506 vdpts = vedo.shapes.Points(inputobj.vertices, r=8, c='k') 2507 if hasattr(inputobj, "vertices_color"): 2508 vcols = (inputobj.vertices_color * 1).astype(np.uint8) 2509 vdpts.pointcolors = vcols 2510 return vdpts 2511 2512 if "path" in inputobj_type: 2513 2514 lines = [] 2515 for e in inputobj.entities: 2516 # print('trimesh entity', e.to_dict()) 2517 l = vedo.shapes.Line(inputobj.vertices[e.points], c="k", lw=2) 2518 lines.append(l) 2519 return vedo.Assembly(lines) 2520 2521 return None 2522 2523 2524def vedo2meshlab(vmesh): 2525 """Convert a `vedo.Mesh` to a Meshlab object.""" 2526 try: 2527 import pymeshlab as mlab 2528 except ModuleNotFoundError: 2529 vedo.logger.error("Need pymeshlab to run:\npip install pymeshlab") 2530 2531 vertex_matrix = vmesh.vertices.astype(np.float64) 2532 2533 try: 2534 face_matrix = np.asarray(vmesh.cells, dtype=np.float64) 2535 except: 2536 print("WARNING: in vedo2meshlab(), need to triangulate mesh first!") 2537 face_matrix = np.array(vmesh.clone().triangulate().cells, dtype=np.float64) 2538 2539 # v_normals_matrix = vmesh.normals(cells=False, recompute=False) 2540 v_normals_matrix = vmesh.vertex_normals 2541 if not v_normals_matrix.shape[0]: 2542 v_normals_matrix = np.empty((0, 3), dtype=np.float64) 2543 2544 # f_normals_matrix = vmesh.normals(cells=True, recompute=False) 2545 f_normals_matrix = vmesh.cell_normals 2546 if not f_normals_matrix.shape[0]: 2547 f_normals_matrix = np.empty((0, 3), dtype=np.float64) 2548 2549 v_color_matrix = vmesh.pointdata["RGBA"] 2550 if v_color_matrix is None: 2551 v_color_matrix = np.empty((0, 4), dtype=np.float64) 2552 else: 2553 v_color_matrix = v_color_matrix.astype(np.float64) / 255 2554 if v_color_matrix.shape[1] == 3: 2555 v_color_matrix = np.c_[ 2556 v_color_matrix, np.ones(v_color_matrix.shape[0], dtype=np.float64) 2557 ] 2558 2559 f_color_matrix = vmesh.celldata["RGBA"] 2560 if f_color_matrix is None: 2561 f_color_matrix = np.empty((0, 4), dtype=np.float64) 2562 else: 2563 f_color_matrix = f_color_matrix.astype(np.float64) / 255 2564 if f_color_matrix.shape[1] == 3: 2565 f_color_matrix = np.c_[ 2566 f_color_matrix, np.ones(f_color_matrix.shape[0], dtype=np.float64) 2567 ] 2568 2569 m = mlab.Mesh( 2570 vertex_matrix=vertex_matrix, 2571 face_matrix=face_matrix, 2572 v_normals_matrix=v_normals_matrix, 2573 f_normals_matrix=f_normals_matrix, 2574 v_color_matrix=v_color_matrix, 2575 f_color_matrix=f_color_matrix, 2576 ) 2577 2578 for k in vmesh.pointdata.keys(): 2579 data = vmesh.pointdata[k] 2580 if data is not None: 2581 if data.ndim == 1: # scalar 2582 m.add_vertex_custom_scalar_attribute(data.astype(np.float64), k) 2583 elif data.ndim == 2: # vectorial data 2584 if "tcoord" not in k.lower() and k not in ["Normals", "TextureCoordinates"]: 2585 m.add_vertex_custom_point_attribute(data.astype(np.float64), k) 2586 2587 for k in vmesh.celldata.keys(): 2588 data = vmesh.celldata[k] 2589 if data is not None: 2590 if data.ndim == 1: # scalar 2591 m.add_face_custom_scalar_attribute(data.astype(np.float64), k) 2592 elif data.ndim == 2 and k != "Normals": # vectorial data 2593 m.add_face_custom_point_attribute(data.astype(np.float64), k) 2594 2595 m.update_bounding_box() 2596 return m 2597 2598 2599def meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()): 2600 """Convert a Meshlab object to `vedo.Mesh`.""" 2601 inputtype = str(type(mmesh)) 2602 2603 if "MeshSet" in inputtype: 2604 mmesh = mmesh.current_mesh() 2605 2606 mpoints, mcells = mmesh.vertex_matrix(), mmesh.face_matrix() 2607 if len(mcells) > 0: 2608 polydata = buildPolyData(mpoints, mcells) 2609 else: 2610 polydata = buildPolyData(mpoints, None) 2611 2612 if mmesh.has_vertex_scalar(): 2613 parr = mmesh.vertex_scalar_array() 2614 parr_vtk = numpy_to_vtk(parr) 2615 parr_vtk.SetName("MeshLabScalars") 2616 polydata.GetPointData().AddArray(parr_vtk) 2617 polydata.GetPointData().SetActiveScalars("MeshLabScalars") 2618 2619 if mmesh.has_face_scalar(): 2620 carr = mmesh.face_scalar_array() 2621 carr_vtk = numpy_to_vtk(carr) 2622 carr_vtk.SetName("MeshLabScalars") 2623 x0, x1 = carr_vtk.GetRange() 2624 polydata.GetCellData().AddArray(carr_vtk) 2625 polydata.GetCellData().SetActiveScalars("MeshLabScalars") 2626 2627 for k in pointdata_keys: 2628 parr = mmesh.vertex_custom_scalar_attribute_array(k) 2629 parr_vtk = numpy_to_vtk(parr) 2630 parr_vtk.SetName(k) 2631 polydata.GetPointData().AddArray(parr_vtk) 2632 polydata.GetPointData().SetActiveScalars(k) 2633 2634 for k in celldata_keys: 2635 carr = mmesh.face_custom_scalar_attribute_array(k) 2636 carr_vtk = numpy_to_vtk(carr) 2637 carr_vtk.SetName(k) 2638 polydata.GetCellData().AddArray(carr_vtk) 2639 polydata.GetCellData().SetActiveScalars(k) 2640 2641 pnorms = mmesh.vertex_normal_matrix() 2642 if len(pnorms) > 0: 2643 polydata.GetPointData().SetNormals(numpy2vtk(pnorms, name="Normals")) 2644 2645 cnorms = mmesh.face_normal_matrix() 2646 if len(cnorms) > 0: 2647 polydata.GetCellData().SetNormals(numpy2vtk(cnorms, name="Normals")) 2648 return vedo.Mesh(polydata) 2649 2650 2651def open3d2vedo(o3d_mesh): 2652 """Convert `open3d.geometry.TriangleMesh` to a `vedo.Mesh`.""" 2653 m = vedo.Mesh([np.array(o3d_mesh.vertices), np.array(o3d_mesh.triangles)]) 2654 # TODO: could also check whether normals and color are present in 2655 # order to port with the above vertices/faces 2656 return m 2657 2658 2659def vedo2open3d(vedo_mesh): 2660 """ 2661 Return an `open3d.geometry.TriangleMesh` version of the current mesh. 2662 """ 2663 try: 2664 import open3d as o3d 2665 except RuntimeError: 2666 vedo.logger.error("Need open3d to run:\npip install open3d") 2667 2668 # create from numpy arrays 2669 o3d_mesh = o3d.geometry.TriangleMesh( 2670 vertices=o3d.utility.Vector3dVector(vedo_mesh.vertices), 2671 triangles=o3d.utility.Vector3iVector(vedo_mesh.cells), 2672 ) 2673 # TODO: need to add some if check here in case color and normals 2674 # info are not existing 2675 # o3d_mesh.vertex_colors = o3d.utility.Vector3dVector(vedo_mesh.pointdata["RGB"]/255) 2676 # o3d_mesh.vertex_normals= o3d.utility.Vector3dVector(vedo_mesh.pointdata["Normals"]) 2677 return o3d_mesh 2678 2679def vedo2madcad(vedo_mesh): 2680 """ 2681 Convert a `vedo.Mesh` to a `madcad.Mesh`. 2682 """ 2683 try: 2684 import madcad 2685 import numbers 2686 except ModuleNotFoundError: 2687 vedo.logger.error("Need madcad to run:\npip install pymadcad") 2688 2689 points = [madcad.vec3(*pt) for pt in vedo_mesh.vertices] 2690 faces = [madcad.vec3(*fc) for fc in vedo_mesh.cells] 2691 2692 options = {} 2693 for key, val in vedo_mesh.pointdata.items(): 2694 vec_type = f"vec{val.shape[-1]}" 2695 is_float = np.issubdtype(val.dtype, np.floating) 2696 madcad_dtype = getattr(madcad, f"f{vec_type}" if is_float else vec_type) 2697 options[key] = [madcad_dtype(v) for v in val] 2698 2699 madcad_mesh = madcad.Mesh(points=points, faces=faces, options=options) 2700 2701 return madcad_mesh 2702 2703 2704def madcad2vedo(madcad_mesh): 2705 """ 2706 Convert a `madcad.Mesh` to a `vedo.Mesh`. 2707 2708 A pointdata or celldata array named "tracks" is added to the output mesh, indicating 2709 the mesh region each point belongs to. 2710 2711 A metadata array named "madcad_groups" is added to the output mesh, indicating 2712 the mesh groups. 2713 2714 See [pymadcad website](https://pymadcad.readthedocs.io/en/latest/index.html) 2715 for more info. 2716 """ 2717 try: 2718 madcad_mesh = madcad_mesh["part"] 2719 except: 2720 pass 2721 2722 madp = [] 2723 for p in madcad_mesh.points: 2724 madp.append([float(p[0]), float(p[1]), float(p[2])]) 2725 madp = np.array(madp) 2726 2727 madf = [] 2728 try: 2729 for f in madcad_mesh.faces: 2730 madf.append([int(f[0]), int(f[1]), int(f[2])]) 2731 madf = np.array(madf).astype(np.uint16) 2732 except AttributeError: 2733 # print("no faces") 2734 pass 2735 2736 made = [] 2737 try: 2738 edges = madcad_mesh.edges 2739 for e in edges: 2740 made.append([int(e[0]), int(e[1])]) 2741 made = np.array(made).astype(np.uint16) 2742 except (AttributeError, TypeError): 2743 # print("no edges") 2744 pass 2745 2746 try: 2747 line = np.array(madcad_mesh.indices).astype(np.uint16) 2748 made.append(line) 2749 except AttributeError: 2750 # print("no indices") 2751 pass 2752 2753 madt = [] 2754 try: 2755 for t in madcad_mesh.tracks: 2756 madt.append(int(t)) 2757 madt = np.array(madt).astype(np.uint16) 2758 except AttributeError: 2759 # print("no tracks") 2760 pass 2761 2762 ############################### 2763 poly = vedo.utils.buildPolyData(madp, madf, made) 2764 if len(madf) == 0 and len(made) == 0: 2765 m = vedo.Points(poly) 2766 else: 2767 m = vedo.Mesh(poly) 2768 2769 if len(madt) == len(madf): 2770 m.celldata["tracks"] = madt 2771 maxt = np.max(madt) 2772 m.mapper.SetScalarRange(0, np.max(madt)) 2773 if maxt==0: m.mapper.SetScalarVisibility(0) 2774 elif len(madt) == len(madp): 2775 m.pointdata["tracks"] = madt 2776 maxt = np.max(madt) 2777 m.mapper.SetScalarRange(0, maxt) 2778 if maxt==0: m.mapper.SetScalarVisibility(0) 2779 2780 try: 2781 m.info["madcad_groups"] = madcad_mesh.groups 2782 except AttributeError: 2783 # print("no groups") 2784 pass 2785 2786 try: 2787 options = dict(madcad_mesh.options) 2788 if "display_wire" in options and options["display_wire"]: 2789 m.lw(1).lc(madcad_mesh.c()) 2790 if "display_faces" in options and not options["display_faces"]: 2791 m.alpha(0.2) 2792 if "color" in options: 2793 m.c(options["color"]) 2794 2795 for key, val in options.items(): 2796 m.pointdata[key] = val 2797 2798 except AttributeError: 2799 # print("no options") 2800 pass 2801 2802 return m 2803 2804 2805def vtk_version_at_least(major, minor=0, build=0) -> bool: 2806 """ 2807 Check the installed VTK version. 2808 2809 Return `True` if the requested VTK version is greater or equal to the actual VTK version. 2810 2811 Arguments: 2812 major : (int) 2813 Major version. 2814 minor : (int) 2815 Minor version. 2816 build : (int) 2817 Build version. 2818 """ 2819 needed_version = 10000000000 * int(major) + 100000000 * int(minor) + int(build) 2820 try: 2821 vtk_version_number = vtki.VTK_VERSION_NUMBER 2822 except AttributeError: # as error: 2823 ver = vtki.vtkVersion() 2824 vtk_version_number = ( 2825 10000000000 * ver.GetVTKMajorVersion() 2826 + 100000000 * ver.GetVTKMinorVersion() 2827 + ver.GetVTKBuildVersion() 2828 ) 2829 return vtk_version_number >= needed_version 2830 2831 2832def ctf2lut(vol, logscale=False): 2833 """Internal use.""" 2834 # build LUT from a color transfer function for tmesh or volume 2835 2836 ctf = vol.properties.GetRGBTransferFunction() 2837 otf = vol.properties.GetScalarOpacity() 2838 x0, x1 = vol.dataset.GetScalarRange() 2839 cols, alphas = [], [] 2840 for x in np.linspace(x0, x1, 256): 2841 cols.append(ctf.GetColor(x)) 2842 alphas.append(otf.GetValue(x)) 2843 2844 if logscale: 2845 lut = vtki.vtkLogLookupTable() 2846 else: 2847 lut = vtki.vtkLookupTable() 2848 2849 lut.SetRange(x0, x1) 2850 lut.SetNumberOfTableValues(len(cols)) 2851 for i, col in enumerate(cols): 2852 r, g, b = col 2853 lut.SetTableValue(i, r, g, b, alphas[i]) 2854 lut.Build() 2855 return lut 2856 2857 2858def get_vtk_name_event(name: str) -> str: 2859 """ 2860 Return the name of a VTK event. 2861 2862 Frequently used events are: 2863 - KeyPress, KeyRelease: listen to keyboard events 2864 - LeftButtonPress, LeftButtonRelease: listen to mouse clicks 2865 - MiddleButtonPress, MiddleButtonRelease 2866 - RightButtonPress, RightButtonRelease 2867 - MouseMove: listen to mouse pointer changing position 2868 - MouseWheelForward, MouseWheelBackward 2869 - Enter, Leave: listen to mouse entering or leaving the window 2870 - Pick, StartPick, EndPick: listen to object picking 2871 - ResetCamera, ResetCameraClippingRange 2872 - Error, Warning, Char, Timer 2873 2874 Check the complete list of events here: 2875 https://vtk.org/doc/nightly/html/classvtkCommand.html 2876 """ 2877 # as vtk names are ugly and difficult to remember: 2878 ln = name.lower() 2879 if "click" in ln or "button" in ln: 2880 event_name = "LeftButtonPress" 2881 if "right" in ln: 2882 event_name = "RightButtonPress" 2883 elif "mid" in ln: 2884 event_name = "MiddleButtonPress" 2885 if "release" in ln: 2886 event_name = event_name.replace("Press", "Release") 2887 else: 2888 event_name = name 2889 if "key" in ln: 2890 if "release" in ln: 2891 event_name = "KeyRelease" 2892 else: 2893 event_name = "KeyPress" 2894 2895 if ("mouse" in ln and "mov" in ln) or "over" in ln: 2896 event_name = "MouseMove" 2897 2898 words = [ 2899 "pick", "timer", "reset", "enter", "leave", "char", 2900 "error", "warning", "start", "end", "wheel", "clipping", 2901 "range", "camera", "render", "interaction", "modified", 2902 ] 2903 for w in words: 2904 if w in ln: 2905 event_name = event_name.replace(w, w.capitalize()) 2906 2907 event_name = event_name.replace("REnd ", "Rend") 2908 event_name = event_name.replace("the ", "") 2909 event_name = event_name.replace(" of ", "").replace(" ", "") 2910 2911 if not event_name.endswith("Event"): 2912 event_name += "Event" 2913 2914 if vtki.vtkCommand.GetEventIdFromString(event_name) == 0: 2915 vedo.printc( 2916 f"Error: '{name}' is not a valid event name.", c='r') 2917 vedo.printc("Check the list of events here:", c='r') 2918 vedo.printc("\thttps://vtk.org/doc/nightly/html/classvtkCommand.html", c='r') 2919 # raise RuntimeError 2920 2921 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.
850def geometry(obj, extent=None) -> "vedo.Mesh": 851 """ 852 Apply the `vtkGeometryFilter` to the input object. 853 This is a general-purpose filter to extract geometry (and associated data) 854 from any type of dataset. 855 This filter also may be used to convert any type of data to polygonal type. 856 The conversion process may be less than satisfactory for some 3D datasets. 857 For example, this filter will extract the outer surface of a volume 858 or structured grid dataset. 859 860 Returns a `vedo.Mesh` object. 861 862 Set `extent` as the `[xmin,xmax, ymin,ymax, zmin,zmax]` bounding box to clip data. 863 """ 864 gf = vtki.new("GeometryFilter") 865 gf.SetInputData(obj) 866 if extent is not None: 867 gf.SetExtent(extent) 868 gf.Update() 869 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:
1038def is_sequence(arg) -> bool: 1039 """Check if the input is iterable.""" 1040 if hasattr(arg, "strip"): 1041 return False 1042 if hasattr(arg, "__getslice__"): 1043 return True 1044 if hasattr(arg, "__iter__"): 1045 return True 1046 return False
Check if the input is iterable.
def
lin_interpolate(x, rangeX, rangeY):
1412def lin_interpolate(x, rangeX, rangeY): 1413 """ 1414 Interpolate linearly the variable `x` in `rangeX` onto the new `rangeY`. 1415 If `x` is a 3D vector the linear weight is the distance to the two 3D `rangeX` vectors. 1416 1417 E.g. if `x` runs in `rangeX=[x0,x1]` and I want it to run in `rangeY=[y0,y1]` then 1418 1419 `y = lin_interpolate(x, rangeX, rangeY)` will interpolate `x` onto `rangeY`. 1420 1421 Examples: 1422 - [lin_interpolate.py](https://github.com/marcomusy/vedo/tree/master/examples/basic/lin_interpolate.py) 1423 1424  1425 """ 1426 if is_sequence(x): 1427 x = np.asarray(x) 1428 x0, x1 = np.asarray(rangeX) 1429 y0, y1 = np.asarray(rangeY) 1430 dx = x1 - x0 1431 dxn = np.linalg.norm(dx) 1432 if not dxn: 1433 return y0 1434 s = np.linalg.norm(x - x0) / dxn 1435 t = np.linalg.norm(x - x1) / dxn 1436 st = s + t 1437 out = y0 * (t / st) + y1 * (s / st) 1438 1439 else: # faster 1440 1441 x0 = rangeX[0] 1442 dx = rangeX[1] - x0 1443 if not dx: 1444 return rangeY[0] 1445 s = (x - x0) / dx 1446 out = rangeY[0] * (1 - s) + rangeY[1] * s 1447 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:
1524def vector(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1525 """ 1526 Return a 3D numpy array representing a vector. 1527 1528 If `y` is `None`, assume input is already in the form `[x,y,z]`. 1529 """ 1530 if y is None: # assume x is already [x,y,z] 1531 return np.asarray(x, dtype=dtype) 1532 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):
1543def mag(v): 1544 """Get the magnitude of a vector or array of vectors.""" 1545 v = np.asarray(v) 1546 if v.ndim == 1: 1547 return np.linalg.norm(v) 1548 return np.linalg.norm(v, axis=1)
Get the magnitude of a vector or array of vectors.
def
mag2(v) -> numpy.ndarray:
1551def mag2(v) -> np.ndarray: 1552 """Get the squared magnitude of a vector or array of vectors.""" 1553 v = np.asarray(v) 1554 if v.ndim == 1: 1555 return np.square(v).sum() 1556 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:
1535def versor(x, y=None, z=0.0, dtype=np.float64) -> np.ndarray: 1536 """Return the unit vector. Input can be a list of vectors.""" 1537 v = vector(x, y, z, dtype) 1538 if isinstance(v[0], np.ndarray): 1539 return np.divide(v, mag(v)[:, None]) 1540 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:
1625def precision(x, p: int, vrange=None, delimiter="e") -> str: 1626 """ 1627 Returns a string representation of `x` formatted to precision `p`. 1628 1629 Set `vrange` to the range in which x exists (to snap x to '0' if below precision). 1630 """ 1631 # Based on the webkit javascript implementation 1632 # `from here <https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp>`_, 1633 # and implemented by `randlet <https://github.com/randlet/to-precision>`_. 1634 # Modified for vedo by M.Musy 2020 1635 1636 if isinstance(x, str): # do nothing 1637 return x 1638 1639 if is_sequence(x): 1640 out = "(" 1641 nn = len(x) - 1 1642 for i, ix in enumerate(x): 1643 1644 try: 1645 if np.isnan(ix): 1646 return "NaN" 1647 except: 1648 # cannot handle list of list 1649 continue 1650 1651 out += precision(ix, p) 1652 if i < nn: 1653 out += ", " 1654 return out + ")" ############ <-- 1655 1656 try: 1657 if np.isnan(x): 1658 return "NaN" 1659 except TypeError: 1660 return "NaN" 1661 1662 x = float(x) 1663 1664 if x == 0.0 or (vrange is not None and abs(x) < vrange / pow(10, p)): 1665 return "0" 1666 1667 out = [] 1668 if x < 0: 1669 out.append("-") 1670 x = -x 1671 1672 e = int(np.log10(x)) 1673 # tens = np.power(10, e - p + 1) 1674 tens = 10 ** (e - p + 1) 1675 n = np.floor(x / tens) 1676 1677 # if n < np.power(10, p - 1): 1678 if n < 10 ** (p - 1): 1679 e = e - 1 1680 # tens = np.power(10, e - p + 1) 1681 tens = 10 ** (e - p + 1) 1682 n = np.floor(x / tens) 1683 1684 if abs((n + 1.0) * tens - x) <= abs(n * tens - x): 1685 n = n + 1 1686 1687 # if n >= np.power(10, p): 1688 if n >= 10 ** p: 1689 n = n / 10.0 1690 e = e + 1 1691 1692 m = "%.*g" % (p, n) 1693 if e < -2 or e >= p: 1694 out.append(m[0]) 1695 if p > 1: 1696 out.append(".") 1697 out.extend(m[1:p]) 1698 out.append(delimiter) 1699 if e > 0: 1700 out.append("+") 1701 out.append(str(e)) 1702 elif e == (p - 1): 1703 out.append(m) 1704 elif e >= 0: 1705 out.append(m[: e + 1]) 1706 if e + 1 < len(m): 1707 out.append(".") 1708 out.extend(m[e + 1 :]) 1709 else: 1710 out.append("0.") 1711 out.extend(["0"] * -(e + 1)) 1712 out.append(m) 1713 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:
1578def round_to_digit(x, p) -> float: 1579 """Round a real number to the specified number of significant digits.""" 1580 if not x: 1581 return 0 1582 r = np.round(x, p - int(np.floor(np.log10(abs(x)))) - 1) 1583 if int(r) == r: 1584 return int(r) 1585 return r
Round a real number to the specified number of significant digits.
def
point_in_triangle(p, p1, p2, p3) -> Optional[bool]:
1121def point_in_triangle(p, p1, p2, p3) -> Union[bool, None]: 1122 """ 1123 Return True if a point is inside (or above/below) 1124 a triangle defined by 3 points in space. 1125 """ 1126 p1 = np.array(p1) 1127 u = p2 - p1 1128 v = p3 - p1 1129 n = np.cross(u, v) 1130 w = p - p1 1131 ln = np.dot(n, n) 1132 if not ln: 1133 return None # degenerate triangle 1134 gamma = (np.dot(np.cross(u, w), n)) / ln 1135 if 0 < gamma < 1: 1136 beta = (np.dot(np.cross(w, v), n)) / ln 1137 if 0 < beta < 1: 1138 alpha = 1 - gamma - beta 1139 if 0 < alpha < 1: 1140 return True 1141 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:
1334def point_line_distance(p, p1, p2) -> float: 1335 """ 1336 Compute the distance of a point to a line (not the segment) 1337 defined by `p1` and `p2`. 1338 """ 1339 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):
1374def closest(point, points, n=1, return_ids=False, use_tree=False): 1375 """ 1376 Returns the distances and the closest point(s) to the given set of points. 1377 Needs `scipy.spatial` library. 1378 1379 Arguments: 1380 n : (int) 1381 the nr of closest points to return 1382 return_ids : (bool) 1383 return the ids instead of the points coordinates 1384 use_tree : (bool) 1385 build a `scipy.spatial.KDTree`. 1386 An already existing one can be passed to avoid rebuilding. 1387 """ 1388 from scipy.spatial import distance, KDTree 1389 1390 points = np.asarray(points) 1391 if n == 1: 1392 dists = distance.cdist([point], points) 1393 closest_idx = np.argmin(dists) 1394 else: 1395 if use_tree: 1396 if isinstance(use_tree, KDTree): # reuse 1397 tree = use_tree 1398 else: 1399 tree = KDTree(points) 1400 dists, closest_idx = tree.query([point], k=n) 1401 closest_idx = closest_idx[0] 1402 else: 1403 dists = distance.cdist([point], points) 1404 closest_idx = np.argsort(dists)[0][:n] 1405 if return_ids: 1406 return dists, closest_idx 1407 else: 1408 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:
1717def grep(filename: str, tag: str, column=None, first_occurrence_only=False) -> list: 1718 """Greps the line in a file that starts with a specific `tag` string inside the file.""" 1719 import re 1720 1721 with open(filename, "r", encoding="UTF-8") as afile: 1722 content = [] 1723 for line in afile: 1724 if re.search(tag, line): 1725 c = line.split() 1726 c[-1] = c[-1].replace("\n", "") 1727 if column is not None: 1728 c = c[column] 1729 content.append(c) 1730 if first_occurrence_only: 1731 break 1732 return content
Greps the line in a file that starts with a specific tag string inside the file.
def
make_bands(inputlist, n):
2015def make_bands(inputlist, n): 2016 """ 2017 Group values of a list into bands of equal value, where 2018 `n` is the number of bands, a positive integer > 2. 2019 2020 Returns a binned list of the same length as the input. 2021 """ 2022 if n < 2: 2023 return inputlist 2024 vmin = np.min(inputlist) 2025 vmax = np.max(inputlist) 2026 bb = np.linspace(vmin, vmax, n, endpoint=0) 2027 dr = bb[1] - bb[0] 2028 bb += dr / 2 2029 tol = dr / 2 * 1.001 2030 newlist = [] 2031 for s in inputlist: 2032 for b in bb: 2033 if abs(s - b) < tol: 2034 newlist.append(b) 2035 break 2036 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:
1588def pack_spheres(bounds, radius) -> np.ndarray: 1589 """ 1590 Packing spheres into a bounding box. 1591 Returns a numpy array of sphere centers. 1592 """ 1593 h = 0.8164965 / 2 1594 d = 0.8660254 1595 a = 0.288675135 1596 1597 if is_sequence(bounds): 1598 x0, x1, y0, y1, z0, z1 = bounds 1599 else: 1600 x0, x1, y0, y1, z0, z1 = bounds.bounds() 1601 1602 x = np.arange(x0, x1, radius) 1603 nul = np.zeros_like(x) 1604 nz = int((z1 - z0) / radius / h / 2 + 1.5) 1605 ny = int((y1 - y0) / radius / d + 1.5) 1606 1607 pts = [] 1608 for iz in range(nz): 1609 z = z0 + nul + iz * h * radius 1610 dx, dy, dz = [radius * 0.5, radius * a, iz * h * radius] 1611 for iy in range(ny): 1612 y = y0 + nul + iy * d * radius 1613 if iy % 2: 1614 xs = x 1615 else: 1616 xs = x + radius * 0.5 1617 if iz % 2: 1618 p = np.c_[xs, y, z] + [dx, dy, dz] 1619 else: 1620 p = np.c_[xs, y, z] + [0, 0, dz] 1621 pts += p.tolist() 1622 return np.array(pts)
Packing spheres into a bounding box. Returns a numpy array of sphere centers.
def
humansort(alist) -> list:
1088def humansort(alist) -> list: 1089 """ 1090 Sort in place a given list the way humans expect. 1091 1092 E.g. `['file11', 'file1'] -> ['file1', 'file11']` 1093 1094 .. warning:: input list is modified in-place by this function. 1095 """ 1096 import re 1097 1098 def alphanum_key(s): 1099 # Turn a string into a list of string and number chunks. 1100 # e.g. "z23a" -> ["z", 23, "a"] 1101 def tryint(s): 1102 if s.isdigit(): 1103 return int(s) 1104 return s 1105 1106 return [tryint(c) for c in re.split("([0-9]+)", s)] 1107 1108 alist.sort(key=alphanum_key) 1109 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:
1790def print_histogram( 1791 data, 1792 bins=10, 1793 height=10, 1794 logscale=False, 1795 minbin=0, 1796 horizontal=True, 1797 char="\U00002589", 1798 c=None, 1799 bold=True, 1800 title="histogram", 1801 spacer="", 1802) -> np.ndarray: 1803 """ 1804 Ascii histogram printing. 1805 1806 Input can be a `vedo.Volume` or `vedo.Mesh`. 1807 Returns the raw data before binning (useful when passing vtk objects). 1808 1809 Arguments: 1810 bins : (int) 1811 number of histogram bins 1812 height : (int) 1813 height of the histogram in character units 1814 logscale : (bool) 1815 use logscale for frequencies 1816 minbin : (int) 1817 ignore bins before minbin 1818 horizontal : (bool) 1819 show histogram horizontally 1820 char : (str) 1821 character to be used 1822 bold : (bool) 1823 use boldface 1824 title : (str) 1825 histogram title 1826 spacer : (str) 1827 horizontal spacer 1828 1829 Example: 1830 ```python 1831 from vedo import print_histogram 1832 import numpy as np 1833 d = np.random.normal(size=1000) 1834 data = print_histogram(d, c='b', logscale=True, title='my scalars') 1835 data = print_histogram(d, c='o') 1836 print(np.mean(data)) # data here is same as d 1837 ``` 1838  1839 """ 1840 # credits: http://pyinsci.blogspot.com/2009/10/ascii-histograms.html 1841 # adapted for vedo by M.Musy, 2019 1842 1843 if not horizontal: # better aspect ratio 1844 bins *= 2 1845 1846 try: 1847 data = vtk2numpy(data.dataset.GetPointData().GetScalars()) 1848 except AttributeError: 1849 # already an array 1850 data = np.asarray(data) 1851 1852 if isinstance(data, vtki.vtkImageData): 1853 dims = data.GetDimensions() 1854 nvx = min(100000, dims[0] * dims[1] * dims[2]) 1855 idxs = np.random.randint(0, min(dims), size=(nvx, 3)) 1856 data = [] 1857 for ix, iy, iz in idxs: 1858 d = data.GetScalarComponentAsFloat(ix, iy, iz, 0) 1859 data.append(d) 1860 data = np.array(data) 1861 1862 elif isinstance(data, vtki.vtkPolyData): 1863 arr = data.GetPointData().GetScalars() 1864 if not arr: 1865 arr = data.GetCellData().GetScalars() 1866 if not arr: 1867 return np.array([]) 1868 data = vtk2numpy(arr) 1869 1870 try: 1871 h = np.histogram(data, bins=bins) 1872 except TypeError as e: 1873 vedo.logger.error(f"cannot compute histogram: {e}") 1874 return np.array([]) 1875 1876 if minbin: 1877 hi = h[0][minbin:-1] 1878 else: 1879 hi = h[0] 1880 1881 if char == "\U00002589" and horizontal: 1882 char = "\U00002586" 1883 1884 title = title.ljust(14) + ":" 1885 entrs = " entries=" + str(len(data)) 1886 if logscale: 1887 h0 = np.log10(hi + 1) 1888 maxh0 = int(max(h0) * 100) / 100 1889 title = title + entrs + " (logscale)" 1890 else: 1891 h0 = hi 1892 maxh0 = max(h0) 1893 title = title + entrs 1894 1895 def _v(): 1896 his = "" 1897 if title: 1898 his += title + "\n" 1899 bars = h0 / maxh0 * height 1900 for l in reversed(range(1, height + 1)): 1901 line = "" 1902 if l == height: 1903 line = "%s " % maxh0 1904 else: 1905 line = " |" + " " * (len(str(maxh0)) - 3) 1906 for c in bars: 1907 if c >= np.ceil(l): 1908 line += char 1909 else: 1910 line += " " 1911 line += "\n" 1912 his += line 1913 his += "%.2f" % h[1][0] + "." * (bins) + "%.2f" % h[1][-1] + "\n" 1914 return his 1915 1916 def _h(): 1917 his = "" 1918 if title: 1919 his += title + "\n" 1920 xl = ["%.2f" % n for n in h[1]] 1921 lxl = [len(l) for l in xl] 1922 bars = h0 / maxh0 * height 1923 his += spacer + " " * int(max(bars) + 2 + max(lxl)) + "%s\n" % maxh0 1924 for i, c in enumerate(bars): 1925 line = xl[i] + " " * int(max(lxl) - lxl[i]) + "| " + char * int(c) + "\n" 1926 his += spacer + line 1927 return his 1928 1929 if horizontal: 1930 height *= 2 1931 vedo.printc(_h(), c=c, bold=bold) 1932 else: 1933 vedo.printc(_v(), c=c, bold=bold) 1934 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:
1991def print_inheritance_tree(C) -> None: 1992 """Prints the inheritance tree of class C.""" 1993 # Adapted from: https://stackoverflow.com/questions/26568976/ 1994 def class_tree(cls): 1995 subc = [class_tree(sub_class) for sub_class in cls.__subclasses__()] 1996 return {cls.__name__: subc} 1997 1998 def print_tree(tree, indent=8, current_ind=0): 1999 for k, v in tree.items(): 2000 if current_ind: 2001 before_dashes = current_ind - indent 2002 m = " " * before_dashes + "└" + "─" * (indent - 1) + " " + k 2003 vedo.printc(m) 2004 else: 2005 vedo.printc(k) 2006 for sub_tree in v: 2007 print_tree(sub_tree, indent=indent, current_ind=current_ind + indent) 2008 2009 if str(C.__class__) != "<class 'type'>": 2010 C = C.__class__ 2011 ct = class_tree(C) 2012 print_tree(ct)
Prints the inheritance tree of class C.
def
camera_from_quaternion( pos, quaternion, distance=10000, ngl_correct=True) -> vtkmodules.vtkRenderingCore.vtkCamera:
2042def camera_from_quaternion(pos, quaternion, distance=10000, ngl_correct=True) -> vtki.vtkCamera: 2043 """ 2044 Define a `vtkCamera` with a particular orientation. 2045 2046 Arguments: 2047 pos: (np.array, list, tuple) 2048 an iterator of length 3 containing the focus point of the camera 2049 quaternion: (np.array, list, tuple) 2050 a len(4) quaternion `(x,y,z,w)` describing the rotation of the camera 2051 such as returned by neuroglancer `x,y,z,w` all in `[0,1]` range 2052 distance: (float) 2053 the desired distance from pos to the camera (default = 10000 nm) 2054 2055 Returns: 2056 `vtki.vtkCamera`, a vtk camera setup according to these rules. 2057 """ 2058 camera = vtki.vtkCamera() 2059 # define the quaternion in vtk, note the swapped order 2060 # w,x,y,z instead of x,y,z,w 2061 quat_vtk = vtki.get_class("Quaternion")( 2062 quaternion[3], quaternion[0], quaternion[1], quaternion[2]) 2063 # use this to define a rotation matrix in x,y,z 2064 # right handed units 2065 M = np.zeros((3, 3), dtype=np.float32) 2066 quat_vtk.ToMatrix3x3(M) 2067 # the default camera orientation is y up 2068 up = [0, 1, 0] 2069 # calculate default camera position is backed off in positive z 2070 pos = [0, 0, distance] 2071 2072 # set the camera rototation by applying the rotation matrix 2073 camera.SetViewUp(*np.dot(M, up)) 2074 # set the camera position by applying the rotation matrix 2075 camera.SetPosition(*np.dot(M, pos)) 2076 if ngl_correct: 2077 # neuroglancer has positive y going down 2078 # so apply these azimuth and roll corrections 2079 # to fix orientatins 2080 camera.Azimuth(-180) 2081 camera.Roll(180) 2082 2083 # shift the camera posiiton and focal position 2084 # to be centered on the desired location 2085 p = camera.GetPosition() 2086 p_new = np.array(p) + pos 2087 camera.SetPosition(*p_new) 2088 camera.SetFocalPoint(*pos) 2089 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=300) -> vtkmodules.vtkRenderingCore.vtkCamera:
2092def camera_from_neuroglancer(state, zoom=300) -> vtki.vtkCamera: 2093 """ 2094 Define a `vtkCamera` from a neuroglancer state dictionary. 2095 2096 Arguments: 2097 state: (dict) 2098 an neuroglancer state dictionary. 2099 zoom: (float) 2100 how much to multiply zoom by to get camera backoff distance 2101 default = 300 > ngl_zoom = 1 > 300 nm backoff distance. 2102 2103 Returns: 2104 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2105 """ 2106 orient = state.get("perspectiveOrientation", [0.0, 0.0, 0.0, 1.0]) 2107 pzoom = state.get("perspectiveZoom", 10.0) 2108 position = state["navigation"]["pose"]["position"] 2109 pos_nm = np.array(position["voxelCoordinates"]) * position["voxelSize"] 2110 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 default = 300 > ngl_zoom = 1 > 300 nm backoff distance.
Returns:
vtki.vtkCamera, a vtk camera setup that matches this state.
def
camera_from_dict(camera, modify_inplace=None) -> vtkmodules.vtkRenderingCore.vtkCamera:
2128def camera_from_dict(camera, modify_inplace=None) -> vtki.vtkCamera: 2129 """ 2130 Generate a `vtkCamera` object from a python dictionary. 2131 2132 Parameters of the camera are: 2133 - `position` or `pos` (3-tuple) 2134 - `focal_point` (3-tuple) 2135 - `viewup` (3-tuple) 2136 - `distance` (float) 2137 - `clipping_range` (2-tuple) 2138 - `parallel_scale` (float) 2139 - `thickness` (float) 2140 - `view_angle` (float) 2141 - `roll` (float) 2142 2143 Exaplanation of the parameters can be found in the 2144 [vtkCamera documentation](https://vtk.org/doc/nightly/html/classvtkCamera.html). 2145 2146 Arguments: 2147 camera: (dict) 2148 a python dictionary containing camera parameters. 2149 modify_inplace: (vtkCamera) 2150 an existing `vtkCamera` object to modify in place. 2151 2152 Returns: 2153 `vtki.vtkCamera`, a vtk camera setup that matches this state. 2154 """ 2155 if modify_inplace: 2156 vcam = modify_inplace 2157 else: 2158 vcam = vtki.vtkCamera() 2159 2160 camera = dict(camera) # make a copy so input is not emptied by pop() 2161 2162 cm_pos = camera.pop("position", camera.pop("pos", None)) 2163 cm_focal_point = camera.pop("focal_point", camera.pop("focalPoint", None)) 2164 cm_viewup = camera.pop("viewup", None) 2165 cm_distance = camera.pop("distance", None) 2166 cm_clipping_range = camera.pop("clipping_range", camera.pop("clippingRange", None)) 2167 cm_parallel_scale = camera.pop("parallel_scale", camera.pop("parallelScale", None)) 2168 cm_thickness = camera.pop("thickness", None) 2169 cm_view_angle = camera.pop("view_angle", camera.pop("viewAngle", None)) 2170 cm_roll = camera.pop("roll", None) 2171 2172 if len(camera.keys()) > 0: 2173 vedo.logger.warning(f"in camera_from_dict, key(s) not recognized: {camera.keys()}") 2174 if cm_pos is not None: vcam.SetPosition(cm_pos) 2175 if cm_focal_point is not None: vcam.SetFocalPoint(cm_focal_point) 2176 if cm_viewup is not None: vcam.SetViewUp(cm_viewup) 2177 if cm_distance is not None: vcam.SetDistance(cm_distance) 2178 if cm_clipping_range is not None: vcam.SetClippingRange(cm_clipping_range) 2179 if cm_parallel_scale is not None: vcam.SetParallelScale(cm_parallel_scale) 2180 if cm_thickness is not None: vcam.SetThickness(cm_thickness) 2181 if cm_view_angle is not None: vcam.SetViewAngle(cm_view_angle) 2182 if cm_roll is not None: vcam.SetRoll(cm_roll) 2183 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:
2185def camera_to_dict(vtkcam) -> dict: 2186 """ 2187 Convert a [vtkCamera](https://vtk.org/doc/nightly/html/classvtkCamera.html) 2188 object into a python dictionary. 2189 2190 Parameters of the camera are: 2191 - `position` (3-tuple) 2192 - `focal_point` (3-tuple) 2193 - `viewup` (3-tuple) 2194 - `distance` (float) 2195 - `clipping_range` (2-tuple) 2196 - `parallel_scale` (float) 2197 - `thickness` (float) 2198 - `view_angle` (float) 2199 - `roll` (float) 2200 2201 Arguments: 2202 vtkcam: (vtkCamera) 2203 a `vtkCamera` object to convert. 2204 """ 2205 cam = dict() 2206 cam["position"] = np.array(vtkcam.GetPosition()) 2207 cam["focal_point"] = np.array(vtkcam.GetFocalPoint()) 2208 cam["viewup"] = np.array(vtkcam.GetViewUp()) 2209 cam["distance"] = vtkcam.GetDistance() 2210 cam["clipping_range"] = np.array(vtkcam.GetClippingRange()) 2211 cam["parallel_scale"] = vtkcam.GetParallelScale() 2212 cam["thickness"] = vtkcam.GetThickness() 2213 cam["view_angle"] = vtkcam.GetViewAngle() 2214 cam["roll"] = vtkcam.GetRoll() 2215 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=(0, 0, 0), up_vector=(0, 1, 0), backoff_vector=(0, 0, 1), backoff=1.0) -> vtkmodules.vtkRenderingCore.vtkCamera:
2113def oriented_camera(center=(0, 0, 0), up_vector=(0, 1, 0), backoff_vector=(0, 0, 1), backoff=1.0) -> vtki.vtkCamera: 2114 """ 2115 Generate a `vtkCamera` pointed at a specific location, 2116 oriented with a given up direction, set to a backoff. 2117 """ 2118 vup = np.array(up_vector) 2119 vup = vup / np.linalg.norm(vup) 2120 pt_backoff = center - backoff * np.array(backoff_vector) 2121 camera = vtki.vtkCamera() 2122 camera.SetFocalPoint(center[0], center[1], center[2]) 2123 camera.SetViewUp(vup[0], vup[1], vup[2]) 2124 camera.SetPosition(pt_backoff[0], pt_backoff[1], pt_backoff[2]) 2125 return camera
Generate a vtkCamera pointed at a specific location,
oriented with a given up direction, set to a backoff.
def
vedo2trimesh(mesh):
2441def vedo2trimesh(mesh): 2442 """ 2443 Convert `vedo.mesh.Mesh` to `Trimesh.Mesh` object. 2444 """ 2445 if is_sequence(mesh): 2446 tms = [] 2447 for a in mesh: 2448 tms.append(vedo2trimesh(a)) 2449 return tms 2450 2451 try: 2452 from trimesh import Trimesh 2453 except ModuleNotFoundError: 2454 vedo.logger.error("Need trimesh to run:\npip install trimesh") 2455 return None 2456 2457 tris = mesh.cells 2458 carr = mesh.celldata["CellIndividualColors"] 2459 ccols = carr 2460 2461 points = mesh.vertices 2462 varr = mesh.pointdata["VertexColors"] 2463 vcols = varr 2464 2465 if len(tris) == 0: 2466 tris = None 2467 2468 return Trimesh(vertices=points, faces=tris, face_colors=ccols, vertex_colors=vcols)
Convert vedo.mesh.Mesh to Trimesh.Mesh object.
def
trimesh2vedo(inputobj):
2471def trimesh2vedo(inputobj): 2472 """ 2473 Convert a `Trimesh` object to `vedo.Mesh` or `vedo.Assembly` object. 2474 """ 2475 if is_sequence(inputobj): 2476 vms = [] 2477 for ob in inputobj: 2478 vms.append(trimesh2vedo(ob)) 2479 return vms 2480 2481 inputobj_type = str(type(inputobj)) 2482 2483 if "Trimesh" in inputobj_type or "primitives" in inputobj_type: 2484 faces = inputobj.faces 2485 poly = buildPolyData(inputobj.vertices, faces) 2486 tact = vedo.Mesh(poly) 2487 if inputobj.visual.kind == "face": 2488 trim_c = inputobj.visual.face_colors 2489 elif inputobj.visual.kind == "texture": 2490 trim_c = inputobj.visual.to_color().vertex_colors 2491 else: 2492 trim_c = inputobj.visual.vertex_colors 2493 2494 if is_sequence(trim_c): 2495 if is_sequence(trim_c[0]): 2496 same_color = len(np.unique(trim_c, axis=0)) < 2 # all vtxs have same color 2497 2498 if same_color: 2499 tact.c(trim_c[0, [0, 1, 2]]).alpha(trim_c[0, 3]) 2500 else: 2501 if inputobj.visual.kind == "face": 2502 tact.cellcolors = trim_c 2503 return tact 2504 2505 if "PointCloud" in inputobj_type: 2506 2507 vdpts = vedo.shapes.Points(inputobj.vertices, r=8, c='k') 2508 if hasattr(inputobj, "vertices_color"): 2509 vcols = (inputobj.vertices_color * 1).astype(np.uint8) 2510 vdpts.pointcolors = vcols 2511 return vdpts 2512 2513 if "path" in inputobj_type: 2514 2515 lines = [] 2516 for e in inputobj.entities: 2517 # print('trimesh entity', e.to_dict()) 2518 l = vedo.shapes.Line(inputobj.vertices[e.points], c="k", lw=2) 2519 lines.append(l) 2520 return vedo.Assembly(lines) 2521 2522 return None
Convert a Trimesh object to vedo.Mesh or vedo.Assembly object.
def
vedo2meshlab(vmesh):
2525def vedo2meshlab(vmesh): 2526 """Convert a `vedo.Mesh` to a Meshlab object.""" 2527 try: 2528 import pymeshlab as mlab 2529 except ModuleNotFoundError: 2530 vedo.logger.error("Need pymeshlab to run:\npip install pymeshlab") 2531 2532 vertex_matrix = vmesh.vertices.astype(np.float64) 2533 2534 try: 2535 face_matrix = np.asarray(vmesh.cells, dtype=np.float64) 2536 except: 2537 print("WARNING: in vedo2meshlab(), need to triangulate mesh first!") 2538 face_matrix = np.array(vmesh.clone().triangulate().cells, dtype=np.float64) 2539 2540 # v_normals_matrix = vmesh.normals(cells=False, recompute=False) 2541 v_normals_matrix = vmesh.vertex_normals 2542 if not v_normals_matrix.shape[0]: 2543 v_normals_matrix = np.empty((0, 3), dtype=np.float64) 2544 2545 # f_normals_matrix = vmesh.normals(cells=True, recompute=False) 2546 f_normals_matrix = vmesh.cell_normals 2547 if not f_normals_matrix.shape[0]: 2548 f_normals_matrix = np.empty((0, 3), dtype=np.float64) 2549 2550 v_color_matrix = vmesh.pointdata["RGBA"] 2551 if v_color_matrix is None: 2552 v_color_matrix = np.empty((0, 4), dtype=np.float64) 2553 else: 2554 v_color_matrix = v_color_matrix.astype(np.float64) / 255 2555 if v_color_matrix.shape[1] == 3: 2556 v_color_matrix = np.c_[ 2557 v_color_matrix, np.ones(v_color_matrix.shape[0], dtype=np.float64) 2558 ] 2559 2560 f_color_matrix = vmesh.celldata["RGBA"] 2561 if f_color_matrix is None: 2562 f_color_matrix = np.empty((0, 4), dtype=np.float64) 2563 else: 2564 f_color_matrix = f_color_matrix.astype(np.float64) / 255 2565 if f_color_matrix.shape[1] == 3: 2566 f_color_matrix = np.c_[ 2567 f_color_matrix, np.ones(f_color_matrix.shape[0], dtype=np.float64) 2568 ] 2569 2570 m = mlab.Mesh( 2571 vertex_matrix=vertex_matrix, 2572 face_matrix=face_matrix, 2573 v_normals_matrix=v_normals_matrix, 2574 f_normals_matrix=f_normals_matrix, 2575 v_color_matrix=v_color_matrix, 2576 f_color_matrix=f_color_matrix, 2577 ) 2578 2579 for k in vmesh.pointdata.keys(): 2580 data = vmesh.pointdata[k] 2581 if data is not None: 2582 if data.ndim == 1: # scalar 2583 m.add_vertex_custom_scalar_attribute(data.astype(np.float64), k) 2584 elif data.ndim == 2: # vectorial data 2585 if "tcoord" not in k.lower() and k not in ["Normals", "TextureCoordinates"]: 2586 m.add_vertex_custom_point_attribute(data.astype(np.float64), k) 2587 2588 for k in vmesh.celldata.keys(): 2589 data = vmesh.celldata[k] 2590 if data is not None: 2591 if data.ndim == 1: # scalar 2592 m.add_face_custom_scalar_attribute(data.astype(np.float64), k) 2593 elif data.ndim == 2 and k != "Normals": # vectorial data 2594 m.add_face_custom_point_attribute(data.astype(np.float64), k) 2595 2596 m.update_bounding_box() 2597 return m
Convert a vedo.Mesh to a Meshlab object.
def
meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()):
2600def meshlab2vedo(mmesh, pointdata_keys=(), celldata_keys=()): 2601 """Convert a Meshlab object to `vedo.Mesh`.""" 2602 inputtype = str(type(mmesh)) 2603 2604 if "MeshSet" in inputtype: 2605 mmesh = mmesh.current_mesh() 2606 2607 mpoints, mcells = mmesh.vertex_matrix(), mmesh.face_matrix() 2608 if len(mcells) > 0: 2609 polydata = buildPolyData(mpoints, mcells) 2610 else: 2611 polydata = buildPolyData(mpoints, None) 2612 2613 if mmesh.has_vertex_scalar(): 2614 parr = mmesh.vertex_scalar_array() 2615 parr_vtk = numpy_to_vtk(parr) 2616 parr_vtk.SetName("MeshLabScalars") 2617 polydata.GetPointData().AddArray(parr_vtk) 2618 polydata.GetPointData().SetActiveScalars("MeshLabScalars") 2619 2620 if mmesh.has_face_scalar(): 2621 carr = mmesh.face_scalar_array() 2622 carr_vtk = numpy_to_vtk(carr) 2623 carr_vtk.SetName("MeshLabScalars") 2624 x0, x1 = carr_vtk.GetRange() 2625 polydata.GetCellData().AddArray(carr_vtk) 2626 polydata.GetCellData().SetActiveScalars("MeshLabScalars") 2627 2628 for k in pointdata_keys: 2629 parr = mmesh.vertex_custom_scalar_attribute_array(k) 2630 parr_vtk = numpy_to_vtk(parr) 2631 parr_vtk.SetName(k) 2632 polydata.GetPointData().AddArray(parr_vtk) 2633 polydata.GetPointData().SetActiveScalars(k) 2634 2635 for k in celldata_keys: 2636 carr = mmesh.face_custom_scalar_attribute_array(k) 2637 carr_vtk = numpy_to_vtk(carr) 2638 carr_vtk.SetName(k) 2639 polydata.GetCellData().AddArray(carr_vtk) 2640 polydata.GetCellData().SetActiveScalars(k) 2641 2642 pnorms = mmesh.vertex_normal_matrix() 2643 if len(pnorms) > 0: 2644 polydata.GetPointData().SetNormals(numpy2vtk(pnorms, name="Normals")) 2645 2646 cnorms = mmesh.face_normal_matrix() 2647 if len(cnorms) > 0: 2648 polydata.GetCellData().SetNormals(numpy2vtk(cnorms, name="Normals")) 2649 return vedo.Mesh(polydata)
Convert a Meshlab object to vedo.Mesh.
def
vedo2open3d(vedo_mesh):
2660def vedo2open3d(vedo_mesh): 2661 """ 2662 Return an `open3d.geometry.TriangleMesh` version of the current mesh. 2663 """ 2664 try: 2665 import open3d as o3d 2666 except RuntimeError: 2667 vedo.logger.error("Need open3d to run:\npip install open3d") 2668 2669 # create from numpy arrays 2670 o3d_mesh = o3d.geometry.TriangleMesh( 2671 vertices=o3d.utility.Vector3dVector(vedo_mesh.vertices), 2672 triangles=o3d.utility.Vector3iVector(vedo_mesh.cells), 2673 ) 2674 # TODO: need to add some if check here in case color and normals 2675 # info are not existing 2676 # o3d_mesh.vertex_colors = o3d.utility.Vector3dVector(vedo_mesh.pointdata["RGB"]/255) 2677 # o3d_mesh.vertex_normals= o3d.utility.Vector3dVector(vedo_mesh.pointdata["Normals"]) 2678 return o3d_mesh
Return an open3d.geometry.TriangleMesh version of the current mesh.
def
open3d2vedo(o3d_mesh):
2652def open3d2vedo(o3d_mesh): 2653 """Convert `open3d.geometry.TriangleMesh` to a `vedo.Mesh`.""" 2654 m = vedo.Mesh([np.array(o3d_mesh.vertices), np.array(o3d_mesh.triangles)]) 2655 # TODO: could also check whether normals and color are present in 2656 # order to port with the above vertices/faces 2657 return m
Convert open3d.geometry.TriangleMesh to a vedo.Mesh.
def
vtk2numpy(varr):
792def vtk2numpy(varr): 793 """Convert a `vtkDataArray`, `vtkIdList` or `vtTransform` into a numpy array.""" 794 if varr is None: 795 return np.array([]) 796 if isinstance(varr, vtki.vtkIdList): 797 return np.array([varr.GetId(i) for i in range(varr.GetNumberOfIds())]) 798 elif isinstance(varr, vtki.vtkBitArray): 799 carr = vtki.vtkCharArray() 800 carr.DeepCopy(varr) 801 varr = carr 802 elif isinstance(varr, vtki.vtkHomogeneousTransform): 803 try: 804 varr = varr.GetMatrix() 805 except AttributeError: 806 pass 807 M = [[varr.GetElement(i, j) for j in range(4)] for i in range(4)] 808 return np.array(M) 809 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 varr = numpy_to_vtkIdTypeArray(arr.astype(np.int64), deep=deep) 780 elif dtype: 781 varr = numpy_to_vtk(arr.astype(dtype), deep=deep) 782 else: 783 # let numpy_to_vtk() decide what is best type based on arr type 784 if arr.dtype == np.bool_: 785 arr = arr.astype(np.uint8) 786 varr = numpy_to_vtk(arr, deep=deep) 787 788 if name: 789 varr.SetName(name) 790 return varr
Convert a numpy array into a vtkDataArray.
Use dtype='id' for vtkIdTypeArray objects.
def
get_uv(p, x, v):
1450def get_uv(p, x, v): 1451 """ 1452 Obtain the texture uv-coords of a point p belonging to a face that has point 1453 coordinates (x0, x1, x2) with the corresponding uv-coordinates v=(v0, v1, v2). 1454 All p and x0,x1,x2 are 3D-vectors, while v are their 2D uv-coordinates. 1455 1456 Example: 1457 ```python 1458 from vedo import * 1459 1460 pic = Image(dataurl+"coloured_cube_faces.jpg") 1461 cb = Mesh(dataurl+"coloured_cube.obj").lighting("off").texture(pic) 1462 1463 cbpts = cb.vertices 1464 faces = cb.cells 1465 uv = cb.pointdata["Material"] 1466 1467 pt = [-0.2, 0.75, 2] 1468 pr = cb.closest_point(pt) 1469 1470 idface = cb.closest_point(pt, return_cell_id=True) 1471 idpts = faces[idface] 1472 uv_face = uv[idpts] 1473 1474 uv_pr = utils.get_uv(pr, cbpts[idpts], uv_face) 1475 print("interpolated uv =", uv_pr) 1476 1477 sx, sy = pic.dimensions() 1478 i_interp_uv = uv_pr * [sy, sx] 1479 ix, iy = i_interp_uv.astype(int) 1480 mpic = pic.tomesh() 1481 rgba = mpic.pointdata["RGBA"].reshape(sy, sx, 3) 1482 print("color =", rgba[ix, iy]) 1483 1484 show( 1485 [[cb, Point(pr), cb.labels("Material")], 1486 [pic, Point(i_interp_uv)]], 1487 N=2, axes=1, sharecam=False, 1488 ).close() 1489 ``` 1490  1491 """ 1492 # Vector vp=p-x0 is representable as alpha*s + beta*t, 1493 # where s = x1-x0 and t = x2-x0, in matrix form 1494 # vp = [alpha, beta] . matrix(s,t) 1495 # M = matrix(s,t) is 2x3 matrix, so (alpha, beta) can be found by 1496 # inverting any of its minor A with non-zero determinant. 1497 # Once found, uv-coords of p are vt0 + alpha (vt1-v0) + beta (vt2-v0) 1498 1499 p = np.asarray(p) 1500 x0, x1, x2 = np.asarray(x)[:3] 1501 vt0, vt1, vt2 = np.asarray(v)[:3] 1502 1503 s = x1 - x0 1504 t = x2 - x0 1505 vs = vt1 - vt0 1506 vt = vt2 - vt0 1507 vp = p - x0 1508 1509 # finding a minor with independent rows 1510 M = np.matrix([s, t]) 1511 mnr = [0, 1] 1512 A = M[:, mnr] 1513 if np.abs(np.linalg.det(A)) < 0.000001: 1514 mnr = [0, 2] 1515 A = M[:, mnr] 1516 if np.abs(np.linalg.det(A)) < 0.000001: 1517 mnr = [1, 2] 1518 A = M[:, mnr] 1519 Ainv = np.linalg.inv(A) 1520 alpha_beta = vp[mnr].dot(Ainv) # [alpha, beta] 1521 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.
