Open main menu

About

Code

File: r1rhor2eff.py

### python imports
import sys
import os
from math import cos, sin

### plotting facility.
import matplotlib.pyplot as plt

# Ordered dictionary
import collections

### relax modules
# Import some tools to loop over the spins.
from pipe_control.mol_res_spin import return_spin, spin_loop
# Import method to calculate the R1_rho offset data
from specific_analyses.relax_disp.disp_data import calc_rotating_frame_params, generate_r20_key, loop_exp_frq

###############

# You have to provide a DPL94 results state file
res_folder = "resultsR1"
res_state = os.path.join(res_folder, "DPL94", "results")
spin_inte = ":52@N"
spin_inte_rep = spin_inte.replace('#', '_').replace(':', '_').replace('@', '_')

# Load the state
state.load(res_state, force=True)

# Show pipes
pipe.display()
pipe.current()

# Get the spin of interest and save it in cdp, to access it after execution of script.
cdp.myspin = return_spin(spin_inte)

# Calculate the offset data
theta_spin_dic, Domega_spin_dic, w_eff_spin_dic, dic_key_list = calc_rotating_frame_params(spin=cdp.myspin, spin_id=spin_inte, verbosity=1)
# Save the data in cdp to access it after execution of script.
cdp.myspin.theta_spin_dic = theta_spin_dic
cdp.myspin.w_eff_spin_dic = w_eff_spin_dic
cdp.myspin.dic_key_list = dic_key_list

# Create spin.r2 keys
r2keys = []
for exp_type, frq in loop_exp_frq():
    r2keys.append(generate_r20_key(exp_type=exp_type, frq=frq) )
# Save the keys
cdp.r2keys = r2keys

x_w_eff = []
x_theta = []
x_disp_point = []
y = []
y_R1_rho = []

cdp.mydic = collections.OrderedDict()

for dic_key in cdp.myspin.dic_key_list:
    R1_rho_prime = cdp.myspin.r2[cdp.r2keys[0]]
    R1 = cdp.myspin.ri_data['R1']
    R1_rho = cdp.myspin.r2eff[dic_key]

    # Store the R1_rho value
    y_R1_rho.append(R1_rho)

    exp_type, frq_1e6, offset, point = dic_key.split("_")

    # Get disp_point, the Spin-lock field strength
    x_disp_point.append(float(point))

    # Get w_eff
    w_eff = cdp.myspin.w_eff_spin_dic[dic_key]
    x_w_eff.append(w_eff)

    # Get theta
    theta = cdp.myspin.theta_spin_dic[dic_key]
    x_theta.append(theta)
    
    # Get Domega
    Domega = Domega_spin_dic[dic_key]

    # Calculate y value
    R1_rho_R2eff = (R1_rho - R1*cos(theta)*cos(theta)) / (sin(theta) * sin(theta))
    y.append(R1_rho_R2eff)

    # Add to dic
    if exp_type not in cdp.mydic:
        cdp.mydic[exp_type] = collections.OrderedDict()
    if frq_1e6 not in cdp.mydic[exp_type]:
        cdp.mydic[exp_type][frq_1e6] = collections.OrderedDict()
    if offset not in cdp.mydic[exp_type][frq_1e6]:
        cdp.mydic[exp_type][frq_1e6][offset] = collections.OrderedDict()
        cdp.mydic[exp_type][frq_1e6][offset]['point'] = []
        cdp.mydic[exp_type][frq_1e6][offset]['theta'] = []
        cdp.mydic[exp_type][frq_1e6][offset]['w_eff'] = []
        cdp.mydic[exp_type][frq_1e6][offset]['R1_rho'] = []
        cdp.mydic[exp_type][frq_1e6][offset]['R1_rho_R2eff'] = []
    cdp.mydic[exp_type][frq_1e6][offset]['point'].append(float(point))
    cdp.mydic[exp_type][frq_1e6][offset]['theta'].append(float(theta))
    cdp.mydic[exp_type][frq_1e6][offset]['w_eff'].append(float(w_eff))
    cdp.mydic[exp_type][frq_1e6][offset]['R1_rho'].append(float(R1_rho))
    cdp.mydic[exp_type][frq_1e6][offset]['R1_rho_R2eff'].append(float(R1_rho_R2eff))
    print cdp.mydic[exp_type][frq_1e6][offset]['R1_rho']

# Print values in dic
for exptype, frq_dic in cdp.mydic.items():
    for frq, offset_dic in frq_dic.items():
        for offset, val_dics in offset_dic.items():
            print exptype, frq, offset, val_dics['R1_rho'], val_dics['theta']

# Modify data
w_eff_div = 10**4
rem_points = 2
x_w_eff_mod = [x/w_eff_div for x in x_w_eff[:-rem_points]]
y_mod = y[:-rem_points]

# Define labels for plotting
plotlabel_R1_rho_R2eff = 'R1_rho_R2eff'
plotlabel_R1_rho = 'R1_rho'
xlabel_w_eff = 'Effective field in rotating frame [%s rad.s^-1]'%(str(w_eff_div))
xlabel_theta = 'Rotating frame tilt angle [rad]'
xlabel_lock = 'Spin-lock field strength [Hz]'
ylabel_R1_rho_R2eff = 'R1_rho_R2eff [rad.s^-1]'
ylabel_R1_rho = 'R1_rho [rad.s^-1]'

# Plot R1_rho_R2eff as function of w_eff
plt.figure()
plt.plot(x_w_eff_mod, y_mod, 'o', label=plotlabel_R1_rho_R2eff)
plt.xlabel(xlabel_w_eff)
plt.ylabel(ylabel_R1_rho_R2eff)
plt.legend(loc='best')
plt.grid(True)
plt.ylim([0,16])
plt.title("%s \n %s as function of %s"%(spin_inte, ylabel_R1_rho_R2eff, xlabel_w_eff))
plt.savefig("matplotlib_%s_%s_w_eff.png"%(spin_inte_rep, plotlabel_R1_rho_R2eff) )

# Plot R1_rho_R2eff as function of theta
plt.figure()
plt.plot(x_theta, y, 'o', label=plotlabel_R1_rho_R2eff)
plt.xlabel(xlabel_theta)
plt.ylabel(ylabel_R1_rho_R2eff)
plt.legend(loc='best')
plt.grid(True)
plt.ylim([0,16])
plt.title("%s \n %s as function of %s"%(spin_inte, ylabel_R1_rho_R2eff, xlabel_theta))
plt.savefig("matplotlib_%s_%s_theta.png"%(spin_inte_rep, plotlabel_R1_rho_R2eff) )

# Plot R1_rho_R2eff as function of disp_point, the Spin-lock field strength
plt.figure()
plt.plot(x_disp_point, y, 'o', label=plotlabel_R1_rho_R2eff)
plt.xlabel(xlabel_lock)
plt.ylabel(ylabel_R1_rho_R2eff)
plt.legend(loc='best')
plt.grid(True)
plt.ylim([0,16])
plt.title("%s \n %s as function of %s"%(spin_inte, ylabel_R1_rho_R2eff, xlabel_lock))
plt.savefig("matplotlib_%s_%s_disp.png"%(spin_inte_rep, plotlabel_R1_rho_R2eff) )

# Plot R1_rho as function of theta.
plt.figure()
plt.plot(x_theta, y_R1_rho, 'o', label=plotlabel_R1_rho)
plt.xlabel(xlabel_theta)
plt.ylabel(ylabel_R1_rho)
plt.legend(loc='best')
plt.grid(True)
plt.ylim([0,16])
plt.title("%s \n %s as function of %s"%(spin_inte, ylabel_R1_rho, xlabel_theta))
plt.savefig("matplotlib_%s_%s_theta.png"%(spin_inte_rep, plotlabel_R1_rho) )

# Plot R1_rho as function of theta.
plt.figure()
for exptype, frq_dic in cdp.mydic.items():
    for frq, offset_dic in frq_dic.items():
        for offset, val_dics in offset_dic.items():
            graphlabel = "%s_%s"%(frq, offset)
            plt.plot(val_dics['theta'], val_dics['R1_rho'], '-o', label=graphlabel)
plt.ylabel(ylabel_R1_rho)
plt.xlabel(xlabel_theta)
plt.legend(loc='best')
plt.grid(True)
plt.ylim([0,16])
plt.title("%s \n %s as function of %s"%(spin_inte, ylabel_R1_rho, xlabel_theta))
plt.savefig("matplotlib_%s_%s_theta_sep.png"%(spin_inte_rep, plotlabel_R1_rho) )

plt.show()

To run

relax -p r1rhor2eff.py

See also