* *
* * This class provides a number of different "entry points" to traverse * the underlying Structure data. The one an application should choose * depeonds mostly on what the application wishes to accomplish. For * example, while a basic sequence viewer might simply walk the raw list * of residues (ie: by calling getResidueCount and getResidue in a loop) * another sequence viewer may want to obtain residues by walking each * chain (ie: by calling getChainCount, plus getChain and chain.getResidue * in a nested loop) so that it knows where the residues of one chain ends * and another begins. Again, its entirely up to the application. *
*
* @author John L. Moreland * @see edu.sdsc.mbt.StructureComponent * @see edu.sdsc.mbt.StructureComponentRegistry */ public class StructureMap { // PROGRAMMER NOTE #1: // Eventually, if/when the Structure class adds "set/add/remove" methods, // this class could be made to listen for StructureComponentEvent messages // and (in response) automatically update the internal map values. // PROGRAMMER NOTE #2: // Conformation interpretation - // a) Although Conformation records contain both start and end chain // identifiers, Conformations never actually cross chain boundaries. // b) Since chain idtentifiers are STRING values, there really is // no practical way to compute and iteratively walk intermediate values. // c) Start and end residue numbers are always constrained within one chain. // d) Residue numbers may be reused (restart at 1) for every chain, // or, may be continuous over the entire structure. So, don't count on // residue numbers for indexing purposes. // // Private state // // Stores a reference to the Structure. private Structure structure; // Primary StructureMap containers. private Vector atoms; // All Atoms in the Stucture. private Vector residues; // All Residues in the Stucture. private Vector fragments; // All Fragments in the Stucture. private Vector chains; // All Chains in the Stucture. private Vector ligands; // Only Ligand Residues. private Vector bonds; // All Bond objects added to this StructureMap. private Hashtable bondUniqueness; // Make sure Bond objects are unique. private Hashtable atomToBonds; // Find all Bonds connected to each Atom. private Hashtable atomToAtoms; //find all Atoms connected to each Atom private Hashtable bondToBonds; //find all Bonds connected to each Bond private Object value = new Object(); // Stores Chain object references by atom.chain_id value. private Hashtable chainById = null; // Stores Residue object references by atom.chain_id+atom_residue_id value. private Hashtable residueByChainAndResidueId = null; private double[] minCoordinate; private double[] maxCoordinate; private double[] center;//geometric center of the structure private HashMap phi = new HashMap();//Residue -> phi angle private HashMap psi = new HashMap();//Residue -> psi angle private Vector endResidues = new Vector();// private HashMap cMap = new HashMap(); private HashMap nMap = new HashMap(); private static final String defaultChainId = "_"; private boolean fillDisorderedGaps = false; // private boolean generateBondsByDistance = false; public boolean bondOrderImported = false;//indicates whether the source had bond order information public boolean bondOrderUserAssigned = false;//indicates whether the user has added bond order information public boolean chargesUserAssigned = false; public int MIN_ATOMS_FOR_FRAGMENTS = 20;//min number of atoms in the structure to run fragment detection private double bondCutoff = 1.9;//structure-specific bond cutoff distance (may be different for salt fragments) public static double hBondCutoff = 1.7;//cutoff for any bond including H /** * Used in the force field calculations */ private double[] gradient;//an array with three elements per atom // // Constructors // /** * Constructs a StructureMap object for a given Structure. */ public StructureMap( Structure structure ) { if ( structure == null ) throw new IllegalArgumentException( "null Structure" ); this.structure = structure; initialize( ); if ( Status.getOutputLevel() >= Status.LEVEL_DUMP ) print( ); // System.exit(1); } /** * Initialize all of the internal StructureMap state. * @param computeSecondaryStructure Should secondary structure be computed * (or loaded as it is defined from the source data)? */ private void initialize( ) { int atomCount = structure.getStructureComponentCount( StructureComponentRegistry.TYPE_ATOM ); int residueCount = structure.getStructureComponentCount( StructureComponentRegistry.TYPE_RESIDUE ); if ( atomCount > 0 ) atoms = new Vector( atomCount ); // All Atoms in the Stucture. else atoms = new Vector( ); // All Atoms in the Stucture. if ( residueCount > 0 ) residues = new Vector( residueCount ); // All Residues in the Stucture. else residues = new Vector( ); // All Residues in the Stucture. bonds = new Vector( ); // All Bond objects added to this StructureMap. bondUniqueness = new Hashtable( ); // Make sure Bond objects are unique. atomToBonds = new Hashtable( ); // Find all Bonds connected to each Atom. atomToAtoms = new Hashtable( ); // Find all Bonds connected to each Atom. bondToBonds = new Hashtable( ); // Find all Bonds connected to each Atom. chains = new Vector( ); // All Chains in the Stucture. ligands = new Vector( ); // Only Ligand Residues. fragments = new Vector( ); // All Fragments in the Stucture. // Stores Chain object references by atom.chain_id value. chainById = new Hashtable( ); // Stores Residue object references by atom.chain_id+atom_residue_id value. // Makes it possible to take an Atom record and find the appropriate residue // in which the Atom should be stored. residueByChainAndResidueId = new Hashtable( ); if ( atomCount > 0 ){ // Its molecule data processAtomRecords( ); //check whether Structure already has Bond records if (structure.getStructureComponentCount(StructureComponentRegistry.TYPE_BOND) > 0){ for (int i = 0; i < structure.getStructureComponentCount(StructureComponentRegistry.TYPE_BOND); i++){ Bond b = (Bond)structure.getStructureComponentByIndex(StructureComponentRegistry.TYPE_BOND, i); b.structure = structure; addBond(b); } bondOrderImported = true; } else{ // System.out.println("generating bonds"); generateBonds( ); } if (bonds.size() < atoms.size()/1.2){ detectBonds(); } else{ detectLigandBonds(); } checkResidueClassification(); //check if it makes sense to derive fragments if (atomCount > MIN_ATOMS_FOR_FRAGMENTS){ generateFragments( ); } extractLigands( ); } else if ( residueCount > 0 ){ // Its sequence data processResidueRecords( ); } else { // Its useless data } // System.out.println("bond count = " + bonds.size()); //set up the gradient array int size = getAtomCount()*3; gradient = new double[size]; //while passing through the residues, determine each atom's backbone assignment for (int i = 0; i < getResidueCount(); i++){ Residue r = getResidue(i); double[] sum = new double[3]; for (int j = 0; j < r.getAtomCount(); j++){ Atom a = r.getAtom(j); sum[0] += a.coordinate[0]; sum[1] += a.coordinate[1]; sum[2] += a.coordinate[2]; } sum[0] /= r.getAtomCount(); sum[1] /= r.getAtomCount(); sum[2] /= r.getAtomCount(); //another pass to determine which atom is the closest to the center double min = 100000; Atom target = null; for (int j = 0; j < r.getAtomCount(); j++){ Atom a = r.getAtom(j); double d = Algebra.distance(a.coordinate, sum); if (d < min){ min = d; target = a; } } if (target == null) target = r.getAtom(0); r.centerAtom = target; if (r.getClassification().equals(Residue.COMPOUND_AMINO_ACID) || r.getClassification().equals(Residue.COMPOUND_NUCLEIC_ACID)){ for (int j = 0; j < r.getAtomCount(); j++){ Atom a = r.getAtom(j); if (ProteinProperties.isProteinBackbone(a) || DNAProperties.isDNABackbone(a)){ a.backbone = true; } } } } } private void checkResidueClassification(){ //at the end, run a check for unnatural aminoacids, which should really still be marked as amino acids, rather than ligands // System.out.println("count = " + residues.size()); for (int i = 0; i < residues.size(); i++){ Residue rr = (Residue)residues.get(i); // System.out.println(i + ": rr = " + rr.getCompoundCode() + ", " + rr.getClassification()); if (rr.getClassification().equals(Residue.COMPOUND_LIGAND)){ // System.out.println(i + ": rr = " + rr.getCompoundCode() + ", " + rr.getClassification()); boolean nBond = false;//there is a N-C bond with a known amino acid boolean nCABond = false;//there is a N-CA bond within this residue boolean cBond = false;//there is a C-N bond with a known amino acid boolean cCABond = false;//there is a C-CA bond Atom ca = null; //check whether it's bonded to other amino acids and has typical backbone atom names for (int j = 0; j < rr.getAtomCount(); j++){ Atom a = rr.getAtom(j); if (a.name.equals("N")){ // System.out.println("Found N"); //check its bonds and find out whether one of them is connected to a CA and the other to a true amino acid Vector bonds = this.getBonds(a); // System.out.println("bonds = " + bonds); if (bonds != null){ for (int k = 0; k < bonds.size(); k++){ Bond b = (Bond)bonds.get(k); if (b.getAtom(0) == a){ Atom other = b.getAtom(1); if (other.residue == rr){ if (other.name.equals("CA")){ nCABond = true; } } else{ // System.out.println("other = " + other.name + ", classification = " + other.residue.getClassification()); if (other.name.equals("C") && other.residue.getClassification().equals(Residue.COMPOUND_AMINO_ACID)){ nBond = true; } } } if (b.getAtom(1) == a){ Atom other = b.getAtom(0); if (other.residue == rr){ if (other.name.equals("CA")){ nCABond = true; } } else{ // System.out.println("other = " + other.name + ", classification = " + other.residue.getClassification()); if (other.name.equals("C") && other.residue.getClassification().equals(Residue.COMPOUND_AMINO_ACID)){ nBond = true; } } } } } } else if (a.name.equals("C")){ // System.out.println("Found C"); Vector bonds = this.getBonds(a); if (bonds != null){ for (int k = 0; k < bonds.size(); k++){ Bond b = (Bond)bonds.get(k); if (b.getAtom(0) == a){ Atom other = b.getAtom(1); if (other.residue == rr){ if (other.name.equals("CA")){ cCABond = true; } } else{ if (other.name.equals("N") && other.residue.getClassification().equals(Residue.COMPOUND_AMINO_ACID)){ cBond = true; } } } if (b.getAtom(1) == a){ Atom other = b.getAtom(0); if (other.residue == rr){ if (other.name.equals("CA")){ cCABond = true; } } else{ if (other.name.equals("N") && other.residue.getClassification().equals(Residue.COMPOUND_AMINO_ACID)){ cBond = true; } } } } } } else if (a.name.equals("CA")){ ca = a; } } // System.out.println("nBond = " + nBond + ", cBond = " + cBond + ", cCABond = " + cCABond + ", ca = " + ca); if (nBond && cBond && nCABond && cCABond && ca != null){ rr.setClassification(Residue.COMPOUND_AMINO_ACID); rr.setAlphaAtom(ca); } } } } // // Initialization helper methods (private) // /** * Processes the atom records for the Structure and builds the hierarchy. *
*
*/ public Chain getChain( Atom atom ) { if ( atom == null ) return null; if ( residueByChainAndResidueId == null ) return null; return (Chain) chainById.get( atom.chain_id ); } /** * Get the Chain object given its ID. *
*/ public Chain getChain( String chainId ) { if ( chainId == null ) return null; if ( chainById == null ) return null; return (Chain) chainById.get( chainId ); } /** * Get a count of bonds contained in this StructureMap. */ public int getBondCount( ) { return bonds.size( ); } /** * Get the Bond at the specified index. *
*/ public Bond getBond( int bondIndex ) { return (Bond) bonds.elementAt( bondIndex ); } /** * Get the index of the specified Bond. *
*/
public int getBondIndex( Bond bond )
{
if ( bond == null ) return -1;
int bondCount = bonds.size( );
if ( bondCount < 1 ) return -1;
Atom atom0 = bond.getAtom( 0 );
if ( atom0 == null ) return -1;
Atom atom1 = bond.getAtom( 1 );
if ( atom1 == null ) return -1;
// Do a binary search to find the Bond index.
int low = 0;
int high = bondCount - 1;
while ( low <= high )
{
int mid = (low + high) / 2;
Bond bond2 = getBond( mid );
if ( bond == bond2 ) return mid;
Atom atom2 = bond2.getAtom( 0 );
Atom atom3 = bond2.getAtom( 1 );
if ( atom0.number < atom2.number )
{
high = mid - 1;
}
else if ( atom0.number == atom2.number )
{
if ( atom1.number < atom3.number )
high = mid - 1;
else
low = mid + 1;
}
else // if ( atom0.number > atom2.number )
{
low = mid + 1;
}
}
// Just in case we didn't find the Bond using the binary search,
// do a linear search of the bonds vector. This would only happen
// if the bond did not exist, or the bonds were not sorted.
for ( int i=0; i
*/
public void addBonds( Vector bondVector )
{
if ( bondVector == null ) return;
int bvCount = bondVector.size( );
for ( int i=0; i
*/
public void removeBond( int bondIndex )
{
removeBond( getBond( bondIndex ) );
}
/**
* Remove all Bond objects from the StructureMap.
*
*/
public void removeAllBonds( )
{
// Since we are taking a short cut below by blowing away the data directly,
// we will eventually have to generate events for Bond in "bonds" here...
bonds.clear( );
atomToBonds.clear( );
bondToBonds.clear();
atomToAtoms.clear();
bondUniqueness.clear( );
}
/**
* Return a Vector of all Bond objects connected to the given Atom object.
* NOTE: The Vector returned here is the internal copy (for speed purposes),
* so DO NOT MODIFY THE VECTOR!
*
*/
public Vector getBonds( Atom atom )
{
return (Vector) atomToBonds.get( atom );
}
/**
* Return a Vector of bonds that are neighbors of the argument
* @param bond
* @return
*/
public Vector getBonds (Bond bond){
return (Vector)bondToBonds.get(bond);
}
/**
* Return a Vector of all Bond objects connected to the given Atom objects.
*
*/
public Vector getBonds( Vector atomVector )
{
Hashtable uniqueBonds = new Hashtable( );
int atomCount = atomVector.size( );
for ( int a=0; a
*/
public Residue getLigandResidue( int ligandIndex )
{
if ( ligands == null ) return null;
return (Residue) ligands.elementAt( ligandIndex );
}
/**
* Return the residue count extracted from the Structure.
*/
public int getResidueCount( )
{
if ( residues == null ) return 0;
return residues.size( );
}
/**
* Return the start Atom index for a given residue index.
*
*/
public Residue getResidue( int residueIndex )
{
if ( residues == null ) return null;
return (Residue) residues.elementAt( residueIndex );
}
/**
* Get the Residue object by its chain and residue id.
*
*/
public Residue getResidue (String chainAndResidue){
return (Residue)residueByChainAndResidueId.get(chainAndResidue);
}
/**
* Get the index of the specified Residue.
*
*/
public int getResidueIndex( Residue residue )
throws IllegalArgumentException
{
if ( residue == null )
throw new IllegalArgumentException( "null residue" );
if ( residues == null )
throw new IllegalArgumentException( "no residues" );
// Do a binary search of the residues vector.
int low = 0;
int high = residues.size() - 1;
while ( low <= high )
{
int mid = (low + high) / 2;
Residue residue2 = (Residue) residues.elementAt( mid );
// If the residue matches, return the index.
if ( residue == residue2 )
{
return mid;
}
else
{
String chainId = residue.getChainId( );
String chainId2 = residue2.getChainId( );
int chainCompare = chainId.compareTo( chainId2 );
if ( chainCompare == 0 )
{
int residueId = residue.getResidueId( );
int residueId2 = residue2.getResidueId( );
if ( residueId == residueId2 )
{
if ( residue == residue2 )
return mid;
else
break; // Give up and do an exhastive linear search
}
else if ( residueId < residueId2 )
{
high = mid - 1;
}
else // ( residueId > residueId2 )
{
low = mid + 1;
}
}
else if ( chainCompare < 0 )
{
high = mid - 1;
}
else // ( chainCompare > 0 )
{
low = mid + 1;
}
}
}
// Try a much more expensive linear search in case the residue
// is out of order.
int residueCount = residues.size();
for ( int i=0; i
*/
public int getFragmentIndex( Fragment fragment )
{
if ( fragment == null )
throw new IllegalArgumentException( "null fragment" );
if ( fragments == null )
throw new IllegalArgumentException( "no fragments!" );
int fragmentCount = fragments.size( );
for ( int f=0; f
*/
public Vector getChildren( Object object )
{
if ( object == null ) return null;
if ( object instanceof Structure )
return getChains( );
if ( ! (object instanceof StructureComponent) )
return null;
StructureComponent sc = (StructureComponent) object;
String type = sc.getStructureComponentType( );
if ( type == StructureComponentRegistry.TYPE_CHAIN )
{
return ((Chain) sc).getFragments( );
}
else if ( type == StructureComponentRegistry.TYPE_FRAGMENT )
{
return ((Fragment) sc).getResidues( );
}
else if ( type == StructureComponentRegistry.TYPE_RESIDUE )
{
return ((Residue) sc).getAtoms();
}
else if ( type == StructureComponentRegistry.TYPE_ATOM )
{
return null;
}
else if ( type == StructureComponentRegistry.TYPE_BOND )
{
// Even though physically Atoms are children, they
// are not concidered children logically.
// Note that Atoms ARE parents to Bonds.
return null;
}
else
{
return null;
}
}
/**
* Return the chains from the Structure.
*/
public Vector getChains( )
{
return chains;
}
/**
* @return Returns the atoms.
*/
public final Vector getAtoms() {
return atoms;
}
/**
* @return Returns the atomToAtoms.
*/
public final Hashtable getAtomToAtoms() {
return atomToAtoms;
}
/**
* @return Returns the atomToBonds.
*/
public final Hashtable getAtomToBonds() {
return atomToBonds;
}
/**
* @return Returns the bondOrderImported.
*/
public final boolean isBondOrderImported() {
return bondOrderImported;
}
/**
* @return Returns the bondOrderUserAssigned.
*/
public final boolean isBondOrderUserAssigned() {
return bondOrderUserAssigned;
}
/**
* @return Returns the bonds.
*/
public final Vector getBonds() {
return bonds;
}
/**
* @return Returns the bondToBonds.
*/
public final Hashtable getBondToBonds() {
return bondToBonds;
}
/**
* @return Returns the bondUniqueness.
*/
public final Hashtable getBondUniqueness() {
return bondUniqueness;
}
/**
* @return Returns the chainById.
*/
public final Hashtable getChainById() {
return chainById;
}
/**
* @return Returns the chargesUserAssigned.
*/
public final boolean isChargesUserAssigned() {
return chargesUserAssigned;
}
/**
* @return Returns the fragments.
*/
public final Vector getFragments() {
return fragments;
}
/**
* @return Returns the ligands.
*/
public final Vector getLigands() {
return ligands;
}
/**
* @return Returns the residueByChainAndResidueId.
*/
public final Hashtable getResidueByChainAndResidueId() {
return residueByChainAndResidueId;
}
/**
* @return Returns the residues.
*/
public final Vector getResidues() {
return residues;
}
/**
* Return the StructureComponent parent objects.
* Atom->Residue->Fragment->Chain->Structure->null, or,
* Bond->{Atom,Atom}->null, or,
* null.
*
*/
public Vector getParents( Object object )
{
if ( object == null )
return null;
if ( object instanceof Structure )
return null;
if ( ! (object instanceof StructureComponent) )
return null;
StructureComponent sc = (StructureComponent) object;
String type = sc.getStructureComponentType( );
Vector parents = null;
if ( type == StructureComponentRegistry.TYPE_ATOM )
{
parents = new Vector( );
parents.add( getResidue( (Atom) sc ) );
}
else if ( type == StructureComponentRegistry.TYPE_BOND )
{
parents = new Vector( );
Bond bond = (Bond) sc;
parents.add( bond.getAtom(0) );
parents.add( bond.getAtom(1) );
}
else if ( type == StructureComponentRegistry.TYPE_RESIDUE )
{
parents = new Vector( );
parents.add( ((Residue)sc).getFragment() );
}
else if ( type == StructureComponentRegistry.TYPE_FRAGMENT )
{
parents = new Vector( );
parents.add( ((Fragment)sc).getChain() );
}
else if ( type == StructureComponentRegistry.TYPE_CHAIN )
{
parents = new Vector( );
parents.add( structure );
}
return parents;
}
/**
* Should we "fill in" the disordered residue gaps with random coil
* when we load fragments?
*
* @see #getFillDisorderedGaps( )
*
*/
public void setFillDisorderedGaps( boolean state )
{
fillDisorderedGaps = state;
}
/**
* Should we "fill in" the disordered residue gaps with random coil
* when we load fragments?
*
* @see #setFillDisorderedGaps( boolean state )
*
*/
public boolean getFillDisorderedGaps( )
{
return fillDisorderedGaps;
}
/**
* This method walks over all bonds in the structure and calculates orders
* based on their geometry
*
*/
public void deriveBondOrders(){
for (int i = 0; i < bonds.size(); i++){
Bond bond = (Bond)bonds.get(i);
Atom a0 = bond.getAtom(0);
Atom a1 = bond.getAtom(1);
if (a0.element == 6 && a1.element == 6){//carbons
//both are carbons
}
else if (a0.element == 6){
double[] center = a0.coordinate;
double[] tip1 = a1.coordinate;
Atom third = null;
Vector atoms = getAtoms(a0);
for (int j = 0; j < atoms.size(); j++){
if ((Atom)atoms.get(j) != a0 && (Atom)atoms.get(j) != a1){
third = (Atom)atoms.get(j);
break;
}
}
if (third == null){
//nothing else is connected to this carbon
//this could be in files with no hydrogens for a methyl group
//if so, go to the other carbon to get the bond order
}
else{
//get the angle
double[] tip2 = third.coordinate;
double angle = Algebra.angle(tip1, center, tip2);
if (angle > 105 && angle <= 112){
//sp3 carbon
bond.setOrder(1.0f);
}
else if (angle > 112 && angle <= 118){
/**
* TODO check here for a possible aromatic case
*/
bond.setOrder(2.0f);
}
else {
bond.setOrder(3.0f);
}
}
}
else if (a1.element == 6){
}
else if (a0.element == 1 || a1.element == 1){
bond.setOrder(1.0f);
}
else{
bond.setOrder(1.0f);
}
}
}
public boolean inRing(Atom atom, int ringSize){
Vector atoms = getAtoms(atom);
visitedAtoms.put(atom, "");
boolean root = false;//root atom was seen as a neighbor
cycles = 1;
for (int i = 0; i < atoms.size(); i++){
Atom a = (Atom)atoms.get(i);
hasNeighbor(a, atom);
visitedAtoms.put(a, "");
}
return false;
}
private int cycles = 0;
private HashMap visitedAtoms = new HashMap();
private boolean hasNeighbor(Atom atom, Atom root){
Vector atoms = getAtoms(atom);
for (int i = 0; i < atoms.size(); i++){
Atom a = (Atom)atoms.get(0);
if (visitedAtoms.containsKey(a))continue;
if (a == root){
return true;
}
else{
return hasNeighbor(a, root);
}
}
return false;
}
/**
* This method removes an atom from the structure
*
*/
public synchronized void removeAtom(Atom atom){
if (!atoms.contains(atom)) return;
Chain chain = getChain( atom );
//get its bonds and remove them as well
Vector bonds = getBonds(atom);
if (bonds != null){
for (int i = bonds.size()-1; i >= 0; i--){
Bond b = (Bond)bonds.get(i);
removeBond(b);
}
}
atomToBonds.remove(atom);
atoms.remove(atom);
atomToAtoms.remove(atom);
//check all values in atomToAtoms to remove references to this Atom
/* Set keys = atomToAtoms.keySet();
Iterator it = keys.iterator();
while (it.hasNext()){
Atom a = (Atom)it.next();
Vector v = (Vector)atomToAtoms.get(a);
if (v.contains(atom)){
v.remove(atom);
}
}
*/
//remove the atom from its residue as well
Residue r = getResidue(atom);
r.getAtoms().remove(atom);
/**
* TODO handling of structure components in the Structure object
*/
//check if there is a residue now with no atoms
if (r.getAtomCount() == 0){
//this residue should be deleted and removed from the chain record
residues.remove(r);
chain.removeResidue(r);
getStructureStyles().registerStructureComponentRemoval(r);
//further check to see if the chain itself is still there
if (chain.getResidueCount() == 0){
//delete the chain
chains.remove(chain);
}
}
}
/**
* This method incorporates a new structure (fragment) into this structure
* by importing all data structures and updating names
* @param map
*/
public void addStructure(StructureMap map){
// System.out.println("this = " + this + ", map = " + map);
Structure s = this.getStructure();
for (int i = 0; i < map.getAtomCount(); i++){
Atom a = map.getAtom(i);
a.number = getNextAtomNumber();
a.structure = s;
}
for (int i = 0; i < map.getBondCount(); i++){
Bond b = map.getBond(i);
b.structure = s;
}
atoms.addAll(map.getAtoms());
bonds.addAll(map.getBonds());
//check the chains
//if a chain has a lettered name, add it as is. otherwise, merge the content with the corresponding chain
//of the parent structure. no need to have multiple chains for merged fragments of small molecules
//no two chains can have the same name. If it's a letter, reassign the name to be the next in alphabet after the last one existing
//in the parent structure. If it's a "_" and there is no "_" in the parent, add it as is. If there are "_" in both
//structures, merge them into one chain and renumber residues
//first, get a list of existing chains in both parts
Vector chains1 = new Vector();
Chain dashChain1 = null;
for (int i = 0; i < this.getChainCount(); i++){
Chain chain = this.getChain(i);
if (chain.getChainId().equals("_")) dashChain1 = chain;
chains1.add(chain.getChainId());
}
for (int i = 0; i < map.getChainCount(); i++){
Chain chain = map.getChain(i);
if (chain.getChainId().equals("_")){
if (dashChain1 == null){
//there is no dash chain in the parent structure, so just add this one
chains.add(chain);
chain.structure = s;
for (int j = 0; j < chain.getResidueCount(); j++){
chain.getResidue(j).structure = s;
residues.add(chain.getResidue(j));
}
chainById.put("_", chain);
residueByChainAndResidueId.putAll(map.getResidueByChainAndResidueId());
}
else{
//content of the new chain must be merged with the parent chain
//the only data that gets added is residue by chain and residue id hash entries
//in addition, the residues get renumbered
chain.structure = s;
Hashtable labels = new Hashtable();
for (int k = 0; k < chain.getResidueCount(); k++){
Residue r = chain.getResidue(k);
int n = getNextResidueNumber();
r.setResidueId(n);
for (int nn = 0; nn < r.getAtomCount(); nn++){
r.getAtom(nn).residue_id = n;
}
r.structure = s;
String label = "_" + n;
labels.put(label, r);
dashChain1.addResidue(r);
residues.add(r);
}
residueByChainAndResidueId.putAll(labels);
}
}
else{
//check whether the new chain has a name that is already taken
String chainId = chain.getChainId();
boolean pass = true;
for (int j = 0; j < this.getChainCount(); j++){
Chain c = this.getChain(j);
if (c.getChainId().equals(chainId)){
pass = false;
break;
}
}
if (pass){
//add the chain as is, with no changes to its name
chains.add(chain);
chain.structure = s;
for (int j = 0; j < chain.getResidueCount(); j++){
chain.getResidue(j).structure = s;
residues.add(chain.getResidue(j));
}
chainById.put(chainId, chain);
residueByChainAndResidueId.putAll(map.getResidueByChainAndResidueId());
}
else{
String id = getNextChainId();
chain.structure = s;
for (int j = 0; j < chain.getResidueCount(); j++){
Residue r = chain.getResidue(j);
r.structure = s;
int ri = r.getResidueId();
residueByChainAndResidueId.put(id+ri, r);
for (int k = 0; k < r.getAtomCount(); k++){
Atom a = r.getAtom(k);
a.chain_id = id;
a.residue_id = ri;
}
residues.add(r);
}
chainById.put(id, chain);
chains.add(chain);
}
}
}
fragments.addAll(map.getFragments());
atomToAtoms.putAll(map.getAtomToAtoms());
atomToBonds.putAll(map.getAtomToBonds());
bondToBonds.putAll(map.getBondToBonds());
bondUniqueness.putAll(map.getBondUniqueness());
ligands.addAll(map.getLigands());
for (int i = 0; i < fragments.size(); i++){
Fragment f = (Fragment)fragments.get(i);
f.structure = s;
}
//recalculate new structure bounds
calculateBounds();
}
public void calculateBounds(){
double maxX=0.0, minX=0.0, maxY=0.0, minY=0.0, maxZ=0.0, minZ=0.0;
for ( int i = 0; i < getAtomCount(); i++ )
{
Atom atom = getAtom(i);
if ( i == 0 )
{
maxX = atom.coordinate[0];
minX = atom.coordinate[0];
maxY = atom.coordinate[1];
minY = atom.coordinate[1];
maxZ = atom.coordinate[2];
minZ = atom.coordinate[2];
}
if ( atom.coordinate[0] >= maxX ) maxX = atom.coordinate[0];
if ( atom.coordinate[0] <= minX ) minX = atom.coordinate[0];
if ( atom.coordinate[1] >= maxY ) maxY = atom.coordinate[1];
if ( atom.coordinate[1] <= minY ) minY = atom.coordinate[1];
if ( atom.coordinate[2] >= maxZ ) maxZ = atom.coordinate[2];
if ( atom.coordinate[2] <= minZ ) minZ = atom.coordinate[2];
}
//calculate the bounds of the structure
maxCoordinate = new double[] {maxX, maxY, maxZ};
minCoordinate = new double[] {minX, minY, minZ};
}
/**
* Get all combination of three atoms connected through the central atom and forming an angle
* Each item is an array of three Atoms
* @return
*/
public Vector getAngleConnections(){
Vector out = new Vector();
//walk through atoms and pick those that have at least two bonds
//for those with more, record all combinations
//IMPORTANT: the central atom should be in the middle of the set, since iteration starts
//at the periphery
for (int i = 0; i < getAtomCount(); i++){
Atom a = getAtom(i);
Vector atoms = getAtoms(a);
int number = atoms.size();
if (number <= 1) continue;
else if (number == 2){
Atom[] set = new Atom[3];
set[0] = (Atom)atoms.get(0);
set[1] = a;
set[2] = (Atom)atoms.get(1);
out.add(set);
}
else{
//there are more than two atoms - record all possible pairwise combinations
//of the peripheral atoms
for (int j = 0; j < number; j++){
for (int k = j+1; k < number; k++){
Atom[] set = new Atom[3];
set[0] = (Atom)atoms.get(j);
set[1] = a;
set[2] = (Atom)atoms.get(k);
out.add(set);
}
}
}
}
return out;
}
public Vector getTorsionConnections(){
Vector out = new Vector();
//loop through all bonds and find those with terminal bonds attached
//record the atoms
for (int i = 0; i < getBondCount(); i++){
Bond b = getBond(i);
Atom a1 = b.getAtom(0);
Atom a2 = b.getAtom(1);
Vector atoms1 = getAtoms(a1);
Vector atoms2 = getAtoms(a2);
if (atoms1.size() <= 1 || atoms2.size() <= 1) continue;//this is a terminal atom, so no torsion at this bond
//walk through the atoms and record them
Atom[] indices1 = new Atom[atoms1.size() - 1];//minus the atom in the parent bond
int counter = 0;
for (int j = 0; j < atoms1.size(); j++){
if (atoms1.get(j) == a2){
continue;
}
indices1[counter++] = (Atom)atoms1.get(j);
}
//walk through the atoms and record their indices
Atom[] indices2 = new Atom[atoms2.size() - 1];//minus the atom in the parent bond
counter = 0;
for (int j = 0; j < atoms2.size(); j++){
if (atoms2.get(j) == a1){
continue;
}
indices2[counter++] = (Atom)atoms2.get(j);
}
//now walk through the collected indices and
for (int n = 0; n < indices1.length; n++){
for (int m = 0; m < indices2.length; m++){
Atom[] set = new Atom[4];
set[0] = indices1[n];
set[1] = a1;
set[2] = a2;
set[3] = indices2[m];
out.add(set);
}
}
}
return out;
}
/**
* Returns a Vector of arrays with pairs of atoms that don't form bonds
* @return
*/
public Vector getNonBondedConnections(){
Vector out = new Vector();
for (int i = 0; i < getAtomCount(); i++){
for (int j = i+1; j < getAtomCount(); j++){
Atom a1 = getAtom(i);
Atom a2 = getAtom(j);
Vector atoms = getAtoms(a1);
if (!atoms.contains(a2)){
Atom[] set = new Atom[2];
set[0] = a1;
set[1] = a2;
out.add(set);
}
}
}
return out;
}
public boolean forceFieldTypesSet = false;
/**
* This method calculates phi-psi angles for the protein structure and saves the values in the
* hashmaps
* If the structure contains no amino acids, false is returned. Otherwise, the result is true
*
*/
public boolean computePhiPsi(){
//first pass: check for amino acid content and save out C and N atoms for each Residue
cMap = new HashMap();
nMap = new HashMap();
boolean peptide = false;
for (int i = 0; i < getResidueCount(); i++){
Residue r = getResidue(i);
if (ProteinProperties.isAminoAcid(r.getCompoundCode())){
peptide = true;
//get the necessary atoms by looping through the list
for (int j = 0; j < r.getAtomCount(); j++){
Atom a = r.getAtom(j);
if (a.name.equals("C")){
cMap.put(r, a);
}
else if (a.name.equals("N")){
nMap.put(r, a);
}
}
}
}
if (!peptide) return false;
//go over all residues, this time calculating the angles
for (int i = 0; i < getResidueCount(); i++){
Residue r = getResidue(i);
if (ProteinProperties.isAminoAcid(r.getCompoundCode())){
Atom CA = null;
Atom localC = null;
Atom localN = null;
for (int j = 0; j < r.getAtomCount(); j++){
Atom a = r.getAtom(j);
if (a.name.equals("CA")){
CA = a;
}
else if (a.name.equals("C")){
localC = a;
}
else if (a.name.equals("N")){
localN = a;
}
}
try{
Atom C = (Atom)cMap.get(getResidue(i-1));
double angle = Algebra.dihedralAngle(C.coordinate, localN.coordinate, CA.coordinate, localC.coordinate);
phi.put(r, new Double(angle));
}
catch (Exception ex){
if (!endResidues.contains(r)) endResidues.add(r);
phi.put(r, new Double(-180));
}
try{
Atom N = (Atom)nMap.get(getResidue(i+1));
double angle = Algebra.dihedralAngle(localN.coordinate, CA.coordinate, localC.coordinate, N.coordinate);
psi.put(r, new Double(angle));
}
catch (Exception ex){
if (!endResidues.contains(r)) endResidues.add(r);
psi.put(r, new Double(180));
}
}
}
return true;
}
public void computePhiPsi(Residue r, HashMap transformed){
//if transformed coordinates are null, use atom records, otherwise, use the hash
int index = getResidueIndex(r);
if (index < 0) return;
if (ProteinProperties.isAminoAcid(r.getCompoundCode())){
//get the necessary atoms by looping through the list
Atom C = null;
Atom N = null;
Atom CA = null;
for (int j = 0; j < r.getAtomCount(); j++){
Atom a = r.getAtom(j);
if (a.name.equals("C")){
C = a;
}
else if (a.name.equals("N")){
N = a;
}
else if (a.name.equals("CA")){
CA = a;
}
}
//get the other atoms from the neighboring residues
try{
Residue last = getResidue(index-1);
double[] coordinateLastC = null;
double[] coordinateN = null;
double[] coordinateCA = null;
double[] coordinateC = null;
if (transformed != null && transformed.containsKey(cMap.get(last))){
coordinateLastC = (double[])transformed.get(cMap.get(last));
}
else{
coordinateLastC = ((Atom)cMap.get(last)).coordinate;
}
if (transformed != null && transformed.containsKey(N)){
coordinateN = (double[])transformed.get(N);
}
else{
coordinateN = N.coordinate;
}
if (transformed != null && transformed.containsKey(CA)){
coordinateCA = (double[])transformed.get(CA);
}
else{
coordinateCA = CA.coordinate;
}
if (transformed != null && transformed.containsKey(C)){
coordinateC = (double[])transformed.get(C);
}
else{
coordinateC = C.coordinate;
}
double angle = Algebra.dihedralAngle(coordinateLastC, coordinateN, coordinateCA, coordinateC);
phi.put(r, new Double(angle));
}
catch (Exception e){
phi.put(r, new Double(-180));
}
try{
Residue next = getResidue(index+1);
double[] coordinateNextN = null;
double[] coordinateN = null;
double[] coordinateCA = null;
double[] coordinateC = null;
if (transformed != null && transformed.containsKey(nMap.get(next))){
coordinateNextN = (double[])transformed.get(nMap.get(next));
}
else{
coordinateNextN = ((Atom)nMap.get(next)).coordinate;
}
if (transformed != null && transformed.containsKey(N)){
coordinateN = (double[])transformed.get(N);
}
else{
coordinateN = N.coordinate;
}
if (transformed != null && transformed.containsKey(CA)){
coordinateCA = (double[])transformed.get(CA);
}
else{
coordinateCA = CA.coordinate;
}
if (transformed != null && transformed.containsKey(C)){
coordinateC = (double[])transformed.get(C);
}
else{
coordinateC = C.coordinate;
}
double angle = Algebra.dihedralAngle(coordinateN, coordinateCA, coordinateC, coordinateNextN);
psi.put(r, new Double(angle));
}
catch (Exception e){
psi.put(r, new Double(180));
}
}
}
public HashMap getPhiMap(){
return phi;
}
public HashMap getPsiMap(){
return psi;
}
public Vector getEndResidues(){
return endResidues;
}
public double[] getPsiPhi(Residue r){
if (phi.containsKey(r) && psi.containsKey(r)){
return new double[]{ ((Double)phi.get(r)).doubleValue(),((Double)psi.get(r)).doubleValue() };
}
return null;
}
}