1 /* $RCSfile$ 2 * $Author: hansonr $ 3 * $Date: 2007-05-18 08:19:45 -0500 (Fri, 18 May 2007) $ 4 * $Revision: 7742 $ 5 6 * 7 * Copyright (C) 2003-2005 The Jmol Development Team 8 * 9 * Contact: jmol-developers@lists.sf.net 10 * 11 * This library is free software; you can redistribute it and/or 12 * modify it under the terms of the GNU Lesser General Public 13 * License as published by the Free Software Foundation; either 14 * version 2.1 of the License, or (at your option) any later version. 15 * 16 * This library is distributed in the hope that it will be useful, 17 * but WITHOUT ANY WARRANTY; without even the implied warranty of 18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 19 * Lesser General Public License for more details. 20 * 21 * You should have received a copy of the GNU Lesser General Public 22 * License along with this library; if not, write to the Free Software 23 * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. 24 */ 25 26 package org.jmol.symmetry; 27 28 import java.util.Hashtable; 29 import java.util.Map; 30 31 import javajs.util.Lst; 32 import javajs.util.M4; 33 import javajs.util.Measure; 34 import javajs.util.P3; 35 import javajs.util.P4; 36 import javajs.util.PT; 37 import javajs.util.Quat; 38 import javajs.util.SB; 39 import javajs.util.T3; 40 import javajs.util.V3; 41 42 import org.jmol.api.SymmetryInterface; 43 import javajs.util.BS; 44 import org.jmol.modelset.Atom; 45 import org.jmol.modelset.ModelSet; 46 import org.jmol.script.T; 47 import org.jmol.util.Escape; 48 import org.jmol.util.Logger; 49 50 /** 51 * A class to handle requests for information about space groups and symmetry 52 * operations. 53 * 54 * Two entry points, both from Symmetry: 55 * 56 * getSymopInfo 57 * 58 * getSpaceGroupInfo 59 * 60 * 61 * 62 */ 63 public class SymmetryDesc { 64 65 private ModelSet modelSet; 66 SymmetryDesc()67 public SymmetryDesc() { 68 // for reflection 69 } 70 set(ModelSet modelSet)71 public SymmetryDesc set(ModelSet modelSet) { 72 this.modelSet = modelSet; 73 return this; 74 } 75 76 private final static String[] keys = { "xyz", "xyzOriginal", "label", 77 null /*draw*/, "fractionalTranslation", "cartesianTranslation", 78 "inversionCenter", null /*point*/, "axisVector", "rotationAngle", 79 "matrix", "unitTranslation", "centeringVector", "timeReversal", "plane", 80 "_type", "id", "cif2", "xyzCanonical" }; 81 82 private final static int KEY_DRAW = 3; 83 private final static int KEY_POINT = 7; 84 85 ////// "public" methods //////// 86 87 /** 88 * 89 * @param iAtom 90 * @param xyz 91 * @param op 92 * @param translation TODO 93 * @param pt 94 * @param pt2 95 * @param id 96 * @param type 97 * @param scaleFactor 98 * @param nth 99 * @param options 0 or T.offset 100 * @return "" or a bitset of matching atoms, or 101 */ getSymopInfo(int iAtom, String xyz, int op, P3 translation, P3 pt, P3 pt2, String id, int type, float scaleFactor, int nth, int options)102 Object getSymopInfo(int iAtom, String xyz, int op, P3 translation, P3 pt, P3 pt2, 103 String id, int type, float scaleFactor, int nth, int options) { 104 if (type == 0) 105 type = getType(id); 106 Object ret = (type == T.atoms ? new BS() : ""); 107 if (iAtom < 0) 108 return ret; 109 110 // get model symmetry 111 112 int iModel = modelSet.at[iAtom].mi; 113 SymmetryInterface uc = modelSet.am[iModel].biosymmetry; 114 if (uc == null && (uc = modelSet.getUnitCell(iModel)) == null) 115 return ret; 116 117 // generally get the result from getSymmetryInfo 118 119 if (type != T.draw || op != Integer.MAX_VALUE) { 120 return getSymmetryInfo((Symmetry) uc, iModel, iAtom, (Symmetry) uc, xyz, 121 op, translation, pt, pt2, id, type, scaleFactor, nth, options); 122 } 123 124 // draw SPACEGROUP 125 126 String s = ""; 127 M4[] ops = uc.getSymmetryOperations(); 128 if (ops != null) { 129 if (id == null) 130 id = "sg"; 131 int n = ops.length; 132 for (op = 0; op < n; op++) 133 s += (String) getSymmetryInfo((Symmetry) uc, iModel, iAtom, 134 (Symmetry) uc, xyz, op, translation, pt, pt2, id + op, T.draw, scaleFactor, nth, options); 135 } 136 return s; 137 } 138 139 @SuppressWarnings("unchecked") getSpaceGroupInfo(Symmetry sym, int modelIndex, String sgName, int symOp, P3 pt1, P3 pt2, String drawID, float scaleFactor, int nth, boolean isFull, boolean isForModel, int options, SymmetryInterface cellInfo)140 Map<String, Object> getSpaceGroupInfo(Symmetry sym, int modelIndex, 141 String sgName, int symOp, P3 pt1, 142 P3 pt2, String drawID, 143 float scaleFactor, int nth, 144 boolean isFull, boolean isForModel, int options, SymmetryInterface cellInfo) { 145 Map<String, Object> info = null; 146 boolean isStandard = (pt1 == null && drawID == null && nth <= 0); 147 boolean isBio = false; 148 String sgNote = null; 149 boolean haveName = (sgName != null && sgName.length() > 0); 150 boolean haveRawName = (haveName && sgName.indexOf("[--]") >= 0); 151 if (isForModel || !haveName) { 152 boolean saveModelInfo = (isStandard && symOp == 0); 153 if (modelIndex < 0) 154 modelIndex = (pt1 instanceof Atom ? ((Atom) pt1).mi 155 : modelSet.vwr.am.cmi); 156 if (modelIndex < 0) 157 sgNote = "no single current model"; 158 else if (cellInfo == null && !(isBio = (cellInfo = modelSet.am[modelIndex].biosymmetry) != null) 159 && (cellInfo = modelSet.getUnitCell(modelIndex)) == null) 160 sgNote = "not applicable"; 161 if (sgNote != null) { 162 info = new Hashtable<String, Object>(); 163 info.put("spaceGroupInfo", ""); 164 info.put("spaceGroupNote", sgNote); 165 info.put("symmetryInfo", ""); 166 } else if (isStandard) { 167 info = (Map<String, Object>) modelSet.getInfo(modelIndex, 168 "spaceGroupInfo"); 169 } 170 // created once 171 if (info != null) 172 return info; 173 174 // full check 175 info = new Hashtable<String, Object>(); 176 sgName = cellInfo.getSpaceGroupName(); 177 SymmetryOperation[] ops = (SymmetryOperation[]) cellInfo 178 .getSymmetryOperations(); 179 SpaceGroup sg = (isBio ? ((Symmetry) cellInfo).spaceGroup : null); 180 String slist = (haveRawName ? "" : null); 181 int opCount = 0; 182 if (ops != null) { 183 if (isBio) 184 sym.spaceGroup = SpaceGroup.getNull(false, false, false); 185 else 186 sym.setSpaceGroup(false); 187 // check to make sure that new group has been created magnetic or not 188 if (ops[0].timeReversal != 0) 189 ((SymmetryOperation) sym.getSpaceGroupOperation(0)).timeReversal = 1; 190 Object[][] infolist = new Object[ops.length][]; 191 String sops = ""; 192 for (int i = 0, nop = 0; i < ops.length && nop != nth; i++) { 193 194 SymmetryOperation op = ops[i]; 195 boolean isNewIncomm = (i == 0 && op.xyz.indexOf("x4") >= 0); 196 int iop = (!isNewIncomm && sym.getSpaceGroupOperation(i) != null ? i : isBio ? sym 197 .addBioMoleculeOperation(sg.finalOperations[i], false) : sym 198 .addSpaceGroupOperation("=" + op.xyz, i + 1)); 199 if (iop < 0) 200 continue; 201 String xyzOriginal = op.xyzOriginal; 202 op = (SymmetryOperation) sym.getSpaceGroupOperation(i); 203 if (op == null) 204 continue; 205 op.xyzOriginal = xyzOriginal; 206 if (op.timeReversal != 0 || op.modDim > 0) 207 isStandard = false; 208 if (slist != null) 209 slist += ";" + op.xyz; 210 Object[] ret = (symOp > 0 && symOp - 1 != iop ? null 211 : createInfoArray(op, cellInfo, pt1, pt2, drawID, scaleFactor, options, false)); 212 if (ret != null) { 213 if (nth > 0 && ++nop != nth) 214 continue; 215 infolist[i] = ret; 216 sops += "\n" + (i + 1) + "\t" + ret[0] + "\t" + ret[2]; 217 opCount++; 218 } 219 } 220 info.put("operations", infolist); 221 info.put("symmetryInfo", (sops.length() == 0 ? "" : sops.substring(1))); 222 } 223 sgNote = (opCount == 0 ? "\n no symmetry operations" 224 : nth <= 0 && symOp <= 0 ? "\n" + opCount 225 + " symmetry operation" + (opCount == 1 ? ":\n" : "s:\n") : ""); 226 if (slist != null) 227 sgName = slist.substring(slist.indexOf(";") + 1); 228 if (saveModelInfo) 229 modelSet.setInfo(modelIndex, "spaceGroupInfo", info); 230 } else { 231 info = new Hashtable<String, Object>(); 232 } 233 info.put("spaceGroupName", sgName); 234 info.put("spaceGroupNote", sgNote == null ? "" : sgNote); 235 Object data; 236 if (isBio) { 237 data = sgName; 238 } else { 239 if (haveName && !haveRawName) 240 sym.setSpaceGroupName(sgName); 241 data = sym.getSpaceGroupInfoObj(sgName, cellInfo, isFull); 242 if (data == null || data.equals("?")) { 243 data = "?"; 244 info.put("spaceGroupNote", "could not identify space group from name: " + sgName 245 + "\nformat: show spacegroup \"2\" or \"P 2c\" " 246 + "or \"C m m m\" or \"x, y, z;-x ,-y, -z\""); 247 } 248 } 249 info.put("spaceGroupInfo", data); 250 return info; 251 } 252 253 //////////// private methods /////////// 254 255 /** Determine the type of this request. 256 * Note that label and xyz will be returned as T.xys and T.label 257 * 258 * @param id 259 * @return a code that identifies this request. 260 */ getType(String id)261 private static int getType(String id) { 262 int type; 263 if (id == null) 264 return T.list; 265 if (id.equalsIgnoreCase("matrix")) 266 return T.matrix4f; 267 if (id.equalsIgnoreCase("description")) 268 return T.label; 269 if (id.equalsIgnoreCase("axispoint")) 270 return T.point; 271 if (id.equalsIgnoreCase("time")) 272 return T.times; 273 if (id.equalsIgnoreCase("info")) 274 return T.array; 275 // center, draw, plane, axis, atom, translation, angle, array, list 276 type = T.getTokFromName(id); 277 if (type != 0) 278 return type; 279 type = getKeyType(id); 280 return (type < 0 ? type : T.all); 281 } 282 getKeyType(String id)283 private static int getKeyType(String id) { 284 if ("type".equals(id)) 285 id = "_type"; 286 for (int type = 0; type < keys.length; type++) 287 if (id.equalsIgnoreCase(keys[type])) 288 return -1 - type; 289 return 0; 290 } 291 nullReturn(int type)292 private static Object nullReturn(int type) { 293 switch(type) { 294 case T.draw: 295 return "draw ID sym_* delete"; 296 case T.full: 297 case T.label: 298 case T.id: 299 case T.xyz: 300 case T.matrix3f: 301 case T.origin: 302 return ""; 303 case T.atoms: 304 return new BS(); 305 default: 306 return null; 307 } 308 } 309 310 /** 311 * Return information about a symmetry operator by type: 312 * 313 * array, angle, axis, center, draw, full, info, label, matrix4f, point, time, 314 * plane, translation, unitcell, xyz, all, 315 * 316 * or a negative number (-length, -1]: 317 * 318 * { "xyz", etc. } 319 * 320 * where "all" is the info array itself, 321 * 322 * @param info 323 * @param type 324 * @return object specified 325 * 326 */ 327 getInfo(Object[] info, int type)328 private static Object getInfo(Object[] info, int type) { 329 if (info.length == 0) 330 return ""; 331 if (type < 0 && -type <= keys.length && -type <= info.length) 332 return info[-1 - type]; 333 switch (type) { 334 case T.all: 335 case T.info: 336 return info; // 'info' is internal only; if user selects "info", it is turned into "array" 337 case T.array: 338 Map<String, Object> lst = new Hashtable<String, Object>(); 339 for (int j = 0, n = info.length; j < n; j++) { 340 String key = (j == KEY_DRAW ? "draw" : j == KEY_POINT ? "axispoint" : keys[j]); 341 if (info[j] != null) 342 lst.put(key, info[j]); 343 } 344 return lst; 345 case T.full: 346 return info[0] + " \t" + info[2]; 347 case T.xyz: 348 return info[0]; 349 case T.origin: 350 return info[1]; 351 default: 352 case T.label: 353 return info[2]; 354 case T.draw: 355 return info[3] + "\nprint " + PT.esc(info[0] + " " + info[2]); 356 case T.fracxyz: 357 return info[4]; // fractional translation "fxyz" 358 case T.translation: 359 return info[5]; // cartesian translation 360 case T.center: 361 return info[6]; 362 case T.point: 363 return info[7]; 364 case T.axis: 365 return info[8]; 366 case T.angle: 367 return info[9]; 368 case T.matrix4f: 369 return info[10]; 370 case T.unitcell: 371 return info[11]; 372 case T.translate: 373 // centering 374 return info[12]; 375 case T.times: 376 return info[13]; 377 case T.plane: 378 return info[14]; 379 case T.type: 380 return info[15]; 381 case T.id: 382 return info[16]; 383 } 384 } 385 386 /** 387 * 388 * @param op 389 * @param uc 390 * @param pta00 391 * optional initial atom point 392 * @param ptTarget 393 * optional target atom point 394 * @param id 395 * @param scaleFactor 396 * scale for rotation vector only 397 * @param options 398 * 0 or T.offset 399 * @param haveTranslation TODO 400 * @return Object[] containing: 401 * 402 * [0] xyz (Jones-Faithful calculated from matrix) 403 * 404 * [1] xyzOriginal (Provided by calling method) 405 * 406 * [2] info ("C2 axis", for example) 407 * 408 * [3] draw commands 409 * 410 * [4] translation vector (fractional) 411 * 412 * [5] translation vector (Cartesian) 413 * 414 * [6] inversion point 415 * 416 * [7] axis point 417 * 418 * [8] axis vector (defines plane if angle = 0 419 * 420 * [9] angle of rotation 421 * 422 * [10] matrix representation 423 * 424 * [11] lattice translation 425 * 426 * [12] centering 427 * 428 * [13] time reversal 429 * 430 * [14] plane 431 * 432 * [15] _type 433 * 434 * [16] id 435 * 436 */ 437 createInfoArray(SymmetryOperation op, SymmetryInterface uc, P3 pta00, P3 ptTarget, String id, float scaleFactor, int options, boolean haveTranslation)438 private Object[] createInfoArray(SymmetryOperation op, SymmetryInterface uc, 439 P3 pta00, P3 ptTarget, String id, 440 float scaleFactor, int options, boolean haveTranslation) { 441 if (!op.isFinalized) 442 op.doFinalize(); 443 boolean isTimeReversed = (op.timeReversal == -1); 444 if (scaleFactor == 0) 445 scaleFactor = 1f; 446 V3 vtemp = new V3(); 447 P3 ptemp = new P3(); 448 P3 ptemp2 = new P3(); 449 P3 pta01 = new P3(); 450 P3 pta02 = new P3(); 451 V3 ftrans = new V3(); 452 V3 vtrans = new V3(); 453 P4 plane = null; 454 455 if (pta00 == null || Float.isNaN(pta00.x)) 456 pta00 = new P3(); 457 if (ptTarget != null) { 458 459 // Check to see that this is the correct operator 460 461 setFractional(uc, pta00, pta01, ptemp); 462 op.rotTrans(pta01); 463 uc.toCartesian(pta01, false); 464 uc.toUnitCell(pta01, ptemp); 465 pta02.setT(ptTarget); 466 uc.toUnitCell(pta02, ptemp); 467 if (pta01.distance(pta02) > 0.1f) 468 return null; 469 470 // Check to see if the two points only differ by 471 // a translation after transformation. 472 // If so, add that difference to the matrix transformation 473 474 setFractional(uc, pta00, pta01, null); 475 op.rotTrans(pta01); 476 setFractional(uc, ptTarget, pta02, null); 477 vtrans.sub2(pta02, pta01); 478 } 479 480 // get the frame vectors and points 481 482 pta01.set(1, 0, 0); 483 pta02.set(0, 1, 0); 484 P3 pta03 = P3.new3(0, 0, 1); 485 pta01.add(pta00); 486 pta02.add(pta00); 487 pta03.add(pta00); 488 489 // target point, rotated, inverted, and translated 490 491 P3 pt0 = rotTransCart(op, uc, pta00, vtrans); 492 P3 pt1 = rotTransCart(op, uc, pta01, vtrans); 493 P3 pt2 = rotTransCart(op, uc, pta02, vtrans); 494 P3 pt3 = rotTransCart(op, uc, pta03, vtrans); 495 496 V3 vt1 = V3.newVsub(pt1, pt0); 497 V3 vt2 = V3.newVsub(pt2, pt0); 498 V3 vt3 = V3.newVsub(pt3, pt0); 499 500 approx(vtrans); 501 502 // check for inversion 503 504 vtemp.cross(vt1, vt2); 505 boolean haveInversion = (vtemp.dot(vt3) < 0); 506 507 // The first trick is to check cross products to see if we still have a 508 // right-hand axis. 509 510 if (haveInversion) { 511 512 // undo inversion for quaternion analysis (requires proper rotations only) 513 514 pt1.sub2(pt0, vt1); 515 pt2.sub2(pt0, vt2); 516 pt3.sub2(pt0, vt3); 517 518 } 519 520 // The second trick is to use quaternions. Each of the three faces of the 521 // frame (xy, yz, and zx) 522 // is checked. The helix() function will return data about the local helical 523 // axis, and the 524 // symop(sym,{0 0 0}) function will return the overall translation. 525 526 T3[] info = Measure.computeHelicalAxis(pta00, pt0, 527 Quat.getQuaternionFrame(pt0, pt1, pt2) 528 .div(Quat.getQuaternionFrame(pta00, pta01, pta02))); 529 // new T3[] { pt_a_prime, n, r, P3.new3(theta, pitch, residuesPerTurn), pt_b_prime }; 530 P3 pa1 = (P3) info[0]; 531 V3 ax1 = (V3) info[1]; 532 int ang1 = (int) Math.abs(PT.approx(((P3) info[3]).x, 1)); 533 float pitch1 = SymmetryOperation.approxF(((P3) info[3]).y); 534 535 if (haveInversion) { 536 537 // redo inversion 538 539 pt1.add2(pt0, vt1); 540 pt2.add2(pt0, vt2); 541 pt3.add2(pt0, vt3); 542 543 } 544 545 V3 trans = V3.newVsub(pt0, pta00); 546 if (trans.length() < 0.1f) 547 trans = null; 548 549 // ////////// determination of type of operation from first principles 550 551 P3 ptinv = null; // inverted point for translucent frame 552 P3 ipt = null; // inversion center 553 P3 ptref = null; // reflection center 554 555 boolean isTranslation = (ang1 == 0); 556 boolean isRotation = !isTranslation; 557 boolean isInversionOnly = false; 558 boolean isMirrorPlane = false; 559 560 if (isRotation || haveInversion) { 561 // the translation will be taken care of other ways 562 trans = null; 563 } 564 565 // handle inversion 566 567 if (haveInversion && isTranslation) { 568 569 // simple inversion operation 570 571 ipt = P3.newP(pta00); 572 ipt.add(pt0); 573 ipt.scale(0.5f); 574 ptinv = pt0; 575 isInversionOnly = true; 576 577 } else if (haveInversion) { 578 579 /* 580 * 581 * We must convert simple rotations to rotation-inversions; 2-screws to 582 * planes and glide planes. 583 * 584 * The idea here is that there is a relationship between the axis for a 585 * simple or screw rotation of an inverted frame and one for a 586 * rotation-inversion. The relationship involves two adjacent equilateral 587 * triangles: 588 * 589 * 590 * o 591 * / \ 592 * / \ i' 593 * / \ 594 * / i \ 595 * A/_________\A' 596 * \ / 597 * \ j / 598 * \ / 599 * \ / 600 * \ / 601 * x 602 * 603 * Points i and j are at the centers of the triangles. Points A and A' are 604 * the frame centers; an operation at point i, j, x, or o is taking A to 605 * A'. Point i is 2/3 of the way from x to o. In addition, point j is half 606 * way between i and x. 607 * 608 * The question is this: Say you have an rotation/inversion taking A to 609 * A'. The relationships are: 610 * 611 * 6-fold screw x for inverted frame corresponds to 6-bar at i for actual 612 * frame 3-fold screw i for inverted frame corresponds to 3-bar at x for 613 * actual frame 614 * 615 * The proof of this follows. Consider point x. Point x can transform A to 616 * A' as a clockwise 6-fold screw axis. So, say we get that result for the 617 * inverted frame. What we need for the real frame is a 6-bar axis 618 * instead. Remember, though, that we inverted the frame at A to get this 619 * result. The real frame isn't inverted. The 6-bar must do that inversion 620 * AND also get the frame to point A' with the same (clockwise) rotation. 621 * The key is to see that there is another axis -- at point i -- that does 622 * the trick. 623 * 624 * Take a look at the angles and distances that arise when you project A 625 * through point i. The result is a frame at i'. Since the distance i-i' 626 * is the same as i-A (and thus i-A') and the angle i'-i-A' is 60 degrees, 627 * point i is also a 6-bar axis transforming A to A'. 628 * 629 * Note that both the 6-fold screw axis at x and the 6-bar axis at i are 630 * both clockwise. 631 * 632 * Similar analysis shows that the 3-fold screw i corresponds to the 3-bar 633 * axis at x. 634 * 635 * So in each case we just calculate the vector i-j or x-o and then factor 636 * appropriately. 637 * 638 * The 4-fold case is simpler -- just a parallelogram. 639 */ 640 641 V3 d = (pitch1 == 0 ? new V3() : ax1); 642 float f = 0; 643 switch (ang1) { 644 case 60: // 6_1 at x to 6-bar at i 645 f = 2f / 3f; 646 break; 647 case 120: // 3_1 at i to 3-bar at x 648 f = 2; 649 break; 650 case 90: // 4_1 to 4-bar at opposite corner 651 f = 1; 652 break; 653 case 180: // 2_1 to mirror plane 654 // C2 with inversion is a mirror plane -- but could have a glide 655 // component. 656 ptref = P3.newP(pta00); 657 ptref.add(d); 658 pa1.scaleAdd2(0.5f, d, pta00); 659 if (ptref.distance(pt0) > 0.1f) { 660 trans = V3.newVsub(pt0, ptref); 661 setFractional(uc, trans, ptemp, null); 662 ftrans.setT(ptemp); 663 } else { 664 trans = null; 665 } 666 isRotation = false; 667 haveInversion = false; 668 isMirrorPlane = true; 669 break; 670 default: 671 haveInversion = false; 672 break; 673 } 674 if (f != 0) { 675 vtemp.sub2(pta00, pa1); 676 vtemp.add(pt0); 677 vtemp.sub(pa1); 678 vtemp.sub(d); 679 vtemp.scale(f); 680 pa1.add(vtemp); 681 ipt = new P3(); 682 ipt.scaleAdd2(0.5f, d, pa1); 683 ptinv = new P3(); 684 ptinv.scaleAdd2(-2, ipt, pta00); 685 ptinv.scale(-1); 686 } 687 688 } else if (trans != null) { 689 690 // get rid of unnecessary translations added to keep most operations 691 // within cell 555 692 693 ptemp.setT(trans); 694 uc.toFractional(ptemp, false); 695 // if (SymmetryOperation.approxF(ptemp.x) == 1) 696 // ptemp.x = 0; 697 // if (SymmetryOperation.approxF(ptemp.y) == 1) 698 // ptemp.y = 0; 699 // if (SymmetryOperation.approxF(ptemp.z) == 1) 700 // ptemp.z = 0; 701 702 ftrans.setT(ptemp); 703 uc.toCartesian(ptemp, false); 704 trans.setT(ptemp); 705 } 706 707 // fix angle based on direction of axis 708 709 int ang = ang1; 710 approx0(ax1); 711 712 if (isRotation) { 713 714 P3 ptr = new P3(); 715 716 vtemp.setT(ax1); 717 718 // draw the lines associated with a rotation 719 720 int ang2 = ang1; 721 if (haveInversion) { 722 ptr.add2(pa1, vtemp); 723 ang2 = Math.round(Measure.computeTorsion(ptinv, pa1, ptr, pt0, true)); 724 } else if (pitch1 == 0) { 725 ptr.setT(pa1); 726 ptemp.scaleAdd2(1, ptr, vtemp); 727 ang2 = Math.round(Measure.computeTorsion(pta00, pa1, ptemp, pt0, true)); 728 } else { 729 ptemp.add2(pa1, vtemp); 730 ptr.scaleAdd2(0.5f, vtemp, pa1); 731 ang2 = Math.round(Measure.computeTorsion(pta00, pa1, ptemp, pt0, true)); 732 } 733 734 if (ang2 != 0) 735 ang1 = ang2; 736 737 if (!haveInversion && pitch1 == 0 && (ax1.z < 0 738 || ax1.z == 0 && (ax1.y < 0 || ax1.y == 0 && ax1.x < 0))) { 739 ax1.scale(-1); 740 ang1 = -ang1; 741 } 742 743 } 744 745 // time to get the description 746 747 String info1 = "identity"; 748 String type = info1; 749 750 if (isInversionOnly) { 751 ptemp.setT(ipt); 752 uc.toFractional(ptemp, false); 753 info1 = "Ci: " + strCoord(op, ptemp, op.isBio); 754 type = "inversion center"; 755 } else if (isRotation) { 756 if (haveInversion) { 757 type = info1 = (360 / ang) + "-bar axis"; 758 } else if (pitch1 != 0) { 759 type = info1 = (360 / ang) + "-fold screw axis"; 760 ptemp.setT(ax1); 761 uc.toFractional(ptemp, false); 762 info1 += "|translation: " + strCoord(op, ptemp, op.isBio); 763 } else { 764 type = info1 = "C" + (360 / ang) + " axis"; 765 } 766 } else if (trans != null) { 767 String s = " " + strCoord(op, ftrans, op.isBio); 768 if (isTranslation) { 769 type = info1 = "translation"; 770 info1 += ":" + s; 771 } else if (isMirrorPlane) { 772 float fx = Math.abs(SymmetryOperation.approxF(ftrans.x)); 773 float fy = Math.abs(SymmetryOperation.approxF(ftrans.y)); 774 float fz = Math.abs(SymmetryOperation.approxF(ftrans.z)); 775 s = " " + strCoord(op, ftrans, op.isBio); 776 // set ITA Table 2.1.2.1 777 if (fx != 0 && fy != 0 && fz != 0) { 778 if (fx == 1 / 4f && fy == 1 / 4f && fz == 1 / 4f) { 779 // diamond 780 info1 = "d-"; 781 } else if (fx == 1 / 2f && fy == 1 / 2f && fz == 1 / 2f) { 782 info1 = "n-"; 783 } else { 784 info1 = "g-"; 785 } 786 } else if (fx != 0 && fy != 0 || fy != 0 && fz != 0 787 || fz != 0 && fx != 0) { 788 // any two 789 if (fx == 1 / 4f && fy == 1 / 4f || fx == 1 / 4f && fz == 1 / 4f 790 || fy == 1 / 4f && fz == 1 / 4f) { 791 info1 = "d-"; 792 } else if (fx == 1 / 2f && fy == 1 / 2f 793 || fx == 1 / 2f && fz == 1 / 2f || fy == 1 / 2f && fz == 1 / 2f) { 794 // making sure here that this is truly a diagonal in the plane, not just 795 // a glide parallel to a face on a diagonal plane! Mois Aroyo 2018 796 if (fx == 0 && ax1.x == 0 || fy == 0 && ax1.y == 0 797 || fz == 0 && ax1.z == 0) { 798 info1 = "g-"; 799 } else { 800 info1 = "n-"; 801 } 802 } else { 803 info1 = "g-"; 804 } 805 } else if (fx != 0) 806 info1 = "a-"; 807 else if (fy != 0) 808 info1 = "b-"; 809 else 810 info1 = "c-"; 811 type = info1 = info1 + "glide plane"; 812 info1 += "|translation:" + s; 813 } 814 } else if (isMirrorPlane) { 815 type = info1 = "mirror plane"; 816 } 817 818 if (haveInversion && !isInversionOnly) { 819 ptemp.setT(ipt); 820 uc.toFractional(ptemp, false); 821 info1 += "|inversion center at " + strCoord(op, ptemp, op.isBio); 822 } 823 824 if (isTimeReversed) { 825 info1 += "|time-reversed"; 826 type += " (time-reversed)"; 827 } 828 829 String cmds = null; 830 String xyzNew = (op.isBio ? op.xyzOriginal 831 : SymmetryOperation.getXYZFromMatrix(op, false, false, false)); 832 833 // check for drawing 834 835 if (id != null) { 836 837 String opType = null; 838 String drawid = "\ndraw ID " + id + "_"; 839 840 // delete previous elements of this user-settable ID 841 842 SB draw1 = new SB(); 843 844 draw1.append(drawid).append("* delete"); 845 // .append( 846 // ("print " + PT.esc( 847 // id + " " + (op.index + 1) + " " + op.fixMagneticXYZ(op, op.xyzOriginal, false) + "|" 848 // + op.fixMagneticXYZ(op, xyzNew, true) + "|" + info1).replace( 849 // '\n', ' '))).append("\n") 850 851 // draw the initial frame 852 853 drawLine(draw1, drawid + "frame1X", 0.15f, pta00, pta01, "red"); 854 drawLine(draw1, drawid + "frame1Y", 0.15f, pta00, pta02, "green"); 855 drawLine(draw1, drawid + "frame1Z", 0.15f, pta00, pta03, "blue"); 856 857 String color; 858 859 if (isRotation) { 860 861 P3 ptr = new P3(); 862 863 color = "red"; 864 865 ang = ang1; 866 float scale = 1f; 867 vtemp.setT(ax1); 868 869 // draw the lines associated with a rotation 870 871 if (haveInversion) { 872 opType = drawid + "rotinv"; 873 ptr.add2(pa1, vtemp); 874 if (pitch1 == 0) { 875 ptr.setT(ipt); 876 vtemp.scale(3 * scaleFactor); 877 ptemp.scaleAdd2(-1, vtemp, pa1); 878 drawVector(draw1, drawid, "rotVector2", "", pa1, ptemp, "red"); 879 } 880 scale = pt0.distance(ptr); 881 draw1.append(drawid).append("rotLine1 ").append(Escape.eP(ptr)) 882 .append(Escape.eP(ptinv)).append(" color red"); 883 draw1.append(drawid).append("rotLine2 ").append(Escape.eP(ptr)) 884 .append(Escape.eP(pt0)).append(" color red"); 885 } else if (pitch1 == 0) { 886 opType = drawid + "rot"; 887 boolean isSpecial = (pta00.distance(pt0) < 0.2f); 888 if (!isSpecial) { 889 draw1.append(drawid).append("rotLine1 ").append(Escape.eP(pta00)) 890 .append(Escape.eP(pa1)).append(" color red"); 891 draw1.append(drawid).append("rotLine2 ").append(Escape.eP(pt0)) 892 .append(Escape.eP(pa1)).append(" color red"); 893 } 894 vtemp.scale(3 * scaleFactor); 895 ptemp.scaleAdd2(-1, vtemp, pa1); 896 drawVector(draw1, drawid, "rotVector2", "", pa1, ptemp, 897 isTimeReversed ? "gray" : "red"); 898 ptr.setT(pa1); 899 if (pitch1 == 0 && pta00.distance(pt0) < 0.2) 900 ptr.scaleAdd2(0.5f, vtemp, ptr); 901 } else { 902 opType = drawid + "screw"; 903 color = "orange"; 904 draw1.append(drawid).append("rotLine1 ").append(Escape.eP(pta00)) 905 .append(Escape.eP(pa1)).append(" color red"); 906 ptemp.add2(pa1, vtemp); 907 draw1.append(drawid).append("rotLine2 ").append(Escape.eP(pt0)) 908 .append(Escape.eP(ptemp)).append(" color red"); 909 ptr.scaleAdd2(0.5f, vtemp, pa1); 910 } 911 912 // draw arc arrow 913 914 ptemp.add2(ptr, vtemp); 915 if (haveInversion && pitch1 != 0) { 916 draw1.append(drawid).append("rotRotLine1").append(Escape.eP(ptr)) 917 .append(Escape.eP(ptinv)).append(" color red"); 918 draw1.append(drawid).append("rotRotLine2").append(Escape.eP(ptr)) 919 .append(Escape.eP(pt0)).append(" color red"); 920 } 921 draw1.append(drawid) 922 .append( 923 "rotRotArrow arrow width 0.1 scale " + PT.escF(scale) + " arc ") 924 .append(Escape.eP(ptr)).append(Escape.eP(ptemp)); 925 ptemp.setT(haveInversion ? ptinv : pta00); 926 if (ptemp.distance(pt0) < 0.1f) 927 ptemp.set((float) Math.random(), (float) Math.random(), 928 (float) Math.random()); 929 draw1.append(Escape.eP(ptemp)); 930 ptemp.set(0, ang, 0); 931 draw1.append(Escape.eP(ptemp)).append(" color red"); 932 933 // draw the main vector 934 935 drawVector(draw1, drawid, "rotVector1", "vector", pa1, vtemp, 936 isTimeReversed ? "Gray" : color); 937 938 } else if (isMirrorPlane) { 939 940 // lavender arrow across plane from pt00 to pt0 941 942 ptemp.sub2(ptref, pta00); 943 if (pta00.distance(ptref) > 0.2) 944 drawVector(draw1, drawid, "planeVector", "vector", pta00, ptemp, 945 isTimeReversed ? "Gray" : "cyan"); 946 947 // faint inverted frame if mirror trans is not null 948 949 opType = drawid + "plane"; 950 if (trans != null) { 951 opType = drawid + "glide"; 952 drawFrameLine("X", ptref, vt1, 0.15f, ptemp, draw1, opType, "red"); 953 drawFrameLine("Y", ptref, vt2, 0.15f, ptemp, draw1, opType, "green"); 954 drawFrameLine("Z", ptref, vt3, 0.15f, ptemp, draw1, opType, "blue"); 955 } 956 957 color = (trans == null ? "green" : "blue"); 958 959 // ok, now HERE's a good trick. We use the Marching Cubes 960 // algorithm to find the intersection points of a plane and the unit 961 // cell. 962 // We expand the unit cell by 5% in all directions just so we are 963 // guaranteed to get cutoffs. 964 965 vtemp.setT(ax1); 966 vtemp.normalize(); 967 // ax + by + cz + d = 0 968 // so if a point is in the plane, then N dot X = -d 969 float w = -vtemp.x * pa1.x - vtemp.y * pa1.y - vtemp.z * pa1.z; 970 plane = P4.new4(vtemp.x, vtemp.y, vtemp.z, w); 971 float margin = (Math.abs(w) < 0.01f && vtemp.x * vtemp.y > 0.4 ? 1.30f 972 : 1.05f); 973 // returns triangles and lines 974 Lst<Object> v = modelSet.vwr.getTriangulator().intersectPlane(plane, 975 uc.getCanonicalCopy(margin, true), 3); 976 if (v != null) 977 for (int i = v.size(); --i >= 0;) { 978 P3[] pts = (P3[]) v.get(i); 979 draw1.append(drawid).append("planep").appendI(i).append(" ") 980 .append(Escape.eP(pts[0])).append(Escape.eP(pts[1])); 981 if (pts.length == 3) 982 draw1.append(Escape.eP(pts[2])); 983 draw1.append(" color translucent ").append(color); 984 } 985 986 // and JUST in case that does not work, at least draw a circle 987 988 if (v == null || v.size() == 0) { 989 ptemp.add2(pa1, ax1); 990 draw1.append(drawid).append("planeCircle scale 2.0 circle ") 991 .append(Escape.eP(pa1)).append(Escape.eP(ptemp)) 992 .append(" color translucent ").append(color).append(" mesh fill"); 993 } 994 } 995 996 if (haveInversion) { 997 opType = drawid + "inv"; 998 draw1.append(drawid).append("invPoint diameter 0.4 ") 999 .append(Escape.eP(ipt)); 1000 ptemp.sub2(ptinv, pta00); 1001 drawVector(draw1, drawid, "invArrow", "vector", pta00, ptemp, 1002 isTimeReversed ? "gray" : "cyan"); 1003 if (!isInversionOnly && options != T.offset) { 1004 // n-bar: draw a faint frame showing the inversion 1005 drawFrameLine("X", ptinv, vt1, 0.15f, ptemp, draw1, opType, "red"); 1006 drawFrameLine("Y", ptinv, vt2, 0.15f, ptemp, draw1, opType, "green"); 1007 drawFrameLine("Z", ptinv, vt3, 0.15f, ptemp, draw1, opType, "blue"); 1008 } 1009 } 1010 1011 // and display translation if still not {0 0 0} 1012 1013 if (trans != null) { 1014 if (ptref == null) 1015 ptref = P3.newP(pta00); 1016 drawVector(draw1, drawid, "transVector", "vector", ptref, trans, 1017 isTimeReversed && !haveInversion && !isMirrorPlane && !isRotation 1018 ? "darkGray" 1019 : "gold"); 1020 } 1021 1022 // draw the final frame just a bit fatter and shorter, in case they 1023 // overlap 1024 1025 ptemp2.setT(pt0); 1026 ptemp.sub2(pt1, pt0); 1027 ptemp.scaleAdd2(0.9f, ptemp, ptemp2); 1028 drawLine(draw1, drawid + "frame2X", 0.2f, ptemp2, ptemp, "red"); 1029 ptemp.sub2(pt2, pt0); 1030 ptemp.scaleAdd2(0.9f, ptemp, ptemp2); 1031 drawLine(draw1, drawid + "frame2Y", 0.2f, ptemp2, ptemp, "green"); 1032 ptemp.sub2(pt3, pt0); 1033 ptemp.scaleAdd2(0.9f, ptemp, ptemp2); 1034 drawLine(draw1, drawid + "frame2Z", 0.2f, ptemp2, ptemp, "purple"); 1035 1036 // color the targeted atoms opaque and add another frame if necessary 1037 1038 draw1.append("\nsym_point = " + Escape.eP(pta00)); 1039 draw1.append("\nvar p0 = " + Escape.eP(ptemp2)); 1040 draw1.append( 1041 "\nvar set2 = within(0.2,p0);if(!set2){set2 = within(0.2,p0.uxyz.xyz)}"); 1042 if (pta00 instanceof Atom) 1043 draw1.append("\n set2 &= {_" + ((Atom) pta00).getElementSymbol() + "}"); 1044 draw1.append("\nsym_target = set2;if (set2) {"); 1045 // if (haveCentering) 1046 // draw1.append(drawid).append( 1047 // "cellOffsetVector arrow @p0 @set2 color grey"); 1048 if (options != T.offset && ptTarget == null && !haveTranslation) { 1049 draw1.append(drawid) 1050 .append("offsetFrameX diameter 0.20 @{set2.xyz} @{set2.xyz + ") 1051 .append(Escape.eP(vt1)).append("*0.9} color red"); 1052 draw1.append(drawid) 1053 .append("offsetFrameY diameter 0.20 @{set2.xyz} @{set2.xyz + ") 1054 .append(Escape.eP(vt2)).append("*0.9} color green"); 1055 draw1.append(drawid) 1056 .append("offsetFrameZ diameter 0.20 @{set2.xyz} @{set2.xyz + ") 1057 .append(Escape.eP(vt3)).append("*0.9} color purple"); 1058 } 1059 draw1.append("\n}\n"); 1060 cmds = draw1.toString(); 1061 if (Logger.debugging) 1062 Logger.info(cmds); 1063 draw1 = null; 1064 drawid = null; 1065 } 1066 if (trans == null) 1067 ftrans = null; 1068 if (isRotation) { 1069 if (haveInversion) { 1070 } else if (pitch1 == 0) { 1071 } else { 1072 // screw 1073 trans = V3.newV(ax1); 1074 ptemp.setT(trans); 1075 uc.toFractional(ptemp, false); 1076 ftrans = V3.newV(ptemp); 1077 } 1078 } 1079 if (isMirrorPlane) { 1080 ang1 = 0; 1081 } 1082 if (haveInversion) { 1083 if (isInversionOnly) { 1084 pa1 = null; 1085 ax1 = null; 1086 trans = null; 1087 ftrans = null; 1088 } 1089 } else if (isTranslation) { 1090 pa1 = null; 1091 ax1 = null; 1092 } 1093 1094 // and display translation if still not {0 0 0} 1095 if (ax1 != null) 1096 ax1.normalize(); 1097 M4 m2 = null; 1098 m2 = M4.newM4(op); 1099 if (vtrans.length() != 0) { 1100 m2.m03 += vtrans.x; 1101 m2.m13 += vtrans.y; 1102 m2.m23 += vtrans.z; 1103 } 1104 // TODO: here we need magnetic 1105 xyzNew = (op.isBio ? m2.toString() 1106 : op.modDim > 0 ? op.xyzOriginal 1107 : SymmetryOperation.getXYZFromMatrix(m2, false, false, false)); 1108 if (op.timeReversal != 0) 1109 xyzNew = op.fixMagneticXYZ(m2, xyzNew, true); 1110 T3 f0 = approx0(ftrans); 1111 T3 cift = null; 1112 int cifi = (op.number < 0 ? 0 : op.number); 1113 if (!xyzNew.equals(op.xyzOriginal)) { 1114 if (op.number > 0) { 1115 M4 orig = SymmetryOperation.getMatrixFromXYZ(op.xyzOriginal); 1116 orig.sub(m2); 1117 cift = new P3(); 1118 orig.getTranslation(cift); 1119 } 1120 } 1121 String cif2 = cifi + (cift == null ? " [0 0 0]" 1122 : " [" + (int) -cift.x + " " + (int) -cift.y + " " + (int) -cift.z + "]"); 1123 return new Object[] { xyzNew, op.xyzOriginal, info1, cmds, f0, 1124 approx0(trans), approx0(ipt), approx0(pa1), 1125 (plane == null ? approx0(ax1) : null), 1126 (ang1 != 0 ? Integer.valueOf(ang1) : null), m2, 1127 (vtrans.lengthSquared() > 0 ? vtrans : null), op.getCentering(), 1128 Integer.valueOf(op.timeReversal), plane, type, 1129 Integer.valueOf(op.number), cif2, op.xyzCanonical }; 1130 } 1131 drawLine(SB s, String id, float diameter, P3 pt0, P3 pt1, String color)1132 private static void drawLine(SB s, String id, float diameter, P3 pt0, P3 pt1, 1133 String color) { 1134 s.append(id).append(" diameter ").appendF(diameter).append(Escape.eP(pt0)) 1135 .append(Escape.eP(pt1)).append(" color ").append(color); 1136 } 1137 drawFrameLine(String xyz, P3 pt, V3 v, float width, P3 ptemp, SB draw1, String key, String color)1138 private static void drawFrameLine(String xyz, P3 pt, V3 v, float width, 1139 P3 ptemp, SB draw1, String key, String color) { 1140 ptemp.setT(pt); 1141 ptemp.add(v); 1142 drawLine(draw1, key + "Pt" + xyz, width, pt, ptemp, "translucent " + color); 1143 } 1144 drawVector(SB draw1, String drawid, String label, String type, T3 pt1, T3 v, String color)1145 private static void drawVector(SB draw1, String drawid, String label, 1146 String type, T3 pt1, T3 v, String color) { 1147 draw1.append(drawid).append(label).append(" diameter 0.1 ").append(type) 1148 .append(Escape.eP(pt1)).append(Escape.eP(v)).append(" color ") 1149 .append(color); 1150 } 1151 1152 /** 1153 * Set pt01 to pt00, possibly adding offset into unit cell 1154 * 1155 * @param uc 1156 * @param pt00 1157 * @param pt01 1158 * @param offset 1159 */ setFractional(SymmetryInterface uc, T3 pt00, P3 pt01, P3 offset)1160 private static void setFractional(SymmetryInterface uc, T3 pt00, P3 pt01, 1161 P3 offset) { 1162 pt01.setT(pt00); 1163 if (offset != null) 1164 uc.toUnitCell(pt01, offset); 1165 uc.toFractional(pt01, false); 1166 } 1167 rotTransCart(SymmetryOperation op, SymmetryInterface uc, P3 pt00, V3 vtrans)1168 private static P3 rotTransCart(SymmetryOperation op, SymmetryInterface uc, 1169 P3 pt00, V3 vtrans) { 1170 P3 p0 = P3.newP(pt00); 1171 uc.toFractional(p0, false); 1172 op.rotTrans(p0); 1173 p0.add(vtrans); 1174 uc.toCartesian(p0, false); 1175 return p0; 1176 } 1177 strCoord(SymmetryOperation op, T3 p, boolean isBio)1178 private static String strCoord(SymmetryOperation op, T3 p, boolean isBio) { 1179 approx0(p); 1180 return (isBio ? p.x + " " + p.y + " " + p.z : op.fcoord2(p)); 1181 } 1182 approx0(T3 pt)1183 private static T3 approx0(T3 pt) { 1184 if (pt != null) { 1185 if (Math.abs(pt.x) < 0.0001f) 1186 pt.x = 0; 1187 if (Math.abs(pt.y) < 0.0001f) 1188 pt.y = 0; 1189 if (Math.abs(pt.z) < 0.0001f) 1190 pt.z = 0; 1191 } 1192 return pt; 1193 } 1194 approx(T3 pt)1195 private static T3 approx(T3 pt) { 1196 if (pt != null) { 1197 pt.x = SymmetryOperation.approxF(pt.x); 1198 pt.y = SymmetryOperation.approxF(pt.y); 1199 pt.z = SymmetryOperation.approxF(pt.z); 1200 } 1201 return pt; 1202 } 1203 1204 /** 1205 * multipurpose function handling a variety of tasks, including: 1206 * 1207 * processing of "lattice", "list", "atom", "point", and some "draw" output 1208 * types 1209 * 1210 * finding the operator in the given space group 1211 * 1212 * creating a temporary space group for an xyz operator 1213 * 1214 * 1215 * @param sym 1216 * @param iModel 1217 * @param iatom 1218 * @param uc 1219 * @param xyz 1220 * @param op 1221 * @param translation [i j k] to be added to operator 1222 * @param pt 1223 * @param pt2 second point or offset 1224 * @param id 1225 * @param type 1226 * @param scaleFactor 1227 * @param nth 1228 * @param options 0 or T.offset 1229 * @return a string or an Object[] containing information 1230 */ getSymmetryInfo(Symmetry sym, int iModel, int iatom, Symmetry uc, String xyz, int op, P3 translation, P3 pt, P3 pt2, String id, int type, float scaleFactor, int nth, int options)1231 private Object getSymmetryInfo(Symmetry sym, int iModel, int iatom, 1232 Symmetry uc, String xyz, int op, P3 translation, 1233 P3 pt, P3 pt2, String id, 1234 int type, float scaleFactor, int nth, int options) { 1235 int returnType = 0; 1236 Object nullRet = nullReturn(type); 1237 switch (type) { 1238 case T.lattice: 1239 return uc.getLatticeType(); 1240 case T.array: 1241 returnType = getType(id); 1242 switch (returnType) { 1243 case T.atoms: 1244 case T.draw: 1245 case T.full: 1246 case T.list: 1247 case T.point: 1248 type = returnType; 1249 break; 1250 default: 1251 returnType = getKeyType(id); 1252 break; 1253 } 1254 break; 1255 } 1256 int iop = op; 1257 P3 offset = (options == T.offset && (type == T.atoms || type == T.point)? pt2 : null); 1258 if (offset != null) 1259 pt2 = null; 1260 Object[] info = null; 1261 String xyzOriginal = null; 1262 if (pt2 == null) { 1263 if (xyz == null) { 1264 SymmetryOperation[] ops = (SymmetryOperation[]) uc 1265 .getSymmetryOperations(); 1266 if (ops == null || op == 0 || Math.abs(op) > ops.length) 1267 return nullRet; 1268 iop = Math.abs(op) - 1; 1269 xyz = (translation == null ? ops[iop].xyz : ops[iop].getxyzTrans(translation)); 1270 xyzOriginal = ops[iop].xyzOriginal; 1271 } else { 1272 iop = op = 0; 1273 } 1274 SymmetryInterface symTemp = new Symmetry(); 1275 symTemp.setSpaceGroup(false); 1276 boolean isBio = uc.isBio(); 1277 int i = (isBio ? symTemp.addBioMoleculeOperation( 1278 uc.spaceGroup.finalOperations[iop], op < 0) : symTemp 1279 .addSpaceGroupOperation((op < 0 ? "!" : "=") + xyz, Math.abs(op))); 1280 1281 if (i < 0) 1282 return nullRet; 1283 SymmetryOperation opTemp = (SymmetryOperation) symTemp 1284 .getSpaceGroupOperation(i); 1285 if (xyzOriginal != null) 1286 opTemp.xyzOriginal = xyzOriginal; 1287 opTemp.number = op; 1288 if (!isBio) 1289 opTemp.getCentering(); 1290 if (pt == null && iatom >= 0) 1291 pt = modelSet.at[iatom]; 1292 if (type == T.point || type == T.atoms) { 1293 if (isBio) 1294 return nullRet; 1295 symTemp.setUnitCell(uc); 1296 pt = P3.newP(pt); 1297 uc.toFractional(pt, false); 1298 if (Float.isNaN(pt.x)) 1299 return nullRet; 1300 P3 sympt = new P3(); 1301 symTemp.newSpaceGroupPoint(i, pt, sympt, 0, 0, 0, null); 1302 1303 if (options == T.offset) { 1304 uc.unitize(sympt);//uc.unitize01(sympt); 1305 if (options == T.offset) 1306 sympt.add(offset); 1307 } 1308 symTemp.toCartesian(sympt, false); 1309 return (type == T.atoms ? getAtom(uc, iModel, iatom, sympt) : sympt); 1310 } 1311 info = createInfoArray(opTemp, uc, pt, null, (id == null ? "sym" : id), 1312 scaleFactor, options, (translation != null)); 1313 if (type == T.array && id != null) { 1314 returnType = getKeyType(id); 1315 } 1316 } else { 1317 // pt1, pt2 1318 String stype = "info"; 1319 boolean asString = false; 1320 switch (type) { 1321 case T.array: // new Jmol 14.29.45 1322 returnType = getKeyType(id); 1323 id = stype = null; 1324 asString = false; 1325 if (nth == 0) 1326 nth = -1; 1327 break; 1328 case T.list: 1329 id = stype = null; 1330 asString = true; 1331 if (nth == 0) 1332 nth = -1; 1333 break; 1334 case T.draw: 1335 if (id == null) 1336 id = "sym"; 1337 stype = "all"; 1338 asString = true; 1339 break; 1340 case T.atoms: 1341 id = stype = null; 1342 //$FALL-THROUGH$ 1343 default: 1344 if (nth == 0) 1345 nth = 1; 1346 } 1347 Object ret1 = getSymopInfoForPoints(sym, iModel, op, translation, pt, pt2, id, 1348 stype, scaleFactor, nth, asString, options); 1349 if (asString) 1350 return ret1; 1351 if (ret1 instanceof String) 1352 return nullRet; // two atoms are not connected, no such oper 1353 info = (Object[]) ret1; 1354 if (type == T.atoms) { 1355 if (!(pt instanceof Atom) && !(pt2 instanceof Atom)) 1356 iatom = -1; 1357 return (info == null ? nullRet : getAtom(uc, iModel, iatom, 1358 (T3) info[7])); 1359 } 1360 } 1361 if (info == null) 1362 return nullRet; 1363 boolean isList = (info.length > 0 && info[0] instanceof Object[]); 1364 if (nth < 0 && op <= 0 && (type == T.array || isList)) { 1365 if (type == T.array && info.length > 0 && !(info[0] instanceof Object[])) 1366 info = new Object[] { info }; 1367 Lst<Object> lst = new Lst<Object>(); 1368 for (int i = 0; i < info.length; i++) 1369 lst.addLast(getInfo((Object[])info[i], returnType < 0 ? returnType : type)); 1370 return lst; 1371 } else if (returnType < 0 && (nth >= 0 || op > 0)) { 1372 type = returnType; 1373 } 1374 if (nth > 0 && isList) 1375 info = (Object[]) ((Object[]) info)[0]; 1376 return getInfo(info, type); 1377 } 1378 getAtom(Symmetry uc, int iModel, int iAtom, T3 sympt)1379 private BS getAtom(Symmetry uc, int iModel, int iAtom, T3 sympt) { 1380 BS bsElement = null; 1381 if (iAtom >= 0) 1382 modelSet.getAtomBitsMDa(T.elemno, 1383 Integer.valueOf(modelSet.at[iAtom].getElementNumber()), 1384 bsElement = new BS()); 1385 BS bsResult = new BS(); 1386 modelSet.getAtomsWithin(0.02f, sympt, bsResult, iModel); 1387 if (bsElement != null) 1388 bsResult.and(bsElement); 1389 if (bsResult.isEmpty()) { 1390 sympt = P3.newP(sympt); 1391 uc.toUnitCell(sympt, null); 1392 uc.toCartesian(sympt, false); 1393 modelSet.getAtomsWithin(0.02f, sympt, bsResult, iModel); 1394 if (bsElement != null) 1395 bsResult.and(bsElement); 1396 } 1397 return bsResult; 1398 } 1399 1400 /** 1401 * get information about a symmetry operation relating two specific points or atoms 1402 * 1403 * @param sym 1404 * @param modelIndex 1405 * @param symOp 1406 * @param translation TODO 1407 * @param pt1 1408 * @param pt2 1409 * @param drawID 1410 * @param stype 1411 * @param scaleFactor 1412 * @param nth 1413 * @param asString 1414 * @param options 0 or T.offset 1415 * @return Object[] or String or Object[Object[]] (nth = 0, "array") 1416 * 1417 */ getSymopInfoForPoints(Symmetry sym, int modelIndex, int symOp, P3 translation, P3 pt1, P3 pt2, String drawID, String stype, float scaleFactor, int nth, boolean asString, int options)1418 private Object getSymopInfoForPoints(Symmetry sym, int modelIndex, int symOp, 1419 P3 translation, P3 pt1, P3 pt2, 1420 String drawID, String stype, float scaleFactor, 1421 int nth, boolean asString, int options) { 1422 Object ret = (asString ? "" : null); 1423 Map<String, Object> sginfo = getSpaceGroupInfo(sym, modelIndex, null, 1424 symOp, pt1, pt2, drawID, scaleFactor, nth, false, true, options, null); 1425 if (sginfo == null) 1426 return ret; 1427 Object[][] infolist = (Object[][]) sginfo.get("operations"); 1428 // at this point, if we have two points, we have a full list of operations, but 1429 // some are null. 1430 if (infolist == null) 1431 return ret; 1432 SB sb = (asString ? new SB() : null); 1433 symOp--; 1434 boolean isAll = (!asString && symOp < 0); 1435 String strOperations = (String) sginfo.get("symmetryInfo"); 1436 boolean labelOnly = "label".equals(stype); 1437 int n = 0; 1438 for (int i = 0; i < infolist.length; i++) { 1439 if (infolist[i] == null || symOp >= 0 && symOp != i) 1440 continue; 1441 if (!asString) { 1442 if (!isAll) 1443 return infolist[i]; 1444 infolist[n++] = infolist[i]; 1445 continue; 1446 } 1447 if (drawID != null) 1448 return ((String) infolist[i][3]) + "\nprint " + PT.esc(strOperations); 1449 if (sb.length() > 0) 1450 sb.appendC('\n'); 1451 if (!labelOnly) { 1452 if (symOp < 0) 1453 sb.appendI(i + 1).appendC('\t'); 1454 sb.append((String) infolist[i][0]).appendC('\t'); //xyz 1455 } 1456 sb.append((String) infolist[i][2]); //desc 1457 } 1458 if (!asString) { 1459 Object[] a = new Object[n]; 1460 for (int i = 0; i < n; i++) 1461 a[i] = infolist[i]; 1462 return a; 1463 } 1464 if (sb.length() == 0) 1465 return (drawID != null ? "draw " + drawID + "* delete" : ret); 1466 return sb.toString(); 1467 } 1468 } 1469