-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathvis03_plot_tangent_product_kymograph.py
65 lines (51 loc) · 1.75 KB
/
vis03_plot_tangent_product_kymograph.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
'''
Result:
- tangents change slowly
- curvature calculation - TODO might be wrong; TODO - fix
'''
import os.path
import numpy as np
from scipy.linalg import norm
import matplotlib.pyplot as plt
input_folder = 'data/rescale'
output_folder = 'res/rescale_new_angles'
iframes_raw = list(range(1, 18, 2)) + [18] + list(range(21, 42, 2))
t = []
t_prod = []
kappalist = []
for iframe_raw in iframes_raw: #[15]: #
xname = os.path.join(output_folder, 'x-{}.dat'.format(iframe_raw))
yname = os.path.join(output_folder, 'y-{}.dat'.format(iframe_raw))
zname = os.path.join(output_folder, 'z-{}.dat'.format(iframe_raw))
xx = np.loadtxt(xname)
yy = np.loadtxt(yname)
zz = np.loadtxt(zname)
rr = np.array([xx, yy, zz]).T
drr = np.diff(rr, axis=0) # OK
tt = np.array([drrr / norm(drrr) for drrr in drr]) # OK
tangent_products = np.sum(tt[:-1] * tt[1:], axis=1)
t_prod.append(tangent_products)
t.append(tt)
# Test: curvature; formula from Taylor expansion of (t_i, t_{i+1}) at i-th point; + Frenet formulas
# Result: doesn't work well; The reason is probably that the formula has O(ds ** (1/2)) accuracy
assert np.std(norm(drr, axis=1)) < 0.0001
ds = np.mean(norm(drr, axis=1))
kappa = (2 * tangent_products) ** (1 / 2) / ds
kappalist.append(kappa)
# tangent
filename = None # os.path.join(output_folder, 'kymograph_psi.png')
plt.title(r'$t_i \cdot t_{i+1}$')
plt.xlabel(r'icoord')
plt.ylabel(r'Frame')
norm = None # quick.MidpointNormalize(vmin=-1, vmax=1, midpoint=0)
cmap = 'Reds'
vals = t_prod
plt.imshow(vals, cmap=cmap, norm=norm, origin='lower')
plt.colorbar()
plt.show()
# Constant along a shape ??
cmap = 'Reds'
vals = kappalist
plt.imshow(vals, cmap=cmap, norm=norm, origin='lower')
plt.colorbar()
plt.show()