User:Remig/plico/plicoCommonNT
< User:Remig | plico
Jump to navigation
Jump to search
This script contains routines used by some other scripts of the Plico suite involved with polynucleotide manipulation. It must be located in the same directory as any script that uses these routines.
Copy and paste the following into a text editor and save in your scripts folder as plicoNTcommon.spt.
# plicoNTcommon - Jmol script by Ron Mignery
# v1.5 beta 8/18/2014 -better A-B conversion
#
# Routines and values common to Plico suite scripts that work with nucleotides
# Must be present in the same directory as other Plico scripts that use them
kNTcommon = 4
kC5O5PO3B = -71.0
kO5PO3C3B = -107.0
kPO3C3C4B = -161.5
kO3C3C4C5B = 140.0
kC3C4C5O5B = 55.65
kC4C5O5PB = 169.0
kO4C4C3C2B = 15.92
kC4O4C1C2B = -41.7
kC4O4C1NxB = -159.0
kC5C4O4C1B = 146.3
kC3C1C2O2B = 120.5
kPuB = 59.0
kPyB = 61.0
kC5O5PO3A = -59.3
kO5PO3C3A = -63.1
kPO3C3C4A = -157.4
kO3C3C4C5A = 75.5
kC3C4C5O5A = 49.55
kC4C5O5PA = 169.2
kO4C4C3C2A = -35.55
kC4O4C1C2A = 3.8
kC4O4C1NxA = -117.4
kC5C4O4C1A = 144.85
kC3C1C2O2A = 116.3
kPuA = 13.5
kPyA = 16.5
gChain1 = "A"
gChain2 = ""
# Select before calling
function force_p_res(cres, iChain) {
var pres = cres-1
var aP = atom_rcn( cres, iChain, "P")
var aO5 = atom_rcn( cres, iChain, "O5\'")
var aC5 = atom_rcn( cres, iChain, "C5\'")
var aC4 = atom_rcn( cres, iChain, "C4\'")
var aOP1 = atom_rcn( cres, iChain, "OP1")
var aOP2 = atom_rcn( cres, iChain, "OP2")
var aO3p = atom_rcn( pres, iChain, "O3\'")
var aC3p = atom_rcn( pres, iChain, "C3\'")
if (aO3p.size > 0) {
var selsave = {selected}
set_distance_atoms(aC5, aO3p, 3.1)
select aO5
var dist = distance(aO3p, aO5)
var widen = (dist < 2.85)
var dir = (widen ? -1 : 1)
var first = true
while (abs(dist-2.85) > kDtolerance) {
rotateSelected @aC4 @aC5 @dir
var newdist = distance(aO3p, aO5)
if (widen ? (newdist < dist) : (newdist > dist)) {
if (first) {
dir = -dir
rotateSelected @aC4 @aC5 @dir
}
else {
break
}
}
dist=newdist
first = false
}
select aP
set_distance_atoms(aO5, aP, 1.73)
set_angle_atoms(aC5, aO5, aP, 110.1)
#set_dihedral_atoms(aC4, aC5, aO5, aP, 150.3)
aOP2.xyz = get_tet_idx(aO3p.atomIndex, aP.atomIndex, aO5.atomIndex, 1.73)
aOP1.xyz = get_tet_idx(aO5.atomIndex, aP.atomIndex, aO3p.atomIndex, 1.73)
var ms = minimizationSilent
set minimizationSilent true
minimize select {(connected(aP) or aP) and not aO3p}
set minimizationSilent ms
select selsave
}
}
function fix_p_res(cres, iChain, force) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var pres = cres-1
var aP = atom_rcn( cres, iChain, "P")
var aO5 = atom_rcn( cres, iChain, "O5\'")
var aC5 = atom_rcn( cres, iChain, "C5\'")
var aC4 = atom_rcn( cres, iChain, "C4\'")
var aOP1 = atom_rcn( cres, iChain, "OP1")
var aOP2 = atom_rcn( cres, iChain, "OP2")
var aO3p = atom_rcn( pres, iChain, "O3\'")
var aC3p = atom_rcn( pres, iChain, "C3\'")
if ((aO3p.size > 0)and (aC4.size > 0)) {
var selsave = {selected}
# If collision
if (force and distance(aC3p, aC5) <= kCtolerance) {
# Push away
select {(resno <= @{aC5.resno}) and (chain=iChain)
and (file=f) and (model=m)}
set_distance_atoms(aC3p, aC5, kCtolerance)
}
select aO5
var dist = distance(aO3p, aO5)
var widen = (dist < 2.85)
var dir = (widen ? -1 : 1)
var first = true
while (abs(dist-2.85) > kDtolerance) {
rotateSelected @aC4 @aC5 @dir
var newdist = distance(aO3p, aO5)
if (widen ? (newdist < dist) : (newdist > dist)) {
if (first) {
dir = -dir
rotateSelected @aC4 @aC5 @dir
}
else {
break
}
}
dist=newdist
first = false
}
if (force and (abs(distance(aO3p, aC5)-4.11) < kDtolerance)) {
# Push away
select {(resno <= @{aC5.resno}) and (chain=iChain)
and (file=f) and (model=m)}
set_distance_atoms(aO3p, aC5, 4.11)
}
select aP
set_distance_atoms(aO5, aP, 1.73)
set_angle_atoms(aC5, aO5, aP, 110.1)
dist = distance(aO3p, aP)
widen = (dist < 1.73)
dir = (widen ? -1 : 1)
first = true
while (abs(dist-1.73) > kDtolerance) {
rotateSelected @aC5 @aO5 @dir
var newdist = distance(aO3p, aP)
if (widen ? (newdist < dist) : (newdist > dist)) {
if (first) {
dir = -dir
rotateSelected @aC5 @aO5 @dir
}
else {
break
}
}
dist=newdist
first = false
}
#set_dihedral_atoms(aC4, aC5, aO5, aP, 150.3)
aOP2.xyz = get_tet_idx(aO3p.atomIndex, aP.atomIndex, aO5.atomIndex, 1.73)
aOP1.xyz = get_tet_idx(aO5.atomIndex, aP.atomIndex, aO3p.atomIndex, 1.73)
if (force) {
var ms = minimizationSilent
set minimizationSilent true
minimize select {(connected(aP) or aP) and not aO3p}
set minimizationSilent ms
}
select selsave
}
}
function get_rotors_res(res) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var rotors = array()
var sRes = res
var mRes = sRes-1
var iChain = {(resno=res) and (atomName="P")
and (file=f) and (model=m)}.chain
var mC4 = atom_rcn( mRes, iChain, "C4\'")
var mC3 = atom_rcn( mRes, iChain, "C3\'")
var mO3 = atom_rcn( mRes, iChain, "O3\'")
var sP = atom_rcn( sRes, iChain, "P" )
var sO5 = atom_rcn( sRes, iChain, "O5\'")
var sC5 = atom_rcn( sRes, iChain, "C5\'")
var sC4 = atom_rcn( sRes, iChain, "C4\'")
var sC3 = atom_rcn( sRes, iChain, "C3\'")
rotors += [mC4.atomIndex, mC3.atomIndex, mO3.atomIndex, sP.atomIndex]
rotors += [mC3.atomIndex, mO3.atomIndex, sP.atomIndex, sO5.atomIndex]
rotors += [mO3.atomIndex, sP.atomIndex, sO5.atomIndex, sC5.atomIndex]
rotors += [sP.atomIndex, sO5.atomIndex, sC5.atomIndex, sC4.atomIndex]
rotors += [sO5.atomIndex, sC5.atomIndex, sC4.atomIndex, sC3.atomIndex]
return rotors
}
function get_nt_chi_rotor_res(res, iChain) {
var rotors = array()
var aO4 = atom_rcn( res, iChain, "O4\'")
var aC1 = atom_rcn( res, iChain, "C1\'")
var isR = ((aC1 and {purine}).size > 0)
var N1or9 = (isR ? "N9" : "N1")
var C6or8 = (isR ? "C8" : "C6")
var aN = atom_rcn(res, iChain, N1or9)
var aC = atom_rcn(res, iChain, C6or8)
rotors = [aO4.atomIndex, aC1.atomIndex, aN.atomIndex, aC.atomIndex]
return rotors
}
function get_nt_ab_rotor_res(res, iChain) {
var rotors = array()
var aC5 = atom_rcn(res, iChain, "C5\'")
var aC4 = atom_rcn(res, iChain, "C4\'")
var aC3 = atom_rcn(res, iChain, "C3\'")
var aO3 = atom_rcn(res, iChain, "O3\'")
rotors = [aO3.atomIndex, aC3.atomIndex, aC4.atomIndex, aC5.atomIndex]
return rotors
}
function gen_nt_rotors(res5, res3, iChain) {
var rotors = array()
for (var i = res5+1; i <= res3; i++) {
rotors += get_rotors_res(i, iChain)
}
return rotors
}
# Moved object must be selected, fixed object not
# as[6] = fixed[1-3] moved[4-6]
# vs[6] = [distance(as[3-4]), angle(as[2-4]),
# dihedral(as[1-4]), angle(as[5-3], dihedral(as[6-3],
# dihedral(as[2-5]
function move_it(as, vs) {
# Distance, angle, dihedral positions atom[4] to a point
set_distance_atoms(as[3], as[4], vs[1])
set_angle_atoms(as[2], as[3], as[4], vs[2])
set_dihedral_atoms(as[1], as[2], as[3], as[4], vs[3])
# Angle and dihedral orients atom[4]'s object
set_angle_atoms(as[3], as[4], as[5], vs[4])
set_dihedral_atoms(as[3], as[4], as[5], as[6], vs[5])
# Dihedral sets TBD
set_dihedral_atoms(as[2], as[3], as[4], as[5], vs[6])
}
function gen_as(res5, res3, iChain, jChain) {
var as = array()
as[1] = atom_rcn(res3, jChain, "C4\'")
as[2] = atom_rcn(res3, jChain, "C1\'")
as[3] = connected(as[2]) and {element="N"}
as[5] = atom_rcn(res5, iChain, "C1\'")
as[6] = atom_rcn(res5, iChain, "C4\'")
as[4] = connected(as[5]) and {element="N"}
return as
}
# Pair res5 on res3 moving res <= res3
function pair_it_res(res5, res3, ares, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var isA = is_form_a(res5, iChain)
var vs = array()
vs[1] = 9.00 # distance iN-jN (1tna[11-24]=9.00)
vs[2] = 124.6 # angle iN-jN-jC1 (1tna[11-24]=99.6)//126.0
vs[3] = (isA ? 160.0 : -140.0) # dihedral iN-jN-jC1-jC4 (1tna[11-24]=160.0)//-142.6
vs[4] = 124.6 # angle iC1-iN-jN (1tna[11-24]=124.6)
vs[5] = (isA ? 160.0 : -140.0) # dihedral iC4-iN-iC1-jN (1tna[11-24]=160.0)//-138.6
vs[6] = (isA ? -5.0 : -1.6) # dihedral iN-iC1-jN-jC1 (1tna[11-24]=-5.0)//-1.6
if (ares < 0) {
select ((resno=res5) and (chain=iChain) and (file=f) and (model=m))
}
else {
select ((resno <= ares) and (chain=iChain) and (file=f) and (model=m))
}
move_it(as, vs)
fix_p_res(res5, iChain, true)
fix_p_res(res3, jChain, true)
}
# Unstack res5 on res3 moving just res5
function single_unstack_res5_on_res3(res5, res3, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var vs = array()
vs[1] = 4.31 # distance iN-jN (1tna[38-39]=4.31)
vs[2] = 99.6 # angle iN-jN-jC1 (1tna[38-39]=99.6)
vs[3] = -61.1 # dihedral iN-jN-jC1-jC4 (1tna[38-39]=-61.1)
vs[4] = 128.8 # angle iC1-iN-jN (1tna[38-39]=128.8)
vs[5] = 58.4 # dihedral iC4-iN-iC1-jN (1tna[38-39]=58.4)
vs[6] = 78.4 # dihedral iN-iC1-jN-jC1 (1tna[38-39]=78.4)
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
move_it(as, vs)
force_p_res(res3, jChain)
move_it(as, vs)
#fix_p_res(res5, iChain, false)
force_p_res(res3, jChain)
}
# Flatstack res5 on res3 moving just res5
function single_flatstack_res5_on_res3(res5, res3, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var vs = array()
vs[1] = 7.00 # distance iN-jN
vs[2] = 89.1 # angle iN-jN-jC1
vs[3] = -49.9 # dihedral iN-jN-jC1-jC4
vs[4] = 83.4 # angle iC1-iN-jN
vs[5] = 125.7 # dihedral iC4-iN-iC1-jN
vs[6] = 5.8 # dihedral iN-iC1-jN-jC1
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
move_it(as, vs)
force_p_res(res3, jChain)
move_it(as, vs)
#fix_p_res(res3, jChain, true)
force_p_res(res3, jChain)
}
# Outstack res5 on res3 moving just res5
function single_outstack_res5_on_res3(res5, res3, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var vs = array()
vs[1] = 8.23 # distance iN-jN
vs[2] = 32.4 # angle iN-jN-jC1
vs[3] = -26.8 # dihedral iN-jN-jC1-jC4
vs[4] = 99.6 # angle iC1-iN-jN
vs[5] = 57.4 # dihedral iC4-iN-iC1-jN
vs[6] = 179.1 # dihedral iN-iC1-jN-jC1
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
move_it(as, vs)
force_p_res(res3, jChain)
move_it(as, vs)
fix_p_res(res3, jChain, true)
}
# Flatstack res3 on res5 moving just res5
function single_flatstack_res3_on_res5(res5, res3, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
vs = array()
vs[1] = 6.00 #4.6# distance iN-jN
vs[2] = 90#75.1 # angle iN-jN-jC1
vs[3] = 90#135.3 # dihedral iN-jN-jC1-jC4
vs[4] = 90#89.9 # angle iC1-iN-jN
vs[5] = -90#-47.3 # dihedral iC4-iN-iC1-jN
vs[6] = 0#1.7 # dihedral iN-iC1-jN-jC1
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
move_it(as, vs)
force_p_res(res5, jChain)
}
# Outstack res3 on res5 moving just res5
function single_outstack_res3_on_res5(res5, res3, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var vs = array()
vs[1] = 8.9 # distance iN-jN
vs[2] = 65.3 # angle iN-jN-jC1
vs[3] = 55.7 # dihedral iN-jN-jC1-jC4
vs[4] = 61.2 # angle iC1-iN-jN
vs[5] = -41.2 # dihedral iC4-iN-iC1-jN
vs[6] = -138.4 # dihedral iN-iC1-jN-jC1
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
move_it(as, vs)
force_p_res(res5, jChain)
move_it(as, vs)
fix_p_res(res5, jChain, true)
}
# Pair A res5 on A res3 Hogsteen (N6-N7)2 moving res5 => res3
function pair_it_h_aa(res5, res3, ares, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = array()
var vs = array()
as[1] = atom_rcn(res3, jChain, "N6")
as[2] = atom_rcn(res3, jChain, "C1\'")
as[3] = atom_rcn(res3, jChain, "N7")
as[4] = atom_rcn(res5, iChain, "N6")
as[5] = atom_rcn(res5, iChain, "C1\'")
as[6] = atom_rcn(res5, iChain, "N7")
var cp = as[5].xyz
select {(resno = res5) and (chain=iChain) and (file=f) and (model=m)
and base} or @{as[5]}
# Set distance of iN6 from jN7 (1tna=2.92)
vs[1] = 2.92
# Set angle of iN6 from jN7 and jC1 (1tna=123.1)
vs[2] = 123.1
# Set dihedral of iN6 from jN7 and jC1 and jN6 (1tna= 154.9)
vs[3] = 154.9
# Set angle of iN7 from iN6 and jC1 (1ana=98.2)
vs[4] = 98.2
# Set dihedral of iC4 from iN and iC1 and jN (1tna=18.2)
vs[5] = 18.2
# Set dihedral of iN7 from iN6 and jN7 and jC1 (1tna=177.6)
vs[6] = 177.6
# Move the base and C1' into final position
move_it(as, vs)
# Mark C1' xyz and move it back to its orginal position
var pt = as[5].xyz
as[5].xyz = cp
# Collect available P rotors
var rotors = gen_nt_rotors(res5, ares, iChain)
# Until there (4 tries)
var cnt = 0
while (distance(as[5], pt) > kDtolerance) {
# Rotate on rotor set to move C1' to its new position
move_atom_nt( as[5].atomIndex, pt, 0, rotors)
# Rotate on anchor chi
select {(resno >= res5) and (resno <= ares and (chain=iChain)
and (file=f) and (model=m)) and not base}
rotate_chi_for_distance_atoms(ares, iChain, as[5], pt, kDtolerance)
cnt++
if (cnt > 3) {
break
}
}
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)
and not base}
set_distance_atoms(pt, as[5] kDtolerance)
force_p_res(ares, iChain)
}
# Pair U res5 on A res3 Hogsteen N3-N7, N6-O2, O4-O1p moving res5 => res3
function pair_it_h_ua(res5, res3, ares, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var as = gen_as(res5, res3, iChain, jChain)
var vs = array()
var cp = as[6].xyz
select {(resno = res5) and (chain=iChain) and (file=f) and (model=m)
and base} or @{as[6]}
# Set distance of iN1or3 from jN1or3 (1tna=5.77)
vs[1] = 5.77
# Set angle of iN1or3 from jN1or3 and jN9or1 (1tna=115.7)
vs[2] = 115.7
# Set dihedral of iN1or3 from jN1or3 and jN9or1 and jC1' (1tna= 19.9)
vs[3] = 19.9
# Set angle of iN9or1 from iN1or3 and jN1or3 (1ana=57.3)
vs[4] = 57.3
# Set dihedral of iC1' from iN9or1 and iN1or3 and jN1or3 (1tna=-177.2)
vs[5] = -177.2
# Set dihedral of iN9or1 from iN1or3 and jN1or3 and jN9or1 (1tna=168.1)
vs[6] = 168.1
# Move the base and C1' into final position
move_it(as, vs)
# Mark C1' xyz and move it back to its orginal position
var pt = as[6].xyz
as[6].xyz = cp
# Collect available P rotors
var rotors = gen_nt_rotors(res5, ares, iChain)
# Until there (4 tries)
var cnt = 0
while (distance(as[5], pt) > kDtolerance) {
# Rotate on rotor set to move C1' to its new position
move_atom_nt( as[5].atomIndex, pt, 0, rotors)
# Rotate on anchor chi
select {(resno >= res5) and (resno <= ares) and (chain=iChain)
and (file=f) and (model=m) and not base}
rotate_chi_for_distance_atoms(ares, iChain, as[6], pt, kDtolerance)
cnt++
if (cnt > 3) {
break
}
}
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)
and not base}
set_distance_atoms(pt, as[6] kDtolerance)
force_p_res(ares, iChain)
}
function level_base(rMove, rFixed, iChain, jChain) {
var selsave = {selected}
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var mC1 = atom_rcn(rMove, iChain, "C1\'")
var mIsR = ({mC1 and purine}.size > 0)
var m9or1 = (mIsR ? "N9" : "N1")
var m6or8 = (mIsR ? "C8" : "C6")
var m4or2 = (mIsR ? "C4" : "C2")
var mN = atom_rcn(rMove, iChain, m9or1)
var mC6or8 = atom_rcn(rMove, iChain, m6or8)
var mC4or2 = atom_rcn(rMove, iChain, m4or2)
var fC1 = atom_rcn(rFixed, jChain, "C1\'")
var fIsR = ({fC1 and purine}.size > 0)
var f9or1 = (fIsR ? "N9" : "N1")
var f6or8 = (fIsR ? "C8" : "C6")
var f4or2 = (fIsR ? "C4" : "C2")
var fN = atom_rcn(rFixed, jChain, f9or1)
var fC6or8 = atom_rcn(rFixed, jChain, f6or8)
var fC4or2 = atom_rcn(rFixed, jChain, f4or2)
var dist = abs(distance(fc4or2, mc4or2) - distance(fc6or8, mc6or8))
var newdist = dist
var dir = 0.1
select {(resno=rMove) and (chain=iChain) and base and (file=f) and (model=m)}
while(newdist > 0.01) {
if (newdist > dist) {
if (dir == 0.1) {
dir = -0.1
}
else {
rotateSelected @mC1 @mN @{-dir}
break
}
}
dist = newdist
rotateSelected @mC1 @mN @{dir}
newdist = abs(distance(fc4or2, mc4or2) - distance(fc6or8, mc6or8))
}
select selsave
}
# Stack res rMove on res rFixed
function base_stack_res( rMove, rFixed, iChain, jChain, sep , ang, single) {
var isA = is_form_a(rMove, iChain)
var is3on5 = (rMove > rFixed)
var j = rFixed
var i = rMove
var as = array()
var vs = array()
as[1] = atom_rcn(j, jChain, "O3\'")
as[2] = atom_rcn(j, jChain, "C5\'")
var Njx = (((as[1] and {purine}).size > 0) ? "N9" : "N1")
as[3] = atom_rcn(j, jChain, Njx)
as[5] = atom_rcn(i, iChain, "C5\'")
as[6] = atom_rcn(i, iChain, "O3\'")
var Nix = (((as[5] and {purine}).size > 0) ? "N9" : "N1")
as[4] = atom_rcn(i, iChain, Nix)
var f = as[1].file
var m = as[1].model
if (single) {
select {(resno = i) and (chain=iChain) and (file=f) and (model=m)}
}
else {
if (is3on5) {
select {(resno >= i) and (chain=iChain) and (file=f) and (model=m)}
}
else {
select {(resno <= i) and (chain=iChain) and (file=f) and (model=m)}
}
}
# Set distance of iNx from jNx (1tna=4.2)
vs[1] = sep
# Set angle Njx Nix C5i (1ana=85.7)
vs[2] = (isA ? (is3on5 ? 112.9 : 85.7) : (is3on5 ? 124.2 : 77.8))
# Set dihedral Njx Nix C5i O3i (1tna=179.9)
vs[3] = (isA ? (is3on5 ? -20 : 179.9) : (is3on5 ? -14.2 : 173.3))
# Set angle C5j Njx Nxif (1tna=112.9)
vs[4] = (isA ? (is3on5 ? 85.7 : 112.9) : (is3on5 ? 77.8 : 124.2))
# Set dihedral O3j C5j Njx Nix (1tna= -20)
vs[5] = (isA ? (is3on5 ? 179.9 : -20) : (is3on5 ? 173.3 : -14.2))
# Set dihedral of C5j Njx Nix C5i (1tna=20)
vs[6] = ang
move_it(as, vs)
level_base(rMove, rFixed, iChain, jChain)
#force_p_res(i, iChain)
fix_p_res(i+1, iChain, true)
}
function single_unpair_yy( rMove, rFixed, iChain, jChain) {
}
function single_unpair_ry( rMove, rFixed, iChain, jChain) {
}
function single_unpair_yr( rMove, rFixed, iChain, jChain) {
}
function single_unpair_rr( res5, res3, ares, iChain, jChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
# Push aside the r stacked stemward on ares
var isYR5 = (res3 = (ares-1))
var as = array()
var vs = array()
if (isYR5) {
select {(resno=res3) and (chain=iChain) and (file=f) and (model=m)}
to_ab_nt_res(res3, -1, jChain, false)
force_p_res(res3, jChain)
as = gen_as(res3, ares, jChain, jChain)
}
else {
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
to_ab_nt_res(res5, -1, iChain, false)
force_p_res(res5, iChain)
as = gen_as(res5, res5-1, iChain, iChain)
}
# Set distance of iN from jN (1ana=6.14)
vs[1] = 6.14
# Set angle of iN from jN and jC1 (1ana=102.5)
vs[2] = 102.5
# Set dihedral of iN from jN and jC1 and jC4 (1ana=-66.2)
vs[3] = -66.2
# Set angle of iC1 from iN nad jN (1ana=72.5)
vs[4] = 72.5
# Set dihedral of iC4 from iN and iC1 and jN (1ana=160.0)
vs[5] = 176.5
# Set dihedral of iN from iC1 and jN and jC1 (1ana=-5.0)
vs[6] = -1.3
move_it(as, vs)
if (isYR5) {
force_p_res(res3+1, jChain)
fix_p_res(res3, jChain, false)
}
else {
force_p_res(res5-1, iChain)
fix_p_res(res5, iChain, false)
}
#########
# Stack the other r on ares
var as = gen_as(res5, ares, iChain, jChain)
select {(resno = res5) and (chain=iChain) and (file=f) and (model=m)}
# Set distance of iN from jN (1ana=9.00)
vs[1] = 6.52
# Set angle of iN from jN and jC1 (1ana=124.6)
vs[2] = 94.0
# Set dihedral of iN from jN and jC1 and jC4 (1ana=160.0)
vs[3] = 176.5
# Set angle of iC1 from iN nad jN (1ana=124.6)
vs[4] = 136.3
# Set dihedral of iC4 from iN and iC1 and jN (1ana=160.0)
vs[5] = -174.7
# Set dihedral of iN from iC1 and jN and jC1 (1ana=-5.0)
vs[6] = -39.7
move_it(as, vs)
fix_p_res(res5, iChain, false)
}
# Rotate rotor set to move target atom to its proper place
function move_atom_nt(targetIdx, targetPt, ares, iRotors) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var pt = targetPt
var rotors = iRotors
var targetNo = {atomIndex=targetIdx}.atomno
var targetRes = {atomIndex=targetIdx}.resno
var iChain = {atomIndex=targetIdx}.chain
gOK = false
var dist = distance(pt, {atomIndex=targetIdx}.xyz)
# If target is a C1' atom, collect its base
var tBase = ({})
var i1 = 0
var i2 = 0
var i3 = 0
var i4 = 0
if ({atomIndex=targetIdx}.atomName == "C1\'") {
tBase = {(resno = targetRes) and base}
}
# For idx number of passes
for (var pass1 = 0; pass1 < 20; pass1++) {
var blocked = ({})
for (var pass2 = 0; pass2 < (rotors.size/4); pass2++) {
var v1 = {atomIndex=targetIdx}.xyz - pt
# Find the most orthgonal unused rotor
var imax = 0
var smax = 0.5
for (var ri = 1; ri < rotors.size; ri += 4) {
i2 = rotors[ri+1]
i3 = rotors[ri+2]
i4 = rotors[ri+3]
if ((i2 != targetIdx) and (i3 != targetIdx) and (i4 != targetIdx)) {
if ({blocked and {atomIndex=i2}}.count == 0) {
v2 = {atomIndex=i3}.xyz - {atomIndex=i2}.xyz
s = sin(abs(angle(v1, {0 0 0}, v2)))
if (s > smax) {
smax = s
imax = ri
}
}
}
}
# If no more rotors, break to next full try
if (imax == 0) {
break
}
i1 = rotors[imax+0]
i2 = rotors[imax+1]
i3 = rotors[imax+2]
i4 = rotors[imax+3]
# Get dihedral of rotor with target point
var dt = (angle({atomIndex=targetIdx}, {atomIndex=i2},
{atomIndex=i3}, pt)/(rotors.size/20))
# Select and rotate
if (ares > targetRes) {
select_3ward_atom({atomIndex=i3}, ares, iChain)
res3 = {atomIndex=i4}.resno
}
else {
select_5ward_atom({atomIndex=i3}, ares, iChain)
res3 = {atomIndex=i1}.resno
}
#***************************************************
rotateSelected {atomIndex=i3} {atomIndex=i2} @{-dt}
# If collisions
var res5 = res3-1
var set3 = {(resno=res3) and (atomName!="P") and (atomName!="OP1")
and (file=f) and (model=m)}
var set5 = {(resno=res5) and (atomName!="P") and (atomName!="OP1")
and (file=f) and (model=m)}
if ((set5 and within(kCtolerance, set3)).size > 0) {
# Binary undo until fixed
while ((abs(dt) > kDtolerance)
and ((set5 and within(kCtolerance, set3)).size > 0)) {
dt /= 2.0
rotateSelected {atomIndex=i3} {atomIndex=i2} @{dt}
}
while ((abs(dt) > kDtolerance)
and ((set5 and within(kCtolerance, set3)).size > 0)) {
dt /= 2.0
rotateSelected {atomIndex=i3} {atomIndex=i2} @{-dt}
}
rotateSelected {atomIndex=i3} {atomIndex=i2} @{dt}
}
# If close enough, stop
dist = distance(pt, {atomIndex=targetIdx})
if (dist < kDtolerance) {
gOK = true
gTargetPt = pt
break
}
# Block rotor
blocked |= {atomIndex=i2}
} # endfor num rotors passes
if (gOK) {
break
}
} # endfor 20 passes
}
# If ares < 0 then adjust iRes only
function to_ab_nt_res(res, ares, iChain, toA) {
var selsave = {selected}
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var aO3 = atom_rcn( res, iChain, "O3\'")
var aC3 = atom_rcn( res, iChain, "C3\'")
var aC4 = atom_rcn( res, iChain, "C4\'")
var aC5 = atom_rcn( res, iChain, "C5\'")
var aC1 = atom_rcn( res, iChain, "C1\'")
var aC2 = atom_rcn( res, iChain, "C2\'")
var aO2 = atom_rcn( res, iChain, "O2\'")
var aO4 = atom_rcn( res, iChain, "O4\'")
if (ares < 0) {
select ((resno=res) and (chain=iChain) and (file=f) and (model=m)
and not aO3 and not aC3 and not aC4)
}
else {
select ((resno >= ares) and (resno <= res) and (chain=iChain)
and (file=f) and (model=m) and not aO3 and not aC3 and not aC4)
}
set_dihedral_atoms(aO3, aC3, aC4, aC5, (toA ? kO3C3C4C5A : kO3C3C4C5B))
# Set chi
var aNx = -1
var aCx = -1
var ang = 0.0
select {(resno=res) and (chain=iChain) and (file=f) and (model=m) and base}
if ((aC1 and {purine}).size > 0) {
aNx = atom_rcn( res, iChain, "N9")
aCx = atom_rcn( res, iChain, "C8")
ang = (toA ? kPuA : kPuB)
}
else {
aNx = atom_rcn(res, iChain, "N1")
aCx = atom_rcn(res, iChain, "C6")
ang = (toA ? kPyA : kPyB)
}
set_dihedral_atoms(aO4, aC1, aNx, aCx, ang)
# Set pucker 3' endo or 2' endo
var pSet = {aC1 or aC2 or aO2}
select pSet or {(resno=res) and (chain=iChain)
and (file=f) and (model=m) and base}
set_dihedral_atoms(aC5, aC4, aO4, aC1, (toA ? kC5C4O4C1A : kC5C4O4C1B))
set_dihedral_atoms(aC4, aO4, aC1, aNx, (toA ? kC4O4C1NxA : kC4O4C1NxB))
set_dihedral_atoms(aC4, aO4, aC1, aC2, (toA ? kC4O4C1C2A : kC4O4C1C2B))
if (aO2.size > 0) {
select aO2 or aC2
ang = (toA ? kC3C1C2O2A : kC3C1C2O2B)
set_dihedral_atoms(aC3, aC1, aC2, aO2, (toA ? kC3C1C2O2A : kC3C1C2O2B))
}
set_distance_atoms(aC3, aC2, 1.52)
set_distance_atoms(aC1, aC2, 1.52)
select selsave
}
function adjust_nts(res5, res3, iChain, toab, a, s) {
var savemt = useMinimizationThread
set useMinimizationThread false
# Collect any pairing
var w = array()
for (var i = res5; i <= res3; i++) {
w = w + [who_pairs(i, iChain)]
}
# Twist and turn
for (var i = res3; i >= res5; i--) {
var j = i-res5+1
if (toab.size > 0) {
to_ab_nt_res(i, -1, iChain, (toab == "A"))
if ((w[j])[1] >= 0) {
to_ab_nt_res((w[j])[1], -1, (w[j])[2], (toab == "A"))
}
}
}
for (var i = res5; i < res3; i++) {
base_stack_res(i, i+1, iChain, iChain, s, a, false)
}
# Restore pairings
for (var i = res3; i >= res5; i--) {
var j = i-res5+1
if ((w[j])[1] >= 0) {
pair_it_res((w[j])[1], i, -1, (w[j])[2], iChain)
}
}
# Clean up
for (var i = res3; i >= res5; i--) {
var j = i-res5+1
fix_p_res(i, iChain, true)
if ((w[j])[1] >= 0) {
fix_p_res((w[j])[1], (w[j])[2], true)
}
}
set useMinimizationThread savemt
}
#########################################################
### STAND ALONE GENERAL PURPOSE FUNCTIONS ###
#########################################################
function is_form_a( iResno, iChain) {
var aO4 = atom_rcn( iResno, iChain, "O4\'")
var aC1 = atom_rcn( iResno, iChain, "C1\'")
var aC2 = atom_rcn( iResno, iChain, "C2\'")
var aC3 = atom_rcn( iResno, iChain, "C3\'")
return (angle(aO4, aC1, aC2, aC3) < 0.0)
}
function is_r_res( iResno, iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
return ({(resno=iResno) and (chain=iChain) and (file=f)
and (model=m) and purine}.size > 0)
}
function repair_p_res(res, iChain) {
var aP = atom_rcn( res, iChain, "P")
var ms = minimizationSilent
set minimizationSilent true
minimize select {connected(aP) or aP}
set minimizationSilent ms
}
function who_pairs(iRes, iChain) {
var aC4or6 = atom_rcn( iRes, iChain, "C4")
var aN1or3 = atom_rcn( iRes, iChain, "N1")
if ({aN1or3 and purine}.size = 0) {
aC4or6 = atom_rcn( iRes, iChain, "C6")
aN1or3 = atom_rcn( iRes, iChain, "N3")
}
var near = within(3.1, aN1or3) and {resno!=iRes} and {element="N"}
for (var i = 1; i <= near.size; i++) {
if (angle(near[i], aN1or3, aC4or6) > 150) {
return [near[i].resno, near[i].chain]
}
}
return [-1, aC4or6.chain]
}
function who_stacks(iRes, iChain) {
var ret = array()
var aNear = ((within(4.0, {(resno=iRes) and base}) )
and {base} and {not resno=iRes})
var done = array()
for (var i = 1; i <= aNear.size; i++) {
var jRes = aNear[i].resno
if (not done.find(jRes)) {
var jChain = aNear[i].chain
var as = gen_as(iRes, jRes, iChain, jChain)
var d = distance({(resno=iRes) and base}, {(resno=jRes) and base})
var a1 = angle(as[2], as[3], as[4])
var a2 = angle(as[5], as[4], as[3])
var dh = angle(as[5], as[4], as[3], as[2])
var bset = ((connected(as[3]) and not as[2])
or (connected(as[4]) and not as[5]))
var a3 = angle(bset[1], bset[2], bset[3])
var a4 = angle(bset[2], bset[3], bset[4])
var isStacked = true
# Bases are parallel as sin(a1) = sin(a2) and sin(a3) = sin(a4)
if (abs(sin(a1)-sin(a2)) > 20) {
isStacked = false
}
if (abs(sin(a3)-sin(a4)) > 20) {
isStacked = false
}
# Bases are stacked as d*sin(a1) < 6.0 and d3 = 0.0
if (d*sin(a1) > 6.2) {
isStacked = false
}
if (abs(dh) > 30) {
#isStacked = false
}
if (isStacked) {
ret += aNear[i].resno
}
done += jRes
}
}
return ret
}
function match_nt(mask, nt) {
var ret = false
switch (mask) {
case "A":
case "U":
case "C":
case "G":
ret = (mask = nt)
break
case "N":
ret = true
break
case "Y":
ret = ((nt=="U") or (nt=="C"))
break
case "R":
ret = ((nt=="A") or (nt=="G"))
break
}
return ret
}
# Calls function match_nt above
function select_seqs(seq, iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
select none
for (var i = {(chain=iChain) and (file=f) and (model=m)}.resno.min;
i <= {(chain=iChain) and (file=f) and (model=m)}.resno.max; i++) {
var j = 1
for (; j <= seq.size; j++) {
var nt = {(chain=iChain) and (file=f) and (model=m) and (resno=@{i+j-1})
and (atomName="C1\'")}.group[1]
if (not match_nt(seq[j], nt)) { # <== external call
break
}
}
if (j > seq.size) {
print format("%s at %d (%s-%s-%s)", seq, i,
gSeq[i-1],
gSeq[i][i+seq.size-1],
gSeq[i+seq.size])
var rset = {(resno=i) and (chain=iChain) and (file=f) and (model=m)}
rset.selected = true
}
}
}
# Calls is_form_a
function select_b_form_nts(iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
select none
for (var i = {(chain=iChain) and (file=f) and (model=m)}.resno.min;
i <= {(chain=iChain) and (file=f) and (model=m)}.resno.max; i++) {
if (not is_form_a(i, iChain)) { # <== external call
print format("Res %d is form B", i)
var rset = {(resno=i) and (chain=iChain) and (file=f) and (model=m)}
rset.selected = true
}
}
}
# From 2LU0
function make_uncg_loop(res5, ares, iChain) {
to_ab_nt_res(res5+1, ares, iChain, false)
to_ab_nt_res(res5+2, ares, iChain, false)
var vs = [
[144.8, 28.8, 138.6, 140.9, 149.6]
[-85.1, -148.4, 56.3, -143.6, 135.7]
[-145.9, -25.6, -93.0, -148.1, 49.6]
[-84.8, 145.9, -7.0, -130.0, -176.2]
[-133.5, -78.0, -60.2, 132.6, 99.0]]
make_tetra_loop(res5, ares, iChain, vs) # <== external call
}
# From 2LU0
function make_gnra_loop(res5, ares, iChain) {
to_ab_nt_res(res5+1, ares, iChain, false)
var vs = [
[147.0, 23.9, 144.6, 144.5, 151.0]
[-118.2, -90.5, 109.0, 166.6, 100.1]
[-155.6, -34.4, -167.7, 115.4, 131.2]
[-136.5, -69.9, -86.8, -170.6, 56.0]
[-147.1, -57.5, -76.8, 147.7, 94.7]]
make_tetra_loop(res5, ares, iChain, vs) # <== external call
}
function make_tetra_loop(res5, ares, iChain, vs) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
for (var i = (res5+4); i >= res5; i--) {
var j = i-1
var as = array()
as += atom_rcn( j, iChain, "C4\'")
as += atom_rcn( j, iChain, "C3\'")
as += atom_rcn( j, iChain, "O3\'")
as += atom_rcn( i, iChain, "P")
as += atom_rcn( i, iChain, "O5\'")
as += atom_rcn( i, iChain, "C5\'")
as += atom_rcn( i, iChain, "C4\'")
as += atom_rcn( i, iChain, "C3\'")
as += atom_rcn( i, iChain, "OP1")
as += atom_rcn( i, iChain, "OP2")
for (var k = 5; k > 0; k--) {
pset = ((k>2) ? (connected(@{as[4]}) or @{as[4]}) : ({}))
select {((resno<i) and (resno>ares) and (chain=iChain)
and (file=f) and (model=m))
or pset}
set_dihedral_atoms(as[k+3], as[k+2], as[k+1], as[k],
(vs[i-res5+1])[k])
}
}
}
function select_3ward_atom(ar3, ares, iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var i = ar3.resno
var aP = atom_rcn( i, iChain, "P")
switch(ar3.atomName) {
case "O3\'" :
select {(resno>i) and (resno<ares) and (chain=iChain)
and (file=f) and (model=m)}
break
case "P" :
select {(resno>=i) and (resno<ares) and (chain=iChain)
and (file=f) and (model=m)}
break
case "O5\'" :
case "C5\'" :
case "C4\'" :
select {(resno>=i) and (resno<ares) and (chain=iChain)
and (file=f) and (model=m) and not (connected(aP) or aP)}
break
case "C3\'" :
var aO3 = atom_rcn( i, iChain, "O3\'")
select {((resno>i) and (resno<ares) and (chain=iChain)
and (file=f) and (model=m)) or aO3}
break
}
}
function select_5ward_atom(ar5, ares, iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var i = ar5.resno
var aP = atom_rcn( i, iChain, "P")
switch(ar5.atomName) {
case "O3\'" :
select {(resno<=i) and (resno>ares) and (chain=iChain)
and (file=f) and (model=m)}
break
case "P" :
case "O5\'" :
case "C5\'" :
select {((resno<i) and (resno>ares) and (chain=iChain)
and (file=f) and (model=m)) or (connected(aP) or aP)}
break
case "C4\'" :
var aC5 = atom_rcn( i, iChain, "C5\'")
select {((resno<i) and (resno>ares) and (chain=iChain)
and (file=f) and (model=m)) or (connected(aP) or aP or aC5)}
break
}
}
function get_atom_rcn( iResno, iChain, iName) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
return {(resno=iResno) and (chain=iChain) and (atomName=iName)
and (file=f) and (model=m)}
}
function plot_ab_chi( iChain) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
select (file=f) and (model=m) and (chain=iChain) and nucleic
for (var i = {(chain=iChain) and (file=f) and (model=m)}.resno.min;
i <= {(chain=iChain) and (file=f) and (model=m)}.resno.max; i++) {
var aO4 = get_atom_rcn(i, iChain, "O4\'")
var aC1 = get_atom_rcn(i, iChain, "C1\'")
var isR = ({aC1 and purine}.size > 0)
var a9or1 = (isR ? "N9" : "N1")
var a6or8 = (isR ? "C8" : "C6")
var aN = get_atom_rcn(i, iChain, a9or1)
var aC = get_atom_rcn(i, iChain, a6or8)
var aO3 = get_atom_rcn(i, iChain, "O3\'")
var aC3 = get_atom_rcn(i, iChain, "C3\'")
var aC4 = get_atom_rcn(i, iChain, "C4\'")
var aC5 = get_atom_rcn(i, iChain, "C5\'")
var chi = angle(aO4, aC1, aN, aC)
{resno=i}.vx = chi
var aorb = angle(aO3, aC3, aC4, aC5)
{resno=i}.vy = aorb
}
select {selected and (atomname="N1")}
plot properties vx vy resno
set echo top left
echo "vx = base chi angle vy = a ==> b form"
}
# end of plicoNTcommon.spt
