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}