Source code

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.quaternary;
022
023import org.biojava.nbio.structure.*;
024
025import java.util.ArrayList;
026import java.util.Arrays;
027import java.util.List;
028
029
030
031
032
033/**
034 * @author Peter Rose
035 *
036 *
037 */
038public class BioAssemblyTools {
039
040
041        /**
042         * Checks if the passed in expression is a unary operator expression
043         * Example: (1,2,3) or (1-60) are unary operator expressions
044         *          (1-60)(61-88) is a binary operator expression, representing
045         *          a cartesian product of the two parenthesised lists
046         *
047         * @param expression
048         * @return true if expression is a unary operator expression
049         */
050        public static boolean isUnaryExpression(String expression) {
051                int first = expression.indexOf("(");
052                int last = expression.lastIndexOf("(");
053                if (first < 0 || last < 0) {
054                        return true;
055                }
056                return ! (first == 0 && last > first);
057        }
058
059        public static List<String> parseUnaryOperatorExpression(String operatorExpression) throws IllegalArgumentException {
060                return parseSubExpression(operatorExpression);
061        }
062
063        private static List<String> parseSubExpression(String expression) throws IllegalArgumentException {
064                // remove parenthesis, if any
065                String tmp = expression.replace("(", "");
066                tmp = tmp.replace(")", "");
067
068                // separate the operators
069                List<String> components = null;
070                try {
071                        components = Arrays.asList(tmp.split(","));
072                } catch (Exception e) {
073                        throw new IllegalArgumentException("Invalid oper_expression: " + expression);
074                }
075
076                // expand ranges if present, i.e. 1-60 -> 1, 2, 3, ..., 60
077                List<String> operators = new ArrayList<String>();
078                for (String component : components) {
079                        if (component.contains("-")) {
080                                operators.addAll(expandRange(component));
081                        } else {
082                                operators.add(component);
083                        }
084                }
085                return operators;
086        }
087
088        /**
089         * Expands a range expression, i.e. (1-6) to a list 1,2,3,4,5,6
090         * @param expression the expression to be expanded
091         * @return list of items in range
092         * @throws IllegalArgumentException
093         */
094        private static List<String> expandRange(String expression) throws IllegalArgumentException {
095                int first = 0;
096                int last = 0;
097                try {
098                        String[] range = expression.split("-");
099                        first = Integer.parseInt(range[0]);
100                        last = Integer.parseInt(range[1]);
101                } catch (Exception e) {
102                        throw new IllegalArgumentException("Invalid range specification in oper_expression: " + expression);
103                }
104
105                List<String> expandedExpression = new ArrayList<String>(last-first+1);
106                for (int i = first; i <= last; i++) {
107                        expandedExpression.add(String.valueOf(i));
108                }
109                return expandedExpression;
110        }
111
112        public static List<OrderedPair<String>> parseBinaryOperatorExpression(String expression)
113                        throws IllegalArgumentException {
114                // split operator expression, i.e. (1,2,3)(4,5) into two subexpressions
115                String[] subExpressions = null;
116                try {
117                        subExpressions = expression.split("\\)\\(");
118                } catch (Exception e) {
119                        throw new IllegalArgumentException("Invalid oper_expression: " + expression);
120                }
121                if (subExpressions.length != 2) {
122                        throw new IllegalArgumentException("Invalid oper_expression: " + expression);
123                }
124                List<String> leftSide = parseSubExpression(subExpressions[0]);
125                List<String> rightSide = parseSubExpression(subExpressions[1]);
126
127                // form the cartesian product of the two lists
128                CartesianProduct<String> product = new CartesianProduct<String>(leftSide, rightSide);
129                return product.getOrderedPairs();
130        }
131
132        public static double[][]  getBiologicalMoleculeBounds(Structure asymStructure,List<BiologicalAssemblyTransformation> transformations) {
133                final double coordinateBounds[][] = new double[2][3];
134                coordinateBounds[0][0] = Double.MAX_VALUE;  // min x
135                coordinateBounds[0][1] = Double.MAX_VALUE;  // min y
136                coordinateBounds[0][2] = Double.MAX_VALUE;  // min z
137                coordinateBounds[1][0] = Double.MIN_VALUE;  // max x
138                coordinateBounds[1][1] = Double.MIN_VALUE;  // max y
139                coordinateBounds[1][2] = Double.MIN_VALUE;  // max z
140
141                double[] transformedCoords = new double[3];
142
143                Atom[] atoms = StructureTools.getAllAtomArray(asymStructure);
144
145                for ( Atom a : atoms) {
146
147                        Chain c = a.getGroup().getChain();
148                        String intChainID = c.getInternalChainID();
149                        if ( intChainID == null)
150                                intChainID = c.getChainID();
151
152                        for (BiologicalAssemblyTransformation m: transformations) {
153                                if ( ! m.getChainId().equals(intChainID))
154                                        continue;
155                                double[] coords = a.getCoords();
156                                transformedCoords[0] = coords[0];
157                                transformedCoords[1] = coords[1];
158                                transformedCoords[2] = coords[2];
159
160                                if (transformedCoords[0] < coordinateBounds[0][0] ) {
161                                        coordinateBounds[0][0] = transformedCoords[0];  // min x
162                                }
163
164                                if (transformedCoords[1] < coordinateBounds[0][1] ) {
165                                        coordinateBounds[0][1] = transformedCoords[1];  // min y
166                                }
167
168                                if (transformedCoords[2] < coordinateBounds[0][2] ) {
169                                        coordinateBounds[0][2] = transformedCoords[2];  // min z
170                                }
171
172                                if (transformedCoords[0] > coordinateBounds[1][0] ) {
173                                        coordinateBounds[1][0] = transformedCoords[0];  // max x
174                                }
175
176                                if (transformedCoords[1] > coordinateBounds[1][1] ) {
177                                        coordinateBounds[1][1] = transformedCoords[1];  // max y
178                                }
179
180                                if (transformedCoords[2] > coordinateBounds[1][2] ) {
181                                        coordinateBounds[1][2] = transformedCoords[2];  // max z
182                                }
183                        }
184                }
185                return coordinateBounds;
186        }
187        public static double[][] getAtomCoordinateBounds(Structure s){
188
189                org.biojava.nbio.structure.Atom[] atoms = StructureTools.getAllAtomArray(s);
190                int atomCount = atoms.length;
191                final double coordinateBounds[][] = new double[2][3];
192                if ( atomCount <= 0 ) {
193                        return coordinateBounds;
194                }
195
196                Atom a = atoms[0];
197
198                coordinateBounds[0][0] = a.getX();  // min x
199                coordinateBounds[0][1] = a.getY();  // min y
200                coordinateBounds[0][2] = a.getZ();  // min z
201                coordinateBounds[1][0] = a.getX();  // max x
202                coordinateBounds[1][1] = a.getY();  // max y
203                coordinateBounds[1][2] = a.getZ();  // max z
204
205                for ( int i=1; i<atomCount; i++ )
206                {
207                        a =atoms[i];
208
209                        if ( a.getX() < coordinateBounds[0][0] ) {
210                                coordinateBounds[0][0] = a.getX();  // min x
211                        }
212
213                        if ( a.getY() < coordinateBounds[0][1] ) {
214                                coordinateBounds[0][1] = a.getY();  // min y
215                        }
216
217                        if ( a.getZ() < coordinateBounds[0][2] ) {
218                                coordinateBounds[0][2] = a.getZ();  // min z
219                        }
220
221                        if ( a.getX() > coordinateBounds[1][0] ) {
222                                coordinateBounds[1][0] = a.getX();  // max x
223                        }
224
225                        if ( a.getY() > coordinateBounds[1][1] ) {
226                                coordinateBounds[1][1] = a.getY();  // max y
227                        }
228
229                        if ( a.getZ() > coordinateBounds[1][2] ) {
230                                coordinateBounds[1][2] = a.getZ();  // max z
231                        }
232                }
233
234                return coordinateBounds;
235        }
236
237        /**
238         * Returns the maximum extend of the structure in the x, y, or z direction.
239         * @param structure
240         * @return maximum extend
241         */
242        public static double getMaximumExtend( final Structure structure ) {
243                double[][] bounds = getAtomCoordinateBounds(structure);
244                double xMax = Math.abs(bounds[0][0] - bounds[1][0]);
245                double yMax = Math.abs(bounds[0][1] - bounds[1][1]);
246                double zMax = Math.abs(bounds[0][2] - bounds[1][2]);
247                return Math.max(xMax, Math.max(yMax, zMax));
248        }
249
250        /**
251         * Returns the maximum extend of the biological molecule in the x, y, or z direction.
252         * @param structure
253         * @return maximum extend
254         */
255        public static double getBiologicalMoleculeMaximumExtend( final Structure structure,List<BiologicalAssemblyTransformation> transformations ) {
256                double[][] bounds = getBiologicalMoleculeBounds(structure, transformations);
257                double xMax = Math.abs(bounds[0][0] - bounds[1][0]);
258                double yMax = Math.abs(bounds[0][1] - bounds[1][1]);
259                double zMax = Math.abs(bounds[0][2] - bounds[1][2]);
260                return Math.max(xMax, Math.max(yMax, zMax));
261        }
262
263        /**
264         * Returns the centroid of the biological molecule.
265         * @param structure
266         * @return centroid
267         * @throws IllegalArgumentException if structure is null
268         */
269
270        public static double[] getBiologicalMoleculeCentroid( final Structure asymUnit,List<BiologicalAssemblyTransformation> transformations ) throws IllegalArgumentException {
271                if ( asymUnit == null ) {
272                        throw new IllegalArgumentException( "null structure" );
273                }
274
275                Atom[] atoms = StructureTools.getAllAtomArray(asymUnit);
276                int atomCount = atoms.length;
277                double[] centroid = new double[3];
278
279                if ( atomCount <= 0 ) {
280                        return centroid;
281                }
282
283                if ( transformations.size() == 0) {
284                        return Calc.getCentroid(atoms).getCoords();
285
286                }
287
288
289
290                int count = 0;
291                double[] transformedCoordinate = new double[3];
292
293                for (int i = 0; i < atomCount; i++)
294                {
295                        Atom atom = atoms[i];
296                        Chain chain = atom.getGroup().getChain();
297                        String intChainID = chain.getInternalChainID();
298                        if ( intChainID == null)
299                                intChainID = chain.getChainID();
300
301
302                        for (BiologicalAssemblyTransformation m: transformations) {
303                                if (!  m.getChainId().equals(intChainID))
304                                        continue;
305
306                                double[] coords = atom.getCoords();
307                                transformedCoordinate[0] = coords[0];
308                                transformedCoordinate[1] = coords[1];
309                                transformedCoordinate[2] = coords[2];
310                                m.transformPoint(transformedCoordinate);
311                                centroid[0] += transformedCoordinate[0];
312                                centroid[1] += transformedCoordinate[1];
313                                centroid[2] += transformedCoordinate[2];
314                                count++;
315                        }
316                }
317
318
319
320                centroid[0] /= count;
321                centroid[1] /= count;
322                centroid[2] /= count;
323
324                return centroid;
325        }
326
327        /** reduce a structure to a Calpha representation only
328         *
329         * @param orig
330         * @return
331         * @deprecated Use the more generic {@link #getReducedStructure(Structure)}
332         */
333        @Deprecated
334        public static Structure getReducedCAStructure(Structure orig){
335                Structure s = new StructureImpl();
336                s.setPDBHeader(orig.getPDBHeader());
337                for ( Chain c : orig.getChains()){
338
339                        Chain c1 = new ChainImpl();
340                        c1.setChainID(c.getChainID());
341                        s.addChain(c1);
342
343                        for (Group g : c.getAtomGroups()){
344
345                                try {
346                                        Atom a = g.getAtom(StructureTools.CA_ATOM_NAME);
347                                        if ( a != null){
348
349                                                Group g1 = (Group)g.clone();
350                                                g1.clearAtoms();
351                                                g1.addAtom(a);
352                                                c1.addGroup(g1);
353
354                                        }
355                                } catch (Exception e){}
356                        }
357
358                }
359                return s;
360        }
361        /** reduce a structure to a single-atom representation (e.g. CA atoms
362         *
363         * @param orig
364         * @return
365         * @since Biojava 4.1.0
366         */
367        public static Structure getReducedStructure(Structure orig){
368                Structure s = new StructureImpl();
369                s.setPDBHeader(orig.getPDBHeader());
370                for ( Chain c : orig.getChains()){
371
372                        Chain c1 = new ChainImpl();
373                        c1.setChainID(c.getChainID());
374                        s.addChain(c1);
375
376                        for (Group g : c.getAtomGroups()){
377
378                                try {
379                                        Atom a = null;
380                                        switch(g.getType()) {
381                                        case AMINOACID:
382                                                a = g.getAtom(StructureTools.CA_ATOM_NAME);
383                                                break;
384                                        case NUCLEOTIDE:
385                                                a = g.getAtom(StructureTools.NUCLEOTIDE_REPRESENTATIVE);
386                                                break;
387                                        default:
388                                                //omit group
389                                        }
390                                        if ( a != null){
391
392                                                Group g1 = (Group)g.clone();
393                                                g1.clearAtoms();
394                                                g1.addAtom(a);
395                                                c1.addGroup(g1);
396
397                                        }
398                                } catch (Exception e){}
399                        }
400
401                }
402                return s;
403        }
404
405}