C++ OBUnitCell::GetA方法代码示例

  void UnitCellExtension::setMolecule(Molecule *molecule)
  {
    m_molecule = molecule;

    if (m_molecule == NULL || m_dialog == NULL)
      return; // nothing we can do

    OBUnitCell *uc = m_molecule->OBUnitCell();
    if (!uc)
      return; // no unit cell

    // We don't want to send signals while we update
    disconnect(m_dialog, SIGNAL(unitCellParametersChanged(double, double, double, double, double, double)),
               this, SLOT(unitCellParametersChanged(double, double, double, double, double, double)));

    m_dialog->aLength(uc->GetA());
    m_dialog->bLength(uc->GetB());
    m_dialog->cLength(uc->GetC());

    m_dialog->alpha(uc->GetAlpha());
    m_dialog->beta(uc->GetBeta());
    m_dialog->gamma(uc->GetGamma());

    // reconnect
    connect(m_dialog, SIGNAL(unitCellParametersChanged(double, double, double, double, double, double)),
            this, SLOT(unitCellParametersChanged(double, double, double, double, double, double)));
  }

  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 FreeFormFractionalFormat::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;

    char buffer[BUFF_SIZE];
    OBUnitCell *uc = NULL;

    ofs << mol.GetTitle() << endl;

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

    vector3 v;
    FOR_ATOMS_OF_MOL(atom, mol)
      {
        v = atom->GetVector();
        if (uc != NULL)
          v = uc->CartesianToFractional(v);

        snprintf(buffer, BUFF_SIZE, "%s %10.5f%10.5f%10.5f",
                 etab.GetSymbol(atom->GetAtomicNum()),
                 v.x(),
                 v.y(),
                 v.z());
        ofs << buffer << endl;
      }

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

    char buffer[BUFF_SIZE];

    if (!mol.HasData(OBGenericDataType::UnitCell))
      {
        snprintf(buffer, BUFF_SIZE,
                 " REFERENCE STRUCTURE = 00000   A,B,C =%8.3f%8.3f%8.3f",
                 1.0,1.0,1.0);
        ofs << buffer << endl;
        snprintf(buffer, BUFF_SIZE,
                 "   ALPHA,BETA,GAMMA =%8.3f%8.3f%8.3f    SPGR =    P1"
                 , 90.0f, 90.0f, 90.0f);
        ofs << buffer << endl;
      }
    else
      {
        OBUnitCell *uc = (OBUnitCell*)mol.GetData(OBGenericDataType::UnitCell);
        snprintf(buffer, BUFF_SIZE,
                 " REFERENCE STRUCTURE = 00000   A,B,C =%8.3f%8.3f%8.3f",
                 uc->GetA(), uc->GetB(), uc->GetC());
        ofs << buffer << endl;
        snprintf(buffer, BUFF_SIZE,
                 "   ALPHA,BETA,GAMMA =%8.3f%8.3f%8.3f    SPGR =    P1",
                 uc->GetAlpha() , uc->GetBeta(), uc->GetGamma());
        ofs << buffer << endl;
      }

    snprintf(buffer, BUFF_SIZE, "%4d   1 %s\n",mol.NumAtoms(), mol.GetTitle());
    ofs << buffer << endl << endl;

    OBAtom *atom,*nbr;
    vector<OBAtom*>::iterator i;
    vector<OBBond*>::iterator j;
    vector<int> vtmp(106,0);
    int bonds;

    for(atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      {
        //assign_pdb_number(pdb_types,atom->GetIdx());
        vtmp[atom->GetAtomicNum()]++;
        snprintf(buffer, BUFF_SIZE, "%4d%2s%-3d  %9.5f %9.5f %9.5f ",
                 atom->GetIdx(),
                 etab.GetSymbol(atom->GetAtomicNum()),
                 vtmp[atom->GetAtomicNum()],
                 atom->x(),
                 atom->y(),
                 atom->z());
        ofs << buffer;
        bonds = 0;
        for (nbr = atom->BeginNbrAtom(j); nbr; nbr = atom->NextNbrAtom(j))
          {
            if (bonds > 8) break;
            snprintf(buffer, BUFF_SIZE, "%4d",nbr->GetIdx());
            ofs << buffer;
            bonds++;
          }
        for (; bonds < 8; bonds ++)
          {
            snprintf(buffer, BUFF_SIZE, "%4d",0);
            ofs << buffer;
          }
        snprintf(buffer, BUFF_SIZE, " %7.3f%4d", atom->GetPartialCharge(), 1);
        ofs << buffer << endl;
      }

    return(true);
  }

  QUndoCommand* UnitCellExtension::performAction(QAction *, GLWidget *widget)
  {
    // FIXME: this is bad mmmkay
    m_widget = widget;

    if (m_molecule == NULL)
      return NULL; // nothing we can do

    OBUnitCell *uc = m_molecule->OBUnitCell();

    if (uc == NULL) {
      // show warning and ask if the user wants to create a unit cell
      // (otherwise this extension isn't very useful)

      QMessageBox::StandardButton ret;
      ret = QMessageBox::warning(qobject_cast<QWidget*>(parent()),
                                 tr("Avogadro"),
                                 tr("This document is currently an isolated molecule.\n\n"
                                    "Do you want to create a crystal unit cell?"),
                                 QMessageBox::Yes
                                 | QMessageBox::No);
      if (ret == QMessageBox::Yes) {
        // Set some initial data (e.g., a box about the size of the molecule)
        // and one unit cell in each direction
        uc = new OBUnitCell;
        double estimatedSize = widget->radius() + 2.0;
        uc->SetData(estimatedSize, estimatedSize, estimatedSize,
                    90.0, 90.0, 90.0);
        m_molecule->setOBUnitCell(uc);

        widget->setUnitCells(1, 1, 1);
        widget->setUnitCellColor(m_color);
      } else { // do nothing -- user picked "Cancel"
        return NULL;
      }

    } // end if (existing unit cell or create a new one)

    // Don't emit signals while we update these
    disconnect(m_dialog, SIGNAL(unitCellDisplayChanged(int, int, int)),
            this, SLOT(unitCellDisplayChanged(int, int, int)));
    disconnect(m_dialog, SIGNAL(unitCellParametersChanged(double, double, double, double, double, double)),
            this, SLOT(unitCellParametersChanged(double, double, double, double, double, double)));

    m_dialog->aCells(widget->aCells());
    m_dialog->bCells(widget->bCells());
    m_dialog->cCells(widget->cCells());

    m_dialog->aLength(uc->GetA());
    m_dialog->bLength(uc->GetB());
    m_dialog->cLength(uc->GetC());

    m_dialog->alpha(uc->GetAlpha());
    m_dialog->beta(uc->GetBeta());
    m_dialog->gamma(uc->GetGamma());

    // OK, now we can handle signal/slots
    connect(m_dialog, SIGNAL(unitCellDisplayChanged(int, int, int)),
            this, SLOT(unitCellDisplayChanged(int, int, int)));
    connect(m_dialog, SIGNAL(unitCellParametersChanged(double, double, double, double, double, double)),
            this, SLOT(unitCellParametersChanged(double, double, double, double, double, double)));

    m_dialog->show();

    return NULL;
  }

  bool PDBFormat::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 type_name[10], padded_name[10];
    char the_res[10];
    char the_chain = ' ';
    const char *element_name;
    int res_num;
    bool het=true;
    int model_num = 0;
    if (!pConv->IsLast() || pConv->GetOutputIndex() > 1)
      { // More than one molecule record
        model_num = pConv->GetOutputIndex(); // MODEL 1-based index
        snprintf(buffer, BUFF_SIZE, "MODEL %8d", model_num);
        ofs << buffer << endl;
      }

    // write back all fields (REMARKS, HELIX, SHEET, SITE, ...)
    bool compndWritten = false;
    bool authorWritten = false;
    std::vector<OBGenericData*> pairData = mol.GetAllData(OBGenericDataType::PairData);
    for (std::vector<OBGenericData*>::iterator data = pairData.begin(); data != pairData.end(); ++data) {
      OBPairData *pd = static_cast<OBPairData*>(*data);
      string attr = pd->GetAttribute();

      // filter to make sure we are writing pdb fields only
      if (attr != "HEADER" && attr != "OBSLTE" && attr != "TITLE" && attr != "SPLIT" &&
          attr != "CAVEAT" && attr != "COMPND" && attr != "SOURCE" && attr != "KEYWDS" &&
          attr != "EXPDTA" && attr != "NUMMDL" && attr != "MDLTYP" && attr != "AUTHOR" &&
          attr != "REVDAT" && attr != "SPRSDE" && attr != "JRNL" && attr != "REMARK" &&
          attr != "DBREF" && attr != "DBREF1" && attr != "DBREF2" && attr != "SEQADV" &&
          attr != "SEQRES" && attr != "MODRES" && attr != "HET" && attr != "HETNAM" &&
          attr != "HETSYN" && attr != "FORMUL" && attr != "HELIX" && attr != "SHEET" &&
          attr != "SSBOND" && attr != "LINK" && attr != "CISPEP" && attr != "SITE" &&
          attr != "ORIGX1" && attr != "ORIGX2" && attr != "ORIGX3" && attr != "SCALE1" &&
          attr != "SCALE2" && attr != "SCALE3" && attr != "MATRIX1" && attr != "MATRIX2" &&
          attr != "MATRIX3" && attr != "MODEL")
        continue;

      if (attr == "COMPND")
        compndWritten = true;
      if (attr == "AUTHOR")
        authorWritten = true;

      // compute spacing needed. HELIX, SITE, HET, ... are trimmed when reading
      int nSpacing = 6 - attr.size();
      for (int i = 0; i < nSpacing; ++i)
        attr += " ";


      std::string lines = pd->GetValue();
      string::size_type last = 0;
      string::size_type pos = lines.find('\n');
      while (last != string::npos) {
        string line = lines.substr(last, pos - last);
        if (pos == string::npos)
          last = string::npos;
        else
          last = pos + 1;
        pos = lines.find('\n', last);

        ofs << attr << line << endl;
      }
    }

    if (!compndWritten) {
      if (strlen(mol.GetTitle()) > 0)
        snprintf(buffer, BUFF_SIZE, "COMPND    %s ",mol.GetTitle());
      else
        snprintf(buffer, BUFF_SIZE, "COMPND    UNNAMED");
      ofs << buffer << endl;
    }

    if (!authorWritten) {
      snprintf(buffer, BUFF_SIZE, "AUTHOR    GENERATED BY OPEN BABEL %s",BABEL_VERSION);
      ofs << buffer << endl;
    }

    // Write CRYST1 record, containing unit cell parameters, space group
    // and Z value (supposed to be 1)
    if (pmol->HasData(OBGenericDataType::UnitCell))
      {
        OBUnitCell *pUC = (OBUnitCell*)pmol->GetData(OBGenericDataType::UnitCell);
        if(pUC->GetSpaceGroup()){
          string tmpHM=pUC->GetSpaceGroup()->GetHMName();
          // Do we have an extended HM symbol, with origin choice as ":1" or ":2" ? If so, remove it.
          size_t n=tmpHM.find(":");
          if(n!=string::npos) tmpHM=tmpHM.substr(0,n);
          snprintf(buffer, BUFF_SIZE,
                   "CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s 1",
                   pUC->GetA(), pUC->GetB(), pUC->GetC(),
//.........这里部分代码省略.........

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

    vector<OBAtom*>::iterator i;
    int max_val;
    OBAtom *atom;
    char buffer[BUFF_SIZE];
    char elmnt_typ[8], dreid_typ[8], atm_sym[16], max_val_str[8];

    mol.Kekulize();

    ofs << "BIOGRF 200\n";
    snprintf(buffer, BUFF_SIZE, "DESCRP %s\n",mol.GetTitle());
    ofs << buffer;
    snprintf(buffer, BUFF_SIZE, "REMARK BGF file created by Open Babel %s\n",BABEL_VERSION);
    ofs << "FORCEFIELD DREIDING  \n";

    // write unit cell if available
    if (mol.HasData(OBGenericDataType::UnitCell))
      {
        OBUnitCell *uc = (OBUnitCell*)mol.GetData(OBGenericDataType::UnitCell);
        // e.g. CRYSTX    49.30287   49.23010   25.45631   90.00008   89.99995   57.10041
        snprintf(buffer, BUFF_SIZE,
                 "CRYSTX%12.5f%12.5f%12.5f%12.5f%12.5f%12.5f",
                 uc->GetA(), uc->GetB(), uc->GetC(),
                 uc->GetAlpha() , uc->GetBeta(), uc->GetGamma());
        ofs << buffer << "\n";
      }

    ofs << "FORMAT ATOM   (a6,1x,i5,1x,a5,1x,a3,1x,a1,1x,a5,3f10.5,1x,a5,i3,i2,1x,f8.5)\n";

    ttab.SetFromType("INT");

    for (atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      {
        strncpy(elmnt_typ,etab.GetSymbol(atom->GetAtomicNum()), 7); // make sure to null-terminate
        elmnt_typ[sizeof(elmnt_typ) - 1] = '0';
        ToUpper(elmnt_typ);

        ttab.SetToType("DRE");
        ttab.Translate(dreid_typ,atom->GetType());
        ttab.SetToType("HAD");
        ttab.Translate(max_val_str,atom->GetType());
        max_val = atoi(max_val_str);
        if (max_val == 0)
          max_val = 1;
        snprintf(atm_sym,16,"%s%d",elmnt_typ,atom->GetIdx());
        snprintf(buffer,BUFF_SIZE,"%6s %5d %-5s %3s %1s %5s%10.5f%10.5f%10.5f %-5s%3d%2d %8.5f\n",
                "HETATM",
                atom->GetIdx(),
                atm_sym,
                "RES",
                "A",
                "444",
                atom->GetX(),
                atom->GetY(),
                atom->GetZ(),
                dreid_typ,
                max_val,
                0,
                atom->GetPartialCharge());
        ofs << buffer;
      }
    ofs<< "FORMAT CONECT (a6,12i6)\n\n";

    OBAtom *nbr;
    vector<OBBond*>::iterator j;
    for (atom = mol.BeginAtom(i);atom;atom = mol.NextAtom(i))
      if (atom->GetValence())
        {
          snprintf(buffer,BUFF_SIZE,"CONECT%6d",atom->GetIdx());
          ofs << buffer;
          for (nbr = atom->BeginNbrAtom(j);nbr;nbr = atom->NextNbrAtom(j))
            {
              snprintf(buffer,BUFF_SIZE,"%6d",nbr->GetIdx());
              ofs << buffer;
            }
          ofs << endl;

          snprintf(buffer,BUFF_SIZE,"ORDER %6d",atom->GetIdx());
          ofs << buffer;
          for (nbr = atom->BeginNbrAtom(j);nbr;nbr = atom->NextNbrAtom(j))
            {
              snprintf(buffer,BUFF_SIZE,"%6d",(*j)->GetBO());
              ofs << buffer;
            }
          ofs << endl;
        }

    ofs << "END" << endl;
    return(true);
  }

  //! \return whether partial charges were successfully assigned to this molecule
  bool EQEqCharges::ComputeCharges(OBMol &mol)
  {
    int i, j, a, c, N = mol.NumAtoms();
    double cellVolume;
    VectorXf chi(N), J(N), b(N), x(N);
    MatrixXf J_ij(N, N), A(N, N);
    OBUnitCell *obuc;
    matrix3x3 unitcell, fourier;
    vector3 dx;
    int numNeighbors[3];
    OBAtom *atom;

    // If parameters have not yet been loaded, do that
    if (!_paramFileLoaded)
    {
      if (ParseParamFile())
      {
        _paramFileLoaded = true;
      } else
      {
        return false;
      }
    }

    // Calculate atomic properties based around their ionic charge
    for (i = 0; i < N; i++)
    {
      atom = mol.GetAtom(i + 1);
      a = atom->GetAtomicNum();
      c = _chargeCenter[a];

      // Fail if ionization data is missing for any atom in the molecule
      if (_ionizations[a][c + 1] == -1 || _ionizations[a][c] == -1 || a > TABLE_OF_ELEMENTS_SIZE)
      {
        obErrorLog.ThrowError(__FUNCTION__, "Insufficient ionization data for atoms in the given molecule. Update `data/eqeqIonizations.txt` with missing information and re-run this function.", obError);
        return false;
      }

      J(i) = _ionizations[a][c + 1] - _ionizations[a][c];
      chi(i) = 0.5 * (_ionizations[a][c + 1] + _ionizations[a][c]) - (a == 1? 0 : c * J(i));
    }

    // If a unit cell is defined, use the periodic Ewald calculation
    if (mol.HasData(OBGenericDataType::UnitCell))
    {
      // Get unit cell and calculate its Fourier transform + volume
      obuc = (OBUnitCell *) mol.GetData(OBGenericDataType::UnitCell);
      unitcell = obuc->GetCellMatrix();
      fourier = (2 * PI * unitcell.inverse()).transpose();
      cellVolume = obuc->GetCellVolume();

      // Get the number of radial unit cells to use in x, y, and z
      numNeighbors[0] = int(ceil(minCellLength / (2.0 * (obuc->GetA())))) - 1;
      numNeighbors[1] = int(ceil(minCellLength / (2.0 * (obuc->GetB())))) - 1;
      numNeighbors[2] = int(ceil(minCellLength / (2.0 * (obuc->GetC())))) - 1;

      for (i = 0; i < N; i++)
      {
        atom = mol.GetAtom(i + 1);
        for (j = 0; j < N; j++)
        {
          dx = atom->GetVector() - (mol.GetAtom(j + 1))->GetVector();
          J_ij(i, j) = GetPeriodicEwaldJij(J(i), J(j), dx, (i == j), unitcell, fourier, cellVolume, numNeighbors);
        }
      }
    // If no unit cell, use the simplified nonperiodic calculation
    } else
    {
      for (i = 0; i < N; i++)
      {
        atom = mol.GetAtom(i + 1);
        for (j = 0; j < N; j++)
        {
          J_ij(i, j) = GetNonperiodicJij(J(i), J(j), atom->GetDistance(j + 1), (i == j));
        }
        return false;
      }
    }

    // Formulate problem as A x = b, where x is the calculated partial charges
    // First equation is a simple overall balance: sum(Q) = 0
    A.row(0) = VectorXf::Ones(N);
    b(0) = 0;

    // Remaining equations are based off of the fact that, at equilibrium, the
    // energy of the system changes equally for a change in any charge:
    //     dE/dQ_1 = dE/dQ_2 = ... = dE/dQ_N
    A.block(1, 0, N - 1, N) = J_ij.block(0, 0, N - 1, N) - J_ij.block(1, 0, N - 1, N);
    b.tail(N - 1) = chi.tail(N - 1) - chi.head(N - 1);

    // The solution is a list of charges in the system
    x = A.colPivHouseholderQr().solve(b);

    // Now we are done calculating, pass all this back to OpenBabel molecule
    mol.SetPartialChargesPerceived();
    OBPairData *dp = new OBPairData;
    dp->SetAttribute("PartialCharges");
    dp->SetValue("EQEq");
    dp->SetOrigin(perceived);
//.........这里部分代码省略.........

  bool PDBFormat::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 type_name[10], padded_name[10];
    char the_res[10];
    char the_chain = ' ';
    const char *element_name;
    int res_num;
    bool het=true;
    int model_num = 0;
    if (!pConv->IsLast() || pConv->GetOutputIndex() > 1)
      { // More than one molecule record
        model_num = pConv->GetOutputIndex(); // MODEL 1-based index
        snprintf(buffer, BUFF_SIZE, "MODEL %8d", model_num);
        ofs << buffer << endl;
      }

    if (strlen(mol.GetTitle()) > 0)
      snprintf(buffer, BUFF_SIZE, "COMPND    %s ",mol.GetTitle());
    else
      snprintf(buffer, BUFF_SIZE, "COMPND    UNNAMED");
    ofs << buffer << endl;

    snprintf(buffer, BUFF_SIZE, "AUTHOR    GENERATED BY OPEN BABEL %s",BABEL_VERSION);
    ofs << buffer << endl;

    // Write CRYST1 record, containing unit cell parameters, space group
    // and Z value (supposed to be 1)
    if (pmol->HasData(OBGenericDataType::UnitCell))
      {
        OBUnitCell *pUC = (OBUnitCell*)pmol->GetData(OBGenericDataType::UnitCell);
	
        snprintf(buffer, BUFF_SIZE,
                 "CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s 1",
                 pUC->GetA(), pUC->GetB(), pUC->GetC(),
                 pUC->GetAlpha(), pUC->GetBeta(), pUC->GetGamma(),
                 pUC->GetSpaceGroup() ?
                 pUC->GetSpaceGroup()->GetHMName().c_str() : "P1");
        ofs << buffer << endl;
      }

    // before we write any records, we should check to see if any coord < -1000
    // which will cause errors in the formatting

    double minX, minY, minZ;
    minX = minY = minZ = -999.0f;
    FOR_ATOMS_OF_MOL(a, mol)
      {
        if (a->GetX() < minX)
          minX = a->GetX();
        if (a->GetY() < minY)
          minY = a->GetY();
        if (a->GetZ() < minZ)
          minZ = a->GetZ();
      }
    vector3 transV = VZero;
    if (minX < -999.0)
      transV.SetX(-1.0*minX - 900.0);
    if (minY < -999.0)
      transV.SetY(-1.0*minY - 900.0);
    if (minZ < -999.0)
      transV.SetZ(-1.0*minZ - 900.0);

    // if minX, minY, or minZ was never changed, shift will be 0.0f
    // otherwise, move enough so that smallest coord is > -999.0f
    mol.Translate(transV);

    OBAtom *atom;
    OBResidue *res;
    for (i = 1; i <= mol.NumAtoms(); i++)
      {
        atom = mol.GetAtom(i);
        strncpy(type_name, etab.GetSymbol(atom->GetAtomicNum()), sizeof(type_name));
        type_name[sizeof(type_name) - 1] = '\0';

        //two char. elements are on position 13 and 14 one char. start at 14
        if (strlen(type_name) > 1)
          type_name[1] = toupper(type_name[1]);
        else
          {
            char tmp[10];
            strncpy(tmp, type_name, 10);
            snprintf(type_name, sizeof(type_name), " %-3s", tmp);
          }

        if ( (res = atom->GetResidue()) != 0 )
          {
            het = res->IsHetAtom(atom);
            snprintf(the_res,4,"%s",(char*)res->GetName().c_str());
            snprintf(type_name,5,"%s",(char*)res->GetAtomID(atom).c_str());
            the_chain = res->GetChain();
//.........这里部分代码省略.........