In [1]:
%qtconsole
In [1]:
#Initialize
import os
import pandas as pd
%pylab inline
In [4]:
#change to working dir
work_dir = "/var/run/media/chen/Simulation/Simulation/Nucleosome/relex_trr_pre"
os.chdir(work_dir)
Render trajectory by different colors according to force value loaded on residue.
After done the partial force calculation, there are three output file.
</font>
0.1: Steps:
0. vmd #start vmd
1. mol load gro *.gro xtc *.xtc #loaded *.gro and correspond *.xtc into vmd
2.1. package require partial_force_render
2.2. partial_force_render "*.csv" #use .csv file to render the whole structure
* for modified version:
2.1.1 package require partial_force_render_threshold
2.1.2 partial_force_render_threshold "*.csv" threshold #only load value > threshold0.2: Display:
Color map correspond force value, together with a scale bar.0.3: Note:
*.csv:
* File is the output from partial_force calculatiion program.
* The first line is the head (time, residues)
* Start from the second line, is the force value use to render trajcetories.
* Notice the script only checks the number of residue corrspond to main molecule.
* User should take the responsibility to load the right file from certain calculation.0.4: Major Patamers:
* color scale method "RGB"
* set color_max 1000
* set color_min 0.0
In [0]:
In [0]:
In [0]:
In [6]:
#residue range:
H3_1 = arange(1,136)
H3_2 = arange(488, 623)
H4_1 = arange(136,238)
H4_2 = arange(623,725)
H2A_1 = arange(238,366)
H2A_2 = arange(725,853)
H2B_1 = arange(366,488)
H2B_2 = arange(853,975)
DNA = arange(975,1307)
DNA_1 = arange(975,1141)
DNA_2 = arange(1141,1307)
#print DNA_2
In [7]:
#
fm = pd.read_csv('force_external_md1000_pbc_mol_skip_10.csv')
cl=fm.columns
In [48]:
print cl[63], cl[83], cl[0], cl[1]
print cl[488 : 623]
print cl[DNA]
In [14]:
dt_63 = fm[cl[63]]
dt_83 = fm[cl[83]]
In [17]:
plot(dt_63)
plot(dt_83)
legend([cl[63],cl[83]])
Out[17]:
In [24]:
dt_570 = fm[cl[570]]
dt_550 = fm[cl[550]]
In [25]:
plot(dt_550)
plot(dt_570)
legend([cl[550],cl[570]])
Out[25]:
In [19]:
time_label = cl[0]
target = cl[1293]
plot(fm[time_label], fm[target]*1.66/1000)
legend(["force loaded on residue :" + target])
xlabel("time[ps]")
ylabel("force[nN]")
Out[19]:
In [18]:
time_label = cl[0]
target = cl[1292]
plot(fm[time_label], fm[target]*1.66/1000)
legend([target])
xlabel("time[ps]")
ylabel("force[nN]")
Out[18]:
In [ ]: