use of ffx.potential.bonded.ResidueEnumerations.AminoAcid3 in project ffx by mjschnie.
the class PhMD method recursiveMap.
/**
* Finds Titration definitions for the given Residue and adds them to the given MultiResidue.
* For three-state transitions, simply populate the enumeration with multiple titrations
* from a shared state and this will include them in MultiResidue construction.
*/
private void recursiveMap(Residue member, MultiResidue multiRes) {
// Map titrations for this member.
Titration[] titrations = Titration.multiLookup(member);
titrationMap.put(member, Arrays.asList(titrations));
// For each titration, check whether it needs added as a MultiResidue option.
for (Titration titration : titrations) {
// Allow manual override of Histidine treatment.
if ((titration.deprotForm == AminoAcid3.HID && config.histidineMode == HistidineMode.HIE_ONLY) || (titration.deprotForm == AminoAcid3.HIE && config.histidineMode == HistidineMode.HID_ONLY)) {
continue;
}
// Find all the choices currently available to this MultiResidue.
List<AminoAcid3> choices = new ArrayList<>();
for (Residue choice : multiRes.getConsideredResidues()) {
choices.add(choice.getAminoAcid3());
}
// If this Titration target is not a choice for the MultiResidue, then add it.
if (!choices.contains(titration.protForm) || !(choices.contains(titration.deprotForm))) {
String targetName = (member.getAminoAcid3() == titration.protForm) ? titration.deprotForm.toString() : titration.protForm.toString();
int resNumber = member.getResidueNumber();
ResidueType resType = member.getResidueType();
Residue newChoice = new Residue(targetName, resNumber, resType);
multiRes.addResidue(newChoice);
titrationMap.put(newChoice, Arrays.asList(Titration.multiLookup(newChoice)));
}
}
}
use of ffx.potential.bonded.ResidueEnumerations.AminoAcid3 in project ffx by mjschnie.
the class RosenbluthChiAllMove method engage_diffs.
private void engage_diffs() {
report.append(String.format(" Rosenbluth CBMC Move: %4d\n", moveNumber));
report.append(String.format(" residue: %s\n", target.toString()));
double origFastEnergy = fastEnergy();
double origSlowEnergy = slowEnergy();
AminoAcid3 name = AminoAcid3.valueOf(target.getName());
double[] chi = RotamerLibrary.measureRotamer(target, false);
// Now we're really just wingin' it: uDep is now "fast energy" (like RESPA does),
// uExt is the slow part, and we operate solely on diffs... so the bko_zero term -> unity.
TrialSet newTrialSet = cheapTorsionSet_diffs(testSetSize, "bkn");
// <-- explore difference between prod W(n) and sum wi
Wn = newTrialSet.sumExtBolt();
if (Wn <= 0) {
StringBuilder sb = new StringBuilder();
sb.append("Numerical instability in CMBC:");
sb.append(" Test set uExt values: ");
for (int i = 0; i < newTrialSet.uExt.length; i++) {
sb.append(String.format("%5.2f, ", newTrialSet.uExt[i]));
}
logger.warning(sb.toString());
}
// Choose a proposal move from amongst this trial set (bn).
double rng = rand.nextDouble(Wn);
double running = 0.0;
for (int j = 0; j < newTrialSet.uExt.length; j++) {
double uExtBolt = FastMath.exp(-beta * newTrialSet.uExt[j]);
running += uExtBolt;
if (rng < running) {
proposedMove = newTrialSet.rotamer[j];
double prob = uExtBolt / Wn * 100;
report.append(String.format(" Chose %d %7.4f\t%7.4f\t %4.1f%%\n", j + 1, newTrialSet.uExt[j], uExtBolt, prob));
break;
}
}
report.append(String.format(" Wn Total: %g\n", Wn));
// Reprise the above procedure for the old configuration.
// The existing conformation forms first member of each test set (wo).
target.revertState(origState);
// original-conf uDep
double ouDep = slowEnergy() - origSlowEnergy;
// original-conf uExt
double ouExt = fastEnergy() - origFastEnergy;
if (Math.abs(ouDep) > tolerance || Math.abs(ouExt) > tolerance) {
report.append(" WAIT WHAT (v2). ouDep/ouExt remainder: " + ouDep + "," + ouExt + "\n");
}
double ouExtBolt = FastMath.exp(-beta * ouExt);
report.append(String.format(" %3s %d: %9.5g %9.5g %9.5g\n", "bko", 0, ouDep, ouExt, ouExtBolt));
writeSnapshot("bko", true);
TrialSet oldTrialSet = cheapTorsionSet_diffs(testSetSize - 1, "bko");
Wo = ouExtBolt + oldTrialSet.sumExtBolt();
report.append(String.format(" Wo Total: %11.4g\n", Wo));
report.append(String.format(" Wn/Wo: %11.4g", Wn / Wo));
target.revertState(origState);
updateAll();
if (verbose) {
logger.info(report.toString());
}
}
use of ffx.potential.bonded.ResidueEnumerations.AminoAcid3 in project ffx by mjschnie.
the class PDBFilter method renameAminoAcidToPDBStandard.
/**
* Renames the Atoms in an amino acid to PDB standard.
*
* @param residue Residue to fix atom names of.
*/
private static void renameAminoAcidToPDBStandard(Residue residue) {
final Atom N = findNitrogenAtom(residue);
AminoAcid3 aa3 = residue.getAminoAcid3();
if (N != null) {
N.setName("N");
Atom CA = getAlphaCarbon(residue, N);
CA.setName("CA");
List<Atom> has = findBondedAtoms(CA, 1);
switch(aa3) {
case NME:
// Do all renaming here then return out of the method.
findBondedAtoms(N, 1).get(0).setName("H");
CA.setName("CH3");
for (int i = 1; i <= 3; i++) {
has.get(i - 1).setName(String.format("H%d", i));
}
return;
case GLY:
has.get(0).setName("HA2");
has.get(1).setName("HA3");
break;
default:
has.get(0).setName("HA");
break;
}
CA.setName("CA");
Atom C = null;
Atom CB = null;
for (Atom carb : findBondedAtoms(CA, 6)) {
// Second check is because of serine/threonine OG bonded straight to CB.
if (hasAttachedAtom(carb, 8) && !hasAttachedAtom(carb, 1)) {
C = carb;
C.setName("C");
} else {
CB = carb;
CB.setName("CB");
}
}
if (C == null) {
throw new IllegalArgumentException(String.format(" Could not find carbonyl carbon for residue %s!", residue));
}
if (CB == null && aa3 != AminoAcid3.GLY) {
throw new IllegalArgumentException(String.format(" Could not find beta carbon for residue %s!", residue));
}
List<Atom> ctermOxygens = findBondedAtoms(C, 8);
switch(ctermOxygens.size()) {
case 1:
ctermOxygens.get(0).setName("O");
break;
case 2:
Atom O = null;
for (Atom oxy : ctermOxygens) {
if (oxy.getBonds().size() == 2) {
O = oxy;
O.setName("O");
findBondedAtoms(O, 1).get(0).setName("HO");
}
}
if (O == null) {
ctermOxygens.get(0).setName("O");
ctermOxygens.get(1).setName("OXT");
}
}
List<Atom> amideProtons = findBondedAtoms(N, 1);
switch(amideProtons.size()) {
case 1:
amideProtons.get(0).setName("H");
break;
default:
// Should catch both N-termini and proline.
for (int i = 1; i <= amideProtons.size(); i++) {
amideProtons.get(i - 1).setName(String.format("H%d", i));
}
break;
}
/**
* All common atoms are now named: N, H[1-3], CA, HA[2-3], CB, C, O[XT], [HO]
*/
renameCommonAminoAcids(residue, aa3, CA, CB);
} else {
switch(aa3) {
case ACE:
{
Atom O = findAtomsOfElement(residue, 8).get(0);
O.setName("O");
Atom C = findBondedAtoms(O, 6).get(0);
C.setName("C");
Atom CH3 = findBondedAtoms(C, 6).get(0);
CH3.setName("CH3");
List<Atom> hydrogens = findBondedAtoms(CH3, 1);
for (int i = 1; i <= 3; i++) {
hydrogens.get(i - 1).setName(String.format("H%d", i));
}
}
break;
default:
throw new IllegalArgumentException(String.format(" Could not find nitrogen atom for residue %s!", residue));
}
}
}
use of ffx.potential.bonded.ResidueEnumerations.AminoAcid3 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;
}
Aggregations