Skip to content

L
libcifpp
Commit eb4f2f77 by maarten
Options

backup 

git-svn-id: svn+ssh://gitlab/srv/svn-repos/pdb-redo/trunk@228 a1961a4f-ab94-4bcc-80e8-33b5a54de466
parent a793eb4c
Showing with 153 additions and 13 deletions
include/cif++/Structure.h
...@@ -164,6 +164,9 @@ class Residue ...@@ -164,6 +164,9 @@ class Residue
const std::string& asymID() const { return mAsymID; } const std::string& asymID() const { return mAsymID; }
int seqID() const { return mSeqID; } int seqID() const { return mSeqID; }
const std::string& altID() const { return mAltID; } const std::string& altID() const { return mAltID; }
// Is this residue a single entity?
bool isEntity() const;
protected: protected:
...@@ -340,12 +343,15 @@ class Structure ...@@ -340,12 +343,15 @@ class Structure
void removeAtom(Atom& a); void removeAtom(Atom& a);
void swapAtoms(Atom& a1, Atom& a2); // swap the labels for these atoms void swapAtoms(Atom& a1, Atom& a2); // swap the labels for these atoms
void moveAtom(Atom& a, Point p); // move atom to a new location void moveAtom(Atom& a, Point p); // move atom to a new location
void changeResidue(Residue& res, const std::string& newCompound,
const std::vector<std::tuple<std::string,std::string>>& remappedAtoms);
private: private:
friend Polymer; friend Polymer;
friend Residue; friend Residue;
cif::Category& category(const char* name) const; cif::Category& category(const char* name) const;
cif::Datablock& datablock() const;
struct StructureImpl* mImpl; struct StructureImpl* mImpl;
}; };
......
src/Cif++.cpp
...@@ -1808,6 +1808,8 @@ void Row::assign(const string& name, const string& value, bool emplacing) ...@@ -1808,6 +1808,8 @@ void Row::assign(const string& name, const string& value, bool emplacing)
if (oldValue != nullptr and value == oldValue) // no need to update if (oldValue != nullptr and value == oldValue) // no need to update
return; return;
string oldStrValue = oldValue ? oldValue : "";
// check the value // check the value
if (col.mValidator) if (col.mValidator)
(*col.mValidator)(value); (*col.mValidator)(value);
...@@ -1883,11 +1885,11 @@ void Row::assign(const string& name, const string& value, bool emplacing) ...@@ -1883,11 +1885,11 @@ void Row::assign(const string& name, const string& value, bool emplacing)
if (childCat == nullptr) if (childCat == nullptr)
continue; continue;
//#if DEBUG #if DEBUG
//cerr << "fixing linked item " << child->mCategory->mName << '.' << child->mTag << endl; cerr << "fixing linked item " << child->mCategory->mName << '.' << child->mTag << endl;
//#endif #endif
auto rows = childCat->find(Key(child->mTag) == oldValue); auto rows = childCat->find(Key(child->mTag) == oldStrValue);
for (auto& cr: rows) for (auto& cr: rows)
cr.assign(child->mTag, value, false); cr.assign(child->mTag, value, false);
} }
...@@ -2008,9 +2010,9 @@ void Row::swap(const string& name, ItemRow* a, ItemRow* b) ...@@ -2008,9 +2010,9 @@ void Row::swap(const string& name, ItemRow* a, ItemRow* b)
if (childCat == nullptr) if (childCat == nullptr)
continue; continue;
//#if DEBUG #if DEBUG
//cerr << "fixing linked item " << child->mCategory->mName << '.' << child->mTag << endl; cerr << "fixing linked item " << child->mCategory->mName << '.' << child->mTag << endl;
//#endif #endif
if (ai != nullptr) if (ai != nullptr)
{ {
auto rows = childCat->find(Key(child->mTag) == string(ai->mText)); auto rows = childCat->find(Key(child->mTag) == string(ai->mText));
......
src/Compound.cpp
...@@ -231,7 +231,7 @@ bool SubStructuresAreIsomeric( ...@@ -231,7 +231,7 @@ bool SubStructuresAreIsomeric(
break; break;
} }
if (result) if (result and na.id != nb.id)
outMapping.emplace_back(na.id, nb.id); outMapping.emplace_back(na.id, nb.id);
return result; return result;
...@@ -368,6 +368,16 @@ string Compound::formula() const ...@@ -368,6 +368,16 @@ string Compound::formula() const
return result; return result;
} }
float Compound::formulaWeight() const
{
float result = 0;
for (auto r: mAtoms)
result += AtomTypeTraits(r.typeSymbol).weight();
return result;
}
int Compound::charge() const int Compound::charge() const
{ {
float result = 0; float result = 0;
...@@ -409,6 +419,8 @@ string Compound::type() const ...@@ -409,6 +419,8 @@ string Compound::type() const
result = "DNA linking"; result = "DNA linking";
else if (cif::iequals(mGroup, "RNA")) else if (cif::iequals(mGroup, "RNA"))
result = "RNA linking"; result = "RNA linking";
else
result = mGroup;
return result; return result;
} }
...@@ -606,9 +618,9 @@ const Compound* CompoundFactory::create(std::string id) ...@@ -606,9 +618,9 @@ const Compound* CompoundFactory::create(std::string id)
row.get("id", "atom_id_centre", "atom_id_1", row.get("id", "atom_id_centre", "atom_id_1",
"atom_id_2", "atom_id_3", "volume_sign"); "atom_id_2", "atom_id_3", "volume_sign");
if (volumeSign == "negativ") if (volumeSign == "negativ" or volumeSign == "negative")
cc.volumeSign = negativ; cc.volumeSign = negativ;
else if (volumeSign == "positiv") else if (volumeSign == "positiv" or volumeSign == "positive")
cc.volumeSign = positiv; cc.volumeSign = positiv;
else if (volumeSign == "both") else if (volumeSign == "both")
cc.volumeSign = both; cc.volumeSign = both;
......
src/Structure.cpp
...@@ -587,6 +587,18 @@ Atom Residue::atomByID(const string& atomID) const ...@@ -587,6 +587,18 @@ Atom Residue::atomByID(const string& atomID) const
throw runtime_error("Atom with atom_id " + atomID + " not found in residue " + mAsymID + ':' + to_string(mSeqID)); throw runtime_error("Atom with atom_id " + atomID + " not found in residue " + mAsymID + ':' + to_string(mSeqID));
} }
// Residue is a single entity if the atoms for the asym with mAsymID is equal
// to the number of atoms in this residue... hope this is correct....
bool Residue::isEntity() const
{
auto& db = mStructure->datablock();
auto a1 = db["atom_site"].find(cif::Key("label_asym_id") == mAsymID);
auto a2 = atoms();
return a1.size() == a2.size();
}
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// monomer // monomer
...@@ -828,12 +840,60 @@ struct StructureImpl ...@@ -828,12 +840,60 @@ struct StructureImpl
void removeAtom(Atom& a); void removeAtom(Atom& a);
void swapAtoms(Atom& a1, Atom& a2); void swapAtoms(Atom& a1, Atom& a2);
void moveAtom(Atom& a, Point p); void moveAtom(Atom& a, Point p);
void changeResidue(Residue& res, const string& newCompound,
const vector<tuple<string,string>>& remappedAtoms);
void insertCompound(const string& compoundID, bool isEntity);
File* mFile; File* mFile;
uint32 mModelNr; uint32 mModelNr;
AtomView mAtoms; AtomView mAtoms;
}; };
void StructureImpl::insertCompound(const string& compoundID, bool isEntity)
{
auto compound = Compound::create(compoundID);
if (compound == nullptr)
throw runtime_error("Trying to insert unknown compound " + compoundID + " (not found in CCP4 monomers lib)");
cif::Datablock& db = *mFile->impl().mDb;
auto& chemComp = db["chem_comp"];
auto r = chemComp.find(cif::Key("id") == compoundID);
if (r.empty())
{
chemComp.emplace({
{ "id", compoundID },
{ "name", compound->name() },
{ "formula", compound->formula() },
{ "type", compound->type() }
});
}
if (isEntity)
{
auto& pdbxEntityNonpoly = db["pdbx_entity_nonpoly"];
if (pdbxEntityNonpoly.find(cif::Key("comp_id") == compoundID).empty())
{
auto& entity = db["entity"];
string entityID = to_string(entity.size() + 1);
entity.emplace({
{ "id", entityID },
{ "type", "non-polymer" },
{ "pdbx_description", compound->name() },
{ "formula_weight", compound->formulaWeight() }
});
pdbxEntityNonpoly.emplace({
{ "entity_id", entityID },
{ "name", compound->name() },
{ "comp_id", compoundID }
});
}
}
}
void StructureImpl::removeAtom(Atom& a) void StructureImpl::removeAtom(Atom& a)
{ {
cif::Datablock& db = *mFile->impl().mDb; cif::Datablock& db = *mFile->impl().mDb;
...@@ -883,6 +943,55 @@ void StructureImpl::moveAtom(Atom& a, Point p) ...@@ -883,6 +943,55 @@ void StructureImpl::moveAtom(Atom& a, Point p)
a.location(p); a.location(p);
} }
void StructureImpl::changeResidue(Residue& res, const string& newCompound,
const vector<tuple<string,string>>& remappedAtoms)
{
cif::Datablock& db = *mFile->impl().mDb;
string entityID;
// First make sure the compound is already known or insert it.
// And if the residue is an entity, we must make sure it exists
insertCompound(newCompound, res.isEntity());
auto& atomSites = db["atom_site"];
auto atoms = res.atoms();
for (auto& a: remappedAtoms)
{
string a1, a2;
tie(a1, a2) = a;
auto i = find_if(atoms.begin(), atoms.end(), [&](const Atom& a) { return a.labelAtomId() == a1; });
if (i == atoms.end())
{
cerr << "Missing atom for atom ID " << a1 << endl;
continue;
}
auto r = atomSites.find(cif::Key("id") == i->id());
if (r.size() != 1)
continue;
if (a1 != a2)
r.front()["label_atom_id"] = a2;
}
for (auto a: atoms)
{
auto r = atomSites.find(cif::Key("id") == a.id());
assert(r.size() == 1);
if (r.size() != 1)
continue;
r.front()["label_comp_id"] = newCompound;
if (not entityID.empty())
r.front()["label_entity_id"] = entityID;
}
}
Structure::Structure(File& f, uint32 modelNr) Structure::Structure(File& f, uint32 modelNr)
: mImpl(new StructureImpl(*this, f, modelNr)) : mImpl(new StructureImpl(*this, f, modelNr))
{ {
...@@ -902,7 +1011,7 @@ AtomView Structure::waters() const ...@@ -902,7 +1011,7 @@ AtomView Structure::waters() const
{ {
AtomView result; AtomView result;
auto& db = *getFile().impl().mDb; auto& db = datablock();
// Get the entity id for water // Get the entity id for water
auto& entityCat = db["entity"]; auto& entityCat = db["entity"];
...@@ -983,7 +1092,7 @@ File& Structure::getFile() const ...@@ -983,7 +1092,7 @@ File& Structure::getFile() const
cif::Category& Structure::category(const char* name) const cif::Category& Structure::category(const char* name) const
{ {
auto& db = *getFile().impl().mDb; auto& db = datablock();
return db[name]; return db[name];
} }
...@@ -1034,7 +1143,7 @@ cif::Category& Structure::category(const char* name) const ...@@ -1034,7 +1143,7 @@ cif::Category& Structure::category(const char* name) const
tuple<string,int,string,string> Structure::MapLabelToPDB( tuple<string,int,string,string> Structure::MapLabelToPDB(
const string& asymId, int seqId, const string& monId) const string& asymId, int seqId, const string& monId)
{ {
auto& db = *getFile().impl().mDb; auto& db = datablock();
tuple<string,int,string,string> result; tuple<string,int,string,string> result;
...@@ -1073,6 +1182,11 @@ tuple<string,int,string,string> Structure::MapLabelToPDB( ...@@ -1073,6 +1182,11 @@ tuple<string,int,string,string> Structure::MapLabelToPDB(
return result; return result;
} }
cif::Datablock& Structure::datablock() const
{
return *mImpl->mFile->impl().mDb;
}
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// actions // actions
...@@ -1091,4 +1205,10 @@ void Structure::moveAtom(Atom& a, Point p) ...@@ -1091,4 +1205,10 @@ void Structure::moveAtom(Atom& a, Point p)
mImpl->moveAtom(a, p); mImpl->moveAtom(a, p);
} }
void Structure::changeResidue(Residue& res, const string& newCompound,
const vector<tuple<string,string>>& remappedAtoms)
{
mImpl->changeResidue(res, newCompound, remappedAtoms);
}
} }
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment