User:Remig/plico/remapNT

From Jmol
< User:Remig‎ | plico
Revision as of 15:23, 16 May 2014 by Remig (talk | contribs) (lc all functions)
Jump to navigation Jump to search

RemapNT allows you to change the chain ID, atom numbers and/or residue numbers of a polynucleotide chain by mouse actions. It also calculates group values [nucleotide names (DU, A, etc.)]. Finally it prints the resultant 1 char string to the console.

When you click on a polynucleotide chain, it gives the current chain ID, residue, residue number and atom number of the most 5'ward atom in that chain. You may then edit any value (except residue). Remap then remaps the entire chain based on those values by conventional increments and identifies each nucleotide residue.

Note that it will also remove all waters, ligands, hydrogens and %B alternates when any chain is updated.

RemapNT is a member of the Plico suite of protein folding tools described here. It may be installed and accessed as a macro with the file:

Title=PLICO Remap Polynucleotide
Script=script <path to your scripts folder>/remapNT.spt;plico_remap_nt

saved as plicoRemapNT.macro in your .jmol/macros directory as described in Macro.

Copy and paste the following to a text editor and save to your scripts directory as remapNT.spt:

#   remapNT - Jmol script by Ron Mignery
#   v1.4 beta    5/16/2014 -lc all functions
#
#   Calculate or change polynucleotide chain, atom number, residue names
#    and/or residue numbers and print the resultant 1 char string
#

function find_5_prime(pIdx) {
    while (pIdx >= 0) {
    
        # Find C3'
        var c3pIdx = -1
        select {atomIndex=pIdx}
        var cSet = {selected}
        for (var i = 1; i <= cSet.size; i++) {
            var ocSet = connected(cSet[i])
            for (var j = 1; j <= ocSet.size; j++) {
                var occSet = connected(ocSet[j]) and not {atomIndex=pIdx}
                for (var k = 1; k <= occSet.size; k++) {
                    if (connected(occSet[k]).size > 2) {
                        c3pIdx = occSet[k].atomIndex
                    }
                }
            }
        }
        
        if (c3pIdx < 0) {
            return pIdx
        }
        
        # Find C4'
        var endIdx = -1
        cSet = connected({atomIndex=c3pIdx}) and not {oxygen}
        for (var i = 1; i <= cSet.size; i++) {
            var ocSet = connected(cSet[i]) and {oxygen}
            for (var j = 1; j <= ocSet.size; j++) {
                if (connected(ocSet[j]).size > 1) {
                    endIdx = cSet[i].atomIndex
                }
            }
        }
        
        # Find C5'
        var c5idx = -1
        cSet = (connected({atomIndex=endIdx}) and not {atomIndex=c3pIdx}
            and not {oxygen})
        if (cSet.size > 0) {
            c5idx = cSet[1].atomIndex
        }
        
        # Find O5'
        var o5idx = -1
        cSet = connected({atomIndex=c5idx}) and {oxygen}
        if (cSet.size > 0) {
            o5idx = cSet[1].atomIndex
        }
        
        if (o5idx < 0) {
            return c5idx
        }
        
        # Find P
        pIdx = -1
        cSet = connected({atomIndex=o5idx}) and {phosphorus}
        if (cSet.size > 0) {
            Pidx =  cSet[1].atomIndex
        }
        
        if (Pidx < 0) {
            return o5idx
        }
    }
    
    return -1    
}

function find_p_idx(idx) {
    select {atomIndex=idx}
    var cSet = {selected}
    while (cSet.size > 0) {
        for (var i = 1; i <= cSet.size; i++) {
            if (cSet[i].element == "P") {
                return cSet[i].atomIndex
            }
        }
        cSet = connected({selected}) and not {selected} and not {hydrogen}
        select {selected} or cset
    }
    return -1
}

# cSet, s, and newGreek are arrays and thus passed by reference
function ring_common(cSet, nIdx, s, newGreek, nextGreek) {
    if (cSet.size > 2) {
        print format("Unrecognized structure with set %s", cSet)
    }
    var oldGreek = 0 + newGreek[1]
    newGreek[1] = nextGreek
    for (var i = 1; i <= cSet.size; i++) {
        var ccSet = connected(cSet[i])
        if (ccSet.size == 1) {
            if (cSet[i].element == ccSet[1].element) {
                s[2] = i
                s[1] = ((i > 1) ? 1 : 2)
                newGreek[i] = 1 + nextGreek
                newGreek[s[1]] = nextGreek
                return cSet[s[1]].atomIndex
            }
            else {
                s[1] = i
                s[2] = ((i > 1) ? 1 : 2)
                newGreek[i] = oldGreek
                newGreek[s[2]] = nextGreek
                return cSet[s[2]].atomIndex
            }
        }
    }
    return cSet[1].atomIndex
}

# Bound to ALT-LEFT-CLICK by plico_remap_nt
function remap_nt_cargo_mb() {
    
    var idx =_atomPicked
    if ({atomIndex=idx}.element == "H") {
        idx = connected({atomIndex=idx})[1].atomIndex
    }
    
    # If P can be found
    var pIdx = find_p_idx(idx)
    var isValid = FALSE
    var newResno = 1
    var newChain = "A"
    var newAtomno = 1
    var t5idx = -1
    if (pIdx >= 0) {
    
        t5idx = find_5_prime(pIdx)
        
        if (t5idx >= 0) {
            c5resno = {atomIndex=t5idx}.resno
            c5chain = {atomIndex=t5idx}.chain
            c5atomno = {atomIndex=t5idx}.atomno
            select {atomIndex=t5idx}
            halo on
            refresh
            
            # Prompt for new designators
            var p = prompt(("Enter 5\'-terminal atom designator as\n"
                + "   <resno>:<chain>#<atomno>"),
                format("%d:%s#%d", c5resno, c5chain, c5atomno))%0
            # If valid
            if (p != "null") {
                var iColon = p.find(":")
                if (iColon > 0) {
                    var iHash = p.find("#")
                    if (iHash > 0) {
                        newResno = 0 + (p[1][iColon-1])
                        newChain = p[iColon+1][iHash-1]
                        newAtomno = 0 + (p[iHash+1][0])
                        if ((newResno > 0)
                            and (newChain.size == 1)
                            and (newAtomno > 0)) {
                                isValid = TRUE
                        }
                    }
                }
                if (not isValid) {
                    prompt ("Entry not valid!")
                }
            }
        }
    }
    
    if (isValid) {
        background ECHO pink
        refresh

        delete {hydrogen}
        delete {hoh}
        delete %B
        delete ligands
        connect
        
        # Build inline pdb file
        ls = "data \"append remap\"\n"    # global PDB atom record
        var rs = ""
        
        select {atomIndex=t5idx}
        var cSet = {selected}
        var nextAtomName = {atomIndex=t5idx}.element
        var newGroup = "UNK"
        var newGreek = array("", "", "", "")
        var nIdx = t5idx
        var c1pIdx = -1
        var o3pIdx = -1
        var stopIdx = -1
        var endIdx = -1
        var isRNA = FALSE
        var first = TRUE
        var psu = FALSE
        while (cSet.size > 0) {
            var s = array(1, 2, 3)
            var iKeep = -1
            var iDrop = -1
            switch( nextAtomName) {
            case "O" :
                newGreek[1] = (first ? "5\'" : "P3")
                nextAtomName = (first ? "C5\'" : "P")
                nIdx = cSet[1].atomIndex
                break
            case "P" :
                newGreek[1] = ""
                nextAtomName = "OP"
                nIdx = cSet[1].atomIndex
                break
            case "OP" :
                var oc5set = ({})
                for (var i = 1; i <= cSet.size; i++) {
                    if (connected(cSet[i]).size > 1) {
                        s[3] = i
                        newGreek[i] = "5\'"
                        oc5set = connected(cSet[i]) and {carbon}
                        nIdx = cSet[i].atomIndex
                    }
                }
                for (var i = 1; i <= cSet.size; i++) {
                    if (i != s[3]) {
                        if (abs(angle(cSet[i], {atomIndex=pIdx}, cSet[s[3]],
                            oc5set[1])) < 90.0) {
                            s[1] = i
                            newGreek[i] = "P1"
                        }
                        else {
                            s[2] = i
                            newGreek[i] = "P2"
                        }
                    }
                }
                #nIdx = pIdx
                nextAtomName = "C5\'"
                break
            case "C5\'" :
                nIdx = cSet[1].atomIndex
                newGreek[1] = "5\'"
                nextAtomName = "C4\'"
                break
            case "C4\'" :
                nIdx = cSet[1].atomIndex
                newGreek[1] = "4\'"
                nextAtomName = "C3\'"
                break
            case "C3\'" :
                for (var i = 1; i <= cSet.size; i++) {
                    if (cSet[i].element == "O") {
                        s[1] = i
                        newGreek[1] = "4\'"
                        cSet[i].selected = FALSE
                        stopIdx = cSet[i].atomIndex
                    }
                    else {
                        s[2] = i
                        newGreek[2] = "3\'"
                        nIdx = cSet[i].atomIndex
                        
                    }
                }
                nextAtomName = "O3\'"
                break
            case "O3\'" :
                for (var i = 1; i <= cSet.size; i++) {
                    if (cSet[i].element == "O") {
                        s[1] = i
                        newGreek[1] = "3\'"
                        o3pIdx = cSet[i].atomIndex
                    }
                    else {
                        s[2] = i
                        newGreek[2] = "2\'"
                        nIdx = cSet[i].atomIndex
                    }
                }
                nextAtomName = "C2\'"
                break
            case "C2\'" :
                pIdx = -1
                for (var i = 1; i <= cSet.size; i++) {
                    if (cSet[i].element == "P") {
                        pIdx = cSet[i].atomIndex
                        cSet = cSet and not cSet[i]
                        break
                    }
                }
                for (var i = 1; i <= cSet.size; i++) {
                    if (cSet[i].element == "O") {
                        s[2] = i
                        newGreek[2] = "2\'"
                        isRNA = TRUE
                    }
                    else {
                        s[1] = i
                        newGreek[1] = "1\'"
                        c1pIdx = cSet[i].atomIndex
                        nIdx = cSet[i].atomIndex
                    }
                }
                nextAtomName = "C1\'"
                break                
            case "C1\'" :
                for (var i = 1; i <= cSet.size; i++) {
                    if (cSet[i].element == "N") {
                        iKeep = i
                        nIdx = cSet[i].atomIndex
                    }
                    else if ((cSet[i].element == "C") and 
                        ((connected(cSet[i]) and {oxygen}).size == 0)) { #PSU
                        psu = TRUE
                        iKeep = i
                        nIdx = cSet[i].atomIndex
                    }
                    else {
                        cSet[i].selected = FALSE
                    }
                }
                cSet = cSet[iKeep]
                var ccSet = connected(cSet[1]) and not {atomIndex=c1pIdx}
                newGreek[1] = "9"
                nextAtomName = "N9u"
                newGroup = "PU"
                for (var j = 1; j <= ccSet.size; j++) {
                    if ((connected(ccSet[j]) and {oxygen}) > 0) {
                        newGreek[1] = "1"
                        nextAtomName = "N1y"
                        newGroup = "PY"
                    }
                }
                break                
            case "N1y" :
                for (var i = 1; i <= cSet.size; i++) {
                    if (connected(cSet[i]) > 2) {
                        iKeep = i
                        nIdx = cSet[i].atomIndex
                    }
                    else {
                        stopIdx = cSet[i].atomIndex
                    }
                }
                cSet = cSet[iKeep]
                newGreek[1] = "2"
                nextAtomName = "C2"
                break                
            case "N9u" :
                # Find N-C-N-C-N  
                for (var i = 1; i <= cSet.size; i++) {
                    n1atom = (connected(cSet[i]) and {nitrogen}
                        and not {atomIndex=nIdx})
                    c2set = connected(n1atom) and {carbon} and not cSet[i]
                    for (var j = 1; j <= c2set.size; j++) {
                        if ((connected(c2set[j]) and {nitrogen}) > 1) {
                            iDrop = i
                        }
                    }
                }
                stopIdx = cSet[iDrop].atomIndex
                cSet = cSet and not cSet[iDrop]
                nIdx = cSet[1].atomIndex
                newGreek[1] = "8"
                nextAtomName = "C8"
                break                
            case "C8" :
                nIdx = ring_common( cSet, nIdx, s, newGreek, "7")
                nextAtomName = "N7"
                break                
            case "N7" :
                nIdx = ring_common( cSet, nIdx, s, newGreek, "5")
                nextAtomName = "C5"
                break                
            case "C5" :
                if (isRNA and (newGroup == "DU ")) {
                    var c5set = {atomIndex=nIdx} or connected({atomIndex=nIdx})
                    if (angle(c5set[1], c5set[2], c5set[3]) < 114.0) {
                        newGroup = "D  "
                    } 
                }                
                nIdx = ring_common( cSet, nIdx, s, newGreek, "6")
                if ((newGroup == "DU ") and (cSet.size > 1)) {
                    newGroup = "DT "
                }
                nextAtomName = "C6"
                break                
            case "C6" :
                if (newGroup == "PU") {
                    nIdx = ring_common( cSet, nIdx, s, newGreek, "1")
                    newGroup = ((cSet[1].element == "O") ? "DG " : "DA ")
                    nextAtomName = "N1"
                }
                else {
                    if (psu) {
                        psu = FALSE
                        newGroup = "DU "
                    }                
                    cSet = ({})
                }
                break                
            case "N1" :
                if (connected({atomIndex=nIdx}).size > 2) { # YG
                    newGroup = "X  "                
                }
                nIdx = ring_common( cSet, nIdx, s, newGreek, "2")
                nextAtomName = "C2"
                break                
            case "C2" :
                nIdx = ring_common( cSet, nIdx, s, newGreek, "3")
                nextAtomName = "N3"
                stopIdx = -1
                break                
            case "N3" :
                nIdx = ring_common( cSet, nIdx, s, newGreek, "4")
                nextAtomName = "C4"
                break                
            case "C4" :
                if (newGroup != "PY") {
                    cSet = ({})
                }
                else {
                    nIdx = ring_common( cSet, nIdx, s, newGreek, "5")
                    newGroup = ((cSet[1].element == "N") ? "DC " : "DU ")
                    nextAtomName = "C5"
                }
                break
            }
            first = FALSE
            
            for (var i = 1; i <= cSet.size; i++) {
                rs += format("ATOM  %5d  %-4sUNK ", newAtomNo,
                    (cSet[s[i]].element + newGreek[s[i]]))
                rs += format("%s%4d    %8.3f", newChain, newResno, cSet[s[i]].x)          
                rs += format("%8.3f%8.3f\n", cSet[s[i]].y, cSet[s[i]].z)
                newAtomno++
            }
            
            cSet = (connected(cSet and not {atomIndex=stopIdx})
                and not cSet and not {atomIndex=stopIdx} and not {atomIndex=endIdx})
            endIdx = nIdx
            
            if (cSet.size == 0) {
                if (isRNA) {
                    newGroup = (newGroup.replace("DA ","A  ").replace("DG ","G  ")
                        .replace("DC ","C  ").replace("DT ","T  ").replace("DU ","U  "))
                }     
                ls += rs.replace("UNK", newGroup)
                rs = ""
                newResno++
                
                if (pIdx >= 0) {
                    cSet = {atomIndex=pIdx}
                    nextAtomName = "P"
                    newGroup = "UNK"
                    greek = ""
                    newGreek[1] = greek
                    c1pIdx = -1
                    stopIdx = o3pIdx
                    endIdx = -1
                    isRNA = FALSE
                }
                else {
                    break
                }
            }
        } # endwhile 

        # Replace chain with new chain
        cset = {atomIndex=idx}
        select cSet
        while (cSet.size > 0) {
            cSet = connected({selected}) and not {selected}
            select {selected} or cSet
        }
        delete {selected}
        
        ls += "end \"append remap\""
        script inline @{ls}
        gEcho += format("|Chain %s has been rebuilt", newChain)
        set echo TOP LEFT
        echo @gEcho
        background ECHO yellow
        refresh
        print_1c_chain( newChain)
    }
}

function print_1c_chain(iChain) {
    var resmin = {chain=iChain}.resno.min
    var resmax = {chain=iChain}.resno.max
    var rchar = (({(resno=resmin) and (chain=iChain)}.group[0].size > 1) ? "" : "R")
    var lcAtoms = (within(3.1, FALSE, {(resno=resmin) and (chain=iChain) and base})
        and not {(resno=resmin) and (chain=iChain)})
    var chain2 = ""
    var schar = "S"
    if (lcAtoms.size > 0) {
        chain2 = lcAtoms[1].chain
        if (((rchar == "R") and (lcAtoms[1].group.size > 1))
            or ((rchar == "") and (lcAtoms[1].group.size == 1))) {
            schar = "M"
        }
        else {
            schar = ""
        }
    }
    var ls = format("%s%s:%s", iChain, chain2, format("%s%s", rchar, schar))
    for (var i = {chain=iChain}.resno.min; i <= {chain=iChain}.resno.max; i++) {
        ls += ({(resno=i) and (chain=iChain)}.group[0])[0]
    }
    print ls
}

# Top level of Remap
function plico_remap_nt() {
    
    # Push selected
    gSelSaves = {selected}
    
    gAppendNew = appendNew
    set appendNew FALSE
    gScheme = defaultColorScheme
    gAltScheme = ((gScheme == "Jmol") ? "Rasmol" : "Jmol")
    set echo TOP LEFT
    background ECHO yellow
    gEcho = "_____REMAP NT_____|ALT-CLICK=select NT chain|DOUBLE-CLICK=exit"
    echo @gEcho
    gChain = ""
    unbind

    set picking ON
    bind "ALT-LEFT-CLICK" "_pickAtom";
    bind "ALT-LEFT-CLICK" "+:remap_nt_cargo_mb";
    bind "DOUBLE" "remap_nt_exit";
}

# Bound to DOUBLE by plicoRemap
function remap_nt_exit() {
    unbind
    halo off 
    echo
    select all
    color {selected} @gScheme
    gBusy = FALSE
    set appendNew gAppendNew
    
    # Pop selected
    select gSelSaves
}

# End of REMAPNT.SPT

Contributors

Remig