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.geometry; 022 023import javax.vecmath.Matrix4d; 024import javax.vecmath.Point3d; 025import javax.vecmath.Vector3d; 026 027import org.biojava.nbio.structure.jama.Matrix; 028 029/** 030 * Utility operations on Point3d. 031 * 032 * @author Aleix Lafita 033 * @since 5.0.0 034 * 035 */ 036public class CalcPoint { 037 038 /** Prevent instantiation */ 039 private CalcPoint() { 040 } 041 042 /** 043 * Center a cloud of points. This means subtracting the {@lin 044 * #centroid(Point3d[])} of the cloud to each point. 045 * 046 * @param x 047 * array of points. Point objects will be modified 048 */ 049 public static void center(Point3d[] x) { 050 Point3d center = centroid(x); 051 center.negate(); 052 translate(new Vector3d(center), x); 053 } 054 055 /** 056 * Calculate the centroid of the point cloud. 057 * 058 * @param x 059 * array of points. Point objects will not be modified 060 * @return centroid as Point3d 061 */ 062 public static Point3d centroid(Point3d[] x) { 063 Point3d center = new Point3d(); 064 for (Point3d p : x) { 065 center.add(p); 066 } 067 center.scale(1.0 / x.length); 068 return center; 069 } 070 071 /** 072 * Transform all points with a 4x4 transformation matrix. 073 * 074 * @param rotTrans 075 * 4x4 transformation matrix 076 * @param x 077 * array of points. Point objects will be modified 078 */ 079 public static void transform(Matrix4d rotTrans, Point3d[] x) { 080 for (Point3d p : x) { 081 rotTrans.transform(p); 082 } 083 } 084 085 /** 086 * Translate all points with a translation vector. 087 * 088 * @param trans 089 * the translation vector to apply 090 * @param x 091 * array of points. Point objects will be modified 092 */ 093 public static void translate(Vector3d trans, Point3d[] x) { 094 for (Point3d p : x) { 095 p.add(trans); 096 } 097 } 098 099 /** 100 * Clone an array of points. 101 * 102 * @param x 103 * original array of points. Point objects will not be modified 104 * @return new array of points, identical clone of x 105 */ 106 public static Point3d[] clonePoint3dArray(Point3d[] x) { 107 Point3d[] clone = new Point3d[x.length]; 108 for (int i = 0; i < x.length; i++) { 109 clone[i] = new Point3d(x[i]); 110 } 111 return clone; 112 } 113 114 /* 115 * Peter can you document this method? TODO 116 * 117 * @param moved 118 * @param fixed 119 * @return 120 */ 121 public static Matrix formMatrix(Point3d[] a, Point3d[] b) { 122 double xx = 0.0, xy = 0.0, xz = 0.0; 123 double yx = 0.0, yy = 0.0, yz = 0.0; 124 double zx = 0.0, zy = 0.0, zz = 0.0; 125 126 for (int i = 0; i < a.length; i++) { 127 xx += a[i].x * b[i].x; 128 xy += a[i].x * b[i].y; 129 xz += a[i].x * b[i].z; 130 yx += a[i].y * b[i].x; 131 yy += a[i].y * b[i].y; 132 yz += a[i].y * b[i].z; 133 zx += a[i].z * b[i].x; 134 zy += a[i].z * b[i].y; 135 zz += a[i].z * b[i].z; 136 } 137 138 double[][] f = new double[4][4]; 139 f[0][0] = xx + yy + zz; 140 f[0][1] = zy - yz; 141 f[1][0] = f[0][1]; 142 f[1][1] = xx - yy - zz; 143 f[0][2] = xz - zx; 144 f[2][0] = f[0][2]; 145 f[1][2] = xy + yx; 146 f[2][1] = f[1][2]; 147 f[2][2] = yy - zz - xx; 148 f[0][3] = yx - xy; 149 f[3][0] = f[0][3]; 150 f[1][3] = zx + xz; 151 f[3][1] = f[1][3]; 152 f[2][3] = yz + zy; 153 f[3][2] = f[2][3]; 154 f[3][3] = zz - xx - yy; 155 156 return new Matrix(f); 157 } 158 159 /** 160 * Returns the TM-Score for two superimposed sets of coordinates Yang Zhang 161 * and Jeffrey Skolnick, PROTEINS: Structure, Function, and Bioinformatics 162 * 57:702–710 (2004) 163 * 164 * @param x 165 * coordinate set 1 166 * @param y 167 * coordinate set 2 168 * @param lengthNative 169 * total length of native sequence 170 * @return 171 */ 172 public static double TMScore(Point3d[] x, Point3d[] y, int lengthNative) { 173 174 if (x.length != y.length) { 175 throw new IllegalArgumentException( 176 "Point arrays are not of the same length."); 177 } 178 179 double d0 = 1.24 * Math.cbrt(x.length - 15.0) - 1.8; 180 double d0Sq = d0 * d0; 181 182 double sum = 0; 183 for (int i = 0; i < x.length; i++) { 184 sum += 1.0 / (1.0 + x[i].distanceSquared(y[i]) / d0Sq); 185 } 186 187 return sum / lengthNative; 188 } 189 190 /* 191 * Needs documentation! 192 * 193 * @param x 194 * 195 * @param y 196 * 197 * @return 198 */ 199 public static double GTSlikeScore(Point3d[] x, Point3d[] y) { 200 201 if (x.length != y.length) { 202 throw new IllegalArgumentException( 203 "Point arrays are not of the same length."); 204 } 205 206 int contacts = 0; 207 208 for (Point3d px : x) { 209 double minDist = Double.MAX_VALUE; 210 211 for (Point3d py : y) { 212 minDist = Math.min(minDist, px.distanceSquared(py)); 213 } 214 215 if (minDist > 64) 216 continue; 217 contacts++; 218 219 if (minDist > 16) 220 continue; 221 contacts++; 222 223 if (minDist > 4) 224 continue; 225 contacts++; 226 227 if (minDist > 1) 228 continue; 229 contacts++; 230 } 231 232 return contacts * 25.0 / x.length; 233 } 234 235 /** 236 * Calculate the RMSD of two point arrays, already superposed. 237 * 238 * @param x 239 * array of points superposed to y 240 * @param y 241 * array of points superposed to x 242 * @return RMSD 243 */ 244 public static double rmsd(Point3d[] x, Point3d[] y) { 245 246 if (x.length != y.length) { 247 throw new IllegalArgumentException( 248 "Point arrays are not of the same length."); 249 } 250 251 double sum = 0.0; 252 for (int i = 0; i < x.length; i++) { 253 sum += x[i].distanceSquared(y[i]); 254 } 255 return Math.sqrt(sum / x.length); 256 } 257 258 /* 259 * Needs documentation! 260 * 261 * @param x 262 * 263 * @param y 264 * 265 * @return 266 */ 267 public static double rmsdMin(Point3d[] x, Point3d[] y) { 268 269 if (x.length != y.length) { 270 throw new IllegalArgumentException( 271 "Point arrays are not of the same length."); 272 } 273 274 double sum = 0.0; 275 for (int i = 0; i < x.length; i++) { 276 double minDist = Double.MAX_VALUE; 277 for (int j = 0; j < y.length; j++) { 278 minDist = Math.min(minDist, x[i].distanceSquared(y[j])); 279 } 280 sum += minDist; 281 } 282 return Math.sqrt(sum / x.length); 283 } 284 285 /* 286 * Needs documentation! 287 * 288 * @param x 289 * 290 * @param y 291 * 292 * @param maxDistance 293 * 294 * @return 295 */ 296 public static int contacts(Point3d[] x, Point3d[] y, double maxDistance) { 297 int contacts = 0; 298 for (int i = 0; i < x.length; i++) { 299 double minDist = Double.MAX_VALUE; 300 for (int j = 0; j < y.length; j++) { 301 minDist = Math.min(minDist, x[i].distanceSquared(y[j])); 302 } 303 if (minDist < maxDistance * maxDistance) { 304 contacts++; 305 } 306 } 307 return contacts; 308 } 309 310}