-
Notifications
You must be signed in to change notification settings - Fork 0
/
Tilting.py
146 lines (114 loc) · 5.47 KB
/
Tilting.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
"""
Tilting.py
Defines functions for optimizing the tilt of a given set of movements with respect
to a CA trace structure defined by XYZ coordinates. The primary function, `optimize_tilt`,
iteratively rotates the movements in 3D space, evaluating each rotation using a cost function
to find the combination of movements that minimizes the cost. Auxiliary functions
`rotate_movements`, `create_cost_matrix`, and `get_cost` support this process by handling
the rotation calculations, creating a cost matrix, and calculating the cost of movements,
respectively.
"""
import numpy as np
def optimize_tilt(ca_trace, movements):
"""
Optimize the tilt of a given set of movements with respect to a CA trace structure
:param ca_trace:
:param movements:
:return: lowest_movements, cost_matrix, lowest_sum, [x, y, z]
"""
ca_trace_np = ca_trace[['X', 'Y', 'Z']].to_numpy()
lowest_sum = get_cost(ca_trace_np, movements)
lowest_movements = movements
x, y, z = 0, 0, 0
for i in np.arange(0, 375, 10):
for j in np.arange(0, 375, 10):
for k in np.arange(0, 375, 10):
if 0 <= i <= 360 and 0 <= j <= 360 and 0 <= k <= 360:
tilt_movements = rotate_movements(movements, i, j, k)
current_sum = get_cost(ca_trace_np, tilt_movements)
if current_sum < lowest_sum:
lowest_sum = current_sum
lowest_movements = tilt_movements
x, y, z = i, j, k
for i in np.arange(x - 10, x + 10):
for j in np.arange(y - 10, y + 10):
for k in np.arange(z - 10, z + 10):
if 0 <= i <= 360 and 0 <= j <= 360 and 0 <= k <= 360:
tilt_movements = rotate_movements(movements, i, j, k)
current_sum = get_cost(ca_trace_np, tilt_movements)
if current_sum < lowest_sum:
lowest_sum = current_sum
lowest_movements = tilt_movements
x, y, z = i, j, k
# Create normalized cost matrix
cost_matrix = create_cost_matrix(ca_trace_np, lowest_movements)
return lowest_movements, cost_matrix, lowest_sum, [x, y, z]
def rotate_movements(movements, degree_x, degree_y, degree_z):
# Convert degrees to radians
angles_rad_x = np.radians(degree_x)
angles_rad_y = np.radians(degree_y)
angles_rad_z = np.radians(degree_z)
# Compute sin and cos for each angle
sin_x, cos_x = np.sin(angles_rad_x), np.cos(angles_rad_x)
sin_y, cos_y = np.sin(angles_rad_y), np.cos(angles_rad_y)
sin_z, cos_z = np.sin(angles_rad_z), np.cos(angles_rad_z)
# Directly compute the combined rotation matrix
rotation_matrix = np.array([
[cos_y * cos_z, cos_z * sin_x * sin_y - cos_x * sin_z, sin_x * sin_z + cos_x * cos_z * sin_y],
[cos_y * sin_z, cos_x * cos_z + sin_x * sin_y * sin_z, cos_x * sin_y * sin_z - cos_z * sin_x],
[-sin_y, cos_y * sin_x, cos_x * cos_y]
])
# Apply the rotation to the movements
changed_movements = np.dot(movements, rotation_matrix)
return changed_movements
def create_cost_matrix(ca_trace, movements):
"""
Create a cost matrix for the given movements
:param ca_trace: NP Array
:param movements:
:return: cost_matrix
"""
# Compute normalized movement vectors
movement_vectors = ca_trace[1:] - ca_trace[:-1]
magnitudes = np.linalg.norm(movement_vectors, axis=1)
# Initialize the cost matrix
cost_matrix = np.zeros((movement_vectors.shape[0], len(movements)), dtype=np.float16)
# Iterate over movements to calculate costs
for i, movement in enumerate(movements):
# Convert movement to float16
movement = np.array(movement, dtype=np.float16)
# Normalize movement vectors and compute the cost
# Handle zero magnitudes to avoid division by zero
norm_movement_vectors = movement_vectors / magnitudes[:, np.newaxis]
norm_movement_vectors[magnitudes == 0] = 0
# Compute the dot product between the normalized movement vectors and the movement
dot_product = np.dot(norm_movement_vectors, movement)
dot_product_clamped = np.clip(dot_product, -1.0, 1.0)
# Compute the arccos of the clamped dot product
cost_matrix[:, i] = np.arccos(dot_product_clamped)
# Find the minimum value in each row
lowest_costs = np.min(cost_matrix, axis=1)
# Subtract the minimum value from each element in the row
cost_matrix = cost_matrix - lowest_costs[:, np.newaxis]
return cost_matrix
def get_cost(ca_trace, movements):
"""
Get the cost of the given movements
:param ca_trace: NP Array
:param movements:
:return: cost
"""
# Compute movement vectors and their magnitudes
movement_vectors = ca_trace[1:] - ca_trace[:-1]
magnitudes = np.linalg.norm(movement_vectors, axis=1)
# Initialize the cost matrix
cost_matrix = np.zeros((movement_vectors.shape[0], len(movements)))
# Iterate over movements to calculate costs
for i, movement in enumerate(movements):
# Normalize movement vectors and compute the cost
# Handle zero magnitudes to avoid division by zero
norm_movement_vectors = movement_vectors / magnitudes[:, np.newaxis]
norm_movement_vectors[magnitudes == 0] = 0
cost_matrix[:, i] = np.linalg.norm(norm_movement_vectors - movement, axis=1)
# Return the sum of minimum costs for each vector
return np.min(cost_matrix, axis=1).sum()