vedo.utils

Utilities submodule.

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

Print the tree of operations.

def show(self, orientation='LR', popup=True):
215    def show(self, orientation="LR", popup=True):
216        """Show the graphviz output for the pipeline of this object"""
217        if not vedo.settings.enable_pipeline:
218            return
219
220        try:
221            from graphviz import Digraph
222        except ImportError:
223            vedo.logger.error("please install graphviz with command\n pip install graphviz")
224            return
225
226        # visualize the entire tree
227        dot = Digraph(
228            node_attr={"fontcolor": "#201010", "fontname": "Helvetica", "fontsize": "12"},
229            edge_attr={"fontname": "Helvetica", "fontsize": "6", "arrowsize": "0.4"},
230        )
231        dot.attr(rankdir=orientation)
232
233        self.counts = 0
234        self._build_tree(dot)
235        self.dot = dot
236
237        home_dir = os.path.expanduser("~")
238        gpath = os.path.join(
239            home_dir, vedo.settings.cache_directory, "vedo", "pipeline_graphviz")
240
241        dot.render(gpath, view=popup)

Show the graphviz output for the pipeline of this object

class ProgressBar:
245class ProgressBar:
246    """
247    Class to print a progress bar.
248    """
249
250    def __init__(
251        self,
252        start,
253        stop,
254        step=1,
255        c=None,
256        bold=True,
257        italic=False,
258        title="",
259        eta=True,
260        delay=-1,
261        width=25,
262        char="\U00002501",
263        char_back="\U00002500",
264    ):
265        """
266        Class to print a progress bar with optional text message.
267
268        Check out also function `progressbar()`.
269
270        Arguments:
271            start : (int)
272                starting value
273            stop : (int)
274                stopping value
275            step : (int)
276                step value
277            c : (str)
278                color in hex format
279            title : (str)
280                title text
281            eta : (bool)
282                estimate time of arrival
283            delay : (float)
284                minimum time before printing anything,
285                if negative use the default value
286                as set in `vedo.settings.progressbar_delay`
287            width : (int)
288                width of the progress bar
289            char : (str)
290                character to use for the progress bar
291            char_back : (str)
292                character to use for the background of the progress bar
293
294        Example:
295            ```python
296            import time
297            from vedo import ProgressBar
298            pb = ProgressBar(0,40, c='r')
299            for i in pb.range():
300                time.sleep(0.1)
301                pb.print()
302            ```
303            ![](https://user-images.githubusercontent.com/32848391/51858823-ed1f4880-2335-11e9-8788-2d102ace2578.png)
304        """
305        self.char = char
306        self.char_back = char_back
307
308        self.title = title + " "
309        if title:
310            self.title = " " + self.title
311
312        if delay < 0:
313            delay = vedo.settings.progressbar_delay
314
315        self.start = start
316        self.stop = stop
317        self.step = step
318
319        self.color = c
320        self.bold = bold
321        self.italic = italic
322        self.width = width
323        self.pbar = ""
324        self.percent = 0.0
325        self.percent_int = 0
326        self.eta = eta
327        self.delay = delay
328
329        self.t0 = time.time()
330        self._remaining = 1e10
331
332        self._update(0)
333
334        self._counts = 0
335        self._oldbar = ""
336        self._lentxt = 0
337        self._range = np.arange(start, stop, step)
338
339    def print(self, txt="", c=None):
340        """Print the progress bar with an optional message."""
341        if not c:
342            c = self.color
343
344        self._update(self._counts + self.step)
345
346        if self.delay:
347            if time.time() - self.t0 < self.delay:
348                return
349
350        if self.pbar != self._oldbar:
351            self._oldbar = self.pbar
352
353            if self.eta and self._counts > 1:
354
355                tdenom = time.time() - self.t0
356                if tdenom:
357                    vel = self._counts / tdenom
358                    self._remaining = (self.stop - self._counts) / vel
359                else:
360                    vel = 1
361                    self._remaining = 0.0
362
363                if self._remaining > 60:
364                    mins = int(self._remaining / 60)
365                    secs = self._remaining - 60 * mins
366                    mins = f"{mins}m"
367                    secs = f"{int(secs + 0.5)}s "
368                else:
369                    mins = ""
370                    secs = f"{int(self._remaining + 0.5)}s "
371
372                vel = round(vel, 1)
373                eta = f"eta: {mins}{secs}({vel} it/s) "
374                if self._remaining < 0.5:
375                    dt = time.time() - self.t0
376                    if dt > 60:
377                        mins = int(dt / 60)
378                        secs = dt - 60 * mins
379                        mins = f"{mins}m"
380                        secs = f"{int(secs + 0.5)}s "
381                    else:
382                        mins = ""
383                        secs = f"{int(dt + 0.5)}s "
384                    eta = f"elapsed: {mins}{secs}({vel} it/s)        "
385                    txt = ""
386            else:
387                eta = ""
388
389            eraser = " " * self._lentxt + "\b" * self._lentxt
390
391            s = f"{self.pbar} {eraser}{eta}{txt}\r"
392            vedo.printc(s, c=c, bold=self.bold, italic=self.italic, end="")
393            if self.percent > 99.999:
394                print("")
395
396            self._lentxt = len(txt)
397
398    def range(self):
399        """Return the range iterator."""
400        return self._range
401
402    def _update(self, counts):
403        if counts < self.start:
404            counts = self.start
405        elif counts > self.stop:
406            counts = self.stop
407        self._counts = counts
408
409        self.percent = (self._counts - self.start) * 100.0
410
411        delta = self.stop - self.start
412        if delta:
413            self.percent /= delta
414        else:
415            self.percent = 0.0
416
417        self.percent_int = int(round(self.percent))
418        af = self.width - 2
419        nh = int(round(self.percent_int / 100 * af))
420        pbar_background = "\x1b[2m" + self.char_back * (af - nh)
421        self.pbar = f"{self.title}{self.char * (nh-1)}{pbar_background}"
422        if self.percent < 100.0:
423            ps = f" {self.percent_int}%"
424        else:
425            ps = ""
426        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='─')
250    def __init__(
251        self,
252        start,
253        stop,
254        step=1,
255        c=None,
256        bold=True,
257        italic=False,
258        title="",
259        eta=True,
260        delay=-1,
261        width=25,
262        char="\U00002501",
263        char_back="\U00002500",
264    ):
265        """
266        Class to print a progress bar with optional text message.
267
268        Check out also function `progressbar()`.
269
270        Arguments:
271            start : (int)
272                starting value
273            stop : (int)
274                stopping value
275            step : (int)
276                step value
277            c : (str)
278                color in hex format
279            title : (str)
280                title text
281            eta : (bool)
282                estimate time of arrival
283            delay : (float)
284                minimum time before printing anything,
285                if negative use the default value
286                as set in `vedo.settings.progressbar_delay`
287            width : (int)
288                width of the progress bar
289            char : (str)
290                character to use for the progress bar
291            char_back : (str)
292                character to use for the background of the progress bar
293
294        Example:
295            ```python
296            import time
297            from vedo import ProgressBar
298            pb = ProgressBar(0,40, c='r')
299            for i in pb.range():
300                time.sleep(0.1)
301                pb.print()
302            ```
303            ![](https://user-images.githubusercontent.com/32848391/51858823-ed1f4880-2335-11e9-8788-2d102ace2578.png)
304        """
305        self.char = char
306        self.char_back = char_back
307
308        self.title = title + " "
309        if title:
310            self.title = " " + self.title
311
312        if delay < 0:
313            delay = vedo.settings.progressbar_delay
314
315        self.start = start
316        self.stop = stop
317        self.step = step
318
319        self.color = c
320        self.bold = bold
321        self.italic = italic
322        self.width = width
323        self.pbar = ""
324        self.percent = 0.0
325        self.percent_int = 0
326        self.eta = eta
327        self.delay = delay
328
329        self.t0 = time.time()
330        self._remaining = 1e10
331
332        self._update(0)
333
334        self._counts = 0
335        self._oldbar = ""
336        self._lentxt = 0
337        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):
339    def print(self, txt="", c=None):
340        """Print the progress bar with an optional message."""
341        if not c:
342            c = self.color
343
344        self._update(self._counts + self.step)
345
346        if self.delay:
347            if time.time() - self.t0 < self.delay:
348                return
349
350        if self.pbar != self._oldbar:
351            self._oldbar = self.pbar
352
353            if self.eta and self._counts > 1:
354
355                tdenom = time.time() - self.t0
356                if tdenom:
357                    vel = self._counts / tdenom
358                    self._remaining = (self.stop - self._counts) / vel
359                else:
360                    vel = 1
361                    self._remaining = 0.0
362
363                if self._remaining > 60:
364                    mins = int(self._remaining / 60)
365                    secs = self._remaining - 60 * mins
366                    mins = f"{mins}m"
367                    secs = f"{int(secs + 0.5)}s "
368                else:
369                    mins = ""
370                    secs = f"{int(self._remaining + 0.5)}s "
371
372                vel = round(vel, 1)
373                eta = f"eta: {mins}{secs}({vel} it/s) "
374                if self._remaining < 0.5:
375                    dt = time.time() - self.t0
376                    if dt > 60:
377                        mins = int(dt / 60)
378                        secs = dt - 60 * mins
379                        mins = f"{mins}m"
380                        secs = f"{int(secs + 0.5)}s "
381                    else:
382                        mins = ""
383                        secs = f"{int(dt + 0.5)}s "
384                    eta = f"elapsed: {mins}{secs}({vel} it/s)        "
385                    txt = ""
386            else:
387                eta = ""
388
389            eraser = " " * self._lentxt + "\b" * self._lentxt
390
391            s = f"{self.pbar} {eraser}{eta}{txt}\r"
392            vedo.printc(s, c=c, bold=self.bold, italic=self.italic, end="")
393            if self.percent > 99.999:
394                print("")
395
396            self._lentxt = len(txt)

Print the progress bar with an optional message.

def range(self):
398    def range(self):
399        """Return the range iterator."""
400        return self._range

Return the range iterator.

def progressbar( iterable, c=None, bold=True, italic=False, title='', eta=True, width=25, delay=-1):
430def progressbar(
431        iterable,
432        c=None, bold=True, italic=False, title="",
433        eta=True, width=25, delay=-1,
434    ):
435    """
436    Function to print a progress bar with optional text message.
437
438    Use delay to set a minimum time before printing anything.
439    If delay is negative, then use the default value
440    as set in `vedo.settings.progressbar_delay`.
441
442    Arguments:
443        start : (int)
444            starting value
445        stop : (int)
446            stopping value
447        step : (int)
448            step value
449        c : (str)
450            color in hex format
451        title : (str)
452            title text
453        eta : (bool)
454            estimate time of arrival
455        delay : (float)
456            minimum time before printing anything,
457            if negative use the default value
458            set in `vedo.settings.progressbar_delay`
459        width : (int)
460            width of the progress bar
461        char : (str)
462            character to use for the progress bar
463        char_back : (str)
464            character to use for the background of the progress bar
465
466    Example:
467        ```python
468        import time
469        for i in progressbar(range(100), c='r'):
470            time.sleep(0.1)
471        ```
472        ![](https://user-images.githubusercontent.com/32848391/51858823-ed1f4880-2335-11e9-8788-2d102ace2578.png)
473    """
474    try:
475        if is_number(iterable):
476            total = int(iterable)
477            iterable = range(total)
478        else:
479            total = len(iterable)
480    except TypeError:
481        iterable = list(iterable)
482        total = len(iterable)
483
484    pb = ProgressBar(
485        0, total, c=c, bold=bold, italic=italic, title=title,
486        eta=eta, delay=delay, width=width,
487    )
488    for item in iterable:
489        pb.print()
490        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:
494class Minimizer:
495    """
496    A function minimizer that uses the Nelder-Mead method.
497
498    The algorithm constructs an n-dimensional simplex in parameter
499    space (i.e. a tetrahedron if the number or parameters is 3)
500    and moves the vertices around parameter space until
501    a local minimum is found. The amoeba method is robust,
502    reasonably efficient, but is not guaranteed to find
503    the global minimum if several local minima exist.
504    
505    Arguments:
506        function : (callable)
507            the function to minimize
508        max_iterations : (int)
509            the maximum number of iterations
510        contraction_ratio : (float)
511            The contraction ratio.
512            The default value of 0.5 gives fast convergence,
513            but larger values such as 0.6 or 0.7 provide greater stability.
514        expansion_ratio : (float)
515            The expansion ratio.
516            The default value is 2.0, which provides rapid expansion.
517            Values between 1.1 and 2.0 are valid.
518        tol : (float)
519            the tolerance for convergence
520    
521    Example:
522        - [nelder-mead.py](https://github.com/marcomusy/vedo/blob/master/examples/others/nelder-mead.py)
523    """
524    def __init__(
525            self,
526            function=None,
527            max_iterations=10000,
528            contraction_ratio=0.5,
529            expansion_ratio=2.0,
530            tol=1e-5,
531        ):
532        self.function = function
533        self.tolerance = tol
534        self.contraction_ratio = contraction_ratio
535        self.expansion_ratio = expansion_ratio
536        self.max_iterations = max_iterations
537        self.minimizer = vtki.new("AmoebaMinimizer")
538        self.minimizer.SetFunction(self._vtkfunc)
539        self.results = {}
540        self.parameters_path = []
541        self.function_path = []
542
543    def _vtkfunc(self):
544        n = self.minimizer.GetNumberOfParameters()
545        ain = [self.minimizer.GetParameterValue(i) for i in range(n)]
546        r = self.function(ain)
547        self.minimizer.SetFunctionValue(r)
548        self.parameters_path.append(ain)
549        self.function_path.append(r)
550        return r
551    
552    def eval(self, parameters=()):
553        """
554        Evaluate the function at the current or given parameters.
555        """
556        if len(parameters) == 0:
557            return self.minimizer.EvaluateFunction()
558        self.set_parameters(parameters)
559        return self.function(parameters)
560    
561    def set_parameter(self, name, value, scale=1.0):
562        """
563        Set the parameter value.
564        The initial amount by which the parameter
565        will be modified during the search for the minimum. 
566        """
567        self.minimizer.SetParameterValue(name, value)
568        self.minimizer.SetParameterScale(name, scale)
569    
570    def set_parameters(self, parameters):
571        """
572        Set the parameters names and values from a dictionary.
573        """
574        for name, value in parameters.items():
575            if len(value) == 2:
576                self.set_parameter(name, value[0], value[1])
577            else:
578                self.set_parameter(name, value)
579    
580    def minimize(self):
581        """
582        Minimize the input function.
583
584        Returns:
585            dict : 
586                the minimization results
587            init_parameters : (dict)
588                the initial parameters
589            parameters : (dict)
590                the final parameters
591            min_value : (float)
592                the minimum value
593            iterations : (int)
594                the number of iterations
595            max_iterations : (int)
596                the maximum number of iterations
597            tolerance : (float)
598                the tolerance for convergence
599            convergence_flag : (int)
600                zero if the tolerance stopping
601                criterion has been met.
602            parameters_path : (np.array)
603                the path of the minimization
604                algorithm in parameter space
605            function_path : (np.array)
606                the path of the minimization
607                algorithm in function space
608            hessian : (np.array)
609                the Hessian matrix of the
610                function at the minimum
611            parameter_errors : (np.array)
612                the errors on the parameters
613        """
614        n = self.minimizer.GetNumberOfParameters()
615        out = [(
616            self.minimizer.GetParameterName(i),
617            (self.minimizer.GetParameterValue(i),
618             self.minimizer.GetParameterScale(i))
619        ) for i in range(n)]
620        self.results["init_parameters"] = dict(out)
621
622        self.minimizer.SetTolerance(self.tolerance)
623        self.minimizer.SetContractionRatio(self.contraction_ratio)
624        self.minimizer.SetExpansionRatio(self.expansion_ratio)
625        self.minimizer.SetMaxIterations(self.max_iterations)
626
627        self.minimizer.Minimize()
628        self.results["convergence_flag"] = not bool(self.minimizer.Iterate())
629
630        out = [(
631            self.minimizer.GetParameterName(i),
632            self.minimizer.GetParameterValue(i),
633        ) for i in range(n)]
634
635        self.results["parameters"] = dict(out)
636        self.results["min_value"] = self.minimizer.GetFunctionValue()
637        self.results["iterations"] = self.minimizer.GetIterations()
638        self.results["max_iterations"] = self.minimizer.GetMaxIterations()
639        self.results["tolerance"] = self.minimizer.GetTolerance()
640        self.results["expansion_ratio"] = self.expansion_ratio
641        self.results["contraction_ratio"] = self.contraction_ratio
642        self.results["parameters_path"] = np.array(self.parameters_path)
643        self.results["function_path"] = np.array(self.function_path)
644        self.results["hessian"] = np.zeros((n,n))
645        self.results["parameter_errors"] = np.zeros(n)
646        return self.results
647
648    def compute_hessian(self, epsilon=0):
649        """
650        Compute the Hessian matrix of `function` at the
651        minimum numerically.
652
653        Arguments:
654            epsilon : (float)
655                Step size used for numerical approximation.
656
657        Returns:
658            array: Hessian matrix of `function` at minimum.
659        """
660        if not epsilon:
661            epsilon = self.tolerance * 10
662        n = self.minimizer.GetNumberOfParameters()
663        x0 = [self.minimizer.GetParameterValue(i) for i in range(n)]
664        hessian = np.zeros((n, n))
665        for i in vedo.progressbar(n, title="Computing Hessian", delay=2):
666            for j in range(n):
667                xijp = np.copy(x0)
668                xijp[i] += epsilon
669                xijp[j] += epsilon
670                xijm = np.copy(x0)
671                xijm[i] += epsilon
672                xijm[j] -= epsilon
673                xjip = np.copy(x0)
674                xjip[i] -= epsilon
675                xjip[j] += epsilon
676                xjim = np.copy(x0)
677                xjim[i] -= epsilon
678                xjim[j] -= epsilon
679                # Second derivative approximation
680                fijp = self.function(xijp)
681                fijm = self.function(xijm)
682                fjip = self.function(xjip)
683                fjim = self.function(xjim)
684                hessian[i, j] = (fijp - fijm - fjip + fjim) / (2 * epsilon**2)
685        self.results["hessian"] = hessian
686        try:
687            ihess = np.linalg.inv(hessian)
688            self.results["parameter_errors"] = np.sqrt(np.diag(ihess))
689        except:
690            vedo.logger.warning("Cannot compute hessian for parameter errors")
691            self.results["parameter_errors"] = np.zeros(n)
692        return hessian
693
694    def __str__(self) -> str:
695        out = vedo.printc(
696            f"vedo.utils.Minimizer at ({hex(id(self))})".ljust(75),
697            bold=True, invert=True, return_string=True,
698        )
699        out += "Function name".ljust(20) + self.function.__name__ + "()\n"
700        out += "-------- parameters initial value -----------\n"
701        out += "Name".ljust(20) + "Value".ljust(20) + "Scale\n"
702        for name, value in self.results["init_parameters"].items():
703            out += name.ljust(20) + str(value[0]).ljust(20) + str(value[1]) + "\n"
704        out += "-------- parameters final value --------------\n"
705        for name, value in self.results["parameters"].items():
706            out += name.ljust(20) + f"{value:.6f}"
707            ierr = list(self.results["parameters"]).index(name)
708            err = self.results["parameter_errors"][ierr]
709            if err:
710                out += f" ± {err:.4f}"
711            out += "\n"
712        out += "Value at minimum".ljust(20)+ f'{self.results["min_value"]}\n'
713        out += "Iterations".ljust(20)      + f'{self.results["iterations"]}\n'
714        out += "Max iterations".ljust(20)  + f'{self.results["max_iterations"]}\n'
715        out += "Convergence flag".ljust(20)+ f'{self.results["convergence_flag"]}\n'
716        out += "Tolerance".ljust(20)       + f'{self.results["tolerance"]}\n'
717        try:
718            arr = np.array2string(
719                self.compute_hessian(),
720                separator=', ', precision=6, suppress_small=True,
721            )
722            out += "Hessian Matrix:\n" + arr
723        except:
724            out += "Hessian Matrix: (not available)"
725        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)
524    def __init__(
525            self,
526            function=None,
527            max_iterations=10000,
528            contraction_ratio=0.5,
529            expansion_ratio=2.0,
530            tol=1e-5,
531        ):
532        self.function = function
533        self.tolerance = tol
534        self.contraction_ratio = contraction_ratio
535        self.expansion_ratio = expansion_ratio
536        self.max_iterations = max_iterations
537        self.minimizer = vtki.new("AmoebaMinimizer")
538        self.minimizer.SetFunction(self._vtkfunc)
539        self.results = {}
540        self.parameters_path = []
541        self.function_path = []
def eval(self, parameters=()):
552    def eval(self, parameters=()):
553        """
554        Evaluate the function at the current or given parameters.
555        """
556        if len(parameters) == 0:
557            return self.minimizer.EvaluateFunction()
558        self.set_parameters(parameters)
559        return self.function(parameters)

Evaluate the function at the current or given parameters.

def set_parameter(self, name, value, scale=1.0):
561    def set_parameter(self, name, value, scale=1.0):
562        """
563        Set the parameter value.
564        The initial amount by which the parameter
565        will be modified during the search for the minimum. 
566        """
567        self.minimizer.SetParameterValue(name, value)
568        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):
570    def set_parameters(self, parameters):
571        """
572        Set the parameters names and values from a dictionary.
573        """
574        for name, value in parameters.items():
575            if len(value) == 2:
576                self.set_parameter(name, value[0], value[1])
577            else:
578                self.set_parameter(name, value)

Set the parameters names and values from a dictionary.

def minimize(self):
580    def minimize(self):
581        """
582        Minimize the input function.
583
584        Returns:
585            dict : 
586                the minimization results
587            init_parameters : (dict)
588                the initial parameters
589            parameters : (dict)
590                the final parameters
591            min_value : (float)
592                the minimum value
593            iterations : (int)
594                the number of iterations
595            max_iterations : (int)
596                the maximum number of iterations
597            tolerance : (float)
598                the tolerance for convergence
599            convergence_flag : (int)
600                zero if the tolerance stopping
601                criterion has been met.
602            parameters_path : (np.array)
603                the path of the minimization
604                algorithm in parameter space
605            function_path : (np.array)
606                the path of the minimization
607                algorithm in function space
608            hessian : (np.array)
609                the Hessian matrix of the
610                function at the minimum
611            parameter_errors : (np.array)
612                the errors on the parameters
613        """
614        n = self.minimizer.GetNumberOfParameters()
615        out = [(
616            self.minimizer.GetParameterName(i),
617            (self.minimizer.GetParameterValue(i),
618             self.minimizer.GetParameterScale(i))
619        ) for i in range(n)]
620        self.results["init_parameters"] = dict(out)
621
622        self.minimizer.SetTolerance(self.tolerance)
623        self.minimizer.SetContractionRatio(self.contraction_ratio)
624        self.minimizer.SetExpansionRatio(self.expansion_ratio)
625        self.minimizer.SetMaxIterations(self.max_iterations)
626
627        self.minimizer.Minimize()
628        self.results["convergence_flag"] = not bool(self.minimizer.Iterate())
629
630        out = [(
631            self.minimizer.GetParameterName(i),
632            self.minimizer.GetParameterValue(i),
633        ) for i in range(n)]
634
635        self.results["parameters"] = dict(out)
636        self.results["min_value"] = self.minimizer.GetFunctionValue()
637        self.results["iterations"] = self.minimizer.GetIterations()
638        self.results["max_iterations"] = self.minimizer.GetMaxIterations()
639        self.results["tolerance"] = self.minimizer.GetTolerance()
640        self.results["expansion_ratio"] = self.expansion_ratio
641        self.results["contraction_ratio"] = self.contraction_ratio
642        self.results["parameters_path"] = np.array(self.parameters_path)
643        self.results["function_path"] = np.array(self.function_path)
644        self.results["hessian"] = np.zeros((n,n))
645        self.results["parameter_errors"] = np.zeros(n)
646        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):
648    def compute_hessian(self, epsilon=0):
649        """
650        Compute the Hessian matrix of `function` at the
651        minimum numerically.
652
653        Arguments:
654            epsilon : (float)
655                Step size used for numerical approximation.
656
657        Returns:
658            array: Hessian matrix of `function` at minimum.
659        """
660        if not epsilon:
661            epsilon = self.tolerance * 10
662        n = self.minimizer.GetNumberOfParameters()
663        x0 = [self.minimizer.GetParameterValue(i) for i in range(n)]
664        hessian = np.zeros((n, n))
665        for i in vedo.progressbar(n, title="Computing Hessian", delay=2):
666            for j in range(n):
667                xijp = np.copy(x0)
668                xijp[i] += epsilon
669                xijp[j] += epsilon
670                xijm = np.copy(x0)
671                xijm[i] += epsilon
672                xijm[j] -= epsilon
673                xjip = np.copy(x0)
674                xjip[i] -= epsilon
675                xjip[j] += epsilon
676                xjim = np.copy(x0)
677                xjim[i] -= epsilon
678                xjim[j] -= epsilon
679                # Second derivative approximation
680                fijp = self.function(xijp)
681                fijm = self.function(xijm)
682                fjip = self.function(xjip)
683                fjim = self.function(xjim)
684                hessian[i, j] = (fijp - fijm - fjip + fjim) / (2 * epsilon**2)
685        self.results["hessian"] = hessian
686        try:
687            ihess = np.linalg.inv(hessian)
688            self.results["parameter_errors"] = np.sqrt(np.diag(ihess))
689        except:
690            vedo.logger.warning("Cannot compute hessian for parameter errors")
691            self.results["parameter_errors"] = np.zeros(n)
692        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 function at minimum.

def geometry(obj, extent=None):
868def geometry(obj, extent=None):
869    """
870    Apply the `vtkGeometryFilter` to the input object.
871    This is a general-purpose filter to extract geometry (and associated data)
872    from any type of dataset.
873    This filter also may be used to convert any type of data to polygonal type.
874    The conversion process may be less than satisfactory for some 3D datasets.
875    For example, this filter will extract the outer surface of a volume
876    or structured grid dataset.
877
878    Returns a `vedo.Mesh` object.
879
880    Set `extent` as the `[xmin,xmax, ymin,ymax, zmin,zmax]` bounding box to clip data.
881    """
882    gf = vtki.new("GeometryFilter")
883    gf.SetInputData(obj)
884    if extent is not None:
885        gf.SetExtent(extent)
886    gf.Update()
887    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):
1090def is_sequence(arg):
1091    """Check if the input is iterable."""
1092    if hasattr(arg, "strip"):
1093        return False
1094    if hasattr(arg, "__getslice__"):
1095        return True
1096    if hasattr(arg, "__iter__"):
1097        return True
1098    return False

Check if the input is iterable.

def lin_interpolate(x, rangeX, rangeY):
1435def lin_interpolate(x, rangeX, rangeY):
1436    """
1437    Interpolate linearly the variable `x` in `rangeX` onto the new `rangeY`.
1438    If `x` is a 3D vector the linear weight is the distance to the two 3D `rangeX` vectors.
1439
1440    E.g. if `x` runs in `rangeX=[x0,x1]` and I want it to run in `rangeY=[y0,y1]` then
1441
1442    `y = lin_interpolate(x, rangeX, rangeY)` will interpolate `x` onto `rangeY`.
1443
1444    Examples:
1445        - [lin_interpolate.py](https://github.com/marcomusy/vedo/tree/master/examples/basic/lin_interpolate.py)
1446
1447            ![](https://vedo.embl.es/images/basic/linInterpolate.png)
1448    """
1449    if is_sequence(x):
1450        x = np.asarray(x)
1451        x0, x1 = np.asarray(rangeX)
1452        y0, y1 = np.asarray(rangeY)
1453        dx = x1 - x0
1454        dxn = np.linalg.norm(dx)
1455        if not dxn:
1456            return y0
1457        s = np.linalg.norm(x - x0) / dxn
1458        t = np.linalg.norm(x - x1) / dxn
1459        st = s + t
1460        out = y0 * (t / st) + y1 * (s / st)
1461
1462    else:  # faster
1463
1464        x0 = rangeX[0]
1465        dx = rangeX[1] - x0
1466        if not dx:
1467            return rangeY[0]
1468        s = (x - x0) / dx
1469        out = rangeY[0] * (1 - s) + rangeY[1] * s
1470    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'>):
1547def vector(x, y=None, z=0.0, dtype=np.float64):
1548    """
1549    Return a 3D numpy array representing a vector.
1550
1551    If `y` is `None`, assume input is already in the form `[x,y,z]`.
1552    """
1553    if y is None:  # assume x is already [x,y,z]
1554        return np.asarray(x, dtype=dtype)
1555    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):
1566def mag(v):
1567    """Get the magnitude of a vector or array of vectors."""
1568    v = np.asarray(v)
1569    if v.ndim == 1:
1570        return np.linalg.norm(v)
1571    return np.linalg.norm(v, axis=1)

Get the magnitude of a vector or array of vectors.

def mag2(v):
1574def mag2(v):
1575    """Get the squared magnitude of a vector or array of vectors."""
1576    v = np.asarray(v)
1577    if v.ndim == 1:
1578        return np.square(v).sum()
1579    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'>):
1558def versor(x, y=None, z=0.0, dtype=np.float64):
1559    """Return the unit vector. Input can be a list of vectors."""
1560    v = vector(x, y, z, dtype)
1561    if isinstance(v[0], np.ndarray):
1562        return np.divide(v, mag(v)[:, None])
1563    return v / mag(v)

Return the unit vector. Input can be a list of vectors.

def precision(x, p, vrange=None, delimiter='e'):
1648def precision(x, p, vrange=None, delimiter="e"):
1649    """
1650    Returns a string representation of `x` formatted to precision `p`.
1651
1652    Set `vrange` to the range in which x exists (to snap x to '0' if below precision).
1653    """
1654    # Based on the webkit javascript implementation
1655    # `from here <https://code.google.com/p/webkit-mirror/source/browse/JavaScriptCore/kjs/number_object.cpp>`_,
1656    # and implemented by `randlet <https://github.com/randlet/to-precision>`_.
1657    # Modified for vedo by M.Musy 2020
1658
1659    if isinstance(x, str):  # do nothing
1660        return x
1661
1662    if is_sequence(x):
1663        out = "("
1664        nn = len(x) - 1
1665        for i, ix in enumerate(x):
1666
1667            try:
1668                if np.isnan(ix):
1669                    return "NaN"
1670            except:
1671                # cannot handle list of list
1672                continue
1673
1674            out += precision(ix, p)
1675            if i < nn:
1676                out += ", "
1677        return out + ")"  ############ <--
1678
1679    try:
1680        if np.isnan(x):
1681            return "NaN"
1682    except TypeError:
1683        return "NaN"
1684
1685    x = float(x)
1686
1687    if x == 0.0 or (vrange is not None and abs(x) < vrange / pow(10, p)):
1688        return "0"
1689
1690    out = []
1691    if x < 0:
1692        out.append("-")
1693        x = -x
1694
1695    e = int(np.log10(x))
1696    # tens = np.power(10, e - p + 1)
1697    tens = 10 ** (e - p + 1)
1698    n = np.floor(x / tens)
1699
1700    # if n < np.power(10, p - 1):
1701    if n < 10 ** (p - 1):
1702        e = e - 1
1703        # tens = np.power(10, e - p + 1)
1704        tens = 10 ** (e - p + 1)
1705        n = np.floor(x / tens)
1706
1707    if abs((n + 1.0) * tens - x) <= abs(n * tens - x):
1708        n = n + 1
1709
1710    # if n >= np.power(10, p):
1711    if n >= 10 ** p:
1712        n = n / 10.0
1713        e = e + 1
1714
1715    m = "%.*g" % (p, n)
1716    if e < -2 or e >= p:
1717        out.append(m[0])
1718        if p > 1:
1719            out.append(".")
1720            out.extend(m[1:p])
1721        out.append(delimiter)
1722        if e > 0:
1723            out.append("+")
1724        out.append(str(e))
1725    elif e == (p - 1):
1726        out.append(m)
1727    elif e >= 0:
1728        out.append(m[: e + 1])
1729        if e + 1 < len(m):
1730            out.append(".")
1731            out.extend(m[e + 1 :])
1732    else:
1733        out.append("0.")
1734        out.extend(["0"] * -(e + 1))
1735        out.append(m)
1736    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):
1601def round_to_digit(x, p):
1602    """Round a real number to the specified number of significant digits."""
1603    if not x:
1604        return 0
1605    r = np.round(x, p - int(np.floor(np.log10(abs(x)))) - 1)
1606    if int(r) == r:
1607        return int(r)
1608    return r

Round a real number to the specified number of significant digits.

def point_in_triangle(p, p1, p2, p3):
1173def point_in_triangle(p, p1, p2, p3):
1174    """
1175    Return True if a point is inside (or above/below)
1176    a triangle defined by 3 points in space.
1177    """
1178    p1 = np.array(p1)
1179    u = p2 - p1
1180    v = p3 - p1
1181    n = np.cross(u, v)
1182    w = p - p1
1183    ln = np.dot(n, n)
1184    if not ln:
1185        return None  # degenerate triangle
1186    gamma = (np.dot(np.cross(u, w), n)) / ln
1187    if 0 < gamma < 1:
1188        beta