Commit 19714ecb by maarten

platonyzer work, adding a symop table generator

git-svn-id: svn+ssh://gitlab/srv/svn-repos/pdb-redo/trunk@483 a1961a4f-ab94-4bcc-80e8-33b5a54de466
parent 9fbd41ae
......@@ -61,7 +61,9 @@ class SymmetryAtomIteratorFactory
// SymmetryAtomIteratorFactory(const Structure& p);
SymmetryAtomIteratorFactory(const Structure& p, const clipper::Spacegroup& spacegroup, const clipper::Cell& cell)
: mD(DistanceMap::CalculateOffsetForCell(p, spacegroup, cell))
, mRtOrth(DistanceMap::AlternativeSites(spacegroup, cell)) {}
, mRtOrth(DistanceMap::AlternativeSites(spacegroup, cell))
, mSpacegroup(spacegroup)
, mCell(cell) {}
SymmetryAtomIteratorFactory(const SymmetryAtomIteratorFactory&) = delete;
SymmetryAtomIteratorFactory& operator=(const SymmetryAtomIteratorFactory&) = delete;
......@@ -156,9 +158,13 @@ class SymmetryAtomIteratorFactory
return SymmetryAtomIteratorRange(*this, a);
}
std::string symop_mmcif(const Atom& a) const;
private:
Point mD; // needed to move atoms to center
std::vector<clipper::RTop_orth> mRtOrth;
clipper::Spacegroup mSpacegroup;
clipper::Cell mCell;
};
}
......@@ -85,6 +85,7 @@ class Atom
Atom symmetryCopy(const Point& d, const clipper::RTop_orth& rt);
bool isSymmetryCopy() const;
std::string symmetry() const;
const clipper::RTop_orth& symop() const;
const Compound& comp() const;
bool isWater() const;
......@@ -112,7 +113,7 @@ class Atom
std::string authAsymId() const;
std::string authSeqId() const;
std::string pdbxAuthInsCode() const;
std::string authAltId() const;
std::string pdbxAuthAltId() const;
std::string labelID() const;// label_comp_id + '_' + label_asym_id + '_' + label_seq_id
std::string pdbID() const; // auth_comp_id + '_' + auth_asym_id + '_' + auth_seq_id + pdbx_PDB_ins_code
......
......@@ -450,4 +450,60 @@ vector<Atom> DistanceMap::near(const Atom& a, float maxDistance) const
return result;
}
// -----------------------------------------------------------------------
string SymmetryAtomIteratorFactory::symop_mmcif(const Atom& a) const
{
string result;
if (not a.isSymmetryCopy())
result = "1_555";
else
{
auto rtop_o = a.symop();
for (int i = 0; i < mSpacegroup.num_symops(); ++i)
{
const auto& symop = mSpacegroup.symop(i);
for (int u: { -1, 0, 1})
for (int v: { -1, 0, 1})
for (int w: { -1, 0, 1})
{
// if (i == 0 and u == 0 and v == 0 and w == 0)
// continue;
auto rtop = clipper::RTop_frac(
symop.rot(), symop.trn() + clipper::Vec3<>(u, v, w)
).rtop_orth(mCell);
if (rtop.rot().equals(rtop_o.rot(), 0.00001) and rtop.trn().equals(rtop_o.trn(), 0.000001))
{
// gotcha
auto rtop_f = rtop.rtop_frac(mCell);
uint32_t t[3] = {
static_cast<uint32_t>(5 + static_cast<int>(rint(rtop_f.trn()[0]))),
static_cast<uint32_t>(5 + static_cast<int>(rint(rtop_f.trn()[1]))),
static_cast<uint32_t>(5 + static_cast<int>(rint(rtop_f.trn()[2])))
};
if (t[0] > 9 or t[1] > 9 or t[2] > 9)
throw runtime_error("Symmetry operation has an out-of-range translation.");
result += to_string(1 + i) + "_"
+ to_string(t[0])
+ to_string(t[1])
+ to_string(t[2]);
return result;
}
}
}
}
return result;
}
}
......@@ -383,6 +383,7 @@ struct AtomImpl
clipper::RTop_orth mRTop;
Point mD;
int32_t mRTix;
};
//Atom::Atom(const File& f, const string& id)
......@@ -530,9 +531,9 @@ string Atom::authAtomId() const
return property<string>("auth_atom_id");
}
string Atom::authAltId() const
string Atom::pdbxAuthAltId() const
{
return property<string>("auth_alt_id");
return property<string>("pdbx_auth_alt_id");
}
string Atom::pdbxAuthInsCode() const
......@@ -589,6 +590,11 @@ string Atom::symmetry() const
return clipper::Symop(mImpl->mRTop).format() + "\n" + mImpl->mRTop.format();
}
const clipper::RTop_orth& Atom::symop() const
{
return mImpl->mRTop;
}
const Compound& Atom::comp() const
{
return mImpl->comp();
......
/*
cd ~/projects/pdb-redo/libcif++/tools/
clang++ -I ~/my-clipper/include -L ~/my-clipper/lib -o symop-map-generator symop-map-generator.cpp -lclipper-core
./symop-map-generator
*/
#include <cassert>
#include <iostream>
#include <fstream>
#include <regex>
#include <cstdlib>
#include <clipper/clipper.h>
using namespace std;
std::regex kNameRx(R"(^(\d+) +(\d+) +(\d+) +(\S+) +(\S+) +(\S+) +'([^']+)'( +'([^']+)')?(?: +!.+)?$)");
class SymopParser
{
public:
SymopParser(const std::string& s)
: m_s(s)
, m_p(s.begin())
, m_e(s.end())
{
m_lookahead = next_char();
}
void parse()
{
parsepart(0);
match(',');
parsepart(1);
match(',');
parsepart(2);
if (m_lookahead != 0)
throw runtime_error("symmetry expression contains more data than expected");
}
private:
char next_char()
{
char result = 0;
while (m_p != m_e)
{
char ch = *m_p++;
if (ch == ' ')
continue;
result = ch;
break;
}
return result;
}
void retract()
{
--m_p;
}
void match(char c)
{
if (m_lookahead != c)
throw runtime_error("Unexpected character " + m_lookahead + " expected " + c);
m_lookahead = next_char();
}
void parsepart(int xyz)
{
if (isdigit(m_lookahead))
{
parserational(xyz);
switch (m_lookahead)
{
case '-':
match('-');
break;
}
}
}
if (m_lookahead == '+')
{
match('+');
}
else if (m_lookahead = '-')
{
}
}
tuple<int,int> parserational();
char m_lookahead;
string m_s;
string::iterator m_p, m_e;
int m_rot[3][3] = {};
int m_trn[3] = {};
};
int main()
{
try
{
const char* CLIBD = getenv("CLIBD");
if (CLIBD == nullptr)
throw runtime_error("CCP4 not sourced");
ifstream file(CLIBD + "/symop.lib"s);
if (not file.is_open())
throw runtime_error("Could not open symop.lib file");
string line;
clipper::Spacegroup spacegroup;
enum class State { Initial, InSpacegroup, Error } state = State::Initial;
while (getline(file, line))
{
if (line.empty())
throw runtime_error("Invalid symop.lib file, contains empty line");
switch (state)
{
case State::Error:
case State::InSpacegroup:
if (line[0] == ' ')
{
if (state == State::Error)
continue;
continue;
}
// fall through
case State::Initial:
{
smatch m;
if (not regex_match(line, m, kNameRx))
{
cerr << line << endl;
throw runtime_error("Name line does not match regular expression");
}
int nr = stoi(m[1].str());
if (nr > 230)
{
state = State::Error;
continue;
}
try
{
spacegroup = clipper::Spacegroup(clipper::Spgr_descr(stoi(m[1].str())));
if (spacegroup.num_symops() != stoi(m[2]))
{
cerr << line << endl
<< "Num of symops do not match: " << spacegroup.num_symops() << " <> " << m[2] << endl;
state = State::Error;
continue;
}
state = State::InSpacegroup;
}
catch (const clipper::Message_fatal& e)
{
cerr << line << endl
<< e.text() << endl;
state = State::Error;
}
break;
}
}
}
}
catch (const exception& ex)
{
cerr << endl
<< "Program terminated due to error:" << endl
<< ex.what() << endl;
}
return 0;
}
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