relax disp.spin lock offset+field figure
From relax wiki
Jump to navigation
Jump to search
Reference to original figure
Comparing to Figure 1 and 10 the reference:
- Palmer, 3rd, A. G. and Massi, F. (2006). Characterization of the dynamics of biomacromolecules using rotating-frame spin relaxation NMR spectroscopy. Chem. Rev., 106(5), 1700-1719. (DOI: 10.1021/cr0404287)
Script to produce figure
Python script: Reproduce figure 1 of Palmer and Massi, 2006 using matplotlib.
#############
# Made by Troels E. Linnet
# PhD student
# Copenhagen University
# SBiNLab, Structural Biology and NMR Laboratory.
# March 2014
#
# Trying to reproduce Figure 1 in:
#
# Palmer, A.G. & Massi, F. (2006). Characterization of the dynamics of biomacromolecules using rotating-frame spin relaxation NMR spectroscopy.
# Chem. Rev. 106, 1700-1719 http://dx.doi.org/10.1021/cr04042875
#
# This script is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#######################
import matplotlib.pyplot as plt
from matplotlib.patches import Circle, Ellipse, PathPatch, FancyArrowPatch
from mpl_toolkits.mplot3d import Axes3D, proj3d
import mpl_toolkits.mplot3d.art3d as art3d
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
import pylab
######## Define helper functions
class Arrow3D(FancyArrowPatch):
def __init__(self, xs, ys, zs, *args, **kwargs):
FancyArrowPatch.__init__(self, (0,0), (0,0), *args, **kwargs)
self._verts3d = xs, ys, zs
def draw(self, renderer):
xs3d, ys3d, zs3d = self._verts3d
xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, renderer.M)
self.set_positions((xs[0],ys[0]),(xs[1],ys[1]))
FancyArrowPatch.draw(self, renderer)
def rotation_matrix(d):
"""
Calculates a rotation matrix given a vector d. The direction of d
corresponds to the rotation axis. The length of d corresponds to
the sin of the angle of rotation.
Variant of: http://mail.scipy.org/pipermail/numpy-discussion/2009-March/040806.html
"""
sin_angle = np.linalg.norm(d)
if sin_angle == 0:
return np.identity(3)
d /= sin_angle
eye = np.eye(3)
ddt = np.outer(d, d)
skew = np.array([[ 0, d[2], -d[1]],
[-d[2], 0, d[0]],
[d[1], -d[0], 0]], dtype=np.float64)
M = ddt + np.sqrt(1 - sin_angle**2) * (eye - ddt) + sin_angle * skew
return M
def pathpatch_2d_to_3d(pathpatch, z = 0, normal = 'z'):
"""
Transforms a 2D Patch to a 3D patch using the given normal vector.
The patch is projected into they XY plane, rotated about the origin
and finally translated by z.
"""
if type(normal) is str: #Translate strings to normal vectors
index = "xyz".index(normal)
normal = np.roll((1,0,0), index)
normal /= np.linalg.norm(normal) #Make sure the vector is normalised
path = pathpatch.get_path() #Get the path and the associated transform
trans = pathpatch.get_patch_transform()
path = trans.transform_path(path) #Apply the transform
pathpatch.__class__ = art3d.PathPatch3D #Change the class
pathpatch._code3d = path.codes #Copy the codes
pathpatch._facecolor3d = pathpatch.get_facecolor #Get the face color
verts = path.vertices #Get the vertices in 2D
d = np.cross(normal, (0, 0, 1)) #Obtain the rotation vector
M = rotation_matrix(d) #Get the rotation matrix
pathpatch._segment3d = np.array([np.dot(M, (x, y, 0)) + (0, 0, z) for x, y in verts])
def pathpatch_translate(pathpatch, delta):
"""
Translates the 3D pathpatch by the amount delta.
"""
pathpatch._segment3d += delta
def update_position(e):
print "From update position"
#Transform co-ordinates to get new 2D projection
tX, tY, _ = proj3d.proj_transform(dataX, dataY, dataZ, ax.get_proj())
for i in range(len(dataX)):
label = labels[i]
label.xy = tX[i],tY[i]
label.update_positions(fig.canvas.renderer)
tX_t, tY_t, _ = proj3d.proj_transform(dataX_t, dataY_t, dataZ_t, ax.get_proj())
for i in range(len(dataX_t)):
label = labels_t[i]
label.xy = tX_t[i],tY_t[i]
label.update_positions(fig.canvas.renderer)
fig.canvas.draw()
return
#########
fig = plt.figure(figsize=(12, 12))
ax=fig.gca(projection='3d')
# Draw a Center point
ax.scatter([0],[0],[0],c="k",s=50)
# Draw circle
circle = Circle((0, 0), 1)
circle.set_color("k")
circle.set_alpha(80)
#circle.set_fill(False)
ax.add_patch(circle)
#art3d.pathpatch_2d_to_3d(circle, z=0, zdir='y')
pathpatch_2d_to_3d(circle, z=0, normal = [1, 0, 1] )
#######
#Input 3D Data for Text
Sx = [1,0,0]
Sy = [0,-1,0]
Sz = [0,0,1]
Szp = [1.1,0,1.1]
Sxp = [0.75,0,-0.75]
w1 = [0.7,0,0]
O1 = [0,0,0.3]
we1 = [w1[0],0,O1[2]]
O_av = [0.0,0,Szp[2]/Szp[0]*w1[0]]
we = [w1[0],0,Szp[2]/Szp[0]*w1[0]]
data = np.array([Sx, Sy, Sz , Szp , Sxp, w1, O1, we1, O_av, we])
textlab = [r"S$_x$", r"S$_y$=S$_y$'", r"S$_z$", r"S$_z$'", r"S$_x$'", r"$\omega_1$", r"$\Omega_1$", r"$\omega_{e1}$", r"$\bar{\Omega}$", r"$\omega_{e}$"]
#Separate into X, Y, Z for greater clarity
dataX = data[:,0]
dataY = data[:,1]
dataZ = data[:,2]
#3D scatter plot
ax.scatter(dataX, dataY, dataZ, marker = 'o', c='k', s=1)
#Transform co-ordinates to get initial 2D projection
tX, tY, _ = proj3d.proj_transform(dataX, dataY, dataZ, ax.get_proj())
#Array of labels
labels = []
#Loop through data points to initially annotate scatter plot
#and populate labels array
for i in range(len(dataX)):
#text='['+str(int(dataX[i]))+','+str(int(dataY[i]))+','+str(int(dataZ[i]))+']'
text=textlab[i]
label = ax.annotate(text,
xycoords='data',
xy = (tX[i], tY[i]), xytext = (+20, -20),
textcoords = 'offset points', ha = 'right', va = 'top', fontsize=12,
bbox = dict(boxstyle = 'round,pad=0.5', fc = 'grey', alpha = 0.8),
arrowprops = dict(arrowstyle = '->', connectionstyle = 'arc3,rad=0'))
labels.append(label)
#######
#######
#Input 3D Data for Ttheta
t1 = [0.5*we1[0],0,0.5*we1[2]]
t = [0.5*we[0],0,0.5*we[2]]
data_t = np.array([t1, t])
textlab_t = [r"$\theta_1$", r"$\theta$"]
#Separate into X, Y, Z for greater clarity
dataX_t = data_t[:,0]
dataY_t = data_t[:,1]
dataZ_t = data_t[:,2]
#3D scatter plot
ax.scatter(dataX_t, dataY_t, dataZ_t, marker = 'o', c='k', s=1)
#Transform co-ordinates to get initial 2D projection
tX_t, tY_t, _ = proj3d.proj_transform(dataX_t, dataY_t, dataZ_t, ax.get_proj())
#Array of labels
labels_t = []
#Loop through data points to initially annotate scatter plot
#and populate labels array
for i in range(len(dataX_t)):
#text='['+str(int(dataX[i]))+','+str(int(dataY[i]))+','+str(int(dataZ[i]))+']'
text=textlab_t[i]
label = ax.annotate(text,
xycoords='data',
xy = (tX[i], tY[i]), xytext = (-100, 50),
textcoords = 'offset points', ha = 'right', va = 'top', fontsize=12,
bbox = dict(boxstyle = 'round,pad=0.5', fc = 'grey', alpha = 0.8),
arrowprops = dict(arrowstyle = '->', connectionstyle = 'angle3,angleA=0,angleB=90'))
labels_t.append(label)
#######
# Make Sx
x, y, z = Sx
Sx = Arrow3D([-1, x],[0, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-|>", color="k")
ax.add_artist(Sx)
# Make Sy
x, y, z = Sy
fSy = Arrow3D([0, x],[1, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-|>", color="k")
ax.add_artist(fSy)
# Make Sz
x, y, z = Sz
fSz = Arrow3D([0, x],[0, y], [-1, z], mutation_scale=20, lw=1, arrowstyle="-|>", color="k")
ax.add_artist(fSz)
# Make Szp
x, y, z = Szp
fSzp = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=3, arrowstyle="-|>", color="k")
ax.add_artist(fSzp)
# Make Sxp
x, y, z = Sxp
fSxp = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=3, arrowstyle="-|>", color="k")
ax.add_artist(fSxp)
# Make w1
x, y, z = w1
fw1 = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-", color="k")
ax.add_artist(fw1)
# Make O1
x, y, z = O1
fO1 = Arrow3D([0, x],[0, y], [0, y], mutation_scale=20, lw=1, arrowstyle="-", color="k")
ax.add_artist(fO1)
# Make we1
x, y, z = we1
fwe1 = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-|>", color="b")
ax.add_artist(fwe1)
# Make O1 -> we1
fO1_we1 = Arrow3D([O1[0], we1[0]],[O1[1], we1[1]], [O1[2], we1[2]], mutation_scale=20, lw=1, arrowstyle="-", color="b")
ax.add_artist(fO1_we1)
# Make O_av
x, y, z = O_av
fO_av = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-", color="k")
ax.add_artist(fO_av)
# Make we
x, y, z = we
fwe = Arrow3D([0, x],[0, y], [0, z], mutation_scale=20, lw=1, arrowstyle="-|>", color="r")
ax.add_artist(fwe)
# Make w1 -> we
fw1_we = Arrow3D([w1[0], we[0]],[w1[1], we[1]], [w1[2], we[2]], mutation_scale=20, lw=1, arrowstyle="-", color="r")
ax.add_artist(fw1_we)
# Make O_ave -> we
fO_av_we = Arrow3D([O_av[0], we[0]],[O_av[1], we[1]], [O_av[2], we[2]], mutation_scale=20, lw=1, arrowstyle="-", color="r")
ax.add_artist(fO_av_we)
ax.set_xlim3d(-1, 1)
ax.set_ylim3d(-1, 1)
ax.set_zlim3d(-1, 1)
fig.canvas.mpl_connect('button_release_event', update_position)
plt.show()