File: Renumber.cpp

package info (click to toggle)
rdkit 201809.1%2Bdfsg-6
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 123,688 kB
  • sloc: cpp: 230,509; python: 70,501; java: 6,329; ansic: 5,427; sql: 1,899; yacc: 1,739; lex: 1,243; makefile: 445; xml: 229; fortran: 183; sh: 123; cs: 93
file content (104 lines) | stat: -rw-r--r-- 3,343 bytes parent folder | download | duplicates (2)
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
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
//
//  Copyright (C) 2013 Greg Landrum
//
//   @@ All Rights Reserved @@
//  This file is part of the RDKit.
//  The contents are covered by the terms of the BSD license
//  which is included in the file license.txt, found at the root
//  of the RDKit source tree.
//
#include <GraphMol/GraphMol.h>
#include <GraphMol/MolOps.h>
#include <RDGeneral/Exceptions.h>
#include <GraphMol/AtomIterators.h>
#include <GraphMol/BondIterators.h>
#include <boost/foreach.hpp>

namespace RDKit {
namespace MolOps {
ROMol *renumberAtoms(const ROMol &mol,
                     const std::vector<unsigned int> &newOrder) {
  unsigned int nAts = mol.getNumAtoms();
  PRECONDITION(newOrder.size() == nAts, "bad newOrder size");

  std::vector<unsigned int> revOrder(nAts);
  for (unsigned int nIdx = 0; nIdx < nAts; ++nIdx) {
    unsigned int oIdx = newOrder[nIdx];
    if (oIdx > nAts) {
      throw ValueErrorException("idx value exceeds numAtoms");
    }
    revOrder[oIdx] = nIdx;
  }

  // ------
  // newOrder[i] : which atom should be in position i of the new mol
  // revOrder[i] : where atom i of the original mol landed in the new mol
  auto *res = new RWMol();

  // copy over the atoms:
  for (unsigned int nIdx = 0; nIdx < nAts; ++nIdx) {
    unsigned int oIdx = newOrder[nIdx];
    const Atom *oAtom = mol.getAtomWithIdx(oIdx);
    Atom *nAtom = oAtom->copy();
    res->addAtom(nAtom, false, true);

    // take care of atom-numbering-dependent properties:
    INT_VECT nAtoms;
    if (nAtom->getPropIfPresent(common_properties::_ringStereoAtoms, nAtoms)) {
      // FIX: ought to be able to avoid this copy.
      BOOST_FOREACH (int &val, nAtoms) {
        if (val < 0) {
          val = -1 * (revOrder[(-val - 1)] + 1);
        } else {
          val = revOrder[val - 1] + 1;
        }
      }
      nAtom->setProp(common_properties::_ringStereoAtoms, nAtoms, true);
    }
  }

  // now the bonds:
  for (ROMol::ConstBondIterator bi = mol.beginBonds(); bi != mol.endBonds();
       ++bi) {
    const Bond *oBond = (*bi);
    Bond *nBond = oBond->copy();
    nBond->setBeginAtomIdx(revOrder[oBond->getBeginAtomIdx()]);
    nBond->setEndAtomIdx(revOrder[oBond->getEndAtomIdx()]);
    res->addBond(nBond, true);
    // take care of atom-numbering-dependent properties:
    BOOST_FOREACH (int &idx, nBond->getStereoAtoms()) { idx = revOrder[idx]; }
  }

  // Conformers:
  for (auto oConf = mol.beginConformers(); oConf != mol.endConformers();
       ++oConf) {
    auto *nConf = new Conformer(nAts);
    for (unsigned int i = 0; i < nAts; ++i) {
      nConf->setAtomPos(i, (*oConf)->getAtomPos(newOrder[i]));
    }
    nConf->setId((*oConf)->getId());
    nConf->set3D((*oConf)->is3D());
    res->addConformer(nConf);
  }

  // update the ring info:
  const RingInfo *oRings = mol.getRingInfo();
  if (oRings && oRings->isInitialized() ) {
    RingInfo *nRings = res->getRingInfo();
    nRings->reset();
    nRings->initialize();
    for (unsigned int i = 0; i < oRings->numRings(); ++i) {
      const INT_VECT &oRing = oRings->atomRings()[i];
      INT_VECT nRing(oRing.size());
      for (unsigned int j = 0; j < oRing.size(); ++j) {
        nRing[j] = revOrder[oRing[j]];
      }
      nRings->addRing(nRing, oRings->bondRings()[i]);
    }
  }

  return dynamic_cast<ROMol *>(res);
}

};  // end of namespace MolOps
};  // end of namespace RDKit