Source code: joelib/util/JHM.java
1 ///////////////////////////////////////////////////////////////////////////////
2 // Filename: $RCSfile: JHM.java,v $
3 // Purpose: Atom representation.
4 // Language: Java
5 // Compiler: JDK 1.4
6 // Authors: Joerg K. Wegner
7 // Version: $Revision: 1.21 $
8 // $Date: 2003/08/22 15:56:21 $
9 // $Author: wegner $
10 // Original Author: ???, OpenEye Scientific Software
11 // Original Version: babel 2.0a1
12 //
13 // Copyright (c) Dept. Computer Architecture, University of Tuebingen, Germany
14 //
15 // This program is free software; you can redistribute it and/or modify
16 // it under the terms of the GNU General Public License as published by
17 // the Free Software Foundation version 2 of the License.
18 //
19 // This program is distributed in the hope that it will be useful,
20 // but WITHOUT ANY WARRANTY; without even the implied warranty of
21 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 // GNU General Public License for more details.
23 ///////////////////////////////////////////////////////////////////////////////
24 package joelib.util;
25
26 import joelib.ext.ExternalHelper;
27
28 import joelib.io.IOType;
29 import joelib.io.JOEFileFormat;
30 import joelib.io.MoleculeFileType;
31 import joelib.io.MoleculeIOException;
32 import joelib.io.PropertyWriter;
33
34 import joelib.math.Matrix3x3;
35 import joelib.math.XYZVector;
36
37 import joelib.molecule.JOEAtom;
38 import joelib.molecule.JOEBond;
39 import joelib.molecule.JOEMol;
40 import joelib.molecule.JOERTree;
41
42 import joelib.smarts.JOESmartsPattern;
43
44 import joelib.util.iterator.AtomIterator;
45 import joelib.util.iterator.BondIterator;
46 import joelib.util.iterator.NbrAtomIterator;
47
48 import joelib.util.types.JOEInternalCoord;
49
50 import wsi.ra.tool.PropertyHolder;
51
52 /*==========================================================================*
53 * IMPORTS
54 *========================================================================== */
55 import java.io.FileOutputStream;
56 import java.io.IOException;
57 import java.io.OutputStream;
58 import java.io.PrintStream;
59
60 import java.util.Properties;
61 import java.util.StringTokenizer;
62 import java.util.Vector;
63
64 import org.apache.log4j.Category;
65
66
67 /*==========================================================================*
68 * CLASS DECLARATION
69 *========================================================================== */
70
71 /**
72 * JOELib hepler methods.
73 *
74 * @author wegnerj
75 * @license GPL
76 * @cvsversion $Revision: 1.21 $, $Date: 2003/08/22 15:56:21 $
77 */
78 public class JHM
79 {
80 //~ Static fields/initializers /////////////////////////////////////////////
81
82 /*-------------------------------------------------------------------------*
83 * public static member variables
84 *------------------------------------------------------------------------- */
85
86 // Obtain a suitable logger.
87 private static Category logger = Category.getInstance("joelib.util.JHM");
88
89 /**
90 * Description of the Field
91 */
92 public final static double RAD_TO_DEG = 180.0f / Math.PI;
93
94 /**
95 * Description of the Field
96 */
97 public final static double DEG_TO_RAD = Math.PI / 180.0f;
98 public static final String eol = System.getProperty("line.separator", "\n");
99
100 /*-------------------------------------------------------------------------*
101 * private member variables
102 *------------------------------------------------------------------------- */
103 private final static int OE_RTREE_CUTOFF = 20;
104 private static JHM jhm;
105
106 //~ Constructors ///////////////////////////////////////////////////////////
107
108 /*-------------------------------------------------------------------------*
109 * constructor
110 *------------------------------------------------------------------------- */
111
112 /**
113 * Constructor for the JHM object
114 */
115 private JHM()
116 {
117 }
118
119 //~ Methods ////////////////////////////////////////////////////////////////
120
121 /**
122 * Gets the currentValence attribute of the JHM class
123 *
124 * @param atom Description of the Parameter
125 * @return The currentValence value
126 */
127 public static int getCurrentValence(JOEAtom atom)
128 {
129 int count = atom.getImplicitValence();
130
131 JOEBond bond;
132 BondIterator bit = atom.bondIterator();
133
134 while (bit.hasNext())
135 {
136 bond = bit.nextBond();
137
138 if (bond.isKDouble())
139 {
140 count++;
141 }
142 else if (bond.isKTriple())
143 {
144 count += 2;
145 }
146 }
147
148 return (count);
149 }
150
151 /*public static int PDBSort(ChainsAtom arg1[], ChainsAtom arg2[])
152 *{
153 *throw Exception("use ChainAtomComparator");
154 *}
155 *static void outputPDBFile(ChainsMolecule mol, OutputStream os)
156 *{
157 *int src,dst;
158 *ChainsAtom atom;
159 *String ptr;
160 *char tmpc[]=new char[1];
161 *int i;
162 *PrintfStream fp = new PrintfStream(os);
163 *for( i=0; i<mol.acount; i++ )
164 *PDBOrder[i] = mol.atom[i];
165 *qsort(PDBOrder,mol.acount, new ChainAtomComparator());
166 *ptr = mol.name;
167 *if( ptr!=null )
168 *{
169 *fp.print("COMPND ");
170 *fp.prinln( ptr.getBytes() );
171 *}
172 *for( i=0; i<mol.acount; i++ )
173 *{
174 *atom = PDBOrder[i];
175 *atom.serno = i+1;
176 *if( atom.hetflag )
177 *{
178 *fp.print( "HETATM" );
179 *}
180 *else
181 *{
182 *fp.print( "ATOM " );
183 *}
184 *fp.printf("%4d ", atom.serno);
185 *if( atom.atomid == -1 )
186 *{
187 *fp.printf("%s ", chainsElem[atom.elem].symbol);
188 *}
189 *else
190 *{
191 *if( atom.elem == 1 )
192 *{
193 *if( atom.hcount )
194 *{
195 *tmpc[0] = atom.hcount+'0';
196 *fp.write(tmpc);
197 *}
198 *else
199 *{
200 *tmpc[0] = ' ';
201 *fp.write(tmpc);
202 *}
203 *fp.printf("H%.2s", chainsAtomName[atom.atomid]+2);
204 *}
205 *else
206 *{
207 *fp.printf( ,"%.4s", chainsAtomName[atom.atomid]);
208 *}
209 *fp.printf( " %s ", chainsResName[atom.resid]);
210 *fp.printf( "%c%4d", atom.chain);
211 *fp.printf( "%4d", atom.resno);
212 *fp.printf( " %8.3lf", atom.x);
213 *fp.printf( "%8.3lf", atom.y);
214 *fp.printf( "%8.3lf", atom.z);
215 *fp.println(" 1.00 0.00");
216 *}
217 *for( i=0; i<mol.bcount; i++ )
218 *{
219 *if( mol.bond[i].flag & BF_DOUBLE )
220 *{
221 *src = mol.atom[mol.bond[i].src].serno;
222 *dst = mol.atom[mol.bond[i].dst].serno;
223 *fp.printf( "CONECT%5d", src);
224 *fp.printf( "%5d", dst);
225 *fp.printf( "%5d\n", dst);
226 *fp.printf( "CONECT%5d", dst);
227 *fp.printf( "%5d", src);
228 *fp.printf( "%5d\n", src);
229 *}
230 *}
231 *fp.println("END ");
232 *} */
233
234 /**
235 * @param dffv {@link java.util.Vector} of <tt>int[1]</tt>
236 * @param mol Description of the Parameter
237 * @param bv Description of the Parameter
238 * @return The dFFVector value
239 */
240 public static boolean getDFFVector(JOEMol mol, Vector dffv, JOEBitVec bv)
241 {
242 dffv.clear();
243 dffv.setSize(mol.numAtoms());
244
245 int dffcount;
246
247 int natom;
248 JOEBitVec used = new JOEBitVec();
249 JOEBitVec curr = new JOEBitVec();
250 JOEBitVec next = new JOEBitVec();
251 JOEAtom atom;
252 JOEAtom atom1;
253 JOEBond bond;
254
255 next.clear();
256
257 AtomIterator ait = mol.atomIterator();
258
259 while (ait.hasNext())
260 {
261 atom = ait.nextAtom();
262
263 if (bv.get(atom.getIdx()))
264 {
265 ((int[]) dffv.get(atom.getIdx() - 1))[0] = 0;
266
267 continue;
268 }
269
270 dffcount = 0;
271 used.clear();
272 curr.clear();
273 used.setBitOn(atom.getIdx());
274 curr.setBitOn(atom.getIdx());
275
276 while (!curr.isEmpty() && ((bv.and(curr)).size() == 0))
277 {
278 next.clear();
279
280 // for(natom=curr.nextSetBit(0); natom>=0; natom=curr.nextSetBit(natom+1))
281 for (natom = curr.nextBit(-1); natom != curr.endBit();
282 natom = curr.nextBit(natom))
283 {
284 atom1 = mol.getAtom(natom);
285
286 BondIterator bit = atom1.bondIterator();
287
288 while (bit.hasNext())
289 {
290 bond = bit.nextBond();
291
292 if (!used.bitIsOn(bond.getNbrAtomIdx(atom1)) &&
293 !curr.bitIsOn(bond.getNbrAtomIdx(atom1)))
294 {
295 if (!(bond.getNbrAtom(atom1)).isHydrogen())
296 {
297 next.setBitOn(bond.getNbrAtomIdx(atom1));
298 }
299 }
300 }
301 }
302
303 used.orSet(next);
304 curr.set(next);
305 dffcount++;
306 }
307
308 dffv.set(atom.getIdx() - 1, new int[]{dffcount});
309 }
310
311 return true;
312 }
313
314 /**
315 * @param v of type <tt>int[1]</tt>
316 * @param value The new intVectorToValue value
317 */
318 public static void setIntVectorToValue(Vector v, int value)
319 {
320 int[] itmp;
321
322 for (int i = 0; i < v.size(); i++)
323 {
324 itmp = (int[]) v.get(i);
325
326 if (itmp != null)
327 {
328 itmp[0] = value;
329 }
330 else
331 {
332 v.set(i, new int[]{value});
333 }
334 }
335 }
336
337 /**
338 * This methodName will rotate the coordinates of 'atoms' such that tor == ang
339 * - atoms in 'tor' should be ordered such that the 3rd atom is the pivot
340 * around which atoms rotate ang is in degrees.
341 *
342 * @param atoms {@link java.util.Vector} of <tt>int[1]</tt>
343 * @param c The new rotorToAngle value
344 * @param ref The new rotorToAngle value
345 * @param ang The new rotorToAngle value
346 */
347 public static void setRotorToAngle(double[] c, JOEAtom[] ref, double ang,
348 Vector atoms)
349 {
350 double v1x;
351 double v1y;
352 double v1z;
353 double v2x;
354 double v2y;
355 double v2z;
356 double v3x;
357 double v3y;
358 double v3z;
359 double c1x;
360 double c1y;
361 double c1z;
362 double c2x;
363 double c2y;
364 double c2z;
365 double c3x;
366 double c3y;
367 double c3z;
368 double c1mag;
369 double c2mag;
370 double radang;
371 double costheta;
372 double[] m = new double[9];
373 double x;
374 double y;
375 double z;
376 double mag;
377 double rotang;
378 double sn;
379 double cs;
380 double t;
381 double tx;
382 double ty;
383 double tz;
384
385 int[] tor = new int[4];
386 tor[0] = ref[0].getCIdx();
387 tor[1] = ref[1].getCIdx();
388 tor[2] = ref[2].getCIdx();
389 tor[3] = ref[3].getCIdx();
390
391 //
392 //calculate the torsion angle
393 //
394 v1x = c[tor[0]] - c[tor[1]];
395 v2x = c[tor[1]] - c[tor[2]];
396 v1y = c[tor[0] + 1] - c[tor[1] + 1];
397 v2y = c[tor[1] + 1] - c[tor[2] + 1];
398 v1z = c[tor[0] + 2] - c[tor[1] + 2];
399 v2z = c[tor[1] + 2] - c[tor[2] + 2];
400 v3x = c[tor[2]] - c[tor[3]];
401 v3y = c[tor[2] + 1] - c[tor[3] + 1];
402 v3z = c[tor[2] + 2] - c[tor[3] + 2];
403
404 c1x = (v1y * v2z) - (v1z * v2y);
405 c2x = (v2y * v3z) - (v2z * v3y);
406 c1y = (-v1x * v2z) + (v1z * v2x);
407 c2y = (-v2x * v3z) + (v2z * v3x);
408 c1z = (v1x * v2y) - (v1y * v2x);
409 c2z = (v2x * v3y) - (v2y * v3x);
410 c3x = (c1y * c2z) - (c1z * c2y);
411 c3y = (-c1x * c2z) + (c1z * c2x);
412 c3z = (c1x * c2y) - (c1y * c2x);
413
414 c1mag = (c1x * c1x) + (c1y * c1y) + (c1z * c1z);
415 c2mag = (c2x * c2x) + (c2y * c2y) + (c2z * c2z);
416
417 if ((c1mag * c2mag) < 0.01)
418 {
419 costheta = 1.0;
420 }
421
422 //avoid div by zero error
423 else
424 {
425 costheta = ((c1x * c2x) + (c1y * c2y) + (c1z * c2z)) / (Math.sqrt(c1mag * c2mag));
426 }
427
428 if (costheta < -0.999999)
429 {
430 costheta = -0.999999f;
431 }
432
433 if (costheta > 0.999999)
434 {
435 costheta = 0.999999f;
436 }
437
438 if (((v2x * c3x) + (v2y * c3y) + (v2z * c3z)) > 0.0)
439 {
440 radang = -Math.acos(costheta);
441 }
442 else
443 {
444 radang = Math.acos(costheta);
445 }
446
447 //
448 // now we have the torsion angle (radang) - set up the rot matrix
449 //
450 //find the difference between current and requested
451 rotang = (DEG_TO_RAD * ang) - radang;
452
453 sn = Math.sin(rotang);
454 cs = Math.cos(rotang);
455 t = 1 - cs;
456
457 //normalize the rotation vector
458 mag = Math.sqrt((v2x * v2x) + (v2y * v2y) + (v2z * v2z));
459 x = v2x / mag;
460 y = v2y / mag;
461 z = v2z / mag;
462
463 //set up the rotation matrix
464 m[0] = (t * x * x) + cs;
465 m[1] = (t * x * y) + (sn * z);
466 m[2] = (t * x * z) - (sn * y);
467 m[3] = (t * x * y) - (sn * z);
468 m[4] = (t * y * y) + cs;
469 m[5] = (t * y * z) + (sn * x);
470 m[6] = (t * x * z) + (sn * y);
471 m[7] = (t * y * z) - (sn * x);
472 m[8] = (t * z * z) + cs;
473
474 //
475 //now the matrix is set - time to rotate the atoms
476 //
477 tx = c[tor[1]];
478 ty = c[tor[1] + 1];
479 tz = c[tor[1] + 2];
480
481 int j;
482
483 for (int i = 0; i < atoms.size(); i++)
484 {
485 int[] atom = (int[]) atoms.get(i);
486 j = (atom[0] - 1) * 3;
487 c[j] -= tx;
488 c[j + 1] -= ty;
489 c[j + 2] -= tz;
490
491 x = (c[j] * m[0]) + (c[j + 1] * m[1]) + (c[j + 2] * m[2]);
492 y = (c[j] * m[3]) + (c[j + 1] * m[4]) + (c[j + 2] * m[5]);
493 z = (c[j] * m[6]) + (c[j + 1] * m[7]) + (c[j + 2] * m[8]);
494 c[j] = x;
495 c[j + 1] = y;
496 c[j + 2] = z;
497 c[j] += tx;
498 c[j + 1] += ty;
499 c[j + 2] += tz;
500 }
501 }
502
503 public static String getTempFileBase()
504 {
505 Properties prop = PropertyHolder.instance().getProperties();
506 String tempDir = "joelib.temporary.directory." +
507 ExternalHelper.getOperationSystemName();
508 String temp = prop.getProperty(tempDir);
509
510 if (temp == null)
511 {
512 logger.error("You must define a temporary directory in '" +
513 tempDir + "'.");
514
515 return null;
516 }
517
518 String fileSeparator = System.getProperty("file.separator");
519
520 return temp + fileSeparator;
521 }
522
523 /**
524 * @param prv {@link java.util.Vector} of <tt>JOERTree</tt>
525 * @param atom Description of the Parameter
526 * @param vt Description of the Parameter
527 * @param bv Description of the Parameter
528 */
529 public static void buildOERTreeVector(JOEAtom atom, JOERTree prv,
530 JOERTree[] vt, JOEBitVec bv)
531 {
532 vt[atom.getIdx()] = new JOERTree(atom, prv);
533
534 int i;
535 JOEAtom nbr;
536 JOEMol mol = atom.getParent();
537 JOEBitVec curr = new JOEBitVec();
538 JOEBitVec next = new JOEBitVec();
539 JOEBitVec used;
540
541 curr.set(atom.getIdx());
542
543 // System.out.println("b20"+bv.get(20));
544 used = JOEBitVec.or(bv, curr);
545
546 // System.out.println("b20"+used.get(20));
547 // System.out.println("first: "+atom.getIdx());
548 int level = 0;
549
550 for (;;)
551 {
552 next.clear();
553
554 // for(i=curr.nextSetBit(0); i>=0; i=curr.nextSetBit(i+1))
555 for (i = curr.nextBit(0); i != bv.endBit(); i = curr.nextBit(i))
556 {
557 // System.out.println("bit "+i);
558 atom = mol.getAtom(i);
559
560 NbrAtomIterator nait = atom.nbrAtomIterator();
561
562 while (nait.hasNext())
563 {
564 nbr = nait.nextNbrAtom();
565
566 // System.out.println("nbr "+nbr.getIdx()+" u "+used.get(nbr.getIdx()));
567 if (!used.get(nbr.getIdx()))
568 {
569 next.set(nbr.getIdx());
570 used.set(nbr.getIdx());
571
572 // System.out.println("tree "+nbr.getIdx());
573 vt[nbr.getIdx()] = new JOERTree(nbr, vt[atom.getIdx()]);
574 }
575 }
576 }
577
578 if (next.size() == 0)
579 {
580 break;
581 }
582
583 curr.set(next);
584 level++;
585
586 if (level > OE_RTREE_CUTOFF)
587 {
588 break;
589 }
590 }
591 }
592
593 /**
594 * Description of the Method
595 *
596 * @param mol Description of the Parameter
597 */
598 public static void correctBadResonanceForm(JOEMol mol)
599 {
600 JOEBond b1;
601 JOEBond b2;
602 JOEBond b3;
603 JOEAtom a1;
604 JOEAtom a2;
605 JOEAtom a3;
606 JOEAtom a4;
607 Vector mlist;
608
609 //Vector of int[]
610 //carboxylic acid
611 JOESmartsPattern acid = new JOESmartsPattern();
612 acid.init("[oD1]c[oD1]");
613
614 //Vector vec;
615 int[] itmp;
616
617 if (acid.match(mol))
618 {
619 mlist = acid.getUMapList();
620
621 for (int i = 0; i < mlist.size(); i++)
622 {
623 itmp = (int[]) mlist.get(i);
624
625 //vec = (Vector)mlist.get(i);
626 //itmp = (int[])vec.get(0); a1 = mol.getAtom( itmp[0] );
627 //itmp = (int[])vec.get(1); a2 = mol.getAtom( itmp[0] );
628 //itmp = (int[])vec.get(2); a3 = mol.getAtom( itmp[0] );
629 a1 = mol.getAtom(itmp[0]);
630 a2 = mol.getAtom(itmp[1]);
631 a3 = mol.getAtom(itmp[2]);
632 b1 = a2.getBond(a1);
633 b2 = a2.getBond(a3);
634
635 if ((b1 == null) || (b2 == null))
636 {
637 continue;
638 }
639
640 b1.setKDouble();
641 b2.setKSingle();
642 }
643 }
644
645 //phosphonic acid
646 JOESmartsPattern phosphate = new JOESmartsPattern();
647 phosphate.init("[p]([oD1])([oD1])([oD1])[#6,#8]");
648
649 //Vector vec;
650 if (phosphate.match(mol))
651 {
652 mlist = phosphate.getUMapList();
653
654 for (int i = 0; i < mlist.size(); i++)
655 {
656 //vec = (Vector)mlist.get(i);
657 //itmp = (int[])vec.get(0); a1 = mol.getAtom( itmp[0] );
658 //itmp = (int[])vec.get(1); a2 = mol.getAtom( itmp[0] );
659 //itmp = (int[])vec.get(2); a3 = mol.getAtom( itmp[0] );
660 //itmp = (int[])vec.get(3); a4 = mol.getAtom( itmp[0] );
661 itmp = (int[]) mlist.get(i);
662 a1 = mol.getAtom(itmp[0]);
663 a2 = mol.getAtom(itmp[1]);
664 a3 = mol.getAtom(itmp[2]);
665 a4 = mol.getAtom(itmp[3]);
666 b1 = a1.getBond(a2);
667 b2 = a1.getBond(a3);
668 b3 = a1.getBond(a4);
669
670 if ((b1 == null) || (b2 == null) || (b3 == null))
671 {
672 continue;
673 }
674
675 b1.setKDouble();
676 b2.setKSingle();
677 b3.setKSingle();
678 }
679 }
680
681 //amidene and guanidine
682 JOESmartsPattern amidene = new JOESmartsPattern();
683 amidene.init("[nD1]c([nD1])*");
684
685 if (amidene.match(mol))
686 {
687 mlist = amidene.getUMapList();
688
689 for (int i = 0; i < mlist.size(); i++)
690 {
691 //vec = (Vector)mlist.get(i);
692 //itmp = (int[])vec.get(0); a1 = mol.getAtom( itmp[0] );
693 //itmp = (int[])vec.get(1); a2 = mol.getAtom( itmp[0] );
694 //itmp = (int[])vec.get(2); a3 = mol.getAtom( itmp[0] );
695 itmp = (int[]) mlist.get(i);
696 a1 = mol.getAtom(itmp[0]);
697 a2 = mol.getAtom(itmp[1]);
698 a3 = mol.getAtom(itmp[2]);
699
700 b1 = a2.getBond(a1);
701 b2 = a2.getBond(a3);
702
703 if ((b1 == null) || (b2 == null))
704 {
705 continue;
706 }
707
708 b1.setKDouble();
709 b2.setKSingle();
710 }
711 }
712 }
713
714 /**
715 * @param va {@link java.util.Vector} of <tt>JOEAtom</tt>
716 * @param mol Description of the Parameter
717 * @param i Description of the Parameter
718 * @param maxv Description of the Parameter
719 * @param secondpass Description of the Parameter
720 * @return Description of the Return Value
721 */
722 public static boolean expandKekule(JOEMol mol, Vector va, AtomIterator ii,
723 int[] maxv, boolean secondpass)
724 {
725 AtomIterator i = (AtomIterator) ii.clone();
726
727 // System.out.println("expKekul:"+i.actualIndex()+" "+i.hasNext());
728 if (!i.hasNext())
729 {
730 //check to see that the ideal valence has been achieved for all atoms
731 JOEAtom atom;
732 AtomIterator ait = new AtomIterator(va);
733
734 while (ait.hasNext())
735 {
736 atom = ait.nextAtom();
737
738 //let erroneously aromatic carboxylates pass
739 if (atom.isOxygen() && (atom.getValence() == 1))
740 {
741 continue;
742 }
743
744 if (getCurrentValence(atom) != maxv[atom.getIdx()])
745 {
746 return false;
747 }
748 }
749
750 return true;
751 }
752
753 //jump to next atom in list if current atom doesn't have any attached
754 //aromatic bonds
755 JOEBond bond;
756 JOEAtom atom = (JOEAtom) i.actual();
757
758 BondIterator bit = atom.bondIterator();
759 boolean done = true;
760
761 while (bit.hasNext())
762 {
763 bond = bit.nextBond();
764
765 if (bond.getBO() == 5)
766 {
767 done = false;
768
769 break;
770 }
771 }
772
773 if (done)
774 {
775 i.setActualIndex(i.actualIndex() + 1);
776
777 return (expandKekule(mol, va, i, maxv, secondpass));
778
779 // i.setActualIndex( i.actualIndex()-1 );
780 }
781
782 //store list of attached aromatic atoms
783 JOEAtom nbr;
784 Vector vb = new Vector();
785
786 //type JOEBond
787 NbrAtomIterator nait = atom.nbrAtomIterator();
788
789 while (nait.hasNext())
790 {
791 nbr = nait.nextNbrAtom();
792 bond = nait.actualBond();
793
794 if (bond.getBO() == 5)
795 {
796 vb.add(bond);
797 bond.setBO(1);
798 bond.setKSingle();
799 }
800 }
801
802 //try setting a double bond
803 if (getCurrentValence(atom) < maxv[atom.getIdx()])
804 {
805 for (int j = 0; j < vb.size(); j++)
806 {
807 bond = (JOEBond) vb.get(j);
808 nbr = bond.getNbrAtom(atom);
809
810 if (getCurrentValence(nbr) <= maxv[nbr.getIdx()])
811 {
812 bond.setKDouble();
813 bond.setBO(2);
814 i.setActualIndex(i.actualIndex() + 1);
815
816 if (expandKekule(mol, va, i, maxv, secondpass))
817 {
818 return true;
819 }
820
821 i.setActualIndex(i.actualIndex() - 1);
822 bond.setKSingle();
823 bond.setBO(1);
824 }
825 }
826
827 if (secondpass && atom.isNitrogen() &&
828 (atom.getFormalCharge() == 0) &&
829 (atom.getImplicitValence() == 2))
830 {
831 atom.incrementImplicitValence();
832 i.setActualIndex(i.actualIndex() + 1);
833
834 if (expandKekule(mol, va, i, maxv, secondpass))
835 {
836 return true;
837 }
838
839 i.setActualIndex(i.actualIndex() - 1);
840 atom.decrementImplicitValence();
841 }
842 }
843 else
844 {
845 //full valence - no double bond to set
846 i.setActualIndex(i.actualIndex() + 1);
847
848 if (expandKekule(mol, va, i, maxv, secondpass))
849 {
850 return true;
851 }
852
853 i.setActualIndex(i.actualIndex() - 1);
854
855 boolean trycharge = false;
856
857 if (secondpass && (atom.getFormalCharge() == 0))
858 {
859 if (atom.isNitrogen() && (atom.getHvyValence() == 3))
860 {
861 trycharge = true;
862 }
863
864 if (atom.isOxygen() && (atom.getHvyValence() == 2))
865 {
866 trycharge = true;
867 }
868
869 if (atom.isSulfur() && (atom.getHvyValence() == 2))
870 {
871 trycharge = true;
872 }
873 }
874
875 if (trycharge)
876 {
877 maxv[atom.getIdx()]++;
878 atom.setFormalCharge(1);
879
880 for (int j = 0; j < vb.size(); j++)
881 {
882 bond = (JOEBond) vb.get(j);
883 nbr = bond.getNbrAtom(atom);
884
885 if (getCurrentValence(nbr) <= maxv[nbr.getIdx()])
886 {
887 bond.setKDouble();
888 bond.setBO(2);
889 i.setActualIndex(i.actualIndex() + 1);
890
891 if (expandKekule(mol, va, i, maxv, secondpass))
892 {
893 return true;
894 }
895
896 i.setActualIndex(i.actualIndex() - 1);
897 bond.setKSingle();
898 bond.setBO(1);
899 }
900 }
901
902 maxv[atom.getIdx()]--;
903 atom.setFormalCharge(0);
904 }
905
906 if (secondpass && atom.isNitrogen() &&
907 (atom.getFormalCharge() == 0) &&
908 (atom.getImplicitValence() == 2))
909 {
910 //try protonating the nitrogen
911 atom.incrementImplicitValence();
912 i.setActualIndex(i.actualIndex() + 1);
913
914 if (expandKekule(mol, va, i, maxv, secondpass))
915 {
916 return (true);
917 }
918
919 i.setActualIndex(i.actualIndex() - 1);
920 atom.decrementImplicitValence();
921 }
922 }
923
924 //failed to find a valid solution - reset attached bonds
925 for (int j = 0; j < vb.size(); j++)
926 {
927 bond = (JOEBond) vb.get(j);
928 bond.setKSingle();
929 bond.setBO(5);
930 }
931
932 return false;
933 }
934
935 /**
936 * @param path {@link java.util.Vector} of <tt>int[1]</tt>
937 * @param mol Description of the Parameter
938 * @param avisit Description of the Parameter
939 * @param bvisit Description of the Parameter
940 * @param natom Description of the Parameter
941 * @param depth Description of the Parameter
942 */
943 public static void findRings(JOEMol mol, Vector path, JOEBitVec avisit,
944 JOEBitVec bvisit, int natom, int depth)
945 {
946 JOEAtom atom;
947 JOEBond bond;
948 int[] itmp;
949
950 if (avisit.get(natom))
951 {
952 int j = depth - 1;
953
954 itmp = (int[]) path.get(j--);
955 bond = mol.getBond(itmp[0]);
956
957 bond.setInRing();
958
959 while (j >= 0)
960 {
961 itmp = (int[]) path.get(j--);
962 bond = mol.getBond(itmp[0]);
963
964 bond.setInRing();
965 (bond.getBeginAtom()).setInRing();
966 (bond.getEndAtom()).setInRing();
967
968 if ((bond.getBeginAtomIdx() == natom) ||
969 (bond.getEndAtomIdx() == natom))
970 {
971 break;
972 }
973 }
974 }
975 else
976 {
977 avisit.setBitOn(natom);
978 atom = mol.getAtom(natom);
979
980 BondIterator bit = atom.bondIterator();
981
982 while (bit.hasNext())
983 {
984 bond = bit.nextBond();
985
986 if (!bvisit.get(bond.getIdx()))
987 {
988 itmp = (int[]) path.get(depth);
989 itmp[0] = bond.getIdx();
990 bvisit.setBitOn(bond.getIdx());
991
992 findRings(mol, path, avisit, bvisit,
993 bond.getNbrAtomIdx(atom), depth + 1);
994 }
995 }
996 }
997 }
998
999 /*-------------------------------------------------------------------------*
1000 * public methods
1001 *------------------------------------------------------------------------- */
1002
1003 /**
1004 * Description of the Method
1005 *
1006 * @return Description of the Return Value
1007 */
1008 public static synchronized JHM instance()
1009 {
1010 if (jhm == null)
1011 {
1012 jhm = new JHM();
1013 }
1014
1015 return jhm;
1016 }
1017
1018 /**
1019 * @param visit {@link java.util.Vector} of <tt>int[0]</tt>
1020 * @param ival {@link java.util.Vector} of <tt>int[0]</tt>
1021 * @param atom Description of the Parameter
1022 * @param depth Description of the Parameter
1023 * @return Description of the Return Value
1024 */
1025 public static boolean kekulePropagate(JOEAtom atom, Vector visit,
1026 Vector ival, int depth)
1027 {
1028 int count = 0;
1029 JOEBond bond;
1030 BondIterator bit = atom.bondIterator();
1031
1032 while (bit.hasNext())
1033 {
1034 bond = bit.nextBond();
1035
1036 if (((int[]) visit.get(bond.getIdx()))[0] == 0)
1037 {
1038 count++;
1039 }
1040 }
1041
1042 if (count == 0)
1043 {
1044 return (valenceSum(atom) == ((int[]) ival.get(atom.getIdx()))[0]);
1045 }
1046
1047 boolean result = true;
1048
1049 JOEAtom nbr;
1050 NbrAtomIterator nait = atom.nbrAtomIterator();
1051
1052 if (valenceSum(atom) >= ((int[]) ival.get(atom.getIdx()))[0])
1053 {
1054 while (nait.hasNext())
1055 {
1056 nbr = nait.nextNbrAtom();
1057 bond = nait.actualBond();
1058
1059 if (nbr.isAromatic() &&
1060 (((int[]) visit.get(bond.getIdx()))[0] == 0))
1061 {
1062 ((int[]) visit.get(bond.getIdx()))[0] = depth;
1063 bond.setKSingle();
1064 result = kekulePropagate(nbr, visit, ival, depth);
1065
1066 if (result)
1067 {
1068 break;
1069 }
1070
1071 // if (!result) break;
1072 }
1073 }
1074 }
1075 else if (count == 1)
1076 {
1077 while (nait.hasNext())
1078 {
1079 nbr = nait.nextNbrAtom();
1080 bond = nait.actualBond();
1081
1082 if (nbr.isAromatic() &&
1083 (((int[]) visit.get(bond.getIdx()))[0] == 0))
1084 {
1085 ((int[]) visit.get(bond.getIdx()))[0] = depth;
1086 bond.setKDouble();
1087 result = kekulePropagate(nbr, visit, ival, depth);
1088
1089 //break;
1090 if (result)
1091 {
1092 break;
1093 }
1094 }
1095 }
1096 }
1097
1098 return result;
1099 }
1100
1101 /**
1102 * Description of the Method
1103 *
1104 * @param mol Description of the Parameter
1105 * @param a Description of the Parameter
1106 * @param b Description of the Parameter
1107 * @param one2one Description of the Parameter
1108 * @return Description of the Return Value
1109 */
1110 public static double minimumPairRMS(JOEMol mol, double[] a, double[] b,
1111 boolean[] one2one)
1112 {
1113 int i;
1114 int j;
1115 int k = 0;
1116 double min;
1117 double tmp;
1118 double d_2 = 0.0;
1119 JOEBitVec bset = new JOEBitVec();
1120 one2one[0] = true;
1121
1122 Vector _atom = new Vector();
1123
1124 // of type JOEAtom
1125 _atom.setSize(mol.numAtoms());
1126
1127 for (i = 0; i < mol.numAtoms(); i++)
1128 {
1129 _atom.set(i, mol.getAtom(i + 1));
1130 }
1131
1132 for (i = 0; i < mol.numAtoms(); i++)
1133 {
1134 min = Double.MAX_VALUE;
1135
1136 for (j = 0; j < mol.numAtoms(); j++)
1137 {
1138 if ((((JOEAtom) _atom.get(i)).getAtomicNum() == ((JOEAtom) _atom.get(
1139 j)).getAtomicNum()) &&
1140 (((JOEAtom) _atom.get(i)).getHyb() == ((JOEAtom) _atom.get(
1141 j)).getHyb()))
1142 {
1143 if (!bset.get(j))
1144 {
1145 double tmp1 = a[3 * i] - b[3 * j];
1146 double tmp2 = a[(3 * i) + 1] - b[(3 * j) + 1];
1147 double tmp3 = a[(3 * i) + 2] - b[(3 * j) + 2];
1148 tmp = (tmp1 * tmp1) + (tmp2 * tmp2) + (tmp3 * tmp3);
1149
1150 if (tmp < min)
1151 {
1152 k = j;
1153 min = tmp;
1154 }
1155 }
1156 }
1157 }
1158
1159 if (i != j)
1160 {
1161 one2one[0] = false;
1162 }
1163
1164 bset.setBitOn(k);
1165 d_2 += min;
1166 }
1167
1168 d_2 /= (double) mol.numAtoms();
1169
1170 return Math.sqrt(d_2);
1171 }
1172
1173 public static synchronized boolean moleculeToFile(String filename,
1174 JOEMol mol, IOType type, boolean writeDescriptors)
1175 {
1176 FileOutputStream out;
1177
1178 try
1179 {
1180 out = new FileOutputStream(filename);
1181 }
1182 catch (Exception ex)
1183 {
1184 ex.printStackTrace();
1185
1186 return false;
1187 }
1188
1189 MoleculeFileType writer = null;
1190
1191 try
1192 {
1193 writer = JOEFileFormat.getMolWriter(out, type);
1194
1195 if (!writer.writeable())
1196 {
1197 logger.warn(type.getRepresentation() + " is not writeable.");
1198
1199 return false;
1200 }
1201
1202 if (!writeDescriptors &&
1203 JOEHelper.hasInterface(writer, "PropertyWriter"))
1204 {
1205 ((PropertyWriter) writer).write(mol, null, false, null);
1206 }
1207 else
1208 {
1209 writer.write(mol, null);
1210 }
1211
1212 // close molecule writer
1213 writer.closeWriter();
1214 out.close();
1215 }
1216 catch (IOException ex)
1217 {
1218 ex.printStackTrace();
1219
1220 return false;
1221 }
1222 catch (MoleculeIOException ex)
1223 {
1224 ex.printStackTrace();
1225
1226 return false;
1227 }
1228 catch (Exception ex)
1229 {
1230 ex.printStackTrace();
1231
1232 return false;
1233 }
1234
1235 return true;
1236 }
1237
1238 /**
1239 * @param visit {@link java.util.Vector} of <tt>int[0]</tt>
1240 * @param mol Description of the Parameter
1241 * @param depth Description of the Parameter
1242 */
1243 public static void resetVisit(JOEMol mol, Vector visit, int depth)
1244 {
1245 JOEBond bond;
1246 BondIterator bit = mol.bondIterator();
1247
1248 while (bit.hasNext())
1249 {
1250 bond = bit.nextBond();
1251
1252 if (bond.isAromatic() &&
1253 (((int[]) visit.get(bond.getIdx()))[0] >= depth))
1254 {
1255 ((int[]) visit.get(bond.getIdx()))[0] = 0;
1256 }
1257 }
1258 }
1259
1260 /**
1261 * @param a Description of the Parameter
1262 * @param b Description of the Parameter
1263 * @param c Description of the Parameter
1264 * @param d Description of the Parameter
1265 * @return Description of the Return Value
1266 */
1267
1268 //public void getChirality(JOEMol mol, Vector chirality)
1269 //{
1270 // int itmp[];
1271 //
1272 // chirality.setSize(mol.numAtoms()+1);
1273 // setIntVectorToValue(chirality, 0);
1274 //
1275 // JOEAtom atom;
1276 // AtomIterator ait = mol.atomIterator();
1277 // while(ait.hasNext())
1278 // {
1279 // atom = ait.nextAtom();
1280 // if (atom.isChiral())
1281 // {
1282 // double sv = calcSignedVolume(mol,atom);
1283 // itmp = (int[])chirality.get(atom.getIdx()-1);
1284 // if (sv < 0.0f) itmp[0] = -1;
1285 // else if (sv > 0.0) itmp[0] = 1;
1286 // }
1287 // }
1288 //}
1289
1290 /**
1291 * Calculate the signed volume for an atom. If the atom has a valence of 3
1292 * the coordinates of an attached hydrogen are calculated
1293 *
1294 * @param a Description of the Parameter
1295 * @param b Description of the Parameter
1296 * @param c Description of the Parameter
1297 * @param d Description of the Parameter
1298 * @return Description of the Return Value
1299 */
1300
1301 //public static double calcSignedVolume(JOEMol mol, JOEAtom atm)
1302 //{
1303 // XYZVector tmp_crd;
1304 // Vector nbr_atms; // of type int[1]
1305 // Vector nbr_crds; // of type XYZVector
1306 // double hbrad = JOEElementTable.instance().correctedBondRad(1,0);
1307 //
1308 // if (atm.getHvyValence() < 3)
1309 // {
1310 // logger.error("Cannot calculate a signed volume for an atom with a heavy atom valence of "+atm.getHvyValence());
1311 // System.exit(0);
1312 // }
1313 //
1314 // // Create a vector with the coordinates of the neighbor atoms
1315 // JOEAtom nbr;
1316 // NbrAtomIterator nait = atm.nbrAtomIterator();
1317 // while(nait.hasNext())
1318 // {
1319 // nbr = nait.nextNbrAtom();
1320 // nbr_atms.add(new int[]{nbr.getIdx()});
1321 // }
1322 //
1323 // // sort the neighbor atoms to insure a consistent ordering
1324 // //QuickInsertSort sorting = new QuickInsertSort();
1325 // //RingSizeComparator ringSizeComparator = new RingSizeComparator();
1326 // //sorting.sort(_rlist, ringSizeComparator);
1327 // sort(nbr_atms.begin(),nbr_atms.end());
1328 // for (int i = 0; i < nbr_atms.size(); i++)
1329 // {
1330 // JOEAtom tmp_atm = mol.getAtom( ((int[])nbr_atms.get(i))[0] );
1331 // nbr_crds.add(tmp_atm.getVector());
1332 // }
1333 //
1334 // // If we have three heavy atoms we need to calculate the position of the fourth
1335 // if (atm.getHvyValence() == 3)
1336 // {
1337 // double bondlen = hbrad+JOEElementTable.instance().correctedBondRad(atm.getAtomicNum(), atm.getHyb());
1338 // atm.getNewBondVector(tmp_crd,bondlen);
1339 // nbr_crds.add(tmp_crd);
1340 // }
1341 //
1342 // return signed_volume( (XYZVector) nbr_crds.get(0),
1343 // (XYZVector) nbr_crds.get(1),
1344 // (XYZVector) nbr_crds.get(2),
1345 // (XYZVector) nbr_crds.get(3));
1346 //}
1347
1348 /**
1349 * Calculate the signed volume for an atom. If the atom has a valence of 3
1350 * the coordinates of an attached hydrogen are calculated
1351 *
1352 * Calculate a signed volume given a set of 4 coordinates.
1353 *
1354 * @param a Description of the Parameter
1355 * @param b Description of the Parameter
1356 * @param c Description of the Parameter
1357 * @param d Description of the Parameter
1358 * @return Description of the Return Value
1359 */
1360 public static double signedVolume(final XYZVector a, final XYZVector b,
1361 final XYZVector c, final XYZVector d)
1362 {
1363 XYZVector A;
1364 XYZVector B;
1365 XYZVector C;
1366 A = XYZVector.sub(b, a);
1367 B = XYZVector.sub(c, a);
1368 C = XYZVector.sub(d, a);
1369
1370 Matrix3x3 m = new Matrix3x3(A, B, C);
1371
1372 return m.determinant();
1373 }
1374
1375 /**
1376 * Description of the Method
1377 *
1378 * @param atom Description of the Parameter
1379 * @return Description of the Return Value
1380 */
1381 public static int valenceSum(JOEAtom atom)
1382 {
1383 int count = atom.getImplicitValence();
1384
1385 JOEBond bond;
1386 BondIterator bit = atom.bondIterator();
1387
1388 while (bit.hasNext())
1389 {
1390 bond = bit.nextBond();
1391
1392 if (bond.isKDouble())
1393 {
1394 count++;
1395 }
1396 }
1397
1398 return count;
1399 }
1400
1401 /**
1402 * Description of the Method
1403 *
1404 * @param a Description of the Parameter
1405 * @param b Description of the Parameter
1406 * @return Description of the Return Value
1407 */
1408 public boolean compareBonds(final JOEBond a, final JOEBond b)
1409 {
1410 if (a.getBeginAtomIdx() == b.getBeginAtomIdx())
1411 {
1412 return (a.getEndAtomIdx() < b.getEndAtomIdx());
1413 }
1414
1415 return (a.getBeginAtomIdx() < b.getBeginAtomIdx());
1416 }
1417
1418 /*void Toupper(string &s)
1419 *{
1420 *unsigned int i;
1421 *for (i = 0;i < s.size();i++)
1422 *s[i] = toupper(s[i]);
1423 *}
1424 *void Tolower(string &s)
1425 *{
1426 *unsigned int i;
1427 *for (i = 0;i < s.size();i++)
1428 *s[i] = tolower(s[i]);
1429 *} */
1430
1431 /**
1432 * Description of the Method
1433 *
1434 * @param mol Description of the Parameter
1435 * @return Description of the Return Value
1436 */
1437 public static int determineFRJ(JOEMol mol)
1438 {
1439 Vector cfl = new Vector();
1440
1441 // of type int[1]-Vector
1442 //find all continuous graphs in the mol area
1443 mol.contiguousFragments(cfl);
1444
1445 // System.out.println("cfl size: "+cfl.size());
1446 if (cfl.size() == 0)
1447 {
1448 return (0);
1449 }
1450
1451 if (cfl.size() == 1)
1452 {
1453 return (mol.numBonds() - mol.numAtoms() + 1);
1454 }
1455
1456 //count up the atoms and bonds belonging to each graph
1457 JOEBond bond = null;
1458
1459 int numatoms;
1460
1461 int numbonds;
1462
1463 int frj = 0;
1464 JOEBitVec frag = new JOEBitVec();
1465
1466 // Vector vtmp;
1467 int[] itmp;
1468
1469 // System.out.println("cfl size: "+cfl.size());
1470 for (int i = 0; i < cfl.size(); i++)
1471 {
1472 //vtmp = (Vector) cfl.get(i);
1473 itmp = (int[]) cfl.get(i);
1474
1475 // for (int ii = 0; ii < vtmp.size(); ii++)
1476 // {
1477 // System.out.print(" "+ ((int[])vtmp.get(ii))[0] );
1478 // }
1479 // System.out.println("");
1480 frag.clear();
1481
1482 //frag.fromVecInt(vtmp);
1483 //numatoms = vtmp.size();
1484 frag.fromIntArray(itmp);
1485 numatoms = itmp.length;
1486 numbonds = 0;
1487
1488 BondIterator bit = mol.bondIterator();
1489
1490 while (bit.hasNext())
1491 {
1492 bond = bit.nextBond();
1493
1494 if (frag.bitIsOn(bond.getBeginAtomIdx()) &&
1495 frag.bitIsOn(bond.getEndAtomIdx()))
1496 {
1497 numbonds++;
1498 }
1499 }
1500
1501 frj += (numbonds - numatoms + 1);
1502
1503 // if(bond!=null)System.out.println("frj: "+frj+" begidx "+bond.getBeginAtomIdx()+" endidx "+bond.getEndAtomIdx());
1504 // else System.out.println("not def");
1505 }
1506
1507 return (frj);
1508 }
1509
1510 /**
1511 * @param vcr {@link java.util.Vector} of <tt>String</tt>
1512 * @param buf Description of the Parameter
1513 * @return Description of the Return Value
1514 */
1515 public static boolean tokenize(Vector vcr, String buf)
1516 {
1517 return tokenize(vcr, buf, " \t\n");
1518 }
1519
1520 /**
1521 * @param vcr {@link java.util.Vector} of <tt>String</tt>
1522 * @param buf Description of the Parameter
1523 * @param delimstr Description of the Parameter
1524 * @return Description of the Return Value
1525 */
1526 public static boolean tokenize(Vector vcr, String buf, String delimstr)
1527 {
1528 vcr.clear();
1529 buf = new String(buf + "\n");
1530
1531 StringTokenizer st = new StringTokenizer(buf, delimstr);
1532
1533 while (st.hasMoreTokens())
1534 {
1535 vcr.add(st.nextToken());
1536 }
1537
1538 return true;
1539 }
1540
1541 /**
1542 * @param vcr {@link java.util.Vector} of <tt>String</tt>
1543 * @param s Description of the Parameter
1544 * @param delimstr Description of the Parameter
1545 * @param limit Description of the Parameter
1546 * @return Description of the Return Value
1547 */
1548 public static boolean tokenize(Vector vcr, String s, String delimstr,
1549 int limit)
1550 {
1551 System.out.println("Warning: tokenize \"" + s + "\"");
1552 vcr.clear();
1553 s = new String(s + "\n");
1554
1555 int endpos = 0;
1556 int matched = 0;
1557
1558 StringTokenizer st = new StringTokenizer(s, delimstr);
1559
1560 while (st.hasMoreTokens())
1561 {
1562 String tmp = st.nextToken();
1563 vcr.add(tmp);
1564
1565 matched++;
1566
1567 if (matched == limit)
1568 {
1569 endpos = s.lastIndexOf(tmp);
1570 vcr.add(s.substring(endpos + tmp.length()));
1571
1572 break;
1573 }
1574 }
1575
1576 return true;
1577 }
1578
1579 //public boolean readMolecule(istream ifs, JOEMol mol, final String filename, String title)
1580 //{
1581 // int i, len, exten;
1582 //
1583 // if (filename==null)
1584 // {
1585 // logger.error("ReadMolecule error: filename is null");
1586 // return false;
1587 // }
1588 //
1589 // exten = filename.indexOf(".");
1590 // final String ext = filename.substring(exten+1);
1591 //
1592 // if ( ext.equals("mol2") ) ) return readMol2(ifs,mol,title);
1593 // else if (ext.equals("sdf") ) return readSDFile(ifs,mol,title);
1594 // else if (ext.equals("pdb") )
1595 // {
1596 // boolean b=readPDB(ifs,mol,title);
1597 // if ( mol.getTitle().length() == 0) mol.setTitle(filename);
1598 // return b;
1599 // }
1600 // else if (ext.equals("box") ) return readBox(ifs,mol,title);
1601 // else if (ext.equals("smi") ) return readSmiles(ifs,mol,title);
1602 // else if (ext.equals("gpMM") || ext.equals("gpQM") ) return readGhemical(ifs,mol,title);
1603 // else
1604 // {
1605 // logger.error("ReadMolecule error: "+filename+" : file extension not recognized");
1606 // return false;
1607 // }
1608 //}
1609 //public boolean writeMolecule(ostream ofs, JOEMol mol, final String filename, String dimension)
1610 //{
1611 // int i, len, exten;
1612 //
1613 // if (filename == null)
1614 // {
1615 // logger.error("WriteMolecule error: filename is null");
1616 // return false;
1617 // }
1618 //
1619 // exten = filename.indexOf(".");
1620 // final String ext = filename.substring(exten+1);
1621 //
1622 // if (ext.equals("mol2") ) return WriteMol2(ofs,mol,dimension);
1623 // else if (ext.equals("sdf") ) return WriteSDFile(ofs,mol,dimension);
1624 // else if (ext.equals("pdb") ) return WriteDelphiPDB(ofs,mol);
1625 // else if (ext.equals("gpMM") ) return WriteGhemical(ofs,mol);
1626 // else
1627 // {
1628 // logger.error("WriteMolecule error: "+filename+" : file extension not recognized");
1629 // return false;
1630 // }
1631 //}
1632
1633 /**
1634 * Description of the Method
1635 *
1636 * @param ofs Description of the Parameter
1637 * @param mol Description of the Parameter
1638 * @return Description of the Return Value
1639 */
1640 public boolean writeTitles(OutputStream ofs, JOEMol mol)
1641 {
1642 PrintStream ps = new PrintStream(ofs);
1643 ps.println(mol.getTitle());
1644
1645 return true;
1646 }
1647
1648 double calcRMS(double[] r, double[] f, int size)
1649 {
1650 int i;
1651 float d2 = 0.0f;
1652
1653 for (i = 0; i < size; i++)
1654 {
1655 d2 += (((r[i * 3] - f[i * 3]) * (r[i * 3] - f[i * 3])) +
1656 ((r[(i * 3) + 1] - f[(i * 3) + 1]) * (r[(i * 3) + 1] -
1657 f[(i * 3) + 1])) +
1658 ((r[(i * 3) + 2] - f[(i * 3) + 2]) * (r[(i * 3) + 2] -
1659 f[(i * 3) + 2])));
1660 }
1661
1662 d2 /= (double) size;
1663
1664 return (Math.sqrt(d2));
1665 }
1666
1667 /**
1668 *
1669 * @param vic of type <tt>JOEInternalCoord</tt>
1670 */
1671 void cartesianToInternal(Vector vic, JOEMol mol)
1672 {
1673 double r;
1674 double sum;
1675 JOEAtom atom;
1676 JOEAtom nbr;
1677 JOEAtom ref;
1678
1679 // vector<OEAtom*>::iterator i,j,m;
1680 //set reference atoms
1681 AtomIterator ait = mol.atomIterator();
1682 AtomIterator aitNbr = mol.atomIterator();
1683 JOEInternalCoord intCoord;
1684 JOEInternalCoord intCoordNbr;
1685
1686 while (ait.hasNext())
1687 {
1688 atom = ait.nextAtom();
1689 intCoord = (JOEInternalCoord) vic.get(atom.getIdx());
1690
1691 if (atom.getIdx() == 1)
1692 {
1693 continue;
1694 }
1695 else if (atom.getIdx() == 2)
1696 {
1697 intCoord._a = mol.getAtom(1);
1698
1699 continue;
1700 }
1701 else if (atom.getIdx() == 3)
1702 {
1703 intCoord._a = mol.getAtom(2);
1704 intCoord._b = mol.getAtom(1);
1705
1706 continue;
1707 }
1708
1709 sum = 1.0E10;
1710 ref = mol.getAtom(1);
1711 aitNbr.reset();
1712
1713 while (aitNbr.hasNext() &&
1714