Skip to content

L
libcifpp
Commit 51aebf84 by maarten
Options

nieuwe symmetrie code 

git-svn-id: svn+ssh://gitlab/srv/svn-repos/pdb-redo/trunk@495 a1961a4f-ab94-4bcc-80e8-33b5a54de466
parent 38121e20
Showing with 115 additions and 60 deletions
include/cif++/Symmetry.h
...@@ -12,8 +12,8 @@ namespace mmcif ...@@ -12,8 +12,8 @@ namespace mmcif
clipper::Coord_orth CalculateOffsetForCell(const Structure& p, const clipper::Spacegroup& spacegroup, const clipper::Cell& cell); clipper::Coord_orth CalculateOffsetForCell(const Structure& p, const clipper::Spacegroup& spacegroup, const clipper::Cell& cell);
std::vector<clipper::RTop_orth> AlternativeSites(const clipper::Spacegroup& spacegroup, const clipper::Cell& cell); std::vector<clipper::RTop_orth> AlternativeSites(const clipper::Spacegroup& spacegroup, const clipper::Cell& cell);
int GetSpacegroupNumber(const std::string& spacegroup); // alternative for clipper's parsing code int GetSpacegroupNumber(std::string spacegroup); // alternative for clipper's parsing code
std::string SpacegroupToxHM(std::string spacegroup); // std::string SpacegroupToHall(std::string spacegroup);
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// To iterate over all near symmetry copies of an atom // To iterate over all near symmetry copies of an atom
...@@ -23,19 +23,14 @@ class SymmetryAtomIteratorFactory ...@@ -23,19 +23,14 @@ class SymmetryAtomIteratorFactory
public: public:
// SymmetryAtomIteratorFactory(const Structure& p); // SymmetryAtomIteratorFactory(const Structure& p);
SymmetryAtomIteratorFactory(const Structure& p, const clipper::Spacegroup& spacegroup, const clipper::Cell& cell) // SymmetryAtomIteratorFactory(const Structure& p, const clipper::Spacegroup& spacegroup, const clipper::Cell& cell)
: mSpacegroupNr(spacegroup.spacegroup_number()) // : mSpacegroupNr(spacegroup.spacegroup_number())
, mSpacegroup(spacegroup) // , mSpacegroup(spacegroup)
, mD(CalculateOffsetForCell(p, spacegroup, cell)) // , mD(CalculateOffsetForCell(p, spacegroup, cell))
, mRtOrth(AlternativeSites(spacegroup, cell)) // , mRtOrth(AlternativeSites(spacegroup, cell))
, mCell(cell) {} // , mCell(cell) {}
SymmetryAtomIteratorFactory(const Structure& p, int spacegroupNr, const clipper::Cell& cell) SymmetryAtomIteratorFactory(const Structure& p, int spacegroupNr, const clipper::Cell& cell);
: mSpacegroupNr(spacegroupNr)
, mSpacegroup(GetCCP4SpacegroupDescr(spacegroupNr))
, mD(CalculateOffsetForCell(p, mSpacegroup, cell))
, mRtOrth(AlternativeSites(mSpacegroup, cell))
, mCell(cell) {}
SymmetryAtomIteratorFactory(const SymmetryAtomIteratorFactory&) = delete; SymmetryAtomIteratorFactory(const SymmetryAtomIteratorFactory&) = delete;
SymmetryAtomIteratorFactory& operator=(const SymmetryAtomIteratorFactory&) = delete; SymmetryAtomIteratorFactory& operator=(const SymmetryAtomIteratorFactory&) = delete;
......
src/Symmetry.cpp
...@@ -206,8 +206,13 @@ int32_t GetRotationalIndexNumber(int spacegroup, const clipper::RTop_frac& rt) ...@@ -206,8 +206,13 @@ int32_t GetRotationalIndexNumber(int spacegroup, const clipper::RTop_frac& rt)
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// And unfortunately, clipper does not know all the spacegroup names mentioned in symop.lib // And unfortunately, clipper does not know all the spacegroup names mentioned in symop.lib
int GetSpacegroupNumber(const std::string& spacegroup) int GetSpacegroupNumber(std::string spacegroup)
{ {
if (spacegroup == "P 21 21 2 A")
spacegroup = "P 21 21 2 (a)";
else if (spacegroup.empty())
throw runtime_error("No spacegroup, cannot continue");
int result = 0; int result = 0;
const size_t N = sizeof(kSpaceGroups) / sizeof(Spacegroup); const size_t N = sizeof(kSpaceGroups) / sizeof(Spacegroup);
...@@ -237,57 +242,98 @@ int GetSpacegroupNumber(const std::string& spacegroup) ...@@ -237,57 +242,98 @@ int GetSpacegroupNumber(const std::string& spacegroup)
// ----------------------------------------------------------------------- // -----------------------------------------------------------------------
std::string SpacegroupToxHM(std::string spacegroup) std::string SpacegroupToHall(std::string spacegroup)
{ {
if (spacegroup == "P 1-") if (spacegroup == "P 21 21 2 A")
spacegroup = "P -1";
else if (spacegroup == "P 21 21 2 A")
spacegroup = "P 21 21 2 (a)"; spacegroup = "P 21 21 2 (a)";
else if (spacegroup.empty()) else if (spacegroup.empty())
throw runtime_error("No spacegroup, cannot continue"); throw runtime_error("No spacegroup, cannot continue");
int nr = GetSpacegroupNumber(spacegroup); string result;
return kSpaceGroups[nr].xHM;
}
clipper::Spgr_descr GetCCP4SpacegroupDescr(int nr) const size_t N = sizeof(kSpaceGroups) / sizeof(Spacegroup);
{
const size_t N = sizeof(kSymopNrTable) / sizeof(SymopDataBlock);
int32_t L = 0, R = static_cast<int32_t>(N - 1); int32_t L = 0, R = static_cast<int32_t>(N - 1);
while (L <= R) while (L <= R)
{ {
int32_t i = (L + R) / 2; int32_t i = (L + R) / 2;
if (kSymopNrTable[i].spacegroupNr < nr)
int d = spacegroup.compare(kSpaceGroups[i].name);
if (d > 0)
L = i + 1; L = i + 1;
else else if (d < 0)
R = i - 1; R = i - 1;
else
{
result = kSpaceGroups[i].Hall;
break;
}
} }
if (L < 0 or L >= N) if (result.empty())
throw runtime_error("Invalid spacegroup number"); throw runtime_error("Spacegroup name " + spacegroup + " was not found in table");
return result;
}
clipper::Spgr_descr::Symop_codes symop_codes; // clipper::Spgr_descr GetCCP4SpacegroupDescr(int nr)
// {
// const size_t N = sizeof(kSymopNrTable) / sizeof(SymopDataBlock);
// int32_t L = 0, R = static_cast<int32_t>(N - 1);
// while (L <= R)
// {
// int32_t i = (L + R) / 2;
// if (kSymopNrTable[i].spacegroupNr < nr)
// L = i + 1;
// else
// R = i - 1;
// }
// if (L < 0 or L >= N)
// throw runtime_error("Invalid spacegroup number");
// clipper::Spgr_descr::Symop_codes symop_codes;
// for (size_t i = L; i < N and kSymopNrTable[i].spacegroupNr == nr; ++i)
// {
// auto& symop = kSymopNrTable[i];
// clipper::ftype m[4][4] = {
// { symop.rt.rot_0_0, symop.rt.rot_0_1, symop.rt.rot_0_2, static_cast<float>(symop.rt.trn_0_0) / symop.rt.trn_0_1 },
// { symop.rt.rot_1_0, symop.rt.rot_1_1, symop.rt.rot_1_2, static_cast<float>(symop.rt.trn_1_0) / symop.rt.trn_1_1 },
// { symop.rt.rot_2_0, symop.rt.rot_2_1, symop.rt.rot_2_2, static_cast<float>(symop.rt.trn_2_0) / symop.rt.trn_2_1 },
// { 0, 0, 0, 1 }
// };
// symop_codes.emplace_back(clipper::Symop(m));
// }
// return clipper::Spgr_descr(symop_codes);
// }
for (size_t i = L; i < N and kSymopNrTable[i].spacegroupNr == nr; ++i) clipper::Spgr_descr GetCCP4SpacegroupDescr(int nr)
{
for (auto& sg: kSpaceGroups)
{ {
auto& symop = kSymopNrTable[i]; if (sg.nr == nr)
return clipper::Spgr_descr(sg.Hall, clipper::Spgr_descr::Hall);
clipper::ftype m[4][4] = {
{ symop.rt.rot_0_0, symop.rt.rot_0_1, symop.rt.rot_0_2, static_cast<float>(symop.rt.trn_0_0) / symop.rt.trn_0_1 },
{ symop.rt.rot_1_0, symop.rt.rot_1_1, symop.rt.rot_1_2, static_cast<float>(symop.rt.trn_1_0) / symop.rt.trn_1_1 },
{ symop.rt.rot_2_0, symop.rt.rot_2_1, symop.rt.rot_2_2, static_cast<float>(symop.rt.trn_2_0) / symop.rt.trn_2_1 },
{ 0, 0, 0, 1 }
};
symop_codes.emplace_back(clipper::Symop(m));
} }
return clipper::Spgr_descr(symop_codes); throw runtime_error("Invalid spacegroup number: " + to_string(nr));
} }
// ----------------------------------------------------------------------- // -----------------------------------------------------------------------
SymmetryAtomIteratorFactory::SymmetryAtomIteratorFactory(const Structure& p, int spacegroupNr, const clipper::Cell& cell)
: mSpacegroupNr(spacegroupNr)
, mSpacegroup(GetCCP4SpacegroupDescr(spacegroupNr))
, mD(CalculateOffsetForCell(p, mSpacegroup, cell))
, mRtOrth(AlternativeSites(mSpacegroup, cell))
, mCell(cell)
{
}
string SymmetryAtomIteratorFactory::symop_mmcif(const Atom& a) const string SymmetryAtomIteratorFactory::symop_mmcif(const Atom& a) const
{ {
string result; string result;
......
tools/symop-map-generator.cpp
...@@ -157,7 +157,7 @@ class SymopParser ...@@ -157,7 +157,7 @@ class SymopParser
int m_trn[3][2] = {}; int m_trn[3][2] = {};
}; };
map<string,string> get_xhm_map(string path) map<string,string> get_Hall_map(string path)
{ {
ifstream file(path); ifstream file(path);
if (not file.is_open()) if (not file.is_open())
...@@ -166,9 +166,10 @@ map<string,string> get_xhm_map(string path) ...@@ -166,9 +166,10 @@ map<string,string> get_xhm_map(string path)
enum class State { skip, spacegroup } state = State::skip; enum class State { skip, spacegroup } state = State::skip;
string line; string line;
string old, xHM; string Hall;
vector<string> old;
const regex rx(R"(^symbol +(xHM|old) +'(.*)'$)"); const regex rx(R"(^symbol +(Hall|old) +'(.+?)'(?: +'(.+?)')?$)");
map<string,string> result; map<string,string> result;
...@@ -186,21 +187,26 @@ map<string,string> get_xhm_map(string path) ...@@ -186,21 +187,26 @@ map<string,string> get_xhm_map(string path)
smatch m; smatch m;
if (regex_match(line, m, rx)) if (regex_match(line, m, rx))
{ {
cerr << "match " << m[1] << " => " << m[2] << endl;
if (m[1] == "old") if (m[1] == "old")
old = m[2]; {
else if (m[1] == "xHM") old.push_back(m[2]);
xHM = m[2]; if (m[3].matched)
old.push_back(m[3]);
}
else
Hall = m[2];
} }
else if (line == "end_spacegroup") else if (line == "end_spacegroup")
{ {
if (not (old.empty() or xHM.empty())) if (not Hall.empty())
result.emplace(old, xHM); {
for (auto& o: old)
result.emplace(o, Hall);
}
state = State::skip; state = State::skip;
old.clear(); old.clear();
xHM.clear(); Hall.clear();
} }
break; break;
} }
...@@ -293,7 +299,7 @@ int main() ...@@ -293,7 +299,7 @@ int main()
// ----------------------------------------------------------------------- // -----------------------------------------------------------------------
auto xhm_map = get_xhm_map(CLIBD + "/syminfo.lib"s); auto Hall_map = get_Hall_map(CLIBD + "/syminfo.lib"s);
// -------------------------------------------------------------------- // --------------------------------------------------------------------
...@@ -307,7 +313,7 @@ int main() ...@@ -307,7 +313,7 @@ int main()
struct Spacegroup struct Spacegroup
{ {
const char* name; const char* name;
const char* xHM; const char* Hall;
int nr; int nr;
} kSpaceGroups[] = } kSpaceGroups[] =
{ {
...@@ -315,11 +321,19 @@ struct Spacegroup ...@@ -315,11 +321,19 @@ struct Spacegroup
for (auto [name, nr]: spacegroups) for (auto [name, nr]: spacegroups)
{ {
string Hall = Hall_map[name];
name = '"' + name + '"' + string(20 - name.length(), ' '); name = '"' + name + '"' + string(20 - name.length(), ' ');
string xHM = xhm_map[name];
xHM = '"' + xHM + '"' + string(20 - xHM.length(), ' ');
cout << "\t{ " << name << ", " << xHM << ", " << nr << " }," << endl; for (string::size_type p = Hall.length(); p > 0; --p)
{
if (Hall[p - 1] == '"')
Hall.insert(p - 1, "\\", 1);
}
Hall = '"' + Hall + '"' + string(40 - Hall.length(), ' ');
cout << "\t{ " << name << ", " << Hall << ", " << nr << " }," << endl;
} }
cout << R"( cout << R"(
......
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