User:Remig/plico/plicoCommon
< User:Remig | plico
Jump to navigation
Jump to search
This script contains routines used by other scripts of the Plico suite. 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 plicoCommon.spt.
# plicoCommon - Jmol script by Ron Mignery
# v1.10 beta 9/11/2015 -use thisModel
#
# Routines and values common to all Plico suite scripts
# Must be present in the same directory as other Plico scripts that use it
kCommon = 6
kDtolerance = 0.1
kAtolerance = 5.0
kCtolerance = 1.85
kMtolerance = 0.8
gMouseX = 0
gMouseY = 0
gMinNo = 1
gMaxNo = 9999
gOK = true # global return value to work around jmol *feature*
gOk2 = true # " "
gScheme = "Jmol"
gAltScheme = "Rasmol"
gBusy = false
gEcho = ""
gMenuMin = false
# Return L tetrahedron point if i1<i2<i3, else R point
function get_tet_idx(i1, i2, i3, dist) {
var v1 = {atomIndex=i3}.xyz - {atomIndex=i2}.xyz
var v2 = {atomIndex=i1}.xyz - {atomIndex=i2}.xyz
var axis = cross(v1, v2)
var pma = ({atomIndex=i1}.xyz + {atomIndex=i3}.xyz)/2
var pmo = {atomIndex=i2}.xyz + {atomIndex=i2}.xyz - pma
var pt = pmo + (axis/axis)
var v = pt - {atomIndex=i2}.xyz
var cdist = distance(pt, {atomIndex=i2})
var factor = (dist/cdist)
var lpt = v * factor
return lpt + {atomIndex=i2}.xyz
}
function get_trigonal_idx(i1, i2, i3, dist) {
var v1 = {atomIndex=i1}.xyz - {atomIndex=i2}.xyz
var v2 = {atomIndex=i3}.xyz - {atomIndex=i2}.xyz
var pt = {atomIndex=i2}.xyz - (v1 + v2)
var v = pt - {atomIndex=i2}.xyz
var cdist = distance(pt, {atomIndex=i2})
var factor = (dist/cdist)
var lpt = (v * factor)
return lpt + {atomIndex=i2}.xyz
}
function abs( x) {
return ((x < 0) ? -x : x)
}
# Selected must include second parameter but not the first
function set_distance_atoms( static, mobile, desired) {
try {
var v = mobile.xyz - static.xyz
var dist = distance(static, mobile)
translateSelected @{((v * (desired/dist)) - v)}
}
catch {
}
}
function set_distance_idx( staticIdx, mobileIdx, desired) {
set_distance_atoms({atomIndex=staticIdx}, {atomIndex=mobileIdx}, desired)
}
# Selected must include third parameter but not the first
function set_angle_atoms( stator, pivot, rotor, toangle) {
try {
var v1=stator.xyz - pivot.xyz
var v2=rotor.xyz - pivot.xyz
var axis = cross(v1, v2) + pivot.xyz
var curangle = angle(stator, pivot, rotor)
rotateselected @axis @pivot @{curangle-toangle}
}
catch {
}
}
function set_angle_idx( statorIdx, pivotIdx, rotorIdx, toangle) {
set_angle_atoms({atomIndex=statorIdx}, {atomIndex=pivotIdx},
{atomIndex=rotorIdx}, toangle)
}
# Selected must include fourth parameter but not the first
function set_dihedral_atoms( stator, pivot1, pivot2, rotor, toangle) {
try {
var curangle = angle(stator, pivot1, pivot2, rotor)
rotateselected {pivot2} {pivot1} @{curangle-toangle}
}
catch {
}
}
function set_dihedral_idx( statorIdx, pivot1idx, pivot2idx, rotorIdx, toangle) {
set_dihedral_atoms({atomIndex = statorIdx}, {atomIndex = pivot1idx},
{atomIndex = pivot2idx}, {atomIndex = rotorIdx}, toangle)
}
function angle_idx_4( a1idx, a2idx, a3idx, a4idx) {
return angle({atomIndex=a1idx}, {atomIndex=a2idx},
{atomIndex=a3idx}, {atomIndex=a4idx})
}
# First and last are BB atoms
# Any side atoms in the range are also selected
function select_nward_idx (firstIdx, lastIdx) {
var firstno = ((firstIdx < 0) ? {atomIndex=lastIdx}.atomno : {atomIndex=firstIdx}.atomno)
var lastno = ((lastIdx < 0) ? firstno : {atomIndex=lastIdx}.atomno)
var iChain = ((firstIdx < 0)
? {atomIndex=lastIdx}.chain : {atomIndex=firstIdx}.chain)
select ((atomno <= firstno) and (atomno >= lastno) and (chain = iChain)
and thisModel)
if ({(atomno=firstno) and (chain=gChain)
and thisModel}.atomName == "C") { # if psi
add_sc_to_select(firstno-1, true, true, iChain)
var a = {(atomno=@{firstno+1}) and (chain=iChain) and thisModel}
a.selected = true # add O
}
if ({(atomno=firstno) and (chain=iChain)
and thisModel}.atomName == "CA") {
add_sc_to_select(firstno, true, false, iChain)
}
if ({(atomno=lastno) and (chain=iChain)
and thisModel}.atomName == "C") { # if psi
add_sc_to_select(lastno-1, false, false, iChain)
}
}
function add_sc_to_select(CAno, isAdd, addOXT, iChain) {
var res = {(atomno=CaNo) and (chain=iChain) and thisModel}.resno
var scset = {(resno=res) and (chain=iChain) and thisModel and sidechain}
if (addOXT) {
scset = scset or get_atom_rcn(res, iChain, "OXT")
}
if (isAdd) {
select add scset
}
else {
select remove scset
}
}
# First and last are BB atoms
# Any side atoms in the range are also selected
function select_cward_idx (firstIdx, lastIdx) {
var firstno = ((firstIdx < 0) ? gMaxNo : {atomIndex=firstIdx}.atomno)
var lastno = ((lastIdx < 0) ? 1 : {atomIndex=lastIdx}.atomno)
var iChain = ((firstIdx < 0)
? {atomIndex=lastIdx}.chain : {atomIndex=firstIdx}.chain)
# If nWard anchor in range, begin selection with it
if ((gNanchorIdx >= 0) and ({atomIndex=gNanchorIdx}.chain == iChain)) {
var aNo = {atomIndex=gNanchorIdx}.atomno
if (aNo > firstNo) {
firstno = aNo
}
}
# If cWard anchor in range, end selection with it
if ((gCanchorIdx >= 0) and ({atomIndex=gCanchorIdx}.chain == iChain)) {
var aNo = {atomIndex=gCanchorIdx}.atomno
if (aNo < lastNo) {
lastno = aNo
}
}
select ((atomno >= firstno) and (atomno <= lastno) and (chain = iChain)
and thisModel)
if ({(atomno=firstno) and (chain=iChain)
and thisModel}.atomName == "C") { # if psi
add_sc_to_select(firstno-1, false, false, iChain)
}
if ({(atomno=lastno) and (chain=iChain)
and thisModel}.atomName == "CA") {
add_sc_to_select(lastno, true, false, iChain)
}
if ({(atomno=lastno) and (chain=iChain)
and thisModel}.atomName == "C") { # if psi
add_sc_to_select(lastno-1, true, true, iChain)
var a = {(atomno=@{lastno+1}) and (chain=iChain) and thisModel}
a.selected = true # add O
}
}
# tug.spt needed to handle collisions
function to_handle_collisions( idx) {
#(Re)load tug functions
script $SCRIPT_PATH$tug.spt
if (kTug < 3) {
prompt ("A newer version of tug.SPT is required")
quit
}
handle_collisions( idx)
}
function count_collisions(rc) {
var cAtoms = ({})
for (var idx = {thisModel}.min.atomIndex;
idx <= {thisModel}.max.atomIndex; idx++) {
if ({atomIndex=idx} and {element!="H"}) {
var lcAtoms = (within(kCtolerance, false, {atomIndex=idx})
and {atomIndex > idx} and {element!="H"}
and not connected({atomIndex=idx}))
if (lcAtoms) {
cAtoms = cAtoms or lcAtoms or {atomIndex=idx}
if (rc == true) {
print format("Collision of idx=%d with %s", idx, lcAtoms)
for (var i = 1; i <= lcAtoms.size; i++) {
measure {atomindex=idx} ({@lcAtoms[i]})
}
}
}
}
}
return cAtoms
}
# Utility
function get_atom_rcn( iResno, iChain, iName) {
return {(resno=iResno) and (chain=iChain) and (atomName=iName)
and thisModel}
}
# Utility
function atom_noc( iNo, iChain) {
return {(atomno=iNo) and (chain=iChain) and thisModel}
}
# A handy debug routine
function hi {
hover "%D %U"
}
function plico_menu_toggle() {
var yr = _height-25
if ((_mouseY > (_height-25)) and (_mouseX < 300)) {
if (gMenuMin) {
echo @gEcho
}
else {
var p = gEcho.find("|")
if (p > 0) {
echo @{gEcho[1][p-1]}
}
else {
echo
}
}
gMenuMin = !gMenuMin
}
}
function plico_prelim(repair, backup) {
# Push selected
gSelSaves = {selected}
select (thisModel)
# Bad idea to proceed when collisions present
while (repair) {
var cc = count_collisions(({})).size
if (cc > 0) {
var p = prompt(format("%d collision%s present!\nProceed anyway?",
cc, ((cc > 1) ? "s" : "")), "OK|Cancel|Repair", true)
if (p == "Cancel") {
quit
}
else if (p == "Repair") {
select (thisModel)
allSet = {selected}
gChain = "XX"
for (var idx = {allSet}.atomIndex.min;
idx <= {allSet}.atomIndex.max; idx++) {
if ({atomIndex=idx}.chain != gChain) {
gChain = {atomIndex=idx}.chain
select {(chain=gChain) and thisModel}
to_handle_collisions( idx)
}
}
}
else {
break
}
}
else {
break
}
} # endwhile
gZoom = script("show zoom")
gRotate = script("show rotation")
if (backup) {
write plico_save.pdb
}
select none
gScheme = defaultColorScheme
gAltScheme = ((gScheme == "jmol") ? "rasmol" : "jmol")
set echo TOP LEFT
background ECHO yellow
gChain = ""
gMenuMin = false
unbind
}
# gPlicoRecord is maintained by the macro plicoRecord
function plico_record(s) {
var g = format("show file \"%s\"", gPlicoRecord)
var ls = script(g)
if (ls.find("FileNotFoundException") > 0) {
ls = ""
}
ls += s
write var ls @gPlicoRecord
}
function plico_exit(undo) {
var p = ""
var done = false
if (gPlico) {
p = prompt(format("Exit %s?", gPlico),
"Yes|No"+(undo ? "|Undo all" : ""), true)
if (p == "Undo all") {
load plico_save.pdb
script inline gZoom
rotate @gRotate
echo Session undone
if (gPlicoRecord != "") {
plico_record("load plico_save.pdb;")
}
reset gPlotEcho
set AnimFrameCallback NONE
}
}
if (p != "No") {
unbind all
halo off
select (thisModel)
halo off
star off
color {selected} @gScheme
draw gSCcircle DELETE
gBusy = false
gEcho = ""
set echo TOP LEFT
background ECHO yellow
echo @gEcho
gMenuMin = false
done = true
# Pop selected
select gSelSaves
}
return done
}
# Prevents minimization from disrespecting planar elements
# and allows calculate hydrogens to place hydrogens correctly
function double_bond_planars(iChain, undo) {
print format("double_bond_planars(%s, %s)", iChain, undo)
for (var i = get_resno_min(iChain); i <= get_resno_max(iChain); i++) {
var aa = {resno=i}.group[0]
var bps = ["O","C"]
switch (aa) {
case 'ARG' :
bps += ["CZ","NH1"]
break
case 'ASN' :
case 'ASP' :
bps += ["OD1","CG"]
break
case 'GLN' :
case 'GLU' :
bps += ["OE1","CD"]
break
case 'HIS' :
bps += ["CG","CD2","CE1","NE2"]
break
case 'PHE' :
case 'TYR' :
bps += ["CG","CD1","CD2","CE2","CE1","CZ"]
break
case 'TRP' :
bps += ["CG","CD1","CD2","CE2","CE3","CZ3","CZ2","CH2"]
break
case 'A' :
case 'G' :
case 'DA' :
case 'DG' :
bps = ["N7","C8","C4","C5","N3","C2", ((aa = "A") ? "N1" : "O6"), "C6"]
break
case 'C' :
case 'U' :
case 'T' :
case 'DC' :
case 'DT' :
case 'DU' :
bps = ["C5","C6","C2","O2","C4","CZ",((aa = "C") ? "N3" : "O4")]
break
}
for (var j = 1; j <= bps.size; j += 2) {
var a1 = get_atom_rcn( i, iChain, bps[j])
var a2 = get_atom_rcn( i, iChain, bps[j+1])
if (x1) throw context t
if (undo) {
connect @a1 @a2 single
}
else {
connect @a1 @a2 double
}
}
var a1 = get_atom_rcn(get_resno_max(iChain), iChain, "O")
var a2 = get_atom_rcn(get_resno_max(iChain), iChain, "C")
if (undo) {
connect @a1 @a2 single
}
else {
connect @a1 @a2 double
}
}
}
function get_resno_max( iChain) {
return {not hoh and not ligand and (chain=iChain) and thisModel}.resno.max
}
function get_resno_min( iChain) {
return {not hoh and not ligand and (chain=iChain) and thisModel}.resno.min
}
function plico_minimize( aset) {
try {
var savemt = useMinimizationThread
set useMinimizationThread false
var ms = minimizationSilent
set minimizationSilent true
minimize @aset
set useMinimizationThread savemt
set minimizationSilent ms
}
catch {
print "minimization error"
}
}
# end of plicocommon.spt
