# Source code for ase.io.pov

```"""
Module for povray file format support.

See http://www.povray.org/ for details on the format.
"""
from collections.abc import Sequence
from subprocess import check_call, DEVNULL
from pathlib import Path

import numpy as np

from ase.io.utils import PlottingVariables
from ase.constraints import FixAtoms
from ase import Atoms

def pa(array):
"""Povray array syntax"""
return '<' + ', '.join(f"{x:>6.2f}" for x in tuple(array)) + '>'

def pc(array):
"""Povray color syntax"""
if isinstance(array, str):
return 'color ' + array
if isinstance(array, float):
return f'rgb <{array:.2f}>*3'.format(array)
L = len(array)
if L > 2 and L < 6:
return f"rgb{'' if L == 3 else 't' if L == 4 else 'ft'} <" +\
', '.join(f"{x:.2f}" for x in tuple(array)) + '>'

"""Get all pairs of bonding atoms

Return all pairs of atoms which are closer than radius times the
sum of their respective covalent radii.  The pairs are returned as
tuples::

(a, b, (i1, i2, i3))

so that atoms a bonds to atom b displaced by the vector::

_     _     _
i c + i c + i c ,
1 1   2 2   3 3

where c1, c2 and c3 are the unit cell vectors and i1, i2, i3 are
integers."""

from ase.neighborlist import NeighborList
nl = NeighborList(cutoffs=cutoffs, self_interaction=False)
nl.update(atoms)
bondpairs = []
for a in range(len(atoms)):
indices, offsets = nl.get_neighbors(a)
bondpairs.extend([(a, a2, offset)
for a2, offset in zip(indices, offsets)])
return bondpairs

def set_high_bondorder_pairs(bondpairs, high_bondorder_pairs=None):
"""Set high bondorder pairs

Modify bondpairs list (from get_bondpairs((atoms)) to include high
bondorder pairs.

Parameters:
-----------
bondpairs: List of pairs, generated from get_bondpairs(atoms)
high_bondorder_pairs: Dictionary of pairs with high bond orders
using the following format:
{ ( a1, b1 ): ( offset1, bond_order1, bond_offset1),
( a2, b2 ): ( offset2, bond_order2, bond_offset2),
...
}
offset, bond_order, bond_offset are optional.
However, if they are provided, the 1st value is
offset, 2nd value is bond_order,
3rd value is bond_offset """

if high_bondorder_pairs is None:
high_bondorder_pairs = dict()
bondpairs_ = []
for pair in bondpairs:
(a, b) = (pair, pair)
if (a, b) in high_bondorder_pairs.keys():
bondpair = [a, b] + [item for item in high_bondorder_pairs[(a, b)]]
bondpairs_.append(bondpair)
elif (b, a) in high_bondorder_pairs.keys():
bondpair = [a, b] + [item for item in high_bondorder_pairs[(b, a)]]
bondpairs_.append(bondpair)
else:
bondpairs_.append(pair)
return bondpairs_

class POVRAY:
# These new styles were an attempt to port the old styles o the correct
# gamma, many had or still have unphysical light properties inorder to
# acheive a certain look.
material_styles_dict = dict(
simple='finish {phong 0.7 ambient 0.4 diffuse 0.55}',
# In general, 'pale' doesn't conserve energy and can look
# strange in many cases.
pale=('finish {ambient 0.9 diffuse 0.30 roughness 0.001 '
'specular 0.2 }'),
intermediate=('finish {ambient 0.4 diffuse 0.6 specular 0.1 '
'roughness 0.04}'),
vmd=(
'finish {ambient 0.2 diffuse 0.80 phong 0.25 phong_size 10.0 '
'specular 0.2 roughness 0.1}'),
jmol=('finish {ambient 0.4 diffuse 0.6 specular 1 roughness 0.001 '
'metallic}'),
ase2=('finish {ambient 0.2 brilliance 3 diffuse 0.6 metallic '
'specular 0.7 roughness 0.04 reflection 0.15}'),
ase3=('finish {ambient 0.4 brilliance 2 diffuse 0.6 metallic '
'specular 1.0 roughness 0.001 reflection 0.0}'),
glass=('finish {ambient 0.4 diffuse 0.35 specular 1.0 '
'roughness 0.001}'),
glass2=('finish {ambient 0.3 diffuse 0.3 specular 1.0 '
'reflection 0.25 roughness 0.001}'),
)

# These styles were made when assumed_gamma was 1.0 which gives poor color
# reproduction, the correct gamma is 2.2 for the sRGB standard.
material_styles_dict_old = dict(
simple='finish {phong 0.7}',
pale=('finish {ambient 0.5 diffuse 0.85 roughness 0.001 '
'specular 0.200 }'),
intermediate=('finish {ambient 0.3 diffuse 0.6 specular 0.1 '
'roughness 0.04}'),
vmd=('finish {ambient 0.0 diffuse 0.65 phong 0.1 phong_size 40.0 '
'specular 0.5 }'),
jmol=('finish {ambient 0.2 diffuse 0.6 specular 1 roughness 0.001 '
'metallic}'),
ase2=('finish {ambient 0.05 brilliance 3 diffuse 0.6 metallic '
'specular 0.7 roughness 0.04 reflection 0.15}'),
ase3=('finish {ambient 0.15 brilliance 2 diffuse 0.6 metallic '
'specular 1.0 roughness 0.001 reflection 0.0}'),
glass=('finish {ambient 0.05 diffuse 0.3 specular 1.0 '
'roughness 0.001}'),
glass2=('finish {ambient 0.01 diffuse 0.3 specular 1.0 '
'reflection 0.25 roughness 0.001}'),
)

def __init__(self, cell, cell_vertices, positions, diameters, colors,
image_width, image_height, constraints=tuple(), isosurfaces=[],
display=False, pause=True, transparent=True, canvas_width=None,
canvas_height=None, camera_dist=50., image_plane=None,
camera_type='orthographic', point_lights=[],
area_light=[(2., 3., 40.), 'White', .7, .7, 3, 3],
background='White', textures=None, transmittances=None,
depth_cueing=False, cue_density=5e-3,
celllinewidth=0.05, bondlinewidth=0.10, bondatoms=[],
exportconstraints=False):
"""
# x, y is the image plane, z is *out* of the screen
cell: ase.cell
cell object
cell_vertices: 2-d numpy array
contains the 8 vertices of the cell, each with three coordinates
positions: 2-d numpy array
number of atoms length array with three coordinates for positions
diameters: 1-d numpy array
diameter of atoms (in order with positions)
colors: list of str
colors of atoms (in order with positions)
image_width: float
image width in pixels
image_height: float
image height in pixels
constraints: Atoms.constraints
constraints to be visualized
isosurfaces: list of POVRAYIsosurface
composite object to write/render POVRAY isosurfaces
display: bool
display while rendering
pause: bool
pause when done rendering (only if display)
transparent: bool
make background transparent
canvas_width: int
width of canvas in pixels
canvas_height: int
height of canvas in pixels
camera_dist: float
distance from camera to front atom
image_plane: float
distance from front atom to image plane
camera_type: str
if 'orthographic' perspective, ultra_wide_angle
point_lights: list of 2-element sequences
like [[loc1, color1], [loc2, color2],...]
area_light: 3-element sequence of location (3-tuple), color (str),
width (float), height (float),
Nlamps_x (int), Nlamps_y (int)
example [(2., 3., 40.), 'White', .7, .7, 3, 3]
background: str
color specification, e.g., 'White'
textures: list of str
length of atoms list of texture names
transmittances: list of floats
length of atoms list of transmittances of the atoms
depth_cueing: bool
whether or not to use depth cueing a.k.a. fog
use with care - adjust the camera_distance to be closer
cue_density: float
if there is depth_cueing, how dense is it (how dense is the fog)
celllinewidth: float
radius of the cylinders representing the cell (Ang.)
bondlinewidth: float
radius of the cylinders representing bonds (Ang.)
bondatoms: list of lists (polymorphic)
[[atom1, atom2], ... ] pairs of bonding atoms
For bond order > 1 = [[atom1, atom2, offset,
bond_order, bond_offset],
... ]
bond_order: 1, 2, 3 for single, double,
and triple bond
bond_offset: vector for shifting bonds from
original position. Coordinates are
in Angstrom unit.
exportconstraints: bool
honour FixAtoms and mark?"""

# attributes from initialization
self.area_light = area_light
self.background = background
self.bondatoms = bondatoms
self.bondlinewidth = bondlinewidth
self.camera_dist = camera_dist
self.camera_type = camera_type
self.celllinewidth = celllinewidth
self.cue_density = cue_density
self.depth_cueing = depth_cueing
self.display = display
self.exportconstraints = exportconstraints
self.isosurfaces = isosurfaces
self.pause = pause
self.point_lights = point_lights
self.textures = textures
self.transmittances = transmittances
self.transparent = transparent

self.image_width = image_width
self.image_height = image_height
self.colors = colors
self.cell = cell
self.diameters = diameters

# calculations based on passed inputs

z0 = positions[:, 2].max()
self.offset = (image_width / 2, image_height / 2, z0)
self.positions = positions - self.offset

if cell_vertices is not None:
self.cell_vertices = cell_vertices - self.offset
self.cell_vertices.shape = (2, 2, 2, 3)
else:
self.cell_vertices = None

ratio = float(self.image_width) / self.image_height
if canvas_width is None:
if canvas_height is None:
self.canvas_width = min(self.image_width * 15, 640)
self.canvas_height = min(self.image_height * 15, 640)
else:
self.canvas_width = canvas_height * ratio
self.canvas_height = canvas_height
elif canvas_height is None:
self.canvas_width = canvas_width
self.canvas_height = self.canvas_width / ratio
else:
raise RuntimeError("Can't set *both* width and height!")

# Distance to image plane from camera
if image_plane is None:
if self.camera_type == 'orthographic':
self.image_plane = 1 - self.camera_dist
else:
self.image_plane = 0
self.image_plane += self.camera_dist

self.constrainatoms = []
for c in constraints:
if isinstance(c, FixAtoms):
# self.constrainatoms.extend(c.index) # is this list-like?
for n, i in enumerate(c.index):
self.constrainatoms += [i]

@classmethod
def from_PlottingVariables(cls, pvars, **kwargs):
cell = pvars.cell
cell_vertices = pvars.cell_vertices
if 'colors' in kwargs.keys():
colors = kwargs.pop('colors')
else:
colors = pvars.colors
diameters = pvars.d
image_height = pvars.h
image_width = pvars.w
positions = pvars.positions
constraints = pvars.constraints
return cls(cell=cell, cell_vertices=cell_vertices, colors=colors,
constraints=constraints, diameters=diameters,
image_height=image_height, image_width=image_width,
positions=positions, **kwargs)

@classmethod
def from_atoms(cls, atoms, **kwargs):
return cls.from_plotting_variables(
PlottingVariables(atoms, scale=1.0), **kwargs)

def write_ini(self, path):
"""Write ini file."""

ini_str = f"""\
Input_File_Name={path.with_suffix('.pov').name}
Output_to_File=True
Output_File_Type=N
Output_Alpha={'on' if self.transparent else 'off'}
; if you adjust Height, and width, you must preserve the ratio
; Width / Height = {self.canvas_width/self.canvas_height:f}
Width={self.canvas_width}
Height={self.canvas_height}
Antialias=True
Antialias_Threshold=0.1
Display={self.display}
Display_Gamma=2.2
Pause_When_Done={self.pause}
Verbose=False
"""
with open(path, 'w') as fd:
fd.write(ini_str)
return path

def write_pov(self, path):
"""Write pov file."""

point_lights = '\n'.join(f"light_source {{{pa(loc)} {pc(rgb)}}}"
for loc, rgb in self.point_lights)

area_light = ''
if self.area_light is not None:
loc, color, width, height, nx, ny = self.area_light
area_light += f"""\nlight_source {{{pa(loc)} {pc(color)}
area_light <{width:.2f}, 0, 0>, <0, {height:.2f}, 0>, {nx:n}, {ny:n}

fog = ''
if self.depth_cueing and (self.cue_density >= 1e-4):
# same way vmd does it
if self.cue_density > 1e4:
# larger does not make any sense
dist = 1e-4
else:
dist = 1. / self.cue_density
fog += f'fog {{fog_type 1 distance {dist:.4f} '\
f'color {pc(self.background)}}}'

mat_style_keys = (f'#declare {k} = {v}'
for k, v in self.material_styles_dict.items())
mat_style_keys = '\n'.join(mat_style_keys)

# Draw unit cell
cell_vertices = ''
if self.cell_vertices is not None:
for c in range(3):
for j in ([0, 0], [1, 0], [1, 1], [0, 1]):
p1 = self.cell_vertices[tuple(j[:c]) + (0,) + tuple(j[c:])]
p2 = self.cell_vertices[tuple(j[:c]) + (1,) + tuple(j[c:])]

distance = np.linalg.norm(p2 - p1)
if distance < 1e-12:
continue

cell_vertices += f'cylinder {{{pa(p1)}, {pa(p2)}, '\
f'Rcell pigment {{Black}}}}\n'
# all strings are f-strings for consistency
cell_vertices = cell_vertices.strip('\n')

# Draw atoms
a = 0
atoms = ''
for loc, dia, col in zip(self.positions, self.diameters, self.colors):
tex = 'ase3'
trans = 0.
if self.textures is not None:
tex = self.textures[a]
if self.transmittances is not None:
trans = self.transmittances[a]
atoms += f'atom({pa(loc)}, {dia/2.:.2f}, {pc(col)}, '\
f'{trans}, {tex}) // #{a:n}\n'
a += 1
atoms = atoms.strip('\n')

# Draw atom bonds
bondatoms = ''
for pair in self.bondatoms:
# Make sure that each pair has 4 componets: a, b, offset,
#                                           bond_order, bond_offset
# a, b: atom index to draw bond
# offset: original meaning to make offset for mid-point.
# bond_oder: if not supplied, set it to 1 (single bond).
#            It can be  1, 2, 3, corresponding to single,
#            double, triple bond
# bond_offset: displacement from original bond position.
#              Default is (bondlinewidth, bondlinewidth, 0)
#              for bond_order > 1.
if len(pair) == 2:
a, b = pair
offset = (0, 0, 0)
bond_order = 1
bond_offset = (0, 0, 0)
elif len(pair) == 3:
a, b, offset = pair
bond_order = 1
bond_offset = (0, 0, 0)
elif len(pair) == 4:
a, b, offset, bond_order = pair
bond_offset = (self.bondlinewidth, self.bondlinewidth, 0)
elif len(pair) > 4:
a, b, offset, bond_order, bond_offset = pair
else:
raise RuntimeError('Each list in bondatom must have at least '
'2 entries. Error at %s' % pair)

if len(offset) != 3:
raise ValueError('offset must have 3 elements. '
'Error at %s' % pair)
if len(bond_offset) != 3:
raise ValueError('bond_offset must have 3 elements. '
'Error at %s' % pair)
if bond_order not in [0, 1, 2, 3]:
raise ValueError('bond_order must be either 0, 1, 2, or 3. '
'Error at %s' % pair)

# Up to here, we should have all a, b, offset, bond_order,
# bond_offset for all bonds.

# Rotate bond_offset so that its direction is 90 deg. off the bond
# Utilize Atoms object to rotate
if bond_order > 1 and np.linalg.norm(bond_offset) > 1.e-9:
tmp_atoms = Atoms('H3')
tmp_atoms.set_cell(self.cell)
tmp_atoms.set_positions([
self.positions[a],
self.positions[b],
self.positions[b] + np.array(bond_offset),
])
tmp_atoms.center()
tmp_atoms.set_angle(0, 1, 2, 90)
bond_offset = tmp_atoms.position - tmp_atoms.position

R = np.dot(offset, self.cell)
mida = 0.5 * (self.positions[a] + self.positions[b] + R)
midb = 0.5 * (self.positions[a] + self.positions[b] - R)
if self.textures is not None:
texa = self.textures[a]
texb = self.textures[b]
else:
texa = texb = 'ase3'

if self.transmittances is not None:
transa = self.transmittances[a]
transb = self.transmittances[b]
else:
transa = transb = 0.

# draw bond, according to its bond_order.
# bond_order == 0: No bond is plotted
# bond_order == 1: use original code
# bond_order == 2: draw two bonds, one is shifted by bond_offset/2,
#                  and another is shifted by -bond_offset/2.
# bond_order == 3: draw two bonds, one is shifted by bond_offset,
#                  and one is shifted by -bond_offset, and the
#                  other has no shift.
# To shift the bond, add the shift to the first two coordinate in
# write statement.

posa = self.positions[a]
posb = self.positions[b]
cola = self.colors[a]
colb = self.colors[b]

if bond_order == 1:
draw_tuples = (
(posa, mida, cola, transa, texa),
(posb, midb, colb, transb, texb))

elif bond_order == 2:
bs = [x / 2 for x in bond_offset]
draw_tuples = (
(posa - bs, mida - bs, cola, transa, texa),
(posb - bs, midb - bs, colb, transb, texb),
(posa + bs, mida + bs, cola, transa, texa),
(posb + bs, midb + bs, colb, transb, texb))

elif bond_order == 3:
bs = bond_offset
draw_tuples = (
(posa, mida, cola, transa, texa),
(posb, midb, colb, transb, texb),
(posa + bs, mida + bs, cola, transa, texa),
(posb + bs, midb + bs, colb, transb, texb),
(posa - bs, mida - bs, cola, transa, texa),
(posb - bs, midb - bs, colb, transb, texb))

bondatoms += ''.join(f'cylinder {{{pa(p)}, '
f'{pa(m)}, Rbond texture{{pigment '
f'{{color {pc(c)} '
f'transmit {tr}}} finish{{{tx}}}}}}}\n'
for p, m, c, tr, tx in
draw_tuples)

bondatoms = bondatoms.strip('\n')

# Draw constraints if requested
constraints = ''
if self.exportconstraints:
for a in self.constrainatoms:
dia = self.diameters[a]
loc = self.positions[a]
trans = 0.0
if self.transmittances is not None:
trans = self.transmittances[a]
constraints += f'constrain({pa(loc)}, {dia/2.:.2f}, Black, '\
f'{trans}, {tex}) // #{a:n} \n'
constraints = constraints.strip('\n')

pov = f"""#version 3.6;
#include "colors.inc"
#include "finish.inc"

global_settings {{assumed_gamma 2.2 max_trace_level 6}}
background {{{pc(self.background)}{' transmit 1.0' if self.transparent else ''}}}
camera {{{self.camera_type}
right -{self.image_width:.2f}*x up {self.image_height:.2f}*y
direction {self.image_plane:.2f}*z
location <0,0,{self.camera_dist:.2f}> look_at <0,0,0>}}
{point_lights}
{area_light if area_light != '' else '// no area light'}
{fog if fog != '' else '// no fog'}
{mat_style_keys}
#declare Rcell = {self.celllinewidth:.3f};
#declare Rbond = {self.bondlinewidth:.3f};

#macro atom(LOC, R, COL, TRANS, FIN)
sphere{{LOC, R texture{{pigment{{color COL transmit TRANS}} finish{{FIN}}}}}}
#end
#macro constrain(LOC, R, COL, TRANS FIN)
union{{torus{{R, Rcell rotate 45*z texture{{pigment{{color COL transmit TRANS}} finish{{FIN}}}}}}
torus{{R, Rcell rotate -45*z texture{{pigment{{color COL transmit TRANS}} finish{{FIN}}}}}}
translate LOC}}
#end

{cell_vertices if cell_vertices != '' else '// no cell vertices'}
{atoms}
{bondatoms}
{constraints if constraints != '' else '// no constraints'}
"""  # noqa: E501

with open(path, 'w') as fd:
fd.write(pov)

return path

def write(self, pov_path):
pov_path = require_pov(pov_path)
ini_path = pov_path.with_suffix('.ini')
self.write_ini(ini_path)
self.write_pov(pov_path)
if self.isosurfaces is not None:
with open(pov_path, 'a') as fd:
for iso in self.isosurfaces:
fd.write(iso.format_mesh())
return POVRAYInputs(ini_path)

def require_pov(path):
path = Path(path)
if path.suffix != '.pov':
raise ValueError(f'Expected .pov path, got {path}')
return path

class POVRAYInputs:
def __init__(self, path):
self.path = path

def render(self, povray_executable='povray', stderr=DEVNULL,
clean_up=False):
cmd = [povray_executable, str(self.path)]

check_call(cmd, stderr=stderr)
png_path = self.path.with_suffix('.png').absolute()
if not png_path.is_file():
raise RuntimeError(f'Povray left no output PNG file "{png_path}"')

if clean_up:

return png_path

class POVRAYIsosurface:
def __init__(self, density_grid, cut_off, cell, cell_origin,
color=(0.85, 0.80, 0.25, 0.2), material='ase3'):
"""
density_grid: 3D float ndarray
A regular grid on that spans the cell. The first dimension
corresponds to the first cell vector and so on.
cut_off: float
The density value of the isosurface.
cell: 2D float ndarray or ASE cell object
The 3 vectors which give the cell's repetition
cell_origin: 4 float tuple
The cell origin as used by POVRAY object
closed_edges: bool
Setting this will fill in isosurface edges at the cell boundaries.
Filling in the edges can help with visualizing
highly porous structures.
Lets you pick the area you want to enclose, i.e., should the denser
or less dense area be filled in.
color: povray color string, float, or float tuple
1 float is interpreted as grey scale, a 3 float tuple is rgb,
4 float tuple is rgbt, and 5 float tuple is rgbft, where
t is transmission fraction and f is filter fraction.
Named Povray colors are set in colors.inc
(http://wiki.povray.org/content/Reference:Colors.inc)
material: string
Can be a finish macro defined by POVRAY.material_styles
or a full Povray material {...} specification. Using a
full material specification willoverride the color parameter.
"""

self.color = color
self.material = material
self.closed_edges = closed_edges
self._cut_off = cut_off

cv = 2 * cut_off
else:
cv = 0

if closed_edges:
shape_old = density_grid.shape
# since well be padding, we need to keep the data at origin
cell_origin += -(1.0 / np.array(shape_old)) @ cell
1,), mode='constant', constant_values=cv)
shape_new = density_grid.shape
s = np.array(shape_new) / np.array(shape_old)
cell = cell @ np.diag(s)

self.cell = cell
self.cell_origin = cell_origin
self.density_grid = density_grid
self.spacing = tuple(1.0 / np.array(self.density_grid.shape))

scaled_verts, faces, normals, values = self.compute_mesh(
self.density_grid,
self.cut_off,
self.spacing,

# The verts are scaled by default, this is the super easy way of
# distributing them in real space but it's easier to do affine
# transformations/rotations on a unit cube so I leave it like that
# verts = scaled_verts.dot(self.cell)
self.verts = scaled_verts
self.faces = faces

@property
def cut_off(self):
return self._cut_off

@cut_off.setter
def cut_off(self, value):
raise Exception("Use the set_cut_off method")

def set_cut_off(self, value):
self._cut_off = value

cv = 2 * self.cut_off
else:
cv = 0

if self.closed_edges:
shape_old = self.density_grid.shape
# since well be padding, we need to keep the data at origin
self.cell_origin += -(1.0 / np.array(shape_old)) @ self.cell
1,), mode='constant', constant_values=cv)
shape_new = self.density_grid.shape
s = np.array(shape_new) / np.array(shape_old)
self.cell = self.cell @ np.diag(s)

self.spacing = tuple(1.0 / np.array(self.density_grid.shape))

scaled_verts, faces, _, _ = self.compute_mesh(
self.density_grid,
self.cut_off,
self.spacing,

self.verts = scaled_verts
self.faces = faces

@classmethod
def from_POVRAY(cls, povray, density_grid, cut_off, **kwargs):
return cls(cell=povray.cell,
cell_origin=povray.cell_vertices[0, 0, 0],
density_grid=density_grid,
cut_off=cut_off, **kwargs)

@staticmethod
def wrapped_triples_section(triple_list,
triple_format="<{:f}, {:f}, {:f}>".format,
triples_per_line=4):

triples = [triple_format(*x) for x in triple_list]
n = len(triples)
s = ''
tpl = triples_per_line
c = 0

while c < n - tpl:
c += tpl
s += '\n     '
s += ', '.join(triples[c - tpl:c])
s += '\n    '
s += ', '.join(triples[c:])
return s

@staticmethod
"""

Import statement is in this method and not file header
since few users will use isosurface rendering.

Returns scaled_verts, faces, normals, values. See skimage docs.

"""

# marching_cubes name was changed in skimage v0.19
try:
# New skimage
from skimage.measure import marching_cubes
except ImportError:
# Old skimage (remove at some point)
from skimage.measure import (
marching_cubes_lewiner as marching_cubes)

return marching_cubes(
density_grid,
level=cut_off,
spacing=spacing,
allow_degenerate=False)

def format_mesh(self):
"""Returns a formatted data output for POVRAY files

Example:
material = '''
material { // This material looks like pink jelly
texture {
pigment { rgbt <0.8, 0.25, 0.25, 0.5> }
finish{ diffuse 0.85 ambient 0.99 brilliance 3 specular 0.5 roughness 0.001
reflection { 0.05, 0.98 fresnel on exponent 1.5 }
conserve_energy
}
}
interior { ior 1.3 }
}
photons {
target
refraction on
reflection on
collect on
}'''
"""  # noqa: E501

if self.material in POVRAY.material_styles_dict:
material = f"""material {{
texture {{
pigment {{ {pc(self.color)} }}
finish {{ {self.material} }}
}}
}}"""
else:
material = self.material

# Start writing the mesh2
vertex_vectors = self.wrapped_triples_section(
triple_list=self.verts,
triple_format="<{:f}, {:f}, {:f}>".format,
triples_per_line=4)

face_indices = self.wrapped_triples_section(
triple_list=self.faces,
triple_format="<{:n}, {:n}, {:n}>".format,
triples_per_line=5)

cell = self.cell
cell_or = self.cell_origin
mesh2 = f"""\n\nmesh2 {{
vertex_vectors {{  {len(self.verts):n},
{vertex_vectors}
}}
face_indices {{ {len(self.faces):n},
{face_indices}
}}
{material if material != '' else '// no material'}
matrix < {cell:f}, {cell:f}, {cell:f},
{cell:f}, {cell:f}, {cell:f},
{cell:f}, {cell:f}, {cell:f},
{cell_or:f}, {cell_or:f}, {cell_or:f}>
}}
"""
return mesh2

def pop_deprecated(dct, name):
import warnings
if name in dct:
del dct[name]
warnings.warn(f'The "{name}" keyword of write_pov() is deprecated '
'and has no effect; this will raise an error in the '
'future.', FutureWarning)

[docs]def write_pov(filename, atoms, *,
povray_settings=None, isosurface_data=None,
**generic_projection_settings):

for name in ['run_povray', 'povray_path', 'stderr', 'extras']:
pop_deprecated(generic_projection_settings, name)

if povray_settings is None:
povray_settings = {}

pvars = PlottingVariables(atoms, scale=1.0, **generic_projection_settings)
pov_obj = POVRAY.from_PlottingVariables(pvars, **povray_settings)

if isosurface_data is None:
isosurface_data = []
elif not isinstance(isosurface_data, Sequence):
isosurface_data = [isosurface_data]

isosurfaces = []
for isodata in isosurface_data:
if isinstance(isodata, POVRAYIsosurface):
iso = isodata
else:
iso = POVRAYIsosurface.from_POVRAY(pov_obj, **isodata)
isosurfaces.append(iso)
pov_obj.isosurfaces = isosurfaces

return pov_obj.write(filename)
```