diff --git a/examples/basic/clustering.py b/examples/basic/clustering.py index f42d61da..d7bce128 100644 --- a/examples/basic/clustering.py +++ b/examples/basic/clustering.py @@ -22,4 +22,4 @@ # find back their identity through clustering cl = cluster(pts, radius=0.1) # returns a vtkAssembly -show(cl, Text(__doc__), axes=1, verbose=0) +show(cl, Text(__doc__), axes=1, viewup='z') diff --git a/examples/basic/mesh_sharemap.py b/examples/basic/mesh_sharemap.py index 93608bfa..11505584 100644 --- a/examples/basic/mesh_sharemap.py +++ b/examples/basic/mesh_sharemap.py @@ -2,7 +2,6 @@ How to share the same color map across different meshes. """ -print(__doc__) from vtkplotter import load, Text, show, datadir @@ -11,13 +10,11 @@ scals = man1.coordinates()[:, 2] * 5 + 27 # pick z coordinates [18->34] man1.pointColors(scals, cmap="jet", vmin=18, vmax=44) -man1.show(N=2, at=0, axes=0, elevation=-80) ##################################### -man2 = load(datadir+"man.vtk").addScalarBar() +man2 = load(datadir+"man.vtk") scals = man2.coordinates()[:, 2] * 5 + 37 # pick z coordinates [28->44] man2.pointColors(scals, cmap="jet", vmin=18, vmax=44) -# man2.show(at=1, interactive=0) -show(man2, Text(__doc__), at=1, interactive=1) +show([[man1, Text(__doc__)], man2], N=2, elevation=-40) diff --git a/examples/notebooks/sphere.ipynb b/examples/notebooks/sphere.ipynb index ffae4e80..a3c49bda 100644 --- a/examples/notebooks/sphere.ipynb +++ b/examples/notebooks/sphere.ipynb @@ -2,13 +2,13 @@ "cells": [ { "cell_type": "code", - "execution_count": 39, + "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { - "model_id": "1c932fd2f6124207afb7e72f40dff52a", + "model_id": "e80883391eaf46d0b90e212fe2a35414", "version_major": 2, "version_minor": 0 }, diff --git a/examples/other/dolfin/curl2d.py b/examples/other/dolfin/curl2d.py new file mode 100644 index 00000000..5505ce81 --- /dev/null +++ b/examples/other/dolfin/curl2d.py @@ -0,0 +1,68 @@ +#Magnetostatic equation that solves for magnetic vector potential, +#where the reluctivity matrix can be defined as inverse of permeability matrix.. +#https://fenicsproject.discourse.group/t/anisotropic-material-definition-and-results-issue/1051 +from dolfin import * +from mshr import * +from scipy import constants +from vtkplotter.dolfin import plot + +domain = Rectangle(Point(-10, -10), Point(10, 10)) +mesh = generate_mesh(domain, 64) + +# function space +V = FunctionSpace(mesh, "P", 1) + +# boundary conditions +walls = "on_boundary && (near(abs(x[0]), 10.0) || near(abs(x[1]), 10.0))" +bc = DirichletBC(V, Constant(0.0), walls) + +tol = 1e-6 + +# Wire +class Omega_0(SubDomain): + def inside(self, x, on_boundary): + return x[0] ** 2 + x[1] ** 2 <= 4 - tol + +# Space +class Omega_1(SubDomain): + def inside(self, x, on_boundary): + return x[0] ** 2 + x[1] ** 2 > 4 + tol + +def curl2D(v): + return as_vector((v.dx(1), -v.dx(0))) + + +materials = MeshFunction("size_t", mesh, mesh.topology().dim()) + +subdomain_0 = Omega_0() +subdomain_1 = Omega_1() +subdomain_0.mark(materials, 0) +subdomain_1.mark(materials, 1) + +dx = Measure("dx", domain=mesh, subdomain_data=materials) + +A_z = Function(V) # magnetic vector potential +v = TestFunction(V) + +J = 5.0e6 +# anisotropic material parameters, reluctivity = 1/constants.mu_0 +reluctivity = as_matrix( + ((1 / (constants.mu_0 * 1000), 0), + (0, 1 / (constants.mu_0 * 1))) +) + +F = inner(reluctivity * curl2D(A_z), curl2D(v)) * dx - J * v * dx(0) +solve(F == 0, A_z, bc) + +W = VectorFunctionSpace(mesh, "P", 1) + +Bx = A_z.dx(1) +By = -A_z.dx(0) +B = project(as_vector((Bx, By)), W) + +plot(B, mode='mesh and arrows', + style=2, + scale=0.01, + lw=0, + warpZfactor=-0.01, + ) diff --git a/examples/other/dolfin/magnetostatics.py b/examples/other/dolfin/magnetostatics.py index 7c587cdc..6baf1894 100644 --- a/examples/other/dolfin/magnetostatics.py +++ b/examples/other/dolfin/magnetostatics.py @@ -88,5 +88,7 @@ def eval_cell(self, values, x, cell): plot(A_z, at=0, N=2, # draw on the first of 2 renderers lw=0, # linewidth of mesh bg='white', - isolines={'n':20, 'lw':1.5, 'c':'black'}) -plot(B, at=1, text=__doc__) # draw on the second renderer + isolines={'n':20, 'lw':1.5, 'c':'black'}, + scalarbar=False, + ) +plot(B, at=1, scalarbar=False, text=__doc__) # draw on the second renderer diff --git a/examples/other/dolfin/run_all.sh b/examples/other/dolfin/run_all.sh index 25aab35f..d4b6aa1c 100755 --- a/examples/other/dolfin/run_all.sh +++ b/examples/other/dolfin/run_all.sh @@ -38,6 +38,8 @@ python elasticbeam.py echo Running magnetostatics.py python magnetostatics.py +echo Running curl2d.py +python curl2d.py ###################################### echo Running ex01_show-mesh.py diff --git a/examples/basic/qt_window.py b/examples/other/qt_window.py similarity index 69% rename from examples/basic/qt_window.py rename to examples/other/qt_window.py index fd57c09d..6daf7ea1 100644 --- a/examples/basic/qt_window.py +++ b/examples/other/qt_window.py @@ -1,12 +1,14 @@ +""" +A sort of minimal example of how to embed a rendering window +into a qt application. +""" +print(__doc__) import sys from PyQt5 import Qt from vtk.qt.QVTKRenderWindowInteractor import QVTKRenderWindowInteractor from vtkplotter import Plotter, Cone -# settings.usingQt = True <-- not needed anymore, automatically triggered by passing a qtWidget to Plotter - - class MainWindow(Qt.QMainWindow): def __init__(self, parent=None): @@ -16,19 +18,19 @@ def __init__(self, parent=None): self.vtkWidget = QVTKRenderWindowInteractor(self.frame) self.vl.addWidget(self.vtkWidget) - vp = Plotter(qtWidget=self.vtkWidget, axes=4, verbose=False) + vp = Plotter(qtWidget=self.vtkWidget, axes=4, bg='white') vp += Cone() - vp.show() # to create renderer and add the actors + vp.show() # create renderer and add the actors # set-up the rest of the Qt window self.frame.setLayout(self.vl) self.setCentralWidget(self.frame) - self.show() # <--- show the Qt Window + self.show() # <--- show the Qt Window if __name__ == "__main__": app = Qt.QApplication(sys.argv) window = MainWindow() - app.exec_() \ No newline at end of file + app.exec_() diff --git a/examples/other/qt_embed.py b/examples/other/qt_window_split.py similarity index 95% rename from examples/other/qt_embed.py rename to examples/other/qt_window_split.py index bf0e4215..ad804acb 100644 --- a/examples/other/qt_embed.py +++ b/examples/other/qt_window_split.py @@ -10,8 +10,6 @@ from vtkplotter import Plotter, Cube, Torus, Cone, settings -settings.usingQt = True - class MainWindow(Qt.QMainWindow): def __init__(self, parent=None): @@ -22,7 +20,7 @@ def __init__(self, parent=None): self.vtkWidget = QVTKRenderWindowInteractor(self.frame) self.vl.addWidget(self.vtkWidget) - vp = Plotter(offscreen=1, interactive=0, axes=2, N=2) + vp = Plotter(qtWidget=self.vtkWidget, axes=2, N=2) cn = Cone() cc = Cube().pos(1, 1, 1).color("pink") diff --git a/examples/other/self_org_maps2d.py b/examples/other/self_org_maps2d.py index 198db320..fa66a5ef 100644 --- a/examples/other/self_org_maps2d.py +++ b/examples/other/self_org_maps2d.py @@ -47,7 +47,8 @@ def learn(self, n_epoch=10000, sigma=(0.25,0.01), lrate=(0.5,0.01)): # Draw network if i>500 and not i%20 or i==n_epoch-1: x, y, z = [self.codebook[:,i].reshape(n,n) for i in range(3)] - grd = Grid(resx=n-1, resy=n-1).wire(False).lw(0.5).bc('lightblue') + grd = Grid(resx=n-1, resy=n-1) + grd.wireframe(False).lw(0.5).bc('lightblue').flat() for i in range(n): for j in range(n): grd.setPoint(i*n+j, (x[i,j], y[i,j], z[i,j])) diff --git a/examples/run_all.sh b/examples/run_all.sh index b62f228d..3ed58787 100755 --- a/examples/run_all.sh +++ b/examples/run_all.sh @@ -352,8 +352,11 @@ python other/icon.py echo Running other/inset.py python other/inset.py -echo Running other/qt_embed.py # needs qt5 -python other/qt_embed.py +echo Running other/qt_window.py # needs qt5 +python other/qt_window.py + +echo Running other/qt_window_split.py # needs qt5 +python other/qt_window_split.py echo Running other/self_org_maps2d.py python other/self_org_maps2d.py diff --git a/examples/simulations/README.md b/examples/simulations/README.md index be2b5f5e..6547b981 100644 --- a/examples/simulations/README.md +++ b/examples/simulations/README.md @@ -25,8 +25,6 @@ python example.py # on mac OSX try 'pythonw' instead | | | | [![mpend](https://user-images.githubusercontent.com/32848391/50738892-db380300-11d8-11e9-807c-fb320c7b7917.gif)](https://github.com/marcomusy/vtkplotter/blob/master/examples/simulations/multiple_pendulum.py)
`multiple_pendulum.py` | Simulation of an elastic multiple pendulums with viscous friction. | | | | -| [![pendulum](https://user-images.githubusercontent.com/32848391/55420020-51e56200-5576-11e9-8513-4a5d93913b17.png)](https://github.com/marcomusy/vtkplotter/blob/master/examples/simulations/pendulum.py)
`pendulum.py` | Visualize the phase space of a simple pendulum (from [3Blue1Brown](https://www.youtube.com/watch?v=p_di4Zn4wz4)). | -| | | | [![hanoi](https://user-images.githubusercontent.com/32848391/56989284-58c1bd80-6b92-11e9-8f82-1ce95813f846.gif)](https://github.com/marcomusy/vtkplotter/blob/master/examples/simulations/hanoi3d.py)
`hanoi3d.py` | Solve the Tower of Hanoi puzzle (contributed by [G. Jacquenot](https://github.com/Gjacquenot)). | | | | | [![airplanes](https://user-images.githubusercontent.com/32848391/57341963-b8910900-713c-11e9-898a-84b6d3712bce.gif)](https://github.com/marcomusy/vtkplotter/blob/master/examples/simulations/airplanes.py)
`airplanes.py` | Two acrobatic planes casting shadows and leaving a trail. | diff --git a/examples/simulations/multiple_pendulum.py b/examples/simulations/multiple_pendulum.py index 04b77fe7..e2590557 100644 --- a/examples/simulations/multiple_pendulum.py +++ b/examples/simulations/multiple_pendulum.py @@ -1,138 +1,110 @@ -""" -A model of an ideal gas with hard-sphere collisions. -""" -## Based on gas.py by Bruce Sherwood for a cube as a container -## Sligthly modified by Andrey Antonov for a torus. -## Adapted by M. Musy for vtkplotter -## relevant points in the code are marked with '### <--' from __future__ import division, print_function -from random import random -from vtkplotter import Plotter, ProgressBar, mag, versor, Text -from vtkplotter import Torus, Sphere +from vtkplotter import Plotter, printc, mag, versor, vector +from vtkplotter import Cylinder, Spring, Box, Sphere import numpy as np -############################################################# -Natoms = 400 # change this to have more or fewer atoms -Nsteps = 350 # nr of steps in the simulation -Matom = 4e-3 / 6e23 # helium mass -Ratom = 0.025 # wildly exaggerated size of helium atom -RingThickness = 0.3 # thickness of the toroid -RingRadius = 1 -k = 1.4e-23 # Boltzmann constant -T = 300 # room temperature -dt = 1.5e-5 -############################################################# - - -def reflection(p, pos): - n = versor(pos) - return np.dot(np.identity(3) - 2 * n * n[:, np.newaxis], p) - - -vp = Plotter(title="gas in toroid", interactive=0, axes=0, bg="w") - -vp += Text(__doc__) -vp += Torus(c="g", r=RingRadius, thickness=RingThickness, alpha=0.1).wire(1) ### <-- - -Atoms = [] -poslist = [] -plist, mlist, rlist = [], [], [] -mass = Matom * Ratom ** 3 / Ratom ** 3 -pavg = np.sqrt(2.0 * mass * 1.5 * k * T) # average kinetic energy p**2/(2mass) = (3/2)kT - -for i in range(Natoms): - alpha = 2 * np.pi * random() - x = RingRadius * np.cos(alpha) * 0.9 - y = RingRadius * np.sin(alpha) * 0.9 - z = 0 - atm = Sphere(pos=(x, y, z), r=Ratom, c=i) - vp += atm - Atoms = Atoms + [atm] ### <-- - theta = np.pi * random() - phi = 2 * np.pi * random() - px = pavg * np.sin(theta) * np.cos(phi) - py = pavg * np.sin(theta) * np.sin(phi) - pz = pavg * np.cos(theta) - poslist.append((x, y, z)) - plist.append((px, py, pz)) - mlist.append(mass) - rlist.append(Ratom) - -pos = np.array(poslist) -poscircle = pos -p = np.array(plist) -m = np.array(mlist) -m.shape = (Natoms, 1) -radius = np.array(rlist) -r = pos - pos[:, np.newaxis] # all pairs of atom-to-atom vectors - -ds = (p / m) * (dt / 2.0) -if "False" not in np.less_equal(mag(ds), radius).tolist(): - pos = pos + (p / mass) * (dt / 2.0) # initial half-step - -pb = ProgressBar(0, Nsteps, c=1) -for i in pb.range(): - - # Update all positions - ds = mag((p / m) * (dt / 2.0)) - if "False" not in np.less_equal(ds, radius).tolist(): - pos = pos + (p / m) * dt - - r = pos - pos[:, np.newaxis] # all pairs of atom-to-atom vectors - rmag = np.sqrt(np.sum(np.square(r), -1)) # atom-to-atom scalar distances - hit = np.less_equal(rmag, radius + radius[:, None]) - np.identity(Natoms) - hitlist = np.sort(np.nonzero(hit.flat)[0]).tolist() # i,j encoded as i*Natoms+j - - # If any collisions took place: - for ij in hitlist: - i, j = divmod(ij, Natoms) # decode atom pair - hitlist.remove(j * Natoms + i) # remove symmetric j,i pair from list - ptot = p[i] + p[j] - mi = m[i, 0] - mj = m[j, 0] - vi = p[i] / mi - vj = p[j] / mj - ri = Ratom - rj = Ratom - a = mag(vj - vi) ** 2 - if a == 0: - continue # exactly same velocities - b = 2 * np.dot(pos[i] - pos[j], vj - vi) - c = mag(pos[i] - pos[j]) ** 2 - (ri + rj) ** 2 - d = b ** 2 - 4.0 * a * c - if d < 0: - continue # something wrong; ignore this rare case - deltat = (-b + np.sqrt(d)) / (2.0 * a) # t-deltat is when they made contact - pos[i] = pos[i] - (p[i] / mi) * deltat # back up to contact configuration - pos[j] = pos[j] - (p[j] / mj) * deltat - mtot = mi + mj - pcmi = p[i] - ptot * mi / mtot # transform momenta to cm frame - pcmj = p[j] - ptot * mj / mtot - rrel = versor(pos[j] - pos[i]) - pcmi = pcmi - 2 * np.dot(pcmi, rrel) * rrel # bounce in cm frame - pcmj = pcmj - 2 * np.dot(pcmj, rrel) * rrel - p[i] = pcmi + ptot * mi / mtot # transform momenta back to lab frame - p[j] = pcmj + ptot * mj / mtot - pos[i] = pos[i] + (p[i] / mi) * deltat # move forward deltat in time - pos[j] = pos[j] + (p[j] / mj) * deltat - - # Bounce off the boundary of the torus - for j in range(Natoms): - poscircle[j] = versor(pos[j]) * RingRadius * [1, 1, 0] - outside = np.greater_equal(mag(poscircle - pos), RingThickness - 2 * Ratom) - - for k in range(len(outside)): - if outside[k] == 1 and np.dot(p[k], pos[k] - poscircle[k]) > 0: - p[k] = reflection(p[k], pos[k] - poscircle[k]) - - # then update positions of display objects - for i in range(Natoms): - Atoms[i].pos(pos[i]) ### <-- - outside = np.greater_equal(mag(pos), RingRadius + RingThickness) - - vp.show() ### <-- - vp.camera.Azimuth(0.5) - vp.camera.Elevation(0.1) - pb.print() - -vp.show(interactive=1) +############## Constants +N = 5 # number of bobs +R = 0.3 # radius of bob (separation between bobs=1) +Ks = 50 # k of springs (masses=1) +g = 9.81 # gravity acceleration +gamma = 0.1 # some friction +Dt = 0.03 # time step + +# Create the initial positions and velocitites (0,0) of the bobs +bob_x = [0] +bob_y = [0] +x_dot = [0] * (N + 1) # velocities +y_dot = [0] * (N + 1) + +for k in range(1, N + 1): + alpha = np.pi / 5 * k / 10 + bob_x.append(bob_x[k - 1] + np.cos(alpha) + np.random.normal(0, 0.1)) + bob_y.append(bob_y[k - 1] + np.sin(alpha) + np.random.normal(0, 0.1)) + +# Create the bobs +vp = Plotter(title="Multiple Pendulum", axes=0, interactive=0) +vp += Box(pos=(0, -5, 0), length=12, width=12, height=0.7, c="k").wire(1) +bob = [vp.add(Sphere(pos=(bob_x[0], bob_y[0], 0), r=R / 2, c="gray"))] +for k in range(1, N + 1): + bob.append(vp.add(Cylinder(pos=(bob_x[k], bob_y[k], 0), r=R, height=0.3, c=k))) + +# Create the springs out of N links +link = [0] * N +for k in range(N): + p0 = bob[k].pos() + p1 = bob[k + 1].pos() + link[k] = vp.add(Spring(p0, p1, thickness=0.015, r=R / 3, c="gray")) + +# Create some auxiliary variables +x_dot_m = [0] * (N + 1) +y_dot_m = [0] * (N + 1) +dij = [0] * (N + 1) # array with distances to previous bob +dij_m = [0] * (N + 1) +for k in range(1, N + 1): + dij[k] = mag([bob_x[k] - bob_x[k - 1], bob_y[k] - bob_y[k - 1]]) + +fctr = lambda x: (x - 1) / x +Dt *= np.sqrt(1 / g) +Dt2 = Dt / 2 # Midpoint time step +DiaSq = (2 * R) ** 2 # Diameter of bob squared + +printc("Press Esc to exit.", c="red") + +while True: + bob_x_m = list(map((lambda x, dx: x + Dt2 * dx), bob_x, x_dot)) # midpoint variables + bob_y_m = list(map((lambda y, dy: y + Dt2 * dy), bob_y, y_dot)) + + for k in range(1, N + 1): + factor = fctr(dij[k]) + x_dot_m[k] = x_dot[k] - Dt2 * (Ks * (bob_x[k] - bob_x[k - 1]) * factor + gamma * x_dot[k]) + y_dot_m[k] = y_dot[k] - Dt2 * ( + Ks * (bob_y[k] - bob_y[k - 1]) * factor + gamma * y_dot[k] + g + ) + + for k in range(1, N): + factor = fctr(dij[k + 1]) + x_dot_m[k] -= Dt2 * Ks * (bob_x[k] - bob_x[k + 1]) * factor + y_dot_m[k] -= Dt2 * Ks * (bob_y[k] - bob_y[k + 1]) * factor + + # Compute the full step variables + bob_x = list(map((lambda x, dx: x + Dt * dx), bob_x, x_dot_m)) + bob_y = list(map((lambda y, dy: y + Dt * dy), bob_y, y_dot_m)) + + for k in range(1, N + 1): + dij[k] = mag([bob_x[k] - bob_x[k - 1], bob_y[k] - bob_y[k - 1]]) + dij_m[k] = mag([bob_x_m[k] - bob_x_m[k - 1], bob_y_m[k] - bob_y_m[k - 1]]) + factor = fctr(dij_m[k]) + x_dot[k] -= Dt * (Ks * (bob_x_m[k] - bob_x_m[k - 1]) * factor + gamma * x_dot_m[k]) + y_dot[k] -= Dt * (Ks * (bob_y_m[k] - bob_y_m[k - 1]) * factor + gamma * y_dot_m[k] + g) + + for k in range(1, N): + factor = fctr(dij_m[k + 1]) + x_dot[k] -= Dt * Ks * (bob_x_m[k] - bob_x_m[k + 1]) * factor + y_dot[k] -= Dt * Ks * (bob_y_m[k] - bob_y_m[k + 1]) * factor + + # Check to see if they are colliding + for i in range(1, N): + for j in range(i + 1, N + 1): + dist2 = (bob_x[i] - bob_x[j]) ** 2 + (bob_y[i] - bob_y[j]) ** 2 + if dist2 < DiaSq: # are colliding + Ddist = np.sqrt(dist2) - 2 * R + tau = versor([bob_x[j] - bob_x[i], bob_y[j] - bob_y[i], 0]) + DR = Ddist / 2 * tau + bob_x[i] += DR[0] # DR.x + bob_y[i] += DR[1] # DR.y + bob_x[j] -= DR[0] # DR.x + bob_y[j] -= DR[1] # DR.y + Vji = vector(x_dot[j] - x_dot[i], y_dot[j] - y_dot[i]) + DV = np.dot(Vji, tau) * tau + x_dot[i] += DV[0] # DV.x + y_dot[i] += DV[1] # DV.y + x_dot[j] -= DV[0] # DV.x + y_dot[j] -= DV[1] # DV.y + + # Update the loations of the bobs and the stretching of the springs + for k in range(1, N + 1): + bob[k].pos([bob_x[k], bob_y[k], 0]) + link[k - 1].stretch(bob[k - 1].pos(), bob[k].pos()) + + vp.show() diff --git a/examples/volumetric/create_imagedata.ipynb b/examples/volumetric/create_imagedata.ipynb new file mode 100644 index 00000000..4d8b13a7 --- /dev/null +++ b/examples/volumetric/create_imagedata.ipynb @@ -0,0 +1,72 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "04b1c7d9346547a196e5868a1707560e", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "Plot(antialias=3, axes=['x', 'y', 'z'], background_color=16777215, camera=[4.5, 4.5, 4.5, 0.0, 0.0, 0.0, 1.0, …" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Use VTK commands to create a vtkImageData object and visualize it\"\"\"\n", + "\n", + "import vtk\n", + "import numpy as np\n", + "from vtkplotter import show\n", + "\n", + "n1, n2, n3 = 5,10,15\n", + "raw = np.random.rand(n1*n2*n3).reshape(n1,n2,n3)\n", + "\n", + "img = vtk.vtkImageData()\n", + "img.SetDimensions(raw.shape)\n", + "\n", + "raw = raw.transpose(2,1,0).flatten()\n", + "arr = vtk.util.numpy_support.numpy_to_vtk(raw, array_type=vtk.VTK_FLOAT)\n", + "img.GetPointData().SetScalars(arr)\n", + "\n", + "show(img)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.6.8" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/examples/volumetric/streamribbons.ipynb b/examples/volumetric/streamribbons.ipynb index 9298ff52..3a59de99 100644 --- a/examples/volumetric/streamribbons.ipynb +++ b/examples/volumetric/streamribbons.ipynb @@ -8,7 +8,7 @@ { "data": { "application/vnd.jupyter.widget-view+json": { - "model_id": "bee02179603045679c2b9ec751fded14", + "model_id": "efb47a2cfe1642d89175808353130121", "version_major": 2, "version_minor": 0 }, diff --git a/examples/volumetric/tensors.py b/examples/volumetric/tensors.py index f540c64e..462abd73 100644 --- a/examples/volumetric/tensors.py +++ b/examples/volumetric/tensors.py @@ -1,3 +1,4 @@ +"""Visualize stress tensors as ellipsoids.""" import vtk from vtkplotter import * @@ -15,8 +16,9 @@ zsl = vol.zSlice(3) # Generate tensor ellipsoids -#tens = Tensors(vol, source='ellipse', scale=10) -tens = Tensors(zsl, source='ellipse', scale=20) +tens1 = Tensors(vol, source='ellipse', scale=10) +tens2 = Tensors(zsl, source='ellipse', scale=20) +t = Text(__doc__, c='k') -#show([vol, [tens, zsl]], N=2, axes=1, viewup='z') -show(vol, tens, zsl, axes=1, bg='w', viewup='z') +show([[vol, t], tens1], N=2, axes=9, bg='w', viewup='z') +show(vol, tens2, zsl, axes=9, bg='w', viewup='z', newPlotter=True) diff --git a/setup.py b/setup.py index 548f85e2..3a51d0d6 100644 --- a/setup.py +++ b/setup.py @@ -50,11 +50,11 @@ # pip install . # ( sudo -H pip install . ) -# cd examples -# ./run_all.sh +# cd examples && ./run_all.sh # check vtkconvert: # vtkconvert vtkplotter/data/290.vtk -to ply +# rm vtkplotter/data/290.ply # check on python2 the same stuff is ok diff --git a/vtkplotter/actors.py b/vtkplotter/actors.py index 189ddafa..d1a5356f 100644 --- a/vtkplotter/actors.py +++ b/vtkplotter/actors.py @@ -873,7 +873,6 @@ def __init__( self._bfprop = None # backface property holder prp = self.GetProperty() - #prp.SetInterpolationToFlat() if settings.renderPointsAsSpheres: if hasattr(prp, 'RenderPointsAsSpheresOn'): @@ -884,8 +883,10 @@ def __init__( prp.RenderLinesAsTubesOn() if c is None: - self.mapper.ScalarVisibilityOn() - prp.SetColor(colors.getColor("gold")) + if self.poly and self.poly.GetPointData().GetScalars(): + self.mapper.ScalarVisibilityOn() + else: + prp.SetColor(colors.getColor("gold")) else: self.mapper.ScalarVisibilityOff() c = colors.getColor(c) @@ -1111,7 +1112,7 @@ def getConnectivity(self): def deletePoints(self, indices): """Delete a list of vertices identified by their index. - + |deleteMeshPoints| |deleteMeshPoints.py|_ """ cellIds = vtk.vtkIdList() @@ -2998,7 +2999,7 @@ def __init__(self, img, for iz in range(nz): vtkimg.SetScalarComponentFromFloat(ix, iy, iz, 0, img[ix, iy, iz]) img = vtkimg - + if origin is not None: img.SetOrigin(origin) if spacing is not None: @@ -3011,7 +3012,7 @@ def __init__(self, img, self.mode(mode) self.color(c).alpha(alpha) - + def N(self): """Retrieve number of volume points. Shortcut for `volume.NPoints()`.""" return self.imagedata().GetNumberOfPoints() @@ -3249,7 +3250,7 @@ def xSlice(self, i): vslice.SetExtent(i,i, 0,ny, 0,nz) vslice.Update() return Actor(vslice.GetOutput()) - + def ySlice(self, j): """Extract the slice at index `j` of volume along y-axis.""" vslice = vtk.vtkImageDataGeometryFilter() @@ -3271,22 +3272,21 @@ def zSlice(self, k): vslice.SetExtent(0,nx, 0,ny, k,k) vslice.Update() return Actor(vslice.GetOutput()) - - - - - - - - - - - - - - - - - - - \ No newline at end of file + + + + + + + + + + + + + + + + + + diff --git a/vtkplotter/addons.py b/vtkplotter/addons.py index 5d0f69ad..1ce30930 100644 --- a/vtkplotter/addons.py +++ b/vtkplotter/addons.py @@ -216,9 +216,9 @@ def _cellColors(scale, scalars, cmap, alpha): tlabs = numpy.linspace(vmin, vmax, num=nlabels, endpoint=True) tacts = [] prec = (vmax - vmin) / abs(vmax + vmin) * 2 - prec = int(3 + abs(numpy.log10(prec + 1))) + prec = int(2 + abs(numpy.log10(prec + 1))) for i, t in enumerate(tlabs): - tx = utils.precision(t, prec) + tx = utils.precision(t, prec, vrange=vmax-vmin) y = -sy / 1.98 + sy * i / (nlabels - 1) a = shapes.Text(tx, pos=[sx * gap, y, 0], s=sy / 50, c=c, alpha=alpha, depth=0) a.PickableOff() @@ -1142,7 +1142,7 @@ def addAxes(axtype=None, c=None): v = (ic/rx, -xLabelOffset, 0) val = v[0]*sizes[0]+min_bns[0] if abs(val)>1 and sizes[0]<1: xLabelPrecision = int(xLabelPrecision-numpy.log10(sizes[0])) - tval = utils.precision(val, xLabelPrecision) + tval = utils.precision(val, xLabelPrecision, vrange=sizes[0]) xlab = shapes.Text(tval, pos=v, s=xLabelSize, justify="center-top", depth=0) if xKeepAspectRatio: xlab.SetScale(x_aspect_ratio_scale) labels.append(xlab.c(xTickColor).lighting(specular=0, ambient=1)) @@ -1153,7 +1153,7 @@ def addAxes(axtype=None, c=None): v = (-yLabelOffset, ic/ry, 0) val = v[1]*sizes[1]+min_bns[2] if abs(val)>1 and sizes[1]<1: yLabelPrecision = int(yLabelPrecision-numpy.log10(sizes[1])) - tval = utils.precision(val, yLabelPrecision) + tval = utils.precision(val, yLabelPrecision, vrange=sizes[1]) ylab = shapes.Text(tval, pos=(0,0,0), s=yLabelSize, justify="center-bottom", depth=0) if yKeepAspectRatio: ylab.SetScale(y_aspect_ratio_scale) ylab.rotateZ(yTitleRotation).pos(v) @@ -1164,7 +1164,7 @@ def addAxes(axtype=None, c=None): for ic in range(1, rz+enableLastLabel): v = (-zLabelOffset, -zLabelOffset, ic/rz) val = v[2]*sizes[2]+min_bns[4] - tval = utils.precision(val, zLabelPrecision) + tval = utils.precision(val, zLabelPrecision, vrange=sizes[2]) if abs(val)>1 and sizes[2]<1: zLabelPrecision = int(zLabelPrecision-numpy.log10(sizes[2])) zlab = shapes.Text(tval, pos=(0,0,0), s=zLabelSize, justify="center-bottom", depth=0) if zKeepAspectRatio: zlab.SetScale(z_aspect_ratio_scale) diff --git a/vtkplotter/analysis.py b/vtkplotter/analysis.py index bd4c4635..8690efc7 100644 --- a/vtkplotter/analysis.py +++ b/vtkplotter/analysis.py @@ -1080,7 +1080,7 @@ def surfaceIntersection(actor1, actor2, tol=1e-06): def probePoints(vol, pts): """ Takes a ``Volume`` and probes its scalars at the specified points in space. - + Note that a mask is also output with valid/invalid points which can be accessed with `actor.scalars()`. """ @@ -1966,7 +1966,7 @@ def interpolateToVolume(actor, kernel='shepard', radius=None, :param list bounds: bounding box of the output Volume object :param list dims: dimensions of the output Volume object :param float nullValue: value to be assigned to invalid points - + |interpolateVolume| |interpolateVolume.py|_ """ output = actor.polydata() @@ -2446,12 +2446,12 @@ def volumeFromMesh(actor, bounds=None, dims=(20,20,20), signed=True, negate=Fals img.SetSpacing(sx, sy, sz) img.SetOrigin(bounds[0], bounds[2], bounds[4]) img.AllocateScalars(vtk.VTK_FLOAT, 1) - + imp = vtk.vtkImplicitPolyDataDistance() imp.SetInput(actor.polydata()) b4 = bounds[4] r2 = range(dims[2]) - + for i in range(dims[0]): x = i*sx+bounds[0] for j in range(dims[1]): @@ -2462,7 +2462,7 @@ def volumeFromMesh(actor, bounds=None, dims=(20,20,20), signed=True, negate=Fals if negate: v = -v else: - v = abs(v) + v = abs(v) img.SetScalarComponentFromFloat(i,j,k, 0, v) return Volume(img) @@ -2473,7 +2473,7 @@ def computeNormalsWithPCA(actor, n=20, orientationPoint=None, negate=False): Basically this estimates a local tangent plane around each sample point p by considering a small neighborhood of points around p, and fitting a plane to the neighborhood (via PCA). - + :param int n: neighborhood size to calculate the normal :param list orientationPoint: adjust the +/- sign of the normals so that the normals all point towards a specified point. If None, perform a traversal @@ -2485,20 +2485,20 @@ def computeNormalsWithPCA(actor, n=20, orientationPoint=None, negate=False): pcan = vtk.vtkPCANormalEstimation() pcan.SetInputData(poly) pcan.SetSampleSize(n) - + if orientationPoint is not None: pcan.SetNormalOrientationToPoint() pcan.SetOrientationPoint(orientationPoint) else: pcan.SetNormalOrientationToGraphTraversal() - + if negate: pcan.FlipNormalsOn() - + pcan.Update() out = pcan.GetOutput() vnorm = out.GetPointData().GetNormals() - + newact = actor.clone() newact.polydata().GetPointData().SetNormals(vnorm) newact.polydata().GetPointData().Modified() @@ -2527,13 +2527,13 @@ def pointDensity(actor, dims=(30,30,30), bounds=None, radius=None, computeGradie img = pdf.GetOutput() vol = Volume(img) return volumeOperation(vol, '/', img.GetScalarRange()[1]) - - - - - - - - + + + + + + + + diff --git a/vtkplotter/colors.py b/vtkplotter/colors.py index cc9cdeea..f4ba9337 100644 --- a/vtkplotter/colors.py +++ b/vtkplotter/colors.py @@ -149,23 +149,23 @@ # available colormap names from matplotlib: # Accent, Accent_r, Blues, Blues_r, BrBG, BrBG_r, BuGn, BuGn_r, BuPu, BuPu_r, CMRmap, -# CMRmap_r, Dark2, Dark2_r, GnBu, GnBu_r, Greens, Greens_r, Greys, Greys_r, OrRd, -# OrRd_r, Oranges, Oranges_r, PRGn, PRGn_r, Paired, Paired_r, Pastel1, Pastel1_r, -# Pastel2, Pastel2_r, PiYG, PiYG_r, PuBu, PuBuGn, PuBuGn_r, PuBu_r, PuOr, PuOr_r, -# PuRd, PuRd_r, Purples, Purples_r, RdBu, RdBu_r, RdGy, RdGy_r, RdPu, RdPu_r, RdYlBu, -# RdYlBu_r, RdYlGn, RdYlGn_r, Reds, Reds_r, Set1, Set1_r, Set2, Set2_r, Set3, Set3_r, -# Spectral, Spectral_r, Wistia, Wistia_r, YlGn, YlGnBu, YlGnBu_r, YlGn_r, YlOrBr, +# CMRmap_r, Dark2, Dark2_r, GnBu, GnBu_r, Greens, Greens_r, Greys, Greys_r, OrRd, +# OrRd_r, Oranges, Oranges_r, PRGn, PRGn_r, Paired, Paired_r, Pastel1, Pastel1_r, +# Pastel2, Pastel2_r, PiYG, PiYG_r, PuBu, PuBuGn, PuBuGn_r, PuBu_r, PuOr, PuOr_r, +# PuRd, PuRd_r, Purples, Purples_r, RdBu, RdBu_r, RdGy, RdGy_r, RdPu, RdPu_r, RdYlBu, +# RdYlBu_r, RdYlGn, RdYlGn_r, Reds, Reds_r, Set1, Set1_r, Set2, Set2_r, Set3, Set3_r, +# Spectral, Spectral_r, Wistia, Wistia_r, YlGn, YlGnBu, YlGnBu_r, YlGn_r, YlOrBr, # YlOrBr_r, YlOrRd, YlOrRd_r, afmhot, afmhot_r, autumn, autumn_r, binary, binary_r -# bone, bone_r, brg, brg_r, bwr, bwr_r, cividis, cividis_r, cool, cool_r, coolwarm, -# coolwarm_r, copper, copper_r, cubehelix, cubehelix_r, flag, flag_r, -# gist_earth, gist_earth_r, gist_gray, gist_gray_r, gist_heat, gist_heat_r, -# gist_ncar, gist_ncar_r, gist_rainbow, gist_rainbow_r, gist_stern, gist_stern_r, -# gist_yarg, gist_yarg_r, gnuplot, gnuplot2, gnuplot2_r, gnuplot_r, gray, gray_r, -# hot, hot_r, hsv, hsv_r, inferno, inferno_r, jet, jet_r, magma, magma_r, -# nipy_spectral, nipy_spectral_r, ocean, ocean_r, pink, pink_r, plasma, plasma_r, -# prism, prism_r, rainbow, rainbow_r, seismic, seismic_r, spring, spring_r, summer, -# summer_r, tab10, tab10_r, tab20, tab20_r, tab20b, tab20b_r, tab20c, tab20c_r, -# terrain, terrain_r, twilight, twilight_r, twilight_shifted, twilight_shifted_r, +# bone, bone_r, brg, brg_r, bwr, bwr_r, cividis, cividis_r, cool, cool_r, coolwarm, +# coolwarm_r, copper, copper_r, cubehelix, cubehelix_r, flag, flag_r, +# gist_earth, gist_earth_r, gist_gray, gist_gray_r, gist_heat, gist_heat_r, +# gist_ncar, gist_ncar_r, gist_rainbow, gist_rainbow_r, gist_stern, gist_stern_r, +# gist_yarg, gist_yarg_r, gnuplot, gnuplot2, gnuplot2_r, gnuplot_r, gray, gray_r, +# hot, hot_r, hsv, hsv_r, inferno, inferno_r, jet, jet_r, magma, magma_r, +# nipy_spectral, nipy_spectral_r, ocean, ocean_r, pink, pink_r, plasma, plasma_r, +# prism, prism_r, rainbow, rainbow_r, seismic, seismic_r, spring, spring_r, summer, +# summer_r, tab10, tab10_r, tab20, tab20_r, tab20b, tab20b_r, tab20c, tab20c_r, +# terrain, terrain_r, twilight, twilight_r, twilight_shifted, twilight_shifted_r, # viridis, viridis_r, winter, winter_r @@ -630,17 +630,19 @@ def printc(*strings, **keys): Available colors are: black, red, green, yellow, blue, magenta, cyan, white. - :param c: foreground color [''] - :param bc: background color [''] - :param hidden: do not show text [False] - :param bold: boldface [True] - :param blink: blinking text [False] - :param underline: underline text [False] - :param dim: make text look dimmer [False] - :param invert: invert background anf forward colors [False] + :param c: foreground color + :param bc: background color + :param bool hidden: do not show text [False] + :param bool bold: boldface [True] + :param bool italic: italic [False] + :param bool blink: blinking text [False] + :param bool underline: underline text [False] + :param bool strike: strike through text [False] + :param bool dim: make text look dimmer [False] + :param bool invert: invert background anf forward colors [False] :param box: print a box with specified text character [''] - :param flush: flush buffer after printing [True] - :param end: end character to be printed [return] + :param bool flush: flush buffer after printing [True] + :param str end: end character to be printed [newline] :Example: .. code-block:: python @@ -674,8 +676,10 @@ def printc(*strings, **keys): bc = keys.pop("bc", None) hidden = keys.pop("hidden", False) bold = keys.pop("bold", True) + italic = keys.pop("italic", False) blink = keys.pop("blink", False) underline = keys.pop("underline", False) + strike = keys.pop("strike", False) dim = keys.pop("dim", False) invert = keys.pop("invert", False) box = keys.pop("box", "") @@ -730,12 +734,16 @@ def printc(*strings, **keys): cseq += "\x1b[" + str(40 + cb) + "m" if underline and not box: special += "\x1b[4m" + if strike and not box: + special += "\x1b[9m" if dim: special += "\x1b[2m" if invert: special += "\x1b[7m" if bold: special += "\x1b[1m" + if italic: + special += "\x1b[3m" if blink: special += "\x1b[5m" diff --git a/vtkplotter/dolfin.py b/vtkplotter/dolfin.py index 2cf27699..e93e7e31 100644 --- a/vtkplotter/dolfin.py +++ b/vtkplotter/dolfin.py @@ -137,6 +137,8 @@ def _inputsort(obj): + import dolfin + u = None mesh = None if not utils.isSequence(obj): @@ -153,44 +155,22 @@ def _inputsort(obj): if "MeshFunction" in inputtype: mesh = ob.mesh() - import dolfin - V = dolfin.FunctionSpace(mesh, "CG", 1) - u = dolfin.Function(V) - #print(mesh.cells()) - #print(len(mesh.cells()), len(mesh.coordinates()), len(ob.array())) - #print(mesh.num_cells()) - #print(u.vector()[:]) - - v2d = dolfin.vertex_to_dof_map(V) - u.vector()[v2d] = ob.array() - -# r = ob.dim() -# if r == 0: -# V = dolfin.FunctionSpace(mesh, "CG", 1) -# elif r == 1: -# V = dolfin.VectorFunctionSpace(mesh, "CG", 1, dim=r) -# else: -# V = dolfin.TensorFunctionSpace(mesh, "CG", 1, shape=(r,r)) -# except: -# printc('~times Sorry could not deal with your MeshFunction', c=1) -# return None -# tdim = mesh.topology().dim() -# d = ob.dim() -# if tdim == 2 and d == 2: -# import matplotlib.tri as tri -# xy = mesh.coordinates() -# mh = buildPolyData(xy, mesh.cells()) -# show(mh) -# print( tri.Triangulation(xy[:, 0], xy[:, 1], mesh.cells()) ) -# exit() - + if ob.dim()>0: + printc('MeshFunction of dim>0 not supported.', c=1) + printc('Try e.g.: MeshFunction("size_t", mesh, 0)', c=1, italic=1) + printc('instead of MeshFunction("size_t", mesh, 1)', c=1, strike=1) + else: + #printc(ob.dim(), mesh.num_cells(), len(mesh.coordinates()), len(ob.array())) + V = dolfin.FunctionSpace(mesh, "CG", 1) + u = dolfin.Function(V) + v2d = dolfin.vertex_to_dof_map(V) + u.vector()[v2d] = ob.array() elif "Function" in inputtype or "Expression" in inputtype: u = ob elif "Mesh" in inputtype: mesh = ob if "str" in inputtype: - import dolfin mesh = dolfin.Mesh(ob) if u and not mesh and hasattr(u, "function_space"): @@ -200,10 +180,6 @@ def _inputsort(obj): if u and not mesh and hasattr(u, "mesh"): mesh = u.mesh() - if not mesh: - printc("~times Error: dolfin mesh is not defined.", c=1) - raise RuntimeError() - #printc('------------------------------------') #printc('mesh.topology dim=', mesh.topology().dim()) #printc('mesh.geometry dim=', mesh.geometry().dim()) @@ -313,17 +289,17 @@ def plot(*inputobj, **options): - distance `(float)`, set the focal point to the specified distance from the camera position. - clippingRange `(float)`, distance of the near and far clipping planes along the direction of projection. - + - parallelScale `(float)`, scaling used for a parallel projection, i.e. the height of the viewport in world-coordinate distances. The default is 1. Note that the "scale" parameter works as an "inverse scale", larger numbers produce smaller images. This method has no effect in perspective projection mode. - + - thickness `(float)`, set the distance between clipping planes. This method adjusts the far clipping plane to be set a distance 'thickness' beyond the near clipping plane. - + - viewAngle `(float)`, the camera view angle, which is the angular height of the camera view measured in degrees. The default angle is 30 degrees. @@ -489,7 +465,7 @@ def plot(*inputobj, **options): if add and settings.plotter_instance: actors = settings.plotter_instance.actors - if 'mesh' in mode or 'color' in mode or 'warp' in mode or 'displac' in mode: + if mesh and ('mesh' in mode or 'color' in mode or 'warp' in mode or 'displac' in mode): if 'warp' in mode: #deprecation printc("~bomb Please use 'displacement' instead of 'warp' in mode!", c=1) @@ -625,9 +601,7 @@ def plot(*inputobj, **options): class MeshActor(Actor): """MeshActor, a vtkActor derived object for dolfin support.""" - def __init__( - self, *inputobj, **options # c="gold", alpha=1, wire=True, bc=None, computeNormals=False - ): + def __init__(self, *inputobj, **options): c = options.pop("c", "gold") alpha = options.pop("alpha", 1) @@ -635,7 +609,9 @@ def __init__( bc = options.pop("bc", None) computeNormals = options.pop("computeNormals", False) - mesh, u = _inputsort(inputobj) + mesh, u = _inputsort(inputobj) + if not mesh: + return poly = vtkio.buildPolyData(mesh) @@ -691,6 +667,8 @@ def MeshPoints(*inputobj, **options): alpha = options.pop("alpha", 1) mesh, u = _inputsort(inputobj) + if not mesh: + return None plist = mesh.coordinates() if u: u_values = np.array([u(p) for p in plist]) @@ -731,6 +709,9 @@ def MeshLines(*inputobj, **options): alpha = options.pop("alpha", 1) mesh, u = _inputsort(inputobj) + if not mesh: + return None + startPoints = mesh.coordinates() u_values = np.array([u(p) for p in mesh.coordinates()]) if not utils.isSequence(u_values[0]): @@ -766,6 +747,9 @@ def MeshArrows(*inputobj, **options): res = options.pop("res", 12) mesh, u = _inputsort(inputobj) + if not mesh: + return None + startPoints = mesh.coordinates() u_values = np.array([u(p) for p in mesh.coordinates()]) if not utils.isSequence(u_values[0]): @@ -787,98 +771,6 @@ def MeshArrows(*inputobj, **options): -#def make_mapping(sub_space, super_space): -# from scipy.spatial import cKDTree -# super_dof_coor = super_space.tabulate_dof_coordinates() -# sub_dof_coor = sub_space.tabulate_dof_coordinates() -# -# tree = cKDTree(super_dof_coor) -# _,mapping = tree.query(sub_dof_coor, k=1) -# return mapping - - - - -#from dolfin import * -#mesh = UnitCubeMesh(10, 10, 10) # this will be the "grandparent" mesh -#mesh.coordinates()[:,0] -= .5 # shift x-coords -#mesh.coordinates()[:,1] -= .5 # shift y-coords -# -#bmesh = BoundaryMesh(mesh, "exterior") # surface boundary mesh -# -## mark the cells on the bottom of the bmesh surface by iterating over bmesh cells. -## Look up ccorresponding facet in mesh and mark if facet normal points in -z direction -#cellmap = bmesh.entity_map(2) -#vertmap = bmesh.entity_map(0) -#pb = MeshFunction("size_t", bmesh, dim =1) -#for c in cells(bmesh): -# if Facet(mesh, cellmap[c.index()]).normal().z() < 0: -# pb[c] = 1 -# -## use the marked bottom of the bmesh to create a Submesh -#submesh = SubMesh(bmesh, pb, 1) # bottom of boundary mesh -# -## FunctionSpaces on main mesh, bmesh, submesh -#V = FunctionSpace(mesh, "CG", 1) # mesh function space -#Vb = FunctionSpace(bmesh, "CG", 1) # surface function space -#Vs = FunctionSpace(submesh, "CG", 1) # submesh function space -# -## mappings we may need: -#m = vertex_to_dof_map(V) -#b = vertex_to_dof_map(Vb) -#s = vertex_to_dof_map(Vs) -# -#mi = dof_to_vertex_map(V) -#bi = dof_to_vertex_map(Vb) -#si = dof_to_vertex_map(Vs) -# -#t = submesh.data().array('parent_vertex_indices', 0) # mapping from submesh back to bmesh -# -## Functions on main mesh, bmesh, and submesh -#u = Function(V) -#ub = Function(Vb) # boundary function -#us = Function(Vs) # surface function -# -## Interpolate the following expr onto u, ub, us -#expr = Expression('sqrt(pow(x[0],2) + pow(x[1], 2))', degree=2) -#u.interpolate(expr) -#ub.interpolate(expr) -#us.interpolate(expr) -# -## Some empty function to test the mappings -#ub_test = Function(Vb) # empty bmesh function -#u_test = Function(V) # empty mesh function -# -## Mapping from submesh to bmesh (works)... Any way to avoid calls to get_local/set_local?? -#us_a = us.vector().get_local() # origin array -#ub_test_a = ub_test.vector().get_local() # destination array -#ub_test_a[b[t]] = us_a[s] # transfer -#ub_test.vector().set_local(ub_test_a) # set destination function values -# -## Mapping from submesh to grandparent mesh (less than satisfying solution to question in fenics forum link) -## Bonus points for getting this kind of map composition to work: -## u_test_a = u_test.vector().get_local() # destination array -## u_test_a[m[b[t]]] = us_a[s] # transfer ( not working ) -## u_test.vector().set_local(u_test_a) -#for Vs_dof, val in enumerate(us.vector().get_local()): -# submesh_vertex = si[Vs_dof] -# boundary_vertex = t[submesh_vertex] -# mesh_vertex = vertmap[int(boundary_vertex)] # np.uint not accepted -# V_dof = m[mesh_vertex] -# u_test.vector()[V_dof] = val -# -#u.rename('u','u') -#ub_test.rename('ub_test','ub_test') -#u_test.rename('u_test','u_test') -#us.rename('us','us') - - - - - - - -# #import dolfin #import dolfin.cpp as cpp #import ufl diff --git a/vtkplotter/plotter.py b/vtkplotter/plotter.py index 7b06dc37..ea8df180 100644 --- a/vtkplotter/plotter.py +++ b/vtkplotter/plotter.py @@ -26,32 +26,7 @@ "closeWindow", "closePlotter", "interactive"] ######################################################################## -def show(*actors, **options -# at=None, -# shape=(1, 1), -# N=None, -# pos=(0, 0), -# size="auto", -# screensize="auto", -# title="", -# bg="blackboard", -# bg2=None, -# axes=4, -# infinity=False, -# verbose=True, -# interactive=None, -# offscreen=False, -# resetcam=True, -# zoom=None, -# viewup="", -# azimuth=0, -# elevation=0, -# roll=0, -# interactorStyle=0, -# newPlotter=False, -# depthpeeling=False, -# q=False, -): +def show(*actors, **options): """ Create on the fly an instance of class ``Plotter`` and show the object(s) provided. @@ -94,17 +69,17 @@ def show(*actors, **options - distance `(float)`, set the focal point to the specified distance from the camera position. - clippingRange `(float)`, distance of the near and far clipping planes along the direction of projection. - + - parallelScale `(float)`, scaling used for a parallel projection, i.e. the height of the viewport in world-coordinate distances. The default is 1. Note that the "scale" parameter works as an "inverse scale", larger numbers produce smaller images. This method has no effect in perspective projection mode. - + - thickness `(float)`, set the distance between clipping planes. This method adjusts the far clipping plane to be set a distance 'thickness' beyond the near clipping plane. - + - viewAngle `(float)`, the camera view angle, which is the angular height of the camera view measured in degrees. The default angle is 30 degrees. @@ -680,7 +655,7 @@ def add(self, actors, render=True): self.actors.append(a) if render and self.renderer: self.renderer.AddActor(a) - if render and self.interactor: + if render and self.interactor: self.interactor.Render() return None else: @@ -706,7 +681,7 @@ def remove(self, actors): del self.actors[i] #################################################### - def load(self, inputobj, c="gold", alpha=1, threshold=False, spacing=(), unpack=True): + def load(self, inputobj, c=None, alpha=1, threshold=False, spacing=(), unpack=True): """ Load Actors and Volumes from file. The output will depend on the file extension. See examples below. @@ -1194,7 +1169,7 @@ def show( :param float azimuth/elevation/roll: move camera accordingly :param str viewup: either ['x', 'y', 'z'] or a vector to set vertical direction :param bool resetcam: re-adjust camera position to fit objects - :param dict camera: Camera parameters can further be specified with a dictionary assigned + :param dict camera: Camera parameters can further be specified with a dictionary assigned to the ``camera`` keyword (E.g. `show(camera={'pos':(1,2,3), 'thickness':1000,})`) - pos, `(list)`, the position of the camera in world coordinates @@ -1203,17 +1178,17 @@ def show( - distance `(float)`, set the focal point to the specified distance from the camera position. - clippingRange `(float)`, distance of the near and far clipping planes along the direction of projection. - + - parallelScale `(float)`, scaling used for a parallel projection, i.e. the height of the viewport in world-coordinate distances. The default is 1. Note that the "scale" parameter works as an "inverse scale", larger numbers produce smaller images. This method has no effect in perspective projection mode. - + - thickness `(float)`, set the distance between clipping planes. This method adjusts the far clipping plane to be set a distance 'thickness' beyond the near clipping plane. - + - viewAngle `(float)`, the camera view angle, which is the angular height of the camera view measured in degrees. The default angle is 30 degrees. @@ -1544,8 +1519,10 @@ def scan(wannabeacts): vtkscals = vtkdata.GetScalars() if vtkscals is not None: + if not vtkscals.GetName(): + vtkscals.SetName('scalars') scals_min, scals_max = ia.mapper.GetScalarRange() - color_attribute=(vtkscals.GetName(), scals_min, scals_max) + color_attribute = (vtkscals.GetName(), scals_min, scals_max) lut = ia.mapper.GetLookupTable() lut.Build() kcmap=[] diff --git a/vtkplotter/shapes.py b/vtkplotter/shapes.py index 7cd22ee9..16c5ca13 100644 --- a/vtkplotter/shapes.py +++ b/vtkplotter/shapes.py @@ -85,19 +85,24 @@ def _colorPoints(plist, cols, r, alpha): vgf.SetInputData(src.GetOutput()) vgf.Update() pd = vgf.GetOutput() + ucols = vtk.vtkUnsignedCharArray() ucols.SetNumberOfComponents(3) ucols.SetName("pointsRGB") + for i in range(len(plist)): c = np.array(colors.getColor(cols[i])) * 255 ucols.InsertNextTuple3(c[0], c[1], c[2]) + + if len(plist[0]) == 2: #make it 3d + plist = np.c_[np.array(plist), np.zeros(len(plist))] + pd.GetPoints().SetData(numpy_to_vtk(plist, deep=True)) pd.GetPointData().SetScalars(ucols) actor = Actor(pd, c, alpha) actor.mapper.ScalarVisibilityOn() actor.GetProperty().SetInterpolationToFlat() actor.GetProperty().SetPointSize(r) - settings.collectable_actors.append(actor) return actor n = len(plist) @@ -134,7 +139,7 @@ def _colorPoints(plist, cols, r, alpha): pd = vtk.vtkPolyData() pd.SetPoints(sourcePoints) pd.SetVerts(sourceVertices) - + if n == 1: # passing just one point pd.GetPoints().SetPoint(0, [0, 0, 0]) else: @@ -173,7 +178,7 @@ def Glyph(actor, glyphObj, orientationArray=None, if c in list(colors._mapscales.cmap_d.keys()): cmap = c c = None - + if tol: actor = actor.clone().clean(tol) poly = actor.polydata() @@ -255,32 +260,32 @@ def Tensors(domain, source='ellipsoid', useEigenValues=True, isSymmetric=True, """Geometric representation of tensors defined on a domain or set of points. Tensors can be scaled and/or rotated according to the source at eache input point. Scaling and rotation is controlled by the eigenvalues/eigenvectors of the symmetrical part - of the tensor as follows: - - For each tensor, the eigenvalues (and associated eigenvectors) are sorted - to determine the major, medium, and minor eigenvalues/eigenvectors. + of the tensor as follows: + + For each tensor, the eigenvalues (and associated eigenvectors) are sorted + to determine the major, medium, and minor eigenvalues/eigenvectors. The eigenvalue decomposition only makes sense for symmetric tensors, hence the need to only consider the symmetric part of the tensor, which is 1/2*(T+T.transposed()). :param str source: preset type of source shape ['ellipsoid', 'cylinder', 'cube' or any specified ``Actor``] - + :param bool useEigenValues: color source glyph using the eigenvalues or by scalars. - + :param bool threeAxes: if `False` scale the source in the x-direction, the medium in the y-direction, - and the minor in the z-direction. Then, the source is rotated so that the glyph's local x-axis lies + and the minor in the z-direction. Then, the source is rotated so that the glyph's local x-axis lies along the major eigenvector, y-axis along the medium eigenvector, and z-axis along the minor. - + If `True` three sources are produced, each of them oriented along an eigenvector and scaled according to the corresponding eigenvector. - + :param bool isSymmetric: If `True` each source glyph is mirrored (2 or 6 glyphs will be produced). The x-axis of the source glyph will correspond to the eigenvector on output. - + :param float length: distance from the origin to the tip of the source glyph along the x-axis - + :param float scale: scaling factor of the source glyph. :param float maxScale: clamp scaling at this factor. - + |tensors| |tensors.py|_ """ if 'ellip' in source: @@ -296,11 +301,11 @@ def Tensors(domain, source='ellipsoid', useEigenValues=True, isSymmetric=True, else: src = source.normalize().polydata(False) src.Update() - + tg = vtk.vtkTensorGlyph() tg.SetInputData(domain.GetMapper().GetInput()) tg.SetSourceData(src.GetOutput()) - + if c is None: tg.ColorGlyphsOn() else: diff --git a/vtkplotter/utils.py b/vtkplotter/utils.py index 18c0f599..9a395068 100644 --- a/vtkplotter/utils.py +++ b/vtkplotter/utils.py @@ -204,7 +204,7 @@ def tryint(s): l.sort(key=alphanum_key) return None # NB: input list is modified - + def lin_interp(x, rangeX, rangeY): """ Interpolate linearly variable x in rangeX onto rangeY. @@ -249,6 +249,8 @@ def precision(x, p, vrange=None): """ Returns a string representation of `x` formatted with precision `p`. + :param float vrange: range in which x exists (to snap it to '0' if below precision). + Based on the webkit javascript implementation taken `from here `_, and implemented by `randlet `_. @@ -257,8 +259,8 @@ def precision(x, p, vrange=None): x = float(x) - if x == 0.0: - return "0." + "0" * (p - 1) + if x == 0.0 or (vrange is not None and abs(x) < vrange/pow(10,p)): + return "0" out = [] if x < 0: @@ -584,14 +586,14 @@ def printvtkactor(actor, tab=""): img = obj.GetMapper().GetInput() colors.printc(" position: ", c="b", bold=1, end="") colors.printc(pos, c="b", bold=0) - + colors.printc(" dimensions: ", c="b", bold=1, end="") colors.printc(img.GetDimensions(), c="b", bold=0) colors.printc(" spacing: ", c="b", bold=1, end="") colors.printc(img.GetSpacing(), c="b", bold=0) colors.printc(" data dimension: ", c="b", bold=1, end="") colors.printc(img.GetDataDimension(), c="b", bold=0) - + colors.printc(" memory size: ", c="b", bold=1, end="") colors.printc(int(img.GetActualMemorySize()/1024), 'Mb', c="b", bold=0) @@ -609,9 +611,9 @@ def printvtkactor(actor, tab=""): colors.printc(" scalar range: ", c="b", bold=1, end="") colors.printc(img.GetScalarRange(), c="b", bold=0) - printHistogram(obj, horizontal=True, - logscale=True, bins=8, height=15, c='b', bold=0) - + printHistogram(obj, horizontal=True, + logscale=True, bins=8, height=15, c='b', bold=0) + elif hasattr(obj, "interactor"): # dumps Plotter info axtype = { 0: "(no axes)", diff --git a/vtkplotter/vtkio.py b/vtkplotter/vtkio.py index c721ef9e..f8fbf6cb 100644 --- a/vtkplotter/vtkio.py +++ b/vtkplotter/vtkio.py @@ -34,7 +34,7 @@ ] -def load(inputobj, c="gold", alpha=1, threshold=False, spacing=(), unpack=True): +def load(inputobj, c=None, alpha=1, threshold=False, spacing=(), unpack=True): """ Load ``Actor`` and ``Volume`` from file. @@ -116,7 +116,7 @@ def load(inputobj, c="gold", alpha=1, threshold=False, spacing=(), unpack=True): if len(spacing) == 3: image.SetSpacing(spacing[0], spacing[1], spacing[2]) if threshold is False: - if c is "gold" and alpha is 1: + if c is None and alpha is 1: c = ['b','lb','lg','y','r'] # good for blackboard background alpha = (0.0, 0.0, 0.2, 0.4, 0.8, 1) #c = ['lb','db','dg','dr'] # good for white backgr @@ -168,7 +168,7 @@ def _load_file(filename, c, alpha, threshold, spacing, unpack): or fl.endswith(".mhd") or fl.endswith(".nrrd") or fl.endswith(".nii"): img = loadImageData(filename, spacing) if threshold is False: - if c is "gold" and alpha is 1: + if c is None and alpha is 1: c = ['b','lb','lg','y','r'] # good for blackboard background alpha = (0.0, 0.0, 0.2, 0.4, 0.8, 1) #c = ['lb','db','dg','dr'] # good for white backgr @@ -216,11 +216,27 @@ def _load_file(filename, c, alpha, threshold, spacing, unpack): elif fl.endswith(".geojson") or fl.endswith(".geojson.gz"): return loadGeoJSON(fl) - + ################################################################# polygonal mesh: else: - if fl.endswith(".vtk"): +# applygeomf = False + if fl.endswith(".vtk"): # legacy readers +# try: +# with open(filename) as myfile: +# head = [next(myfile) for x in range(3)] +# head = ''.join(head).lower() +# if ' unstructured_grid' in head: +# reader = vtk.vtkUnstructuredGridReader() +# applygeomf = True +# elif ' structured_grid' in head: +# reader = vtk.vtkStructuredGridReader() +# applygeomf = True +# else: +# reader = vtk.vtkPolyDataReader() +# except: + reader = vtk.vtkPolyDataReader() + elif fl.endswith(".ply"): reader = vtk.vtkPLYReader() elif fl.endswith(".obj"): @@ -237,6 +253,8 @@ def _load_file(filename, c, alpha, threshold, spacing, unpack): reader = vtk.vtkXMLStructuredGridReader() elif fl.endswith(".vtu"): reader = vtk.vtkXMLUnstructuredGridReader() + elif fl.endswith(".vtr"): + reader = vtk.vtkXMLRectilinearGridReader() elif fl.endswith(".txt"): reader = vtk.vtkParticleReader() # (format is x, y, z, scalar) elif fl.endswith(".xyz"): @@ -251,7 +269,7 @@ def _load_file(filename, c, alpha, threshold, spacing, unpack): reader = vtk.vtkXMLPImageDataReader() else: reader = vtk.vtkDataReader() - + reader.SetFileName(filename) if hasattr(reader, 'ReadAllScalarsOn'): reader.ReadAllScalarsOn() @@ -261,12 +279,6 @@ def _load_file(filename, c, alpha, threshold, spacing, unpack): reader.Update() poly = reader.GetOutput() - if fl.endswith(".vts") or fl.endswith(".vtu"): # un/structured grid - gf = vtk.vtkGeometryFilter() - gf.SetInputData(poly) - gf.Update() - poly = gf.GetOutput() - if not poly: colors.printc("~noentry Unable to load", filename, c=1) return None @@ -436,7 +448,7 @@ def loadGeoJSON(filename): jr.SetFileName(filename) jr.Update() return Actor(jr.GetOutput()) - + def loadDolfin(filename): """Reads a `Fenics/Dolfin` file format. Return an ``Actor(vtkActor)`` object."""