RDKit
Open-source cheminformatics and machine learning.
new_canon.h
Go to the documentation of this file.
1 //
2 // Copyright (C) 2014 Greg Landrum
3 // Adapted from pseudo-code from Roger Sayle
4 //
5 // @@ All Rights Reserved @@
6 // This file is part of the RDKit.
7 // The contents are covered by the terms of the BSD license
8 // which is included in the file license.txt, found at the root
9 // of the RDKit source tree.
10 //
11 
12 #include <RDGeneral/export.h>
13 #include <RDGeneral/hanoiSort.h>
14 #include <GraphMol/ROMol.h>
15 #include <GraphMol/RingInfo.h>
17 #include <cstdint>
18 #include <boost/foreach.hpp>
19 #include <boost/dynamic_bitset.hpp>
21 #include <cstring>
22 #include <iostream>
23 #include <cassert>
24 #include <cstring>
25 #include <vector>
26 
27 //#define VERBOSE_CANON 1
28 
29 namespace RDKit {
30 namespace Canon {
31 
34  unsigned int bondStereo;
35  unsigned int nbrSymClass;
36  unsigned int nbrIdx;
38  : bondType(Bond::UNSPECIFIED),
39  bondStereo(static_cast<unsigned int>(Bond::STEREONONE)),
40  nbrSymClass(0),
41  nbrIdx(0){};
42  bondholder(Bond::BondType bt, Bond::BondStereo bs, unsigned int ni,
43  unsigned int nsc)
44  : bondType(bt),
45  bondStereo(static_cast<unsigned int>(bs)),
46  nbrSymClass(nsc),
47  nbrIdx(ni){};
48  bondholder(Bond::BondType bt, unsigned int bs, unsigned int ni,
49  unsigned int nsc)
50  : bondType(bt), bondStereo(bs), nbrSymClass(nsc), nbrIdx(ni){};
51  bool operator<(const bondholder &o) const {
52  if (bondType != o.bondType) return bondType < o.bondType;
53  if (bondStereo != o.bondStereo) return bondStereo < o.bondStereo;
54  return nbrSymClass < o.nbrSymClass;
55  }
56  static bool greater(const bondholder &lhs, const bondholder &rhs) {
57  if (lhs.bondType != rhs.bondType) return lhs.bondType > rhs.bondType;
58  if (lhs.bondStereo != rhs.bondStereo)
59  return lhs.bondStereo > rhs.bondStereo;
60  return lhs.nbrSymClass > rhs.nbrSymClass;
61  }
62 
63  static int compare(const bondholder &x, const bondholder &y,
64  unsigned int div = 1) {
65  if (x.bondType < y.bondType)
66  return -1;
67  else if (x.bondType > y.bondType)
68  return 1;
69  if (x.bondStereo < y.bondStereo)
70  return -1;
71  else if (x.bondStereo > y.bondStereo)
72  return 1;
73  return x.nbrSymClass / div - y.nbrSymClass / div;
74  }
75 };
76 
78  public:
79  const Atom *atom;
80  int index;
81  unsigned int degree;
82  unsigned int totalNumHs;
83  bool hasRingNbr;
85  int *nbrIds;
86  const std::string *p_symbol; // if provided, this is used to order atoms
87  std::vector<int> neighborNum;
88  std::vector<int> revistedNeighbors;
89  std::vector<bondholder> bonds;
90 
92  : atom(NULL),
93  index(-1),
94  degree(0),
95  totalNumHs(0),
96  hasRingNbr(false),
97  isRingStereoAtom(false),
98  nbrIds(NULL),
99  p_symbol(NULL){};
100 
102  free(nbrIds);
103  }
104 };
105 
107  canon_atom *atoms, std::vector<bondholder> &nbrs);
108 
110  canon_atom *atoms, std::vector<bondholder> &nbrs, unsigned int atomIdx,
111  std::vector<std::pair<unsigned int, unsigned int>> &result);
112 
113 /*
114  * Different types of atom compare functions:
115  *
116  * - SpecialChiralityAtomCompareFunctor: Allows canonizing molecules exhibiting
117  *dependent chirality
118  * - SpecialSymmetryAtomCompareFunctor: Very specialized, allows canonizing
119  *highly symmetrical graphs/molecules
120  * - AtomCompareFunctor: Basic atom compare function which also allows to
121  *include neighbors within the ranking
122  */
123 
125  public:
127  const ROMol *dp_mol;
128  const boost::dynamic_bitset<> *dp_atomsInPlay, *dp_bondsInPlay;
129 
131  : dp_atoms(NULL),
132  dp_mol(NULL),
133  dp_atomsInPlay(NULL),
134  dp_bondsInPlay(NULL){};
136  Canon::canon_atom *atoms, const ROMol &m,
137  const boost::dynamic_bitset<> *atomsInPlay = NULL,
138  const boost::dynamic_bitset<> *bondsInPlay = NULL)
139  : dp_atoms(atoms),
140  dp_mol(&m),
141  dp_atomsInPlay(atomsInPlay),
142  dp_bondsInPlay(bondsInPlay){};
143  int operator()(int i, int j) const {
144  PRECONDITION(dp_atoms, "no atoms");
145  PRECONDITION(dp_mol, "no molecule");
146  PRECONDITION(i != j, "bad call");
147  if (dp_atomsInPlay && !((*dp_atomsInPlay)[i] || (*dp_atomsInPlay)[j])) {
148  return 0;
149  }
150 
151  if ((dp_atomsInPlay && (*dp_atomsInPlay)[i]) || !dp_atomsInPlay) {
152  updateAtomNeighborIndex(dp_atoms, dp_atoms[i].bonds);
153  }
154  if ((dp_atomsInPlay && (*dp_atomsInPlay)[j]) || !dp_atomsInPlay) {
155  updateAtomNeighborIndex(dp_atoms, dp_atoms[j].bonds);
156  }
157  for (unsigned int ii = 0;
158  ii < dp_atoms[i].bonds.size() && ii < dp_atoms[j].bonds.size(); ++ii) {
159  int cmp =
160  bondholder::compare(dp_atoms[i].bonds[ii], dp_atoms[j].bonds[ii]);
161  if (cmp) return cmp;
162  }
163 
164  std::vector<std::pair<unsigned int, unsigned int>> swapsi;
165  std::vector<std::pair<unsigned int, unsigned int>> swapsj;
166  if ((dp_atomsInPlay && (*dp_atomsInPlay)[i]) || !dp_atomsInPlay) {
167  updateAtomNeighborNumSwaps(dp_atoms, dp_atoms[i].bonds, i, swapsi);
168  }
169  if ((dp_atomsInPlay && (*dp_atomsInPlay)[j]) || !dp_atomsInPlay) {
170  updateAtomNeighborNumSwaps(dp_atoms, dp_atoms[j].bonds, j, swapsj);
171  }
172  for (unsigned int ii = 0; ii < swapsi.size() && ii < swapsj.size(); ++ii) {
173  int cmp = swapsi[ii].second - swapsj[ii].second;
174  if (cmp) return cmp;
175  }
176  return 0;
177  }
178 };
179 
181  public:
183  const ROMol *dp_mol;
184  const boost::dynamic_bitset<> *dp_atomsInPlay, *dp_bondsInPlay;
185 
187  : dp_atoms(NULL),
188  dp_mol(NULL),
189  dp_atomsInPlay(NULL),
190  dp_bondsInPlay(NULL){};
192  Canon::canon_atom *atoms, const ROMol &m,
193  const boost::dynamic_bitset<> *atomsInPlay = NULL,
194  const boost::dynamic_bitset<> *bondsInPlay = NULL)
195  : dp_atoms(atoms),
196  dp_mol(&m),
197  dp_atomsInPlay(atomsInPlay),
198  dp_bondsInPlay(bondsInPlay){};
199  int operator()(int i, int j) const {
200  PRECONDITION(dp_atoms, "no atoms");
201  PRECONDITION(dp_mol, "no molecule");
202  PRECONDITION(i != j, "bad call");
203  if (dp_atomsInPlay && !((*dp_atomsInPlay)[i] || (*dp_atomsInPlay)[j])) {
204  return 0;
205  }
206 
207  if (dp_atoms[i].neighborNum < dp_atoms[j].neighborNum) {
208  return -1;
209  } else if (dp_atoms[i].neighborNum > dp_atoms[j].neighborNum) {
210  return 1;
211  }
212 
213  if (dp_atoms[i].revistedNeighbors < dp_atoms[j].revistedNeighbors) {
214  return -1;
215  } else if (dp_atoms[i].revistedNeighbors > dp_atoms[j].revistedNeighbors) {
216  return 1;
217  }
218 
219  if ((dp_atomsInPlay && (*dp_atomsInPlay)[i]) || !dp_atomsInPlay) {
220  updateAtomNeighborIndex(dp_atoms, dp_atoms[i].bonds);
221  }
222  if ((dp_atomsInPlay && (*dp_atomsInPlay)[j]) || !dp_atomsInPlay) {
223  updateAtomNeighborIndex(dp_atoms, dp_atoms[j].bonds);
224  }
225  for (unsigned int ii = 0;
226  ii < dp_atoms[i].bonds.size() && ii < dp_atoms[j].bonds.size(); ++ii) {
227  int cmp =
228  bondholder::compare(dp_atoms[i].bonds[ii], dp_atoms[j].bonds[ii]);
229  if (cmp) return cmp;
230  }
231 
232  if (dp_atoms[i].bonds.size() < dp_atoms[j].bonds.size()) {
233  return -1;
234  } else if (dp_atoms[i].bonds.size() > dp_atoms[j].bonds.size()) {
235  return 1;
236  }
237  return 0;
238  }
239 };
240 
242  unsigned int getAtomRingNbrCode(unsigned int i) const {
243  if (!dp_atoms[i].hasRingNbr) return 0;
244 
245  int *nbrs = dp_atoms[i].nbrIds;
246  unsigned int code = 0;
247  for (unsigned j = 0; j < dp_atoms[i].degree; ++j) {
248  if (dp_atoms[nbrs[j]].isRingStereoAtom) {
249  code += dp_atoms[nbrs[j]].index * 10000 + 1; // j;
250  }
251  }
252  return code;
253  }
254 
255  int basecomp(int i, int j) const {
256  PRECONDITION(dp_atoms, "no atoms");
257  unsigned int ivi, ivj;
258 
259  // always start with the current class:
260  ivi = dp_atoms[i].index;
261  ivj = dp_atoms[j].index;
262  if (ivi < ivj)
263  return -1;
264  else if (ivi > ivj)
265  return 1;
266 
267  // use the atom-mapping numbers if they were assigned
268  /* boost::any_cast FILED here:
269  std::string molAtomMapNumber_i="";
270  std::string molAtomMapNumber_j="";
271  */
272  int molAtomMapNumber_i = 0;
273  int molAtomMapNumber_j = 0;
274  dp_atoms[i].atom->getPropIfPresent(common_properties::molAtomMapNumber,
275  molAtomMapNumber_i);
276  dp_atoms[j].atom->getPropIfPresent(common_properties::molAtomMapNumber,
277  molAtomMapNumber_j);
278  if (molAtomMapNumber_i < molAtomMapNumber_j)
279  return -1;
280  else if (molAtomMapNumber_i > molAtomMapNumber_j)
281  return 1;
282 
283  // start by comparing degree
284  ivi = dp_atoms[i].degree;
285  ivj = dp_atoms[j].degree;
286  if (ivi < ivj)
287  return -1;
288  else if (ivi > ivj)
289  return 1;
290 
291  if (dp_atoms[i].p_symbol && dp_atoms[j].p_symbol) {
292  if (*(dp_atoms[i].p_symbol) < *(dp_atoms[j].p_symbol))
293  return -1;
294  else if (*(dp_atoms[i].p_symbol) > *(dp_atoms[j].p_symbol))
295  return 1;
296  else
297  return 0;
298  }
299 
300  // move onto atomic number
301  ivi = dp_atoms[i].atom->getAtomicNum();
302  ivj = dp_atoms[j].atom->getAtomicNum();
303  if (ivi < ivj)
304  return -1;
305  else if (ivi > ivj)
306  return 1;
307 
308  // isotopes if we're using them
309  if (df_useIsotopes) {
310  ivi = dp_atoms[i].atom->getIsotope();
311  ivj = dp_atoms[j].atom->getIsotope();
312  if (ivi < ivj)
313  return -1;
314  else if (ivi > ivj)
315  return 1;
316  }
317 
318  // nHs
319  ivi = dp_atoms[i].totalNumHs;
320  ivj = dp_atoms[j].totalNumHs;
321  if (ivi < ivj)
322  return -1;
323  else if (ivi > ivj)
324  return 1;
325 
326  // charge
327  ivi = dp_atoms[i].atom->getFormalCharge();
328  ivj = dp_atoms[j].atom->getFormalCharge();
329  if (ivi < ivj)
330  return -1;
331  else if (ivi > ivj)
332  return 1;
333 
334  // chirality if we're using it
335  if (df_useChirality) {
336  // first atom stereochem:
337  ivi = 0;
338  ivj = 0;
339  std::string cipCode;
340  if (dp_atoms[i].atom->getPropIfPresent(common_properties::_CIPCode,
341  cipCode)) {
342  ivi = cipCode == "R" ? 2 : 1;
343  }
344  if (dp_atoms[j].atom->getPropIfPresent(common_properties::_CIPCode,
345  cipCode)) {
346  ivj = cipCode == "R" ? 2 : 1;
347  }
348  if (ivi < ivj)
349  return -1;
350  else if (ivi > ivj)
351  return 1;
352 
353  // can't actually use values here, because they are arbitrary
354  ivi = dp_atoms[i].atom->getChiralTag() != 0;
355  ivj = dp_atoms[j].atom->getChiralTag() != 0;
356  if (ivi < ivj)
357  return -1;
358  else if (ivi > ivj)
359  return 1;
360  }
361 
362  if (df_useChiralityRings) {
363  // ring stereochemistry
364  ivi = getAtomRingNbrCode(i);
365  ivj = getAtomRingNbrCode(j);
366  if (ivi < ivj)
367  return -1;
368  else if (ivi > ivj)
369  return 1;
370  // bond stereo is taken care of in the neighborhood comparison
371  }
372  return 0;
373  }
374 
375  public:
377  const ROMol *dp_mol;
378  const boost::dynamic_bitset<> *dp_atomsInPlay, *dp_bondsInPlay;
383 
385  : dp_atoms(NULL),
386  dp_mol(NULL),
387  dp_atomsInPlay(NULL),
388  dp_bondsInPlay(NULL),
389  df_useNbrs(false),
390  df_useIsotopes(true),
391  df_useChirality(true),
392  df_useChiralityRings(true){};
394  const boost::dynamic_bitset<> *atomsInPlay = NULL,
395  const boost::dynamic_bitset<> *bondsInPlay = NULL)
396  : dp_atoms(atoms),
397  dp_mol(&m),
398  dp_atomsInPlay(atomsInPlay),
399  dp_bondsInPlay(bondsInPlay),
400  df_useNbrs(false),
401  df_useIsotopes(true),
402  df_useChirality(true),
403  df_useChiralityRings(true){};
404  int operator()(int i, int j) const {
405  PRECONDITION(dp_atoms, "no atoms");
406  PRECONDITION(dp_mol, "no molecule");
407  PRECONDITION(i != j, "bad call");
408  if (dp_atomsInPlay && !((*dp_atomsInPlay)[i] || (*dp_atomsInPlay)[j])) {
409  return 0;
410  }
411 
412  int v = basecomp(i, j);
413  if (v) {
414  return v;
415  }
416 
417  if (df_useNbrs) {
418  if ((dp_atomsInPlay && (*dp_atomsInPlay)[i]) || !dp_atomsInPlay) {
419  updateAtomNeighborIndex(dp_atoms, dp_atoms[i].bonds);
420  }
421  if ((dp_atomsInPlay && (*dp_atomsInPlay)[j]) || !dp_atomsInPlay) {
422  updateAtomNeighborIndex(dp_atoms, dp_atoms[j].bonds);
423  }
424 
425  for (unsigned int ii = 0;
426  ii < dp_atoms[i].bonds.size() && ii < dp_atoms[j].bonds.size();
427  ++ii) {
428  int cmp =
429  bondholder::compare(dp_atoms[i].bonds[ii], dp_atoms[j].bonds[ii]);
430  if (cmp) return cmp;
431  }
432 
433  if (dp_atoms[i].bonds.size() < dp_atoms[j].bonds.size()) {
434  return -1;
435  } else if (dp_atoms[i].bonds.size() > dp_atoms[j].bonds.size()) {
436  return 1;
437  }
438  }
439  return 0;
440  }
441 };
442 
443 /*
444  * A compare function to discriminate chiral atoms, similar to the CIP rules.
445  * This functionality is currently not used.
446  */
447 
448 const unsigned int ATNUM_CLASS_OFFSET = 10000;
450  void getAtomNeighborhood(std::vector<bondholder> &nbrs) const {
451  for (unsigned j = 0; j < nbrs.size(); ++j) {
452  unsigned int nbrIdx = nbrs[j].nbrIdx;
453  if (nbrIdx == ATNUM_CLASS_OFFSET) {
454  // Ignore the Hs
455  continue;
456  }
457  const Atom *nbr = dp_atoms[nbrIdx].atom;
458  nbrs[j].nbrSymClass =
459  nbr->getAtomicNum() * ATNUM_CLASS_OFFSET + dp_atoms[nbrIdx].index + 1;
460  }
461  std::sort(nbrs.begin(), nbrs.end(), bondholder::greater);
462  // FIX: don't want to be doing this long-term
463  }
464 
465  int basecomp(int i, int j) const {
466  PRECONDITION(dp_atoms, "no atoms");
467  unsigned int ivi, ivj;
468 
469  // always start with the current class:
470  ivi = dp_atoms[i].index;
471  ivj = dp_atoms[j].index;
472  if (ivi < ivj)
473  return -1;
474  else if (ivi > ivj)
475  return 1;
476 
477  // move onto atomic number
478  ivi = dp_atoms[i].atom->getAtomicNum();
479  ivj = dp_atoms[j].atom->getAtomicNum();
480  if (ivi < ivj)
481  return -1;
482  else if (ivi > ivj)
483  return 1;
484 
485  // isotopes:
486  ivi = dp_atoms[i].atom->getIsotope();
487  ivj = dp_atoms[j].atom->getIsotope();
488  if (ivi < ivj)
489  return -1;
490  else if (ivi > ivj)
491  return 1;
492 
493  // atom stereochem:
494  ivi = 0;
495  ivj = 0;
496  std::string cipCode;
497  if (dp_atoms[i].atom->getPropIfPresent(common_properties::_CIPCode,
498  cipCode)) {
499  ivi = cipCode == "R" ? 2 : 1;
500  }
501  if (dp_atoms[j].atom->getPropIfPresent(common_properties::_CIPCode,
502  cipCode)) {
503  ivj = cipCode == "R" ? 2 : 1;
504  }
505  if (ivi < ivj)
506  return -1;
507  else if (ivi > ivj)
508  return 1;
509 
510  // bond stereo is taken care of in the neighborhood comparison
511  return 0;
512  }
513 
514  public:
516  const ROMol *dp_mol;
519  : dp_atoms(NULL), dp_mol(NULL), df_useNbrs(false){};
521  : dp_atoms(atoms), dp_mol(&m), df_useNbrs(false){};
522  int operator()(int i, int j) const {
523  PRECONDITION(dp_atoms, "no atoms");
524  PRECONDITION(dp_mol, "no molecule");
525  PRECONDITION(i != j, "bad call");
526  int v = basecomp(i, j);
527  if (v) return v;
528 
529  if (df_useNbrs) {
530  getAtomNeighborhood(dp_atoms[i].bonds);
531  getAtomNeighborhood(dp_atoms[j].bonds);
532 
533  // we do two passes through the neighbor lists. The first just uses the
534  // atomic numbers (by passing the optional 10000 to bondholder::compare),
535  // the second takes the already-computed index into account
536  for (unsigned int ii = 0;
537  ii < dp_atoms[i].bonds.size() && ii < dp_atoms[j].bonds.size();
538  ++ii) {
539  int cmp = bondholder::compare(
540  dp_atoms[i].bonds[ii], dp_atoms[j].bonds[ii], ATNUM_CLASS_OFFSET);
541  if (cmp) return cmp;
542  }
543  for (unsigned int ii = 0;
544  ii < dp_atoms[i].bonds.size() && ii < dp_atoms[j].bonds.size();
545  ++ii) {
546  int cmp =
547  bondholder::compare(dp_atoms[i].bonds[ii], dp_atoms[j].bonds[ii]);
548  if (cmp) return cmp;
549  }
550  if (dp_atoms[i].bonds.size() < dp_atoms[j].bonds.size()) {
551  return -1;
552  } else if (dp_atoms[i].bonds.size() > dp_atoms[j].bonds.size()) {
553  return 1;
554  }
555  }
556  return 0;
557  }
558 };
559 
560 /*
561  * Basic canonicalization function to organize the partitions which will be
562  * sorted next.
563  * */
564 
565 template <typename CompareFunc>
566 void RefinePartitions(const ROMol &mol, canon_atom *atoms, CompareFunc compar,
567  int mode, int *order, int *count, int &activeset,
568  int *next, int *changed, char *touchedPartitions) {
569  unsigned int nAtoms = mol.getNumAtoms();
570  int partition;
571  int symclass = 0;
572  int *start;
573  int offset;
574  int index;
575  int len;
576  int i;
577  // std::vector<char> touchedPartitions(mol.getNumAtoms(),0);
578 
579  // std::cerr<<"&&&&&&&&&&&&&&&& RP"<<std::endl;
580  while (activeset != -1) {
581  // std::cerr<<"ITER: "<<activeset<<" next: "<<next[activeset]<<std::endl;
582  // std::cerr<<" next: ";
583  // for(unsigned int ii=0;ii<nAtoms;++ii){
584  // std::cerr<<ii<<":"<<next[ii]<<" ";
585  // }
586  // std::cerr<<std::endl;
587  // for(unsigned int ii=0;ii<nAtoms;++ii){
588  // std::cerr<<order[ii]<<" count: "<<count[order[ii]]<<" index:
589  // "<<atoms[order[ii]].index<<std::endl;
590  // }
591 
592  partition = activeset;
593  activeset = next[partition];
594  next[partition] = -2;
595 
596  len = count[partition];
597  offset = atoms[partition].index;
598  start = order + offset;
599  // std::cerr<<"\n\n**************************************************************"<<std::endl;
600  // std::cerr<<" sort - class:"<<atoms[partition].index<<" len: "<<len<<":";
601  // for(unsigned int ii=0;ii<len;++ii){
602  // std::cerr<<" "<<order[offset+ii]+1;
603  // }
604  // std::cerr<<std::endl;
605  // for(unsigned int ii=0;ii<nAtoms;++ii){
606  // std::cerr<<order[ii]+1<<" count: "<<count[order[ii]]<<" index:
607  // "<<atoms[order[ii]].index<<std::endl;
608  // }
609  hanoisort(start, len, count, changed, compar);
610  // std::cerr<<"*_*_*_*_*_*_*_*_*_*_*_*_*_*_*_*"<<std::endl;
611  // std::cerr<<" result:";
612  // for(unsigned int ii=0;ii<nAtoms;++ii){
613  // std::cerr<<order[ii]+1<<" count: "<<count[order[ii]]<<" index:
614  // "<<atoms[order[ii]].index<<std::endl;
615  // }
616  for (int k = 0; k < len; ++k) {
617  changed[start[k]] = 0;
618  }
619 
620  index = start[0];
621  // std::cerr<<" len:"<<len<<" index:"<<index<<"
622  // count:"<<count[index]<<std::endl;
623  for (i = count[index]; i < len; i++) {
624  index = start[i];
625  if (count[index]) symclass = offset + i;
626  atoms[index].index = symclass;
627  // std::cerr<<" "<<index+1<<"("<<symclass<<")";
628  // if(mode && (activeset<0 || count[index]>count[activeset]) ){
629  // activeset=index;
630  //}
631  for (unsigned j = 0; j < atoms[index].degree; ++j) {
632  changed[atoms[index].nbrIds[j]] = 1;
633  }
634  }
635  // std::cerr<<std::endl;
636 
637  if (mode) {
638  index = start[0];
639  for (i = count[index]; i < len; i++) {
640  index = start[i];
641  for (unsigned j = 0; j < atoms[index].degree; ++j) {
642  unsigned int nbor = atoms[index].nbrIds[j];
643  touchedPartitions[atoms[nbor].index] = 1;
644  }
645  }
646  for (unsigned int ii = 0; ii < nAtoms; ++ii) {
647  if (touchedPartitions[ii]) {
648  partition = order[ii];
649  if ((count[partition] > 1) && (next[partition] == -2)) {
650  next[partition] = activeset;
651  activeset = partition;
652  }
653  touchedPartitions[ii] = 0;
654  }
655  }
656  }
657  }
658 } // end of RefinePartitions()
659 
660 template <typename CompareFunc>
661 void BreakTies(const ROMol &mol, canon_atom *atoms, CompareFunc compar,
662  int mode, int *order, int *count, int &activeset, int *next,
663  int *changed, char *touchedPartitions) {
664  unsigned int nAtoms = mol.getNumAtoms();
665  int partition;
666  int offset;
667  int index;
668  int len;
669  int oldPart = 0;
670 
671  for (unsigned int i = 0; i < nAtoms; i++) {
672  partition = order[i];
673  oldPart = atoms[partition].index;
674  while (count[partition] > 1) {
675  len = count[partition];
676  offset = atoms[partition].index + len - 1;
677  index = order[offset];
678  atoms[index].index = offset;
679  count[partition] = len - 1;
680  count[index] = 1;
681 
682  // test for ions, water molecules with no
683  if (atoms[index].degree < 1) {
684  continue;
685  }
686  for (unsigned j = 0; j < atoms[index].degree; ++j) {
687  unsigned int nbor = atoms[index].nbrIds[j];
688  touchedPartitions[atoms[nbor].index] = 1;
689  changed[nbor] = 1;
690  }
691 
692  for (unsigned int ii = 0; ii < nAtoms; ++ii) {
693  if (touchedPartitions[ii]) {
694  int npart = order[ii];
695  if ((count[npart] > 1) && (next[npart] == -2)) {
696  next[npart] = activeset;
697  activeset = npart;
698  }
699  touchedPartitions[ii] = 0;
700  }
701  }
702  RefinePartitions(mol, atoms, compar, mode, order, count, activeset, next,
703  changed, touchedPartitions);
704  }
705  // not sure if this works each time
706  if (atoms[partition].index != oldPart) {
707  i -= 1;
708  }
709  }
710 } // end of BreakTies()
711 
712 RDKIT_GRAPHMOL_EXPORT void CreateSinglePartition(unsigned int nAtoms,
713  int *order, int *count,
714  canon_atom *atoms);
715 
716 RDKIT_GRAPHMOL_EXPORT void ActivatePartitions(unsigned int nAtoms, int *order,
717  int *count, int &activeset,
718  int *next, int *changed);
719 
721  std::vector<unsigned int> &res,
722  bool breakTies = true,
723  bool includeChirality = true,
724  bool includeIsotopes = true);
725 
727  const ROMol &mol, std::vector<unsigned int> &res,
728  const boost::dynamic_bitset<> &atomsInPlay,
729  const boost::dynamic_bitset<> &bondsInPlay,
730  const std::vector<std::string> *atomSymbols = NULL, bool breakTies = true,
731  bool includeChirality = true, bool includeIsotopes = true);
732 
734  std::vector<unsigned int> &res);
735 
737  std::vector<Canon::canon_atom> &atoms,
738  bool includeChirality = true);
739 
740 } // namespace Canon
741 } // namespace RDKit
bool operator<(const bondholder &o) const
Definition: new_canon.h:51
RDKIT_RDGENERAL_EXPORT const std::string molAtomMapNumber
unsigned int degree
Definition: new_canon.h:81
const boost::dynamic_bitset * dp_bondsInPlay
Definition: new_canon.h:128
int operator()(int i, int j) const
Definition: new_canon.h:522
RDKIT_GRAPHMOL_EXPORT void updateAtomNeighborNumSwaps(canon_atom *atoms, std::vector< bondholder > &nbrs, unsigned int atomIdx, std::vector< std::pair< unsigned int, unsigned int >> &result)
Bond::BondType bondType
Definition: new_canon.h:33
unsigned int getNumAtoms(bool onlyExplicit=1) const
returns our number of atoms
unsigned int nbrIdx
Definition: new_canon.h:36
std::vector< int > neighborNum
Definition: new_canon.h:87
Defines the primary molecule class ROMol as well as associated typedefs.
RDKIT_RDGENERAL_EXPORT const std::string _CIPCode
void hanoisort(int *base, int nel, int *count, int *changed, CompareFunc compar)
Definition: hanoiSort.h:141
BondStereo
the nature of the bond&#39;s stereochem (for cis/trans)
Definition: Bond.h:95
AtomCompareFunctor(Canon::canon_atom *atoms, const ROMol &m, const boost::dynamic_bitset<> *atomsInPlay=NULL, const boost::dynamic_bitset<> *bondsInPlay=NULL)
Definition: new_canon.h:393
bondholder(Bond::BondType bt, Bond::BondStereo bs, unsigned int ni, unsigned int nsc)
Definition: new_canon.h:42
BondType
the type of Bond
Definition: Bond.h:56
RDKIT_GRAPHMOL_EXPORT void updateAtomNeighborIndex(canon_atom *atoms, std::vector< bondholder > &nbrs)
int operator()(int i, int j) const
Definition: new_canon.h:404
const boost::dynamic_bitset * dp_bondsInPlay
Definition: new_canon.h:378
static bool greater(const bondholder &lhs, const bondholder &rhs)
Definition: new_canon.h:56
#define RDKIT_GRAPHMOL_EXPORT
Definition: export.h:307
ChiralAtomCompareFunctor(Canon::canon_atom *atoms, const ROMol &m)
Definition: new_canon.h:520
RDKIT_GRAPHMOL_EXPORT void rankMolAtoms(const ROMol &mol, std::vector< unsigned int > &res, bool breakTies=true, bool includeChirality=true, bool includeIsotopes=true)
const unsigned int ATNUM_CLASS_OFFSET
Definition: new_canon.h:448
RDKIT_GRAPHMOL_EXPORT void ActivatePartitions(unsigned int nAtoms, int *order, int *count, int &activeset, int *next, int *changed)
SpecialSymmetryAtomCompareFunctor(Canon::canon_atom *atoms, const ROMol &m, const boost::dynamic_bitset<> *atomsInPlay=NULL, const boost::dynamic_bitset<> *bondsInPlay=NULL)
Definition: new_canon.h:191
RDKIT_GRAPHMOL_EXPORT void CreateSinglePartition(unsigned int nAtoms, int *order, int *count, canon_atom *atoms)
std::vector< bondholder > bonds
Definition: new_canon.h:89
RDKIT_GRAPHMOL_EXPORT void rankFragmentAtoms(const ROMol &mol, std::vector< unsigned int > &res, const boost::dynamic_bitset<> &atomsInPlay, const boost::dynamic_bitset<> &bondsInPlay, const std::vector< std::string > *atomSymbols=NULL, bool breakTies=true, bool includeChirality=true, bool includeIsotopes=true)
Std stuff.
Definition: Atom.h:30
unsigned int nbrSymClass
Definition: new_canon.h:35
int getAtomicNum() const
returns our atomic number
Definition: Atom.h:115
void BreakTies(const ROMol &mol, canon_atom *atoms, CompareFunc compar, int mode, int *order, int *count, int &activeset, int *next, int *changed, char *touchedPartitions)
Definition: new_canon.h:661
const boost::dynamic_bitset * dp_bondsInPlay
Definition: new_canon.h:184
std::vector< int > revistedNeighbors
Definition: new_canon.h:88
class for representing a bond
Definition: Bond.h:47
static int compare(const bondholder &x, const bondholder &y, unsigned int div=1)
Definition: new_canon.h:63
void RefinePartitions(const ROMol &mol, canon_atom *atoms, CompareFunc compar, int mode, int *order, int *count, int &activeset, int *next, int *changed, char *touchedPartitions)
Definition: new_canon.h:566
#define PRECONDITION(expr, mess)
Definition: Invariant.h:108
RDKIT_GRAPHMOL_EXPORT void initCanonAtoms(const ROMol &mol, std::vector< Canon::canon_atom > &atoms, bool includeChirality=true)
unsigned int totalNumHs
Definition: new_canon.h:82
RDKIT_GRAPHMOL_EXPORT void chiralRankMolAtoms(const ROMol &mol, std::vector< unsigned int > &res)
bondholder(Bond::BondType bt, unsigned int bs, unsigned int ni, unsigned int nsc)
Definition: new_canon.h:48
SpecialChiralityAtomCompareFunctor(Canon::canon_atom *atoms, const ROMol &m, const boost::dynamic_bitset<> *atomsInPlay=NULL, const boost::dynamic_bitset<> *bondsInPlay=NULL)
Definition: new_canon.h:135
const std::string * p_symbol
Definition: new_canon.h:86
unsigned int bondStereo
Definition: new_canon.h:34
The class for representing atoms.
Definition: Atom.h:69
Canon::canon_atom * dp_atoms
Definition: new_canon.h:376