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 * Created on Apr 7, 2010
021 * Author: Andreas Prlic
022 *
023 */
024
025package org.biojava.nbio.structure.align.util;
026
027
028import javax.vecmath.Matrix4d;
029
030import org.biojava.nbio.structure.Atom;
031import org.biojava.nbio.structure.Calc;
032import org.biojava.nbio.structure.StructureException;
033import org.biojava.nbio.structure.align.model.AFPChain;
034import org.biojava.nbio.structure.geometry.SuperPositions;
035import org.slf4j.Logger;
036import org.slf4j.LoggerFactory;
037
038public class AFPChainScorer {
039
040        private static final Logger logger = LoggerFactory.getLogger(AFPChainScorer.class);
041
042
043        public  static double getTMScore(AFPChain align, Atom[] ca1, Atom[] ca2) throws StructureException
044        {
045                return getTMScore(align, ca1, ca2, true);
046        }
047
048        public  static double getTMScore(AFPChain align, Atom[] ca1, Atom[] ca2, boolean normalizeMin) throws StructureException
049        {
050                if ( align.getNrEQR() == 0)
051                        return -1;
052
053
054                // Create new arrays for the subset of atoms in the alignment.
055                Atom[] ca1aligned = new Atom[align.getOptLength()];
056                Atom[] ca2aligned = new Atom[align.getOptLength()];
057                int pos=0;
058                int[] blockLens = align.getOptLen();
059                int[][][] optAln = align.getOptAln();
060                assert(align.getBlockNum() <= optAln.length);
061
062                for(int block=0;block< align.getBlockNum();block++) {
063
064                        if ( ! ( blockLens[block] <= optAln[block][0].length)) {
065                                logger.warn("AFPChainScorer getTMScore: errors reconstructing alignment block [" + block + "]. Length is " + blockLens[block] + " but should be <=" + optAln[block][0].length);
066                        }
067
068                        for(int i=0;i<blockLens[block];i++) {
069                                int pos1 = optAln[block][0][i];
070                                int pos2 = optAln[block][1][i];
071                                Atom a1 = ca1[pos1];
072                                Atom a2 = (Atom) ca2[pos2].clone();
073
074                                ca1aligned[pos] = a1;
075                                ca2aligned[pos] = a2;
076                                pos++;
077                        }
078                }
079
080                // this can happen when we load an old XML serialization which did not support modern ChemComp representation of modified residues.
081                if ( pos != align.getOptLength()){
082                        logger.warn("AFPChainScorer getTMScore: Problems reconstructing alignment! nr of loaded atoms is " + pos + " but should be " + align.getOptLength());
083                        // we need to resize the array, because we allocated too many atoms earlier on.
084                        ca1aligned = (Atom[]) resizeArray(ca1aligned, pos);
085                        ca2aligned = (Atom[]) resizeArray(ca2aligned, pos);
086                }
087                //Superimpose
088                Matrix4d trans = SuperPositions.superpose(Calc.atomsToPoints(ca1aligned),
089                                Calc.atomsToPoints(ca2aligned));
090
091                Calc.transform(ca2aligned, trans);
092
093                return Calc.getTMScore(ca1aligned, ca2aligned, ca1.length, ca2.length, normalizeMin);
094        }
095
096        /**
097         * Reallocates an array with a new size, and copies the contents
098         * of the old array to the new array.
099         * @param oldArray  the old array, to be reallocated.
100         * @param newSize   the new array size.
101         * @return          A new array with the same contents.
102         */
103        private static Object resizeArray (Object oldArray, int newSize) {
104                int oldSize = java.lang.reflect.Array.getLength(oldArray);
105                @SuppressWarnings("rawtypes")
106                Class elementType = oldArray.getClass().getComponentType();
107                Object newArray = java.lang.reflect.Array.newInstance(
108                                elementType,newSize);
109                int preserveLength = Math.min(oldSize,newSize);
110                if (preserveLength > 0)
111                        System.arraycopy (oldArray,0,newArray,0,preserveLength);
112                return newArray; }
113
114}