Examples with Atom ffx.potential.bonded.Atom used on opensource projects

Example 81 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class PDBFilter method assignAtomTypes.

/**
 * Assign force field atoms types to common chemistries using "biotype"
 * records.
 */
private void assignAtomTypes() {
    /**
     * Create a list to store bonds defined by PDB atom names.
     */
    bondList = new ArrayList<>();
    /**
     * To Do: Look for cyclic peptides and disulfides.
     */
    Polymer[] polymers = activeMolecularAssembly.getChains();
    /**
     * Loop over chains.
     */
    if (polymers != null) {
        logger.info(format("\n Assigning atom types for %d chains.", polymers.length));
        for (Polymer polymer : polymers) {
            List<Residue> residues = polymer.getResidues();
            int numberOfResidues = residues.size();
            /**
             * Check if all residues are known amino acids.
             */
            boolean isProtein = true;
            for (int residueNumber = 0; residueNumber < numberOfResidues; residueNumber++) {
                Residue residue = residues.get(residueNumber);
                String name = residue.getName().toUpperCase();
                boolean aa = false;
                for (AminoAcid3 amino : aminoAcidList) {
                    if (amino.toString().equalsIgnoreCase(name)) {
                        aa = true;
                        checkHydrogenAtomNames(residue);
                        break;
                    }
                }
                // Check for a patch.
                if (!aa) {
                    logger.info(" Checking for non-standard amino acid patch " + name);
                    HashMap<String, AtomType> types = forceField.getAtomTypes(name);
                    if (types.isEmpty()) {
                        isProtein = false;
                        break;
                    } else {
                        logger.info(" Patch found for non-standard amino acid " + name);
                    }
                }
            }
            /**
             * If all the residues in this chain have known amino acids
             * names, then attempt to assign atom types.
             */
            if (isProtein) {
                try {
                    logger.info(format(" Amino acid chain %s", polymer.getName()));
                    double dist = properties.getDouble("chainbreak", 3.0);
                    // Detect main chain breaks!
                    List<List<Residue>> subChains = findChainBreaks(residues, dist);
                    for (List<Residue> subChain : subChains) {
                        assignAminoAcidAtomTypes(subChain);
                    }
                } catch (MissingHeavyAtomException missingHeavyAtomException) {
                    logger.severe(missingHeavyAtomException.toString());
                } catch (MissingAtomTypeException missingAtomTypeException) {
                    logger.severe(missingAtomTypeException.toString());
                }
                continue;
            }
            /**
             * Check if all residues have known nucleic acids names.
             */
            boolean isNucleicAcid = true;
            for (int residueNumber = 0; residueNumber < numberOfResidues; residueNumber++) {
                Residue residue = residues.get(residueNumber);
                String name = residue.getName().toUpperCase();
                /**
                 * Convert 1 and 2-character nucleic acid names to
                 * 3-character names.
                 */
                if (name.length() == 1) {
                    if (name.equals("A")) {
                        name = NucleicAcid3.ADE.toString();
                    } else if (name.equals("C")) {
                        name = NucleicAcid3.CYT.toString();
                    } else if (name.equals("G")) {
                        name = NucleicAcid3.GUA.toString();
                    } else if (name.equals("T")) {
                        name = NucleicAcid3.THY.toString();
                    } else if (name.equals("U")) {
                        name = NucleicAcid3.URI.toString();
                    }
                } else if (name.length() == 2) {
                    if (name.equals("YG")) {
                        name = NucleicAcid3.YYG.toString();
                    }
                }
                residue.setName(name);
                NucleicAcid3 nucleicAcid = null;
                for (NucleicAcid3 nucleic : nucleicAcidList) {
                    String nuc3 = nucleic.toString();
                    nuc3 = nuc3.substring(nuc3.length() - 3);
                    if (nuc3.equalsIgnoreCase(name)) {
                        nucleicAcid = nucleic;
                        break;
                    }
                }
                if (nucleicAcid == null) {
                    logger.info(format("Nucleic acid was not recognized %s.", name));
                    isNucleicAcid = false;
                    break;
                }
            }
            /**
             * If all the residues in this chain have known nucleic acids
             * names, then attempt to assign atom types.
             */
            if (isNucleicAcid) {
                try {
                    logger.info(format(" Nucleic acid chain %s", polymer.getName()));
                    assignNucleicAcidAtomTypes(residues, forceField, bondList);
                } catch (MissingHeavyAtomException | MissingAtomTypeException e) {
                    logger.severe(e.toString());
                }
            }
        }
    }
    // Assign ion atom types.
    ArrayList<MSNode> ions = activeMolecularAssembly.getIons();
    if (ions != null && ions.size() > 0) {
        logger.info(format(" Assigning atom types for %d ions.", ions.size()));
        for (MSNode m : ions) {
            Molecule ion = (Molecule) m;
            String name = ion.getResidueName().toUpperCase();
            HetAtoms hetatm = HetAtoms.valueOf(name);
            Atom atom = ion.getAtomList().get(0);
            if (ion.getAtomList().size() != 1) {
                logger.severe(format(" Check residue %s of chain %s.", ion.toString(), ion.getChainID()));
            }
            try {
                switch(hetatm) {
                    case NA:
                        atom.setAtomType(findAtomType(2003));
                        break;
                    case K:
                        atom.setAtomType(findAtomType(2004));
                        break;
                    case MG:
                    case MG2:
                        atom.setAtomType(findAtomType(2005));
                        break;
                    case CA:
                    case CA2:
                        atom.setAtomType(findAtomType(2006));
                        break;
                    case CL:
                        atom.setAtomType(findAtomType(2007));
                        break;
                    case ZN:
                    case ZN2:
                        atom.setAtomType(findAtomType(2008));
                        break;
                    case BR:
                        atom.setAtomType(findAtomType(2009));
                        break;
                    default:
                        logger.severe(format(" Check residue %s of chain %s.", ion.toString(), ion.getChainID()));
                }
            } catch (Exception e) {
                String message = "Error assigning atom types.";
                logger.log(Level.SEVERE, message, e);
            }
        }
    }
    // Assign water atom types.
    ArrayList<MSNode> water = activeMolecularAssembly.getWaters();
    if (water != null && water.size() > 0) {
        logger.info(format(" Assigning atom types for %d waters.", water.size()));
        for (MSNode m : water) {
            Molecule wat = (Molecule) m;
            try {
                Atom O = buildHeavy(wat, "O", null, 2001);
                Atom H1 = buildHydrogen(wat, "H1", O, 0.96e0, null, 109.5e0, null, 120.0e0, 0, 2002);
                H1.setHetero(true);
                Atom H2 = buildHydrogen(wat, "H2", O, 0.96e0, H1, 109.5e0, null, 120.0e0, 0, 2002);
                H2.setHetero(true);
            } catch (Exception e) {
                String message = "Error assigning atom types to a water.";
                logger.log(Level.SEVERE, message, e);
            }
        }
    }
    // Assign small molecule atom types.
    ArrayList<Molecule> molecules = activeMolecularAssembly.getMolecules();
    for (MSNode m : molecules) {
        Molecule molecule = (Molecule) m;
        String moleculeName = molecule.getResidueName();
        logger.info(" Attempting to patch " + moleculeName);
        ArrayList<Atom> moleculeAtoms = molecule.getAtomList();
        boolean patched = true;
        HashMap<String, AtomType> types = forceField.getAtomTypes(moleculeName);
        /**
         * Assign atom types for all known atoms.
         */
        for (Atom atom : moleculeAtoms) {
            String atomName = atom.getName().toUpperCase();
            AtomType atomType = types.get(atomName);
            if (atomType == null) {
                logger.info(" No atom type was found for " + atomName + " of " + moleculeName + ".");
                patched = false;
                break;
            } else {
                logger.fine(" " + atom.toString() + " -> " + atomType.toString());
                atom.setAtomType(atomType);
                types.remove(atomName);
            }
        }
        /**
         * Create missing hydrogen atoms. Check for missing heavy atoms.
         */
        if (patched && !types.isEmpty()) {
            for (AtomType type : types.values()) {
                if (type.atomicNumber != 1) {
                    logger.info(" Missing heavy atom " + type.name);
                    patched = false;
                    break;
                }
            }
        }
        // Create bonds between known atoms.
        if (patched) {
            for (Atom atom : moleculeAtoms) {
                String atomName = atom.getName();
                String[] bonds = forceField.getBonds(moleculeName, atomName);
                if (bonds != null) {
                    for (String name : bonds) {
                        Atom atom2 = molecule.getAtom(name);
                        if (atom2 != null && !atom.isBonded(atom2)) {
                            buildBond(atom, atom2);
                        }
                    }
                }
            }
        }
        // Create missing hydrogen atoms.
        if (patched && !types.isEmpty()) {
            // Create a hashmap of the molecule's atoms
            HashMap<String, Atom> atomMap = new HashMap<String, Atom>();
            for (Atom atom : moleculeAtoms) {
                atomMap.put(atom.getName().toUpperCase(), atom);
            }
            for (String atomName : types.keySet()) {
                AtomType type = types.get(atomName);
                String[] bonds = forceField.getBonds(moleculeName, atomName.toUpperCase());
                if (bonds == null || bonds.length != 1) {
                    patched = false;
                    logger.info(" Check biotype for hydrogen " + type.name + ".");
                    break;
                }
                // Get the heavy atom the hydrogen is bonded to.
                Atom ia = atomMap.get(bonds[0].toUpperCase());
                Atom hydrogen = new Atom(0, atomName, ia.getAltLoc(), new double[3], ia.getResidueName(), ia.getResidueNumber(), ia.getChainID(), ia.getOccupancy(), ia.getTempFactor(), ia.getSegID());
                logger.fine(" Created hydrogen " + atomName + ".");
                hydrogen.setAtomType(type);
                hydrogen.setHetero(true);
                molecule.addMSNode(hydrogen);
                int valence = ia.getAtomType().valence;
                List<Bond> aBonds = ia.getBonds();
                int numBonds = aBonds.size();
                /**
                 * Try to find the following configuration: ib-ia-ic
                 */
                Atom ib = null;
                Atom ic = null;
                Atom id = null;
                if (numBonds > 0) {
                    Bond bond = aBonds.get(0);
                    ib = bond.get1_2(ia);
                }
                if (numBonds > 1) {
                    Bond bond = aBonds.get(1);
                    ic = bond.get1_2(ia);
                }
                if (numBonds > 2) {
                    Bond bond = aBonds.get(2);
                    id = bond.get1_2(ia);
                }
                /**
                 * Building the hydrogens depends on hybridization and the
                 * locations of other bonded atoms.
                 */
                logger.fine(" Bonding " + atomName + " to " + ia.getName() + " (" + numBonds + " of " + valence + ").");
                switch(valence) {
                    case 4:
                        switch(numBonds) {
                            case 3:
                                // Find the average coordinates of atoms ib, ic and id.
                                double[] b = ib.getXYZ(null);
                                double[] c = ib.getXYZ(null);
                                double[] d = ib.getXYZ(null);
                                double[] a = new double[3];
                                a[0] = (b[0] + c[0] + d[0]) / 3.0;
                                a[1] = (b[1] + c[1] + d[1]) / 3.0;
                                a[2] = (b[2] + c[2] + d[2]) / 3.0;
                                // Place the hydrogen at chiral position #1.
                                intxyz(hydrogen, ia, 1.0, ib, 109.5, ic, 109.5, 1);
                                double[] e1 = new double[3];
                                hydrogen.getXYZ(e1);
                                double[] ret = new double[3];
                                diff(a, e1, ret);
                                double l1 = r(ret);
                                // Place the hydrogen at chiral position #2.
                                intxyz(hydrogen, ia, 1.0, ib, 109.5, ic, 109.5, -1);
                                double[] e2 = new double[3];
                                hydrogen.getXYZ(e2);
                                diff(a, e2, ret);
                                double l2 = r(ret);
                                // Revert to #1 if it is farther from the average.
                                if (l1 > l2) {
                                    hydrogen.setXYZ(e1);
                                }
                                break;
                            case 2:
                                intxyz(hydrogen, ia, 1.0, ib, 109.5, ic, 109.5, 0);
                                break;
                            case 1:
                                intxyz(hydrogen, ia, 1.0, ib, 109.5, null, 0.0, 0);
                                break;
                            case 0:
                                intxyz(hydrogen, ia, 1.0, null, 0.0, null, 0.0, 0);
                                break;
                            default:
                                logger.info(" Check biotype for hydrogen " + atomName + ".");
                                patched = false;
                        }
                        break;
                    case 3:
                        switch(numBonds) {
                            case 2:
                                intxyz(hydrogen, ia, 1.0, ib, 120.0, ic, 0.0, 0);
                                break;
                            case 1:
                                intxyz(hydrogen, ia, 1.0, ib, 120.0, null, 0.0, 0);
                                break;
                            case 0:
                                intxyz(hydrogen, ia, 1.0, null, 0.0, null, 0.0, 0);
                                break;
                            default:
                                logger.info(" Check biotype for hydrogen " + atomName + ".");
                                patched = false;
                        }
                        break;
                    case 2:
                        switch(numBonds) {
                            case 1:
                                intxyz(hydrogen, ia, 1.0, ib, 120.0, null, 0.0, 0);
                                break;
                            case 0:
                                intxyz(hydrogen, ia, 1.0, null, 0.0, null, 0.0, 0);
                                break;
                            default:
                                logger.info(" Check biotype for hydrogen " + atomName + ".");
                                patched = false;
                        }
                        break;
                    case 1:
                        switch(numBonds) {
                            case 0:
                                intxyz(hydrogen, ia, 1.0, null, 0.0, null, 0.0, 0);
                                break;
                            default:
                                logger.info(" Check biotype for hydrogen " + atomName + ".");
                                patched = false;
                        }
                        break;
                    default:
                        logger.info(" Check biotype for hydrogen " + atomName + ".");
                        patched = false;
                }
                if (!patched) {
                    break;
                } else {
                    buildBond(ia, hydrogen);
                }
            }
        }
        if (!patched) {
            logger.log(Level.WARNING, format(" Deleting unrecognized molecule %s.", m.toString()));
            activeMolecularAssembly.deleteMolecule((Molecule) m);
        } else {
            logger.info(" Patch for " + moleculeName + " succeeded.");
        }
    }
    resolvePolymerLinks(molecules);
}

Example 82 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class PDBFilter method checkHydrogenAtomNames.

/**
 * Ensures proper naming of hydrogens according to latest PDB format.
 * Presently mostly guesses at which hydrogens to re-assign, which may cause
 * chirality errors for prochiral hydrogens. If necessary, we will implement
 * more specific mapping.
 *
 * @param residue
 */
private void checkHydrogenAtomNames(Residue residue) {
    switch(fileStandard) {
        case VERSION3_3:
            return;
        case VERSION3_2:
        default:
            break;
    }
    // May have to get position.
    String residueType = residue.getName().toUpperCase();
    ArrayList<Atom> resAtoms = residue.getAtomList();
    for (Atom atom : resAtoms) {
        if (atom == null) {
            continue;
        }
        String atomName = atom.getName().toUpperCase();
        // Handles situations such as 1H where it should be NA_H1, etc.
        if (atomName.contains("H")) {
            try {
                String firstChar = atomName.substring(0, 1);
                Integer.parseInt(firstChar);
                atomName = atomName.substring(1);
                atomName = atomName.concat(firstChar);
                atom.setName(atomName);
            } catch (NumberFormatException e) {
            // Do nothing.
            }
        }
    }
    // Ensures proper hydrogen assignment; for example, Gln should have HB2,
    // HB3 instead of HB1, HB2.
    ArrayList<Atom> betas;
    ArrayList<Atom> gammas;
    ArrayList<Atom> deltas;
    ArrayList<Atom> epsilons;
    ArrayList<Atom> zetas;
    String atomName;
    Atom OH;
    Atom HH;
    Atom HG;
    Atom HD2;
    switch(getAminoAcid(residueType)) {
        case GLY:
            ArrayList<Atom> alphas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HA")) {
                    alphas.add(atom);
                }
            }
            renameGlycineAlphaHydrogens(residue, alphas);
            break;
        case ALA:
            // No known errors with alanine
            break;
        case VAL:
            // No known errors with valine
            break;
        case LEU:
            betas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HB")) {
                    betas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case ILE:
            ArrayList<Atom> ileAtoms = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HG1")) {
                    ileAtoms.add(atom);
                }
            }
            renameIsoleucineHydrogens(residue, ileAtoms);
            break;
        case SER:
            betas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HB")) {
                    betas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case THR:
            Atom HG1 = (Atom) residue.getAtomNode("HG1");
            if (HG1 == null) {
                for (Atom atom : resAtoms) {
                    atomName = atom.getName().toUpperCase();
                    // Length < 4 avoids bringing in HG21, HG22, or HG23.
                    if (atomName.length() < 4 && atomName.contains("HG")) {
                        atom.setName("HG1");
                        break;
                    }
                }
            }
            break;
        case CYS:
            betas = new ArrayList<>();
            HG = (Atom) residue.getAtomNode("HG");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (HG == null && atomName.contains("HG")) {
                    HG = atom;
                    HG.setName("HG");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case CYX:
            // I pray this is never important, because I don't have an example CYX to work from.
            break;
        case CYD:
            betas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HB")) {
                    betas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case PRO:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            deltas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            renameDeltaHydrogens(residue, deltas, 23);
            break;
        case PHE:
            betas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            Atom HZ = (Atom) residue.getAtomNode("HZ");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (HZ == null && atomName.contains("HZ")) {
                    HZ = atom;
                    HZ.setName("HZ");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameDeltaHydrogens(residue, deltas, 12);
            renameEpsilonHydrogens(residue, epsilons, 12);
            break;
        case TYR:
            betas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            HH = (Atom) residue.getAtomNode("HH");
            OH = (Atom) residue.getAtomNode("OH");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (HH == null && atomName.contains("HH")) {
                    HH = atom;
                    HH.setName("HH");
                } else if (OH == null && atomName.contains("O") && atomName.contains("H")) {
                    OH = atom;
                    OH.setName("OH");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameDeltaHydrogens(residue, deltas, 12);
            renameEpsilonHydrogens(residue, epsilons, 12);
            break;
        case TYD:
            betas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            OH = (Atom) residue.getAtomNode("OH");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (OH == null && atomName.contains("O") && atomName.contains("H")) {
                    OH = atom;
                    OH.setName("OH");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameDeltaHydrogens(residue, deltas, 12);
            renameEpsilonHydrogens(residue, epsilons, 12);
            break;
        case TRP:
            betas = new ArrayList<>();
            epsilons = new ArrayList<>();
            zetas = new ArrayList<>();
            Atom HD1 = (Atom) residue.getAtomNode("HD1");
            Atom HH2 = (Atom) residue.getAtomNode("HH2");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (atomName.contains("HZ")) {
                    zetas.add(atom);
                } else if (HD1 == null && atomName.contains("HD")) {
                    HD1 = atom;
                    HD1.setName("HD1");
                } else if (HH2 == null && atomName.contains("HH")) {
                    HH2 = atom;
                    HH2.setName("HH2");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameEpsilonHydrogens(residue, epsilons, 13);
            renameZetaHydrogens(residue, zetas, 23);
            break;
        case HIS:
            betas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameDeltaHydrogens(residue, deltas, 12);
            renameEpsilonHydrogens(residue, epsilons, 12);
            break;
        case HID:
            betas = new ArrayList<>();
            deltas = new ArrayList<>();
            Atom HE1 = (Atom) residue.getAtomNode("HE1");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (HE1 == null && atomName.contains("HE")) {
                    HE1 = atom;
                    HE1.setName("HE1");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameDeltaHydrogens(residue, deltas, 12);
            break;
        case HIE:
            betas = new ArrayList<>();
            epsilons = new ArrayList<>();
            HD2 = (Atom) residue.getAtomNode("HD2");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (HD2 == null && atomName.contains("HD")) {
                    HD2 = atom;
                    HD2.setName("HD2");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameEpsilonHydrogens(residue, epsilons, 12);
            break;
        case ASP:
            betas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                if (atom.getName().toUpperCase().contains("HB")) {
                    betas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case ASH:
            betas = new ArrayList<>();
            HD2 = (Atom) residue.getAtomNode("HD2");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (HD2 == null && atomName.contains("HD")) {
                    HD2 = atom;
                    HD2.setName("HD2");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            break;
        case ASN:
            betas = new ArrayList<>();
            ArrayList<Atom> HD2s = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HD")) {
                    HD2s.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameAsparagineHydrogens(residue, HD2s);
            break;
        case GLU:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            break;
        case GLH:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            Atom HE2 = (Atom) residue.getAtomNode("HE2");
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (HE2 == null && atomName.contains("HE")) {
                    HE2 = atom;
                    HE2.setName("HE2");
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            break;
        case GLN:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            epsilons = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            renameGlutamineHydrogens(residue, epsilons);
            break;
        case MET:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            // Epsilons should not break, as they are 1-3.
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            break;
        case LYS:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            // Zetas are 1-3, should not break.
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            renameDeltaHydrogens(residue, deltas, 23);
            renameEpsilonHydrogens(residue, epsilons, 23);
            break;
        case LYD:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            deltas = new ArrayList<>();
            epsilons = new ArrayList<>();
            zetas = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (atomName.contains("HE")) {
                    epsilons.add(atom);
                } else if (atomName.contains("HZ")) {
                    zetas.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            renameDeltaHydrogens(residue, deltas, 23);
            renameEpsilonHydrogens(residue, epsilons, 23);
            renameZetaHydrogens(residue, zetas, 12);
            break;
        case ARG:
            betas = new ArrayList<>();
            gammas = new ArrayList<>();
            deltas = new ArrayList<>();
            Atom HE = (Atom) residue.getAtomNode("HE");
            ArrayList<Atom> HHn = new ArrayList<>();
            for (Atom atom : resAtoms) {
                atomName = atom.getName().toUpperCase();
                if (atomName.contains("HB")) {
                    betas.add(atom);
                } else if (atomName.contains("HG")) {
                    gammas.add(atom);
                } else if (atomName.contains("HD")) {
                    deltas.add(atom);
                } else if (HE == null && atomName.contains("HE")) {
                    HE = atom;
                    HE.setName("HE");
                } else if (atomName.contains("HH")) {
                    HHn.add(atom);
                }
            }
            renameBetaHydrogens(residue, betas, 23);
            renameGammaHydrogens(residue, gammas, 23);
            renameDeltaHydrogens(residue, deltas, 23);
            renameArginineHydrogens(residue, HHn);
            break;
        case ORN:
        case AIB:
        case PCA:
        case UNK:
        default:
            // labeled under older PDB standards.
            break;
    }
}

Example 83 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class SystemFilter method applyAtomProperties.

/**
 * Automatically sets atom-specific flags, particularly nouse and inactive,
 * and apply harmonic restraints. Intended to be called at the end of
 * readFile() implementations.
 */
public void applyAtomProperties() {
    /**
     * What may be a more elegant implementation is to make readFile() a
     * public concrete, but skeletal method, and then have readFile() call a
     * protected abstract readFile method for each implementation.
     */
    Atom[] molaAtoms = activeMolecularAssembly.getAtomArray();
    int nmolaAtoms = molaAtoms.length;
    String[] nouseKeys = properties.getStringArray("nouse");
    for (String nouseKey : nouseKeys) {
        /*try {
                int[] nouseRange = parseAtNumArg("nouse", nouseKey, nmolaAtoms);
                logger.log(Level.INFO, String.format(" Atoms %d-%d set to be not "
                        + "used", nouseRange[0] + 1, nouseRange[1] + 1));
                for (int i = nouseRange[0]; i <= nouseRange[1]; i++) {
                    molaAtoms[i].setUse(false);
                }
            } catch (IllegalArgumentException ex) {
                logger.log(Level.INFO, ex.getLocalizedMessage());
            }*/
        String[] toks = nouseKey.split("\\s+");
        for (String tok : toks) {
            try {
                int[] nouseRange = parseAtNumArg("restraint", tok, nmolaAtoms);
                logger.info(String.format(" Setting atoms %d-%d to not be used", nouseRange[0] + 1, nouseRange[1] + 1));
                for (int j = nouseRange[0]; j <= nouseRange[1]; j++) {
                    molaAtoms[j].setUse(false);
                }
            } catch (IllegalArgumentException ex) {
                boolean atomFound = false;
                for (Atom atom : molaAtoms) {
                    if (atom.getName().equalsIgnoreCase(tok)) {
                        atomFound = true;
                        atom.setUse(false);
                    }
                }
                if (atomFound) {
                    logger.info(String.format(" Setting atoms with name %s to not be used", tok));
                } else {
                    logger.log(Level.INFO, ex.getLocalizedMessage());
                }
            }
        }
    }
    String[] inactiveKeys = properties.getStringArray("inactive");
    for (String inactiveKey : inactiveKeys) {
        try {
            int[] inactiveRange = parseAtNumArg("inactive", inactiveKey, nmolaAtoms);
            logger.log(Level.INFO, String.format(" Atoms %d-%d set to be not " + "used", inactiveRange[0] + 1, inactiveRange[1] + 1));
            for (int i = inactiveRange[0]; i <= inactiveRange[1]; i++) {
                molaAtoms[i].setActive(false);
            }
        } catch (IllegalArgumentException ex) {
            logger.log(Level.INFO, ex.getLocalizedMessage());
        }
    /*Matcher m = intrangePattern.matcher(inactiveKey);
            if (m.matches()) {
                int start = Integer.parseInt(m.group(1)) - 1;
                int end = Integer.parseInt(m.group(2)) - 1;
                if (start > end) {
                    logger.log(Level.INFO, String.format(" Input %s not valid; "
                            + "start > end", inactiveKey));
                } else if (start < 0) {
                    logger.log(Level.INFO, String.format(" Input %s not valid; "
                            + "atoms should be indexed starting from 1", inactiveKey));
                } else {
                    logger.log(Level.INFO, String.format(" Atoms %s set to be "
                            + "inactive", inactiveKey));
                    for (int i = start; i <= end; i++) {
                        if (i >= nmolaAtoms) {
                            logger.log(Level.INFO, String.format(" Atom index %d is "
                                    + "out of bounds for molecular assembly of "
                                    + "length %d", i + 1, nmolaAtoms));
                            break;
                        }
                        molaAtoms[i].setActive(false);
                    }
                }
            } else {
                try {
                    int atNum = Integer.parseUnsignedInt(inactiveKey) - 1;
                    if (atNum >= nmolaAtoms) {
                        logger.log(Level.INFO, String.format(" Atom index %d is "
                                + "out of bounds for molecular assembly of "
                                + "length %d", atNum + 1, nmolaAtoms));
                    } else if (atNum < 0) {
                        logger.log(Level.INFO, String.format(" Input %s not valid; "
                                + "atoms should be indexed starting from 1", inactiveKey));
                    } else {
                        logger.log(Level.INFO, String.format(" Atom %s set to be "
                                + "inactive", inactiveKey));
                        molaAtoms[atNum].setActive(false);
                    }
                } catch (NumberFormatException ex) {
                    logger.log(Level.INFO, String.format(" inactive key %s cannot "
                            + "be interpreted as an atom number or range of atom "
                            + "numbers.", inactiveKey));
                }
            }*/
    }
    coordRestraints = new ArrayList<>();
    String[] cRestraintStrings = properties.getStringArray("restraint");
    for (String coordRestraint : cRestraintStrings) {
        String[] toks = coordRestraint.split("\\s+");
        double forceconst;
        try {
            forceconst = Double.parseDouble(toks[0]);
        } catch (NumberFormatException ex) {
            logger.log(Level.INFO, " First argument to coordinate restraint must be a positive force constant; discarding coordinate restraint.");
            continue;
        }
        if (forceconst < 0) {
            logger.log(Level.INFO, " Force constants must be positive. Discarding coordinate restraint.");
            continue;
        }
        logger.info(String.format(" Adding lambda-disabled coordinate restraint " + "with force constant %10.4f kcal/mol/A", forceconst));
        Set<Atom> restraintAtoms = new HashSet<>();
        for (int i = 1; i < toks.length; i++) {
            try {
                int[] restrRange = parseAtNumArg("restraint", toks[i], nmolaAtoms);
                logger.info(String.format(" Adding atoms %d-%d to restraint", restrRange[0] + 1, restrRange[1] + 1));
                for (int j = restrRange[0]; j <= restrRange[1]; j++) {
                    restraintAtoms.add(molaAtoms[j]);
                }
            } catch (IllegalArgumentException ex) {
                boolean atomFound = false;
                for (Atom atom : molaAtoms) {
                    if (atom.getName().equalsIgnoreCase(toks[i])) {
                        atomFound = true;
                        restraintAtoms.add(atom);
                    }
                }
                if (atomFound) {
                    logger.info(String.format(" Added atoms with name %s to restraint", toks[i]));
                } else {
                    logger.log(Level.INFO, ex.getLocalizedMessage());
                }
            }
        }
        if (!restraintAtoms.isEmpty()) {
            Atom[] ats = restraintAtoms.toArray(new Atom[restraintAtoms.size()]);
            coordRestraints.add(new CoordRestraint(ats, forceField, false, forceconst));
        } else {
            logger.warning(String.format(" Empty or unparseable restraint argument %s", coordRestraint));
        }
    }
    String[] lamRestraintStrings = properties.getStringArray("lamrestraint");
    for (String coordRestraint : lamRestraintStrings) {
        String[] toks = coordRestraint.split("\\s+");
        double forceconst = Double.parseDouble(toks[0]);
        logger.info(String.format(" Adding lambda-enabled coordinate restraint " + "with force constant %10.4f kcal/mol/A", forceconst));
        Set<Atom> restraintAtoms = new HashSet<>();
        for (int i = 1; i < toks.length; i++) {
            try {
                int[] restrRange = parseAtNumArg("restraint", toks[i], nmolaAtoms);
                logger.info(String.format(" Adding atoms %d-%d to restraint", restrRange[0] + 1, restrRange[1] + 1));
                for (int j = restrRange[0]; j <= restrRange[1]; j++) {
                    restraintAtoms.add(molaAtoms[j]);
                }
            } catch (IllegalArgumentException ex) {
                boolean atomFound = false;
                for (Atom atom : molaAtoms) {
                    if (atom.getName().equalsIgnoreCase(toks[i])) {
                        atomFound = true;
                        restraintAtoms.add(atom);
                    }
                }
                if (atomFound) {
                    logger.info(String.format(" Added atoms with name %s to restraint", toks[i]));
                } else {
                    logger.log(Level.INFO, String.format(" Restraint input %s " + "could not be parsed as a numerical range or " + "an atom type present in assembly", toks[i]));
                }
            }
        }
        if (!restraintAtoms.isEmpty()) {
            Atom[] ats = restraintAtoms.toArray(new Atom[restraintAtoms.size()]);
            coordRestraints.add(new CoordRestraint(ats, forceField, true, forceconst));
        } else {
            logger.warning(String.format(" Empty or unparseable restraint argument %s", coordRestraint));
        }
    }
    String[] xyzRestStrings = properties.getStringArray("xyzRestraint");
    // Variables to let sequential and otherwise identical xyzRestraint strings to be globbed into one restraint.
    List<Atom> restraintAts = new ArrayList<>();
    List<double[]> coords = new ArrayList<>();
    double lastForceConst = 0;
    boolean lastUseLam = false;
    for (String xR : xyzRestStrings) {
        String[] toks = xR.split("\\s+");
        int nToks = toks.length;
        if (nToks != 6) {
            logger.info(" XYZ restraint rejected: must have force constant, lambda boolean (true/false), 3 coordinates, and an atom number");
            logger.info(" For a coordinate restraint centered on original coordinates, use restraint or lamrestraint keys.");
            logger.info(String.format(" Rejected restraint string: %s", xR));
        } else {
            try {
                double forceConst = Double.parseDouble(toks[0]);
                boolean useLambda = Boolean.parseBoolean(toks[1]);
                if (forceConst != lastForceConst || useLambda != lastUseLam) {
                    int nAts = restraintAts.size();
                    if (nAts != 0) {
                        double[][] restXYZ = new double[3][nAts];
                        Atom[] atArr = restraintAts.toArray(new Atom[nAts]);
                        for (int i = 0; i < 3; i++) {
                            for (int j = 0; j < nAts; j++) {
                                restXYZ[i][j] = coords.get(j)[i];
                            }
                        }
                        CoordRestraint thePin = new CoordRestraint(atArr, forceField, lastUseLam, lastForceConst);
                        thePin.setCoordinatePin(restXYZ);
                        thePin.setIgnoreHydrogen(false);
                        coordRestraints.add(thePin);
                    }
                    restraintAts = new ArrayList<>();
                    coords = new ArrayList<>();
                    lastForceConst = forceConst;
                    lastUseLam = useLambda;
                }
                double[] atXYZ = new double[3];
                for (int i = 0; i < 3; i++) {
                    atXYZ[i] = Double.parseDouble(toks[i + 2]);
                }
                int atNum = Integer.parseInt(toks[5]) - 1;
                restraintAts.add(molaAtoms[atNum]);
                coords.add(atXYZ);
            } catch (Exception ex) {
                logger.info(String.format(" Exception parsing xyzRestraint %s: %s", xR, ex.toString()));
            }
        }
    }
    int nAts = restraintAts.size();
    if (nAts != 0) {
        double[][] restXYZ = new double[3][nAts];
        Atom[] atArr = restraintAts.toArray(new Atom[nAts]);
        for (int i = 0; i < 3; i++) {
            for (int j = 0; j < nAts; j++) {
                restXYZ[i][j] = coords.get(j)[i];
            }
        }
        CoordRestraint thePin = new CoordRestraint(atArr, forceField, lastUseLam, lastForceConst);
        thePin.setCoordinatePin(restXYZ);
        thePin.setIgnoreHydrogen(false);
        coordRestraints.add(thePin);
    }
    String[] noElStrings = properties.getStringArray("noElectro");
    for (String noE : noElStrings) {
        String[] toks = noE.split("\\s+");
        for (String tok : toks) {
            try {
                int[] noERange = parseAtNumArg("noElectro", tok, nmolaAtoms);
                for (int i = noERange[0]; i <= noERange[1]; i++) {
                    molaAtoms[i].setElectrostatics(false);
                }
                logger.log(Level.INFO, String.format(" Disabled electrostatics " + "for atoms %d-%d", noERange[0] + 1, noERange[1] + 1));
            } catch (IllegalArgumentException ex) {
                boolean atomFound = false;
                for (Atom atom : molaAtoms) {
                    if (atom.getName().equalsIgnoreCase(tok)) {
                        atomFound = true;
                        atom.setElectrostatics(false);
                    }
                }
                if (atomFound) {
                    logger.info(String.format(" Disabled electrostatics for atoms with name %s", tok));
                } else {
                    logger.log(Level.INFO, String.format(" No electrostatics " + "input %s could not be parsed as a numerical " + "range or atom type present in assembly", tok));
                }
            }
        }
    }
}

Example 84 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class XYZFilter method readOnto.

/**
 * <p>
 * readOnto</p>
 *
 * @param newFile a {@link java.io.File} object.
 * @param oldSystem a {@link ffx.potential.MolecularAssembly} object.
 * @return a boolean.
 */
public static boolean readOnto(File newFile, MolecularAssembly oldSystem) {
    if (newFile == null || !newFile.exists() || oldSystem == null) {
        return false;
    }
    try {
        FileReader fr = new FileReader(newFile);
        BufferedReader br = new BufferedReader(fr);
        String data = br.readLine();
        if (data == null) {
            return false;
        }
        String[] tokens = data.trim().split(" +");
        int num_atoms = Integer.parseInt(tokens[0]);
        if (num_atoms != oldSystem.getAtomList().size()) {
            return false;
        }
        br.mark(10000);
        data = br.readLine();
        if (!readPBC(data, oldSystem)) {
            br.reset();
        }
        double[][] d = new double[num_atoms][3];
        for (int i = 0; i < num_atoms; i++) {
            if (!br.ready()) {
                return false;
            }
            data = br.readLine();
            if (data == null) {
                logger.warning(format(" Check atom %d.", (i + 1)));
                return false;
            }
            tokens = data.trim().split(" +");
            if (tokens == null || tokens.length < 6) {
                logger.warning(format(" Check atom %d.", (i + 1)));
                return false;
            }
            d[i][0] = Double.parseDouble(tokens[2]);
            d[i][1] = Double.parseDouble(tokens[3]);
            d[i][2] = Double.parseDouble(tokens[4]);
        }
        ArrayList<Atom> atoms = oldSystem.getAtomList();
        for (Atom a : atoms) {
            int index = a.getIndex() - 1;
            a.setXYZ(d[index]);
        }
        oldSystem.center();
        oldSystem.setFile(newFile);
        br.close();
        fr.close();
        return true;
    } catch (Exception e) {
        return false;
    }
}

Example 85 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class XYZFilter method readNext.

/**
 * Reads the next snap-shot of an archive into the activeMolecularAssembly.
 * After calling this function, a BufferedReader will remain open until the
 * <code>close</code> method is called.
 *
 * @return true if successful.
 */
@Override
public boolean readNext(boolean resetPosition) {
    try {
        String data;
        Atom[] atoms = activeMolecularAssembly.getAtomArray();
        int nSystem = atoms.length;
        if (bufferedReader == null || resetPosition) /* || !bin.ready()*/
        {
            bufferedReader = new BufferedReader(new FileReader(currentFile));
            // Read past the first N + 1 lines that begin with an integer.
            for (int i = 0; i < nSystem + 1; i++) {
                data = bufferedReader.readLine();
                while (!firstTokenIsInteger(data)) {
                    data = bufferedReader.readLine();
                }
            }
            snapShot = 1;
        }
        snapShot++;
        data = bufferedReader.readLine();
        // Read past blank lines
        while (data != null && data.trim().equals("")) {
            data = bufferedReader.readLine();
        }
        if (data == null) {
            return false;
        }
        logger.info(String.format(" Attempting to read snapshot %d.", snapShot));
        try {
            int nArchive = Integer.parseInt(data.trim().split(" +")[0]);
            if (nArchive != nSystem) {
                String message = String.format("Number of atoms mismatch (Archive: %d, System: %d).", nArchive, nSystem);
                if (dieOnMissingAtom) {
                    logger.severe(message);
                }
                logger.warning(message);
                return false;
            }
        } catch (NumberFormatException e) {
            logger.warning(e.toString());
            return false;
        }
        // The header line is reasonable. Check for periodic box dimensions.
        bufferedReader.mark(10000);
        data = bufferedReader.readLine();
        if (!readPBC(data, activeMolecularAssembly)) {
            bufferedReader.reset();
        }
        for (int i = 0; i < nSystem; i++) {
            data = bufferedReader.readLine();
            // Read past blank lines
            while (data != null && data.trim().equals("")) {
                data = bufferedReader.readLine();
            }
            String[] tokens = data.trim().split(" +");
            if (tokens == null || tokens.length < 6) {
                String message = String.format("Check atom %d in %s.", (i + 1), currentFile.getName());
                logger.warning(message);
                return false;
            }
            double x = Double.parseDouble(tokens[2]);
            double y = Double.parseDouble(tokens[3]);
            double z = Double.parseDouble(tokens[4]);
            int xyzIndex = atoms[i].getIndex();
            if (xyzIndex != i + 1) {
                String message = String.format("Archive atom index %d being read onto system atom index %d.", i + 1, xyzIndex);
                logger.warning(message);
            }
            atoms[i].moveTo(x, y, z);
        }
        return true;
    } catch (FileNotFoundException e) {
        String message = String.format("Exception opening file %s.", currentFile);
        logger.log(Level.WARNING, message, e);
    } catch (IOException e) {
        String message = String.format("Exception reading from file %s.", currentFile);
        logger.log(Level.WARNING, message, e);
    }
    return false;
}