Coverage for surface.py : 100%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1r"""
2.. role:: raw-math(raw) :format: latex html
4--------------------
5Surface
6--------------------
8Calculation of curvature requires a surface of reference. In MembraneCurvature,
9the surface of reference is defined by the `z` position of the `atoms` in `AtomGroup`.
12Functions
13---------
16"""
18import numpy as np
19import warnings
20import MDAnalysis
21import logging
23MDAnalysis.start_logging()
24logger = logging.getLogger("MDAnalysis.MDAKit.membrane_curvature")
27class WarnOnce:
28 """
29 Class to warn atoms out of grid boundaries only once with full message.
30 After the first ocurrance, message will be generic.
31 """
32 def __init__(self, msg, msg_multiple) -> None:
33 self.msg = msg
34 self.warned = False
35 self.warned_multiple = False
36 self.msg_multiple = msg_multiple
38 def warn(self, *args):
39 if not self.warned:
40 self.warned = True
41 warnings.warn(self.msg.format(*args))
42 logger.warning(self.msg.format(*args))
43 elif not self.warned_multiple:
44 warnings.warn(self.msg_multiple)
47def derive_surface(atoms, n_cells_x, n_cells_y, max_width_x, max_width_y):
48 """
49 Derive surface from `atom` positions in `AtomGroup`.
51 Parameters
52 ----------
53 atoms: AtomGroup.
54 AtomGroup of reference selection to define the surface
55 of the membrane.
56 n_cells_x: int.
57 number of cells in the grid of size `max_width_x`, `x` axis.
58 n_cells_y: int.
59 number of cells in the grid of size `max_width_y`, `y` axis.
60 max_width_x: float.
61 Maximum width of simulation box in x axis. (Determined by simulation box dimensions)
62 max_width_y: float.
63 Maximum width of simulation box in y axis. (Determined by simulation box dimensions)
65 Returns
66 -------
67 z_coordinates: numpy.ndarray
68 Average z-coordinate values. Return Numpy array of floats of
69 shape `(n_cells_x, n_cells_y)`.
71 """
72 coordinates = atoms.positions
73 return get_z_surface(coordinates, n_x_bins=n_cells_x, n_y_bins=n_cells_y,
74 x_range=(0, max_width_x), y_range=(0, max_width_y))
77# messages for warnings, negative and positive coordinates.
78msg_exceeds = "More than one atom exceed boundaries of grid."
79negative_coord_warning = WarnOnce("Atom with negative coordinates falls "
80 "outside grid boundaries. Element "
81 "({},{}) in grid can't be assigned."
82 " Skipping atom.", msg_exceeds)
83positive_coord_warning = WarnOnce("Atom coordinates exceed size of grid "
84 "and element ({},{}) can't be assigned. "
85 "Maximum (x,y) coordinates must be < ({}, {}). "
86 "Skipping atom.", msg_exceeds)
89def get_z_surface(coordinates, n_x_bins=10, n_y_bins=10, x_range=(0, 100), y_range=(0, 100)):
90 """
91 Derive surface from distribution of z coordinates in grid.
93 Parameters
94 ----------
95 coordinates : numpy.ndarray
96 Coordinates of AtomGroup. Numpy array of shape=(n_atoms, 3).
97 n_x_bins : int.
98 Number of bins in grid in the `x` dimension.
99 n_y_bins : int.
100 Number of bins in grid in the `y` dimension.
101 x_range : tuple of (float, float)
102 Range of coordinates (min, max) in the `x` dimension with shape=(2,).
103 y_range : tuple of (float, float)
104 Range of coordinates (min, max) in the `y` dimension with shape=(2,).
106 Returns
107 -------
108 z_surface: np.ndarray
109 Surface derived from set of coordinates in grid of `x_range, y_range` dimensions.
110 Returns Numpy array of floats of shape (`n_x_bins`, `n_y_bins`)
112 """
114 grid_z_coordinates = np.zeros((n_x_bins, n_y_bins))
115 grid_norm_unit = np.zeros((n_x_bins, n_y_bins))
117 x_factor = n_x_bins / (x_range[1] - x_range[0])
118 y_factor = n_y_bins / (y_range[1] - y_range[0])
120 x_coords, y_coords, z_coords = coordinates.T
122 cell_x_floor = np.floor(x_coords * x_factor).astype(int)
123 cell_y_floor = np.floor(y_coords * y_factor).astype(int)
125 for l, m, z in zip(cell_x_floor, cell_y_floor, z_coords):
127 try:
128 # negative coordinates
129 if l < 0 or m < 0:
130 negative_coord_warning.warn(l, m)
131 continue
133 grid_z_coordinates[l, m] += z
134 grid_norm_unit[l, m] += 1
136 # too large positive coordinates
137 except IndexError:
138 positive_coord_warning.warn(l, m, x_range[1], y_range[1])
140 z_surface = normalized_grid(grid_z_coordinates, grid_norm_unit)
142 return z_surface
145def normalized_grid(grid_z_coordinates, grid_norm_unit):
146 """
147 Calculates average `z` coordinates in unit cell.
149 Parameters
150 ----------
152 z_ref: np.array
153 Empty array of `(l,m)`
154 grid_z_coordinates: np.array
155 Array of size `(l,m)` with `z` coordinates stored in unit cell.
156 grid_norm_unit: np.array
157 Array of size `(l,m)` with number of atoms in unit cell.
159 Returns
160 -------
161 z_surface: np.ndarray
162 Normalized `z` coordinates in grid.
163 Returns Numpy array of floats of shape (`n_x_bins`, `n_y_bins`)
165 """
167 grid_norm_unit = np.where(grid_norm_unit > 0, grid_norm_unit, np.nan)
168 z_normalized = grid_z_coordinates / grid_norm_unit
170 return z_normalized