Search in sources :

Example 1 with Group

use of org.biojava.bio.structure.Group in project ffx by mjschnie.

the class PDBFileMatcher method compareAtoms.

/**
 * Mildly robust atom comparison tool which attempts to account for things
 * like missing hydrogens causing inequal PDB serial numbering. Compares on
 * element, group PDB code, atom serial (if false, also checks group
 * number), and atom name (currently only accepts H-H1 mismatches).
 *
 * @param a1
 * @param a2
 * @return if considered equivalent.
 */
private boolean compareAtoms(Atom a1, Atom a2) {
    if (!a1.getElement().equals(a2.getElement())) {
        return false;
    }
    Group group1 = a1.getGroup();
    Group group2 = a2.getGroup();
    if (!group1.getPDBName().equalsIgnoreCase(group2.getPDBName())) {
        return false;
    }
    if (a1.getPDBserial() != a2.getPDBserial()) {
        // Second check accounts for situations where one structure may be H-trimmed.
        if (!group1.getResidueNumber().equals(group2.getResidueNumber())) {
            return false;
        }
    }
    String a1Name = a1.getName();
    String a2Name = a2.getName();
    if (!a1Name.equalsIgnoreCase(a2Name)) {
        // Courtesy of MSMBuilder, I have to keep track of Hs which should be H1.
        if (!((a1Name.equalsIgnoreCase("H") && a2Name.equalsIgnoreCase("H1")) || (a2Name.equalsIgnoreCase("H") && a1Name.equalsIgnoreCase("H1")))) {
            return false;
        }
    }
    return true;
}
Also used : Group(org.biojava.bio.structure.Group)

Example 2 with Group

use of org.biojava.bio.structure.Group in project ffx by mjschnie.

the class PDBFileMatcher method getMatchingAtom.

/**
 * Finds atom1's match in structure2; uses atom1's residue number, atom
 * type, and PDB serial number as a guess; if robustMatch is set true, will
 * search all Atoms in structure2.
 *
 * @param atom1 An Atom
 * @param structure2 A Structure to search
 * @param searchAll Whether to search all Atoms in structure2.
 * @return atom1's match in structure2
 * @throws IllegalArgumentException If no match can be found.
 */
private Atom getMatchingAtom(Atom atom1, Structure structure2, boolean searchAll) throws IllegalArgumentException {
    ResidueNumber res1 = atom1.getGroup().getResidueNumber();
    String chainID = res1.getChainId();
    Atom atom2 = null;
    try {
        Chain chain2 = structure2.getChainByPDB(chainID);
        Group group2 = chain2.getGroupByPDB(res1);
        try {
            atom2 = group2.getAtom(atom1.getName());
        } catch (StructureException ex) {
            if (atom1.getName().equalsIgnoreCase("H")) {
                atom2 = group2.getAtom("H1");
            } else if (atom1.getName().equalsIgnoreCase("H1")) {
                atom2 = group2.getAtom("H");
            }
            atom2 = group2.getAtom(atom1.getPDBserial());
        }
    } catch (StructureException ex) {
        if (!searchAll) {
            throw new IllegalArgumentException("Matching atom not found.");
        }
        for (Chain chain : structure2.getChains()) {
            for (Group group : chain.getAtomGroups()) {
                for (Atom atom : group.getAtoms()) {
                    if (compareAtoms(atom1, atom)) {
                        return atom2;
                    }
                }
            }
        }
    }
    if (atom2 != null && compareAtoms(atom1, atom2)) {
        return atom2;
    }
    throw new IllegalArgumentException("Matching atom not found.");
}
Also used : Chain(org.biojava.bio.structure.Chain) Group(org.biojava.bio.structure.Group) StructureException(org.biojava.bio.structure.StructureException) ResidueNumber(org.biojava.bio.structure.ResidueNumber) Atom(org.biojava.bio.structure.Atom)

Example 3 with Group

use of org.biojava.bio.structure.Group in project ffx by mjschnie.

the class PDBFileMatcher method calculateRMSD.

private double calculateRMSD(FileFilePair matchFile, Structure sourceStructure, StructurePairAligner aligner) throws StructureException {
    Structure matchStructure = matchFile.getStructure();
    if (superpose) {
        double rmsd = Double.MAX_VALUE;
        aligner.align(matchStructure, matchStructure);
        AlternativeAlignment[] alignments = aligner.getAlignments();
        for (AlternativeAlignment alignment : alignments) {
            double alignRMSD = alignment.getRmsd();
            rmsd = alignRMSD < rmsd ? alignRMSD : rmsd;
        }
        return rmsd;
    }
    Atom[] matchArray;
    Atom[] sourceArray;
    switch(atomsUsed) {
        case 1:
            String[] atomNames = { "CA", "N1", "N9" };
            matchArray = StructureTools.getAtomArray(matchStructure, atomNames);
            sourceArray = StructureTools.getAtomArray(sourceStructure, atomNames);
            break;
        case 2:
            List<Atom> matchAtoms = new ArrayList<>();
            List<Chain> matchChains = matchStructure.getChains();
            for (Chain chain : matchChains) {
                List<Group> matchGroups = chain.getAtomGroups(GroupType.AMINOACID);
                matchGroups.addAll(chain.getAtomGroups(GroupType.NUCLEOTIDE));
                for (Group group : matchGroups) {
                    matchAtoms.addAll(group.getAtoms());
                }
            }
            matchArray = matchAtoms.toArray(new Atom[matchAtoms.size()]);
            List<Atom> sourceAtoms = new ArrayList<>();
            List<Chain> sourceChains = sourceStructure.getChains();
            for (Chain chain : sourceChains) {
                List<Group> sourceGroups = chain.getAtomGroups(GroupType.AMINOACID);
                sourceGroups.addAll(chain.getAtomGroups(GroupType.NUCLEOTIDE));
                for (Group group : sourceGroups) {
                    sourceAtoms.addAll(group.getAtoms());
                }
            }
            sourceArray = sourceAtoms.toArray(new Atom[sourceAtoms.size()]);
            break;
        case 3:
            matchAtoms = new ArrayList<>();
            matchChains = matchStructure.getChains();
            for (Chain chain : matchChains) {
                List<Group> matchGroups = chain.getAtomGroups();
                for (Group group : matchGroups) {
                    if (!group.isWater()) {
                        matchAtoms.addAll(group.getAtoms());
                    }
                }
            }
            matchArray = matchAtoms.toArray(new Atom[matchAtoms.size()]);
            sourceAtoms = new ArrayList<>();
            sourceChains = sourceStructure.getChains();
            for (Chain chain : sourceChains) {
                List<Group> matchGroups = chain.getAtomGroups();
                for (Group group : matchGroups) {
                    if (!group.isWater()) {
                        sourceAtoms.addAll(group.getAtoms());
                    }
                }
            }
            sourceArray = sourceAtoms.toArray(new Atom[sourceAtoms.size()]);
            break;
        case 4:
            matchArray = StructureTools.getAllAtomArray(matchStructure);
            sourceArray = StructureTools.getAllAtomArray(sourceStructure);
            break;
        default:
            throw new IllegalArgumentException("atomsUsed is not 1-4: this has not been properly checked at an earlier stage.");
    }
    /*List<Atom> matchAtoms = new ArrayList<>();
         Structure matchStructure = matchFile.getStructure();
         List<Chain> matchChains = matchStructure.getChains();
         if (atomsUsed == 4) {
         matchAtoms = StructureTools.getAllAtomArray(matchStructure);
         }
         for (Chain chain : matchChains) {
         List<Group> matchGroups = chain.getSeqResGroups();
         for (Group group : matchGroups) {
         String groupType = group.getType();
         if (groupType.equalsIgnoreCase("HETATOM")) {
         if (atomsUsed < 3) {
         continue;
         } else if (atomsUsed < 4 && group.isWater()) {
         continue;
         }
         }
         if (atomsUsed != 1) {
         matchAtoms.addAll(group.getAtoms());
         } else {
         try {
         matchAtoms.add(getReferenceAtom(group));
         } catch (StructureException ex) {
         String refAtomType = (groupType.equalsIgnoreCase("AminoAcid") ? "CA" : "N1/9");
         // refAtomType should be binary between amino acid and nucleic acid, as HETATOM is eliminated.
         logger.info(String.format(" Reference atom %s could not be found for group %s",
         refAtomType, group.toString()));
         }
         }
         }
         }
         Atom[] match = new Atom[matchAtoms.size()];
         matchAtoms.toArray(match);
         matchAtoms.clear();

         List<Atom> sourceAtoms = new ArrayList<>();
         List<Chain> sourceChains = sourceStructure.getChains();
         for (Chain chain : sourceChains) {
         List<Group> sourceGroups = chain.getSeqResGroups();
         for (Group group : sourceGroups) {
         String groupType = group.getType();
         if (groupType.equalsIgnoreCase("HETATOM")) {
         if (atomsUsed < 3) {
         continue;
         } else if (atomsUsed < 4 && group.isWater()) {
         continue;
         }
         }
         if (atomsUsed != 1) {
         sourceAtoms.addAll(group.getAtoms());
         } else {
         try {
         sourceAtoms.add(getReferenceAtom(group));
         } catch (StructureException ex) {
         logger.info(String.format(" Reference atom could not be found for group %s", group.toString()));
         }
         }
         }
         }
         Atom[] sourceArray = new Atom[sourceAtoms.size()];
         sourceAtoms.toArray(sourceArray);
         sourceAtoms.clear();*/
    return simpleRMSD(sourceArray, matchArray);
}
Also used : Chain(org.biojava.bio.structure.Chain) Group(org.biojava.bio.structure.Group) ArrayList(java.util.ArrayList) Atom(org.biojava.bio.structure.Atom) AlternativeAlignment(org.biojava.bio.structure.align.pairwise.AlternativeAlignment) Structure(org.biojava.bio.structure.Structure)

Example 4 with Group

use of org.biojava.bio.structure.Group in project ffx by mjschnie.

the class BiojavaFilter method convert.

/**
 * {@inheritDoc}
 *
 * Parse the Biojava Structure
 *
 * @return true if the structure is successfully converted
 */
@Override
public boolean convert() {
    // First atom is #1, to match xyz file format
    int xyzIndex = 1;
    setConverted(false);
    systems.add(activeMolecularAssembly);
    if (mutate) {
        List<Character> chainIDs = new ArrayList<>();
        for (Chain chain : structure.getChains()) {
            chainIDs.add(chain.getChainID().charAt(0));
        }
        if (!chainIDs.contains(mutateChainID)) {
            if (chainIDs.size() == 1) {
                logger.warning(String.format(" Chain ID %c for " + "mutation not found: only one chain %c " + "found.", mutateChainID, chainIDs.get(0)));
                mutateChainID = chainIDs.get(0);
            } else {
                logger.warning(String.format(" Chain ID %c for " + "mutation not found: mutation will not " + "proceed.", mutateChainID));
            }
        }
    }
    /**
     * Echo the alternate location being parsed.
     */
    if (currentAltLoc == 'A') {
        logger.info(String.format(" Reading %s", structure.getName()));
    } else {
        logger.info(String.format(" Reading %s alternate location %s", structure.getName(), currentAltLoc));
    }
    org.biojava.bio.structure.Atom[] bjAtoms = StructureTools.getAllAtomArray(structure);
    int nAtoms = bjAtoms.length;
    Atom[] ffxAtoms = new Atom[nAtoms];
    /**
     * Reset the current chain and segID.
     */
    currentChainID = null;
    currentSegID = null;
    PDBCrystallographicInfo cInfo = structure.getCrystallographicInfo();
    if (cInfo.isCrystallographic()) {
        // I do not think we need to check if it already has these properties,
        // but it can be done.
        properties.addProperty("a-axis", cInfo.getA());
        properties.addProperty("b-axis", cInfo.getB());
        properties.addProperty("c-axis", cInfo.getC());
        properties.addProperty("alpha", cInfo.getAlpha());
        properties.addProperty("beta", cInfo.getBeta());
        properties.addProperty("gamma", cInfo.getGamma());
        properties.addProperty("spacegroup", SpaceGroup.pdb2ShortName(cInfo.getSpaceGroup()));
    }
    for (org.biojava.bio.structure.Atom atom : bjAtoms) {
        String name = atom.getName().toUpperCase().trim();
        double[] xyz = new double[3];
        xyz[0] = atom.getX();
        xyz[1] = atom.getY();
        xyz[2] = atom.getZ();
        char altLoc = atom.getAltLoc();
        if (!altLocs.contains(altLoc)) {
            altLocs.add(altLoc);
        }
        if (altLoc != ' ' && altLoc != 'A' && altLoc != currentAltLoc) {
            break;
        }
        if (name.contains("1H") || name.toUpperCase().contains("2H") || name.toUpperCase().contains("3H")) {
            // VERSION3_2 is presently just a placeholder for "anything non-standard".
            fileStandard = VERSION3_2;
        }
        Group group = atom.getGroup();
        ResidueNumber resnum = group.getResidueNumber();
        int resSeq = resnum.getSeqNum();
        String resName = group.getPDBName().trim().toUpperCase();
        Chain chain = group.getChain();
        char chainID = chain.getChainID().charAt(0);
        String segID = getSegID(chainID);
        boolean printAtom = false;
        if (mutate && chainID == mutateChainID && mutateResID == resSeq) {
            if (name.equals("N") || name.equals("C") || name.equals("O") || name.equals("CA")) {
                printAtom = true;
                name = mutateToResname;
            } else {
                logger.info(String.format(" Deleting atom %s of %s %d", name, resName, resSeq));
                break;
            }
        }
        Atom newAtom = new Atom(0, name, altLoc, xyz, resName, resSeq, chainID, atom.getOccupancy(), atom.getTempFactor(), segID);
        /* Biojava sets at least some capping groups, and possibly nonstandard
             amino acids to be heteroatoms. */
        boolean hetatm = true;
        for (AminoAcid3 aa3Name : AminoAcid3.values()) {
            if (aa3Name.name().equals(resName)) {
                hetatm = false;
                break;
            }
        }
        newAtom.setHetero(hetatm);
        // Look Ma, Biojava doesn't care about anisou!
        Atom returnedAtom = (Atom) activeMolecularAssembly.addMSNode(newAtom);
        if (returnedAtom != newAtom) {
            atoms.put(atom.getPDBserial(), returnedAtom);
            if (logger.isLoggable(Level.FINE)) {
                logger.fine(String.format("%s has been retained over\n%s", returnedAtom.toString(), newAtom.toString()));
            }
        } else {
            atoms.put(atom.getPDBserial(), newAtom);
            if (newAtom.getIndex() == 0) {
                newAtom.setXyzIndex(xyzIndex++);
            }
            if (printAtom) {
                logger.info(newAtom.toString());
            }
        }
    }
    List<Bond> ssBondList = new ArrayList<>();
    for (SSBond ssBond : structure.getSSBonds()) {
        Polymer c1 = activeMolecularAssembly.getChain(ssBond.getChainID1());
        Polymer c2 = activeMolecularAssembly.getChain(ssBond.getChainID2());
        int rn1;
        int rn2;
        try {
            rn1 = Integer.parseInt(ssBond.getResnum1());
            rn2 = Integer.parseInt(ssBond.getResnum1());
        } catch (NumberFormatException ex) {
            logger.warning(String.format(" Could not parse SSbond %d", ssBond.getSerNum()));
            continue;
        }
        Residue r1 = c1.getResidue(rn1);
        Residue r2 = c2.getResidue(rn2);
        List<Atom> atoms1 = r1.getAtomList();
        List<Atom> atoms2 = r2.getAtomList();
        Atom SG1 = null;
        Atom SG2 = null;
        for (Atom atom : atoms1) {
            if (atom.getName().equalsIgnoreCase("SG")) {
                SG1 = atom;
                break;
            }
        }
        for (Atom atom : atoms2) {
            if (atom.getName().equalsIgnoreCase("SG")) {
                SG2 = atom;
                break;
            }
        }
        if (SG1 == null) {
            logger.warning(String.format(" SG atom 1 of SS-bond %s is null", ssBond.toString()));
        }
        if (SG2 == null) {
            logger.warning(String.format(" SG atom 2 of SS-bond %s is null", ssBond.toString()));
        }
        if (SG1 == null || SG2 == null) {
            continue;
        }
        double d = VectorMath.dist(SG1.getXYZ(null), SG2.getXYZ(null));
        if (d < 3.0) {
            r1.setName("CYX");
            r2.setName("CYX");
            for (Atom atom : atoms1) {
                atom.setResName("CYX");
            }
            for (Atom atom : atoms2) {
                atom.setResName("CYX");
            }
            Bond bond = new Bond(SG1, SG2);
            ssBondList.add(bond);
        } else {
            String message = String.format("Ignoring [%s]\n due to distance %8.3f A.", ssBond.toString(), d);
            logger.log(Level.WARNING, message);
        }
    }
    int pdbAtoms = activeMolecularAssembly.getAtomArray().length;
    assignAtomTypes();
    StringBuilder sb = new StringBuilder(" Disulfide Bonds:");
    for (Bond bond : ssBondList) {
        Atom a1 = bond.getAtom(0);
        Atom a2 = bond.getAtom(1);
        int[] c = new int[2];
        c[0] = a1.getAtomType().atomClass;
        c[1] = a2.getAtomType().atomClass;
        String key = ffx.potential.parameters.BondType.sortKey(c);
        ffx.potential.parameters.BondType bondType = forceField.getBondType(key);
        if (bondType == null) {
            logger.severe(String.format("No BondType for key: %s\n %s\n %s", key, a1, a2));
        } else {
            bond.setBondType(bondType);
        }
        double d = VectorMath.dist(a1.getXYZ(null), a2.getXYZ(null));
        Polymer c1 = activeMolecularAssembly.getChain(a1.getSegID());
        Polymer c2 = activeMolecularAssembly.getChain(a2.getSegID());
        Residue r1 = c1.getResidue(a1.getResidueNumber());
        Residue r2 = c2.getResidue(a2.getResidueNumber());
        sb.append(String.format("\n S-S distance of %6.2f for %s and %s.", d, r1.toString(), r2.toString()));
        bondList.add(bond);
    }
    if (ssBondList.size() > 0) {
        logger.info(sb.toString());
    }
    /**
     * Finally, re-number the atoms if missing atoms were created.
     */
    if (pdbAtoms != activeMolecularAssembly.getAtomArray().length) {
        numberAtoms();
    }
    return true;
}
Also used : Chain(org.biojava.bio.structure.Chain) Group(org.biojava.bio.structure.Group) MSGroup(ffx.potential.bonded.MSGroup) SpaceGroup(ffx.crystal.SpaceGroup) ArrayList(java.util.ArrayList) PDBCrystallographicInfo(org.biojava.bio.structure.PDBCrystallographicInfo) SSBond(org.biojava.bio.structure.SSBond) AminoAcid3(ffx.potential.bonded.ResidueEnumerations.AminoAcid3) Polymer(ffx.potential.bonded.Polymer) Atom(ffx.potential.bonded.Atom) Residue(ffx.potential.bonded.Residue) ResidueNumber(org.biojava.bio.structure.ResidueNumber) Bond(ffx.potential.bonded.Bond) SSBond(org.biojava.bio.structure.SSBond)

Aggregations

Group (org.biojava.bio.structure.Group)4 Chain (org.biojava.bio.structure.Chain)3 ArrayList (java.util.ArrayList)2 Atom (org.biojava.bio.structure.Atom)2 ResidueNumber (org.biojava.bio.structure.ResidueNumber)2 SpaceGroup (ffx.crystal.SpaceGroup)1 Atom (ffx.potential.bonded.Atom)1 Bond (ffx.potential.bonded.Bond)1 MSGroup (ffx.potential.bonded.MSGroup)1 Polymer (ffx.potential.bonded.Polymer)1 Residue (ffx.potential.bonded.Residue)1 AminoAcid3 (ffx.potential.bonded.ResidueEnumerations.AminoAcid3)1 PDBCrystallographicInfo (org.biojava.bio.structure.PDBCrystallographicInfo)1 SSBond (org.biojava.bio.structure.SSBond)1 Structure (org.biojava.bio.structure.Structure)1 StructureException (org.biojava.bio.structure.StructureException)1 AlternativeAlignment (org.biojava.bio.structure.align.pairwise.AlternativeAlignment)1