C++ OBMol::GetSSSR方法代码示例

///////////////////////////////////////////////////////////////////////////////
//! \return the number of size of the set of smallest rings (SSSR)
int nrings(OBMol &mol)
{
  int nr;
  vector<OBRing*> vr;
  
  vr = mol.GetSSSR();
  nr = vr.size();
  return (nr);
}

  // Adapted from OBAtom::IsInRingSize()
  OBRing* OBBond::FindSmallestRing() const
  {
    vector<OBRing*> rlist;
    vector<OBRing*>::iterator i;

    OBMol *mol = (OBMol*)((OBBond*)this)->GetParent();

    rlist = mol->GetSSSR();
    OBRing* result = (OBRing*) NULL;
    size_t min_size = UINT_MAX;
    for (i = rlist.begin();i != rlist.end();++i) {
      if ((*i)->IsMember((OBBond*)this) && (*i)->Size() < min_size) {
        min_size = (*i)->Size();
        result = *i;
      }
    }
    return result;
  }

  void OBRingTyper::AssignTypes(OBMol &mol)
  {
    if (!_init)
      Init();

    obErrorLog.ThrowError(__FUNCTION__,
                          "Ran OBRing::AssignTypes", obAuditMsg);

    mol.SetRingTypesPerceived();

    vector<vector<int> >::iterator j2;
    vector<pair<OBSmartsPattern*,string> >::iterator i2;

    vector<OBRing*>::iterator i;
    vector<int>::iterator j;
    vector<OBRing*> rlist = mol.GetSSSR();

    int member_count;
    for (i2 = _ringtyp.begin();i2 != _ringtyp.end();++i2) { // for each ring type
      if (i2->first->Match(mol)) {
        _mlist = i2->first->GetMapList();
        for (j2 = _mlist.begin();j2 != _mlist.end();++j2) { // for each found match

          for (i = rlist.begin();i != rlist.end();++i) { // for each ring
            member_count = 0;
            
            for(j = j2->begin(); j != j2->end(); ++j) { // for each atom in the match
              if ((*i)->IsMember(mol.GetAtom(*j))) 
                member_count++;
            }
	    
            if ((*i)->Size() == member_count)
              (*i)->SetType(i2->second);
          }
        }
      }
    }
  }

void ring_test()
{
  ostringstream os;

#ifdef TESTDATADIR
  string testdatadir = TESTDATADIR;
  string results_file = testdatadir + "ringresults.txt";
  string smilestypes_file = testdatadir + "attype.00.smi";
#else
  string results_file = "files/ringresults.txt";
  string smilestypes_file = "files/attype.00.smi";
#endif

  cout << "# Testing ring perception..." << endl;

  std::ifstream mifs;
  os << "Bail out! Cannot read test data " << smilestypes_file.c_str();
  BOOST_REQUIRE_MESSAGE( SafeOpen(mifs, smilestypes_file.c_str()), os.str().c_str() );

  std::ifstream rifs;
  os.str("");
  os << "Bail out! Cannot read test data " << results_file.c_str();
  BOOST_REQUIRE_MESSAGE( SafeOpen(rifs, results_file.c_str()), os.str().c_str() );

  unsigned int size;
  OBBond *bond;
  OBAtom *atom;
  int count;
  char buffer[BUFF_SIZE];
  vector<string> vs;
  vector<OBRing*> vr;
  vector<bool> vb;
  vector<int> vi;
  OBMol mol;
  vector<string>::iterator i;
  vector<OBBond*>::iterator j;
  vector<OBAtom*>::iterator k;
  vector<OBRing*>::iterator m;
  OBConversion conv(&mifs, &cout);
  unsigned int currentTest = 0;

  BOOST_REQUIRE_MESSAGE( conv.SetInAndOutFormats("SMI","SMI"), "Bail out! SMILES format is not loaded" );

  for (;mifs;)
    {
      mol.Clear();
      conv.Read(&mol);
      if (mol.Empty())
        continue;
      
      BOOST_REQUIRE_MESSAGE( rifs.getline(buffer,BUFF_SIZE), "Bail out! error reading reference data" );

      vb.clear();
      vb.resize(mol.NumBonds(),false);
      //check ring bonds
      tokenize(vs,buffer);
      for (i = vs.begin();i != vs.end();i++)
        vb[atoi(i->c_str())] = true;

      for (bond = mol.BeginBond(j);bond;bond = mol.NextBond(j))
        {
          os.str("");
	  os << "ring bond data different than reference # Molecule: " << mol.GetTitle();
          BOOST_CHECK_MESSAGE( vb[bond->GetIdx()] == bond->IsInRing(), os.str().c_str() );
        }

      vr = mol.GetSSSR();
      BOOST_REQUIRE_MESSAGE( rifs.getline(buffer,BUFF_SIZE), "Bail out! error reading reference data" );
      sscanf(buffer,"%d",&size);

      os.str("");
      os << "SSSR size different than reference # Molecule: " << mol.GetTitle();
      BOOST_CHECK_MESSAGE( vr.size() == size, os.str().c_str() ); //check SSSR size

      BOOST_REQUIRE_MESSAGE( rifs.getline(buffer,BUFF_SIZE), "Bail out! error reading reference data" );

      tokenize(vs,buffer);
      i = vs.begin();
      for (atom = mol.BeginAtom(k);atom;atom = mol.NextAtom(k))
        {
          os.str("");
          os << "error in SSSR count # Molecule: " << mol.GetTitle();
          BOOST_CHECK_MESSAGE( i != vs.end(), os.str().c_str() ); //check SSSR size

          count = 0;
          for (m = vr.begin();m != vr.end();m++)
            if ((*m)->_pathset[atom->GetIdx()])
              count++;

	  os.str("");
          os << "ring membership test failed # Molecule: " << mol.GetTitle();
          BOOST_CHECK_MESSAGE( atoi(i->c_str()) == count, os.str().c_str() ); 

          i++;
        }
    }

}

  /**
   * \brief Select the root atoms for traversing atoms in rings.
   *
   * Picking only the begin atom of a closure bond can cause
   * difficulties when the selected atom is an inner atom
   * with three neighbour ring atoms. Why ? Because this atom
   * can get trapped by the other atoms when determining aromaticity,
   * because a simple visited flag is used in the
   * OBAromaticTyper::TraverseCycle() method.
   *
   * Ported from JOELib, copyright Joerg Wegner, 2003 under the GPL version 2
   * Improved by Fabian (fab5) in 2009 -- PR#2889708
   *
   * @param mol the molecule
   * @param avoidInnerRingAtoms inner closure ring atoms with more than 2 neighbours will be avoided
   *
   */
  void OBAromaticTyper::SelectRootAtoms(OBMol &mol, bool avoidInnerRingAtoms)
  {
    OBBond *bond;
    OBAtom *atom, *nbr, *nbr2;
    OBRing *ring;
    //    vector<OBAtom*>::iterator i;
    vector<OBBond*>::iterator j, l, nbr2Iter;
    vector<OBRing*> sssRings = mol.GetSSSR();
    vector<OBRing*>::iterator k;
    
    int rootAtom;
    int ringNbrs;
    int heavyNbrs;
    int newRoot = -1;
    vector<int> tmpRootAtoms;
    vector<int> tmp;

    vector<OBBond*> cbonds;
    vector< vector<OBRing*> > ringAtoms; // store ring pointers on an atom basis
    
    //generate list of closure bonds 
    for (bond = mol.BeginBond(j);bond;bond = mol.NextBond(j))
      {
        if( bond->IsClosure() )
          {
            cbonds.push_back(bond);
            if(avoidInnerRingAtoms)
              {
                tmpRootAtoms.push_back(bond->GetBeginAtomIdx());
              }
          }
      }
    
    if(avoidInnerRingAtoms)
      {
        //for every atom fill vector with ring pointer it's associated with
        ringAtoms.resize(mol.NumAtoms()+1);
        for (k = sssRings.begin();k != sssRings.end();++k)
          {
            tmp = (*k)->_path;
            for (unsigned int j (0),j_end(tmp.size()); j < j_end; ++j)
              {
                ringAtoms[tmp[j]].push_back(*k);
              }
          }
      }
    
    
    //loop over closure bonds
    for(OBBondIterator bd(cbonds.begin()),bd_end(cbonds.end());bd!=bd_end;++bd)
      {
        bond = *bd;

        // BASIC APPROACH
        // pick beginning atom at closure bond
        // this is really ready, isn't it ! ;-)
        rootAtom = bond->GetBeginAtomIdx();
        _root[rootAtom] = true;
	
        // EXTENDED APPROACH
        if (avoidInnerRingAtoms)
          {
            // count the number of neighbor ring atoms
            atom = mol.GetAtom(rootAtom);
            ringNbrs = heavyNbrs = 0;
	    
            for (nbr = atom->BeginNbrAtom(l);nbr;nbr = atom->NextNbrAtom(l))
              {
                // we can get this from atom->GetHvyValence()
                // but we need to find neighbors in rings too
                // so let's save some time
                if (!nbr->IsHydrogen())
                  {
                    heavyNbrs++;
                    if (nbr->IsInRing())
                      ringNbrs++;
                  }
		
                // if this atom has more than 2 neighbor ring atoms
                // we could get trapped later when traversing cycles
                // which can cause aromaticity false detection
                newRoot = -1;
		
//.........这里部分代码省略.........