C++ OBBond::GetIdx方法代码示例

  vector<OBBond*> OBResidue::GetBonds(bool exterior) const
  {
    OBAtom         *atom;
    vector<OBBond*> bonds;
    OBBitVec        idxs;
    unsigned int    sz;

    sz = (unsigned int) _atoms.size();
    for ( unsigned int i = 0 ; i < sz ; ++i )
      {
        atom = _atoms[i];
        OBBond *bond;
        vector<OBBond*>::iterator b;
        for (bond = atom->BeginBond(b) ; bond ; bond = atom->NextBond(b))
          {
            if (!idxs.BitIsOn(bond->GetIdx()))
              {
                if (!exterior)
                  {
                    if (bond->GetNbrAtom(atom)->GetResidue() == this)
                      bonds.push_back(&(*bond));
                  }
                else
                  bonds.push_back(&(*bond));

                idxs.SetBitOn(bond->GetIdx());
              }
          }
      }

    return bonds;
  }

  bool OBBond::IsClosure()
  {
    OBMol *mol = (OBMol*)GetParent();
    if (!mol)
      return(false);
    if (mol->HasClosureBondsPerceived())
      return(HasFlag(OB_CLOSURE_BOND));

    mol->SetClosureBondsPerceived();

    obErrorLog.ThrowError(__FUNCTION__,
                          "Ran OpenBabel::PerceiveClosureBonds", obAuditMsg);

    OBBond *bond;
    OBAtom *atom,*nbr;
    OBBitVec uatoms,ubonds;
    vector<OBAtom*> curr,next;
    vector<OBAtom*>::iterator i;
    vector<OBBond*>::iterator j;

    uatoms.Resize(mol->NumAtoms()+1);
    ubonds.Resize(mol->NumAtoms()+1);

    for (;static_cast<unsigned int>(uatoms.CountBits()) < mol->NumAtoms();)
      {
        if (curr.empty())
          for (atom = mol->BeginAtom(i);atom;atom = mol->NextAtom(i))
            if (!uatoms[atom->GetIdx()])
              {
                uatoms |= atom->GetIdx();
                curr.push_back(atom);
                break;
              }

        for (;!curr.empty();)
          {
            for (i = curr.begin();i != curr.end();++i)
              for (nbr = ((OBAtom*)*i)->BeginNbrAtom(j);nbr;nbr = ((OBAtom*)*i)->NextNbrAtom(j))
                if (!uatoms[nbr->GetIdx()])
                  {
                    uatoms |= nbr->GetIdx();
                    ubonds |= (*j)->GetIdx();
                    next.push_back(nbr);
                  }

            curr = next;
            next.clear();
          }
      }

    for (bond = mol->BeginBond(j);bond;bond = mol->NextBond(j))
      if (!ubonds[bond->GetIdx()])
        bond->SetClosure();

    return(HasFlag(OB_CLOSURE_BOND));
  }

void GenerateRingReference()
{
  std::ifstream ifs;
  if (!SafeOpen(ifs,"attype.00.smi")) return;

  std::ofstream ofs;
  if (!SafeOpen(ofs,"ringresults.txt")) return;

  int count;
  OBAtom *atom;
  OBBond *bond;
  char buffer[BUFF_SIZE];
  vector<OBRing*> vr;
  vector<OBEdgeBase*>::iterator i;
  vector<OBNodeBase*>::iterator j;
  vector<OBRing*>::iterator k;
  OBMol mol(SMI,SMI);
  OBFileFormat ff;

  for (;ifs;)
    {
      mol.Clear();
      ff.ReadMolecule(ifs, mol);
      if (mol.Empty()) continue;
      
      //write out ring bonds
      for (bond = mol.BeginBond(i);bond;bond = mol.NextBond(i))
	if (bond->IsInRing())
	  {
	    sprintf(buffer,"%3d",bond->GetIdx());
	    ofs << buffer;
	  }
      ofs << endl;

      vr = mol.GetSSSR();
      //write the total number of rings
      ofs << vr.size() << endl;

      //write the number of rings that each atom is a member of
      for (atom = mol.BeginAtom(j);atom;atom = mol.NextAtom(j))
	{
	  count = 0;
	  for (k = vr.begin();k != vr.end();k++)
	    if ((*k)->_pathset[atom->GetIdx()])
	      count++;

	  sprintf(buffer,"%3d",count);
	  ofs << buffer;
	}
      ofs << endl;

    }

  ThrowError("Ring perception test results written successfully");
}

  static void FindRings(OBMol &mol,vector<int> &path,OBBitVec &avisit,
                        OBBitVec &bvisit, int natom,int depth )
  {
    OBAtom *atom;
    OBBond *bond;
    vector<OBBond*>::iterator k;

    // don't return if all atoms are visited
    // (For example, some atoms are in multiple rings!) -GRH
      
    if (avisit[natom])
      {
        int j = depth-1;
        bond=mol.GetBond(path[j--]);
        bond->SetInRing();
        while( j >= 0 )
          {
            bond=mol.GetBond(path[j--]);
            bond->SetInRing();
            (bond->GetBeginAtom())->SetInRing();
            (bond->GetEndAtom())->SetInRing();
            if(bond->GetBeginAtomIdx()==static_cast<unsigned int>(natom) || bond->
               GetEndAtomIdx()==static_cast<unsigned int>(natom))
              break;
          }
      }
    else
      {
        avisit.SetBitOn(natom);
        atom = mol.GetAtom(natom);
        for(bond = atom->BeginBond(k);bond;bond=atom->NextBond(k))
          if( !bvisit[bond->GetIdx()])
            {
              path[depth] = bond->GetIdx();
              bvisit.SetBitOn(bond->GetIdx());
              FindRings(mol,path,avisit,bvisit,bond->GetNbrAtomIdx(atom),
                        depth+1);
            }
      }
  }

bool WriteAlchemy(ostream &ofs,OBMol &mol)
{
  unsigned int i;
  char buffer[BUFF_SIZE];
  char bond_string[10];
  
  snprintf(buffer, BUFF_SIZE, "%5d ATOMS, %5d BONDS,     0 CHARGES",
	  mol.NumAtoms(),
	  mol.NumBonds());
  ofs << buffer << endl;
  ttab.SetFromType("INT"); ttab.SetToType("ALC");

  OBAtom *atom;
  string str,str1;
  for(i = 1;i <= mol.NumAtoms(); i++)
  {
    atom = mol.GetAtom(i);
    str = atom->GetType();
    ttab.Translate(str1,str);
    snprintf(buffer, BUFF_SIZE, "%5d %-6s%8.4f %8.4f %8.4f     0.0000",
	    i,
	    (char*)str1.c_str(),
	    atom->GetX(),
	    atom->GetY(),
	    atom->GetZ());
    ofs << buffer << endl;
  }

  OBBond *bond;
  vector<OBEdgeBase*>::iterator j;

  for (bond = mol.BeginBond(j);bond;bond = mol.NextBond(j))
  {
    switch(bond->GetBO())
    {
    case 1 :  strcpy(bond_string,"SINGLE");  break;
    case 2 :  strcpy(bond_string,"DOUBLE");  break;
    case 3 :  strcpy(bond_string,"TRIPLE");  break;
    case 5 :  strcpy(bond_string,"AROMATIC"); break;
    default : strcpy(bond_string,"SINGLE");
    }
    snprintf(buffer, BUFF_SIZE, "%5d  %4d  %4d  %s",
	    bond->GetIdx()+1,
	    bond->GetBeginAtomIdx(),
	    bond->GetEndAtomIdx(),
	    bond_string);
    ofs << buffer << endl;
  }
  return(true);
}

void testOBRotorGetSet()
{
  OBBond bond;
  OBRotor rotor;
  
  // SetBond/GetBond
  rotor.SetBond(&bond);
  OB_ASSERT(rotor.GetBond()->GetIdx() == bond.GetIdx());
  // SetIdx/GetIdx
  rotor.SetIdx(45);
  OB_ASSERT(rotor.GetIdx() == 45);
  // SetDihedralAtoms/GetDihedralAtoms
  int ref_ptr[4] = {1, 2, 3, 4};
  rotor.SetDihedralAtoms(ref_ptr);
  rotor.GetDihedralAtoms(ref_ptr);
  OB_ASSERT(ref_ptr[0] == 1);
  OB_ASSERT(ref_ptr[1] == 2);
  OB_ASSERT(ref_ptr[2] == 3);
  OB_ASSERT(ref_ptr[3] == 4);
  std::vector<int> ref_vector = rotor.GetDihedralAtoms();
  OB_ASSERT(ref_vector[0] == 1);
  OB_ASSERT(ref_vector[1] == 2);
  OB_ASSERT(ref_vector[2] == 3);
  OB_ASSERT(ref_vector[3] == 4);
  rotor.SetDihedralAtoms(ref_vector);
  ref_vector = rotor.GetDihedralAtoms();
  OB_ASSERT(ref_vector[0] == 1);
  OB_ASSERT(ref_vector[1] == 2);
  OB_ASSERT(ref_vector[2] == 3);
  OB_ASSERT(ref_vector[3] == 4);
  // SetTorsionValues/GetTorsionValues/Size
  std::vector<double> angles;
  angles.push_back(0.0);
  angles.push_back(3.1415);
  rotor.SetTorsionValues(angles);
  OB_ASSERT(rotor.GetTorsionValues().size() == 2);
  OB_ASSERT(rotor.Size() == 2);
}

 /* A recursive O(N) traversal of the molecule */
 static int FindRings(OBAtom *atom, int *avisit, unsigned char *bvisit,
                      unsigned int &frj, int depth)
 {
   OBBond *bond;
   int result = -1;
   vector<OBBond*>::iterator k;
   for(bond = atom->BeginBond(k);bond;bond=atom->NextBond(k)) {
     unsigned int bidx = bond->GetIdx();
     if (bvisit[bidx] == 0) {
       bvisit[bidx] = 1;
       OBAtom *nbor = bond->GetNbrAtom(atom);
       unsigned int nidx = nbor->GetIdx();
       int nvisit = avisit[nidx];
       if (nvisit == 0) {
         avisit[nidx] = depth+1;
         nvisit = FindRings(nbor,avisit,bvisit,frj,depth+1);
         if (nvisit > 0) {
           if (nvisit <= depth) {
             bond->SetInRing();
             if (result < 0 || nvisit < result)
               result = nvisit;
           }
         }
       } else {
         if (result < 0 || nvisit < result)
           result = nvisit;
         bond->SetClosure();
         bond->SetInRing();
         frj++;
       }
     }
   }
   if (result > 0 && result <= depth)
     atom->SetInRing();
   return result;
 }

  bool MOL2Format::WriteMolecule(OBBase* pOb, OBConversion* pConv)
  {
    OBMol* pmol = dynamic_cast<OBMol*>(pOb);
    if(pmol==NULL)
      return false;

    //Define some references so we can use the old parameter names
    ostream &ofs = *pConv->GetOutStream();
    OBMol &mol = *pmol;
    bool ligandsOnly = pConv->IsOption("l", OBConversion::OUTOPTIONS)!=NULL;
 
    //The old code follows....
    string str,str1;
    char buffer[BUFF_SIZE],label[BUFF_SIZE];
    char rnum[BUFF_SIZE],rlabel[BUFF_SIZE];

    ofs << "@<TRIPOS>MOLECULE" << endl;
    str = mol.GetTitle();
    if (str.empty())
      ofs << "*****" << endl;
    else
      ofs << str << endl;

    snprintf(buffer, BUFF_SIZE," %d %d 0 0 0", mol.NumAtoms(),mol.NumBonds());
    ofs << buffer << endl;
    ofs << "SMALL" << endl;

    OBPairData *dp = (OBPairData*)mol.GetData("PartialCharges");
    if (dp != NULL) {
        // Tripos spec says:
        // NO_CHARGES, DEL_RE, GASTEIGER, GAST_HUCK, HUCKEL, PULLMAN, 
        // GAUSS80_CHARGES, AMPAC_CHARGES, MULLIKEN_CHARGES, DICT_ CHARGES,
        // MMFF94_CHARGES, USER_CHARGES
      if (dp->GetValue() == "Mulliken")
        ofs << "MULLIKEN_CHARGES" << endl;
      else // should pick from the Tripos types
        ofs << "GASTEIGER" << endl;
    }
    else { // No idea what these charges are... all our code sets "PartialCharges"
        ofs << "GASTEIGER" << endl;
    }

    ofs << "Energy = " << mol.GetEnergy() << endl;

    if (mol.HasData(OBGenericDataType::CommentData))
      {
        OBCommentData *cd = (OBCommentData*)mol.GetData(OBGenericDataType::CommentData);
        ofs << cd->GetData();
      }

    ofs << endl;
    ofs << "@<TRIPOS>ATOM" << endl;

    OBAtom *atom;
    OBResidue *res;

    vector<OBAtom*>::iterator i;
    vector<int> labelcount;
    labelcount.resize( etab.GetNumberOfElements() );

    ttab.SetFromType("INT");
    ttab.SetToType("SYB");

    for (atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      {

        //
        //  Use sequentially numbered atom names if no residues
        //

        snprintf(label,BUFF_SIZE, "%s%d",
                 etab.GetSymbol(atom->GetAtomicNum()),
                 ++labelcount[atom->GetAtomicNum()]);
        strcpy(rlabel,"<1>");
        strcpy(rnum,"1");

        str = atom->GetType();
        ttab.Translate(str1,str);

        //
        //  Use original atom names if there are residues
        //

        if (!ligandsOnly && (res = atom->GetResidue()) )
          {
            // use original atom names defined by residue
            snprintf(label,BUFF_SIZE,"%s",(char*)res->GetAtomID(atom).c_str());
            // make sure that residue name includes its number
            snprintf(rlabel,BUFF_SIZE,"%s%d",res->GetName().c_str(), res->GetNum());
            snprintf(rnum,BUFF_SIZE,"%d",res->GetNum());
          }

        snprintf(buffer,BUFF_SIZE,"%7d%1s%-6s%12.4f%10.4f%10.4f%1s%-5s%4s%1s %-8s%10.4f",
                 atom->GetIdx(),"",label,
                 atom->GetX(),atom->GetY(),atom->GetZ(),
                 "",str1.c_str(),
                 rnum,"",rlabel,
                 atom->GetPartialCharge());
        ofs << buffer << endl;
      }
//.........这里部分代码省略.........

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++;
        }
    }

}

  bool AlchemyFormat::WriteMolecule(OBBase* pOb, OBConversion* pConv)
  {
    OBMol* pmol = dynamic_cast<OBMol*>(pOb);
    if(pmol==NULL)
      return false;

    //Define some references so we can use the old parameter names
    ostream &ofs = *pConv->GetOutStream();
    OBMol &mol = *pmol;

    unsigned int i;
    char buffer[BUFF_SIZE];
    char bond_string[10];

    snprintf(buffer, BUFF_SIZE, "%5d ATOMS, %5d BONDS,     0 CHARGES",
             mol.NumAtoms(),
             mol.NumBonds());
    ofs << buffer << endl;

    OBAtom *atom;
    string str,str1;
    for(i = 1;i <= mol.NumAtoms(); i++)
      {
        atom = mol.GetAtom(i);
        str = atom->GetType();
        ttab.SetFromType("INT");
        ttab.SetToType("ALC");
        ttab.Translate(str1,str);
        snprintf(buffer, BUFF_SIZE, "%5d %-6s%8.4f %8.4f %8.4f     0.0000",
                 i,
                 (char*)str1.c_str(),
                 atom->GetX(),
                 atom->GetY(),
                 atom->GetZ());
        ofs << buffer << endl;
      }

    OBBond *bond;
    vector<OBBond*>::iterator j;

    for (bond = mol.BeginBond(j);bond;bond = mol.NextBond(j))
      {
        switch(bond->GetBO())
          {
          case 1 :
            strcpy(bond_string,"SINGLE");
            break;
          case 2 :
            strcpy(bond_string,"DOUBLE");
            break;
          case 3 :
            strcpy(bond_string,"TRIPLE");
            break;
          case 5 :
            strcpy(bond_string,"AROMATIC");
            break;
          default :
            strcpy(bond_string,"SINGLE");
          }
        snprintf(buffer, BUFF_SIZE, "%5d  %4d  %4d  %s",
                 bond->GetIdx()+1,
                 bond->GetBeginAtomIdx(),
                 bond->GetEndAtomIdx(),
                 bond_string);
        ofs << buffer << endl;
      }
    return(true);
  }

int main(int argc,char *argv[])
{
  if (argc != 1) {
    if (strncmp(argv[1], "-g", 2))
      {
	cout << "Usage: ringtest" << endl;
	cout << "   Tests Open Babel ring perception testing." << endl;
	return 0;
      } else {
	GenerateRingReference();
	return 0;
      }
  }

  cout << endl << "Testing RINGS..." << endl;
  
#ifdef TESTDATADIR
  string testdatadir = TESTDATADIR;
  string results_file = testdatadir + "ringresults.txt";
  string smilestypes_file = testdatadir + "attype.00.smi";
#else
  string results_file = "ringresults.txt";
  string smilestypes_file = "attype.00.smi";
#endif

  std::ifstream mifs;
  if (!SafeOpen(mifs, (char*)smilestypes_file.c_str())) {
    return -1; // test failed
  }

  std::ifstream rifs;
  if (!SafeOpen(rifs, (char*)results_file.c_str())) {
    return -1; // test failed
  }

  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(SMI,SMI);
  vector<string>::iterator i;
  vector<OBEdgeBase*>::iterator j;
  vector<OBNodeBase*>::iterator k;
  vector<OBRing*>::iterator m;
  OBFileFormat ff;

  for (;mifs;)
    {
      mol.Clear();
      ff.ReadMolecule(mifs, mol);
      if (mol.Empty()) continue;
      if (!rifs.getline(buffer,BUFF_SIZE))
	{
	  ThrowError("error reading reference data");
	  return -1; // test failed
	}

      vb.clear();
      vb.resize(mol.NumBonds(),false);
      //check ring bonds
      tokenize(vs,buffer);
      for (i = vs.begin();i != vs.end();i++)
	  vb[atoi((char*)i->c_str())] = true;
      
      for (bond = mol.BeginBond(j);bond;bond = mol.NextBond(j))
	if (vb[bond->GetIdx()] != bond->IsInRing())
	  {
	    ThrowError("ring bond data different than reference");
	    ThrowError((char*)mol.GetTitle());
	    return -1; // test failed
	  }

      vr = mol.GetSSSR();
      if (!rifs.getline(buffer,BUFF_SIZE))
	{
	  ThrowError("error reading reference data");
	  return -1; // test failed
	}
      sscanf(buffer,"%d",&size);
      if (vr.size() != size) //check SSSR size
	{
	  ThrowError("SSSR size different than reference");
	  ThrowError((char*)mol.GetTitle());
	  return -1; // test failed
	}

      if (!rifs.getline(buffer,BUFF_SIZE))
	{
	  ThrowError("error reading reference data");
	  return -1; // test failed
	}

      tokenize(vs,buffer);
      i = vs.begin();
      for (atom = mol.BeginAtom(k);atom;atom = mol.NextAtom(k))
//.........这里部分代码省略.........