001/*
002 *                    BioJava development code
003 *
004 * This code may be freely distributed and modified under the
005 * terms of the GNU Lesser General Public Licence.  This should
006 * be distributed with the code.  If you do not have a copy,
007 * see:
008 *
009 *      http://www.gnu.org/copyleft/lesser.html
010 *
011 * Copyright for this code is held jointly by the individual
012 * authors.  These should be listed in @author doc comments.
013 *
014 * For more information on the BioJava project and its aims,
015 * or to join the biojava-l mailing list, visit the home page
016 * at:
017 *
018 *      http://www.biojava.org/
019 *
020 */
021package org.biojava.nbio.structure.symmetry.core;
022
023import javax.vecmath.Point3d;
024
025import org.jgrapht.UndirectedGraph;
026import org.jgrapht.graph.DefaultEdge;
027import org.jgrapht.graph.SimpleGraph;
028import org.slf4j.Logger;
029import org.slf4j.LoggerFactory;
030
031import java.util.List;
032
033public class SubunitGraph {
034
035        private static final Logger logger = LoggerFactory
036                        .getLogger(SubunitGraph.class);
037
038        private static final double DISTANCE_CUTOFF = 8;
039        private static final int MIN_CONTACTS = 10;
040        private List<Point3d[]> caCoords = null;
041
042        public SubunitGraph(List<Point3d[]> caCoords) {
043                this.caCoords = caCoords;
044        }
045
046        public UndirectedGraph<Integer, DefaultEdge> getProteinGraph() {
047                int n = caCoords.size();
048
049                // add vertex for each chain center
050                UndirectedGraph<Integer, DefaultEdge> graph = new SimpleGraph<Integer, DefaultEdge>(DefaultEdge.class);
051                for (int i = 0; i < n; i++) {
052                        graph.addVertex(i);
053                }
054
055                // add edges if there are 10 or more contact of Calpha atoms
056                for (int i = 0; i < n - 1; i++) {
057                        for (int j = i + 1; j < n; j++) {
058                                int numContacts = calcContactNumber(caCoords.get(i), caCoords.get(j));
059                                logger.debug("Calpha contacts between subunits {},{}: {}", i, j, numContacts);
060                                if (numContacts >= MIN_CONTACTS) {
061                                        graph.addEdge(i, j);
062                                }
063                        }
064                }
065
066                return graph;
067        }
068
069        private int calcContactNumber(Point3d[] a, Point3d[] b) {
070                double distCutoffSq = DISTANCE_CUTOFF * DISTANCE_CUTOFF;
071                int contacts = 0;
072                for (Point3d pa : a) {
073                        for (Point3d pb : b) {
074                                if (pa.distanceSquared(pb) < distCutoffSq) {
075                                        contacts++;
076                                }
077                        }
078                }
079                return contacts;
080        }
081}