Example 26 with Bond
use of ffx.potential.bonded.Bond in project ffx by mjschnie.
the class PDBFilter method readFile.
/**
* {@inheritDoc}
*
* Parse the PDB File
*
* @return true if the file is read successfully.
*/
@Override
public boolean readFile() {
// First atom is #1, to match xyz file format
int xyzIndex = 1;
setFileRead(false);
systems.add(activeMolecularAssembly);
List<String> conects = new ArrayList<>();
List<String> links = new ArrayList<>();
List<String> ssbonds = new ArrayList<>();
List<String> structs = new ArrayList<>();
BufferedReader br = null;
try {
for (File file : files) {
currentFile = file;
if (mutate) {
List<Character> chainIDs = new ArrayList<>();
try (BufferedReader cr = new BufferedReader(new FileReader(file))) {
String line = cr.readLine();
while (line != null) {
// Replace all tabs w/ 4x spaces
line = line.replaceAll("\t", " ");
String identity = line;
if (line.length() > 6) {
identity = line.substring(0, 6);
}
identity = identity.trim().toUpperCase();
Record record;
try {
record = Record.valueOf(identity);
} catch (Exception e) {
/**
* Continue until the record is recognized.
*/
line = cr.readLine();
continue;
}
switch(record) {
case ANISOU:
case HETATM:
case ATOM:
char c22 = line.charAt(21);
boolean idFound = false;
for (Character chainID : chainIDs) {
if (c22 == chainID) {
idFound = true;
break;
}
}
if (!idFound) {
chainIDs.add(c22);
}
break;
}
line = cr.readLine();
}
for (Mutation mtn : mutations) {
if (!chainIDs.contains(mtn.chainChar)) {
if (chainIDs.size() == 1) {
logger.warning(String.format(" Chain ID %c for " + "mutation not found: only one chain %c " + "found.", mtn.chainChar, chainIDs.get(0)));
mtn = new Mutation(mtn.resID, chainIDs.get(0), mtn.resName);
} else {
logger.warning(String.format(" Chain ID %c for " + "mutation not found: mutation will not " + "proceed.", mtn.chainChar));
}
}
}
} catch (IOException ex) {
logger.finest(String.format(" Exception %s in parsing file to find chain IDs", ex.toString()));
}
}
/**
* Check that the current file exists and that we can read it.
*/
if (currentFile == null || !currentFile.exists() || !currentFile.canRead()) {
return false;
}
/**
* Open the current file for parsing.
*/
FileReader fr = new FileReader(currentFile);
br = new BufferedReader(fr);
/**
* Echo the alternate location being parsed.
*/
if (currentAltLoc == 'A') {
logger.info(format(" Reading %s", currentFile.getName()));
} else {
logger.info(format(" Reading %s alternate location %s", currentFile.getName(), currentAltLoc));
}
/**
* Reset the current chain and segID.
*/
currentChainID = null;
currentSegID = null;
boolean containsInsCode = false;
/**
* Read the first line of the file.
*/
String line = br.readLine();
/**
* Parse until END or ENDMDL is found, or to the end of the
* file.
*/
while (line != null) {
// Replace all tabs w/ 4x spaces
line = line.replaceAll("\t", " ");
String identity = line;
if (line.length() > 6) {
identity = line.substring(0, 6);
}
identity = identity.trim().toUpperCase();
Record record;
try {
record = Record.valueOf(identity);
} catch (Exception e) {
/**
* Continue until the record is recognized.
*/
line = br.readLine();
continue;
}
/**
* Switch on the known record.
*/
switch(record) {
case ENDMDL:
case END:
/**
* Setting "line" to null will exit the loop.
*/
line = null;
continue;
case DBREF:
// =============================================================================
// 1 - 6 Record name "DBREF "
// 8 - 11 IDcode idCode ID code of this entry.
// 13 Character chainID Chain identifier.
// 15 - 18 Integer seqBegin Initial sequence number of the
// PDB sequence segment.
// 19 AChar insertBegin Initial insertion code of the
// PDB sequence segment.
// 21 - 24 Integer seqEnd Ending sequence number of the
// PDB sequence segment.
// 25 AChar insertEnd Ending insertion code of the
// PDB sequence segment.
// 27 - 32 LString database Sequence database name.
// 34 - 41 LString dbAccession Sequence database accession code.
// 43 - 54 LString dbIdCode Sequence database identification code.
// 56 - 60 Integer dbseqBegin Initial sequence number of the
// database seqment.
// 61 AChar idbnsBeg Insertion code of initial residue of the
// segment, if PDB is the reference.
// 63 - 67 Integer dbseqEnd Ending sequence number of the
// database segment.
// 68 AChar dbinsEnd Insertion code of the ending residue of
// the segment, if PDB is the reference.
// =============================================================================
Character chainID = line.substring(12, 13).toUpperCase().charAt(0);
int seqBegin = Integer.parseInt(line.substring(14, 18).trim());
int seqEnd = Integer.parseInt(line.substring(20, 24).trim());
int[] seqRange = dbref.get(chainID);
if (seqRange == null) {
seqRange = new int[2];
dbref.put(chainID, seqRange);
}
seqRange[0] = seqBegin;
seqRange[1] = seqEnd;
break;
case SEQRES:
// =============================================================================
// 1 - 6 Record name "SEQRES"
// 8 - 10 Integer serNum Serial number of the SEQRES record for the
// current chain. Starts at 1 and increments
// by one each line. Reset to 1 for each chain.
// 12 Character chainID Chain identifier. This may be any single
// legal character, including a blank which is
// is used if there is only one chain.
// 14 - 17 Integer numRes Number of residues in the chain.
// This value is repeated on every record.
// 20 - 22 Residue name resName Residue name.
// 24 - 26 Residue name resName Residue name.
// 28 - 30 Residue name resName Residue name.
// 32 - 34 Residue name resName Residue name.
// 36 - 38 Residue name resName Residue name.
// 40 - 42 Residue name resName Residue name.
// 44 - 46 Residue name resName Residue name.
// 48 - 50 Residue name resName Residue name.
// 52 - 54 Residue name resName Residue name.
// 56 - 58 Residue name resName Residue name.
// 60 - 62 Residue name resName Residue name.
// 64 - 66 Residue name resName Residue name.
// 68 - 70 Residue name resName Residue name.
// =============================================================================
activeMolecularAssembly.addHeaderLine(line);
chainID = line.substring(11, 12).toUpperCase().charAt(0);
int serNum = Integer.parseInt(line.substring(7, 10).trim());
String[] chain = seqres.get(chainID);
int numRes = Integer.parseInt(line.substring(13, 17).trim());
if (chain == null) {
chain = new String[numRes];
seqres.put(chainID, chain);
}
int resID = (serNum - 1) * 13;
int end = line.length();
for (int start = 19; start + 3 <= end; start += 4) {
String res = line.substring(start, start + 3).trim();
if (res == null || res.length() < 1) {
break;
}
chain[resID++] = res;
}
break;
case MODRES:
String modResName = line.substring(12, 15).trim();
String stdName = line.substring(24, 27).trim();
modres.put(modResName.toUpperCase(), stdName.toUpperCase());
activeMolecularAssembly.addHeaderLine(line);
// =============================================================================
break;
case ANISOU:
// =============================================================================
// 1 - 6 Record name "ANISOU"
// 7 - 11 Integer serial Atom serial number.
// 13 - 16 Atom name Atom name.
// 17 Character altLoc Alternate location indicator
// 18 - 20 Residue name resName Residue name.
// 22 Character chainID Chain identifier.
// 23 - 26 Integer resSeq Residue sequence number.
// 27 AChar iCode Insertion code.
// 29 - 35 Integer u[0][0] U(1,1)
// 36 - 42 Integer u[1][1] U(2,2)
// 43 - 49 Integer u[2][2] U(3,3)
// 50 - 56 Integer u[0][1] U(1,2)
// 57 - 63 Integer u[0][2] U(1,3)
// 64 - 70 Integer u[1][2] U(2,3)
// 77 - 78 LString(2) element Element symbol, right-justified.
// 79 - 80 LString(2) charge Charge on the atom.
// =============================================================================
boolean deleteAnisou = properties.getBoolean("delete-anisou", false);
if (deleteAnisou) {
break;
}
Integer serial = Hybrid36.decode(5, line.substring(6, 11));
Character altLoc = line.substring(16, 17).toUpperCase().charAt(0);
if (!altLocs.contains(altLoc)) {
altLocs.add(altLoc);
}
if (!altLoc.equals(' ') && !altLoc.equals('A') && !altLoc.equals(currentAltLoc)) {
break;
}
double[] adp = new double[6];
adp[0] = new Integer(line.substring(28, 35).trim()) * 1.0e-4;
adp[1] = new Integer(line.substring(35, 42).trim()) * 1.0e-4;
adp[2] = new Integer(line.substring(42, 49).trim()) * 1.0e-4;
adp[3] = new Integer(line.substring(49, 56).trim()) * 1.0e-4;
adp[4] = new Integer(line.substring(56, 63).trim()) * 1.0e-4;
adp[5] = new Integer(line.substring(63, 70).trim()) * 1.0e-4;
if (atoms.containsKey(serial)) {
Atom a = atoms.get(serial);
a.setAltLoc(altLoc);
a.setAnisou(adp);
} else {
logger.info(format(" No ATOM record for ANISOU serial number %d has been found.", serial));
logger.info(format(" This ANISOU record will be ignored:\n %s", line));
}
break;
case ATOM:
// =============================================================================
// 1 - 6 Record name "ATOM "
// 7 - 11 Integer serial Atom serial number.
// 13 - 16 Atom name Atom name.
// 17 Character altLoc Alternate location indicator.
// 18 - 20 Residue name resName Residue name.
// 22 Character chainID Chain identifier.
// 23 - 26 Integer resSeq Residue sequence number.
// 27 AChar iCode Code for insertion of residues.
// 31 - 38 Real(8.3) x Orthogonal coordinates for X in Angstroms.
// 39 - 46 Real(8.3) y Orthogonal coordinates for Y in Angstroms.
// 47 - 54 Real(8.3) z Orthogonal coordinates for Z in Angstroms.
// 55 - 60 Real(6.2) occupancy Occupancy.
// 61 - 66 Real(6.2) tempFactor Temperature factor.
// 77 - 78 LString(2) element Element symbol, right-justified.
// 79 - 80 LString(2) charge Charge on the atom.
// =============================================================================
String name;
String resName;
String segID;
int resSeq;
boolean printAtom;
double[] d;
double occupancy;
double tempFactor;
Atom newAtom;
Atom returnedAtom;
// If it's a misnamed water, it will fall through to HETATM.
if (!line.substring(17, 20).trim().equals("HOH")) {
serial = Hybrid36.decode(5, line.substring(6, 11));
name = line.substring(12, 16).trim();
if (name.toUpperCase().contains("1H") || name.toUpperCase().contains("2H") || name.toUpperCase().contains("3H")) {
// VERSION3_2 is presently just a placeholder for "anything non-standard".
fileStandard = VERSION3_2;
}
altLoc = line.substring(16, 17).toUpperCase().charAt(0);
if (!altLocs.contains(altLoc)) {
altLocs.add(altLoc);
}
if (!altLoc.equals(' ') && !altLoc.equals('A') && !altLoc.equals(currentAltLoc)) {
break;
}
resName = line.substring(17, 20).trim();
chainID = line.substring(21, 22).charAt(0);
segID = getSegID(chainID);
resSeq = Hybrid36.decode(4, line.substring(22, 26));
char insertionCode = line.charAt(26);
if (insertionCode != ' ' && !containsInsCode) {
containsInsCode = true;
logger.warning(" FFX support for files with " + "insertion codes is experimental. " + "Residues will be renumbered to " + "eliminate insertion codes (52A " + "becomes 53, 53 becomes 54, etc)");
}
int offset = inscodeCount.getOrDefault(chainID, 0);
String pdbResNum = String.format("%c%d%c", chainID, resSeq, insertionCode);
if (!pdbToNewResMap.containsKey(pdbResNum)) {
if (insertionCode != ' ') {
++offset;
inscodeCount.put(chainID, offset);
}
resSeq += offset;
if (offset != 0) {
logger.info(String.format(" Chain %c " + "residue %s-%s renumbered to %c %s-%d", chainID, pdbResNum.substring(1).trim(), resName, chainID, resName, resSeq));
}
String newNum = String.format("%c%d", chainID, resSeq);
pdbToNewResMap.put(pdbResNum, newNum);
} else {
resSeq += offset;
}
printAtom = false;
if (mutate) {
boolean doBreak = false;
for (Mutation mtn : mutations) {
if (chainID == mtn.chainChar && resSeq == mtn.resID) {
String atomName = name.toUpperCase();
/*if (atomName.equals("N") || atomName.equals("C")
|| atomName.equals("O") || atomName.equals("CA")) {*/
if (backboneNames.contains(atomName)) {
printAtom = true;
resName = mtn.resName;
} else {
logger.info(String.format(" Deleting atom %s of %s %d", atomName, resName, resSeq));
doBreak = true;
break;
}
}
}
if (doBreak) {
break;
}
}
d = new double[3];
d[0] = new Double(line.substring(30, 38).trim());
d[1] = new Double(line.substring(38, 46).trim());
d[2] = new Double(line.substring(46, 54).trim());
occupancy = 1.0;
tempFactor = 1.0;
try {
occupancy = new Double(line.substring(54, 60).trim());
tempFactor = new Double(line.substring(60, 66).trim());
} catch (NumberFormatException | StringIndexOutOfBoundsException e) {
// Use default values.
if (printMissingFields) {
logger.info(" Missing occupancy and b-factors set to 1.0.");
printMissingFields = false;
} else if (logger.isLoggable(Level.FINE)) {
logger.fine(" Missing occupancy and b-factors set to 1.0.");
}
}
newAtom = new Atom(0, name, altLoc, d, resName, resSeq, chainID, occupancy, tempFactor, segID);
// Check if this is a modified residue.
if (modres.containsKey(resName.toUpperCase())) {
newAtom.setModRes(true);
}
returnedAtom = (Atom) activeMolecularAssembly.addMSNode(newAtom);
if (returnedAtom != newAtom) {
// A previously added atom has been retained.
atoms.put(serial, returnedAtom);
if (logger.isLoggable(Level.FINE)) {
logger.fine(returnedAtom + " has been retained over\n" + newAtom);
}
} else {
// The new atom has been added.
atoms.put(serial, newAtom);
// Check if the newAtom took the xyzIndex of a previous alternate conformer.
if (newAtom.getIndex() == 0) {
newAtom.setXyzIndex(xyzIndex++);
}
if (printAtom) {
logger.info(newAtom.toString());
}
}
break;
}
case HETATM:
// =============================================================================
// 1 - 6 Record name "HETATM"
// 7 - 11 Integer serial Atom serial number.
// 13 - 16 Atom name Atom name.
// 17 Character altLoc Alternate location indicator.
// 18 - 20 Residue name resName Residue name.
// 22 Character chainID Chain identifier.
// 23 - 26 Integer resSeq Residue sequence number.
// 27 AChar iCode Code for insertion of residues.
// 31 - 38 Real(8.3) x Orthogonal coordinates for X.
// 39 - 46 Real(8.3) y Orthogonal coordinates for Y.
// 47 - 54 Real(8.3) z Orthogonal coordinates for Z.
// 55 - 60 Real(6.2) occupancy Occupancy.
// 61 - 66 Real(6.2) tempFactor Temperature factor.
// 77 - 78 LString(2) element Element symbol; right-justified.
// 79 - 80 LString(2) charge Charge on the atom.
// =============================================================================
serial = Hybrid36.decode(5, line.substring(6, 11));
name = line.substring(12, 16).trim();
altLoc = line.substring(16, 17).toUpperCase().charAt(0);
if (!altLocs.contains(altLoc)) {
altLocs.add(altLoc);
}
if (!altLoc.equals(' ') && !altLoc.equals('A') && !altLoc.equals(currentAltLoc)) {
break;
}
resName = line.substring(17, 20).trim();
chainID = line.substring(21, 22).charAt(0);
segID = getSegID(chainID);
resSeq = Hybrid36.decode(4, line.substring(22, 26));
char insertionCode = line.charAt(26);
if (insertionCode != ' ' && !containsInsCode) {
containsInsCode = true;
logger.warning(" FFX support for files with " + "insertion codes is experimental. " + "Residues will be renumbered to " + "eliminate insertion codes (52A " + "becomes 53, 53 becomes 54, etc)");
}
int offset = inscodeCount.getOrDefault(chainID, 0);
String pdbResNum = String.format("%c%d%c", chainID, resSeq, insertionCode);
if (!pdbToNewResMap.containsKey(pdbResNum)) {
if (insertionCode != ' ') {
++offset;
inscodeCount.put(chainID, offset);
}
resSeq += offset;
if (offset != 0) {
logger.info(String.format(" Chain %c " + "molecule %s-%s renumbered to %c %s-%d", chainID, pdbResNum.substring(1).trim(), resName, chainID, resName, resSeq));
}
String newNum = String.format("%c%d", chainID, resSeq);
pdbToNewResMap.put(pdbResNum, newNum);
} else {
resSeq += offset;
}
d = new double[3];
d[0] = new Double(line.substring(30, 38).trim());
d[1] = new Double(line.substring(38, 46).trim());
d[2] = new Double(line.substring(46, 54).trim());
occupancy = 1.0;
tempFactor = 1.0;
try {
occupancy = new Double(line.substring(54, 60).trim());
tempFactor = new Double(line.substring(60, 66).trim());
} catch (NumberFormatException | StringIndexOutOfBoundsException e) {
// Use default values.
if (printMissingFields) {
logger.info(" Missing occupancy and b-factors set to 1.0.");
printMissingFields = false;
} else if (logger.isLoggable(Level.FINE)) {
logger.fine(" Missing occupancy and b-factors set to 1.0.");
}
}
newAtom = new Atom(0, name, altLoc, d, resName, resSeq, chainID, occupancy, tempFactor, segID);
newAtom.setHetero(true);
// Check if this is a modified residue.
if (modres.containsKey(resName.toUpperCase())) {
newAtom.setModRes(true);
}
returnedAtom = (Atom) activeMolecularAssembly.addMSNode(newAtom);
if (returnedAtom != newAtom) {
// A previously added atom has been retained.
atoms.put(serial, returnedAtom);
if (logger.isLoggable(Level.FINE)) {
logger.fine(returnedAtom + " has been retained over\n" + newAtom);
}
} else {
// The new atom has been added.
atoms.put(serial, newAtom);
newAtom.setXyzIndex(xyzIndex++);
}
break;
case CRYST1:
// =============================================================================
// The CRYST1 record presents the unit cell parameters, space group, and Z
// value. If the structure was not determined by crystallographic means, CRYST1
// simply provides the unitary values, with an appropriate REMARK.
//
// 7 - 15 Real(9.3) a a (Angstroms).
// 16 - 24 Real(9.3) b b (Angstroms).
// 25 - 33 Real(9.3) c c (Angstroms).
// 34 - 40 Real(7.2) alpha alpha (degrees).
// 41 - 47 Real(7.2) beta beta (degrees).
// 48 - 54 Real(7.2) gamma gamma (degrees).
// 56 - 66 LString sGroup Space group.
// 67 - 70 Integer z Z value.
// =============================================================================
double aaxis = new Double(line.substring(6, 15).trim());
double baxis = new Double(line.substring(15, 24).trim());
double caxis = new Double(line.substring(24, 33).trim());
double alpha = new Double(line.substring(33, 40).trim());
double beta = new Double(line.substring(40, 47).trim());
double gamma = new Double(line.substring(47, 54).trim());
int limit = 66;
if (line.length() < 66) {
limit = line.length();
}
String sg = line.substring(55, limit).trim();
properties.addProperty("a-axis", aaxis);
properties.addProperty("b-axis", baxis);
properties.addProperty("c-axis", caxis);
properties.addProperty("alpha", alpha);
properties.addProperty("beta", beta);
properties.addProperty("gamma", gamma);
properties.addProperty("spacegroup", SpaceGroup.pdb2ShortName(sg));
break;
case CONECT:
// =============================================================================
// 7 - 11 Integer serial Atom serial number
// 12 - 16 Integer serial Serial number of bonded atom
// 17 - 21 Integer serial Serial number of bonded atom
// 22 - 26 Integer serial Serial number of bonded atom
// 27 - 31 Integer serial Serial number of bonded atom
//
// CONECT records involving atoms for which the coordinates are not present
// in the entry (e.g., symmetry-generated) are not given.
// CONECT records involving atoms for which the coordinates are missing due
// to disorder, are also not provided.
// =============================================================================
conects.add(line);
break;
case LINK:
// =============================================================================
// The LINK records specify connectivity between residues that is not implied by
// the primary structure. Connectivity is expressed in terms of the atom names.
// They also include the distance associated with the each linkage following the
// symmetry operations at the end of each record.
// 13 - 16 Atom name1 Atom name.
// 17 Character altLoc1 Alternate location indicator.
// 18 - 20 Residue name resName1 Residue name.
// 22 Character chainID1 Chain identifier.
// 23 - 26 Integer resSeq1 Residue sequence number.
// 27 AChar iCode1 Insertion code.
// 43 - 46 Atom name2 Atom name.
// 47 Character altLoc2 Alternate location indicator.
// 48 - 50 Residue name resName2 Residue name.
// 52 Character chainID2 Chain identifier.
// 53 - 56 Integer resSeq2 Residue sequence number.
// 57 AChar iCode2 Insertion code.
// 60 - 65 SymOP sym1 Symmetry operator atom 1.
// 67 - 72 SymOP sym2 Symmetry operator atom 2.
// 74 – 78 Real(5.2) Length Link distance
// =============================================================================
Character a1 = line.charAt(16);
Character a2 = line.charAt(46);
if (a1 != a2) {
logger.info(format(" Ignoring LINK record as alternate locations do not match\n %s.", line));
break;
}
if (currentAltLoc == 'A') {
if ((a1 == ' ' || a1 == 'A') && (a2 == ' ' || a2 == 'A')) {
links.add(line);
}
} else if (a1 == currentAltLoc && a2 == currentAltLoc) {
links.add(line);
}
break;
case SSBOND:
// =============================================================================
if (currentAltLoc == 'A') {
ssbonds.add(line);
}
break;
case HELIX:
// =============================================================================
case SHEET:
// =============================================================================
// SHEET records are used to identify the position of sheets in the molecule.
// Sheets are both named and numbered. The residues where the sheet begins and
// ends are noted.
//
// 8 - 10 Integer strand Strand number which starts at 1 for each
// strand within a sheet and increases by one.
// 12 - 14 LString(3) sheetID Sheet identifier.
// 15 - 16 Integer numStrands Number of strands in sheet.
// 18 - 20 Residue name initResName Residue name of initial residue.
// 22 Character initChainID Chain identifier of initial residue
// in strand.
// 23 - 26 Integer initSeqNum Sequence number of initial residue
// in strand.
// 27 AChar initICode Insertion code of initial residue
// in strand.
// 29 - 31 Residue name endResName Residue name of terminal residue.
// 33 Character endChainID Chain identifier of terminal residue.
// 34 - 37 Integer endSeqNum Sequence number of terminal residue.
// 38 AChar endICode Insertion code of terminal residue.
// 39 - 40 Integer sense Sense of strand with respect to previous
// strand in the sheet. 0 if first strand,
// 1 if parallel,and -1 if anti-parallel.
// 42 - 45 Atom curAtom Registration. Atom name in current strand.
// 46 - 48 Residue name curResName Registration. Residue name in current strand
// 50 Character curChainId Registration. Chain identifier in
// current strand.
// 51 - 54 Integer curResSeq Registration. Residue sequence number
// in current strand.
// 55 AChar curICode Registration. Insertion code in
// current strand.
// 57 - 60 Atom prevAtom Registration. Atom name in previous strand.
// 61 - 63 Residue name prevResName Registration. Residue name in
// previous strand.
// 65 Character prevChainId Registration. Chain identifier in
// previous strand.
// 66 - 69 Integer prevResSeq Registration. Residue sequence number
// in previous strand.
// 70 AChar prevICode Registration. Insertion code in
// previous strand.
// =============================================================================
structs.add(line);
break;
// Currently, no handling in initial read.
case MODEL:
default:
break;
}
line = br.readLine();
}
br.close();
}
xyzIndex--;
setFileRead(true);
} catch (IOException e) {
logger.exiting(PDBFilter.class.getName(), "readFile", e);
return false;
}
/**
* Locate disulfide bonds; bond parameters are assigned below.
*/
List<Bond> ssBondList = locateDisulfideBonds(ssbonds);
/**
* Record the number of atoms read in from the PDB file before applying
* algorithms that may build new atoms.
*/
int pdbAtoms = activeMolecularAssembly.getAtomArray().length;
/**
* Build missing backbone atoms in loops.
*/
buildMissingResidues(xyzIndex);
/**
* Assign atom types. Missing side-chains atoms and missing hydrogens
* will be built in.
*/
assignAtomTypes();
/**
* Assign disulfide bonds parameters and log their creation.
*/
buildDisulfideBonds(ssBondList);
/**
* Finally, re-number the atoms if missing atoms were created.
*/
if (pdbAtoms != activeMolecularAssembly.getAtomArray().length) {
numberAtoms();
}
return true;
}
Also used :
ArrayList(java.util.ArrayList)
FileReader(java.io.FileReader)
IOException(java.io.IOException)
MissingHeavyAtomException(ffx.potential.bonded.BondedUtils.MissingHeavyAtomException)
IOException(java.io.IOException)
MissingAtomTypeException(ffx.potential.bonded.BondedUtils.MissingAtomTypeException)
Atom(ffx.potential.bonded.Atom)
BufferedReader(java.io.BufferedReader)
Bond(ffx.potential.bonded.Bond)
File(java.io.File)
Example 27 with Bond
use of ffx.potential.bonded.Bond in project ffx by mjschnie.
the class PDBFilter method buildDisulfideBonds.
/**
* Assign parameters to disulfide bonds.
*
* @param ssBondList
*/
private void buildDisulfideBonds(List<Bond> ssBondList) {
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 = BondType.sortKey(c);
BondType bondType = forceField.getBondType(key);
if (bondType == null) {
logger.severe(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(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());
}
}Example 28 with Bond
use of ffx.potential.bonded.Bond in project ffx by mjschnie.
the class PDBFilter method locateDisulfideBonds.
/**
* Locate disulfide bonds based on SSBOND records.
*
* @param ssbonds
*
* @return An ArrayList of Bond instances for SS-Bonds.
*/
private List<Bond> locateDisulfideBonds(List<String> ssbonds) {
List<Bond> ssBondList = new ArrayList<>();
for (String ssbond : ssbonds) {
// =============================================================================
try {
char c1ch = ssbond.charAt(15);
char c2ch = ssbond.charAt(29);
Polymer c1 = activeMolecularAssembly.getChain(String.format("%c", c1ch));
Polymer c2 = activeMolecularAssembly.getChain(String.format("%c", c2ch));
String origResNum1 = ssbond.substring(17, 21).trim();
char insChar1 = ssbond.charAt(21);
String origResNum2 = ssbond.substring(31, 35).trim();
char insChar2 = ssbond.charAt(35);
String pdbResNum1 = String.format("%c%s%c", c1ch, origResNum1, insChar1);
String pdbResNum2 = String.format("%c%s%c", c2ch, origResNum2, insChar2);
String resnum1 = pdbToNewResMap.get(pdbResNum1);
String resnum2 = pdbToNewResMap.get(pdbResNum2);
if (resnum1 == null) {
logger.warning(String.format(" Could not find residue %s for SS-bond %s", pdbResNum1, ssbond));
continue;
}
if (resnum2 == null) {
logger.warning(String.format(" Could not find residue %s for SS-bond %s", pdbResNum2, ssbond));
continue;
}
Residue r1 = c1.getResidue(Integer.parseInt(resnum1.substring(1)));
Residue r2 = c2.getResidue(Integer.parseInt(resnum2.substring(1)));
/*Residue r1 = c1.getResidue(Hybrid36.decode(4, ssbond.substring(17, 21)));
Residue r2 = c2.getResidue(Hybrid36.decode(4, ssbond.substring(31, 35)));*/
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));
}
if (SG2 == null) {
logger.warning(String.format(" SG atom 2 of SS-bond %s is null", ssbond));
}
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 = format("Ignoring [%s]\n due to distance %8.3f A.", ssbond, d);
logger.log(Level.WARNING, message);
}
} catch (Exception e) {
String message = format("Ignoring [%s]", ssbond);
logger.log(Level.WARNING, message, e);
}
}
return ssBondList;
}
Also used :
Residue(ffx.potential.bonded.Residue)
ArrayList(java.util.ArrayList)
Polymer(ffx.potential.bonded.Polymer)
Bond(ffx.potential.bonded.Bond)
Atom(ffx.potential.bonded.Atom)
MissingHeavyAtomException(ffx.potential.bonded.BondedUtils.MissingHeavyAtomException)
IOException(java.io.IOException)
MissingAtomTypeException(ffx.potential.bonded.BondedUtils.MissingAtomTypeException)
Example 29 with Bond
use of ffx.potential.bonded.Bond in project ffx by mjschnie.
the class XYZFilter method writeFileAsP1.
/**
* <p>
* writeFileAsP1</p>
*
* @param saveFile a {@link java.io.File} object.
* @param append a boolean.
* @param crystal a {@link ffx.crystal.Crystal} object.
* @return a boolean.
*/
public boolean writeFileAsP1(File saveFile, boolean append, Crystal crystal) {
if (saveFile == null) {
return false;
}
try {
File newFile = saveFile;
if (!append) {
newFile = version(saveFile);
}
activeMolecularAssembly.setFile(newFile);
activeMolecularAssembly.setName(newFile.getName());
FileWriter fw = new FileWriter(newFile, append);
BufferedWriter bw = new BufferedWriter(fw);
int nSymm = crystal.spaceGroup.symOps.size();
// XYZ File First Line
int numberOfAtoms = activeMolecularAssembly.getAtomList().size() * nSymm;
String output = format("%7d %s\n", numberOfAtoms, activeMolecularAssembly.toString());
bw.write(output);
if (!crystal.aperiodic()) {
Crystal uc = crystal.getUnitCell();
String params = String.format("%14.8f%14.8f%14.8f%14.8f%14.8f%14.8f\n", uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma);
bw.write(params);
}
Atom a2;
StringBuilder line;
StringBuilder[] lines = new StringBuilder[numberOfAtoms];
// XYZ File Atom Lines
Atom[] atoms = activeMolecularAssembly.getAtomArray();
double[] xyz = new double[3];
for (int iSym = 0; iSym < nSymm; iSym++) {
SymOp symOp = crystal.spaceGroup.getSymOp(iSym);
int indexOffset = iSym * atoms.length;
for (Atom a : atoms) {
int index = a.getIndex() + indexOffset;
String id = a.getAtomType().name;
if (vdwH) {
a.getRedXYZ(xyz);
} else {
xyz[0] = a.getX();
xyz[1] = a.getY();
xyz[2] = a.getZ();
}
crystal.applySymOp(xyz, xyz, symOp);
int type = a.getType();
line = new StringBuilder(format("%7d %3s%14.8f%14.8f%14.8f%6d", index, id, xyz[0], xyz[1], xyz[2], type));
for (Bond b : a.getBonds()) {
a2 = b.get1_2(a);
line.append(format("%8d", a2.getIndex() + indexOffset));
}
lines[index - 1] = line.append("\n");
}
}
try {
for (int i = 0; i < numberOfAtoms; i++) {
bw.write(lines[i].toString());
}
} catch (IOException e) {
String message = format(" Their was an unexpected error writing to %s.", getActiveMolecularSystem().toString());
logger.log(Level.WARNING, message, e);
return false;
}
bw.close();
fw.close();
} catch (IOException e) {
String message = format(" Their was an unexpected error writing to %s.", getActiveMolecularSystem().toString());
logger.log(Level.WARNING, message, e);
return false;
}
return true;
}
Also used :
SymOp(ffx.crystal.SymOp)
FileWriter(java.io.FileWriter)
IOException(java.io.IOException)
Atom(ffx.potential.bonded.Atom)
BufferedWriter(java.io.BufferedWriter)
Bond(ffx.potential.bonded.Bond)
File(java.io.File)
Crystal(ffx.crystal.Crystal)
Example 30 with Bond
use of ffx.potential.bonded.Bond in project ffx by mjschnie.
the class XYZFilter method readFile.
/**
* {@inheritDoc}
*
* Parse the XYZ File
*/
@Override
public boolean readFile() {
File xyzFile = activeMolecularAssembly.getFile();
if (forceField == null) {
logger.warning(format(" No force field is associated with %s.", xyzFile.toString()));
return false;
}
try {
FileReader fr = new FileReader(xyzFile);
BufferedReader br = new BufferedReader(fr);
String data = br.readLine();
// Read blank lines at the top of the file
while (data != null && data.trim().equals("")) {
data = br.readLine();
}
if (data == null) {
return false;
}
String[] tokens = data.trim().split(" +", 2);
int numberOfAtoms = Integer.parseInt(tokens[0]);
if (numberOfAtoms < 1) {
return false;
}
if (tokens.length == 2) {
getActiveMolecularSystem().setName(tokens[1]);
}
logger.info(format("\n Opening %s with %d atoms\n", xyzFile.getName(), numberOfAtoms));
// The header line is reasonable. Check for periodic box dimensions.
br.mark(10000);
data = br.readLine();
if (!readPBC(data, activeMolecularAssembly)) {
br.reset();
}
// Prepare to parse atom lines.
HashMap<Integer, Integer> labelHash = new HashMap<>();
int[] label = new int[numberOfAtoms];
int[][] bonds = new int[numberOfAtoms][8];
double[][] d = new double[numberOfAtoms][3];
boolean renumber = false;
atomList = new ArrayList<>();
// Loop over the expected number of atoms.
for (int i = 0; i < numberOfAtoms; i++) {
if (!br.ready()) {
return false;
}
data = br.readLine();
if (data == null) {
logger.warning(format(" Check atom %d in %s.", (i + 1), activeMolecularAssembly.getFile().getName()));
return false;
}
tokens = data.trim().split(" +");
if (tokens == null || tokens.length < 6) {
logger.warning(format(" Check atom %d in %s.", (i + 1), activeMolecularAssembly.getFile().getName()));
return false;
}
// Valid number of tokens, so try to parse this line.
label[i] = Integer.parseInt(tokens[0]);
// Check for valid atom numbering, or flag for re-numbering.
if (label[i] != i + 1) {
renumber = true;
}
String atomName = tokens[1];
d[i][0] = Double.parseDouble(tokens[2]);
d[i][1] = Double.parseDouble(tokens[3]);
d[i][2] = Double.parseDouble(tokens[4]);
int type = Integer.parseInt(tokens[5]);
AtomType atomType = forceField.getAtomType(Integer.toString(type));
if (atomType == null) {
StringBuilder message = new StringBuilder("Check atom type ");
message.append(type).append(" for Atom ").append(i + 1);
message.append(" in ").append(activeMolecularAssembly.getFile().getName());
logger.warning(message.toString());
return false;
}
Atom a = new Atom(i + 1, atomName, atomType, d[i]);
atomList.add(a);
// Bond Data
int numberOfBonds = tokens.length - 6;
for (int b = 0; b < 8; b++) {
if (b < numberOfBonds) {
int bond = Integer.parseInt(tokens[6 + b]);
bonds[i][b] = bond;
} else {
bonds[i][b] = 0;
}
}
}
// Check if this is an archive.
if (br.ready()) {
// Read past blank lines between archive files
data = br.readLine().trim();
while (data != null && data.equals("")) {
data = br.readLine().trim();
}
if (data != null) {
tokens = data.split(" +", 2);
if (tokens != null && tokens.length > 0) {
try {
int archiveNumberOfAtoms = Integer.parseInt(tokens[0]);
if (archiveNumberOfAtoms == numberOfAtoms) {
setType(FileType.ARC);
}
} catch (NumberFormatException e) {
tokens = null;
}
}
}
}
br.close();
fr.close();
// Try to renumber
if (renumber) {
for (int i = 0; i < numberOfAtoms; i++) {
if (labelHash.containsKey(label[i])) {
logger.warning(format(" Two atoms have the same index: %d.", label[i]));
return false;
}
labelHash.put(label[i], i + 1);
}
for (int i = 0; i < numberOfAtoms; i++) {
int j = -1;
while (j < 3 && bonds[i][++j] > 0) {
bonds[i][j] = labelHash.get(bonds[i][j]);
}
}
}
bondList = new ArrayList<>();
int[] c = new int[2];
for (int i = 1; i <= numberOfAtoms; i++) {
int a1 = i;
int j = -1;
while (j < 7 && bonds[i - 1][++j] > 0) {
int a2 = bonds[i - 1][j];
if (a1 < a2) {
if (a1 > numberOfAtoms || a1 < 1 || a2 > numberOfAtoms || a2 < 1) {
logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
return false;
}
// Check for bidirectional connection
boolean bidirectional = false;
int k = -1;
while (k < 7 && bonds[a2 - 1][++k] > 0) {
int a3 = bonds[a2 - 1][k];
if (a3 == a1) {
bidirectional = true;
break;
}
}
if (!bidirectional) {
logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
return false;
}
Atom atom1 = atomList.get(a1 - 1);
Atom atom2 = atomList.get(a2 - 1);
if (atom1 == null || atom2 == null) {
logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
return false;
}
Bond bond = new Bond(atom1, atom2);
c[0] = atom1.getAtomType().atomClass;
c[1] = atom2.getAtomType().atomClass;
String key = BondType.sortKey(c);
BondType bondType = forceField.getBondType(key);
if (bondType == null) {
logger.severe(format(" No BondType for key %s", key));
} else {
bond.setBondType(bondType);
}
bondList.add(bond);
}
}
}
return true;
} catch (IOException e) {
logger.severe(e.toString());
}
return false;
}
Also used :
BondType(ffx.potential.parameters.BondType)
HashMap(java.util.HashMap)
IOException(java.io.IOException)
Atom(ffx.potential.bonded.Atom)
BufferedReader(java.io.BufferedReader)
FileReader(java.io.FileReader)
AtomType(ffx.potential.parameters.AtomType)
Bond(ffx.potential.bonded.Bond)
File(java.io.File)
Aggregations
Bond (ffx.potential.bonded.Bond)38
Atom (ffx.potential.bonded.Atom)37
IOException (java.io.IOException)12
ArrayList (java.util.ArrayList)12
Polymer (ffx.potential.bonded.Polymer)10
Residue (ffx.potential.bonded.Residue)10
File (java.io.File)10
Angle (ffx.potential.bonded.Angle)9
Crystal (ffx.crystal.Crystal)8
MissingAtomTypeException (ffx.potential.bonded.BondedUtils.MissingAtomTypeException)8
MissingHeavyAtomException (ffx.potential.bonded.BondedUtils.MissingHeavyAtomException)8
MSNode (ffx.potential.bonded.MSNode)7
Molecule (ffx.potential.bonded.Molecule)7
AtomType (ffx.potential.parameters.AtomType)7
BufferedWriter (java.io.BufferedWriter)6
FileWriter (java.io.FileWriter)6
RestraintBond (ffx.potential.bonded.RestraintBond)5
Torsion (ffx.potential.bonded.Torsion)5
CoordRestraint (ffx.potential.nonbonded.CoordRestraint)5
MolecularAssembly (ffx.potential.MolecularAssembly)4
