本文整理汇总了C++中Topology::AddTopAtom方法的典型用法代码示例。如果您正苦于以下问题:C++ Topology::AddTopAtom方法的具体用法?C++ Topology::AddTopAtom怎么用?C++ Topology::AddTopAtom使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类Topology
的用法示例。
在下文中一共展示了Topology::AddTopAtom方法的7个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。
示例1: ReadParm
int Parm_PDB::ReadParm(std::string const& fname, Topology &TopIn) {
PDBfile infile;
double XYZ[3];
int current_res = 0;
int last_res = -1;
if (infile.OpenRead(fname)) return 1;
// Loop over PDB records
while ( infile.NextLine() != 0 ) {
if (infile.IsPDBatomKeyword()) {
// If this is an ATOM / HETATM keyword, add to topology
infile.pdb_XYZ(XYZ);
NameType pdbresname = infile.pdb_Residue( current_res );
TopIn.AddTopAtom(infile.pdb_Atom(), pdbresname, current_res, last_res, XYZ);
} else if (infile.IsPDB_TER() || infile.IsPDB_END()) {
// Indicate end of molecule for TER/END. Finish if END.
TopIn.StartNewMol();
if (infile.IsPDB_END()) break;
}
}
// If Topology name not set with TITLE etc, use base filename.
// TODO: Read in title.
std::string pdbtitle;
TopIn.SetParmName( pdbtitle, infile.Filename().Base() );
infile.CloseFile();
return 0;
}
示例2: mprintf
/** Read file as a Tripos Mol2 file. */
int Parm_Mol2::ReadParm(FileName const& fname, Topology &parmOut) {
Mol2File infile;
if (infile.OpenRead(fname)) return 1;
mprintf(" Reading Mol2 file %s as topology file.\n",infile.Filename().base());
// Get @<TRIPOS>MOLECULE information
if (infile.ReadMolecule()) return 1;
parmOut.SetParmName( infile.Mol2Title(), infile.Filename() );
// Get @<TRIPOS>ATOM information
if (infile.ScanTo( Mol2File::ATOM)) return 1;
double XYZ[3];
for (int atom=0; atom < infile.Mol2Natoms(); atom++) {
if ( infile.Mol2XYZ(XYZ) ) return 1;
parmOut.AddTopAtom( infile.Mol2Atom(), infile.Mol2Residue(), XYZ );
}
// Get @<TRIPOS>BOND information [optional]
int at1 = 0;
int at2 = 0;
if (infile.ScanTo(Mol2File::BOND)==0) {
for (int bond=0; bond < infile.Mol2Nbonds(); bond++) {
if (infile.Mol2Bond(at1, at2)) return 1;
// mol2 atom #s start from 1
parmOut.AddBond(at1-1, at2-1);
}
needsBondSearch_ = false;
} else {
mprintf(" Mol2 file does not contain bond information.\n");
needsBondSearch_ = true;
}
// No box
parmOut.SetParmBox( Box() );
mprintf(" Mol2 contains %i atoms, %i residues,\n", parmOut.Natom(),parmOut.Nres());
//mprintf(" %i bonds to H, %i other bonds.\n", parmOut.NbondsWithH,parmOut.NbondsWithoutH);
infile.CloseFile();
return 0;
}
示例3: ReadParm
/** Read file as a Tinker file. */
int Parm_Tinker::ReadParm(FileName const& fname, Topology &parmOut) {
TinkerFile infile;
infile.SetTinkerName( fname );
if (infile.OpenTinker()) return 1;
mprintf("\tReading Tinker file %s as topology file.\n",infile.Filename().base());
// Allocate memory for coordinates.
double* Coords = new double[ infile.TinkerNatom() * 3 ];
std::vector<int> Bonds;
std::vector<Atom> Atoms = infile.ReadTinkerAtoms(Coords, Bonds);
if (Atoms.empty()) return 1;
// Use up to first 3 chars of title as residue name.
std::string resname;
for (std::string::const_iterator c = infile.TinkerTitle().begin();
c != infile.TinkerTitle().end(); ++c)
resname += *c;
if (resname.size() > 3) resname.resize(3);
Residue tinker_res( resname, 0, ' ', ' ' );
// Put atoms into topology
const double* XYZ = Coords;
for (std::vector<Atom>::const_iterator atom = Atoms.begin();
atom != Atoms.end();
++atom, XYZ += 3)
parmOut.AddTopAtom( *atom, tinker_res, XYZ );
delete[] Coords;
// Add bond information
for (std::vector<int>::const_iterator bond = Bonds.begin();
bond != Bonds.end(); bond += 2)
parmOut.AddBond( *bond, *(bond+1) );
// Try to set up residue info based on bonds.
if (parmOut.Setup_NoResInfo()) return 1;
// Set topology box info.
parmOut.SetParmBox( infile.TinkerBox() );
parmOut.SetParmName( infile.TinkerTitle(), infile.Filename() );
mprintf("\tTinker file contains %i atoms, %i residues,\n", parmOut.Natom(),parmOut.Nres());
//mprintf(" %i bonds to H, %i other bonds.\n", parmOut.NbondsWithH,parmOut.NbondsWithoutH);
infile.CloseFile();
return 0;
}
示例4: Analyze
//.........这里部分代码省略.........
for (std::vector<DataSet_1D*>::iterator ds = histdata_.begin();
ds != histdata_.end(); ++ds, ++dim, ++bOff)
{
double dval = (*ds)->Dval( n );
// Check if data is out of bounds for this dimension.
if (dval > dim->Max() || dval < dim->Min()) {
index = -1L;
break;
}
// Calculate index for this particular dimension (idx)
long int idx = (long int)((dval - dim->Min()) / dim->Step());
if (debug_>1) mprintf(" [%s:%f (%li)],", dim->label(), dval, idx);
// Calculate overall index in Bins, offset has already been calcd.
index += (idx * (*bOff));
}
// If index was successfully calculated, populate bin
if (index > -1L && index < (long int)Bins_.size()) {
if (debug_ > 1) mprintf(" |index=%li",index);
if (calcAMD_)
Bins_[index] += exp( amddata_->Dval(n) );
else
Bins_[index]++;
} else {
mprintf("\tWarning: Frame %zu Coordinates out of bounds (%li)\n", n+1, index);
}
if (debug_>1) mprintf("}\n");
}
// Calc free energy if requested
if (calcFreeE_) CalcFreeE();
// Normalize if requested
if (normalize_ != NO_NORM) Normalize();
if (nativeOut_) {
// Use Histogram built-in output
PrintBins();
} else {
// Using DataFileList framework, set-up labels etc.
if (N_dimensions_ == 1) {
DataSet_double& dds = static_cast<DataSet_double&>( *hist_ );
// Since Allocate1D only reserves data, use assignment op.
dds = Bins_;
hist_->SetDim(Dimension::X, dimensions_[0]);
} else if (N_dimensions_ == 2) {
DataSet_MatrixDbl& mds = static_cast<DataSet_MatrixDbl&>( *hist_ );
mds.Allocate2D( dimensions_[0].Bins(), dimensions_[1].Bins() );
std::copy( Bins_.begin(), Bins_.end(), mds.begin() );
hist_->SetDim(Dimension::X, dimensions_[0]);
hist_->SetDim(Dimension::Y, dimensions_[1]);
outfile_->ProcessArgs("noxcol usemap nolabels");
} else if (N_dimensions_ == 3) {
DataSet_GridFlt& gds = static_cast<DataSet_GridFlt&>( *hist_ );
//gds.Allocate3D( dimensions_[0].Bins(), dimensions_[1].Bins(), dimensions_[2].Bins() );
gds.Allocate_N_O_D( dimensions_[0].Bins(), dimensions_[1].Bins(), dimensions_[2].Bins(),
Vec3(dimensions_[0].Min(), dimensions_[1].Min(), dimensions_[2].Min()),
Vec3(dimensions_[0].Step(), dimensions_[1].Step(), dimensions_[2].Step())
);
//std::copy( Bins_.begin(), Bins_.end(), gds.begin() );
// FIXME: Copy will not work since in grids data is ordered with Z
// changing fastest. Should the ordering in grid be changed?
size_t idx = 0;
for (size_t z = 0; z < gds.NZ(); z++)
for (size_t y = 0; y < gds.NY(); y++)
for (size_t x = 0; x < gds.NX(); x++)
gds.SetElement( x, y, z, (float)Bins_[idx++] );
hist_->SetDim(Dimension::X, dimensions_[0]);
hist_->SetDim(Dimension::Y, dimensions_[1]);
hist_->SetDim(Dimension::Z, dimensions_[2]);
outfile_->ProcessArgs("noxcol usemap nolabels");
// Create pseudo-topology/trajectory
if (!traj3dName_.empty()) {
Topology pseudo;
pseudo.AddTopAtom(Atom("H3D", 0), Residue("H3D", 1, ' ', ' '));
pseudo.CommonSetup();
if (!parmoutName_.empty()) {
ParmFile pfile;
if (pfile.WriteTopology( pseudo, parmoutName_, ParmFile::UNKNOWN_PARM, 0 ))
mprinterr("Error: Could not write pseudo topology to '%s'\n", parmoutName_.c_str());
}
Trajout_Single out;
if (out.PrepareTrajWrite(traj3dName_, ArgList(), &pseudo, CoordinateInfo(),
Ndata, traj3dFmt_) == 0)
{
Frame outFrame(1);
for (size_t i = 0; i < Ndata; ++i) {
outFrame.ClearAtoms();
outFrame.AddVec3( Vec3(histdata_[0]->Dval(i),
histdata_[1]->Dval(i),
histdata_[2]->Dval(i)) );
out.WriteSingle(i, outFrame);
}
out.EndTraj();
} else
mprinterr("Error: Could not set up '%s' for write.\n", traj3dName_.c_str());
}
}
}
return Analysis::OK;
}
示例5: ReadParm
// Parm_CIF::ReadParm()
int Parm_CIF::ReadParm(FileName const& fname, Topology &TopIn) {
CIFfile infile;
CIFfile::DataBlock::data_it line;
if (infile.Read( fname, debug_ )) return 1;
CIFfile::DataBlock const& block = infile.GetDataBlock("_atom_site");
if (block.empty()) {
mprinterr("Error: CIF data block '_atom_site' not found.\n");
return 1;
}
// Does this CIF contain multiple models?
int Nmodels = 0;
int model_col = block.ColumnIndex("pdbx_PDB_model_num");
if (model_col != -1) {
line = block.end();
--line;
Nmodels = convertToInteger( (*line)[model_col] );
if (Nmodels > 1)
mprintf("Warning: CIF '%s' contains %i models. Using first model for topology.\n",
fname.full(), Nmodels);
}
// Get essential columns
int COL[NENTRY];
for (int i = 0; i < (int)NENTRY; i++) {
COL[i] = block.ColumnIndex(Entries[i]);
if (COL[i] == -1) {
mprinterr("Error: In CIF file '%s' could not find entry '%s' in block '%s'\n",
fname.full(), Entries[i], block.Header().c_str());
return 1;
}
if (debug_>0) mprintf("DEBUG: '%s' column = %i\n", Entries[i], COL[i]);
}
// Get optional columns
int occ_col = block.ColumnIndex("occupancy");
int bfac_col = block.ColumnIndex("B_iso_or_equiv");
int icode_col = block.ColumnIndex("pdbx_PDB_ins_code");
int altloc_col = block.ColumnIndex("label_alt_id");
std::vector<AtomExtra> extra;
// Loop over all atom sites
int current_res = 0;
double XYZ[3];
double occupancy = 1.0;
double bfactor = 0.0;
char altloc = ' ';
char icode;
icode = ' ';
Frame Coords;
for (line = block.begin(); line != block.end(); ++line) {
// If more than 1 model check if we are done.
if (Nmodels > 1) {
if ( convertToInteger( (*line)[model_col] ) > 1 )
break;
}
if (occ_col != -1) occupancy = convertToDouble( (*line)[ occ_col ] );
if (bfac_col != -1) bfactor = convertToDouble( (*line)[ bfac_col ] );
if (altloc_col != -1) altloc = (*line)[ altloc_col ][0];
// '.' altloc means blank?
if (altloc == '.') altloc = ' ';
extra.push_back( AtomExtra(occupancy, bfactor, altloc) );
if (icode_col != -1) {
icode = (*line)[ icode_col ][0];
// '?' icode means blank
if (icode == '?') icode = ' ';
}
XYZ[0] = convertToDouble( (*line)[ COL[X] ] );
XYZ[1] = convertToDouble( (*line)[ COL[Y] ] );
XYZ[2] = convertToDouble( (*line)[ COL[Z] ] );
NameType currentResName( (*line)[ COL[RNAME] ] );
// It seems that in some CIF files, there doesnt have to be a residue
// number. Check if residue name has changed.
if ( (*line)[ COL[RNUM] ][0] == '.' ) {
Topology::res_iterator lastResidue = TopIn.ResEnd();
--lastResidue;
if ( currentResName != (*lastResidue).Name() )
current_res = TopIn.Nres() + 1;
} else
current_res = convertToInteger( (*line)[ COL[RNUM] ] );
TopIn.AddTopAtom( Atom((*line)[ COL[ANAME] ], " "),
Residue(currentResName, current_res, icode,
(*line)[ COL[CHAINID] ][0]) );
Coords.AddXYZ( XYZ );
}
if (TopIn.SetExtraAtomInfo( 0, extra )) return 1;
// Search for bonds // FIXME nobondsearch?
BondSearch( TopIn, Coords, Offset_, debug_ );
// Get title.
CIFfile::DataBlock const& entryblock = infile.GetDataBlock("_entry");
std::string ciftitle;
if (!entryblock.empty())
ciftitle = entryblock.Data("id");
TopIn.SetParmName( ciftitle, infile.CIFname() );
// Get unit cell parameters if present.
CIFfile::DataBlock const& cellblock = infile.GetDataBlock("_cell");
if (!cellblock.empty()) {
double cif_box[6];
cif_box[0] = convertToDouble( cellblock.Data("length_a") );
cif_box[1] = convertToDouble( cellblock.Data("length_b") );
cif_box[2] = convertToDouble( cellblock.Data("length_c") );
//.........这里部分代码省略.........
示例6: SequenceAlign
//.........这里部分代码省略.........
else
QmaskExp.append( "," + integerToString(Smap[sres]+qmaskoffset) );
}
if (State.Debug() > 0) mprintf(" %3s\n", qres);
}
mprintf("Smask: %s\n", SmaskExp.c_str());
mprintf("Qmask: %s\n", QmaskExp.c_str());
// Check that query residues match reference.
for (unsigned int sres = 0; sres != Sbjct.size(); sres++) {
int qres = Smap[sres];
if (qres != -1) {
if (Query[qres] != qref.Parm().Res(qres).SingleCharName()) {
mprintf("Warning: Potential residue mismatch: Query %s reference %s\n",
Residue::ConvertResName(Query[qres]), qref.Parm().Res(qres).c_str());
}
}
}
// Build subject using coordinate from reference.
//AtomMask sMask; // Contain atoms that should be in sTop
Topology sTop;
Frame sFrame;
Iarray placeHolder; // Atom indices of placeholder residues.
for (unsigned int sres = 0; sres != Sbjct.size(); sres++) {
int qres = Smap[sres];
NameType SresName( Residue::ConvertResName(Sbjct[sres]) );
if (qres != -1) {
Residue const& QR = qref.Parm().Res(qres);
Residue SR(SresName, sres+1, ' ', QR.ChainID());
if (Query[qres] == Sbjct[sres]) { // Exact match. All non-H atoms.
for (int qat = QR.FirstAtom(); qat != QR.LastAtom(); qat++)
{
if (qref.Parm()[qat].Element() != Atom::HYDROGEN)
sTop.AddTopAtom( qref.Parm()[qat], SR );
sFrame.AddXYZ( qref.Coord().XYZ(qat) );
//sMask.AddAtom(qat);
}
} else { // Partial match. Copy only backbone and CB.
for (int qat = QR.FirstAtom(); qat != QR.LastAtom(); qat++)
{
if ( qref.Parm()[qat].Name().Match("N" ) ||
qref.Parm()[qat].Name().Match("CA") ||
qref.Parm()[qat].Name().Match("CB") ||
qref.Parm()[qat].Name().Match("C" ) ||
qref.Parm()[qat].Name().Match("O" ) )
{
sTop.AddTopAtom( qref.Parm()[qat], SR );
sFrame.AddXYZ( qref.Coord().XYZ(qat) );
}
}
}
} else {
// Residue in query does not exist for subject. Just put placeholder CA for now.
Vec3 Zero(0.0);
placeHolder.push_back( sTop.Natom() );
sTop.AddTopAtom( Atom("CA", "C "), Residue(SresName, sres+1, ' ', ' ') );
sFrame.AddXYZ( Zero.Dptr() );
}
}
//sTop.PrintAtomInfo("*");
mprintf("\tPlaceholder residue indices:");
for (Iarray::const_iterator p = placeHolder.begin(); p != placeHolder.end(); ++p)
mprintf(" %i", *p + 1);
mprintf("\n");
// Try to give placeholders more reasonable coordinates.
if (!placeHolder.empty()) {
示例7: ReadParm
/** Open the Charmm PSF file specified by filename and set up topology data.
* Mask selection requires natom, nres, names, resnames, resnums.
*/
int Parm_CharmmPsf::ReadParm(FileName const& fname, Topology &parmOut) {
const size_t TAGSIZE = 10;
char tag[TAGSIZE];
tag[0]='\0';
CpptrajFile infile;
if (infile.OpenRead(fname)) return 1;
mprintf(" Reading Charmm PSF file %s as topology file.\n",infile.Filename().base());
// Read the first line, should contain PSF...
const char* buffer = 0;
if ( (buffer=infile.NextLine()) == 0 ) return 1;
// Advance to <ntitle> !NTITLE
int ntitle = FindTag(tag, "!NTITLE", 7, infile);
// Only read in 1st title. Skip any asterisks.
std::string psftitle;
if (ntitle > 0) {
buffer = infile.NextLine();
const char* ptr = buffer;
while (*ptr != '\0' && (*ptr == ' ' || *ptr == '*')) ++ptr;
psftitle.assign( ptr );
}
parmOut.SetParmName( NoTrailingWhitespace(psftitle), infile.Filename() );
// Advance to <natom> !NATOM
int natom = FindTag(tag, "!NATOM", 6, infile);
if (debug_>0) mprintf("\tPSF: !NATOM tag found, natom=%i\n", natom);
// If no atoms, probably issue with PSF file
if (natom < 1) {
mprinterr("Error: No atoms in PSF file.\n");
return 1;
}
// Read the next natom lines
int psfresnum = 0;
char psfresname[6];
char psfname[6];
char psftype[6];
double psfcharge;
double psfmass;
for (int atom=0; atom < natom; atom++) {
if ( (buffer=infile.NextLine()) == 0 ) {
mprinterr("Error: ReadParmPSF(): Reading atom %i\n",atom+1);
return 1;
}
// Read line
// ATOM# SEGID RES# RES ATNAME ATTYPE CHRG MASS (REST OF COLUMNS ARE LIKELY FOR CMAP AND CHEQ)
sscanf(buffer,"%*i %*s %i %s %s %s %lf %lf",&psfresnum, psfresname,
psfname, psftype, &psfcharge, &psfmass);
parmOut.AddTopAtom( Atom( psfname, psfcharge, psfmass, psftype),
Residue( psfresname, psfresnum, ' ', ' '), 0 );
} // END loop over atoms
// Advance to <nbond> !NBOND
int bondatoms[9];
int nbond = FindTag(tag, "!NBOND", 6, infile);
if (nbond > 0) {
if (debug_>0) mprintf("\tPSF: !NBOND tag found, nbond=%i\n", nbond);
int nlines = nbond / 4;
if ( (nbond % 4) != 0) nlines++;
for (int bondline=0; bondline < nlines; bondline++) {
if ( (buffer=infile.NextLine()) == 0 ) {
mprinterr("Error: ReadParmPSF(): Reading bond line %i\n",bondline+1);
return 1;
}
// Each line has 4 pairs of atom numbers
int nbondsread = sscanf(buffer,"%i %i %i %i %i %i %i %i",bondatoms,bondatoms+1,
bondatoms+2,bondatoms+3, bondatoms+4,bondatoms+5,
bondatoms+6,bondatoms+7);
// NOTE: Charmm atom nums start from 1
for (int bondidx=0; bondidx < nbondsread; bondidx+=2)
parmOut.AddBond(bondatoms[bondidx]-1, bondatoms[bondidx+1]-1);
}
} else
mprintf("Warning: PSF has no bonds.\n");
// Advance to <nangles> !NTHETA
int nangle = FindTag(tag, "!NTHETA", 7, infile);
if (nangle > 0) {
if (debug_>0) mprintf("\tPSF: !NTHETA tag found, nangle=%i\n", nangle);
int nlines = nangle / 3;
if ( (nangle % 3) != 0) nlines++;
for (int angleline=0; angleline < nlines; angleline++) {
if ( (buffer=infile.NextLine()) == 0) {
mprinterr("Error: Reading angle line %i\n", angleline+1);
return 1;
}
// Each line has 3 groups of 3 atom numbers
int nanglesread = sscanf(buffer,"%i %i %i %i %i %i %i %i %i",bondatoms,bondatoms+1,
bondatoms+2,bondatoms+3, bondatoms+4,bondatoms+5,
bondatoms+6,bondatoms+7, bondatoms+8);
for (int angleidx=0; angleidx < nanglesread; angleidx += 3)
parmOut.AddAngle( bondatoms[angleidx ]-1,
bondatoms[angleidx+1]-1,
bondatoms[angleidx+2]-1 );
}
} else
mprintf("Warning: PSF has no angles.\n");
// Advance to <ndihedrals> !NPHI
int ndihedral = FindTag(tag, "!NPHI", 5, infile);
if (ndihedral > 0) {
if (debug_>0) mprintf("\tPSF: !NPHI tag found, ndihedral=%i\n", ndihedral);
//.........这里部分代码省略.........