001/* This class is based on the original FATCAT implementation by 002 * <pre> 003 * Yuzhen Ye & Adam Godzik (2003) 004 * Flexible structure alignment by chaining aligned fragment pairs allowing twists. 005 * Bioinformatics vol.19 suppl. 2. ii246-ii255. 006 * https://www.ncbi.nlm.nih.gov/pubmed/14534198 007 * </pre> 008 * 009 * Thanks to Yuzhen Ye and A. Godzik for granting permission to freely use and redistribute this code. 010 * 011 * This code may be freely distributed and modified under the 012 * terms of the GNU Lesser General Public Licence. This should 013 * be distributed with the code. If you do not have a copy, 014 * see: 015 * 016 * http://www.gnu.org/copyleft/lesser.html 017 * 018 * Copyright for this code is held jointly by the individual 019 * authors. These should be listed in @author doc comments. 020 * 021 * 022 * Created on Jun 17, 2009 023 * Created by Andreas Prlic - RCSB PDB 024 * 025 */ 026 027package org.biojava.nbio.structure.align.fatcat.calc; 028 029import org.biojava.nbio.structure.Atom; 030import org.biojava.nbio.structure.Group; 031import org.biojava.nbio.structure.StructureException; 032import org.biojava.nbio.structure.StructureTools; 033import org.biojava.nbio.structure.align.AFPTwister; 034import org.biojava.nbio.structure.align.fatcat.FatCat; 035import org.biojava.nbio.structure.align.model.AFP; 036import org.biojava.nbio.structure.align.model.AFPChain; 037import org.biojava.nbio.structure.align.util.AFPAlignmentDisplay; 038 039import java.util.List; 040 041 042/** A class that does calculations on an AFPChain 043 * 044 * @author Andreas Prlic 045 * 046 */ 047public class FatCatAligner 048{ 049 050 public static final boolean debug = false; 051 public static final boolean printTimeStamps = false; 052 053 AFPChain afpChain ; 054 Group[] twistedGroups; 055 056 057 058 public AFPChain getAfpChain() 059 { 060 return afpChain; 061 } 062 063 064 public Group[] getTwistedGroups() 065 { 066 067 return twistedGroups; 068 } 069 070 public void setTwistedGroups(Group[] twistedGroups) 071 { 072 this.twistedGroups = twistedGroups; 073 } 074 075 076 public void align(Atom[] ca1, Atom[] ca2, boolean doRigid, FatCatParameters params) throws StructureException{ 077 078 long tstart = System.currentTimeMillis(); 079 080 afpChain = new AFPChain(FatCat.algorithmName); 081 afpChain.setCa1Length(ca1.length); 082 afpChain.setCa2Length(ca2.length); 083 084 085 086 AFPCalculator.extractAFPChains(params, afpChain,ca1, ca2); 087 088 long cend = System.currentTimeMillis(); 089 090 if (printTimeStamps) 091 System.out.println("calculation took:" + (cend - tstart) + " ms."); 092 093 094 AFPCalculator.sortAfps(afpChain,ca1,ca2); 095 096 if ( printTimeStamps) { 097 long send = System.currentTimeMillis(); 098 099 100 System.out.println("sorting took:" + (send - cend) + " ms."); 101 } 102 103 if ( doRigid) 104 this.twistedGroups = rChainAfp(params, afpChain,ca1,ca2); 105 106 else { 107 this.twistedGroups = chainAfp(params,afpChain,ca1,ca2); 108 } 109 110 // long start = System.currentTimeMillis(); 111 long end = System.currentTimeMillis(); 112 afpChain.setCalculationTime(end-tstart); 113 if ( printTimeStamps) 114 System.out.println("TOTAL calc time: " + (end -tstart) / 1000.0); 115 116 } 117 118 119 120 121 /** runs rigid chaining process 122 * 123 */ 124 private static Group[] rChainAfp(FatCatParameters params, AFPChain afpChain, Atom[] ca1, Atom[] ca2) throws StructureException{ 125 params.setMaxTra(0); 126 afpChain.setMaxTra(0); 127 return chainAfp(params,afpChain,ca1,ca2); 128 } 129 130 /** 131 * run AFP chaining allowing up to maxTra flexible regions. 132 * Input is original coordinates. 133 * 134 */ 135 private static Group[] chainAfp(FatCatParameters params,AFPChain afpChain, Atom[] ca1, Atom[] ca2) throws StructureException{ 136 137 // we don;t want to rotate input atoms, do we? 138 Atom[] ca2clone = StructureTools.cloneAtomArray(ca2); 139 140 List<AFP> afpSet = afpChain.getAfpSet(); 141 142 if (debug) 143 System.out.println("entering chainAfp"); 144 int afpNum = afpSet.size(); 145 146 if ( afpNum < 1) 147 return new Group[0]; 148 149 long bgtime = System.currentTimeMillis(); 150 if(debug) { 151 System.out.println(String.format("total AFP %d\n", afpNum)); 152 } 153 154 //run AFP chaining 155 156 AFPChainer.doChainAfp(params,afpChain ,ca1,ca2); 157 158 int afpChainLen = afpChain.getAfpChainLen(); 159 160 if(afpChainLen < 1) { 161 162 afpChain.setShortAlign(true); 163 return new Group[0]; 164 } //very short alignment 165 166 long chaintime = System.currentTimeMillis(); 167 if(debug) { 168 169 System.out.println("Afp chaining: time " + (chaintime-bgtime)); 170 } 171 172 // do post processing 173 174 AFPPostProcessor.postProcess(params, afpChain,ca1,ca2); 175 176 // Optimize the final alignment 177 178 AFPOptimizer.optimizeAln(params, afpChain,ca1,ca2); 179 180 AFPOptimizer.blockInfo( afpChain); 181 182 AFPOptimizer.updateScore(params,afpChain); 183 184 AFPAlignmentDisplay.getAlign(afpChain,ca1,ca2); 185 186 Group[] twistedPDB = AFPTwister.twistPDB(afpChain, ca1, ca2clone); 187 188 SigEva sig = new SigEva(); 189 double probability = sig.calSigAll(params, afpChain); 190 afpChain.setProbability(probability); 191 double normAlignScore = sig.calNS(params,afpChain); 192 afpChain.setNormAlignScore(normAlignScore); 193 194 /* 195 196 SIGEVA sig; 197 probability = sig.calSigAll(maxTra, sparse, pro1Len, pro2Len, alignScore, totalRmsdOpt, optLength, blockNum - 1); 198 normAlignScore = sig.calNS(pro1Len, pro2Len, alignScore, totalRmsdOpt, optLength, blockNum - 1); 199 200 */ 201 202 //if(maxTra == 0) probability = sig.calSigRigid(pro1Len, pro2Len, alignScore, totalRmsdOpt, optLength); 203 //else probability = sig.calSigFlexi(pro1Len, pro2Len, alignScore, totalRmsdOpt, optLength, blockNum - 1); 204 205 if(debug) { 206 207 long nowtime = System.currentTimeMillis(); 208 long diff = nowtime - chaintime; 209 System.out.println("Alignment optimization: time "+ diff); 210 211 System.out.println("score: " + afpChain.getAlignScore()); 212 System.out.println("opt length: " + afpChain.getOptLength()); 213 System.out.println("opt rmsd: "+ afpChain.getTotalRmsdOpt()); 214 215 } 216 return twistedPDB; 217 218 } 219 220}