User:Remig/plico/tug

Tug allows the user to pull or push by mouse actions to move or rotate one part of a polypeptide against the rest by rotation on its psi and phi bonds with collision detection and restriction.

Tug is a member of the Plico suite of protein folding tools described in User:Remig/plico . It may be installed and accessed as a macro with the file:

Title=PLICO Tug
Script=script <path to your scripts folder>/tug.spt;plicotug

saved as plicotug.macro in your .jmol/macros folder as described in Macro.

Copy and paste the following into a text editor and save in your scripts folder as tug.spt.

#   tug - Jmol script by Ron Mignery
#   v1.0 beta    2/4/2014 for Jmol 13.4
#
#   These functions support protein manipulation in Jmol
var kTug = 1
var kDtolerance = 0.2
var kAtolerance = 5.0
var kCtolerance = 1.85
var kMtolerance = 0.8
var kMinRama = 2
var kCollisionLimit = 200
var gCanchorIdx = -1
var gCanchorNo = -1
var gNanchorIdx = -1
var gNanchorNo = -1
var gNanchorPidx = -1
var gCanchorXyz = {0 0 0}
var gNanchorXyz = {0 0 0}
var gNanchorPxyz = {0 0 0}
var gCcargoIdx = -1
var gNcargoIdx = -1
var gCcargoXyz = {0 0 0}
var gNcargoXyz = {0 0 0}
var gCcargoNo = -1
var gNcargoNo = -1
var gDestAtomIdx = -1
var g1pivotIdx = -1
var g2pivotIdx = -1
var gSelSaves = ({}) 
var gCrotors = array()
var gNrotors = array() 
var gMouseX = 0
var gMouseY = 0
var gAc = ({})
var gChain = ""
var gMinNo = 1    
var gMaxNo = 9999    
var gScheme = "Jmol"
var gAltScheme = "Rasmol"
var gCargoAtoms = ({})
var gSeed = 23423.52
var gBusy = FALSE
var gSCidx = -1
var gSCcircle = -1
var gSCpt = {0 0 0}
var gOK = TRUE # global return value to work around jmol *feature*
var gOk2 = TRUE # "    "
var gTargetPt = {0 0 0}
var gNewDrag = FALSE
var gEcho = ""
var gCountdown = 0
var gZoom = ""
var gRotate = ""

# Return L tetrahedron point if i1<i2<i3, else R point
function getTet(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 getTrigonal(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) {
    if (x < 0) {
        x = -x
    }
    return x
}

function getCAmNo (iNo) {
    while ((iNo > 0) and ({(atomno=iNo) and (chain=gChain)}.atomName != "CA")) {
        iNo--
    }
    return iNo
}

function getCAmIdx (idx) {
    var no = {atomIndex=idx and (chain=gChain)}.atomno
    no = getCAmNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getCApNo (iNo) {
    while ((iNo < gMaxNo) and ({(atomno=iNo) and (chain=gChain)}.atomName != "CA")) {
        iNo++
    }
    return iNo
}

function getCpIdx (idx) {
    var no = {atomIndex=idx and (chain=gChain)}.atomno
    no = getCpNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getCpNo (iNo) {
    while ((iNo < gMaxNo) and ({(atomno=iNo) and (chain=gChain)}.atomName != "C")) {
        iNo++
    }
    return iNo
}

function getCApIdx (idx) {
    var no = {atomIndex=idx and (chain=gChain)}.atomno
    no = getCApNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getCmNo (iNo) {
    while ((iNo > 0) and ({(atomno=iNo) and (chain=gChain)}.atomName != "C")) {
        iNo--
    }
    return iNo
}

function getCmIdx (idx) {
    var no = {atomIndex=idx and (chain=gChain)}.atomno
    no = getCmNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getNmNo (iNo) {
    while ((iNo > 0) and ({(atomno=iNo) and (chain=gChain)}.atomName != "N")) {
        iNo--
    }
    return iNo
}

function getNmIdx (idx) {
    var no = {atomIndex=idx and (chain=gChain)}.atomno
    no = getNmNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getNpNo (iNo) {
    while ((iNo < gMaxNo) and ({(atomno=iNo) and (chain=gChain)}.atomName != "N")) {
        iNo++
    }
    return iNo
}

function getNpIdx (idx) {
    var no = {atomIndex=idx}.atomno
    no = getNpNo( no)
    return ({(atomno=no) and (chain=gChain)}.atomIndex)
}

function getCBidx (BBidx) {
    var no = {atomIndex=BBidx}.atomno
    var i = 1
    for (; i < 5; i++) {
        if ({(atomno=@{no+i}) and (chain=gChain)}.atomName == "CB") {
            break
        }
    }
    return {(atomno=@{no+i}) and (chain=gChain)}.atomIndex
}

function getOidx (BBidx) {
    var no = {atomIndex=BBidx}.atomno
    var i = 1
    for (; i < 4; i++) {
        if ({(atomno=@{no+i}) and (chain=gChain)}.atomName == "O") {
            break
        }
    }
    return {(atomno=@{no+i}) and (chain=gChain)}.atomIndex
}

function getNwardBBno (iNo) {
    while ((iNo >= 0) and (
        ({(atomno=iNo) and (chain=gChain)}.atomName != "N")
        and ({(atomno=iNo) and (chain=gChain)}.atomName != "C")
        and ({(atomno=iNo) and (chain=gChain)}.atomName != "CA"))) {
        iNo--
    }
    return iNo
}

function getNwardBBidx (idx) {
    var no = {atomIndex=idx}.atomno - 1
    no = getNwardBBno( no)
    return ((no >= 0) ? ({(atomno=no) and (chain=gChain)}.atomIndex) : -1)
}

function getCwardBBno (iNo) {
    while ((iNo < gMaxNo) and (
        ({(atomno=iNo) and (chain=gChain)}.atomName != "N")
        and ({(atomno=iNo) and (chain=gChain)}.atomName != "C")
        and ({(atomno=iNo) and (chain=gChain)}.atomName != "CA"))) {
        iNo++
    }
    return iNo
}

function getCwardBBidx (idx) {
    var no = {atomIndex=idx}.atomno + 1
    no = getCwardBBno( no)
    return ((no >= 0) ? ({(atomno=no) and (chain=gChain)}.atomIndex) : -1)
}

function getScBBidx (idx) {
    var no = {atomIndex=idx}.atomno
    for (; no > 0; no--) {
        if ({(atomno=no) and (chain=gChain)}.atomName == "CA") {
            break
        }
        else if ({(atomno=no) and (chain=gChain)}.atomName == "C") {
            break
        }
        else if ({(atomno=no) and (chain=gChain)}.atomName == "N") {
            break
        }
        else if ({(atomno=no) and (chain=gChain)}.atomName == "CB") {
            no -= 3
            break
        }
    }
    return {(atomno=no) and (chain=gChain)}.atomIndex
}

function isBBidx(aIdx) {
    var ret = FALSE
    switch({atomIndex=aIdx}.atomName) {
    case "N":
    case "CA":
    case "C":
        ret = TRUE
        break
    }
    return ret
}

function isSCidx(aIdx) {

    var ret = FALSE
    if (isBBidx(aIDx) == FALSE) {
        
        ret = TRUE
        switch({atomIndex=aIdx}.atomName) {
        case "O":
        case "CB":
            ret = FALSE
            break
        }
    }
    return ret
}

function addSideChainToSelection(CAno, isAdd, addOXT) {
    var iNo = CAno+3
    while ({(atomno=iNo) and (chain=gChain)}.resno == {(atomno=CAno) and (chain=gChain)}.resno) {
        {(atomno=iNo) and (chain=gChain)}.selected = isAdd
        if ({(atomno=iNo) and (chain=gChain)}.atomName == "OXT") {
            {(atomno=iNo) and (chain=gChain)}.selected = addOXT
        }
        iNo++    
    }
}

function selectAddSideChain(fromIdx) {
    var iNo = {atomIndex=fromIdx}.atomno
    select none
    while ({(atomno=iNo) and (chain=gChain)}.atomName != "N") {
        {(atomno=iNo) and (chain=gChain)}.selected = TRUE
        iNo++  
        if (iNo > gMaxNo) {
            break
        }  
    }
}

# First and last are BB atoms
# Any side atoms in the range are also selected
function selectNward (firstIdx, lastIdx) {
    var firstno = ((firstIdx < 0) ? {atomIndex=lastIdx}.atomno : {atomIndex=firstIdx}.atomno)
    var lastno = ((lastIdx < 0) ? firstno : {atomIndex=lastIdx}.atomno)
  
    select (atomno <= firstno) and (atomno >= lastno)
    
    if ({(atomno=firstno) and (chain=gChain)}.atomName == "C") { # if psi
        addSideChainToSelection(firstno-1, TRUE, TRUE)
        {(atomno=@{firstno+1}) and (chain=gChain)}.selected = TRUE # add O
    }
    if ({(atomno=firstno) and (chain=gChain)}.atomName == "CA") {
        addSideChainToSelection(firstno, TRUE, FALSE)
    }
    if ({(atomno=lastno) and (chain=gChain)}.atomName == "C") { # if psi
        addSideChainToSelection(lastno-1, FALSE, FALSE)
    }
}

# First and last are BB atoms
# Any side atoms in the range are also selected
function selectCward (firstIdx, lastIdx) {
    var firstno = ((firstIdx < 0) ? gMaxNo : {atomIndex=firstIdx}.atomno)
    var lastno = ((lastIdx < 0) ? 1 : {atomIndex=lastIdx}.atomno)
    
    # If nWard anchor in range, begin selection with it
    if (gNanchorIdx >= 0) {
        var aNo = {atomIndex=gNanchorIdx}.atomno
        if (aNo > firstNo) {
            firstno = aNo
        }
    }
    
    # If cWard anchor in range, end selection with it
    if (gCanchorIdx >= 0) {
        var aNo = {atomIndex=gCanchorIdx}.atomno
        if (aNo < lastNo) {
            lastno = aNo
        }
    }
    
    select (atomno >= firstno) and (atomno <= lastno)
    
    if ({(atomno=firstno) and (chain=gChain)}.atomName == "C") { # if psi
        addSideChainToSelection(firstno-1, FALSE, FALSE)
    }
    if ({(atomno=lastno) and (chain=gChain)}.atomName == "CA") {
        addSideChainToSelection(lastno, TRUE, FALSE)
    }
    if ({(atomno=lastno) and (chain=gChain)}.atomName == "C") { # if psi
        addSideChainToSelection(lastno-1, TRUE, TRUE)
        {(atomno=@{lastno+1}) and (chain=gChain)}.selected = TRUE # add O
    }
}

# Selected must include second parameter but not the first
function setDistanceIdx (staticIdx, mobileIdx, desired) {
    try {
        var v = {atomIndex=mobileIdx}.xyz - {atomIndex=staticIdx}.xyz
        var dist = distance({atomIndex=staticIdx}, {atomIndex=mobileIdx})
        translateSelected @{((v * (desired/dist)) - v)}
    }
    catch {
    }
}


# Selected must include third parameter but not the first
function setAngleIdx (statorIdx, pivotIdx, rotorIdx, toangle) {
    try {
        var v1={atomIndex=statorIdx}.xyz - {atomIndex=pivotIdx}.xyz
        var v2={atomIndex=rotorIdx}.xyz - {atomIndex=pivotIdx}.xyz
        var axis = cross(v1, v2) + {atomIndex=pivotIdx}.xyz
        var curangle =  angle({atomIndex=statorIdx},
            {atomIndex=pivotIdx}, {atomIndex=rotorIdx})
        rotateselected @axis {atomIndex = pivotIdx} @{curangle-toangle}
    }
    catch {
    }
}

# Selected must include fourth parameter but not the first
function setDihedralIdx (stator1idx, stator2idx, rotor1idx, rotor2idx, toangle) {
    try {
        var a1={atomIndex = stator1idx}.xyz
        var a2={atomIndex = stator2idx}.xyz
        var a3={atomIndex = rotor1idx}.xyz
        var a4={atomIndex = rotor2idx}.xyz
        var curangle =  angle(a1, a2, a3, a4)
        rotateselected {a3} {a2} @{curangle-toangle}
    }
    catch {
    }
}

function countCollisions(rc) {
    var lcAtoms = ({})
    for (var idx = {*}.min.atomIndex; idx <= {*}.max.atomIndex; idx++) {
        lcAtoms = lcAtoms or (within(kCtolerance, FALSE, {atomIndex=idx})
            and {atomIndex > idx}
            and not {rc}
            and not connected({atomIndex=idx}))
    }
    return lcAtoms.size
}

# Resolve collisions
function handleCollisions( nWard, targetIdx) {

    # For all selected atoms
    for (var iNo = {selected}.min.atomno; iNo <= {selected}.max.atomno; iNo++) {
        var idx = {(atomno=iNo) and (chain=gchain)}.atomIndex
        if ({atomindex=idx}.selected) {
       
            # Collect local colliders
            var lcAtoms = (within(kCtolerance, FALSE, {atomIndex=idx})
                and not {atomIndex=idx}
                and not {gOkCollide}
                and not connected({atomIndex=idx}))
            
            if (lcAtoms.size > 0) {
            
                # Ignore kinked BB
                if (isBBidx(idx) and angle({atomIndex=@{getCwardBBidx(idx)}},
                    {atomIndex=idx} , {atomIndex=@{getNwardBBidx(idx)}}) < 100) {
                    continue
                }
            
                #gCountdown--
                if (gCountdown < 0) {
                    timedOut("Too many collsions")
                }    
                
                # For all local water colliders, delete
                var recollect = FALSE
                for (var c = 1; c <= lcAtoms.size; c++ ) {
                    if (lcAtoms[c].group = "HOH") {
                        delete {atomIndex=@{lcAtoms[c].atomIndex}}
                        recollect = TRUE
                    }
                }
                
                # Recollect local colliders if needed
                if (recollect) {
                    lcAtoms = (within(kCtolerance, FALSE, {atomIndex=idx})
                        and not {atomIndex=idx}
                        and not {gOkCollide}
                        and not connected({atomIndex=idx}))
                    recollect = FALSE
                }
                    
                # For all local colliders
                for (var c = 1; c <= lcAtoms.size; c++ ) {
                
                    # If it is with side chain not proline, fix it
                    var cidx = lcAtoms[c].atomIndex
                    if (isSCidx(cidx) and ({atomIndex=cidx}.group != "PRO")) {
                        fixSCcollision2(cidx)
                        recollect = TRUE
                        
                        # If not fixed, exit fail
                        if (gOk2 == FALSE) {
                            return # early exit (break n jmol bug)
                        }
                    }
                }
                
                # Recollect local colliders if needed
                if (recollect) {
                    lcAtoms = (within(kCtolerance, FALSE, {atomIndex=idx})
                        and not {atomIndex=idx}
                        and not {gOkCollide}
                        and not connected({atomIndex=idx}))
                }
                    
                # For all local colliders
                for (var c = 1; c <= lcAtoms.size; c++ ) {
                
                    # If it is with O, counter-rotate
                    if ({atomIndex=@{lcAtoms[c].atomIndex}}.atomName = "O") {
                        counterRotate2(idx, lcAtoms[c].atomIndex, targetIdx)
                        
                        # If not fixed, exit fail
                        if (gOk2 == FALSE) {
                            return # early exit (break n jmol bug)
                        }
                    }
                }
            }
        }
    } # endfor iNo
}

# Rotate rotor set to move target atom to its proper place
function tugTrackIdx(targetIdx, targetPt, nWard, cDetect) {
    gOK = FALSE
    var pt = targetPt
    var dist = distance(pt, {atomIndex=targetIdx}.xyz)

    var rotors = (nWard ? gNrotors : gCrotors)
    # For a 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 i = 1; i < rotors.size; i += 4) {
                var i2 = rotors[i+1]
                var i3 = rotors[i+2]
                var i4 = rotors[i+3]
                if ((i2 != targetIdx) and (i3 != targetIdx) and (i4 != targetIdx)) {
                    if ({blocked and {atomIndex=i2}}.count == 0) {
                        var v2 = {atomIndex=i3}.xyz - {atomIndex=i2}.xyz
                        
                        var s = sin(abs(angle(v1, {0 0 0}, v2)))
                        if (s > smax) {
                            smax = s
                            imax = i
                        }
                    }
                }
            }
            
            # If no more rotors, break to next full try
            if (imax == 0) {
               break
            }
            var i1 = rotors[imax+0]
            var i2 = rotors[imax+1]
            var i3 = rotors[imax+2]
            var i4 = rotors[imax+3]
            
            # Get dihedral of rotor with target point
            var dt = angle({atomIndex=targetIdx}, {atomIndex=i2}, {atomIndex=i3}, pt)
            var dh = angle({atomIndex=i1}, {atomIndex=i2}, {atomIndex=i3}, {atomIndex=i4})
            if (dh == "NaN") {
                dh = -50
            }
            var psi = dh + dt
            var phi = dh + dt
            
            # Compute resultant psi and phi
            #  and select from target atom to first half of rotor
            var movePt = FALSE
            if (nWard) {
                if ({atomIndex=i2}.atomName="CA") {
                    psi = angle({atomIndex=@{getCwardBBidx(i1)}}, {atomIndex=i1},
                        {atomIndex=i2}, {atomIndex=i3}) + dt
                }
                else {
                    phi = angle({atomIndex=i1}, {atomIndex=i2},
                        {atomIndex=i3}, {atomIndex=@{getNwardBBidx(i3)}}) + dt
                }
                
                if ({atomIndex=i2}.atomno > {atomIndex=targetIdx}.atomno) {
                    movePt = TRUE
                    selectNward(i3, getCwardBBidx(targetIdx))
                    {atomIndex=targetIdx}.selected = TRUE
                }
                else {
                    selectCward(i2, targetIdx)
                }
            }
            else {
                if (({atomIndex=i2}.atomName="CA")) {
                    phi = angle({atomIndex=@{getNwardBBidx(i1)}}, {atomIndex=i1},
                        {atomIndex=i2}, {atomIndex=i3}) + dt
                }
                else {
                    psi = angle({atomIndex=i2}, {atomIndex=i3},
                        {atomIndex=i4}, {atomIndex=@{getCwardBBidx(i4)}}) + dt
                }
                
                if ({atomIndex=i2}.atomno < {atomIndex=targetIdx}.atomno) {
                    movePt = TRUE
                    selectCward(i3, getNwardBBidx(targetIdx))
                    {atomIndex=targetIdx}.selected = TRUE
                }
                else {
                    selectNward(i2, targetIdx)
                }
            }
            
            # If rotation within ramachandran limits
            #var ridx = 1 + (36*(((-psi\10)%180)+18)+(((phi\10)%180)+18))
            if ((abs(dt) >= 0.1) and 
                (({atomIndex=i2}.group=="GLY") or (phi < 0))) {#kRama[ridx] >= kMinRama))) {

                # If moving target point, put the target atom there
                var cp = {atomIndex=targetIdx}.xyz
                if (movePt) {
                    dt = -dt
                    {atomIndex=targetIdx}.xyz = pt
                }
                
                # Rotate to minimize vector ====================
                rotateSelected {atomIndex=i2} {atomIndex=i3} @dt

                # If collision checking
                if (cDetect) {
                
                    # If collision, back off by eighths
                    var wasCollision = FALSE
                    for (var ci = 0; ci < 4; ci++) {
                        if (ci < 3) {
                            dt /= 2
                        }                        
                        handleCollisions( nWard, targetIdx)
                        if (gOk2==FALSE) {
                            wasCollision = TRUE
                            rotateSelected {atomIndex=i2} {atomIndex=i3} @{-dt}
                        }
                        else if (wasCollision) {
                            if (ci <3) {
                                rotateSelected {atomIndex=i2} {atomIndex=i3} @{dt}
                            }
                        }
                        else {
                            break
                        }
                        
                        if (dt < 0.01) {
                            break
                        }
                    } # endfor
                }
                
                # If moving target point, put the target atom back
                if (movePt) {
                    pt = {atomIndex=targetIdx}.xyz
                    {atomIndex=targetIdx}.xyz = cp
                }
                
            }
            
            # If close enough, stop
            if (distance(pt, {atomIndex=targetIdx}) < kDtolerance) {
                gOK = TRUE
                gTargetPt = pt
                break
            }
                
            # Block rotor
            blocked |= {atomIndex=i2}
            
        }   # endfor num rotors passes
        
        if (gOK) {
            break
        }
    }   # endfor 10 passes
}

# Counter rotate bonds on either side of a BB O
function docounterRotate2(caPhiIDx, nIdx, cIdx, oIdx, caPsiIdx, dir, nWard) {

    # Rotate phi
    {atomIndex=nIdx}.selected = nWard
    {atomIndex=cIdx}.selected = nWard
    {atomIndex=oIdx}.selected = nWard
    rotateSelected {atomIndex=caPsiIdx} {atomIndex=cIdx} @{dir}

    # Counter-rotate psi
    {atomIndex=nIdx}.selected = not nWard
    {atomIndex=cIdx}.selected = not nWard
    {atomIndex=oIdx}.selected = not nWard
    rotateSelected {atomIndex=nIdx} {atomIndex=caPhiIdx} @{-dir}
}
        
function counterRotate2(aIdx, bIdx, targetIdx) {
    var selsave = {selected}
    gOk2 = TRUE
    var oIdx = aIdx
    var xIdx = bIdx
    if ({atomIndex=bIdx}.atomName=="O") {
        oIdx = bIdx
        xIdx = aIdx
    }
    var cIdx = getScBBidx(oIdx)
    var nIdx = getCwardBBidx(cIdx)
    var caPhiIdx = getCwardBBidx(nIdx)
    var caPsiIdx = getNwardBBidx(cIdx)
    
    var nWard = ({atomIndex=oIdx}.atomno < {atomIndex=targetIdx}.atomno)
    if (nWard) {
        selectCward(cIdx, targetIdx)
    }
    else {
        selectNward(nIdx, targetIdx)
    }
     
    # Until all collisions cancelled
    var dir = 5
    var ang = angle({atomIndex=xIdx}, {atomIndex=oIdx}, {atomIndex=cIdx})
    var tcount = 0
    while (within(kCtolerance, FALSE, {atomIndex=oIdx})
            and not {atomIndex=oIdx} and not connected({atomIndex=oIdx})
            and not {gOkCollide} > 0) {
        
        # Counter-rotate
        docounterRotate2(caPhiIDx, nIdx, cIdx, oIdx, caPsiIdx, dir, nWard)
        var newang = angle({atomIndex=xIdx}, {atomIndex=oIdx}, {atomIndex=cIdx})
        
        # If wrong direction once, undo and reverse
        if (newang > ang) {
            docounterRotate2(caPhiIDx, nIdx, cIdx, oIdx, caPsiIdx, -dir, nWard)
            
            # If first time, continue in opposite direction
            dir *= -1
            if (dir < 0) {
                continue
            }
        }

        # If no go, undo and exit
        tcount++
        if (tcount > (360/abs(dir))) {
            gOk2 = FALSE
            break
        }
        
    } # endwhile
    select selsave    
}

# Repair proline
function repairProline(idx) {
    var cbidx = getCBidx(idx)
    var cbno = {atomIndex=cbidx}.atomno
    var cgidx = {(atomno=@{cbno+1}) and (chain=gChain)}.atomIndex
    var cdidx = {(atomno=@{cbno+2}) and (chain=gChain)}.atomIndex
    var caidx = {(atomno=@{cbno-3}) and (chain=gChain)}.atomIndex
    var nidx = {(atomno=@{cbno-4}) and (chain=gChain)}.atomIndex
    
    select {atomIndex=cbidx}
    setAngleIdx(nidx, caidx, cbidx, 109.5)
    
    select {atomIndex=cdidx}
    setDistanceIdx(nidx, cdidx, 1.47)
    setAngleIdx(caidx, nidx, cdidx, 102.7)
    setDihedralIdx(cbidx, caidx, nidx, cdidx, 16.2)
    
    select {atomIndex=cgidx}
    setDistanceIdx(cdidx, cgidx, 1.51)
    setAngleIdx(nidx, cdidx, cgidx, 106.4)
    setDihedralIdx(caidx, nidx, cdidx, cgidx, 16.2)
}

# Repair side chain
function plicoRepairSideChain(targetIdx, nWard) {

    var idx = (nWard ? getCwardBBidx(targetIdx) : getNwardBBidx(targetIdx))

    if (({atomIndex=targetIdx}.atomName == "CA")
        and ({atomIndex=targetIdx}.group != "GLY")) {
        var cbidx = getCBidx(targetIdx)
        select none
        selectAddSideChain(cbidx)
        setAngleIdx(idx, targetIdx, cbidx, 110.0)
        setDistanceIdx(targetIdx, cbidx, 1.5)
        if ({atomIndex=targetIdx}.group != "PRO") {
            var colliders = (within(kCtolerance, FALSE, {selected})
                and not {atomIndex=targetIdx} and not {selected})
            if (colliders.size > 0) {
                if ({atomIndex=targetIdx}.group != "ALA") {
                    fixSCcollision2(cbidx)
                }
            }
        }
        else {
            if (nWard) {
            }
            else {
                setDihedralIdx(getNwardBBidx(idx), idx, targetIdx, cbidx, 174.2)
            }
        }
    }
    
    else if ({atomIndex=targetIdx}.atomName == "C") {
        var oidx = getOidx(targetIdx)
        select {atomIndex=oidx}
        setAngleIdx(idx, targetIdx, oidx, 120.0)
        setDistanceIdx(targetIdx, oidx, 1.21)
        if (nWard) {
            setDihedralIdx(getCwardBBidx(idx), idx, targetIdx, oidx, 0.0)
        }
        if ({atomIndex=idx}.group == "PRO") {
            repairProline(idx)
            var dNo = {atomIndex=targetIdx}.atomno + 4
            var dIdx = {(atomno=dNO) and (chain=gChain)}.atomIndex
            var colliders = (within(kCtolerance, FALSE, {atomIndex=dIdx})
                and not connected({atomIndex=dIdx})
                and not {atomIndex=dIdx})
            for (var i = 1; i <= colliders.size; i++) {
                if (colliders[i].atomName == "O") {
                    counterRotate2(dIdx, colliders[i].atomIndex, targetIdx)
                }
            }
        }
    }
}

# Rebuild Cward rotors set
function tugTrackC() {

    # For all bb atoms cWard of cargo
    var targetIdx = gCcargoIdx
    var okCount = 0
        
    # Allow collisions with cargo
    gOkCollide = gCargoAtoms
    var tcount = 0                
    while (targetIdx != gCanchorIdx) {
    
        # Step to next atom
        targetIdx = getCwardBBidx(targetIdx)
        
        # No collision with cargo allowed after two atoms placed
        if (tcount == 2) {
           gOkCollide = ({})
        }
        tcount++
        
        # Compute targets desired coords
        var c1idx = getCwardBBidx(targetIdx)
        var n1idx = getNwardBBidx(targetIdx)
        var n2idx = getNwardBBidx(n1Idx)
        var n3idx = getNwardBBidx(n2Idx)
        var pt = {0 0 0}
        if ({atomIndex=targetIdx}.atomName == "N") {
            var oidx = getOidx(n1idx)
            select {atomIndex=oidx}
            
            # Desired target location is trigonal O       
            setDistanceIdx(n1idx, oidx, 1.5)
            pt = getTrigonal(n2idx, n1idx, oidx, 1.37)
            setDistanceIdx(n1idx, oidx, 1.21)
        }
        else if (({atomIndex=targetIdx}.atomName == "C")
            and ({atomIndex=targetIdx}.group != "GLY")) {
            
            # Desired target location is tetragonal CB
            var cbidx = getCBidx(n1idx)
            pt = getTet(n2idx, n1idx, cbidx, 1.5)
        }
        else { # CA (or GLY C)
            
            # Save current target coords
            var cp = {atomIndex=targetIdx}.xyz
    
            # Set target atom at desired distance and angle
            select {atomIndex=targetIdx}
            setDistanceIdx(n1idx, targetIdx, 1.5)
            setAngleIdx(n2idx, n1idx, targetIdx, 110.0)
            if ({atomIndex=targetIdx}.atomName == "CA") {
                setDihedralIdx(n3idx, n2idx, n1idx, targetIdx, 180)
            }
            
            # Record and restore target
            pt = {atomIndex=targetIdx}.xyz
            {atomIndex=targetIdx}.xyz = cp
        }

        # If target not at desired location        
        if (distance(pt, {atomIndex=targetIdx}) > kDtolerance) {
            okCount = 0
            gTargetPt = pt
            var xcount = 0
            gOK = FALSE
            while ((xcount < 20) and (gOK == FALSE)) {
            
                # Rotate on cWard rotor set to move it there
                tugTrackIdx(targetIdx, pt, FALSE, FALSE)
                    #(distance(pt, {atomIndex=targetIdx}) > kMtolerance))
                xcount++
            }
        }
        else {
            gOK = TRUE
            okCount++
        }    
        
        # If successful
        if (gOK == TRUE) {
        
            # Adust any side atoms
            plicoRepairSideChain(targetIdx, FALSE)
        }
        
        # Else fail
        else {
            break
        }
        
        # If no movement in 4 tries, we are done
        if (okCount > 3) {
            break
        }
    } # endwhile (targetIdx != gCanchorIdx) {
}

# Rebuild Nward rotors set
function tugTrackN() {

    gOK = TRUE
    
    # For all bb atoms nWard of cargo
    var targetIdx = gNcargoIdx
    var okCount = 0
    
    # Allow collisions with cargo
    gOkCollide = gCargoAtoms
    var tcount = 0                
    while (targetIdx != gNanchorIdx) {

        # Step to next atom
        targetIdx = getNwardBBidx(targetIdx)
        
        # No collision with cargo allowed after two atoms placed
        if (tcount == 2) {
           gOkCollide = ({})
        }
        tcount++
        
        # Compute targets desired coords
        var n1idx = getNwardBBidx(targetIdx)
        var c1idx = getCwardBBidx(targetIdx)
        var c2idx = getCwardBBidx(c1idx)
        var c3idx = getCwardBBidx(c2idx)
        var pt = {0 0 0}
        if ({atomIndex=targetIdx}.atomName == "CA") {
        
            # Desired target location is trigonal O       
            var oidx = getOidx(c1idx)
            select {atomIndex=oidx}
            setDistanceIdx(c1idx, oidx, 1.39)
            pt = getTrigonal(c2idx, c1idx, oidx, 1.41)
            setDistanceIdx(c1idx, oidx, 1.21)
        }
        else if (({atomIndex=targetIdx}.atomName == "N")
            and ({atomIndex=targetIdx}.group != "GLY")) {
            
            # Desired target location is r-tetragonal CB       
            var cbidx = getCBidx(c1idx)
            pt = getTet(cbidx, c1idx, c2idx, 1.5)
        }
        else { # C
        
            # Save current target coords
            var cp = {atomIndex=targetIdx}.xyz
            
            # Set target atom at desired distance and angle
            select {atomIndex=targetIdx}
            setDistanceIdx(c1idx, targetIdx, 1.37)
            setAngleIdx(c2idx, c1idx, targetIdx, 110.0)
                
            if ({atomIndex=targetIdx}.group == "PRO") {
                setDihedralIdx(c3idx, c2idx, c1idx, targetIdx, -57.0)
            }
        
            # Record and restore target
            pt = {atomIndex=targetIdx}.xyz
            {atomIndex=targetIdx}.xyz = cp
        }
        
        # If target not at desired location        
        if (distance(pt, {atomIndex=targetIdx}) > kDtolerance) {
            var okCount = 0
            gTargetPt = pt
            var xcount = 0
            gOK = FALSE
            while ((xcount < 20) and (gOK == FALSE)) {
            
                # Rotate on cWard rotor set to move it there
                tugTrackIdx(targetIdx, pt, TRUE, FALSE)
                    #(distance(pt, {atomIndex=targetIdx}) > kMtolerance))
                xcount++
            }
        }
        else {
            gOK = TRUE
            okCount++
        }    

        # If sucessful
        if (gOK == TRUE) {
        
            # Adust any side atoms
            plicoRepairSideChain(targetIdx, TRUE)
        }
        
        # Else fail
        else {
            break
        }
        
        # If no movement in 4 tries, we are done
        if (okCount > 3) {
            break
        }
        
    }   # endwhile (targetIdx != gNanchorIdx) {

}

# gPlicoRecord is maintained by the macro pilcoRecord
function plicoRecord(s) {
    var g = format("show file \"%s\"", gPlicoRecord)
    var ls = script(g)
    if (ls.find("FileNotFoundException")) {
        ls = ""
    }
    ls += s
    write var ls @gPlicoRecord
}

# gPlicoRecord is maintained by the macro pilcoRecord
function translateSelectedRecord(pt) {
    if (gPlicoRecord != "") {
        plicoRecord(format("select %s;translateSelected %s;", {selected}, pt))
    }
    translateSelected @pt
}

# gPlicoRecord is maintained by the macro pilcoRecord
function rotateSelectedRecord(g1pivotIdx, axis, a) {
    if (gPlicoRecord != "") {
        plicoRecord(format("select %s;", {selected}))
        plicoRecord(format("rotateSelected {atomIndex=%d} @%s @%s;",
            g1pivotIdx, axis, a))
    }
    rotateSelected {atomIndex=g1pivotIdx} @axis @a
}

function collectSCrotors(no) {
    var scBondIdxs = array()
    for (var iNo = no; iNo >= 0; iNo--) {
        var ile = 0
        switch ({(atomno=iNo) and (chain=gChain)}.atomName) {
        case "CA" :
            return scBondIdxs # Early exit since break 1 appears broken
            #break 1
        case "CZ" :
            if ({(atomno=iNo) and (chain=gChain)}.group == "TYR") {
                break
            }
        case "CE" :
            if ({(atomno=iNo) and (chain=gChain)}.group == "MET") {
                break
            }
        case "CG1" :
            if ({(atomno=iNo) and (chain=gChain)}.group == "VAL") {
                break
            }
            if ({(atomno=iNo) and (chain=gChain)}.group == "ILE") {
                ile = 1
            }
        case "NE" :
        case "CD" :
        case "SD" :
        case "CG" :
        case "CB" :
            scBondIdxs += {(atomno=@{iNo+1+ile}) and (chain=gChain)}.atomIndex
            scBondIdxs += {(atomno=@{iNo+0}) and (chain=gChain)}.atomIndex
            if ({(atomno=iNo) and (chain=gChain)}.atomName%2 == "CG") {
                scBondIdxs += {(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex
                scBondIdxs += {(atomno=@{iNo-4}) and (chain=gChain)}.atomIndex
            }
            else if ({(atomno=iNo) and (chain=gChain)}.atomName == "CB") {
                scBondIdxs += {(atomno=@{iNo-3}) and (chain=gChain)}.atomIndex
                scBondIdxs += {(atomno=@{iNo-4}) and (chain=gChain)}.atomIndex
            }
            else {
                scBondIdxs += {(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex
                scBondIdxs += {(atomno=@{iNo-2}) and (chain=gChain)}.atomIndex
            }
            break
        }
    
    }
    
    return scBondIdxs 
}

# Drag Side Chain
function dragSC() {

    var iNo = {atomIndex=gSCidx}.atomno
    
    if ({atomIndex=gSCidx}.group != "PRO") {
    
        var scBondIdxs = collectSCrotors( iNo)
        var numChi = scBondIdxs.size / 4
        var dist = distance({atomIndex=gSCidx}.xyz, gSCpt)
        # For all rotor combinations
        var dh = array()
        for (var i = 0; i < numChi; i++) {
            dh += angle({atomIndex=@{scBondIdxs[4+(4*i)]}},
                {atomIndex=@{scBondIdxs[3+(4*i)]}},
                {atomIndex=@{scBondIdxs[2+(4*i)]}},
                {atomIndex=@{scBondIdxs[1+(4*i)]}})
        }
        for (var i = 0; i < numChi; i++) {
            var rot = -120
            for (var j = 0; j < 6; j++) {
                rot += 60*j
                selectAddSideChain(scBondIdxs[1+(4*i)])
                setDihedralIdx(scBondIdxs[4+(4*i)], scBondIdxs[3+(4*i)],
                    scBondIdxs[2+(4*i)], scBondIdxs[1+(4*i)], rot)
                var newDist = distance({atomIndex=gSCidx}.xyz, gSCpt)
                
                # Find the best
                if (newDist < dist) {
                    dist = newDist
                    for (var k = 0; k < numChi; k++) {
                        dh[k+1] = angle({atomIndex=@{scBondIdxs[4+(4*k)]}},
                        {atomIndex=@{scBondIdxs[3+(4*k)]}},
                        {atomIndex=@{scBondIdxs[2+(4*k)]}},
                        {atomIndex=@{scBondIdxs[1+(4*k)]}})
                    }
                }
            }
        }
        
        # Now set the best
        for (var i = 0; i < numChi; i++) {
            selectAddSideChain(scBondIdxs[1+(4*i)])
            setDihedralIdx(scBondIdxs[4+(4*i)], scBondIdxs[3+(4*i)],
                scBondIdxs[2+(4*i)], scBondIdxs[1+(4*i)], dh[i+1])
        }
    }
    else { # PRO - toggle between puckers up and down
        var icd = {(atomno=@{iNo+1}) and (chain=gChain)}.atomIndex
        var icb = {(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex
        var ica = {(atomno=@{iNo-4}) and (chain=gChain)}.atomIndex
        var in = {(atomno=@{iNo-5}) and (chain=gChain)}.atomIndex
        select {atomIndex=gSCidx}
       
        if (angle({atomIndex=ica}, {atomIndex=in},
            {atomIndex=icd}, {atomIndex=gSCidx}) < -10.0) {
            setDihedralIdx(ica, in, icd, gSCidx, 8.7)
            setAngleIdx(in, icd, gSCidx, 110.0)
            setDistanceIdx(icd, gSCidx, 1.5)
        }
        else {
            setDihedralIdx(ica, in, icd, gSCidx, -29.5)
            setAngleIdx(in, icd, gSCidx, 108.8)
            setDistanceIdx(icd, gSCidx, 1.5)
        }
    }
    
    draw gSCcircle CIRCLE {atomIndex=gSCidx} MESH NOFILL
    gSCpt = {atomIndex=gSCidx}.xyz
}
 
# Fix side chain collisions
function fixSCcollision2(idx) {
    gOk2 = FALSE
    var iNo = {atomIndex=idx}.atomno
    var resno = {(atomno=iNo) and (chain=gChain)}.resno
    
    # Get SC terminus
    while (resno == {(atomno=iNo) and (chain=gChain)}.resno) {
        iNo++
    }
    iNo--
    
    var sc = array()
    var iBno = iNo
    while ({(atomno=iBno) and (chain=gChain)}.atomName != "CB") {
        sc += {(atomno=iBno) and (chain=gChain)}
        iBno--
    }
    var cbidx = {(atomno=iBno) and (chain=gChain)}.atomIndex
    
    var scBondIdxs = collectSCrotors( iNo)
    var numChi = scBondIdxs.size / 4
    # For all rotor combinations
    for (var i = 0; i < numChi; i++) {
        var rot = -120
        for (var j = 0; j < 6; j++) {
            rot += 60
            selectAddSideChain(scBondIdxs[1+(4*i)])
            #setDihedralIdx(scBondIdxs[4+(4*i)], scBondIdxs[3+(4*i)],
            #    scBondIdxs[2+(4*i)], scBondIdxs[1+(4*i)], rot)
            setDihedralIdx(scBondIdxs[1+(4*i)], scBondIdxs[2+(4*i)],
                scBondIdxs[3+(4*i)], scBondIdxs[4+(4*i)], rot)
            # If no collision, exit
            colliders = (within(kCtolerance, FALSE, {sc})
                and not {atomIndex=cbidx} and not {sc})
                        
            # If it is with water, delete the water
            for (var c = 1; c < colliders.size; c++ ) {
                if (colliders[c].group = "HOH") {
                    delete {atomIndex=@{colliders[c].atomIndex}}
                    colliders = {colliders and not @{colliders[c]}}
                }
            }
                
            if (colliders.size == 0) {
                gOk2 = TRUE
                return # Early exit since break 1 appears broken
                #break 1
            }
            
        }
    }
}

function isMovableSideChain(aIdx) {

    var ret = (({atomIndex=aIdx}.group != "PRO")
        or ({atomIndex=aIdx}.atomName == "CG"))
    switch({atomIndex=aIdx}.atomName) {
    case "N":
    case "CA":
    case "C":
    case "CB":
    case "O":
    case "O4\'":
        ret = FALSE
        break
    }
    return ret
}

# gFreeze is maintained by the script plicoFreeze that allows the user
# to inhibit rotation on selected rotors
function isRotorAvailable(idx) {
    return (gFreeze.find(idx) == 0)
}

function collectBBrotors(nWard) {
    var anchorNo = (nWard
        ? ((gNanchorIdx >= 0) ? {atomIndex=gNanchorIdx}.atomno : gMinNo)
        : ((gCanchorIdx >= 0) ? {atomIndex=gCanchorIdx}.atomno : gMaxNo))
    var cargoNo = (nWard
        ? ((gNcargoIdx >= 0) ? {atomIndex=gNcargoIdx}.atomno
        : {atomIndex=gCcargoIdx}.atomno)
        : {atomIndex=gCcargoIdx}.atomno)
    var rotors = array()
    if (cargoNo < anchorNo) {
    
        # If cWard cargo is C, include its psi
        if ({(atomno=cargoNo) and (chain=gChain)}.atomName == "C") {
            #cargoNo--
        }
        
        for (var iNo = cargoNo; iNo <= anchorNo; iNo++) {
            if ({(atomno=iNo) and (chain=gChain)}.atomName == "CA") {
                if (isRotorAvailable(iNo)) {
                    if (({(atomno=iNo) and (chain=gChain)}.group != "PRO") and (iNo > cargoNo)) { # phi
                        rotors += [{(atomno=@{getCmNo(iNo-1)}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex]
                        rotors += [{(atomno=@{iNo}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo+1}) and (chain=gChain)}.atomIndex]
                    }
                    if (iNo != (anchorNo-1)) { # psi
                        rotors += [{(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo}) and (chain=gChain)}.atomIndex]
                        rotors += [{(atomno=@{iNo+1}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{getNpNo(iNo+2)}) and (chain=gChain)}.atomIndex]
                    }
                }
            }
        }
    }
    else {
    
        # If nWard cargo is C, include its phi
        if ({(atomno=cargoNo) and (chain=gChain)}.atomName == "N") {
            #cargoNo++
        }
        
        for (var iNo = cargoNo; iNo >= anchorNo; iNo--) {
            if ({(atomno=iNo) and (chain=gChain)}.atomName == "CA") {
                if (isRotorAvailable(iNo)) {
                    if ((iNo != (anchorNo-1)) and (iNo < cargoNo)) { # psi
                        rotors += [{(atomno=@{getNpNo(iNo+2)}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo+1}) and (chain=gChain)}.atomIndex]
                        rotors += [{(atomno=@{iNo}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex]
                    }
                    if ({(atomno=iNo) and (chain=gChain)}.group != "PRO") { # phi
                        rotors += [{(atomno=@{iNo+1}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{iNo}) and (chain=gChain)}.atomIndex]
                        rotors += [{(atomno=@{iNo-1}) and (chain=gChain)}.atomIndex,
                            {(atomno=@{getCmNo(iNo-1)}) and (chain=gChain)}.atomIndex]
                    }
                }
            }
        }
    }
    
    if (nWard) {
        gNrotors = rotors
        if (FALSE) {#DEBUG
            print ("gNrotors atomnos")#DEBUG
            for (var i = 1; i < gNrotors.size; i += 4) {#DEBUG
                print format("%d %d %s", {atomIndex=@{gNrotors[i]}}.atomno,
                    {atomIndex=@{gNrotors[i+1]}}.atomno, format("%d %d",
                    {atomIndex=@{gNrotors[i+2]}}.atomno,
                    {atomIndex=@{gNrotors[i+3]}}.atomno))
            }
        }#DEBUG
    }
    else {
        gCrotors = rotors
        if (FALSE) {#DEBUG
            print ("gCrotors atomnos")#DEBUG
            for (var i = 1; i < gCrotors.size; i += 4) {#DEBUG
                print format("%d %d %s", {atomIndex=@{gCrotors[i]}}.atomno,
                    {atomIndex=@{gCrotors[i+1]}}.atomno, format("%d %d",
                    {atomIndex=@{gCrotors[i+2]}}.atomno,
                    {atomIndex=@{gCrotors[i+3]}}.atomno))
            }
        }#DEBUG
    }
}

function collectRotors() {
    collectBBrotors(FALSE)
    collectBBrotors(TRUE)
}

function tugSideChain(pt) {
                
    # If destination atom defined
    if (gDestAtomIdx >= 0) {
        v = {atomIndex=gDestAtomIdx}.xyz - {atomIndex=gSCidx}.xyz
        if (abs(angle({atomIndex=gDestAtomIdx}.xyz, {0 0 0}, pt)) < 90) {
            pt = -v/20.0
        }
        else {
            pt = v/20.0
        }
    }
    gSCpt += pt
    draw arrow {atomIndex=gSCidx} @gSCpt 
}

function setColors() {
    select all
    color {selected} @gScheme
    color {atomIndex=gCanchorIdx} @gAltScheme
    color {atomIndex=gNanchorIdx} @gAltScheme
    color {atomIndex=g1pivotIdx} green
    color {atomIndex=g2pivotIdx} green
    color @gCargoAtoms @gAltScheme
    select {(atomIndex=gCcargoIdx) or (atomIndex=gNcargoIdx)}
    halo on
    select {atomIndex=gDestAtomIdx}
    star on
    select none
}

function clearAtomIdxs() {
    gCcargoIdx = -1
    gNcargoIdx = -1
    gCanchorIdx = -1
    gNanchorIdx = -1
    g1pivotIdx = -1
    g2pivotIdx = -1
    gDestAtomIdx = -1
    gSCidx = -1
}

function timedOut (s) {
    timeout ID"tug" OFF
    prompt(s)
    gBusy = FALSE
    background ECHO yellow
    restore state gState
    select gCargoAtoms
    refresh
    quit
}

function recordDrag() {
    var ls = format("select %s;", {selected})
    ls += format("gCanchorIdx = %d;", gCanchorIdx)
    ls += format("gCanchorNo = %d;", gCanchorNo)
    ls += format("gNanchorIdx = %d;", gNanchorIdx)
    ls += format("gNanchorNo = %d;", gNanchorNo)
    ls += format("gNanchorPidx = %d;", gNanchorPidx)
    ls += format("gCanchorXyz = %s;", gCanchorXyz)
    ls += format("gNanchorXyz = %s;", gNanchorXyz)
    ls += format("gNanchorPxyz = %s;", gNanchorPxyz)
    ls += format("gCcargoIdx = %d;", gCcargoIdx)
    ls += format("gNcargoIdx = %d;", gNcargoIdx)
    ls += format("gCcargoXyz = %s;", gCcargoXyz)
    ls += format("gNcargoXyz = %s;", gNcargoXyz)
    ls += format("gCcargoNo = %d;", gCcargoNo)
    ls += format("gNcargoNo = %d;", gNcargoNo)
    ls += format("gDestAtomIdx = %d;", gDestAtomIdx)
    ls += format("g1pivotIdx = %d;", g1pivotIdx)
    ls += format("g2pivotIdx = %d;", g2pivotIdx)
    ls += format("gAc = %s;", gAc)
    ls += format("gChain = \"%s\";", gChain)
    ls += format("gMinNo = %d;", gMinNo)
    ls += format("gMaxNo = %d;", gMaxNo)    
    ls += format("gCargoAtoms = %s;", gCargoAtoms)
    ls += format("gSCidx = %d;", gSCidx)
    ls += format("gSCcircle = %d;", gSCcircle)
    ls += format("gSCpt = %s;", gSCpt)
    ls += "collectRotors();"
    ls += "tugDragDoneMB();"
    plicoRecord(ls)
}

# Bound to LEFT-UP by tugEnableDrag
function tugDragDoneMB() {
    if (gPlicoRecord != "") {
        recordDrag()
    }
 
    # Move by rotation on rotor sets, smallest first
    gBusy = TRUE
    background ECHO pink
    refresh

    # If side chain mode
    if (gSCidx >= 0) {
        dragSC()
    }
    
    # Else
    else {
        gOK = TRUE
        timeout ID"tug" 20.0 "timedOut(\"Tug timed out\")"
        gCountdown = kCollisionLimit
        if ((gCrotors.size < gNrotors.size) or (gNanchorIdx < 0)) {
            if (gCrotors.size > 4) {
                tugTrackC()  # PULLSTRING MODEL
            }
            if (gOK and (gNrotors.size > 4)) {
                tugTrackN()  # PULLSTRING MODEL
            }
        }
        else {
            if (gNrotors.size > 4) {
                tugTrackN()  # PULLSTRING MODEL
            }
            if (gOK and (gCrotors.size > 4)) {
                tugTrackC()  # PULLSTRING MODEL
            }
        }
        timeout ID"tug" OFF
        
        # If anchor angles acute, fail
        if (gOK == TRUE) {
            if (gCanchorIdx >= 0) {
                var ic = getCwardBBidx(gCanchorIdx)
                var in = getNwardBBidx(gCanchorIdx)
                if ((ic >= 0) and
                    angle({atomIndex=ic}, {atomIndex=gCanchorIdx}, {atomIndex=in})
                    < 100.0) {
                    gOK = FALSE
                }
            }
            if (gNanchorIdx >= 0) {
                var ic = getCwardBBidx(gNanchorIdx)
                var in = getNwardBBidx(gNanchorIdx)
                if ((in >= 0) and
                    angle({atomIndex=ic}, {atomIndex=gNanchorIdx}, {atomIndex=in})
                    < 100.0) {
                    gOK = FALSE
                }
            }
        }
        
        # If too far
        if (gOK == FALSE) {
            timedOut("TUG TOO FAR!")
        }
        
        # Else OK
        else {

            select ((atomno >= gNanchorNo) and (atomno <= gCanchorNo)
                and (chain = gChain))
            var idx = {atomno=@{{chain=gChain}.min.atomno}}.atomIndex

            
            var i = 0
            for (; i < 3; i++) {
                handleCollisions(FALSE, idx)
                if (countCollisions() == 0) {
                    break
                }
            }
            if (i == 3) {
                var p = prompt("Unable to handle all collisions!")
                restore state gState
            }                    
        }
    }
    select {gCargoAtoms}
    gBusy = FALSE
    background ECHO yellow
    refresh
}

# Bound to ALT-LEFT-DRAG by tugEnableDrag
function tugDragMB() {
    if (gBusy == FALSE) {
        var dx = (20.0 * (_mouseX - gMouseX))/_width
        var dy = (20.0 * (_mouseY - gMouseY))/_height
        var q = quaternion(script("show rotation"))
        var ptd = {@dx @dy 0}
        var pt = (!q)%ptd
    
        if (distance(pt,  {0 0 0}) > 0.004) {

            # If sidechain mode
            if (gSCidx >= 0) {
                tugSideChain(pt)
            }

            # Else            
            else {
            
                # If new drag
                if (gNewDrag) {
                    gNewDrag = FALSE
                    save state gState
                }
                
                #mp = ({atomIndex=gNcargoIdx}.xyz + {atomIndex=gCcargoIdx}.xyz )/2.0
                
                # If destination atom defined
                if (gDestAtomIdx >= 0) {
                    var v = {atomIndex=gDestAtomIdx}.xyz - {selected}.xyz
                    if (abs(angle({atomIndex=gDestAtomIdx}.xyz, {0 0 0}, pt)) < 90) {
                        pt = -v/20.0
                    }
                    else {
                        pt = v/20.0
                    }
                }
                
                # Move the cargo
                select {gCargoAtoms}
                
                # If pivots defined, rotate it
                if (g1pivotIdx >= 0) {
                
                    # If two pivots
                    var axis = {0 0 0}
                    if (g2pivotIdx >= 0) {
                        axis = {atomIndex=g2pivotIdx}
                    }
                    
                    # Else
                    else {
                        axis = cross(pt, {0 0 0}) + {atomIndex=g1pivotIdx}.xyz
                    }
                    
                    var a = ((abs(angle(pt + {atomIndex=g1pivotIdx}.xyz,
                        {atomIndex=g1pivotIdx}.xyz, axis)) < 90) ? -2 : 2)
                    #rotateSelected {atomIndex=g1pivotIdx} @axis @a
                    rotateSelectedRecord(g1pivotIdx, axis, a)
        
                }
                
                # Else translate it
                else {
                    #translateSelected @pt
                    translateSelectedRecord(pt)
                }
                
                # If collisions
                var aset = {((atomno < gCanchorNo) and (atomno > gNanchorNo)
                    and (chain=gChain))}   
                var ac = (within(kCtolerance, FALSE, {aset}) and not {aset}
                    and not connected(aset))
                if (ac.size > 0) {
                    # Resolve them
                    for (var i = 1; i <= ac.size;  i++) {
                        select ac[i]
                        handleCollisions()
                    }
                    
                    # If unable
                    if (gOk2 == FALSE) {
                    
                        # Back off
                        background ECHO pink
                        delay 1
                        if (g1pivotIdx >= 0) {
                            rotateSelectedRecord(g1pivotIdx, axis, -a)
                        }
                        else {
                            translateSelectedRecord(-pt)
                        }
                        background ECHO yellow
                    }
                }
                #tugDragDoneMB()
            }
            
            gMouseX = _mouseX
            gMouseY = _mouseY
        }
        select {gCargoAtoms}
    }
}

# Bound to ALT-LEFT-DOWN by tugEnableDrag
function tugMarkMB() {
    gMouseX = _mouseX
    gMouseY = _mouseY
    gNcargoXyz = {atomIndex=gNcargoIdx}.xyz
    gCcargoXyz = {atomIndex=gCcargoIdx}.xyz
    gNewDrag = TRUE
}

# Called by tugCargoMB and by tugAnchorMB or tugPivotMB when cargo exists
function tugEnableDrag() {
    gEcho = "ALT-CLICK=set cargo range|ALT-DRAG=move|SHIFT=set anchors|ALT-CTRL=set pivots|ALT-SHIFT=set dest atom |DOUBLE-CLICK=exit"
    echo @gEcho
    
    # Allow atoms to be dragged
    bind "ALT-LEFT-DOWN" "tugMarkMB";
    bind "ALT-LEFT-UP" "tugDragDoneMB";#ONCEMODE
    bind "ALT-LEFT-DRAG" "tugDragMB";
}

# Bound to SHIFT-LEFT-CLICK by tugCargoMB
function tugAnchorMB() {
    if ({atomIndex=_atomPicked}.chain == gChain) {
        var aPidx = getScBBidx( _atomPicked)
        
        var pno = {atomIndex=aPidx}.atomno
        if (pno > {atomIndex=gCcargoIdx}.atomno) {
            color {atomIndex=gCanchorIdx} @gScheme
            if (gCanchorIdx == aPidx) {
                gCanchorIdx = -1
                gCanchorNo = gMaxNo + 1
            }
            else {
                gCanchorIdx = aPidx
                gCanchorNo = {atomIndex=gCanchorIdx}.atomno
                gCanchorXyz = {atomIndex=gCanchorIdx}.xyz
                color {atomIndex=gCanchorIdx} @gAltScheme
            }
            collectBBrotors(FALSE)
        }
        if (pno < {atomIndex=gNcargoIdx}.atomno) {
            color {atomIndex=gNanchorIdx} @gScheme
            if (gNanchorIdx == aPidx) {
                gNanchorIdx = -1
                gNanchorNo = gMinNo - 1
            }
            else {
                gNanchorIdx = aPidx
                gNanchorNo = {atomIndex=gNanchorIdx}.atomno
                gNanchorPidx = getCwardBBidx( aPidx)
                gNanchorXyz = {atomIndex=gNanchorIdx}.xyz
                gNanchorPxyz = {atomIndex=gNanchorPidx}.xyz
                color {atomIndex=gNanchorIdx} @gAltScheme
            }
            collectBBrotors(TRUE)
        }
    }
    
    # Get connectors between fixed and moving part
    var aset = {((atomno < gCanchorNo) and (atomno > gNanchorNo)
        and (chain=gChain))}   
    gAc = (within(kCtolerance, FALSE, {aset}) and not {aset})
    select {gCargoAtoms}
}

# Bound to ALT-SHIFT-LEFT-CLICK by tugCargoMB
function tugDestAtomMB() {
    var aOk = TRUE
    if ({atomIndex=_atomPicked}.chain == gChain) {

        var pno = {atomIndex=_atomPicked}.atomno
        if ((pno <= {atomIndex=gCcargoIdx}.atomno) and (pno >= {atomIndex=gNcargoIdx}.atomno)) {
            aOk = FALSE
        }
    }
    if (aOk) {
        select {atomIndex=gDestAtomIdx}
        star off
        if (gDestAtomIdx == _atomPicked) {
            gDestAtomIdx = -1
        }
        else {
            gDestAtomIdx = _atomPicked
            select {atomIndex=gDestAtomIdx}
            star on
        }
        select {gCargoAtoms}
    }
}

# Bound to CTRL-LEFT-CLICK by tugCargoMB
function tugPivotMB() {
    if (g1pivotIdx == _atomPicked) {
        color {atomIndex=g1pivotIdx} @gScheme
        if (g2pivotIdx >= 0) {
            g1pivotIdx = g2pivotIdx
            g2pivotIdx = -1
        }
        else {
            g1pivotIdx = -1
        }
    }
    else if (g2pivotIdx == _atomPicked) {
        color {atomIndex=g2pivotIdx} @gScheme
        g2pivotIdx = -1
    }
    else if (g1pivotIdx >= 0) {
        if (g2pivotIdx >= 0) {
            color {atomIndex=g2pivotIdx} @gScheme
        }
            
        g2pivotIdx = _atomPicked
        color {atomIndex=g2pivotIdx} green
    }
    else {
        g1pivotIdx = _atomPicked
        color {atomIndex=g1pivotIdx} green
    }
    select {gCargoAtoms}
}

# Bound to ALT-LEFT-CLICK by plicoTug
function tugCargoMB() {
    
    if (gChain != {atomIndex=_atomPicked}.chain) {
        clearAtomIdxs()
        setColors()
        gChain = {atomIndex=_atomPicked}.chain
    }
    
    # If movable side chain atom picked
    if (isMovableSideChain( _atomPicked)) {
        if (gSCidx >= 0) {
            draw gSCcircle DELETE
        }
        if (gSCidx == _atomPicked) {
            gSCidx = -1
        }        
        else {
            gSCidx = _atomPicked
            gSCpt = {atomIndex=gSCidx}.xyz
            draw gSCcircle CIRCLE {atomIndex=gSCidx} MESH NOFILL
        }
    }
    else if ((gChain == "") or ({atomIndex=_atomPicked}.chain == gChain)) {
        gMinNo = {chain=gChain}.atomno.min    
        gMaxNo = {chain=gChain}.atomno.max
        if (gNanchorIdx < 0) {
            gNanchorNo = gMinNo - 1
        }
        if (gCanchorIdx < 0) {
            gCanchorNo = gMaxNo + 1    
        }
        var aPidx = getScBBidx( _atomPicked)
        
        gSCidx = -1
        draw gSCcircle DELETE

        # If existing cWard cargo picked
        if (gCcargoIdx == aPidx) {
            
            # Clear the highlight
            select {atomIndex=gCcargoIdx}
            halo off
            
            # If nWard cargo exists, mark it as the cWard cargo
            if (gNcargoIdx != gCcargoIdx) {
                gCcargoIdx = getCpIdx(gNcargoIdx)
                gCcargoNo = {atomIndex=gCcargoIdx}.atomno
            }
            else {
                gCcargoIdx = -1
                gNcargoIdx = -1
            }
        }
        else if (gNcargoIdx == aPidx) {
            select {atomIndex=gNcargoIdx}
            halo off
            gNcargoIdx = getNmIdx(gCcargoIdx)
            gNcargoNo = {atomIndex=gNcargoIdx}.atomno
        }
        else if (gCcargoIdx >= 0) {
            var no = {atomIndex=aPidx}.atomno
        
            # If pick is nWard of it
            if (no < {atomIndex=gCcargoIdx}.atomno) {
            
                # If exists, clear its highlight
                if (gNcargoIdx != gCcargoIdx) {
                    select {atomIndex=gNcargoIdx}
                    halo off
                }
                
                # Set new nWard cargo and highlight it
                gNcargoIdx = getNmIdx(aPidx)
                gNcargoNo = {atomIndex=gNcargoIdx}.atomno
            }
            
            # Else cWard
            else {
            
                # Clear its old highlight
                select {atomIndex=gCcargoIdx}
                if (gNcargoIdx != gCcargoIdx) {
                    halo off
                }
             
                # Set new cWard cargo and highlight
                gCcargoIdx = getCpIdx(aPidx)
                gCcargoNo = {atomIndex=gCcargoIdx}.atomno
            }
        }
        
        # Else no cWard cargo
        else {
        
            # Set new cWard cargo and highlight
            gCcargoIdx = getCpIdx(aPidx)
            gCcargoNo = {atomIndex=gCcargoIdx}.atomno
            gNcargoIdx = getNmIdx(gCcargoIdx)
            gNcargoNo = {atomIndex=gNcargoIdx}.atomno
            
            # Set default cWard anchor at cWard N
            var iNo = gMaxNo
            for (; iNo > 0; iNo--) {
                if ({(atomno=iNo) and (chain=gChain)}.atomName == "N") {
                    break;
                }
            }
            gCanchorIdx = {(atomno=iNo) and (chain=gChain)}.atomIndex
            gCanchorNo = {atomIndex=gCanchorIdx}.atomno
            gCanchorXyz = {atomIndex=gCanchorIdx}.xyz
        }
        
        # If any anchor now inside cargo cluster, kill it
        if ({atomIndex=gCanchorIdx}.atomno <= {atomIndex=gCcargoIdx}.atomno) {
            gCanchorIdx = -1
            gCanchorNo = gMaxNo + 1
        }
        if ({atomIndex=gNanchorIdx}.atomno >= {atomIndex=gNcargoIdx}.atomno) {
            gNanchorIdx = -1
            gNanchorNo = gMinNo - 1
        }
        
        # Highlight cargo cluster
        selectNward(gCcargoIdx, gNcargoIdx)
        gCargoAtoms = {selected}
        setColors()
                             
        # Collect the rotor sets
        collectRotors()
    
        # Get connectors between fixed and moving part
        var aset = {((atomno < gCanchorNo) and (atomno > gNanchorNo)
            and (chain=gChain))}   
        gAc = (within(kCtolerance, FALSE, {aset}) and not {aset})
    }
    
    # Enable dragging
    tugEnableDrag()

    # Bind other keys
    bind "SHIFT-LEFT-CLICK" "_pickAtom";
    bind "SHIFT-LEFT-CLICK" "+:tugAnchorMB";
    
    bind "ALT-CTRL-LEFT-CLICK" "_pickAtom";
    bind "ALT-CTRL-LEFT-CLICK" "+:tugPivotMB";
    
    bind "ALT-SHIFT-LEFT-CLICK" "_pickAtom";
    bind "ALT-SHIFT-LEFT-CLICK" "+:tugDestAtomMB";
    
    select {gCargoAtoms}
}

# Top level of Tug
function plicoTug() {
    set allowModelKit TRUE
    set allowRotateSelected TRUE
    set allowMoveAtoms TRUE
    
    # Push selected
    gSelSaves = {selected}
    select all
    gZoom = script("show zoom")
    gRotate = script("show rotation")
    write tugsave.pdb
    select none
    
    gScheme = defaultColorScheme
    gAltScheme = ((gScheme == "Jmol") ? "Rasmol" : "Jmol")
    set echo TOP LEFT
    background ECHO yellow
    gEcho = "ALT-CLICK=set cargo range"
    echo @gEcho
    gCrotors = array()
    gNrotors = array()
    clearAtomIdxs()
    gChain = ""
    unbind

    set picking ON
    bind "ALT-LEFT-CLICK" "_pickAtom";
    bind "ALT-LEFT-CLICK" "+:tugCargoMB";
    bind "DOUBLE" "tugExit";
}

# Bound to DOUBLE by plicoTug
function tugExit() {
    var p = prompt("Exit tug?", "Yes|No|Undo", TRUE)
    if (p == "Undo") {
        load tugsave.pdb
        script inline gZoom
        rotate @gRotate
        echo Tug session undone
        if (gPlicoRecord != "") {
            plicoRecord("load tugsave.pdb;")
        }
    }
    if (p != "No") {
        unbind
        halo off 
        set allowMoveAtoms FALSE
        echo
        select all
        halo off
        star off
        color {selected} @gScheme
        draw gSCcircle DELETE
        gBusy = FALSE
        background ECHO yellow
        
        # Pop selected
        select gSelSaves
    }
}
# End of TUG.SPT

Contributors

Remig