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