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.4 beta 7/8/2014 -handle multiple frames
#
# 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 it
kNTcommon = 4
kC5O5PO3B = -71.0
kO5PO3C3B = -107.0
kPO3C3C4B = -161.5
kO3C3C4C5B = 125.5
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) {
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)
minimize select {connected(aP) or aP}
}
}
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)) {
# 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) {
minimize select {connected(aP) or aP}
}
}
}
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)
}
# Stack res rMove on res rFixed
function base_stack_res( rMove, rFixed, iChain, jChain, sep , ang) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
var isA = is_form_a(rMove, iChain)
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, m)
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, m)
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 ? 85.7 : 77.8)
# Set dihedral Njx Nix C5i O3i (1tna=179.9)
vs[3] = (isA ? 179.9 : 173.3)
# Set angle C5j Njx Nxif (1tna=112.9)
vs[4] = (isA ? 112.9 : 124.2)
# Set dihedral O3j C5j Njx Nix (1tna= -20)
vs[5] = (isA ? -20 : -14.2)
# Set dihedral of C5j Njx Nix C5i (1tna=20)
vs[6] = ang
move_it(as, vs)
#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) {
to_ab_nt_res(res3, res3-1, jChain, FALSE)
force_p_res(res3, jChain)
as = gen_as(res3, ares, jChain, jChain)
select {(resno=res3) and (chain=iChain) and (file=f) and (model=m)}
}
else {
to_ab_nt_res(res5, res5-1, iChain, FALSE)
force_p_res(res5, iChain)
as = gen_as(res5, res5-1, iChain, iChain)
select {(resno=res5) and (chain=iChain) and (file=f) and (model=m)}
}
# 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] = 9.00
# Set angle of iN from jN and jC1 (1ana=124.6)
vs[2] = 124.6
# Set dihedral of iN from jN and jC1 and jC4 (1ana=160.0)
vs[3] = 160.0
# Set angle of iC1 from iN nad jN (1ana=124.6)
vs[4] = 124.6
# Set dihedral of iC4 from iN and iC1 and jN (1ana=160.0)
vs[5] = 160.0
# Set dihedral of iN from iC1 and jN and jC1 (1ana=-5.0)
vs[6] = -5.0
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 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 (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 : kO3C3C4C5A))
# 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)
}
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--) {
if (toab.size > 0) {
to_ab_nt_res(i, -1, iChain, (toab == "A"))
if ((w[i])[1] >= 0) {
to_ab_nt_res((w[i])[1], -1, (w[i])[2], (toab == "A"))
}
}
if (i < res3) {
base_stack_res(i, i+1, iChain, iChain, s, a)
}
}
# Restore pairings
for (var i = res3; i >= res5; i--) {
if ((w[i])[1] >= 0) {
pair_it_res((w[i])[1], i, -1, (w[i])[2], iChain)
}
}
# Clean up
for (var i = res3; i >= res5; i--) {
fix_p_res(i, iChain, TRUE)
if ((w[i])[1] >= 0) {
fix_p_res((w[i])[1], (w[i])[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")
minimize select {connected(aP) or aP}
}
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, m, 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, m, vs) # <== external call
}
function make_tetra_loop(res5, ares, iChain, m, vs) {
var f = (_frameID/1000000)
var m = (_frameID%1000000)
for (var i = (res5+4); i >= ires5; 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
}
}
# end of plicoNTcommon.spt
