Cleanup: simplify file name incrementing logic
[blender-addons.git] / io_mesh_atomic / pdb_import.py
blobdc687c602cec8690f8796c9e5279901cd9a08c18
1 # ##### BEGIN GPL LICENSE BLOCK #####
3 # This program is free software; you can redistribute it and/or
4 # modify it under the terms of the GNU General Public License
5 # as published by the Free Software Foundation; either version 2
6 # of the License, or (at your option) any later version.
8 # This program is distributed in the hope that it will be useful,
9 # but WITHOUT ANY WARRANTY; without even the implied warranty of
10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 # GNU General Public License for more details.
13 # You should have received a copy of the GNU General Public License
14 # along with this program; if not, write to the Free Software Foundation,
15 # Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17 # ##### END GPL LICENSE BLOCK #####
19 import os
20 import bpy
21 import bmesh
22 from math import pi, cos, sin, sqrt, ceil
23 from mathutils import Vector, Matrix
24 from copy import copy
26 # -----------------------------------------------------------------------------
27 # Atom, stick and element data
30 # This is a list that contains some data of all possible elements. The structure
31 # is as follows:
33 # 1, "Hydrogen", "H", [0.0,0.0,1.0], 0.32, 0.32, 0.32 , -1 , 1.54 means
35 # No., name, short name, color, radius (used), radius (covalent), radius (atomic),
37 # charge state 1, radius (ionic) 1, charge state 2, radius (ionic) 2, ... all
38 # charge states for any atom are listed, if existing.
39 # The list is fixed and cannot be changed ... (see below)
41 ELEMENTS_DEFAULT = (
42 ( 1, "Hydrogen", "H", ( 1.0, 1.0, 1.0, 1.0), 0.32, 0.32, 0.79 , -1 , 1.54 ),
43 ( 2, "Helium", "He", ( 0.85, 1.0, 1.0, 1.0), 0.93, 0.93, 0.49 ),
44 ( 3, "Lithium", "Li", ( 0.8, 0.50, 1.0, 1.0), 1.23, 1.23, 2.05 , 1 , 0.68 ),
45 ( 4, "Beryllium", "Be", ( 0.76, 1.0, 0.0, 1.0), 0.90, 0.90, 1.40 , 1 , 0.44 , 2 , 0.35 ),
46 ( 5, "Boron", "B", ( 1.0, 0.70, 0.70, 1.0), 0.82, 0.82, 1.17 , 1 , 0.35 , 3 , 0.23 ),
47 ( 6, "Carbon", "C", ( 0.56, 0.56, 0.56, 1.0), 0.77, 0.77, 0.91 , -4 , 2.60 , 4 , 0.16 ),
48 ( 7, "Nitrogen", "N", ( 0.18, 0.31, 0.97, 1.0), 0.75, 0.75, 0.75 , -3 , 1.71 , 1 , 0.25 , 3 , 0.16 , 5 , 0.13 ),
49 ( 8, "Oxygen", "O", ( 1.0, 0.05, 0.05, 1.0), 0.73, 0.73, 0.65 , -2 , 1.32 , -1 , 1.76 , 1 , 0.22 , 6 , 0.09 ),
50 ( 9, "Fluorine", "F", ( 0.56, 0.87, 0.31, 1.0), 0.72, 0.72, 0.57 , -1 , 1.33 , 7 , 0.08 ),
51 (10, "Neon", "Ne", ( 0.70, 0.89, 0.96, 1.0), 0.71, 0.71, 0.51 , 1 , 1.12 ),
52 (11, "Sodium", "Na", ( 0.67, 0.36, 0.94, 1.0), 1.54, 1.54, 2.23 , 1 , 0.97 ),
53 (12, "Magnesium", "Mg", ( 0.54, 1.0, 0.0, 1.0), 1.36, 1.36, 1.72 , 1 , 0.82 , 2 , 0.66 ),
54 (13, "Aluminium", "Al", ( 0.74, 0.65, 0.65, 1.0), 1.18, 1.18, 1.82 , 3 , 0.51 ),
55 (14, "Silicon", "Si", ( 0.94, 0.78, 0.62, 1.0), 1.11, 1.11, 1.46 , -4 , 2.71 , -1 , 3.84 , 1 , 0.65 , 4 , 0.42 ),
56 (15, "Phosphorus", "P", ( 1.0, 0.50, 0.0, 1.0), 1.06, 1.06, 1.23 , -3 , 2.12 , 3 , 0.44 , 5 , 0.35 ),
57 (16, "Sulfur", "S", ( 1.0, 1.0, 0.18, 1.0), 1.02, 1.02, 1.09 , -2 , 1.84 , 2 , 2.19 , 4 , 0.37 , 6 , 0.30 ),
58 (17, "Chlorine", "Cl", ( 0.12, 0.94, 0.12, 1.0), 0.99, 0.99, 0.97 , -1 , 1.81 , 5 , 0.34 , 7 , 0.27 ),
59 (18, "Argon", "Ar", ( 0.50, 0.81, 0.89, 1.0), 0.98, 0.98, 0.88 , 1 , 1.54 ),
60 (19, "Potassium", "K", ( 0.56, 0.25, 0.83, 1.0), 2.03, 2.03, 2.77 , 1 , 0.81 ),
61 (20, "Calcium", "Ca", ( 0.23, 1.0, 0.0, 1.0), 1.74, 1.74, 2.23 , 1 , 1.18 , 2 , 0.99 ),
62 (21, "Scandium", "Sc", ( 0.90, 0.90, 0.90, 1.0), 1.44, 1.44, 2.09 , 3 , 0.73 ),
63 (22, "Titanium", "Ti", ( 0.74, 0.76, 0.78, 1.0), 1.32, 1.32, 2.00 , 1 , 0.96 , 2 , 0.94 , 3 , 0.76 , 4 , 0.68 ),
64 (23, "Vanadium", "V", ( 0.65, 0.65, 0.67, 1.0), 1.22, 1.22, 1.92 , 2 , 0.88 , 3 , 0.74 , 4 , 0.63 , 5 , 0.59 ),
65 (24, "Chromium", "Cr", ( 0.54, 0.6, 0.78, 1.0), 1.18, 1.18, 1.85 , 1 , 0.81 , 2 , 0.89 , 3 , 0.63 , 6 , 0.52 ),
66 (25, "Manganese", "Mn", ( 0.61, 0.47, 0.78, 1.0), 1.17, 1.17, 1.79 , 2 , 0.80 , 3 , 0.66 , 4 , 0.60 , 7 , 0.46 ),
67 (26, "Iron", "Fe", ( 0.87, 0.4, 0.2, 1.0), 1.17, 1.17, 1.72 , 2 , 0.74 , 3 , 0.64 ),
68 (27, "Cobalt", "Co", ( 0.94, 0.56, 0.62, 1.0), 1.16, 1.16, 1.67 , 2 , 0.72 , 3 , 0.63 ),
69 (28, "Nickel", "Ni", ( 0.31, 0.81, 0.31, 1.0), 1.15, 1.15, 1.62 , 2 , 0.69 ),
70 (29, "Copper", "Cu", ( 0.78, 0.50, 0.2, 1.0), 1.17, 1.17, 1.57 , 1 , 0.96 , 2 , 0.72 ),
71 (30, "Zinc", "Zn", ( 0.49, 0.50, 0.69, 1.0), 1.25, 1.25, 1.53 , 1 , 0.88 , 2 , 0.74 ),
72 (31, "Gallium", "Ga", ( 0.76, 0.56, 0.56, 1.0), 1.26, 1.26, 1.81 , 1 , 0.81 , 3 , 0.62 ),
73 (32, "Germanium", "Ge", ( 0.4, 0.56, 0.56, 1.0), 1.22, 1.22, 1.52 , -4 , 2.72 , 2 , 0.73 , 4 , 0.53 ),
74 (33, "Arsenic", "As", ( 0.74, 0.50, 0.89, 1.0), 1.20, 1.20, 1.33 , -3 , 2.22 , 3 , 0.58 , 5 , 0.46 ),
75 (34, "Selenium", "Se", ( 1.0, 0.63, 0.0, 1.0), 1.16, 1.16, 1.22 , -2 , 1.91 , -1 , 2.32 , 1 , 0.66 , 4 , 0.50 , 6 , 0.42 ),
76 (35, "Bromine", "Br", ( 0.65, 0.16, 0.16, 1.0), 1.14, 1.14, 1.12 , -1 , 1.96 , 5 , 0.47 , 7 , 0.39 ),
77 (36, "Krypton", "Kr", ( 0.36, 0.72, 0.81, 1.0), 1.31, 1.31, 1.24 ),
78 (37, "Rubidium", "Rb", ( 0.43, 0.18, 0.69, 1.0), 2.16, 2.16, 2.98 , 1 , 1.47 ),
79 (38, "Strontium", "Sr", ( 0.0, 1.0, 0.0, 1.0), 1.91, 1.91, 2.45 , 2 , 1.12 ),
80 (39, "Yttrium", "Y", ( 0.58, 1.0, 1.0, 1.0), 1.62, 1.62, 2.27 , 3 , 0.89 ),
81 (40, "Zirconium", "Zr", ( 0.58, 0.87, 0.87, 1.0), 1.45, 1.45, 2.16 , 1 , 1.09 , 4 , 0.79 ),
82 (41, "Niobium", "Nb", ( 0.45, 0.76, 0.78, 1.0), 1.34, 1.34, 2.08 , 1 , 1.00 , 4 , 0.74 , 5 , 0.69 ),
83 (42, "Molybdenum", "Mo", ( 0.32, 0.70, 0.70, 1.0), 1.30, 1.30, 2.01 , 1 , 0.93 , 4 , 0.70 , 6 , 0.62 ),
84 (43, "Technetium", "Tc", ( 0.23, 0.61, 0.61, 1.0), 1.27, 1.27, 1.95 , 7 , 0.97 ),
85 (44, "Ruthenium", "Ru", ( 0.14, 0.56, 0.56, 1.0), 1.25, 1.25, 1.89 , 4 , 0.67 ),
86 (45, "Rhodium", "Rh", ( 0.03, 0.49, 0.54, 1.0), 1.25, 1.25, 1.83 , 3 , 0.68 ),
87 (46, "Palladium", "Pd", ( 0.0, 0.41, 0.52, 1.0), 1.28, 1.28, 1.79 , 2 , 0.80 , 4 , 0.65 ),
88 (47, "Silver", "Ag", ( 0.75, 0.75, 0.75, 1.0), 1.34, 1.34, 1.75 , 1 , 1.26 , 2 , 0.89 ),
89 (48, "Cadmium", "Cd", ( 1.0, 0.85, 0.56, 1.0), 1.48, 1.48, 1.71 , 1 , 1.14 , 2 , 0.97 ),
90 (49, "Indium", "In", ( 0.65, 0.45, 0.45, 1.0), 1.44, 1.44, 2.00 , 3 , 0.81 ),
91 (50, "Tin", "Sn", ( 0.4, 0.50, 0.50, 1.0), 1.41, 1.41, 1.72 , -4 , 2.94 , -1 , 3.70 , 2 , 0.93 , 4 , 0.71 ),
92 (51, "Antimony", "Sb", ( 0.61, 0.38, 0.70, 1.0), 1.40, 1.40, 1.53 , -3 , 2.45 , 3 , 0.76 , 5 , 0.62 ),
93 (52, "Tellurium", "Te", ( 0.83, 0.47, 0.0, 1.0), 1.36, 1.36, 1.42 , -2 , 2.11 , -1 , 2.50 , 1 , 0.82 , 4 , 0.70 , 6 , 0.56 ),
94 (53, "Iodine", "I", ( 0.58, 0.0, 0.58, 1.0), 1.33, 1.33, 1.32 , -1 , 2.20 , 5 , 0.62 , 7 , 0.50 ),
95 (54, "Xenon", "Xe", ( 0.25, 0.61, 0.69, 1.0), 1.31, 1.31, 1.24 ),
96 (55, "Caesium", "Cs", ( 0.34, 0.09, 0.56, 1.0), 2.35, 2.35, 3.35 , 1 , 1.67 ),
97 (56, "Barium", "Ba", ( 0.0, 0.78, 0.0, 1.0), 1.98, 1.98, 2.78 , 1 , 1.53 , 2 , 1.34 ),
98 (57, "Lanthanum", "La", ( 0.43, 0.83, 1.0, 1.0), 1.69, 1.69, 2.74 , 1 , 1.39 , 3 , 1.06 ),
99 (58, "Cerium", "Ce", ( 1.0, 1.0, 0.78, 1.0), 1.65, 1.65, 2.70 , 1 , 1.27 , 3 , 1.03 , 4 , 0.92 ),
100 (59, "Praseodymium", "Pr", ( 0.85, 1.0, 0.78, 1.0), 1.65, 1.65, 2.67 , 3 , 1.01 , 4 , 0.90 ),
101 (60, "Neodymium", "Nd", ( 0.78, 1.0, 0.78, 1.0), 1.64, 1.64, 2.64 , 3 , 0.99 ),
102 (61, "Promethium", "Pm", ( 0.63, 1.0, 0.78, 1.0), 1.63, 1.63, 2.62 , 3 , 0.97 ),
103 (62, "Samarium", "Sm", ( 0.56, 1.0, 0.78, 1.0), 1.62, 1.62, 2.59 , 3 , 0.96 ),
104 (63, "Europium", "Eu", ( 0.38, 1.0, 0.78, 1.0), 1.85, 1.85, 2.56 , 2 , 1.09 , 3 , 0.95 ),
105 (64, "Gadolinium", "Gd", ( 0.27, 1.0, 0.78, 1.0), 1.61, 1.61, 2.54 , 3 , 0.93 ),
106 (65, "Terbium", "Tb", ( 0.18, 1.0, 0.78, 1.0), 1.59, 1.59, 2.51 , 3 , 0.92 , 4 , 0.84 ),
107 (66, "Dysprosium", "Dy", ( 0.12, 1.0, 0.78, 1.0), 1.59, 1.59, 2.49 , 3 , 0.90 ),
108 (67, "Holmium", "Ho", ( 0.0, 1.0, 0.61, 1.0), 1.58, 1.58, 2.47 , 3 , 0.89 ),
109 (68, "Erbium", "Er", ( 0.0, 0.90, 0.45, 1.0), 1.57, 1.57, 2.45 , 3 , 0.88 ),
110 (69, "Thulium", "Tm", ( 0.0, 0.83, 0.32, 1.0), 1.56, 1.56, 2.42 , 3 , 0.87 ),
111 (70, "Ytterbium", "Yb", ( 0.0, 0.74, 0.21, 1.0), 1.74, 1.74, 2.40 , 2 , 0.93 , 3 , 0.85 ),
112 (71, "Lutetium", "Lu", ( 0.0, 0.67, 0.14, 1.0), 1.56, 1.56, 2.25 , 3 , 0.85 ),
113 (72, "Hafnium", "Hf", ( 0.30, 0.76, 1.0, 1.0), 1.44, 1.44, 2.16 , 4 , 0.78 ),
114 (73, "Tantalum", "Ta", ( 0.30, 0.65, 1.0, 1.0), 1.34, 1.34, 2.09 , 5 , 0.68 ),
115 (74, "Tungsten", "W", ( 0.12, 0.58, 0.83, 1.0), 1.30, 1.30, 2.02 , 4 , 0.70 , 6 , 0.62 ),
116 (75, "Rhenium", "Re", ( 0.14, 0.49, 0.67, 1.0), 1.28, 1.28, 1.97 , 4 , 0.72 , 7 , 0.56 ),
117 (76, "Osmium", "Os", ( 0.14, 0.4, 0.58, 1.0), 1.26, 1.26, 1.92 , 4 , 0.88 , 6 , 0.69 ),
118 (77, "Iridium", "Ir", ( 0.09, 0.32, 0.52, 1.0), 1.27, 1.27, 1.87 , 4 , 0.68 ),
119 (78, "Platinum", "Pt", ( 0.81, 0.81, 0.87, 1.0), 1.30, 1.30, 1.83 , 2 , 0.80 , 4 , 0.65 ),
120 (79, "Gold", "Au", ( 1.0, 0.81, 0.13, 1.0), 1.34, 1.34, 1.79 , 1 , 1.37 , 3 , 0.85 ),
121 (80, "Mercury", "Hg", ( 0.72, 0.72, 0.81, 1.0), 1.49, 1.49, 1.76 , 1 , 1.27 , 2 , 1.10 ),
122 (81, "Thallium", "Tl", ( 0.65, 0.32, 0.30, 1.0), 1.48, 1.48, 2.08 , 1 , 1.47 , 3 , 0.95 ),
123 (82, "Lead", "Pb", ( 0.34, 0.34, 0.38, 1.0), 1.47, 1.47, 1.81 , 2 , 1.20 , 4 , 0.84 ),
124 (83, "Bismuth", "Bi", ( 0.61, 0.30, 0.70, 1.0), 1.46, 1.46, 1.63 , 1 , 0.98 , 3 , 0.96 , 5 , 0.74 ),
125 (84, "Polonium", "Po", ( 0.67, 0.36, 0.0, 1.0), 1.46, 1.46, 1.53 , 6 , 0.67 ),
126 (85, "Astatine", "At", ( 0.45, 0.30, 0.27, 1.0), 1.45, 1.45, 1.43 , -3 , 2.22 , 3 , 0.85 , 5 , 0.46 ),
127 (86, "Radon", "Rn", ( 0.25, 0.50, 0.58, 1.0), 1.00, 1.00, 1.34 ),
128 (87, "Francium", "Fr", ( 0.25, 0.0, 0.4, 1.0), 1.00, 1.00, 1.00 , 1 , 1.80 ),
129 (88, "Radium", "Ra", ( 0.0, 0.49, 0.0, 1.0), 1.00, 1.00, 1.00 , 2 , 1.43 ),
130 (89, "Actinium", "Ac", ( 0.43, 0.67, 0.98, 1.0), 1.00, 1.00, 1.00 , 3 , 1.18 ),
131 (90, "Thorium", "Th", ( 0.0, 0.72, 1.0, 1.0), 1.65, 1.65, 1.00 , 4 , 1.02 ),
132 (91, "Protactinium", "Pa", ( 0.0, 0.63, 1.0, 1.0), 1.00, 1.00, 1.00 , 3 , 1.13 , 4 , 0.98 , 5 , 0.89 ),
133 (92, "Uranium", "U", ( 0.0, 0.56, 1.0, 1.0), 1.42, 1.42, 1.00 , 4 , 0.97 , 6 , 0.80 ),
134 (93, "Neptunium", "Np", ( 0.0, 0.50, 1.0, 1.0), 1.00, 1.00, 1.00 , 3 , 1.10 , 4 , 0.95 , 7 , 0.71 ),
135 (94, "Plutonium", "Pu", ( 0.0, 0.41, 1.0, 1.0), 1.00, 1.00, 1.00 , 3 , 1.08 , 4 , 0.93 ),
136 (95, "Americium", "Am", ( 0.32, 0.36, 0.94, 1.0), 1.00, 1.00, 1.00 , 3 , 1.07 , 4 , 0.92 ),
137 (96, "Curium", "Cm", ( 0.47, 0.36, 0.89, 1.0), 1.00, 1.00, 1.00 ),
138 (97, "Berkelium", "Bk", ( 0.54, 0.30, 0.89, 1.0), 1.00, 1.00, 1.00 ),
139 (98, "Californium", "Cf", ( 0.63, 0.21, 0.83, 1.0), 1.00, 1.00, 1.00 ),
140 (99, "Einsteinium", "Es", ( 0.70, 0.12, 0.83, 1.0), 1.00, 1.00, 1.00 ),
141 (100, "Fermium", "Fm", ( 0.70, 0.12, 0.72, 1.0), 1.00, 1.00, 1.00 ),
142 (101, "Mendelevium", "Md", ( 0.70, 0.05, 0.65, 1.0), 1.00, 1.00, 1.00 ),
143 (102, "Nobelium", "No", ( 0.74, 0.05, 0.52, 1.0), 1.00, 1.00, 1.00 ),
144 (103, "Lawrencium", "Lr", ( 0.78, 0.0, 0.4, 1.0), 1.00, 1.00, 1.00 ),
145 (104, "Vacancy", "Vac", ( 0.5, 0.5, 0.5, 1.0), 1.00, 1.00, 1.00),
146 (105, "Default", "Default", ( 1.0, 1.0, 1.0, 1.0), 1.00, 1.00, 1.00),
147 (106, "Stick", "Stick", ( 0.5, 0.5, 0.5, 1.0), 1.00, 1.00, 1.00),
150 # This list here contains all data of the elements and will be used during
151 # runtime. It is a list of classes.
152 # During executing Atomic Blender, the list will be initialized with the fixed
153 # data from above via the class structure below (ElementProp). We
154 # have then one fixed list (above), which will never be changed, and a list of
155 # classes with same data. The latter can be modified via loading a separate
156 # custom data file.
157 ELEMENTS = []
159 # This is the class, which stores the properties for one element.
160 class ElementProp(object):
161 __slots__ = ('number', 'name', 'short_name', 'color', 'radii', 'radii_ionic')
162 def __init__(self, number, name, short_name, color, radii, radii_ionic):
163 self.number = number
164 self.name = name
165 self.short_name = short_name
166 self.color = color
167 self.radii = radii
168 self.radii_ionic = radii_ionic
170 # This is the class, which stores the properties of one atom.
171 class AtomProp(object):
172 __slots__ = ('element', 'name', 'location', 'radius', 'color', 'material')
173 def __init__(self, element, name, location, radius, color, material):
174 self.element = element
175 self.name = name
176 self.location = location
177 self.radius = radius
178 self.color = color
179 self.material = material
181 # This is the class, which stores the two atoms of one stick.
182 class StickProp(object):
183 __slots__ = ('atom1', 'atom2', 'number', 'dist')
184 def __init__(self, atom1, atom2, number, dist):
185 self.atom1 = atom1
186 self.atom2 = atom2
187 self.number = number
188 self.dist = dist
190 # -----------------------------------------------------------------------------
191 # Some basic routines
194 # The function, which reads all necessary properties of the elements.
195 def read_elements():
197 del ELEMENTS[:]
199 for item in ELEMENTS_DEFAULT:
201 # All three radii into a list
202 radii = [item[4],item[5],item[6]]
203 # The handling of the ionic radii will be done later. So far, it is an
204 # empty list.
205 radii_ionic = []
207 li = ElementProp(item[0],item[1],item[2],item[3],
208 radii,radii_ionic)
209 ELEMENTS.append(li)
212 # The function, which reads the x,y,z positions of all atoms in a PDB
213 # file.
215 # filepath_pdb: path to pdb file
216 # radiustype : '0' default
217 # '1' atomic radii
218 # '2' van der Waals
219 def read_pdb_file(filepath_pdb, radiustype):
221 # The list of all atoms as read from the PDB file.
222 all_atoms = []
224 # Open the pdb file ...
225 filepath_pdb_p = open(filepath_pdb, "r")
227 #Go to the line, in which "ATOM" or "HETATM" appears.
228 for line in filepath_pdb_p:
229 split_list = line.split(' ')
230 if "ATOM" in split_list[0]:
231 break
232 if "HETATM" in split_list[0]:
233 break
235 j = 0
236 # This is in fact an endless 'while loop', ...
237 while j > -1:
239 # ... the loop is broken here (EOF) ...
240 if line == "":
241 break
243 # If there is a "TER" we need to put empty entries into the lists
244 # in order to not destroy the order of atom numbers and same numbers
245 # used for sticks. "TER? What is that?" TER indicates the end of a
246 # list of ATOM/HETATM records for a chain.
247 if "TER" in line:
248 short_name = "TER"
249 name = "TER"
250 radius = 0.0
251 # 2019-03-14, New
252 color = [0,0,0, 0]
253 location = Vector((0,0,0))
254 # Append the TER into the list. Material remains empty so far.
255 all_atoms.append(AtomProp(short_name,
256 name,
257 location,
258 radius,
259 color,[]))
261 # If 'ATOM or 'HETATM' appears in the line then do ...
262 elif "ATOM" in line or "HETATM" in line:
264 # What follows is due to deviations which appear from PDB to
265 # PDB file. It is very special!
267 # PLEASE, DO NOT CHANGE! ............................... from here
268 if line[12:13] == " " or line[12:13].isdigit() == True:
269 short_name = line[13:14]
270 if line[14:15].islower() == True:
271 short_name = short_name + line[14:15]
272 elif line[12:13].isupper() == True:
273 short_name = line[12:13]
274 if line[13:14].isalpha() == True:
275 short_name = short_name + line[13:14]
276 else:
277 print("Atomic Blender: Strange error in PDB file.\n"
278 "Look for element names at positions 13-16 and 78-79.\n")
279 return -1
281 if len(line) >= 78:
283 if line[76:77] == " ":
284 short_name2 = line[76:77]
285 else:
286 short_name2 = line[76:78]
288 if short_name2.isalpha() == True:
289 FOUND = False
290 for element in ELEMENTS:
291 if str.upper(short_name2) == str.upper(element.short_name):
292 FOUND = True
293 break
294 if FOUND == False:
295 short_name = short_name2
297 # ....................................................... to here.
299 # Go through all elements and find the element of the current atom.
300 FLAG_FOUND = False
301 for element in ELEMENTS:
302 if str.upper(short_name) == str.upper(element.short_name):
303 # Give the atom its proper names, color and radius:
304 short_name = str.upper(element.short_name)
305 name = element.name
306 # int(radiustype) => type of radius:
307 # pre-defined (0), atomic (1) or van der Waals (2)
308 radius = float(element.radii[int(radiustype)])
309 color = element.color
310 FLAG_FOUND = True
311 break
313 # Is it a vacancy or an 'unknown atom' ?
314 if FLAG_FOUND == False:
315 # Give this atom also a name. If it is an 'X' then it is a
316 # vacancy. Otherwise ...
317 if "X" in short_name:
318 short_name = "VAC"
319 name = "Vacancy"
320 radius = float(ELEMENTS[-3].radii[int(radiustype)])
321 color = ELEMENTS[-3].color
322 # ... take what is written in the PDB file. These are somewhat
323 # unknown atoms. This should never happen, the element list is
324 # almost complete. However, we do this due to security reasons.
325 else:
326 short_name = str.upper(short_name)
327 name = str.upper(short_name)
328 radius = float(ELEMENTS[-2].radii[int(radiustype)])
329 color = ELEMENTS[-2].color
331 # x,y and z are at fixed positions in the PDB file.
332 x = float(line[30:38].rsplit()[0])
333 y = float(line[38:46].rsplit()[0])
334 z = float(line[46:55].rsplit()[0])
336 location = Vector((x,y,z))
338 j += 1
340 # Append the atom to the list. Material remains empty so far.
341 all_atoms.append(AtomProp(short_name,
342 name,
343 location,
344 radius,
345 color,[]))
347 line = filepath_pdb_p.readline()
348 line = line[:-1]
350 filepath_pdb_p.close()
351 # From above it can be clearly seen that j is now the number of all atoms.
352 Number_of_total_atoms = j
354 return (Number_of_total_atoms, all_atoms)
357 # The function, which reads the sticks in a PDB file.
358 def read_pdb_file_sticks(filepath_pdb, use_sticks_bonds, all_atoms):
360 # The list of all sticks.
361 all_sticks = []
363 # Open the PDB file.
364 filepath_pdb_p = open(filepath_pdb, "r")
366 line = filepath_pdb_p.readline()
367 split_list = line.split(' ')
369 # Go to the first entry
370 if "CONECT" not in split_list[0]:
371 for line in filepath_pdb_p:
372 split_list = line.split(' ')
373 if "CONECT" in split_list[0]:
374 break
376 Number_of_sticks = 0
377 sticks_double = 0
378 j = 0
379 # This is in fact an endless while loop, ...
380 while j > -1:
382 # ... which is broken here (EOF) ...
383 if line == "":
384 break
385 # ... or here, when no 'CONECT' appears anymore.
386 if "CONECT" not in line:
387 break
389 # Note 2019-03-16: in a PDB file the identifier for sticks is called
390 # 'CONECT' and NOT 'CONNECT'! Please leave this as is, otherwise the
391 # sticks are NOT correctly loaded.
393 # The strings of the atom numbers do have a clear position in the file
394 # (From 7 to 12, from 13 to 18 and so on.) and one needs to consider
395 # this. One could also use the split function but then one gets into
396 # trouble if there are lots of atoms: For instance, it may happen that
397 # one has
398 # CONECT 11111 22244444
400 # In Fact it means that atom No. 11111 has a connection with atom
401 # No. 222 but also with atom No. 44444. The split function would give
402 # me only two numbers (11111 and 22244444), which is wrong.
404 # Cut spaces from the right and 'CONECT' at the beginning
405 line = line.rstrip()
406 line = line[6:]
407 # Amount of loops
408 length = len(line)
409 loops = int(length/5)
411 # List of atoms
412 atom_list = []
413 for i in range(loops):
414 number = line[5*i:5*(i+1)].rsplit()
415 if number != []:
416 if number[0].isdigit() == True:
417 atom_number = int(number[0])
418 atom_list.append(atom_number)
420 # The first atom is connected with all the others in the list.
421 atom1 = atom_list[0]
423 # For all the other atoms in the list do:
424 for atom2 in atom_list[1:]:
426 if use_sticks_bonds == True:
427 number = atom_list[1:].count(atom2)
429 if number == 2 or number == 3:
430 basis_list = list(set(atom_list[1:]))
432 if len(basis_list) > 1:
433 basis1 = (all_atoms[atom1-1].location
434 - all_atoms[basis_list[0]-1].location)
435 basis2 = (all_atoms[atom1-1].location
436 - all_atoms[basis_list[1]-1].location)
437 plane_n = basis1.cross(basis2)
439 dist_n = (all_atoms[atom1-1].location
440 - all_atoms[atom2-1].location)
441 dist_n = dist_n.cross(plane_n)
442 dist_n = dist_n / dist_n.length
443 else:
444 dist_n = (all_atoms[atom1-1].location
445 - all_atoms[atom2-1].location)
446 dist_n = Vector((dist_n[1],-dist_n[0],0))
447 dist_n = dist_n / dist_n.length
449 elif number > 3:
450 number = 1
451 dist_n = None
452 else:
453 dist_n = None
454 else:
455 number = 1
456 dist_n = None
458 # Note that in a PDB file, sticks of one atom pair can appear a
459 # couple of times. (Only god knows why ...)
460 # So, does a stick between the considered atoms already exist?
461 FLAG_BAR = False
462 for k in range(Number_of_sticks):
463 if ((all_sticks[k].atom1 == atom1 and all_sticks[k].atom2 == atom2) or
464 (all_sticks[k].atom2 == atom1 and all_sticks[k].atom1 == atom2)):
465 sticks_double += 1
466 # If yes, then FLAG on 'True'.
467 FLAG_BAR = True
468 break
470 # If the stick is not yet registered (FLAG_BAR == False), then
471 # register it!
472 if FLAG_BAR == False:
473 all_sticks.append(StickProp(atom1,atom2,number,dist_n))
474 Number_of_sticks += 1
475 j += 1
477 line = filepath_pdb_p.readline()
478 line = line.rstrip()
480 filepath_pdb_p.close()
482 return all_sticks
485 # Function, which produces a cylinder. All is somewhat easy to understand.
486 def build_stick(radius, length, sectors, element_name):
488 dphi = 2.0 * pi/(float(sectors)-1)
490 # Vertices
491 vertices_top = [Vector((0,0,length / 2.0))]
492 vertices_bottom = [Vector((0,0,-length / 2.0))]
493 vertices = []
494 for i in range(sectors-1):
495 x = radius * cos( dphi * i )
496 y = radius * sin( dphi * i )
497 z = length / 2.0
498 vertex = Vector((x,y,z))
499 vertices_top.append(vertex)
500 z = -length / 2.0
501 vertex = Vector((x,y,z))
502 vertices_bottom.append(vertex)
503 vertices = vertices_top + vertices_bottom
505 # Side facets (Cylinder)
506 faces1 = []
507 for i in range(sectors-1):
508 if i == sectors-2:
509 faces1.append( [i+1, 1, 1+sectors, i+1+sectors] )
510 else:
511 faces1.append( [i+1, i+2, i+2+sectors, i+1+sectors] )
513 # Top facets
514 faces2 = []
515 for i in range(sectors-1):
516 if i == sectors-2:
517 face_top = [0,sectors-1,1]
518 face_bottom = [sectors,2*sectors-1,sectors+1]
519 else:
520 face_top = [0]
521 face_bottom = [sectors]
522 for j in range(2):
523 face_top.append(i+j+1)
524 face_bottom.append(i+j+1+sectors)
525 faces2.append(face_top)
526 faces2.append(face_bottom)
528 # Build the mesh, Cylinder
529 cylinder = bpy.data.meshes.new(element_name+"_sticks_cylinder")
530 cylinder.from_pydata(vertices, [], faces1)
531 cylinder.update()
532 new_cylinder = bpy.data.objects.new(element_name+"_sticks_cylinder", cylinder)
533 # Attention: the linking will be done a few moments later, after this
534 # is done definition.
536 # Build the mesh, Cups
537 cups = bpy.data.meshes.new(element_name+"_sticks_cup")
538 cups.from_pydata(vertices, [], faces2)
539 cups.update()
540 new_cups = bpy.data.objects.new(element_name+"_sticks_cup", cups)
541 # Attention: the linking will be done a few moments later, after this
542 # is done definition.
544 return (new_cylinder, new_cups)
547 # Rotate an object.
548 def rotate_object(rot_mat, obj):
550 bpy.ops.object.select_all(action='DESELECT')
551 obj.select_set(True)
553 # Decompose world_matrix's components, and from them assemble 4x4 matrices.
554 orig_loc, orig_rot, orig_scale = obj.matrix_world.decompose()
556 orig_loc_mat = Matrix.Translation(orig_loc)
557 orig_rot_mat = orig_rot.to_matrix().to_4x4()
558 orig_scale_mat = (Matrix.Scale(orig_scale[0],4,(1,0,0)) @
559 Matrix.Scale(orig_scale[1],4,(0,1,0)) @
560 Matrix.Scale(orig_scale[2],4,(0,0,1)))
562 # Assemble the new matrix.
563 obj.matrix_world = orig_loc_mat @ rot_mat @ orig_rot_mat @ orig_scale_mat
566 # Function, which puts a camera and light source into the 3D scene
567 def camera_light_source(use_camera,
568 use_light,
569 object_center_vec,
570 object_size):
572 camera_factor = 15.0
574 # If chosen a camera is put into the scene.
575 if use_camera == True:
577 # Assume that the object is put into the global origin. Then, the
578 # camera is moved in x and z direction, not in y. The object has its
579 # size at distance sqrt(object_size) from the origin. So, move the
580 # camera by this distance times a factor of camera_factor in x and z.
581 # Then add x, y and z of the origin of the object.
582 object_camera_vec = Vector((sqrt(object_size) * camera_factor,
583 0.0,
584 sqrt(object_size) * camera_factor))
585 camera_xyz_vec = object_center_vec + object_camera_vec
587 # Create the camera
588 camera_data = bpy.data.cameras.new("A_camera")
589 camera_data.lens = 45
590 camera_data.clip_end = 500.0
591 camera = bpy.data.objects.new("A_camera", camera_data)
592 camera.location = camera_xyz_vec
593 bpy.context.collection.objects.link(camera)
595 # Here the camera is rotated such it looks towards the center of
596 # the object. The [0.0, 0.0, 1.0] vector along the z axis
597 z_axis_vec = Vector((0.0, 0.0, 1.0))
598 # The angle between the last two vectors
599 angle = object_camera_vec.angle(z_axis_vec, 0)
600 # The cross-product of z_axis_vec and object_camera_vec
601 axis_vec = z_axis_vec.cross(object_camera_vec)
602 # Rotate 'axis_vec' by 'angle' and convert this to euler parameters.
603 # 4 is the size of the matrix.
604 camera.rotation_euler = Matrix.Rotation(angle, 4, axis_vec).to_euler()
606 # Rotate the camera around its axis by 90° such that we have a nice
607 # camera position and view onto the object.
608 bpy.ops.object.select_all(action='DESELECT')
609 camera.select_set(True)
611 # Rotate the camera around its axis 'object_camera_vec' by 90° such
612 # that we have a nice camera view onto the object.
613 matrix_rotation = Matrix.Rotation(90/360*2*pi, 4, object_camera_vec)
614 rotate_object(matrix_rotation, camera)
616 # Here a lamp is put into the scene, if chosen.
617 if use_light == True:
619 # This is the distance from the object measured in terms of %
620 # of the camera distance. It is set onto 50% (1/2) distance.
621 light_dl = sqrt(object_size) * 15 * 0.5
622 # This is a factor to which extend the lamp shall go to the right
623 # (from the camera point of view).
624 light_dy_right = light_dl * (3.0/4.0)
626 # Create x, y and z for the lamp.
627 object_light_vec = Vector((light_dl,light_dy_right,light_dl))
628 light_xyz_vec = object_center_vec + object_light_vec
630 # Create the lamp
631 light_data = bpy.data.lights.new(name="A_light", type="SUN")
632 light_data.distance = 500.0
633 light_data.energy = 3.0
634 lamp = bpy.data.objects.new("A_light", light_data)
635 lamp.location = light_xyz_vec
636 bpy.context.collection.objects.link(lamp)
638 # Some settings for the World: a bit ambient occlusion
639 bpy.context.scene.world.light_settings.use_ambient_occlusion = True
640 bpy.context.scene.world.light_settings.ao_factor = 0.1
641 # Some properties for cycles
642 lamp.data.use_nodes = True
643 lmp_P_BSDF = lamp.data.node_tree.nodes['Emission']
644 lmp_P_BSDF.inputs['Strength'].default_value = 5
647 # Function, which draws the atoms of one type (balls). This is one
648 # dupliverts structure then.
649 # Return: the dupliverts structure
650 def draw_atoms_one_type(draw_all_atoms_type,
651 Ball_type,
652 Ball_azimuth,
653 Ball_zenith,
654 Ball_radius_factor,
655 object_center_vec,
656 collection_molecule):
658 # Create the vertices composed of the coordinates of all atoms of one type
659 atom_vertices = []
660 for atom in draw_all_atoms_type:
661 # In fact, the object is created in the World's origin.
662 # This is why 'object_center_vec' is subtracted. At the end
663 # the whole object is translated back to 'object_center_vec'.
664 atom_vertices.append(atom[2] - object_center_vec)
666 # IMPORTANT: First, we create a collection of the element, which contains
667 # the atoms (balls + mesh) AND the sticks! The definition dealing with the
668 # sticks will put the sticks inside this collection later on.
669 coll_element_name = atom[0] # the element name
670 # Create the new collection and ...
671 coll_element = bpy.data.collections.new(coll_element_name)
672 # ... link it to the collection, which contains all parts of the
673 # molecule.
674 collection_molecule.children.link(coll_element)
676 # Now, create a collection for the atoms, which includes the representative
677 # ball and the mesh.
678 coll_atom_name = atom[0] + "_atom"
679 # Create the new collection and ...
680 coll_atom = bpy.data.collections.new(coll_atom_name)
681 # ... link it to the collection, which contains all parts of the
682 # element (ball and mesh).
683 coll_element.children.link(coll_atom)
685 # Build the mesh
686 atom_mesh = bpy.data.meshes.new("Mesh_"+atom[0])
687 atom_mesh.from_pydata(atom_vertices, [], [])
688 atom_mesh.update()
689 new_atom_mesh = bpy.data.objects.new(atom[0] + "_mesh", atom_mesh)
691 # Link active object to the new collection
692 coll_atom.objects.link(new_atom_mesh)
694 # Now, build a representative sphere (atom).
695 if atom[0] == "Vacancy":
696 bpy.ops.mesh.primitive_cube_add(
697 align='WORLD', enter_editmode=False,
698 location=(0.0, 0.0, 0.0),
699 rotation=(0.0, 0.0, 0.0))
700 else:
701 # NURBS balls
702 if Ball_type == "0":
703 bpy.ops.surface.primitive_nurbs_surface_sphere_add(
704 align='WORLD', enter_editmode=False,
705 location=(0,0,0), rotation=(0.0, 0.0, 0.0))
706 # UV balls
707 elif Ball_type == "1":
708 bpy.ops.mesh.primitive_uv_sphere_add(
709 segments=Ball_azimuth, ring_count=Ball_zenith,
710 align='WORLD', enter_editmode=False,
711 location=(0,0,0), rotation=(0, 0, 0))
712 # Meta balls
713 elif Ball_type == "2":
714 bpy.ops.object.metaball_add(type='BALL', align='WORLD',
715 enter_editmode=False, location=(0, 0, 0),
716 rotation=(0, 0, 0))
718 ball = bpy.context.view_layer.objects.active
719 # Hide this ball because its appearance has no meaning. It is just the
720 # representative ball. The ball is visible at the vertices of the mesh.
721 # Rememmber, this is a dupliverts construct!
722 # However, hiding does not work with meta balls!
723 if Ball_type == "0" or Ball_type == "1":
724 ball.hide_set(True)
725 # Scale up/down the ball radius.
726 ball.scale = (atom[3]*Ball_radius_factor,) * 3
728 if atom[0] == "Vacancy":
729 ball.name = atom[0] + "_cube"
730 else:
731 ball.name = atom[0] + "_ball"
733 ball.active_material = atom[1]
734 ball.parent = new_atom_mesh
735 new_atom_mesh.instance_type = 'VERTS'
736 # The object is back translated to 'object_center_vec'.
737 new_atom_mesh.location = object_center_vec
739 # Note the collection where the ball was placed into.
740 coll_all = ball.users_collection
741 if len(coll_all) > 0:
742 coll_past = coll_all[0]
743 else:
744 coll_past = bpy.context.scene.collection
746 # Put the atom into the new collection 'atom' and ...
747 coll_atom.objects.link(ball)
748 # ... unlink the atom from the other collection.
749 coll_past.objects.unlink(ball)
751 return new_atom_mesh, coll_element
754 # Function, which draws the sticks with help of the dupliverts technique.
755 # Return: list of dupliverts structures.
756 def draw_sticks_dupliverts(all_atoms,
757 atom_all_types_list,
758 center,
759 all_sticks,
760 Stick_diameter,
761 Stick_sectors,
762 Stick_unit,
763 Stick_dist,
764 use_sticks_smooth,
765 use_sticks_color,
766 list_coll_elements):
768 dl = Stick_unit
770 if use_sticks_color == False:
771 stick_material = bpy.data.materials.new(ELEMENTS[-1].name)
772 stick_material.diffuse_color = ELEMENTS[-1].color
774 # Sort the sticks and put them into a new list such that ...
775 sticks_all_lists = []
776 if use_sticks_color == True:
777 for atom_type in atom_all_types_list:
778 if atom_type[0] == "TER":
779 continue
780 sticks_list = []
781 for stick in all_sticks:
782 for repeat in range(stick.number):
784 atom1 = copy(all_atoms[stick.atom1-1].location)-center
785 atom2 = copy(all_atoms[stick.atom2-1].location)-center
787 dist = Stick_diameter * Stick_dist
789 if stick.number == 2:
790 if repeat == 0:
791 atom1 += (stick.dist * dist)
792 atom2 += (stick.dist * dist)
793 if repeat == 1:
794 atom1 -= (stick.dist * dist)
795 atom2 -= (stick.dist * dist)
797 if stick.number == 3:
798 if repeat == 0:
799 atom1 += (stick.dist * dist)
800 atom2 += (stick.dist * dist)
801 if repeat == 2:
802 atom1 -= (stick.dist * dist)
803 atom2 -= (stick.dist * dist)
805 dv = atom1 - atom2
806 n = dv / dv.length
807 if atom_type[0] == all_atoms[stick.atom1-1].name:
808 location = atom1
809 name = "_" + all_atoms[stick.atom1-1].name
810 material = all_atoms[stick.atom1-1].material
811 sticks_list.append([name, location, dv, material])
812 if atom_type[0] == all_atoms[stick.atom2-1].name:
813 location = atom1 - n * dl * int(ceil(dv.length / (2.0 * dl)))
814 name = "_" + all_atoms[stick.atom2-1].name
815 material = all_atoms[stick.atom2-1].material
816 sticks_list.append([name, location, dv, material])
818 if sticks_list != []:
819 sticks_all_lists.append(sticks_list)
820 else:
821 sticks_list = []
822 for stick in all_sticks:
824 if stick.number > 3:
825 stick.number = 1
827 for repeat in range(stick.number):
829 atom1 = copy(all_atoms[stick.atom1-1].location)-center
830 atom2 = copy(all_atoms[stick.atom2-1].location)-center
832 dist = Stick_diameter * Stick_dist
834 if stick.number == 2:
835 if repeat == 0:
836 atom1 += (stick.dist * dist)
837 atom2 += (stick.dist * dist)
838 if repeat == 1:
839 atom1 -= (stick.dist * dist)
840 atom2 -= (stick.dist * dist)
841 if stick.number == 3:
842 if repeat == 0:
843 atom1 += (stick.dist * dist)
844 atom2 += (stick.dist * dist)
845 if repeat == 2:
846 atom1 -= (stick.dist * dist)
847 atom2 -= (stick.dist * dist)
849 dv = atom1 - atom2
850 n = dv / dv.length
851 location = atom1
852 material = stick_material
853 sticks_list.append(["", location, dv, material])
855 sticks_all_lists.append(sticks_list)
857 atom_object_list = []
858 # ... the sticks in the list can be drawn:
859 for stick_list in sticks_all_lists:
860 vertices = []
861 faces = []
862 i = 0
864 # What follows is school mathematics! :-)
865 for stick in stick_list:
867 dv = stick[2]
868 v1 = stick[1]
869 n = dv / dv.length
870 gamma = -n.dot(v1)
871 b = v1 + gamma * n
872 n_b = b / b.length
874 if use_sticks_color == True:
875 loops = int(ceil(dv.length / (2.0 * dl)))
876 else:
877 loops = int(ceil(dv.length / dl))
879 for j in range(loops):
881 g = v1 - n * dl / 2.0 - n * dl * j
882 p1 = g + n_b * Stick_diameter
883 p2 = g - n_b * Stick_diameter
884 p3 = g - n_b.cross(n) * Stick_diameter
885 p4 = g + n_b.cross(n) * Stick_diameter
887 vertices.append(p1)
888 vertices.append(p2)
889 vertices.append(p3)
890 vertices.append(p4)
891 faces.append((i*4+0,i*4+2,i*4+1,i*4+3))
892 i += 1
894 # Create a collection for the sticks, which includes the representative
895 # cylinders, cups and the mesh.
896 coll_name = stick[0][1:] + "_sticks"
897 # Create the collection and ...
898 coll = bpy.data.collections.new(coll_name)
899 # ... link it to the collection, which contains all parts of the
900 # element. 'stick[0][1:]' contains the name of the element!
901 for coll_element_from_list in list_coll_elements:
902 if stick[0][1:] in coll_element_from_list.name:
903 break
904 coll_element_from_list.children.link(coll)
906 # Build the mesh.
907 mesh = bpy.data.meshes.new("Sticks_"+stick[0][1:])
908 mesh.from_pydata(vertices, [], faces)
909 mesh.update()
910 new_mesh = bpy.data.objects.new(stick[0][1:]+"_sticks_mesh", mesh)
911 # Link active object to the new collection
912 coll.objects.link(new_mesh)
914 # Build the object.
915 # Get the cylinder from the 'build_stick' function.
916 object_stick = build_stick(Stick_diameter,
918 Stick_sectors,
919 stick[0][1:])
920 # Link active object to the new collection
921 coll.objects.link(object_stick[0])
922 coll.objects.link(object_stick[1])
924 # Hide these objects because their appearance has no meaning. They are
925 # just the representative objects. The cylinder and cups are visible at
926 # the vertices of the mesh. Rememmber, this is a dupliverts construct!
927 object_stick[0].hide_set(True)
928 object_stick[1].hide_set(True)
930 stick_cylinder = object_stick[0]
931 stick_cylinder.active_material = stick[3]
932 stick_cups = object_stick[1]
933 stick_cups.active_material = stick[3]
935 # Smooth the cylinders.
936 if use_sticks_smooth == True:
937 bpy.ops.object.select_all(action='DESELECT')
938 stick_cylinder.select_set(True)
939 stick_cups.select_set(True)
940 bpy.ops.object.shade_smooth()
942 # Parenting the mesh to the cylinder.
943 stick_cylinder.parent = new_mesh
944 stick_cups.parent = new_mesh
945 new_mesh.instance_type = 'FACES'
946 new_mesh.location = center
947 atom_object_list.append(new_mesh)
949 # Return the list of dupliverts structures.
950 return atom_object_list
953 # Function, which draws the sticks with help of the skin and subdivision
954 # modifiers.
955 def draw_sticks_skin(all_atoms,
956 all_sticks,
957 Stick_diameter,
958 use_sticks_smooth,
959 sticks_subdiv_view,
960 sticks_subdiv_render,
961 coll_molecule):
963 # These counters are for the edges, in the shape [i,i+1].
964 i = 0
966 # This is the list of vertices, containing the atom position
967 # (vectors)).
968 stick_vertices = []
969 # This is the 'same' list, which contains not vector position of
970 # the atoms but their numbers. It is used to handle the edges.
971 stick_vertices_nr = []
972 # This is the list of edges.
973 stick_edges = []
975 # Go through the list of all sticks. For each stick do:
976 for stick in all_sticks:
978 # Each stick has two atoms = two vertices.
981 [ 0,1 , 3,4 , 0,8 , 7,3]
982 [[0,1], [2,3], [4,5], [6,7]]
984 [ 0,1 , 3,4 , x,8 , 7,x] x:deleted
985 [[0,1], [2,3], [0,5], [6,2]]
988 # Check, if the vertex (atom) is already in the vertex list.
989 # edge: [s1,s2]
990 FLAG_s1 = False
991 s1 = 0
992 for stick2 in stick_vertices_nr:
993 if stick2 == stick.atom1-1:
994 FLAG_s1 = True
995 break
996 s1 += 1
997 FLAG_s2 = False
998 s2 = 0
999 for stick2 in stick_vertices_nr:
1000 if stick2 == stick.atom2-1:
1001 FLAG_s2 = True
1002 break
1003 s2 += 1
1005 # If the vertex (atom) is not yet in the vertex list:
1006 # append the number of atom and the vertex to the two lists.
1007 # For the first atom:
1008 if FLAG_s1 == False:
1009 atom1 = copy(all_atoms[stick.atom1-1].location)
1010 stick_vertices.append(atom1)
1011 stick_vertices_nr.append(stick.atom1-1)
1012 # For the second atom:
1013 if FLAG_s2 == False:
1014 atom2 = copy(all_atoms[stick.atom2-1].location)
1015 stick_vertices.append(atom2)
1016 stick_vertices_nr.append(stick.atom2-1)
1018 # Build the edges:
1020 # If both vertices (atoms) were not in the lists, then
1021 # the edge is simply [i,i+1]. These are two new vertices
1022 # (atoms), so increase i by 2.
1023 if FLAG_s1 == False and FLAG_s2 == False:
1024 stick_edges.append([i,i+1])
1025 i += 2
1026 # Both vertices (atoms) were already in the list, so then
1027 # use the vertices (atoms), which already exist. They are
1028 # at positions s1 and s2.
1029 if FLAG_s1 == True and FLAG_s2 == True:
1030 stick_edges.append([s1,s2])
1031 # The following two if cases describe the situation that
1032 # only one vertex (atom) was in the list. Since only ONE
1033 # new vertex was added, increase i by one.
1034 if FLAG_s1 == True and FLAG_s2 == False:
1035 stick_edges.append([s1,i])
1036 i += 1
1037 if FLAG_s1 == False and FLAG_s2 == True:
1038 stick_edges.append([i,s2])
1039 i += 1
1041 # Build the mesh of the sticks
1042 stick_mesh = bpy.data.meshes.new("Mesh_sticks")
1043 stick_mesh.from_pydata(stick_vertices, stick_edges, [])
1044 stick_mesh.update()
1045 new_stick_mesh = bpy.data.objects.new("Sticks", stick_mesh)
1046 # Link the active mesh to the molecule collection
1047 coll_molecule.objects.link(new_stick_mesh)
1049 # Apply the skin modifier.
1050 new_stick_mesh.modifiers.new(name="Sticks_skin", type='SKIN')
1051 # Smooth the skin surface if this option has been chosen.
1052 new_stick_mesh.modifiers[0].use_smooth_shade = use_sticks_smooth
1053 # Apply the Subdivision modifier.
1054 new_stick_mesh.modifiers.new(name="Sticks_subsurf", type='SUBSURF')
1055 # Options: choose the levels
1056 new_stick_mesh.modifiers[1].levels = sticks_subdiv_view
1057 new_stick_mesh.modifiers[1].render_levels = sticks_subdiv_render
1059 stick_material = bpy.data.materials.new(ELEMENTS[-1].name)
1060 stick_material.diffuse_color = ELEMENTS[-1].color
1061 new_stick_mesh.active_material = stick_material
1063 # This is for putting the radius of the sticks onto
1064 # the desired value 'Stick_diameter'
1065 bpy.context.view_layer.objects.active = new_stick_mesh
1066 # EDIT mode
1067 bpy.ops.object.mode_set(mode='EDIT', toggle=False)
1068 bm = bmesh.from_edit_mesh(new_stick_mesh.data)
1069 bpy.ops.mesh.select_all(action='DESELECT')
1071 # Select all vertices
1072 for v in bm.verts:
1073 v.select = True
1075 # This is somewhat a factor for the radius.
1076 r_f = 4.0
1077 # Apply operator 'skin_resize'.
1078 bpy.ops.transform.skin_resize(
1079 value=(
1080 Stick_diameter * r_f,
1081 Stick_diameter * r_f,
1082 Stick_diameter * r_f,
1084 constraint_axis=(False, False, False),
1085 orient_type='GLOBAL',
1086 mirror=False,
1087 use_proportional_edit=False,
1088 snap=False,
1089 snap_target='CLOSEST',
1090 snap_point=(0, 0, 0),
1091 snap_align=False,
1092 snap_normal=(0, 0, 0),
1093 release_confirm=False,
1095 # Back to the OBJECT mode.
1096 bpy.ops.object.mode_set(mode='OBJECT', toggle=False)
1098 return new_stick_mesh
1101 # Draw the sticks the normal way: connect the atoms by simple cylinders.
1102 # Two options: 1. single cylinders parented to an empty
1103 # 2. one single mesh object
1104 def draw_sticks_normal(all_atoms,
1105 all_sticks,
1106 center,
1107 Stick_diameter,
1108 Stick_sectors,
1109 use_sticks_smooth,
1110 use_sticks_one_object,
1111 use_sticks_one_object_nr,
1112 coll_molecule):
1114 stick_material = bpy.data.materials.new(ELEMENTS[-1].name)
1115 stick_material.diffuse_color = ELEMENTS[-1].color
1117 up_axis = Vector([0.0, 0.0, 1.0])
1119 # For all sticks, do ...
1120 list_group = []
1121 list_group_sub = []
1122 counter = 0
1123 for stick in all_sticks:
1125 # The vectors of the two atoms
1126 atom1 = all_atoms[stick.atom1-1].location-center
1127 atom2 = all_atoms[stick.atom2-1].location-center
1128 # Location
1129 location = (atom1 + atom2) * 0.5
1130 # The difference of both vectors
1131 v = (atom2 - atom1)
1132 # Angle with respect to the z-axis
1133 angle = v.angle(up_axis, 0)
1134 # Cross-product between v and the z-axis vector. It is the
1135 # vector of rotation.
1136 axis = up_axis.cross(v)
1137 # Calculate Euler angles
1138 euler = Matrix.Rotation(angle, 4, axis).to_euler()
1139 # Create stick
1140 stick = bpy.ops.mesh.primitive_cylinder_add(vertices=Stick_sectors,
1141 radius=Stick_diameter,
1142 depth=v.length,
1143 end_fill_type='NGON',
1144 align='WORLD',
1145 enter_editmode=False,
1146 location=location,
1147 rotation=(0, 0, 0))
1148 # Put the stick into the scene ...
1149 stick = bpy.context.view_layer.objects.active
1150 # ... and rotate the stick.
1151 stick.rotation_euler = euler
1152 # ... and name
1153 stick.name = "Stick_Cylinder"
1154 counter += 1
1156 # Smooth the cylinder.
1157 if use_sticks_smooth == True:
1158 bpy.ops.object.select_all(action='DESELECT')
1159 stick.select_set(True)
1160 bpy.ops.object.shade_smooth()
1162 list_group_sub.append(stick)
1164 if use_sticks_one_object == True:
1165 if counter == use_sticks_one_object_nr:
1166 bpy.ops.object.select_all(action='DESELECT')
1167 for stick in list_group_sub:
1168 stick.select_set(True)
1169 bpy.ops.object.join()
1170 list_group.append(bpy.context.view_layer.objects.active)
1171 bpy.ops.object.select_all(action='DESELECT')
1172 list_group_sub = []
1173 counter = 0
1174 else:
1175 # Material ...
1176 stick.active_material = stick_material
1178 if use_sticks_one_object == True:
1179 bpy.ops.object.select_all(action='DESELECT')
1180 for stick in list_group_sub:
1181 stick.select_set(True)
1182 bpy.ops.object.join()
1183 list_group.append(bpy.context.view_layer.objects.active)
1184 bpy.ops.object.select_all(action='DESELECT')
1186 for group in list_group:
1187 group.select_set(True)
1188 bpy.ops.object.join()
1189 bpy.ops.object.origin_set(type='ORIGIN_GEOMETRY',
1190 center='MEDIAN')
1191 sticks = bpy.context.view_layer.objects.active
1192 sticks.active_material = stick_material
1194 sticks.location += center
1196 # Collections
1197 # ===========
1198 # Note the collection where the sticks were placed into.
1199 coll_all = sticks.users_collection
1200 if len(coll_all) > 0:
1201 coll_past = coll_all[0]
1202 else:
1203 coll_past = bpy.context.scene.collection
1205 # Link the sticks with the collection of the molecule ...
1206 coll_molecule.objects.link(sticks)
1207 # ... and unlink them from the collection it has been before.
1208 coll_past.objects.unlink(sticks)
1210 return sticks
1211 else:
1212 # Here we use an empty ...
1213 bpy.ops.object.empty_add(type='ARROWS',
1214 align='WORLD',
1215 location=(0, 0, 0),
1216 rotation=(0, 0, 0))
1217 sticks_empty = bpy.context.view_layer.objects.active
1218 sticks_empty.name = "A_sticks_empty"
1219 # ... that is parent to all sticks. With this, we can better move
1220 # all sticks if necessary.
1221 for stick in list_group_sub:
1222 stick.parent = sticks_empty
1224 sticks_empty.location += center
1226 # Collections
1227 # ===========
1228 # Create a collection that will contain all sticks + the empty and ...
1229 coll = bpy.data.collections.new("Sticks")
1230 # ... link it to the collection, which contains all parts of the
1231 # molecule.
1232 coll_molecule.children.link(coll)
1233 # Now, create a collection that only contains the sticks and ...
1234 coll_cylinder = bpy.data.collections.new("Sticks_cylinders")
1235 # ... link it to the collection, which contains the sticks and empty.
1236 coll.children.link(coll_cylinder)
1238 # Note the collection where the empty was placed into, ...
1239 coll_all = sticks_empty.users_collection
1240 if len(coll_all) > 0:
1241 coll_past = coll_all[0]
1242 else:
1243 coll_past = bpy.context.scene.collection
1244 # ... link the empty with the new collection ...
1245 coll.objects.link(sticks_empty)
1246 # ... and unlink it from the old collection where it has been before.
1247 coll_past.objects.unlink(sticks_empty)
1249 # Note the collection where the cylinders were placed into, ...
1250 coll_all = list_group_sub[0].users_collection
1251 if len(coll_all) > 0:
1252 coll_past = coll_all[0]
1253 else:
1254 coll_past = bpy.context.scene.collection
1256 for stick in list_group_sub:
1257 # ... link each stick with the new collection ...
1258 coll_cylinder.objects.link(stick)
1259 # ... and unlink it from the old collection.
1260 coll_past.objects.unlink(stick)
1262 return sticks_empty
1265 # -----------------------------------------------------------------------------
1266 # The main routine
1268 def import_pdb(Ball_type,
1269 Ball_azimuth,
1270 Ball_zenith,
1271 Ball_radius_factor,
1272 radiustype,
1273 Ball_distance_factor,
1274 use_sticks,
1275 use_sticks_type,
1276 sticks_subdiv_view,
1277 sticks_subdiv_render,
1278 use_sticks_color,
1279 use_sticks_smooth,
1280 use_sticks_bonds,
1281 use_sticks_one_object,
1282 use_sticks_one_object_nr,
1283 Stick_unit, Stick_dist,
1284 Stick_sectors,
1285 Stick_diameter,
1286 put_to_center,
1287 use_camera,
1288 use_light,
1289 filepath_pdb):
1291 # List of materials
1292 atom_material_list = []
1294 # A list of ALL objects which are loaded (needed for selecting the loaded
1295 # structure.
1296 atom_object_list = []
1298 # ------------------------------------------------------------------------
1299 # INITIALIZE THE ELEMENT LIST
1301 read_elements()
1303 # ------------------------------------------------------------------------
1304 # READING DATA OF ATOMS
1306 (Number_of_total_atoms, all_atoms) = read_pdb_file(filepath_pdb, radiustype)
1308 # ------------------------------------------------------------------------
1309 # MATERIAL PROPERTIES FOR ATOMS
1311 # The list that contains info about all types of atoms is created
1312 # here. It is used for building the material properties for
1313 # instance (see below).
1314 atom_all_types_list = []
1316 for atom in all_atoms:
1317 FLAG_FOUND = False
1318 for atom_type in atom_all_types_list:
1319 # If the atom name is already in the list, FLAG on 'True'.
1320 if atom_type[0] == atom.name:
1321 FLAG_FOUND = True
1322 break
1323 # No name in the current list has been found? => New entry.
1324 if FLAG_FOUND == False:
1325 # Stored are: Atom label (e.g. 'Na'), the corresponding atom
1326 # name (e.g. 'Sodium') and its color.
1327 atom_all_types_list.append([atom.name, atom.element, atom.color])
1329 # The list of materials is built.
1330 # Note that all atoms of one type (e.g. all hydrogens) get only ONE
1331 # material! This is good because then, by activating one atom in the
1332 # Blender scene and changing the color of this atom, one changes the color
1333 # of ALL atoms of the same type at the same time.
1335 # Create first a new list of materials for each type of atom
1336 # (e.g. hydrogen)
1337 for atom_type in atom_all_types_list:
1338 material = bpy.data.materials.new(atom_type[1])
1339 material.name = atom_type[0]
1340 material.diffuse_color = atom_type[2]
1341 atom_material_list.append(material)
1343 # Now, we go through all atoms and give them a material. For all atoms ...
1344 for atom in all_atoms:
1345 # ... and all materials ...
1346 for material in atom_material_list:
1347 # ... select the correct material for the current atom via
1348 # comparison of names ...
1349 if atom.name in material.name:
1350 # ... and give the atom its material properties.
1351 # However, before we check, if it is a vacancy, because then it
1352 # gets some additional preparation. The vacancy is represented
1353 # by a transparent cube.
1354 if atom.name == "Vacancy":
1355 # Some properties for eevee.
1356 material.metallic = 0.8
1357 material.specular_intensity = 0.5
1358 material.roughness = 0.3
1359 material.blend_method = 'OPAQUE'
1360 material.show_transparent_back = False
1361 # Some properties for cycles
1362 material.use_nodes = True
1363 mat_P_BSDF = material.node_tree.nodes['Principled BSDF']
1364 mat_P_BSDF.inputs['Metallic'].default_value = 0.1
1365 mat_P_BSDF.inputs['Roughness'].default_value = 0.2
1366 mat_P_BSDF.inputs['Transmission'].default_value = 0.97
1367 mat_P_BSDF.inputs['IOR'].default_value = 0.8
1368 # The atom gets its properties.
1369 atom.material = material
1371 # ------------------------------------------------------------------------
1372 # READING DATA OF STICKS
1374 all_sticks = read_pdb_file_sticks(filepath_pdb,
1375 use_sticks_bonds,
1376 all_atoms)
1378 # So far, all atoms, sticks and materials have been registered.
1381 # ------------------------------------------------------------------------
1382 # TRANSLATION OF THE STRUCTURE TO THE ORIGIN
1384 # It may happen that the structure in a PDB file already has an offset
1385 # If chosen, the structure is first put into the center of the scene
1386 # (the offset is subtracted).
1388 if put_to_center == True:
1389 sum_vec = Vector((0.0,0.0,0.0))
1390 # Sum of all atom coordinates
1391 sum_vec = sum([atom.location for atom in all_atoms], sum_vec)
1392 # Then the average is taken
1393 sum_vec = sum_vec / Number_of_total_atoms
1394 # After, for each atom the center of gravity is subtracted
1395 for atom in all_atoms:
1396 atom.location -= sum_vec
1398 # ------------------------------------------------------------------------
1399 # SCALING
1401 # Take all atoms and adjust their radii and scale the distances.
1402 for atom in all_atoms:
1403 atom.location *= Ball_distance_factor
1405 # ------------------------------------------------------------------------
1406 # DETERMINATION OF SOME GEOMETRIC PROPERTIES
1408 # In the following, some geometric properties of the whole object are
1409 # determined: center, size, etc.
1410 sum_vec = Vector((0.0,0.0,0.0))
1412 # First the center is determined. All coordinates are summed up ...
1413 sum_vec = sum([atom.location for atom in all_atoms], sum_vec)
1415 # ... and the average is taken. This gives the center of the object.
1416 object_center_vec = sum_vec / Number_of_total_atoms
1418 # Now, we determine the size.The farthest atom from the object center is
1419 # taken as a measure. The size is used to place well the camera and light
1420 # into the scene.
1421 object_size_vec = [atom.location - object_center_vec for atom in all_atoms]
1422 object_size = max(object_size_vec).length
1424 # ------------------------------------------------------------------------
1425 # SORTING THE ATOMS
1427 # Lists of atoms of one type are created. Example:
1428 # draw_all_atoms = [ data_hydrogen,data_carbon,data_nitrogen ]
1429 # data_hydrogen = [["Hydrogen", Material_Hydrogen, Vector((x,y,z)), 109], ...]
1431 # Go through the list which contains all types of atoms. It is the list,
1432 # which has been created on the top during reading the PDB file.
1433 # Example: atom_all_types_list = ["hydrogen", "carbon", ...]
1434 draw_all_atoms = []
1435 for atom_type in atom_all_types_list:
1437 # Don't draw 'TER atoms'.
1438 if atom_type[0] == "TER":
1439 continue
1441 # This is the draw list, which contains all atoms of one type (e.g.
1442 # all hydrogens) ...
1443 draw_all_atoms_type = []
1445 # Go through all atoms ...
1446 for atom in all_atoms:
1447 # ... select the atoms of the considered type via comparison ...
1448 if atom.name == atom_type[0]:
1449 # ... and append them to the list 'draw_all_atoms_type'.
1450 draw_all_atoms_type.append([atom.name,
1451 atom.material,
1452 atom.location,
1453 atom.radius])
1455 # Now append the atom list to the list of all types of atoms
1456 draw_all_atoms.append(draw_all_atoms_type)
1458 # ------------------------------------------------------------------------
1459 # COLLECTION
1461 # Before we start to draw the atoms and sticks, we first create a
1462 # collection for the molecule. All atoms (balls) and sticks (cylinders)
1463 # are put into this collection.
1464 coll_molecule_name = os.path.basename(filepath_pdb)
1465 scene = bpy.context.scene
1466 coll_molecule = bpy.data.collections.new(coll_molecule_name)
1467 scene.collection.children.link(coll_molecule)
1469 # ------------------------------------------------------------------------
1470 # DRAWING THE ATOMS
1472 bpy.ops.object.select_all(action='DESELECT')
1474 list_coll_elements = []
1475 # For each list of atoms of ONE type (e.g. Hydrogen)
1476 for draw_all_atoms_type in draw_all_atoms:
1478 atom_mesh, coll_element = draw_atoms_one_type(draw_all_atoms_type,
1479 Ball_type,
1480 Ball_azimuth,
1481 Ball_zenith,
1482 Ball_radius_factor,
1483 object_center_vec,
1484 coll_molecule)
1485 atom_object_list.append(atom_mesh)
1486 list_coll_elements.append(coll_element)
1488 # ------------------------------------------------------------------------
1489 # DRAWING THE STICKS: cylinders in a dupliverts structure
1491 if use_sticks == True and use_sticks_type == '0' and all_sticks != []:
1493 sticks = draw_sticks_dupliverts(all_atoms,
1494 atom_all_types_list,
1495 object_center_vec,
1496 all_sticks,
1497 Stick_diameter,
1498 Stick_sectors,
1499 Stick_unit,
1500 Stick_dist,
1501 use_sticks_smooth,
1502 use_sticks_color,
1503 list_coll_elements)
1504 for stick in sticks:
1505 atom_object_list.append(stick)
1507 # ------------------------------------------------------------------------
1508 # DRAWING THE STICKS: skin and subdivision modifier
1510 if use_sticks == True and use_sticks_type == '1' and all_sticks != []:
1512 sticks = draw_sticks_skin(all_atoms,
1513 all_sticks,
1514 Stick_diameter,
1515 use_sticks_smooth,
1516 sticks_subdiv_view,
1517 sticks_subdiv_render,
1518 coll_molecule)
1519 atom_object_list.append(sticks)
1521 # ------------------------------------------------------------------------
1522 # DRAWING THE STICKS: normal cylinders
1524 if use_sticks == True and use_sticks_type == '2' and all_sticks != []:
1526 sticks = draw_sticks_normal(all_atoms,
1527 all_sticks,
1528 object_center_vec,
1529 Stick_diameter,
1530 Stick_sectors,
1531 use_sticks_smooth,
1532 use_sticks_one_object,
1533 use_sticks_one_object_nr,
1534 coll_molecule)
1535 atom_object_list.append(sticks)
1537 # ------------------------------------------------------------------------
1538 # CAMERA and LIGHT SOURCES
1540 camera_light_source(use_camera,
1541 use_light,
1542 object_center_vec,
1543 object_size)
1545 # ------------------------------------------------------------------------
1546 # SELECT ALL LOADED OBJECTS
1547 bpy.ops.object.select_all(action='DESELECT')
1548 obj = None
1549 for obj in atom_object_list:
1550 obj.select_set(True)
1552 # activate the last selected object
1553 if obj:
1554 bpy.context.view_layer.objects.active = obj