use of ffx.potential.parameters.AtomType in project ffx by mjschnie.
the class SimulationFilter method readFile.
/**
* {@inheritDoc}
*/
@Override
public boolean readFile() {
// information
for (int i = 0; i < system.numatoms; i++) {
AtomType atomType = atomTypes.get(system.types[i]);
if (atomType == null) {
atomType = new AtomType(system.types[i], -1, system.name[i], system.story[i], system.atomic[i], system.mass[i], 0);
atomTypes.put(system.types[i], atomType);
}
}
atomList = new ArrayList<Atom>();
Vector<Integer> bonds1 = new Vector<Integer>();
Vector<Integer> bonds2 = new Vector<Integer>();
double[] d = new double[3];
int[] b = new int[4];
for (int i = 0; i < system.numatoms; i++) {
d[0] = system.coordinates[0][i];
d[1] = system.coordinates[1][i];
d[2] = system.coordinates[2][i];
String s = new String("" + system.types[i]);
AtomType atomType = atomTypes.get(s);
Atom a = new Atom(i + 1, new String("" + atomType.type), atomType, d);
atomList.add(a);
int b1 = i + 1;
b[0] = system.connectivity[0][i];
b[1] = system.connectivity[1][i];
b[2] = system.connectivity[2][i];
b[3] = system.connectivity[3][i];
int j = 0;
while (j < 4 && b[j] != 0) {
int b2 = b[j];
bonds1.add(b1);
bonds2.add(b2);
j++;
}
}
bondList = new ArrayList<Bond>();
for (int i = 0; i < bonds1.size(); i++) {
int a1 = bonds1.get(i);
int a2 = bonds2.get(i);
if (a1 < a2) {
Atom atom1 = atomList.get(a1 - 1);
Atom atom2 = atomList.get(a2 - 1);
bondList.add(new Bond(atom1, atom2));
}
}
setFileRead(true);
return true;
}
use of ffx.potential.parameters.AtomType in project ffx by mjschnie.
the class AminoAcidUtils method assignAminoAcidSideChain.
/**
* Assign atom types to a single amino acid side chain.
*
* @param position The position of this amino acid in the chain.
* @param aminoAcid The amino acid to use.
* @param residue The residue node.
* @param CA The C-alpha carbon of this residue.
* @param N The peptide nitrogen of this residue.
* @param C The peptide carbonyl carbon.
* @param forceField
* @param bondList
* @throws ffx.potential.bonded.BondedUtils.MissingHeavyAtomException this
* exception is thrown if when heavy is atom is missing that cannot be
* built.
*/
public static void assignAminoAcidSideChain(ResiduePosition position, AminoAcid3 aminoAcid, Residue residue, Atom CA, Atom N, Atom C, ForceField forceField, ArrayList<Bond> bondList) throws MissingHeavyAtomException {
int k = AA_CB[aminoAcid.ordinal()];
switch(aminoAcid) {
case GLY:
switch(position) {
case FIRST_RESIDUE:
k = AA_HA[0][k];
break;
case LAST_RESIDUE:
k = AA_HA[2][k];
break;
default:
k = AA_HA[1][k];
}
buildHydrogen(residue, "HA3", CA, 1.10, N, 109.5, C, 109.5, 1, k, forceField, bondList);
break;
case ALA:
buildAlanine(residue, CA, N, C, k, forceField, bondList);
break;
case VAL:
buildValine(residue, CA, N, C, k, forceField, bondList);
break;
case LEU:
buildLeucine(residue, CA, N, C, k, forceField, bondList);
break;
case ILE:
buildIsoleucine(residue, CA, N, C, k, forceField, bondList);
break;
case SER:
buildSerine(residue, CA, N, C, k, forceField, bondList);
break;
case THR:
buildThreonine(residue, CA, N, C, k, forceField, bondList);
break;
case CYS:
buildCysteine(residue, CA, N, C, k, forceField, bondList);
break;
case CYX:
buildCystine(residue, CA, N, C, k, forceField, bondList);
break;
case CYD:
buildDeprotonatedCysteine(residue, CA, N, C, k, forceField, bondList);
break;
case PRO:
buildProline(residue, CA, N, C, k, position, forceField, bondList);
break;
case PHE:
buildPhenylalanine(residue, CA, N, C, k, forceField, bondList);
break;
case TYR:
buildTyrosine(residue, CA, N, C, k, forceField, bondList);
break;
case TYD:
buildDeprotonatedTyrosine(residue, CA, N, C, k, forceField, bondList);
break;
case TRP:
buildTryptophan(residue, CA, N, C, k, forceField, bondList);
break;
case HIS:
buildHistidine(residue, CA, N, C, k, forceField, bondList);
break;
case HID:
buildNeutralHistidineD(residue, CA, N, C, k, forceField, bondList);
break;
case HIE:
buildNeutralHistidineE(residue, CA, N, C, k, forceField, bondList);
break;
case ASP:
buildAspartate(residue, CA, N, C, k, forceField, bondList);
break;
case ASH:
buildNeutralAsparticAcid(residue, CA, N, C, k, forceField, bondList);
break;
case ASN:
buildAsparagine(residue, CA, N, C, k, forceField, bondList);
break;
case GLU:
buildGlutamate(residue, CA, N, C, k, forceField, bondList);
break;
case GLH:
buildNeutralGlutamicAcid(residue, CA, N, C, k, forceField, bondList);
break;
case GLN:
buildGlutamine(residue, CA, N, C, k, forceField, bondList);
break;
case MET:
buildMethionine(residue, CA, N, C, k, forceField, bondList);
break;
case LYS:
buildLysine(residue, CA, N, C, k, forceField, bondList);
break;
case LYD:
buildDeprotonatedLysine(residue, CA, N, C, k, forceField, bondList);
break;
case ARG:
buildArginine(residue, CA, N, C, k, forceField, bondList);
break;
case ORN:
buildOrnithine(residue, CA, N, C, k, forceField, bondList);
break;
case AIB:
buildAIB(residue, CA, N, C, k, forceField, bondList);
break;
case PCA:
buildPCA(residue, CA, N, C, k, forceField, bondList);
break;
case UNK:
String residueName = residue.getName();
logger.log(Level.INFO, " Patching side-chain {0}", residueName);
HashMap<String, AtomType> types = forceField.getAtomTypes(residueName);
if (!types.isEmpty()) {
boolean patched = true;
ArrayList<Atom> residueAtoms = residue.getAtomList();
// Assign atom types for side-chain atoms.
for (Atom atom : residueAtoms) {
String atomName = atom.getName().toUpperCase();
AtomType type = atom.getAtomType();
if (type == null) {
type = types.get(atomName);
atom.setAtomType(type);
types.remove(atomName);
}
}
// Create bonds between known atoms.
if (patched) {
for (Atom atom : residueAtoms) {
String atomName = atom.getName();
String[] bonds = forceField.getBonds(residueName, atomName);
if (bonds != null) {
for (String name : bonds) {
Atom atom2 = (Atom) residue.getAtomNode(name);
if (atom2 != null && !atom.isBonded(atom2)) {
buildBond(atom, atom2, forceField, bondList);
}
}
}
}
}
// Create missing hydrogen atoms.
if (patched && !types.isEmpty()) {
// Create a hashmap of the molecule's atoms
HashMap<String, Atom> atomMap = new HashMap<>();
for (Atom atom : residueAtoms) {
atomMap.put(atom.getName().toUpperCase(), atom);
}
for (String atomName : types.keySet()) {
AtomType type = types.get(atomName);
String[] bonds = forceField.getBonds(residueName, atomName.toUpperCase());
if (bonds == null || bonds.length != 1) {
patched = false;
logger.log(Level.INFO, " Check biotype for hydrogen {0}.", 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.log(Level.FINE, " Created hydrogen {0}.", atomName);
hydrogen.setAtomType(type);
hydrogen.setHetero(true);
residue.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.log(Level.FINE, " Bonding {0} to {1} ({2} of {3}).", new Object[] { atomName, ia.getName(), numBonds, 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, 0);
double[] e1 = hydrogen.getXYZ(null);
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 = hydrogen.getXYZ(null);
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.log(Level.INFO, " Check biotype for hydrogen {0}.", atomName);
patched = false;
}
break;
case 3:
switch(numBonds) {
case 2:
intxyz(hydrogen, ia, 1.0, ib, 120.0, ic, 180.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.log(Level.INFO, " Check biotype for hydrogen {0}.", 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.log(Level.INFO, " Check biotype for hydrogen {0}.", 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.log(Level.INFO, " Check biotype for hydrogen {0}.", atomName);
patched = false;
}
break;
default:
logger.log(Level.INFO, " Check biotype for hydrogen {0}.", atomName);
patched = false;
}
if (!patched) {
break;
} else {
buildBond(ia, hydrogen, forceField, bondList);
}
}
}
if (!patched) {
logger.log(Level.SEVERE, format(" Could not patch %s.", residueName));
} else {
logger.log(Level.INFO, " Patch for {0} succeeded.", residueName);
residueAtoms = residue.getAtomList();
// Assign atom types for side-chain atoms.
double charge = 0.0;
for (Atom atom : residueAtoms) {
logger.info(atom.toString() + " -> " + atom.getAtomType().toString());
}
}
} else {
switch(position) {
case FIRST_RESIDUE:
buildHydrogen(residue, "HA2", CA, 1.10e0, N, 109.5e0, C, 109.5e0, 1, 355, forceField, bondList);
break;
case LAST_RESIDUE:
buildHydrogen(residue, "HA2", CA, 1.10e0, N, 109.5e0, C, 109.5e0, 1, 506, forceField, bondList);
break;
default:
buildHydrogen(residue, "HA2", CA, 1.10e0, N, 109.5e0, C, 109.5e0, 1, 6, forceField, bondList);
}
}
break;
}
}
use of ffx.potential.parameters.AtomType in project ffx by mjschnie.
the class AminoAcidUtils method assignAminoAcidAtomTypes.
public static void assignAminoAcidAtomTypes(Residue residue, Residue previousResidue, Residue nextResidue, ForceField forceField, ArrayList<Bond> bondList) throws MissingHeavyAtomException, MissingAtomTypeException {
String residueName = residue.getName().toUpperCase();
int j = 1;
ResiduePosition position = MIDDLE_RESIDUE;
if (previousResidue == null) {
j = 0;
position = FIRST_RESIDUE;
} else if (nextResidue == null) {
j = 2;
position = LAST_RESIDUE;
/**
* If the last residue only contains a nitrogen turn it into an NH2
* group.
*/
Atom N = (Atom) residue.getAtomNode("N");
if (residue.getAtomNodeList().size() == 1 && N != null) {
residueName = "NH2".intern();
residue.setName(residueName);
}
}
AminoAcid3 aminoAcid = getAminoAcid(residueName);
int aminoAcidNumber = getAminoAcidNumber(residueName);
/**
* Non-standard Amino Acid; use ALA backbone types.
*/
boolean nonStandard = false;
if (aminoAcid == AminoAcid3.UNK) {
aminoAcidNumber = getAminoAcidNumber("ALA");
nonStandard = true;
}
/**
* Only the last residue in a chain should have an OXT/OT2 atom.
*/
if (nextResidue != null) {
removeOXT_OT2(residue);
}
/**
* Only the first nitrogen should have H1, H2 and H3 atoms, unless it's
* an NME cap.
*/
if (previousResidue != null) {
removeH1_H2_H3(aminoAcid, residue);
}
/**
* Check for missing heavy atoms. This check ignores special terminating
* groups like FOR, NH2, etc.
*/
if (!nonStandard) {
try {
checkForMissingHeavyAtoms(aminoAcidNumber, aminoAcid, position, residue);
} catch (BondedUtils.MissingHeavyAtomException e) {
logger.log(Level.INFO, " {0} could not be parsed.", residue.toString());
logger.warning("MissingHeavyAtomException incoming from 194.");
throw e;
}
}
Atom pC = null;
Atom pCA = null;
if (previousResidue != null) {
pC = (Atom) previousResidue.getAtomNode("C");
pCA = (Atom) previousResidue.getAtomNode("CA");
}
/**
* Backbone heavy atoms.
*/
Atom N = (Atom) residue.getAtomNode("N");
if (N != null) {
N.setAtomType(BondedUtils.findAtomType(AA_N[j][aminoAcidNumber], forceField));
if (position != FIRST_RESIDUE) {
buildBond(pC, N, forceField, bondList);
}
}
Atom CA = null;
Atom C = null;
Atom O = null;
if (!(position == LAST_RESIDUE && aminoAcid == AminoAcid3.NH2)) {
if (aminoAcid == AminoAcid3.ACE || aminoAcid == AminoAcid3.NME) {
CA = buildHeavy(residue, "CH3", N, AA_CA[j][aminoAcidNumber], forceField, bondList);
} else {
CA = buildHeavy(residue, "CA", N, AA_CA[j][aminoAcidNumber], forceField, bondList);
}
if (!(position == LAST_RESIDUE && aminoAcid == AminoAcid3.NME)) {
C = buildHeavy(residue, "C", CA, AA_C[j][aminoAcidNumber], forceField, bondList);
O = (Atom) residue.getAtomNode("O");
if (O == null) {
O = (Atom) residue.getAtomNode("OT1");
}
AtomType atomType = findAtomType(AA_O[j][aminoAcidNumber], forceField);
if (O == null) {
MissingHeavyAtomException missingHeavyAtom = new MissingHeavyAtomException("O", atomType, C);
logger.warning(" MissingHeavyAtomException incoming from 234.");
throw missingHeavyAtom;
}
O.setAtomType(atomType);
buildBond(C, O, forceField, bondList);
}
}
/**
* Nitrogen hydrogen atoms.
*/
AtomType atomType = findAtomType(AA_HN[j][aminoAcidNumber], forceField);
switch(position) {
case FIRST_RESIDUE:
switch(aminoAcid) {
case PRO:
buildHydrogenAtom(residue, "H2", N, 1.02, CA, 109.5, C, 0.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H3", N, 1.02, CA, 109.5, C, -120.0, 0, atomType, forceField, bondList);
break;
case PCA:
buildHydrogenAtom(residue, "H", N, 1.02, CA, 109.5, C, -60.0, 0, atomType, forceField, bondList);
break;
case ACE:
break;
default:
buildHydrogenAtom(residue, "H1", N, 1.02, CA, 109.5, C, 180.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H2", N, 1.02, CA, 109.5, C, 60.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H3", N, 1.02, CA, 109.5, C, -60.0, 0, atomType, forceField, bondList);
}
break;
case LAST_RESIDUE:
switch(aminoAcid) {
case NH2:
buildHydrogenAtom(residue, "H1", N, 1.02, pC, 119.0, pCA, 0.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H2", N, 1.02, pC, 119.0, pCA, 180.0, 0, atomType, forceField, bondList);
break;
case NME:
buildHydrogenAtom(residue, "H", N, 1.02, pC, 118.0, CA, 121.0, 1, atomType, forceField, bondList);
break;
default:
buildHydrogenAtom(residue, "H", N, 1.02, pC, 119.0, CA, 119.0, 1, atomType, forceField, bondList);
}
break;
default:
// Mid-chain nitrogen hydrogen.
buildHydrogenAtom(residue, "H", N, 1.02, pC, 119.0, CA, 119.0, 1, atomType, forceField, bondList);
}
/**
* C-alpha hydrogen atoms.
*/
String haName = "HA";
if (aminoAcid == AminoAcid3.GLY) {
haName = "HA2";
}
atomType = findAtomType(AA_HA[j][aminoAcidNumber], forceField);
switch(position) {
case FIRST_RESIDUE:
switch(aminoAcid) {
case FOR:
buildHydrogenAtom(residue, "H", C, 1.12, O, 0.0, null, 0.0, 0, atomType, forceField, bondList);
break;
case ACE:
buildHydrogenAtom(residue, "H1", CA, 1.10, C, 109.5, O, 180.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H2", CA, 1.10, C, 109.5, O, 60.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H3", CA, 1.10, C, 109.5, O, -60.0, 0, atomType, forceField, bondList);
break;
default:
buildHydrogenAtom(residue, haName, CA, 1.10, N, 109.5, C, 109.5, -1, atomType, forceField, bondList);
break;
}
break;
case LAST_RESIDUE:
switch(aminoAcid) {
case NME:
buildHydrogenAtom(residue, "H1", CA, 1.10, N, 109.5, pC, 180.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H2", CA, 1.10, N, 109.5, pC, 60.0, 0, atomType, forceField, bondList);
buildHydrogenAtom(residue, "H3", CA, 1.10, N, 109.5, pC, -60.0, 0, atomType, forceField, bondList);
break;
default:
buildHydrogenAtom(residue, haName, CA, 1.10, N, 109.5, C, 109.5, -1, atomType, forceField, bondList);
}
break;
default:
buildHydrogenAtom(residue, haName, CA, 1.10, N, 109.5, C, 109.0, -1, atomType, forceField, bondList);
}
/**
* Build the amino acid side chain.
*/
assignAminoAcidSideChain(position, aminoAcid, residue, CA, N, C, forceField, bondList);
/**
* Build the terminal oxygen if the residue is not NH2 or NME.
*/
if (position == LAST_RESIDUE && !(aminoAcid == AminoAcid3.NH2 || aminoAcid == AminoAcid3.NME)) {
atomType = findAtomType(AA_O[2][aminoAcidNumber], forceField);
Atom OXT = (Atom) residue.getAtomNode("OXT");
if (OXT == null) {
OXT = (Atom) residue.getAtomNode("OT2");
if (OXT != null) {
OXT.setName("OXT");
}
}
if (OXT == null) {
String resName = C.getResidueName();
int resSeq = C.getResidueNumber();
Character chainID = C.getChainID();
Character altLoc = C.getAltLoc();
String segID = C.getSegID();
double occupancy = C.getOccupancy();
double tempFactor = C.getTempFactor();
OXT = new Atom(0, "OXT", altLoc, new double[3], resName, resSeq, chainID, occupancy, tempFactor, segID);
OXT.setAtomType(atomType);
residue.addMSNode(OXT);
intxyz(OXT, C, 1.25, CA, 117.0, O, 126.0, 1);
} else {
OXT.setAtomType(atomType);
}
buildBond(C, OXT, forceField, bondList);
}
/**
* Do some checks on the current residue to make sure all atoms have
* been assigned an atom type.
*/
List<Atom> resAtoms = residue.getAtomList();
for (Atom atom : resAtoms) {
atomType = atom.getAtomType();
if (atomType == null) {
/**
* Sometimes, with deuterons, a proton has been constructed in
* its place, so we have a "dummy" deuteron still hanging
* around.
*/
String protonEq = atom.getName().replaceFirst("D", "H");
Atom correspH = (Atom) residue.getAtomNode(protonEq);
if (correspH == null || correspH.getAtomType() == null) {
MissingAtomTypeException missingAtomTypeException = new MissingAtomTypeException(residue, atom);
logger.warning("MissingAtomTypeException incoming from 393.");
throw missingAtomTypeException;
} else {
correspH.setName(atom.getName());
atom.removeFromParent();
atom = correspH;
atomType = atom.getAtomType();
}
}
int numberOfBonds = atom.getNumBonds();
if (numberOfBonds != atomType.valence) {
if (atom == C && numberOfBonds == atomType.valence - 1 && position != LAST_RESIDUE) {
continue;
}
logger.warning(format(" An atom for residue %s has the wrong number of bonds:\n %s", residueName, atom.toString()));
logger.warning(format(" Expected: %d Actual: %d.", atomType.valence, numberOfBonds));
}
}
}
use of ffx.potential.parameters.AtomType in project ffx by mjschnie.
the class BondedUtils method buildHeavy.
public static Atom buildHeavy(MSGroup residue, String atomName, Atom bondedTo, int key, ForceField forceField, ArrayList<Bond> bondList) throws MissingHeavyAtomException {
Atom atom = (Atom) residue.getAtomNode(atomName);
AtomType atomType = findAtomType(key, forceField);
if (atom == null) {
MissingHeavyAtomException missingHeavyAtom = new MissingHeavyAtomException(atomName, atomType, bondedTo);
throw missingHeavyAtom;
}
atom.setAtomType(atomType);
if (bondedTo != null) {
buildBond(atom, bondedTo, forceField, bondList);
}
return atom;
}
use of ffx.potential.parameters.AtomType in project ffx by mjschnie.
the class RosenbluthChiAllMove method torsionSampler.
private void torsionSampler(List<Torsion> allTors) {
// Collects data for plot of uTors vs. theta.
// This is for finding the appropriate offset to add to each uTors such that the minimum is zero.
// double origChi[] = RotamerLibrary.measureRotamer(target, false);
updateAll();
double[] origChi = measureLysine(target, false);
StringBuilder sb = new StringBuilder();
sb.append(String.format(" Torsion Sampling for %s\n", target.getName()));
sb.append(String.format(" origChi:"));
for (int k = 0; k < origChi.length; k++) {
sb.append(String.format(" %6.2f", origChi[k]));
}
sb.append("\n");
for (int k = 0; k < origChi.length; k++) {
Torsion tors = allTors.get(k);
AtomType type1 = tors.getAtomArray()[0].getAtomType();
AtomType type2 = tors.getAtomArray()[1].getAtomType();
AtomType type3 = tors.getAtomArray()[2].getAtomType();
AtomType type4 = tors.getAtomArray()[3].getAtomType();
sb.append(String.format(" %d: \"(%3d %3d %3s) (%3d %3d %3s) (%3d %3d %3s) (%3d %3d %3s)\"\n", k, type1.type, type1.atomClass, type1.name, type2.type, type2.atomClass, type2.name, type3.type, type3.atomClass, type3.name, type4.type, type4.atomClass, type4.name));
}
logger.info(sb.toString());
sb = new StringBuilder();
String type = "combinatorial";
if (type.equals("combinatorial")) {
sb.append(String.format(" Resi chi0 chi1 chi2 chi3 List<uTors> uTorsSum uTotal\n"));
logger.info(sb.toString());
sb = new StringBuilder();
boolean doChi0 = false, doChi1 = false;
boolean doChi2 = true, doChi3 = true;
double increment = 1.0;
for (double chi0 = -180.0; chi0 < +180.0; chi0 += increment) {
for (double chi1 = -180.0; chi1 <= +180.0; chi1 += increment) {
for (double chi2 = -180.0; chi2 <= +180.0; chi2 += increment) {
for (double chi3 = -180.0; chi3 <= +180.0; chi3 += increment) {
sb.append(String.format(" %3s", target.getName()));
double[] newChi = new double[4];
if (doChi0) {
newChi[0] = chi0;
} else {
newChi[0] = origChi[0];
}
if (doChi1) {
newChi[1] = chi1;
} else {
newChi[1] = origChi[1];
}
if (doChi2) {
newChi[2] = chi2;
} else {
newChi[2] = origChi[2];
}
if (doChi3) {
newChi[3] = chi3;
} else {
newChi[3] = origChi[3];
}
Torsion chi0Tors = allTors.get(0);
Torsion chi1Tors = allTors.get(1);
Torsion chi2Tors = allTors.get(2);
Torsion chi3Tors = allTors.get(3);
Rotamer newState = createRotamer(target, newChi);
RotamerLibrary.applyRotamer(target, newState);
for (int wut = 0; wut < origChi.length; wut++) {
sb.append(String.format(" %3.0f", newChi[wut]));
}
writeSnapshot(String.format("%.0f", chi1), false);
double uTorsSum = 0.0;
for (int k = 0; k < allTors.size(); k++) {
double uTors = allTors.get(k).energy(false);
sb.append(String.format(" %5.2f", uTors));
uTorsSum += uTors;
}
sb.append(String.format(" %5.2f", uTorsSum));
double totalE = 1000.0;
try {
totalE = totalEnergy();
} catch (Exception ex) {
}
if (totalE > 1000.0) {
totalE = 1000.0;
}
sb.append(String.format(" %10.4f", totalE));
sb.append(String.format("\n"));
logger.info(sb.toString());
sb = new StringBuilder();
if (!doChi3) {
break;
}
}
if (!doChi2) {
break;
}
}
if (!doChi1) {
break;
}
}
if (!doChi0) {
break;
}
}
} else {
sb.append(String.format(" Resi Theta ( chi ) List<uTors,uTotal> \n"));
logger.info(sb.toString());
sb = new StringBuilder();
for (double testTheta = -180.0; testTheta <= +180.0; testTheta += 0.01) {
sb.append(String.format(" %3s %6.2f", target.getName(), testTheta));
for (int k = 0; k < origChi.length; k++) {
double[] newChi = new double[origChi.length];
System.arraycopy(origChi, 0, newChi, 0, origChi.length);
Torsion tors = allTors.get(k);
newChi[k] = testTheta;
Rotamer newState = createRotamer(target, newChi);
RotamerLibrary.applyRotamer(target, newState);
sb.append(" (");
for (int wut = 0; wut < origChi.length; wut++) {
sb.append(String.format(" %6.2f", newChi[wut]));
}
sb.append(" )");
writeSnapshot(String.format("%.0f", testTheta), false);
double uTors = allTors.get(k).energy(false);
double totalE = 0.0;
try {
totalE = totalEnergy();
} catch (Exception ex) {
}
sb.append(String.format(" %9.6f %9.6f", uTors, totalE));
}
sb.append(String.format("\n"));
logger.info(sb.toString());
sb = new StringBuilder();
}
}
if (System.getProperty("cbmc-torsionSampler") != null) {
try {
File output = new File(System.getProperty("cbmc-torsionSampler"));
BufferedWriter bw = new BufferedWriter(new FileWriter(output));
bw.write(sb.toString());
bw.close();
} catch (IOException ex) {
}
}
System.exit(0);
}
Aggregations