1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.potential.commands.test;
39
40 import ffx.potential.MolecularAssembly;
41 import ffx.potential.Utilities;
42 import ffx.potential.bonded.Angle;
43 import ffx.potential.bonded.Atom;
44 import ffx.potential.bonded.Bond;
45 import ffx.potential.bonded.MSGroup;
46 import ffx.potential.bonded.MSNode;
47 import ffx.potential.bonded.Molecule;
48 import ffx.potential.bonded.Polymer;
49 import ffx.potential.bonded.Residue;
50 import ffx.potential.cli.PotentialCommand;
51 import ffx.potential.parameters.ForceField;
52 import ffx.potential.parsers.CIFFilter;
53 import ffx.potential.parsers.PDBFilter;
54 import ffx.potential.parsers.SystemFilter;
55 import ffx.potential.parsers.XYZFilter;
56 import ffx.utilities.FFXBinding;
57 import org.apache.commons.configuration2.CompositeConfiguration;
58 import org.openscience.cdk.AtomContainer;
59 import org.openscience.cdk.config.AtomTypeFactory;
60 import org.openscience.cdk.graph.rebond.RebondTool;
61 import org.openscience.cdk.interfaces.IAtom;
62 import org.openscience.cdk.interfaces.IBond;
63 import org.openscience.cdk.isomorphism.AtomMatcher;
64 import org.openscience.cdk.isomorphism.BondMatcher;
65 import org.openscience.cdk.isomorphism.Pattern;
66 import org.openscience.cdk.isomorphism.VentoFoggia;
67 import picocli.CommandLine.Command;
68 import picocli.CommandLine.Option;
69 import picocli.CommandLine.Parameters;
70
71 import javax.vecmath.Point3d;
72 import java.io.BufferedWriter;
73 import java.io.File;
74 import java.io.FileWriter;
75 import java.util.ArrayList;
76 import java.util.Arrays;
77 import java.util.Comparator;
78 import java.util.List;
79 import java.util.logging.Level;
80
81 import static ffx.numerics.math.DoubleMath.dihedralAngle;
82 import static ffx.potential.bonded.BondedUtils.intxyz;
83 import static ffx.potential.utils.Superpose.applyRotation;
84 import static ffx.potential.utils.Superpose.applyTranslation;
85 import static ffx.potential.utils.Superpose.calculateRotation;
86 import static ffx.potential.utils.Superpose.calculateTranslation;
87 import static ffx.potential.utils.Superpose.rmsd;
88 import static ffx.potential.utils.Superpose.superpose;
89 import static java.lang.Double.MAX_VALUE;
90 import static java.lang.String.format;
91 import static org.apache.commons.io.FilenameUtils.getFullPath;
92 import static org.apache.commons.io.FilenameUtils.getName;
93 import static org.apache.commons.io.FilenameUtils.removeExtension;
94 import static org.apache.commons.math3.util.FastMath.sqrt;
95 import static org.openscience.cdk.tools.periodictable.PeriodicTable.getSymbol;
96
97
98
99
100
101 @Command(description = " Reorder the atoms of an XYZ file based on atomic weight.",
102 name = "test.ReorderAtoms")
103 public class ReorderAtoms extends PotentialCommand {
104
105
106
107
108 @Option(names = {"-a", "--atomType"}, paramLabel = "false", defaultValue = "false",
109 description = "Change atom types of first file to match second.")
110 private static boolean atomType = false;
111
112
113
114
115 @Option(names = {"--sw", "--swap"}, paramLabel = "false", defaultValue = "false",
116 description = "Switch positions of equivalent atoms (untested for >2 atoms).")
117 private static boolean swap = false;
118
119
120
121
122 @Option(names = {"-m", "--match"}, paramLabel = "false", defaultValue = "false",
123 description = "Match atom indices based on structure.")
124 private static boolean match = false;
125
126
127
128
129 @Option(names = {"-r", "--reorder"}, paramLabel = "false", defaultValue = "false",
130 description = "Reorder atoms based on environment.")
131 private static boolean reorder = false;
132
133
134
135
136 @Parameters(arity = "1..*", paramLabel = "files",
137 description = "Atomic coordinate files to reorder in XYZ/ARC format. If 2, change file 1 to match file 2.")
138 private List<String> filenames;
139
140
141
142
143 private double bondTolerance = 0.2;
144
145
146
147
148 public ReorderAtoms() {
149 super();
150 }
151
152
153
154
155 public ReorderAtoms(FFXBinding binding) {
156 super(binding);
157 }
158
159
160
161
162 public ReorderAtoms(String[] args) {
163 super(args);
164 }
165
166 @Override
167 public ReorderAtoms run() {
168 logger.warning(" The ReorderAtoms command is under development.\n " +
169 " Please verify the output.");
170
171 System.setProperty("vdwterm", "false");
172 if (!init()) {
173 return null;
174 }
175
176 if (reorder || swap) {
177 if (match || atomType) {
178 logger.warning(" Currently each flag (i.e., atomType, match, or reorder must be run individually. ");
179 }
180 for (String filename : filenames) {
181 potentialFunctions.openAll(filename);
182 SystemFilter systemFilter = potentialFunctions.getFilter();
183 int numModels = systemFilter.countNumModels();
184
185 File saveLocation = new File(filename);
186 File saveFile = potentialFunctions.versionFile(saveLocation);
187 for (int i = 0; i < numModels; i++) {
188 MolecularAssembly assembly = systemFilter.getActiveMolecularSystem();
189 Atom[] atoms = assembly.getAtomArray();
190 List<Atom> atomList0 = assembly.getAtomList();
191 int nAtoms = atoms.length;
192 MolecularAssembly newAssembly = new MolecularAssembly(assembly.getName());
193 List<Bond> bondList = assembly.getBondList();
194
195
196 int atomListSize = atomList0.size();
197 for (int j = 0; j < atomListSize; j++) {
198 for (int k = 0; k < atomListSize; k++) {
199 if (j != k) {
200 if (atomList0.get(j).getIndex() != atomList0.get(k).getIndex()) {
201 AtomComparator ac = new AtomComparator();
202 if (logger.isLoggable(Level.FINER)) {
203 logger.finer(" Compare atom: " + atomList0.get(j) + " and " + atomList0.get(k));
204 }
205 int value = ac.compare(atomList0.get(j), atomList0.get(k));
206 if (value == 0) {
207 logger.info(" This " + atomList0.get(j) + " and " + atomList0.get(k) +
208 " are \"equivalent\".");
209 if (swap) {
210 Atom temp = atomList0.get(j);
211 atomList0.set(j, atomList0.get(k));
212 atomList0.set(k, temp);
213 }
214 } else if (value < 0 && reorder) {
215 if (logger.isLoggable(Level.FINER)) {
216 logger.finer(" Swapped " + atomList0.get(j) + " and " + atomList0.get(k) + ".");
217 }
218 Atom temp = atomList0.get(j);
219 atomList0.set(j, atomList0.get(k));
220 atomList0.set(k, temp);
221 }
222 }
223 }
224 }
225 }
226
227 int[] originalOrder = new int[nAtoms];
228 for (int j = 0; j < nAtoms; j++) {
229 originalOrder[j] = atoms[j].getIndex();
230 if (logger.isLoggable(Level.FINE)) {
231 logger.fine(format(" Original index for atom %3d: %3d New index: %3d", j, originalOrder[j],
232 atomList0.get(j).getIndex()));
233 }
234 }
235
236
237 ArrayList<Atom> atomList = new ArrayList<>();
238 int atomIndex = 1;
239 for (int j = 0; j < nAtoms; j++) {
240 Atom a = atomList0.get(j);
241 double[] xyz = new double[]{a.getX(), a.getY(), a.getZ()};
242 Atom atom = new Atom(atomIndex++, a.getName(), a.getAtomType(), xyz);
243 atomList.add(atom);
244 }
245
246
247 for (Bond bond : bondList) {
248 Atom a1 = bond.getAtom(0);
249 Atom a2 = bond.getAtom(1);
250
251 int index1 = -1;
252 int index2 = -1;
253 int counter = 0;
254 for (Atom atom : atomList0) {
255 int index = atom.getIndex();
256 if (index == a1.getIndex()) {
257 index1 = counter;
258 } else if (index == a2.getIndex()) {
259 index2 = counter;
260 }
261 counter++;
262 }
263 if (index1 == -1 || index2 == -1) {
264 logger.severe(format(" Index1 (%3d) and Index2 (%3d) suggests error for bond: %s", index1, index2, bond));
265 }
266 Atom newA1 = atomList.get(index1);
267 Atom newA2 = atomList.get(index2);
268 Bond b = new Bond(newA1, newA2);
269 b.setBondType(bond.getBondType());
270 }
271
272
273 newAssembly.setForceField(assembly.getForceField());
274 newAssembly.setPotential(assembly.getPotentialEnergy());
275 newAssembly.setCrystal(assembly.getCrystal());
276 newAssembly.setFile(assembly.getFile());
277 Utilities.biochemistry(newAssembly, atomList);
278
279 XYZFilter filter = new XYZFilter(saveFile, newAssembly, assembly.getForceField(),
280 assembly.getProperties());
281 if (numModels > 1) {
282 filter.writeFile(saveFile, true);
283 } else {
284 filter.writeFile(saveFile, false);
285 }
286 systemFilter.readNext(false, false);
287
288 if (numModels > 1) {
289 logger.info(format(" Saved structure %d to " + saveFile.getName(), i + 1));
290 } else {
291 logger.info(format(" Saved to " + saveFile.getName()));
292 }
293 }
294 }
295 } else if (match && filenames != null && filenames.size() > 1) {
296 potentialFunctions.openAll(filenames.get(0));
297 SystemFilter systemFilter1 = potentialFunctions.getFilter();
298 potentialFunctions.openAll(filenames.get(1));
299 SystemFilter systemFilter2 = potentialFunctions.getFilter();
300 logger.info(" Matching atom order between files.");
301 MolecularAssembly assembly1 = systemFilter1.getActiveMolecularSystem();
302 MolecularAssembly assembly2 = systemFilter2.getActiveMolecularSystem();
303
304 List<MSNode> molecules1 = assembly1.getMolecules();
305 List<MSNode> molecules2 = assembly2.getMolecules();
306 int numMol2 = molecules2.size();
307 boolean[] molDone = new boolean[numMol2];
308 int totalNumAtoms = 0;
309
310 ArrayList<Atom> atomList = new ArrayList<>();
311 int atomIndex = 1;
312 File saveLocation = new File(filenames.get(0));
313 File saveFile = potentialFunctions.versionFile(saveLocation);
314 MolecularAssembly newAssembly = new MolecularAssembly(assembly1.getName());
315 newAssembly.setForceField(assembly1.getForceField());
316 newAssembly.setPotential(assembly1.getPotentialEnergy());
317 newAssembly.setCrystal(assembly1.getCrystal());
318 newAssembly.setFile(assembly1.getFile());
319
320 for (MSNode mol1 : molecules1) {
321 List<Atom> atoms1List = mol1.getAtomList();
322 int numAtoms = atoms1List.size();
323 logger.info(format(" Molecule1 with %3d atoms: %s", numAtoms, mol1.toString()));
324
325
326 List<List<Integer>> equivalents = new ArrayList<>();
327 boolean[] equivalent = new boolean[numAtoms];
328 int numEquiv = 0;
329 for (int i = 0; i < numAtoms; i++) {
330 for (int j = 1; j < numAtoms; j++) {
331 if (i == j) {
332 continue;
333 }
334 AtomComparator ac = new AtomComparator();
335 int value = ac.compare(atoms1List.get(i), atoms1List.get(j));
336 if (value == 0) {
337 equivalent[i] = true;
338 equivalent[j] = true;
339 int iIndex = -1;
340 int jIndex = -1;
341 for (int k = 0; k < numEquiv; k++) {
342 List<Integer> list = equivalents.get(k);
343 if (list.contains(i)) {
344 iIndex = k;
345 }
346 if (list.contains(j)) {
347 jIndex = k;
348 }
349 if (iIndex >= 0 && jIndex >= 0) {
350 break;
351 }
352 }
353 if (iIndex == -1 && jIndex == -1) {
354 List<Integer> ij = new ArrayList<>();
355 ij.add(i);
356 ij.add(j);
357 equivalents.add(ij);
358 numEquiv++;
359 } else if (iIndex >= 0 && jIndex == -1) {
360 equivalents.get(iIndex).add(j);
361 } else if (iIndex == -1 && jIndex >= 0) {
362 equivalents.get(jIndex).add(i);
363 }
364 }
365 }
366 }
367 int numUnique = 0;
368 for (boolean value : equivalent) {
369 if (!value) {
370 numUnique++;
371 }
372 }
373
374
375 double[] xyz1 = new double[numAtoms * 3];
376 double[] mass = new double[numAtoms];
377 double[] compCoords1 = new double[numUnique * 3];
378 double[] compMass = new double[numUnique];
379 int index = 0;
380 int indexUnique = 0;
381 for (int i = 0; i < numAtoms; i++) {
382 Atom a = atoms1List.get(i);
383 if (!equivalent[i]) {
384 compCoords1[indexUnique * 3] = a.getX();
385 compCoords1[indexUnique * 3 + 1] = a.getY();
386 compCoords1[indexUnique * 3 + 2] = a.getZ();
387 compMass[indexUnique++] = a.getMass();
388 }
389 xyz1[index * 3] = a.getX();
390 xyz1[index * 3 + 1] = a.getY();
391 xyz1[index * 3 + 2] = a.getZ();
392 mass[index++] = a.getMass();
393 }
394
395 for (int l = 0; l < numMol2; l++) {
396 MSNode mol2 = molecules2.get(l);
397 List<Atom> atoms2List = mol2.getAtomList();
398 int numAtoms2 = atoms2List.size();
399 logger.info(format(" Molecule2 with %3d atoms: %s", numAtoms2, mol2.toString()));
400 if (numAtoms != numAtoms2 || molDone[l]) {
401 continue;
402 }
403 double[] xyz2 = new double[numAtoms2 * 3];
404 double[] compCoords2 = new double[numUnique * 3];
405 index = 0;
406 indexUnique = 0;
407 for (int i = 0; i < numAtoms; i++) {
408 Atom a = atoms2List.get(i);
409 if (!equivalent[i]) {
410 compCoords2[indexUnique * 3] = a.getX();
411 compCoords2[indexUnique * 3 + 1] = a.getY();
412 compCoords2[indexUnique++ * 3 + 2] = a.getZ();
413 }
414 xyz2[index * 3] = a.getX();
415 xyz2[index * 3 + 1] = a.getY();
416 xyz2[index++ * 3 + 2] = a.getZ();
417 }
418
419
420 if (compCoords1.length == 0 || compCoords2.length == 0) {
421 logger.info(format(" Comparison Coordinates invalid (1: %3d 2: %3d) %3d of %3d unique atoms.", compCoords1.length, compCoords2.length, numUnique, numAtoms));
422 compCoords1 = xyz1;
423 compCoords2 = xyz2;
424 compMass = mass;
425 }
426
427 double[] translate1 = calculateTranslation(compCoords1, compMass);
428 double[] translate2 = calculateTranslation(compCoords2, compMass);
429 applyTranslation(compCoords1, translate1);
430 applyTranslation(compCoords2, translate2);
431
432 applyTranslation(xyz1, translate1);
433 applyTranslation(xyz2, translate2);
434 double[][] rotation = calculateRotation(compCoords1, compCoords2, compMass);
435
436 applyRotation(compCoords2, rotation);
437 applyRotation(xyz2, rotation);
438 double compRMSD = superpose(compCoords1, compCoords2, compMass);
439 logger.info(format(" Unique atom RMSD: %9.4f A", compRMSD));
440 if (compRMSD > 1.0) {
441 logger.warning(format(" Value of %9.3f A may lead to insufficient overlap. Double check produced ordering.", compRMSD));
442 }
443
444
445 for (int i = 0; i < numAtoms; i++) {
446 if (equivalent[i]) {
447 int ind = i * 3;
448 double[] coord1 = new double[]{xyz1[ind], xyz1[ind + 1], xyz1[ind + 2]};
449 double[] coord2 = new double[]{xyz2[ind], xyz2[ind + 1], xyz2[ind + 2]};
450 double value = rmsd(coord1, coord2, new double[]{mass[i]});
451 logger.info(format(" Atom %2d distance of %9.3f A", i, value));
452 for (List<Integer> list : equivalents) {
453 if (list.contains(i)) {
454 int minIndex = -1;
455 double minValue = MAX_VALUE;
456 for (Integer aIndex : list) {
457 if (i != aIndex) {
458 int ind2 = aIndex * 3;
459 double[] tempCoord = new double[]{xyz1[ind2], xyz1[ind2 + 1], xyz1[ind2 + 2]};
460 double value2 = rmsd(tempCoord, coord2, new double[]{mass[i]});
461 if (value2 < minValue) {
462 minIndex = aIndex;
463 minValue = value2;
464 }
465 }
466 }
467 if (minValue < value) {
468
469 Atom temp = atoms1List.get(i);
470 atoms1List.set(i, atoms1List.get(minIndex));
471 atoms1List.set(minIndex, temp);
472
473 double[] tempCoord = new double[]{xyz1[ind], xyz1[ind + 1], xyz1[ind + 2]};
474 int ind2 = minIndex * 3;
475 xyz1[ind] = xyz1[ind2];
476 xyz1[ind + 1] = xyz1[ind2 + 1];
477 xyz1[ind + 2] = xyz1[ind2 + 2];
478 xyz1[ind2] = tempCoord[0];
479 xyz1[ind2 + 1] = tempCoord[1];
480 xyz1[ind2 + 2] = tempCoord[2];
481 }
482 value = rmsd(new double[]{xyz1[ind], xyz1[ind + 1], xyz1[ind + 2]},
483 new double[]{xyz2[ind], xyz2[ind + 1], xyz2[ind + 2]}, new double[]{mass[i]});
484 logger.info(format(" Atom %2d distance updated to %9.3f A", i, value));
485 break;
486 }
487 }
488 }
489 }
490
491 compRMSD = superpose(xyz1, xyz2, mass);
492 logger.info(format(" Final RMSD: %9.4f A", compRMSD));
493
494
495 for (int j = 0; j < numAtoms; j++) {
496 Atom a = atoms1List.get(j);
497 double[] xyz = new double[]{a.getX(), a.getY(), a.getZ()};
498 Atom atom = new Atom(atomIndex++, a.getName(), a.getAtomType(), xyz);
499 atomList.add(atom);
500 }
501
502
503 List<Bond> bondList = mol1.getBondList();
504 for (Bond bond : bondList) {
505 Atom a1 = bond.getAtom(0);
506 Atom a2 = bond.getAtom(1);
507
508 int index1 = -1;
509 int index2 = -1;
510 int counter = 0;
511 for (Atom atom : atoms1List) {
512 int aIndex = atom.getIndex();
513 if (aIndex == a1.getIndex()) {
514 index1 = counter;
515 } else if (aIndex == a2.getIndex()) {
516 index2 = counter;
517 }
518 counter++;
519 }
520 if (index1 == -1 || index2 == -1) {
521 logger.severe(format(" Index1 (%3d) and Index2 (%3d) suggests error for bond: %s", index1, index2, bond));
522 }
523 Atom newA1 = atomList.get(index1 + totalNumAtoms);
524 Atom newA2 = atomList.get(index2 + totalNumAtoms);
525 logger.info(format(" index1 %3d Index2 %3d", index1, index2));
526 if (!newA1.isBonded(newA2)) {
527 Bond b = new Bond(newA1, newA2);
528 b.setBondType(bond.getBondType());
529 }
530 }
531
532 molDone[l] = true;
533 totalNumAtoms += numAtoms;
534 }
535 if (!Arrays.asList(molDone).contains(false)) {
536 break;
537 }
538 }
539
540 Utilities.biochemistry(newAssembly, atomList);
541 XYZFilter filter = new XYZFilter(saveFile, newAssembly, newAssembly.getForceField(),
542 newAssembly.getProperties());
543 filter.writeFile(saveFile, true);
544 logger.info(format(" Saved to " + saveFile.getName()));
545 } else if (atomType && filenames != null && filenames.size() > 1) {
546 potentialFunctions.openAll(filenames.get(0));
547 SystemFilter systemFilter1 = potentialFunctions.getFilter();
548 potentialFunctions.openAll(filenames.get(1));
549 SystemFilter systemFilter2 = potentialFunctions.getFilter();
550 logger.info(" Changing atom types in file 1 to match file 2.");
551
552 MolecularAssembly assembly1 = systemFilter1.getActiveMolecularSystem();
553 do {
554 MolecularAssembly assembly2 = systemFilter2.getActiveMolecularSystem();
555 Atom[] atoms1 = assembly1.getAtomArray();
556 Atom[] atoms2 = assembly2.getAtomArray();
557 ArrayList<ArrayList<Atom>> xyzatoms = new ArrayList<>();
558 MolecularAssembly outputAssembly = new MolecularAssembly(assembly1.getName());
559
560 int numHydrogen = 0;
561 for (Atom atom : atoms2) {
562 if (atom.isHydrogen()) {
563 numHydrogen++;
564 }
565 }
566
567 List<MSNode> entitiesInput = assembly2.getAllBondedEntities();
568 int[] molIndices = assembly2.getMoleculeNumbers();
569 int numEntitiesInput = entitiesInput.size();
570 int[] shiftIndices = new int[numEntitiesInput];
571 int count = 0;
572 for (int i = 0; i < numEntitiesInput; i++) {
573 for (int ind : molIndices) {
574 if (ind == i) {
575 count++;
576 }
577 }
578 shiftIndices[i] = count;
579 }
580 if (logger.isLoggable(Level.FINE)) {
581 logger.fine(format(" Number of entities in input: %d", numEntitiesInput));
582 for (MSNode entity : entitiesInput) {
583 logger.fine(format(" Entity: %s", entity.getName()));
584 int size = entity.getAtomList().size();
585 logger.fine(format(" Entity Size: %3d", size));
586 if (size > 0) {
587 logger.fine(format(" Entity First Atom: %s", entity.getAtomList().get(0).toString()));
588 } else {
589 logger.warning(" Entity did not contain atoms.");
590 }
591 }
592 }
593
594 for (MSNode mol : entitiesInput) {
595
596 xyzatoms.add((ArrayList<Atom>) mol.getAtomList());
597 }
598 int atomIndex = 1;
599
600 logger.info(" Num Entities Input: " + numEntitiesInput);
601 for (int i = 0; i < numEntitiesInput; i++) {
602 MSNode mol = entitiesInput.get(i);
603 int numInputMolAtoms = xyzatoms.get(i).size();
604 int numMolHydrogen = 0;
605 for (Atom atom : xyzatoms.get(i)) {
606 if (atom.isHydrogen()) {
607 numMolHydrogen++;
608 }
609 }
610 if (logger.isLoggable(Level.FINE)) {
611 logger.fine(format(" Current entity number of atoms: %d (%d + %dH)", numInputMolAtoms,
612 numInputMolAtoms - numMolHydrogen, numMolHydrogen));
613 }
614
615
616 AtomContainer xyzCDKAtoms = new AtomContainer();
617 for (Atom atom : xyzatoms.get(i)) {
618 String atomName = getSymbol(atom.getAtomType().atomicNumber);
619 xyzCDKAtoms.addAtom(new org.openscience.cdk.Atom(atomName, new Point3d(atom.getXYZ(null))));
620 }
621
622 int zPrime;
623 int nAtoms = atoms1.length;
624 int nInputAtoms = atoms2.length;
625 if (nAtoms % nInputAtoms == 0) {
626 zPrime = (nAtoms / nInputAtoms);
627 } else if (nAtoms % (nInputAtoms - numHydrogen) == 0) {
628 zPrime = (nAtoms / (nInputAtoms - numHydrogen));
629 } else {
630 zPrime = 1;
631 }
632 logger.info(" Original ZPrime: " + zPrime);
633
634 List<Bond> bonds = mol.getBondList();
635 IBond.Order order = IBond.Order.SINGLE;
636 int xyzBonds = bonds.size();
637 if (xyzBonds == 0) {
638 logger.warning(" No bonds detected in input structure. Please check input.\n " +
639 "If correct, separate non-bonded entities into multiple CIFs.");
640 }
641 int shiftIndex;
642 for (Bond xyzBond : bonds) {
643 int atom0Index = xyzBond.getAtom(0).getXyzIndex();
644 int atom1Index = xyzBond.getAtom(1).getXyzIndex();
645 shiftIndex = (molIndices[atom0Index] == 0) ? 0 : shiftIndices[molIndices[atom0Index] - 1];
646 if (logger.isLoggable(Level.FINER)) {
647 logger.finer(format(" Bonded atom 1: %d, Bonded atom 2: %d", atom0Index, atom1Index));
648 logger.finer(format(" Atom0: %3d Atom1: %3d Shift: %3d final0,1: %3d, %3d MolIndex: %3d ShiftIndex: %3d",
649 atom0Index, atom1Index, shiftIndex, atom0Index - shiftIndex - 1, atom1Index - shiftIndex - 1,
650 molIndices[atom0Index], shiftIndices[molIndices[atom0Index]]));
651 }
652 xyzCDKAtoms.addBond(atom0Index - shiftIndex - 1, atom1Index - shiftIndex - 1, order);
653 }
654
655
656 AtomTypeFactory factory = AtomTypeFactory.getInstance(
657 "org/openscience/cdk/config/data/jmol_atomtypes.txt", xyzCDKAtoms.getBuilder());
658
659
660 for (IAtom atom : xyzCDKAtoms.atoms()) {
661 CIFFilter.setAtomTypes(factory, atom);
662 try {
663 factory.configure(atom);
664 } catch (Exception ex) {
665 logger.info(" Failed to configure atoms from CIF.\n" + ex + "\n" + Arrays.toString(
666 ex.getStackTrace()));
667 }
668 }
669
670 ArrayList<ArrayList<Integer>> zindices = new ArrayList<>();
671 int counter = 0;
672
673 int cifBonds = CIFFilter.bondAtoms(atoms1, bondTolerance);
674 if (logger.isLoggable(Level.FINE)) {
675 logger.fine(
676 format(" Created %d bonds between input atoms (%d in input).", cifBonds, xyzBonds));
677 }
678 List<Atom> atomPool = new ArrayList<>(Arrays.asList(atoms1));
679
680 try {
681 while (!atomPool.isEmpty()) {
682 ArrayList<Atom> molecule = new ArrayList<>();
683 CIFFilter.collectAtoms(atomPool.get(0), molecule);
684 if (logger.isLoggable(Level.FINER)) {
685 logger.finer(format(" Molecule (%d) Size: %d", counter, molecule.size()));
686 }
687 ArrayList<Integer> indices = new ArrayList<>();
688 while (!molecule.isEmpty()) {
689 Atom atom = molecule.remove(0);
690 indices.add(atom.getIndex());
691 atomPool.remove(atom);
692 }
693
694 if (logger.isLoggable(Level.FINER)) {
695 logger.finer(format(
696 " Molecule %d: %d atoms are ready and %d remain (%d atoms in input, %d atoms in CIF). ",
697 counter + 1, indices.size(), atomPool.size(), nInputAtoms, nAtoms));
698 }
699 zindices.add(indices);
700 counter++;
701 }
702 } catch (Exception e) {
703 logger.severe(" Failed to separate copies within the asymmetric unit." + e + "\n" + Utilities.stackTraceToString(e));
704 }
705 zPrime = zindices.size();
706 logger.info(" Zprime: " + zPrime);
707
708 AtomContainer[] cifCDKAtomsArr = new AtomContainer[zPrime];
709 for (int j = 0; j < zPrime; j++) {
710 if (zPrime > 1) {
711 logger.info(format("\n Attempting entity %d of %d", j + 1, zPrime));
712 }
713 ArrayList<Integer> currentList = zindices.get(j);
714 int cifMolAtoms = currentList.size();
715 if (logger.isLoggable(Level.FINE)) {
716 logger.fine(format(" CIF atoms in current: %d", cifMolAtoms));
717 }
718
719 if (cifMolAtoms % numInputMolAtoms == 0 || cifMolAtoms % (numInputMolAtoms - numMolHydrogen) == 0) {
720 cifCDKAtomsArr[j] = new AtomContainer();
721 for (Integer integer : currentList) {
722 String atomName = getSymbol(atoms1[integer - 1].getAtomType().atomicNumber);
723 cifCDKAtomsArr[j].addAtom(new org.openscience.cdk.Atom(atomName, new Point3d(atoms1[integer - 1].getXYZ(null))));
724 }
725 AtomContainer nullCDKAtoms = new AtomContainer();
726
727 for (IAtom atom : cifCDKAtomsArr[j].atoms()) {
728 if (atom.toString() == null) {
729 nullCDKAtoms.addAtom(atom);
730 }
731 }
732 cifCDKAtomsArr[j].remove(nullCDKAtoms);
733
734
735 factory = AtomTypeFactory.getInstance(
736 "org/openscience/cdk/config/data/jmol_atomtypes.txt", cifCDKAtomsArr[j].getBuilder());
737 for (IAtom atom : cifCDKAtomsArr[j].atoms()) {
738 CIFFilter.setAtomTypes(factory, atom);
739 try {
740 factory.configure(atom);
741 } catch (Exception ex) {
742 logger.info(" Failed to configure CIF atoms.\n" + ex + "\n" + Utilities.stackTraceToString(ex));
743 }
744 }
745
746 RebondTool rebonder = new RebondTool(CIFFilter.MAX_COVALENT_RADIUS, CIFFilter.MIN_BOND_DISTANCE,
747 bondTolerance);
748 try {
749 rebonder.rebond(cifCDKAtomsArr[j]);
750 } catch (Exception ex) {
751 logger.info("Failed to rebond CIF atoms.\n" + ex + "\n" + Utilities.stackTraceToString(ex));
752 }
753
754 int cifMolBonds = cifCDKAtomsArr[j].getBondCount();
755 if (logger.isLoggable(Level.FINE)) {
756 logger.fine(format(" Number of CIF bonds: %d (%d in input)", cifMolBonds, xyzBonds));
757 }
758
759
760 if (cifMolBonds != 0 && cifMolBonds % xyzBonds == 0) {
761 Pattern pattern = VentoFoggia.findIdentical(
762 xyzCDKAtoms, AtomMatcher.forElement(), BondMatcher.forAny());
763 int[] p = pattern.match(cifCDKAtomsArr[j]);
764 int pLength = p.length;
765 if (p != null && pLength == numInputMolAtoms) {
766
767 for (int k = 0; k < pLength; k++) {
768 if (logger.isLoggable(Level.FINEST)) {
769 logger.finest(
770 format(" %d input %s -> CIF %s", k, xyzCDKAtoms.getAtom(k).getSymbol(),
771 cifCDKAtomsArr[j].getAtom(p[k]).getSymbol()));
772 }
773 Point3d point3d = cifCDKAtomsArr[j].getAtom(p[k]).getPoint3d();
774 xyzatoms.get(i).get(k).setXYZ(new double[]{point3d.x, point3d.y, point3d.z});
775 }
776 } else {
777 if (logger.isLoggable(Level.FINE)) {
778 logger.fine(
779 format(" Atoms from CIF (%d) and input (%d) structures don't match.", p.length,
780 nAtoms));
781 }
782 continue;
783 }
784 } else if ((xyzBonds - numMolHydrogen) == 0 || cifMolBonds % ((xyzBonds - numMolHydrogen)) == 0) {
785
786 if (logger.isLoggable(Level.FINE)) {
787 logger.info(" CIF may contain implicit hydrogen -- attempting to patch.");
788 }
789
790 Pattern pattern = VentoFoggia.findSubstructure(
791 cifCDKAtomsArr[j], AtomMatcher.forElement(), BondMatcher.forAny());
792 int[] p = pattern.match(xyzCDKAtoms);
793 int pLength = p.length;
794 if (p != null && pLength == numInputMolAtoms - numMolHydrogen) {
795
796 for (int k = 0; k < pLength; k++) {
797 Point3d point3d = cifCDKAtomsArr[j].getAtom(k).getPoint3d();
798 xyzatoms.get(i).get(p[k]).setXYZ(new double[]{point3d.x, point3d.y, point3d.z});
799 }
800
801 Atom lastKnownAtom1 = null;
802 int knownCount = 0;
803 for (Atom hydrogen : xyzatoms.get(i)) {
804 if (hydrogen.isHydrogen()) {
805 Bond bond0 = hydrogen.getBonds().get(0);
806 Atom atom1 = bond0.get1_2(hydrogen);
807 double angle0_2 = 0;
808 List<Angle> anglesList = hydrogen.getAngles();
809 Atom atom2 = null;
810 switch (anglesList.size()) {
811 case 0 ->
812
813
814 hydrogen.moveTo(new double[]{atom1.getX() - 0.6, atom1.getY() - 0.6,
815 atom1.getZ() - 0.6});
816 case 1 -> {
817
818
819 for (Angle angle : anglesList) {
820 atom2 = angle.get1_3(hydrogen);
821 if (atom2 != null) {
822 angle0_2 = angle.angleType.angle[0];
823 }
824 if (angle0_2 != 0) {
825
826 break;
827 }
828 }
829 assert atom2 != null;
830 List<Bond> bonds2 = atom2.getBonds();
831 Atom atom3 =
832 (bonds2.size() > 1 && atom1 == bonds2.get(0).get1_2(atom2)) ? bonds2.get(
833 1).get1_2(atom2) : bonds2.get(0).get1_2(atom2);
834 double diAng = dihedralAngle(hydrogen.getXYZ(null), atom1.getXYZ(null),
835 atom2.getXYZ(null), atom3.getXYZ(null));
836 if (atom1 != atom3) {
837 intxyz(hydrogen, atom1, bond0.bondType.distance, atom2, angle0_2, atom3,
838 Math.toDegrees(diAng), 0);
839 } else {
840
841
842 double[] coord = new double[]{atom2.getX(), atom2.getY(), atom3.getZ()};
843 double mag = sqrt(
844 coord[0] * coord[0] + coord[1] * coord[1] + coord[2] * coord[2]);
845 coord[0] /= mag;
846 coord[1] /= mag;
847 coord[2] /= mag;
848 hydrogen.moveTo(atom1.getX() - coord[0], atom1.getY() - coord[1],
849 atom1.getZ() - coord[2]);
850 }
851 }
852 default -> {
853
854 Atom atom2B = null;
855 double angle0_2B = 0;
856 Atom proposedAtom;
857 double proposedAngle = 0;
858 int chiral = 1;
859 for (Angle angle : anglesList) {
860 proposedAtom = angle.get1_3(hydrogen);
861 if (proposedAtom != null && !proposedAtom.isHydrogen()) {
862 proposedAngle = angle.angleType.angle[0];
863 }
864 if (proposedAngle != 0) {
865
866 if (angle0_2 != 0) {
867 atom2B = proposedAtom;
868 angle0_2B = proposedAngle;
869 } else {
870 atom2 = proposedAtom;
871 angle0_2 = proposedAngle;
872 }
873 proposedAngle = 0.0;
874 }
875 if (angle0_2 != 0 && angle0_2B != 0) {
876 break;
877 }
878 }
879 if (lastKnownAtom1 == null || lastKnownAtom1 != atom1) {
880 lastKnownAtom1 = atom1;
881 knownCount = 0;
882 } else {
883 knownCount++;
884 }
885 if (angle0_2B == 0.0) {
886
887 chiral = 0;
888 assert atom2 != null;
889 List<Bond> bonds2 = atom2.getBonds();
890 atom2B =
891 (atom1 == bonds2.get(0).get1_2(atom2)) ? bonds2.get(1).get1_2(atom2)
892 : bonds2.get(0).get1_2(atom2);
893
894 angle0_2B = 180.0 - 120.0 * knownCount;
895 } else if (anglesList.size() == 2) {
896
897 chiral = 3;
898 } else if (knownCount == 1) {
899
900 chiral = -1;
901 }
902
903 intxyz(hydrogen, atom1, bond0.bondType.distance, atom2, angle0_2, atom2B,
904 angle0_2B, chiral);
905 }
906 }
907 }
908 }
909 } else {
910 logger.info(format(" Substructure match failed (%d vs expected %d)", pLength, numInputMolAtoms - numMolHydrogen));
911 }
912 }
913
914
915 MSGroup molecule;
916 if (mol instanceof Polymer) {
917 molecule = new Polymer(((Polymer) mol).getChainID(), mol.getName(), true);
918 } else {
919 molecule = new Molecule(mol.getName());
920 }
921 List<Atom> newAtomList = new ArrayList<>();
922 for (Atom atom : xyzatoms.get(i)) {
923 Atom molAtom;
924 if (molecule instanceof Polymer) {
925 molAtom = new Atom(atomIndex++, atom.getName(), atom.getAltLoc(), atom.getXYZ(null),
926 atom.getResidueName(), atom.getResidueNumber(), atom.getChainID(), atom.getOccupancy(),
927 atom.getTempFactor(), atom.getSegID());
928 molAtom.setAtomType(atom.getAtomType());
929 } else {
930 molAtom = new Atom(atomIndex++, atom.getName(), atom.getAtomType(), atom.getXYZ(null));
931 if (atom.getResidueName() != null) {
932 molAtom.setResName(atom.getResidueName());
933 }
934 }
935 newAtomList.add(molAtom);
936 }
937 List<Bond> bondList = mol.getBondList();
938 for (Bond bond : bondList) {
939 Atom a1 = bond.getAtom(0);
940 Atom a2 = bond.getAtom(1);
941 if (logger.isLoggable(Level.FINE)) {
942 logger.fine(format(" Bonded atom 1: %d, Bonded atom 2: %d", a1.getXyzIndex(), a2.getXyzIndex()));
943 }
944 int molIndex = a1.getXyzIndex();
945 shiftIndex = (molIndices[molIndex] == 0) ? 0 : shiftIndices[molIndices[molIndex] - 1];
946 Atom newA1 = newAtomList.get(a1.getIndex() - shiftIndex - 1);
947 Atom newA2 = newAtomList.get(a2.getIndex() - shiftIndex - 1);
948 Bond bond2 = new Bond(newA1, newA2);
949 bond2.setBondType(bond.getBondType());
950 }
951 Residue res = null;
952 for (Atom atom : newAtomList) {
953 if (molecule instanceof Polymer) {
954 Polymer polyMol = (Polymer) mol;
955 if (res == null) {
956 res = new Residue(atom.getResidueName(), atom.getResidueNumber(), polyMol.getResidue(atom.getResidueNumber()).getResidueType());
957 } else if (res.getResidueNumber() != atom.getResidueNumber()) {
958 if (logger.isLoggable(Level.FINER)) {
959 logger.finer(" Added Residue " + res.getResidueNumber() + " " + res.getName());
960 }
961 molecule.addMSNode(res);
962 res = new Residue(atom.getResidueName(), atom.getResidueNumber(), polyMol.getResidue(atom.getResidueNumber()).getResidueType());
963 }
964 res.addMSNode(atom);
965 } else {
966 molecule.addMSNode(atom);
967 }
968 }
969 if (molecule instanceof Polymer && res != null) {
970 if (logger.isLoggable(Level.FINER)) {
971 logger.finer(" Added Final Residue " + res.getResidueNumber() + " " + res.getName());
972 }
973 molecule.addMSNode(res);
974 }
975 outputAssembly.addMSNode(molecule);
976 } else {
977 if (logger.isLoggable(Level.INFO)) {
978 logger.info(
979 format(" Number of atoms in CIF (%d) molecule do not match input (%d + %dH = %d).",
980 cifMolAtoms, nInputAtoms - numMolHydrogen, numMolHydrogen, nInputAtoms));
981 }
982 }
983 }
984 }
985 if (logger.isLoggable(Level.FINE)) {
986 logger.fine(format("\n Output Assembly Atoms: %d", outputAssembly.getAtomList().size()));
987 }
988 outputAssembly.setPotential(assembly2.getPotentialEnergy());
989 if (assembly1.getCrystal() != null) {
990 outputAssembly.setCrystal(assembly1.getCrystal());
991 }
992 outputAssembly.setForceField(assembly2.getForceField());
993 outputAssembly.setFile(assembly1.getFile());
994 outputAssembly.setName(assembly1.getName());
995
996 if (outputAssembly.getAtomList().size() < 1) {
997 logger.info(" Atom types could not be matched. File could not be written.");
998 } else if (!writeOutputFile(outputAssembly, systemFilter1)) {
999 logger.info(" Output assembly file could not be written.");
1000 }
1001
1002 outputAssembly.destroy();
1003 } while (systemFilter1.readNext());
1004 } else {
1005 logger.info(helpString());
1006 }
1007 return this;
1008 }
1009
1010
1011
1012
1013
1014
1015
1016
1017 private static boolean writeOutputFile(MolecularAssembly ma, SystemFilter systemFilter) {
1018 File file = systemFilter.getFiles().get(0);
1019 String dir = getFullPath(file.getAbsolutePath()) + File.separator;
1020 String fileName = removeExtension(getName(file.getAbsolutePath()));
1021 String spacegroup;
1022 if (ma.getCrystal() != null) {
1023 spacegroup = ma.getCrystal().getUnitCell().spaceGroup.shortName;
1024 } else {
1025 spacegroup = null;
1026 }
1027 List<MSNode> entities = ma.getAllBondedEntities();
1028 File saveFile;
1029 if (systemFilter instanceof PDBFilter) {
1030
1031 if (entities.size() > 1) {
1032 fileName += "_z" + entities.size();
1033 }
1034 saveFile = new File(dir + fileName + ".pdb");
1035 } else {
1036 if (systemFilter.countNumModels() > 1) {
1037 saveFile = new File(dir + fileName + "_reorder.arc");
1038 } else {
1039
1040 saveFile = XYZFilter.version(new File(dir + fileName + ".xyz"));
1041 }
1042 }
1043 ForceField ff = ma.getForceField();
1044 CompositeConfiguration properties = ma.getProperties();
1045 if (systemFilter instanceof PDBFilter) {
1046 PDBFilter pdbFilter = new PDBFilter(saveFile, ma, ff, properties);
1047 pdbFilter.writeFile(saveFile, true);
1048 } else {
1049 XYZFilter xyzFilter = new XYZFilter(saveFile, ma, ff, properties);
1050 xyzFilter.writeFile(saveFile, true);
1051 }
1052 logger.info("\n Saved output file: " + saveFile.getAbsolutePath());
1053 File propertyFile = new File(saveFile.getParentFile(), removeExtension(saveFile.getName()) + ".properties");
1054 if (!propertyFile.exists()) {
1055 try {
1056 FileWriter fw = new FileWriter(propertyFile, false);
1057 BufferedWriter bw = new BufferedWriter(fw);
1058 ForceField forceField = ma.getForceField();
1059 String parameters = forceField.getString("parameters", "none");
1060 if (parameters != null && !parameters.equalsIgnoreCase("none")) {
1061 bw.write(format("parameters %s\n", parameters));
1062 }
1063 String forceFieldProperty = forceField.getString("forcefield", "none");
1064 if (forceFieldProperty != null && !forceFieldProperty.equalsIgnoreCase("none")) {
1065 bw.write(format("forcefield %s\n", forceFieldProperty));
1066 }
1067 String patch = forceField.getString("patch", "none");
1068 if (patch != null && !patch.equalsIgnoreCase("none")) {
1069 bw.write(format("patch %s\n", patch));
1070 }
1071 bw.write(format("spacegroup %s\n", spacegroup));
1072 if (entities.size() > 1) {
1073 bw.write("intermolecular-softcore true\n");
1074 }
1075 logger.info("\n Saved properties file: " + propertyFile.getAbsolutePath() + "\n");
1076 bw.close();
1077 } catch (Exception ex) {
1078 logger.info("Failed to write files.\n" + Utilities.stackTraceToString(ex));
1079 }
1080 } else {
1081 logger.info("\n Property file already exists: " + propertyFile.getAbsolutePath() + "\n");
1082 }
1083 return true;
1084 }
1085
1086
1087
1088
1089 private final class AtomComparator implements Comparator<Atom> {
1090 @Override
1091 public int compare(Atom a1, Atom a2) {
1092 int comp = Integer.compare(a1.getMoleculeNumber(), a2.getMoleculeNumber());
1093 if (logger.isLoggable(Level.FINER) && comp != 0) {
1094 logger.finer(format(" Different Molecule Number (%d vs %d)", a1.getMoleculeNumber(), a2.getMoleculeNumber()));
1095 }
1096 if (comp == 0) {
1097 comp = Integer.compare(a1.getResidueNumber(), a2.getResidueNumber());
1098 if (logger.isLoggable(Level.FINER) && comp != 0) {
1099 logger.finer(format(" Different Residue Numbers %d vs %d", a1.getResidueNumber(), a2.getResidueNumber()));
1100 }
1101 if (comp == 0) {
1102
1103 comp = -Double.compare(a1.getMass(), a2.getMass());
1104 if (logger.isLoggable(Level.FINER) && comp != 0) {
1105 logger.finer(format(" Different Masses %6.3f %6.3f", a1.getMass(), a2.getMass()));
1106 }
1107 if (comp == 0) {
1108 comp = Integer.compare(a1.getAtomType().type, a2.getAtomType().type);
1109 if (logger.isLoggable(Level.FINER) && comp != 0) {
1110 logger.finer(format(" Different Atom Types (%d vs %d)", a1.getAtomType().type, a2.getAtomType().type));
1111 }
1112 if (comp == 0) {
1113 comp = Integer.compare(a1.getBonds().size(), a2.getBonds().size());
1114 if (logger.isLoggable(Level.FINER) && comp != 0) {
1115 logger.finer(format(" Different Bond Sizes (%d vs %d)", a1.getBonds().size(), a2.getBonds().size()));
1116 }
1117 if (comp == 0) {
1118 comp = Integer.compare(a1.get12List().size(), a2.get12List().size());
1119 if (logger.isLoggable(Level.FINER) && comp != 0) {
1120 logger.finer(format(" Different 1-2 Size (%d vs %d).", a1.get12List().size(), a2.get12List().size()));
1121 }
1122 if (comp == 0) {
1123 for (int i = 0; i < a1.get12List().size(); i++) {
1124 Atom a12 = a1.get12List().get(i);
1125 Atom a22 = a2.get12List().get(i);
1126 comp = shallowCompare(a12, a22);
1127 if (comp != 0) {
1128 break;
1129 }
1130 }
1131 if (logger.isLoggable(Level.FINER) && comp != 0) {
1132 logger.finer(" Different 1-2 Shallow.");
1133 }
1134 if (comp == 0) {
1135 comp = Integer.compare(a1.get13List().size(), a2.get13List().size());
1136 if (logger.isLoggable(Level.FINER) && comp != 0) {
1137 logger.finer(" Different 1-3 Size.");
1138 }
1139 if (comp == 0) {
1140 for (int i = 0; i < a1.get13List().size(); i++) {
1141 Atom a13 = a1.get13List().get(i);
1142 Atom a23 = a2.get13List().get(i);
1143 comp = shallowCompare(a13, a23);
1144 if (comp != 0) {
1145 break;
1146 }
1147 }
1148 if (logger.isLoggable(Level.FINER) && comp != 0) {
1149 logger.finer(" Different 1-3 Shallow.");
1150 }
1151 if (comp == 0) {
1152 comp = Integer.compare(a1.get14List().size(), a2.get14List().size());
1153 if (logger.isLoggable(Level.FINER) && comp != 0) {
1154 logger.finer(" Different 1-4 Size.");
1155 }
1156 if (comp == 0) {
1157 for (int i = 0; i < a1.get14List().size(); i++) {
1158 Atom a14 = a1.get14List().get(i);
1159 Atom a24 = a2.get14List().get(i);
1160 comp = shallowCompare(a14, a24);
1161 if (comp != 0) {
1162 break;
1163 }
1164 }
1165 if (logger.isLoggable(Level.FINER) && comp != 0) {
1166 logger.finer(" Different 1-4 Shallow.");
1167 }
1168 if (comp == 0) {
1169 comp = Integer.compare(a1.get15List().size(), a2.get15List().size());
1170 if (logger.isLoggable(Level.FINER) && comp != 0) {
1171 logger.finer(" Different 1-5 Size.");
1172 }
1173 if (comp == 0) {
1174 for (int i = 0; i < a1.get15List().size(); i++) {
1175 Atom a15 = a1.get15List().get(i);
1176 Atom a25 = a2.get15List().get(i);
1177 comp = shallowCompare(a15, a25);
1178 if (comp != 0) {
1179 break;
1180 }
1181 }
1182 if (logger.isLoggable(Level.FINER) && comp != 0) {
1183 logger.finer(" Different 1-5 Shallow.");
1184 }
1185 }
1186 }
1187 }
1188 }
1189 }
1190 }
1191 }
1192 }
1193 }
1194 }
1195 }
1196 }
1197 return comp;
1198 }
1199
1200
1201
1202
1203
1204
1205
1206
1207 private int shallowCompare(Atom a1, Atom a2) {
1208 int comp = Integer.compare(a1.getBonds().size(), a2.getBonds().size());
1209 if (logger.isLoggable(Level.FINER) && comp != 0) {
1210 logger.finer(" Different Shallow Bond Size.");
1211 }
1212 if (comp == 0) {
1213 comp = Integer.compare(a1.getAtomType().type, a2.getAtomType().type);
1214 if (logger.isLoggable(Level.FINER) && comp != 0) {
1215 logger.finer(" Different Shallow Atom Type.");
1216 }
1217 if (comp == 0) {
1218 comp = Integer.compare(a1.get12List().size(), a2.get12List().size());
1219 if (logger.isLoggable(Level.FINER) && comp != 0) {
1220 logger.finer(" Different Shallow 1-2 Size.");
1221 }
1222 if (comp == 0) {
1223 comp = Integer.compare(a1.get13List().size(), a2.get13List().size());
1224 if (logger.isLoggable(Level.FINER) && comp != 0) {
1225 logger.finer(" Different Shallow 1-3 Size.");
1226 }
1227 if (comp == 0) {
1228 comp = Integer.compare(a1.get14List().size(), a2.get14List().size());
1229 if (logger.isLoggable(Level.FINER) && comp != 0) {
1230 logger.finer(" Different Shallow 1-4 Size.");
1231 }
1232 if (comp == 0) {
1233 comp = Integer.compare(a1.get15List().size(), a2.get15List().size());
1234 if (logger.isLoggable(Level.FINER) && comp != 0) {
1235 logger.finer(" Different Shallow 1-5 Size.");
1236 }
1237 }
1238 }
1239 }
1240 }
1241 }
1242 return comp;
1243 }
1244 }
1245 }