C++ OBAtom::x方法代码示例

bool WriteGromos96(ostream &ofs,OBMol &mol,double fac)
{ 
  char type_name[10];
  char res_name[10],padded_name[10];
  char buffer[BUFF_SIZE];
  int res_num;

  sprintf(buffer,"#GENERATED BY OPEN BABEL %s",BABEL_VERSION);
  ofs << buffer << endl;

  /* GROMOS wants a TITLE block, so let's write one*/
  sprintf(buffer,"TITLE\n%s\nEND",mol.GetTitle());
  ofs << buffer << endl;
  ofs << "POSITION" << endl;

  OBAtom *atom;
  OBResidue *res;
  vector<OBNodeBase*>::iterator i;

  for(atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
    {
      if (res = atom->GetResidue())
	{
	  strcpy(res_name,(char*)res->GetName().c_str());
	  strcpy(type_name,(char*)res->GetAtomID(atom).c_str());
	  res_num = res->GetNum();
	}
      else
	{
	  strcpy(type_name,etab.GetSymbol(atom->GetAtomicNum()));
	  strcpy(res_name,"UNK");
	  sprintf(padded_name,"%2s",type_name);
	  strcpy(type_name,padded_name);
	  res_num = 1;
	}
      
      sprintf(buffer,"%5d %5s %5s %6d %15.5f %15.5f %15.5f",
	      res_num,res_name,type_name,atom->GetIdx(),
	      atom->x()*fac,atom->y()*fac,atom->z()*fac);
      ofs << buffer << endl;

      if (!(atom->GetIdx()%10))
      {
	sprintf(buffer,"# %d",atom->GetIdx());
	ofs << buffer << endl;
      }
    }

  ofs << "END" << endl;

  return(true);
}

void CSRFormat::WriteCSRCoords(ostream &ofs,OBMol &mol)
{
    int the_size,jconf;
    double x,y,z,energy;
    char title[100];
    char *tag;

    the_size = sizeof(int) + sizeof(double) + (80 * sizeof(char));

    jconf = 1;
    energy = -2.584565;

    snprintf(title, 80, "%s:%d",mol.GetTitle(),MolCount);
    tag = PadString(title,80);

    WriteSize(the_size,ofs);
    ofs.write((char*)&jconf,sizeof(int));
    ofs.write((char*)&energy,sizeof(double));
    ofs.write(tag,80*sizeof(char));
    WriteSize(the_size,ofs);

    WriteSize(mol.NumAtoms()*sizeof(double),ofs);

    OBAtom *atom;
    vector<OBAtom*>::iterator i;
    for (atom = mol.BeginAtom(i); atom; atom = mol.NextAtom(i))
    {
        x = atom->x();
        ofs.write((char*)&x,sizeof(double));
    }
    WriteSize(mol.NumAtoms()*sizeof(double),ofs);

    WriteSize(mol.NumAtoms()*sizeof(double),ofs);
    for (atom = mol.BeginAtom(i); atom; atom = mol.NextAtom(i))
    {
        y = atom->y();
        ofs.write((char*)&y,sizeof(double));
    }
    WriteSize(mol.NumAtoms()*sizeof(double),ofs);

    WriteSize(mol.NumAtoms()*sizeof(double),ofs);
    for (atom = mol.BeginAtom(i); atom; atom = mol.NextAtom(i))
    {
        z = atom->z();
        ofs.write((char*)&z,sizeof(double));
    }
    WriteSize(mol.NumAtoms()*sizeof(double),ofs);

    delete [] tag;
}

  bool CacaoFormat::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;

    OBAtom *atom;
    char buffer[BUFF_SIZE];
    vector<OBAtom*>::iterator i;

    snprintf(buffer, BUFF_SIZE, "%s\n",mol.GetTitle());
    ofs << buffer;
    snprintf(buffer, BUFF_SIZE, "%3d   DIST  0  0  0\n",mol.NumAtoms());
    ofs << buffer;

    if (!mol.HasData(OBGenericDataType::UnitCell))
      ofs << "CELL 1.,1.,1.,90.,90.,90.\n";
    else
      {
        OBUnitCell *uc = (OBUnitCell*)mol.GetData(OBGenericDataType::UnitCell);
        snprintf(buffer, BUFF_SIZE, "CELL %f,%f,%f,%f,%f,%f\n",
                 uc->GetA(), uc->GetB(), uc->GetC(),
                 uc->GetAlpha(), uc->GetBeta(), uc->GetGamma());
        ofs << buffer;
      }

    for (atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      {
        snprintf(buffer,BUFF_SIZE,"%2s %7.4f, %7.4f, %7.4f\n",
                 etab.GetSymbol(atom->GetAtomicNum()),
                 atom->x(),
                 atom->y(),
                 atom->z());
        ofs << buffer;
      }

    return(true);
  }

bool WriteFeat(ostream &ofs,OBMol &mol)
{ 
  char buffer[BUFF_SIZE];
  
  ofs << mol.NumAtoms() << endl;
  ofs << mol.GetTitle() << endl;

  OBAtom *atom;
  vector<OBNodeBase*>::iterator i;
  for(atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
  {
    sprintf(buffer,"%-3s %8.5f  %8.5f  %8.5f ",
	    etab.GetSymbol(atom->GetAtomicNum()),
	    atom->x(),
	    atom->y(),
	    atom->z());
    ofs << buffer << endl;
  }

  return(true);
}

  bool GROMOS96Format::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;
    double fac = pConv->IsOption("n") ? 0.1 : 1.0; //new framework

    char type_name[16];
    char res_name[16];
    char buffer[BUFF_SIZE];
    string res_num;

    snprintf(buffer, BUFF_SIZE, "#GENERATED BY OPEN BABEL %s\n",BABEL_VERSION);
    ofs << buffer;

    /* GROMOS wants a TITLE block, so let's write one*/
    ofs << "TITLE\n" << mol.GetTitle() << "\nEND\n";
    ofs << "POSITION\n";

    OBAtom *atom;
    OBResidue *res;
    vector<OBAtom*>::iterator i;

    for(atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      {
        if ( (res = atom->GetResidue()) )
          {
            // 16 = sizeof(res_name) and sizeof(type_name)
            strncpy(res_name,(char*)res->GetName().c_str(), 16);
            res_name[15] = '\0';
            strncpy(type_name,(char*)res->GetAtomID(atom).c_str(), 16);
            type_name[15] = '\0';
            res_num = res->GetNumString();
          }
        else
          {
            strncpy(type_name,OBElements::GetSymbol(atom->GetAtomicNum()), 16);
            strcpy(res_name,"UNK");
            res_num = "1";
          }

        snprintf(buffer, BUFF_SIZE, "%5s %5s %5s %6d %15.5f %15.5f %15.5f\n",
                res_num.c_str(),res_name,type_name,atom->GetIdx(),
                atom->x()*fac,atom->y()*fac,atom->z()*fac);
        ofs << buffer;

        if (!(atom->GetIdx()%10))
          {
            snprintf(buffer, BUFF_SIZE, "# %d\n",atom->GetIdx());
            ofs << buffer;
          }
      }

    ofs << "END\n";

    return(true);
  }

bool BoxFormat::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;

    //margin hardwired in new framework. Also was in old fileformat
    double margin=1.0;

    char buffer[BUFF_SIZE];
    vector3 vcenter,vmin,vmax,vmid,vdim;

    OBAtom *atom;
    vector<OBAtom*>::iterator i;
    vmax.Set(-10E10,-10E10,-10E10);
    vmin.Set( 10E10, 10E10, 10E10);

    for (atom = mol.BeginAtom(i); atom; atom = mol.NextAtom(i))
    {
        vcenter += atom->GetVector();
        if (atom->x() < vmin.x())
            vmin.SetX(atom->x());
        if (atom->y() < vmin.y())
            vmin.SetY(atom->y());
        if (atom->z() < vmin.z())
            vmin.SetZ(atom->z());

        if (atom->x() > vmax.x())
            vmax.SetX(atom->x());
        if (atom->y() > vmax.y())
            vmax.SetY(atom->y());
        if (atom->z() > vmax.z())
            vmax.SetZ(atom->z());
    }
    vcenter /= (double)mol.NumAtoms();

    vector3 vmarg(margin,margin,margin);
    vmin -= vmarg;
    vmax += vmarg;
    vdim = vmax - vmin;
    vmid = vmin+vmax;
    vmid /= 2.0;

    ofs << "HEADER    CORNERS OF BOX" << endl;
    snprintf(buffer, BUFF_SIZE, "REMARK    CENTER (X Y Z)      %10.3f %10.3f %10.3f",
             vmid.x(),vmid.y(),vmid.z());
    ofs << buffer << endl;
    snprintf(buffer, BUFF_SIZE, "REMARK    DIMENSIONS (X Y Z)  %10.3f %10.3f %10.3f",
             vdim.x(),vdim.y(),vdim.z());
    ofs << buffer << endl;
    vdim /= 2.0;

    vector3 vtmp;
    int j;
    for (j = 1; j <= 8; j++)
    {
        switch(j)
        {
        case 1:
            vtmp = vmid-vdim;
            break;
        case 2:
            vtmp.SetX(vmid.x()+vdim.x());
            break;
        case 3:
            vtmp.SetZ(vmid.z()+vdim.z());
            break;
        case 4:
            vtmp.SetX(vmid.x()-vdim.x());
            break;
        case 5:
            vtmp = vmid-vdim;
            vtmp.SetY(vmid.y()+vdim.y());
            break;
        case 6:
            vtmp = vmid+vdim;
            vtmp.SetZ(vmid.z()-vdim.z());
            break;
        case 7:
            vtmp = vmid+vdim;
            break;
        case 8:
            vtmp.SetX(vmid.x()-vdim.x());
            break;
        }
        snprintf(buffer, BUFF_SIZE, "ATOM      %d  DUA BOX     1    %8.3f%8.3f%8.3f",
                 j,vtmp.x(),vtmp.y(),vtmp.z());
        ofs << buffer << endl;
    }

    ofs << "CONECT    1    2    4    5" << endl;
    ofs << "CONECT    2    1    3    6" << endl;
    ofs << "CONECT    3    2    4    7" << endl;
    ofs << "CONECT    4    1    3    8" << endl;
    ofs << "CONECT    5    1    6    8" << endl;
    ofs << "CONECT    6    2    5    7" << endl;
//.........这里部分代码省略.........

int main(int argc,char *argv[])
{
  // turn off slow sync with C-style output (we don't use it anyway).
  std::ios::sync_with_stdio(false);

  OBConversion conv;
  OBFormat *inFormat, *canFormat;
  OBMol mol;
  ifstream ifs;
  vector<OBMol> fragments;
  unsigned int fragmentCount = 0; // track how many in library -- give a running count
  map<string, int> index; // index of cansmi
  string currentCAN;
  unsigned int size;
  OBAtom *atom;
  OBBond *bond;
  bool nonRingAtoms, nonRingBonds;
  char buffer[BUFF_SIZE];

  canFormat = conv.FindFormat("can");
  conv.SetOutFormat(canFormat);

  if (argc < 2)
    {
      cout << "Usage: obfragment <file>" << endl;
      return(-1);
    }

  for (int i = 1; i < argc; i++) {
    cerr << " Reading file " << argv[i] << endl;

    inFormat = conv.FormatFromExt(argv[i]);
    if(inFormat==NULL || !conv.SetInFormat(inFormat))
      {
        cerr << " Cannot read file format for " << argv[i] << endl;
        continue; // try next file
      }
    
    ifs.open(argv[i]);
    
    if (!ifs)
      {
        cerr << "Cannot read input file: " << argv[i] << endl;
        continue;
      }
    
    
    while(ifs.peek() != EOF && ifs.good())
      {
        conv.Read(&mol, &ifs);
        if (!mol.Has3D()) continue; // invalid coordinates!
        mol.DeleteHydrogens(); // remove these before we do anything else
        
        do {
          nonRingAtoms = false;
          size = mol.NumAtoms();
          for (unsigned int i = 1; i <= size; ++i)
            {
              atom = mol.GetAtom(i);
              if (!atom->IsInRing()) {
                mol.DeleteAtom(atom);
                nonRingAtoms = true;
                break; // don't know how many atoms there are
              } 
              // Previously, we changed atoms to carbon here.
              // Now we perform this alchemy in terms of string-rewriting
              // once the canonical SMILES is generated
            }
        } while (nonRingAtoms);
        
        if (mol.NumAtoms() < 3)
          continue;
        
        if (mol.NumBonds() == 0)
          continue;
        
        do {
          nonRingBonds = false;
          size = mol.NumBonds();
          for (unsigned int i = 0; i < size; ++i)
            {
              bond = mol.GetBond(i);
              if (!bond->IsInRing()) {
                mol.DeleteBond(bond);
                nonRingBonds = true;
                break; // don't know how many bonds there are
              }
            }        
        } while (nonRingBonds);

        fragments = mol.Separate();
        for (unsigned int i = 0; i < fragments.size(); ++i)
          {
            if (fragments[i].NumAtoms() < 3) // too small to care
              continue;
              
            currentCAN = conv.WriteString(&fragments[i], true);
            currentCAN = RewriteSMILES(currentCAN); // change elements to "a/A" for compression
            if (index.find(currentCAN) != index.end()) { // already got this
              index[currentCAN] += 1; // add to the count for bookkeeping
//.........这里部分代码省略.........

//computes a set of solitary grid points that represent the interaction between
//this ligand and the provided receptor in some way
void OBAMolecule::computeInteractionGridPoints(OBAMolecule& receptor,
		MGrid& grid, double interactionDist,
		double maxClusterDist,
		unsigned minClusterPoints, double interactionPointRadius)
{
	grid.clear();
	//first construct a bounding box for the ligand while assembling a
	//vector of atomic coordinates
	BoundingBox ligandBox;
	vector<AtomPoint> points;
	points.reserve(mol.NumAtoms());
	for (OBAtomIterator aitr = mol.BeginAtoms(); aitr != mol.EndAtoms(); ++aitr)
	{
		OBAtom* atom = *aitr;
		points.push_back(AtomPoint(atom->x(), atom->y(), atom->z()));
		ligandBox.update(atom->x(), atom->y(), atom->z());
	}
	ligandBox.extend(interactionDist);

	//then identify all coordinates that are interacting
	double idistSq = interactionDist * interactionDist;
	OBMol& rmol = receptor.getMol();
	for (OBAtomIterator aitr = rmol.BeginAtoms(); aitr != rmol.EndAtoms();
			++aitr)
	{
		OBAtom* a = *aitr;
		if (ligandBox.contains(a->x(), a->y(), a->z()))
		{
			for (unsigned i = 0, n = points.size(); i < n; i++)
			{
				if (points[i].distSq(a->x(), a->y(), a->z()) <= idistSq)
				{
					points[i].interactingCnt++;
				}
			}
		}
	}

	//prune out non-interacting poitns
	vector<AtomPoint> tmp;
	tmp.reserve(points.size());
	for (unsigned i = 0, n = points.size(); i < n; i++)
	{
		if (points[i].interactingCnt > 0)
			tmp.push_back(points[i]);
	}
	points.swap(tmp);
	tmp.clear();

	//cluster these coordinates
	vector<vector<unsigned> > clusters;
	clusterPoints(points, maxClusterDist, clusters);

	//make the cluster centers the interaction grid points
	for (unsigned i = 0, n = clusters.size(); i < n; i++)
	{
		double xtot = 0, ytot = 0, ztot = 0;
		unsigned npts = clusters[i].size();

		if (npts >= minClusterPoints)
		{
			for (unsigned j = 0; j < npts; j++)
			{
				xtot += points[clusters[i][j]].x;
				ytot += points[clusters[i][j]].y;
				ztot += points[clusters[i][j]].z;
			}
			double xave = xtot / (double) npts;
			double yave = ytot / (double) npts;
			double zave = ztot / (double) npts;

			grid.setPoint(xave, yave, zave);
			if(interactionPointRadius > 0)
			{
				grid.markXYZSphere(xave,yave,zave,interactionPointRadius);
			}
		}
	}
}