Compare commits

...

10 Commits

Author SHA1 Message Date
Maarten L. Hekkelman
46ea0ca930 Clear pdbx_nonpoly_scheme before filling it in reconstruction 2026-04-14 15:49:20 +02:00
Maarten L. Hekkelman
0da8ba85bf right aligning values in cif 2026-04-14 11:22:26 +02:00
Maarten L. Hekkelman
5d3734e2a5 Merge branch 'develop' of github.com:PDB-REDO/libcifpp into develop 2026-04-14 09:35:46 +02:00
Maarten L. Hekkelman
e692ed6c87 Remove runtime lib statements for msvc 2026-04-14 09:35:38 +02:00
Maarten L. Hekkelman
da8e81f871 Merge branch 'develop' of github.com:PDB-REDO/libcifpp into develop 2026-04-14 09:23:51 +02:00
Maarten L. Hekkelman
0aef601bfd Remove spherical_dots, fix for finding fast_float, I hope 2026-04-14 09:23:39 +02:00
Maarten L. Hekkelman
39644c1bd1 Update validation code, added pdbx_item_enumeration check and test case-sensitive when required. 2026-04-13 10:38:43 +02:00
Maarten L. Hekkelman
421758a01c automatic include generation is evil 2026-04-02 09:36:34 +02:00
Maarten L. Hekkelman
a3a0e5cc70 Fix symmetry operations 2026-03-31 15:48:21 +02:00
Maarten L. Hekkelman
9268aa3cfc Remove debug statement 2026-03-31 08:41:04 +02:00
16 changed files with 3810 additions and 522 deletions

View File

@@ -34,7 +34,7 @@ endif()
# set the project name
project(
libcifpp
VERSION 10.0.2
VERSION 10.0.3
LANGUAGES CXX C)
list(PREPEND CMAKE_MODULE_PATH "${CMAKE_CURRENT_SOURCE_DIR}/cmake")
@@ -138,13 +138,6 @@ if(MSVC)
# make msvc standards compliant...
add_compile_options(/permissive- /bigobj)
add_link_options(/NODEFAULTLIB:library)
# This is dubious...
if(BUILD_SHARED_LIBS)
set(CMAKE_MSVC_RUNTIME_LIBRARY "MultiThreaded$<$<CONFIG:Debug>:Debug>DLL")
else()
set(CMAKE_MSVC_RUNTIME_LIBRARY "MultiThreaded$<$<CONFIG:Debug>:Debug>")
endif()
endif()
# Libraries
@@ -367,8 +360,9 @@ else()
endif()
if(NOT STD_CHARCONV_COMPILING)
target_include_directories(cifpp PRIVATE ${FastFloat_INCLUDE_DIRS})
target_compile_definitions(cifpp PRIVATE USE_FAST_FLOAT)
get_target_property(FF_INC_DIR FastFloat::fast_float INTERFACE_INCLUDE_DIRECTORIES)
target_include_directories(cifpp PRIVATE ${FF_INC_DIR})
target_compile_definitions(cifpp PRIVATE USE_FAST_FLOAT)
endif()
if(CMAKE_CXX_COMPILER_ID STREQUAL "AppleClang")

View File

@@ -1,5 +1,11 @@
Version 10.0.3
- Clear pdbx_nonpoly_scheme before filling it in reconstruction
Version 10.0.2
- Fixed regression in reconstruction introduced in 10.0.1
- Fixed symmetry operations
- Added validation for pdbx_item_enumeration as well as
case-sensitive checks for enumerations when needed.
Version 10.0.1
- Fixed some regressions, like assigning to items.

View File

@@ -902,74 +902,4 @@ quaternion align_points(const std::vector<point> &a, const std::vector<point> &b
/// \brief The RMSd for the points in \a a and \a b
double RMSd(const std::vector<point> &a, const std::vector<point> &b);
// --------------------------------------------------------------------
/**
* @brief Helper class to generate evenly divided points on a sphere
*
* We use a fibonacci sphere to calculate even distribution of the dots
*
* @tparam N The number of points on the sphere is 2 * N + 1
*/
template <int N>
class spherical_dots
{
public:
/// \brief the number of points
constexpr static int P = 2 * N * 1;
/// \brief the *weight* of the fibonacci sphere
constexpr static double W = (4 * std::numbers::pi) / P;
/// \brief the internal storage type
using array_type = typename std::array<point, P>;
/// \brief iterator type
using iterator = typename array_type::const_iterator;
/// \brief singleton instance
static spherical_dots &instance()
{
static spherical_dots sInstance;
return sInstance;
}
/// \brief The number of points
[[nodiscard]] std::size_t size() const { return P; }
/// \brief Access a point by index
const point operator[](uint32_t inIx) const { return m_points[inIx]; }
/// \brief iterator pointing to the first point
[[nodiscard]] iterator begin() const { return m_points.begin(); }
/// \brief iterator pointing past the last point
[[nodiscard]] iterator end() const { return m_points.end(); }
/// \brief return the *weight*,
[[nodiscard]] double weight() const { return W; }
spherical_dots()
{
const double
kGoldenRatio = std::numbers::phi;
auto p = m_points.begin();
for (int32_t i = -N; i <= N; ++i)
{
double lat = std::asin((2.0 * i) / P);
double lon = std::fmod(i, kGoldenRatio) * 2 * std::numbers::pi / kGoldenRatio;
p->m_x = std::sin(lon) * std::cos(lat);
p->m_y = std::cos(lon) * std::cos(lat);
p->m_z = std::sin(lat);
++p;
}
}
private:
array_type m_points;
};
} // namespace cif

View File

@@ -196,7 +196,7 @@ class validation_exception : public std::runtime_error
public:
// Constructors
/// @cond
validation_exception(validation_error err)
: validation_exception(make_error_code(err))
{
@@ -337,7 +337,7 @@ struct item_validator
std::string m_item_name; ///< The item name
bool m_mandatory; ///< Flag indicating this item is mandatory
const type_validator *m_type; ///< The type for this item
cif::iset m_enums; ///< If filled, the set of allowed values
std::set<std::string> m_enums; ///< If filled, the set of allowed values
std::string m_default; ///< If filled, a default value for this item
std::string m_category; ///< The category this item_validator belongs to
std::vector<item_alias> m_aliases; ///< The aliases for this item
@@ -502,15 +502,29 @@ class validator
/// @brief Bottleneck function to report an error in validation
void report_error(std::error_code ec, bool fatal = true) const;
/// @brief Bottleneck function to report an error in validation
void report_error(validation_error err, std::string value, std::string_view category,
std::string_view item, bool fatal = true) const
{
report_error(make_error_code(err), value, category, item, fatal);
}
/// @brief Bottleneck function to report an error in validation
void report_error(validation_error err, std::string_view category,
std::string_view item, bool fatal = true) const
{
report_error(make_error_code(err), category, item, fatal);
report_error(make_error_code(err), "", category, item, fatal);
}
/// @brief Bottleneck function to report an error in validation
void report_error(std::error_code ec, std::string_view category,
std::string_view item, bool fatal = true) const
{
report_error(ec, "", category, item, fatal);
}
/// @brief Bottleneck function to report an error in validation
void report_error(std::error_code ec, std::string value, std::string_view category,
std::string_view item, bool fatal = true) const;
/// @brief Write out the audit_conform data for this validator

File diff suppressed because it is too large Load Diff

View File

@@ -749,16 +749,38 @@ void category::set_validator(const validator *v, datablock &db)
bool number = type->m_primitive_type == DDL_PrimitiveType::Numb;
if (number)
continue;
for (auto row = m_head; row != nullptr; row = row->m_next)
{
if (cix >= row->size() or row->operator[](cix).empty())
continue;
for (auto row = m_head; row != nullptr; row = row->m_next)
{
if (cix >= row->size() or row->operator[](cix).empty())
continue;
item_value &v = row->operator[](cix);
if (v.is_number())
v = v.str();
item_value &v = row->operator[](cix);
if (not v.is_number())
{
// Try cast the value to a number and throw in case of failure
try
{
v.cast_to_int();
}
catch (...)
{
v.cast_to_float();
}
}
}
}
else
{
for (auto row = m_head; row != nullptr; row = row->m_next)
{
if (cix >= row->size() or row->operator[](cix).empty())
continue;
item_value &v = row->operator[](cix);
if (v.is_number())
v = v.str();
}
}
}
@@ -890,11 +912,13 @@ bool category::is_valid() const
seen = true;
std::error_code ec;
iv->validate_value(*ri->get(cix), ec);
auto &v = *ri->get(cix);
iv->validate_value(v, ec);
if (ec != std::errc{})
{
m_validator->report_error(ec, m_name, m_items[cix].m_name, false);
m_validator->report_error(ec, v.str(), m_name, m_items[cix].m_name, false);
continue;
}
}
@@ -1450,7 +1474,12 @@ void category::update_value(const std::vector<row_handle> &rows, std::string_vie
std::error_code ec;
col.m_validator->validate_value(value, ec);
if (ec)
if (ec == cif::make_error_code(cif::validation_error::value_is_not_in_enumeration_list))
{
if (cif::VERBOSE >= 0)
m_validator->report_error(ec, m_name, item_name, false);
}
else if (ec)
throw validation_exception(ec, m_name, item_name);
}
}
@@ -1577,7 +1606,18 @@ void category::update_value(row *row, uint16_t item, item_value value, bool upda
// check the value
if (col.m_validator and validate)
col.m_validator->validate_value(value);
{
std::error_code ec;
col.m_validator->validate_value(value, ec);
if (ec == cif::make_error_code(cif::validation_error::value_is_not_in_enumeration_list))
{
if (cif::VERBOSE >= 0)
m_validator->report_error(ec, m_name, m_items[item].m_name, false);
}
else if (ec)
throw validation_exception(ec, m_name, m_items[item].m_name);
}
// If the item is part of the Key for this category, remove it from the index
// before updating
@@ -1789,6 +1829,12 @@ category::iterator category::insert_impl(const_iterator pos, row *n)
}
else if (ec == cif::make_error_code(cif::validation_error::value_is_not_a_char_string))
v->cast_to_string();
else if (ec == cif::make_error_code(cif::validation_error::value_is_not_in_enumeration_list))
{
if (cif::VERBOSE >= 0)
m_validator->report_error(ec, m_name, m_items[ix].m_name, false);
continue;
}
else
throw std::system_error(ec, "Attempt to insert invalid value");
@@ -2158,7 +2204,7 @@ void category::write_cif(std::ostream &os, const std::vector<uint16_t> &order, b
offset = 0;
}
offset = detail::write_value(os, s, offset, w, /* right_aligned[cix] */ iv->is_number());
offset = detail::write_value(os, s, offset, w, right_aligned[cix]/* iv->is_number() */);
if (offset > 132)
{

View File

@@ -25,6 +25,8 @@
*/
#include "cif++/cif++.hpp"
#include "cif++/text.hpp"
#include "cif++/validate.hpp"
#include <cstddef>
#include <exception>
@@ -103,7 +105,19 @@ class dictionary_parser : public parser
error("Undefined category '" + iv.first);
for (auto &v : iv.second)
{
// enums, make lower case if needed
auto tv = v.m_type;
if (tv and tv->m_primitive_type == DDL_PrimitiveType::UChar)
{
std::set<std::string> es;
for (auto &e : v.m_enums)
es.emplace(cif::to_lower_copy(e));
std::swap(es, v.m_enums);
}
const_cast<category_validator *>(cv)->add_item_validator(std::move(v));
}
}
// check all item validators for having a typeValidator
@@ -275,9 +289,11 @@ class dictionary_parser : public parser
if (typeCode.has_value())
tv = m_validator.get_validator_for_type(*typeCode);
iset ess;
std::set<std::string> ess;
for (auto e : dict["item_enumeration"])
ess.insert(e["value"].get<std::string>());
for (auto e : dict["pdbx_item_enumeration"])
ess.insert(e["value"].get<std::string>());
std::string defaultValue;
if (auto &cat = dict["item_default"]; not cat.empty())

View File

@@ -2178,8 +2178,8 @@ std::string structure::create_non_poly(const std::string &entity_id, const std::
{ "label_comp_id", comp_id },
{ "label_asym_id", asym_id },
{ "label_entity_id", entity_id },
{ "label_seq_id", nullptr },
{ "pdbx_PDB_ins_code", "" },
{ "label_seq_id", cif::item_value_type::INAPPLICABLE },
{ "pdbx_PDB_ins_code", cif::item_value_type::MISSING },
{ "Cartn_x", atom.get_property_value("Cartn_x") },
{ "Cartn_y", atom.get_property_value("Cartn_y") },
{ "Cartn_z", atom.get_property_value("Cartn_z") },
@@ -2190,7 +2190,7 @@ std::string structure::create_non_poly(const std::string &entity_id, const std::
{ "auth_comp_id", comp_id },
{ "auth_asym_id", asym_id },
{ "auth_atom_id", atom.get_property_value("label_atom_id") },
{ "pdbx_PDB_model_num", 1 } });
{ "pdbx_PDB_model_num", m_model_nr } });
auto &newAtom = emplace_atom(std::make_shared<atom::atom_impl>(m_db, atom_id));
res.add_atom(newAtom);
@@ -2208,7 +2208,7 @@ std::string structure::create_non_poly(const std::string &entity_id, const std::
{ "pdb_mon_id", comp_id },
{ "auth_mon_id", comp_id },
{ "pdb_strand_id", asym_id },
{ "pdb_ins_code", nullptr },
{ "pdb_ins_code", cif::item_value_type::MISSING },
});
return asym_id;
@@ -2244,12 +2244,14 @@ std::string structure::create_non_poly(const std::string &entity_id, std::vector
atom.set_value_if_empty({ "group_PDB", "HETATM" });
atom.set_value_if_empty({ "label_comp_id", comp_id });
atom.set_value_if_empty({ "label_seq_id", nullptr });
atom.set_value_if_empty({ "label_seq_id", cif::item_value_type::INAPPLICABLE });
atom.set_value_if_empty({ "auth_comp_id", comp_id });
atom.set_value_if_empty({ "auth_seq_id", "1" });
atom.set_value_if_empty({ "pdbx_PDB_model_num", 1 });
atom.set_value_if_empty({ "label_alt_id", "" });
atom.set_value_if_empty({ "pdbx_PDB_model_num", m_model_nr });
atom.set_value_if_empty({ "label_alt_id", cif::item_value_type::INAPPLICABLE });
atom.set_value_if_empty({ "occupancy", { 1.0, 2 } });
atom.set_value_if_empty({ "pdbx_formal_charge", cif::item_value_type::MISSING });
atom.set_value_if_empty({ "pdbx_tls_group_id", cif::item_value_type::MISSING });
auto row = atom_site.emplace(atom.begin(), atom.end());
@@ -2299,7 +2301,7 @@ std::string structure::create_non_poly(const std::string &compound_id, bool skip
{ "Cartn_x", ax },
{ "Cartn_y", ay },
{ "Cartn_z", az },
{ "B_iso_or_equiv", 30.00 } });
{ "B_iso_or_equiv", { 30.00, 2 } } });
}
return create_non_poly(create_non_poly_entity(compound_id), atoms);

View File

@@ -2193,7 +2193,7 @@ void PDBFileParser::ParseRemarks()
if (std::regex_match(line, m, rx))
{
models[0] = std::stoi(m[1].str());
models[1] = stoi(m[2].str());
models[1] = std::stoi(m[2].str());
}
else
headerSeen = cif::contains(line, "RES C SSSEQI");

View File

@@ -1340,6 +1340,8 @@ void createPdbxNonpolyScheme(datablock &db)
auto &pdbx_nonpoly_scheme = db["pdbx_nonpoly_scheme"];
auto &atom_site = db["atom_site"];
pdbx_nonpoly_scheme.clear();
for (const auto &[entity_id, comp_id] : pdbx_entity_nonpoly.rows<std::string, std::string>("entity_id", "comp_id"))
{
for (int ndb_nr = 1; auto row : atom_site.find("label_entity_id"_key == entity_id and "label_comp_id"_key == comp_id))

View File

@@ -24,6 +24,7 @@
* SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
#include "cif++/symmetry.hpp"
#include "cif++/cif++.hpp"
#include "symop_table_data.hpp"
@@ -292,12 +293,7 @@ point spacegroup::operator()(const point &pt, const cell &c, sym_op symop) const
t.m_translation.m_y += symop.m_tb - 5;
t.m_translation.m_z += symop.m_tc - 5;
auto fpt = fractional(pt, c);
auto o = offsetToOriginFractional(fpt);
auto spt = t(fpt + o) - o;
return orthogonal(spt, c);
return orthogonal(t(fractional(pt, c)), c);
}
point spacegroup::inverse(const point &pt, const cell &c, sym_op symop) const
@@ -311,13 +307,8 @@ point spacegroup::inverse(const point &pt, const cell &c, sym_op symop) const
t.m_translation.m_y += symop.m_tb - 5;
t.m_translation.m_z += symop.m_tc - 5;
auto fpt = fractional(pt, c);
auto o = offsetToOriginFractional(fpt);
auto it = cif::inverse(t);
auto spt = it(fpt + o) - o;
return orthogonal(spt, c);
return orthogonal(it(fractional(pt, c)), c);
}
// --------------------------------------------------------------------
@@ -462,9 +453,9 @@ std::tuple<float, point, sym_op> crystal::closest_symmetry_copy(point a, point b
a = orthogonal(fa, m_cell);
for (std::size_t i = 0; i < m_spacegroup.size(); ++i)
for (uint8_t i = 0; std::cmp_less(i, m_spacegroup.size()); ++i)
{
sym_op s(static_cast<uint8_t>(i + 1));
sym_op s(i + 1);
auto &t = m_spacegroup[i];
auto fsb = t(fb);

View File

@@ -27,6 +27,7 @@
#include "cif++/validate.hpp"
#include "cif++/cif++.hpp"
#include "cif++/text.hpp"
#include <cassert>
#include <charconv>
@@ -245,8 +246,17 @@ bool item_validator::validate_value(const item_value &value, std::error_code &ec
ec = make_error_code(validation_error::value_does_not_match_rx);
}
if (ec == std::errc{} and not m_enums.empty() and m_enums.count(value.str()) == 0)
ec = make_error_code(validation_error::value_is_not_in_enumeration_list);
if (ec == std::errc{} and not m_enums.empty())
{
bool valid =
m_type->m_primitive_type == DDL_PrimitiveType::UChar ? //
m_enums.contains(cif::to_lower_copy(value.sv()))
: //
m_enums.contains(std::string{ value.sv() });
if (not valid)
ec = make_error_code(validation_error::value_is_not_in_enumeration_list);
}
}
}
}
@@ -436,8 +446,8 @@ void validator::report_error(std::error_code ec, bool fatal) const
std::cerr << ec.message() << '\n';
}
void validator::report_error(std::error_code ec, std::string_view category,
std::string_view item, bool fatal) const
void validator::report_error(std::error_code ec, std::string value,
std::string_view category, std::string_view item, bool fatal) const
{
if (m_strict or fatal)
{
@@ -448,10 +458,19 @@ void validator::report_error(std::error_code ec, std::string_view category,
}
if (VERBOSE > 0)
std::cerr << ec.message()
<< "; category: " << std::quoted(category)
<< " item: " << std::quoted(item)
<< '\n';
{
if (value.empty())
std::cerr << ec.message()
<< "; category: " << std::quoted(category)
<< " item: " << std::quoted(item)
<< '\n';
else
std::cerr << ec.message()
<< "; value: " << std::quoted(value)
<< "; category: " << std::quoted(category)
<< " item: " << std::quoted(item)
<< '\n';
}
}
void validator::fill_audit_conform(category &audit_conform) const

BIN
test/476d.cif.gz Normal file

Binary file not shown.

View File

@@ -443,8 +443,6 @@ TEST_CASE("test_alternates_1")
const std::filesystem::path example(gTestDir / ".." / "examples" / "1cbs.cif.gz");
cif::file file(example.string());
auto &db = file.front();
cif::mm::structure s(file);
for (auto atom : s.atoms())

View File

@@ -41,10 +41,6 @@ TEST_CASE("reconstruct")
for (std::filesystem::directory_iterator i(gTestDir / "reconstruct"); i != std::filesystem::directory_iterator{}; ++i)
{
if (i->path().filename().string() != "9d3j_fixDMC.cif")
continue;
std::cout << i->path() << '\n';
if (i->path().extension() == ".pdb")

View File

@@ -25,6 +25,7 @@
*/
#include "cif++/cif++.hpp"
#include "cif++/symmetry.hpp"
#include "test-main.hpp"
#include <catch2/catch_test_macros.hpp>
@@ -35,7 +36,7 @@
# pragma warning(disable : 4127) // conditional expression is constant
#endif
// #include <Eigen/Eigenvalues>
#include <Eigen/Eigen>
// --------------------------------------------------------------------
@@ -229,102 +230,31 @@ TEST_CASE("dh_q_1")
}
}
/* TEST_CASE("m2q_0a")
{
for (std::size_t i = 0; i < cif::kSymopNrTableSize; ++i)
{
auto d = cif::kSymopNrTable[i].symop().data();
Eigen::Matrix3f rot;
rot << static_cast<float>(d[0]), static_cast<float>(d[1]), static_cast<float>(d[2]), static_cast<float>(d[3]), static_cast<float>(d[4]), static_cast<float>(d[5]), static_cast<float>(d[6]), static_cast<float>(d[7]), static_cast<float>(d[8]);
// check to see if this matrix contains a true rotation
if (rot * rot.transpose() != Eigen::Matrix3f::Identity() or rot.determinant() != 1)
continue;
Eigen::Quaternionf qe(rot);
auto q = normalize(cif::quaternion{ qe.w(), qe.x(), qe.y(), qe.z() });
cif::point p1{ 1, 1, 1 };
cif::point p2 = p1;
p2.rotate(q);
cif::matrix3x3<float> rot_c({ static_cast<float>(d[0]),
static_cast<float>(d[1]),
static_cast<float>(d[2]),
static_cast<float>(d[3]),
static_cast<float>(d[4]),
static_cast<float>(d[5]),
static_cast<float>(d[6]),
static_cast<float>(d[7]),
static_cast<float>(d[8]) });
cif::point p3 = rot_c * p1;
CHECK_THAT(p2.m_x, Catch::Matchers::WithinRel(p3.m_x, 0.01f));
CHECK_THAT(p2.m_y, Catch::Matchers::WithinRel(p3.m_y, 0.01f));
CHECK_THAT(p2.m_z, Catch::Matchers::WithinRel(p3.m_z, 0.01f));
}
}
*/
// "TEST_CASE(m2q_1")
// TEST_CASE("m2q_1")
// {
// for (std::size_t i = 0; i < cif::kSymopNrTableSize; ++i)
// {
// auto d = cif::kSymopNrTable[i].symop().data();
// cif::matrix3x3<float> rot;
// float Qxx = rot(0, 0) = d[0];
// float Qxy = rot(0, 1) = d[1];
// float Qxz = rot(0, 2) = d[2];
// float Qyx = rot(1, 0) = d[3];
// float Qyy = rot(1, 1) = d[4];
// float Qyz = rot(1, 2) = d[5];
// float Qzx = rot(2, 0) = d[6];
// float Qzy = rot(2, 1) = d[7];
// float Qzz = rot(2, 2) = d[8];
// Eigen::Matrix3f rot;
// cif::matrix4x4<float> m({
// Qxx - Qyy - Qzz, Qyx + Qxy, Qzx + Qxz, Qzy - Qyz,
// Qyx + Qxy, Qyy - Qxx - Qzz, Qzy + Qyz, Qxz - Qzx,
// Qzx + Qxz, Qzy + Qyz, Qzz - Qxx - Qyy, Qyx - Qxy,
// Qzy - Qyz, Qxz - Qzx, Qyx - Qxy, Qxx + Qyy + Qzz
// });
// rot << d[0], d[1], d[2], d[3], d[4], d[5], d[6], d[7], d[8];
// auto &&[ev, em] = cif::eigen(m * (1/3.0f), false);
// std::size_t bestJ = 0;
// float bestEV = -1;
// for (std::size_t j = 0; j < 4; ++j)
// if (rot * rot.transpose() == Eigen::Matrix3f::Identity() and rot.determinant() == 1)
// {
// if (bestEV < ev[j])
// {
// bestEV = ev[j];
// bestJ = j;
// }
// Eigen::Quaternionf qe(rot);
// auto q = normalize(cif::quaternion{ qe.w(), qe.x(), qe.y(), qe.z() });
// cif::point p1{ 1, 1, 1 };
// cif::point p2 = p1;
// p2.rotate(q);
// auto p3 = rot * Eigen::Vector3f{ p1.m_x, p1.m_y, p1.m_z };
// CHECK_THAT(p2.m_x, Catch::Matchers::WithinRel(p3[0], 0.01f));
// CHECK_THAT(p2.m_y, Catch::Matchers::WithinRel(p3[1], 0.01f));
// CHECK_THAT(p2.m_z, Catch::Matchers::WithinRel(p3[2], 0.01f));
// }
// if (std::abs(bestEV - 1) > 0.01)
// continue; // not a rotation matrix
// auto q = normalize(cif::quaternion{
// static_cast<float>(em(bestJ, 3)),
// static_cast<float>(em(bestJ, 0)),
// static_cast<float>(em(bestJ, 1)),
// static_cast<float>(em(bestJ, 2)) });
// cif::point p1{ 1, 1, 1 };
// cif::point p2 = p1;
// p2.rotate(q);
// cif::point p3 = rot * p1;
// REQUIRE(p2.m_x == p3.m_x);
// REQUIRE(p2.m_y == p3.m_y);
// REQUIRE(p2.m_z == p3.m_z);
// }
// }
@@ -394,35 +324,6 @@ TEST_CASE("symm_4")
// --------------------------------------------------------------------
TEST_CASE("symm_4wvp_1")
{
using namespace cif::literals;
cif::file f(gTestDir / "4wvp.cif.gz");
f.front().set_validator(cif::validator_factory::instance().get("mmcif_pdbx.dic"));
auto &db = f.front();
cif::mm::structure s(db);
cif::crystal c(db);
cif::point p{ -78.722f, 98.528f, 11.994f };
auto a = s.get_residue("A", 10, "").get_atom_by_atom_id("O");
auto sp1 = c.symmetry_copy(a.get_location(), "2_565"_symop);
CHECK_THAT(sp1.m_x, Catch::Matchers::WithinAbs(p.m_x, 0.5f));
CHECK_THAT(sp1.m_y, Catch::Matchers::WithinAbs(p.m_y, 0.5f));
CHECK_THAT(sp1.m_z, Catch::Matchers::WithinAbs(p.m_z, 0.5f));
const auto &[d, sp2, so] = c.closest_symmetry_copy(p, a.get_location());
REQUIRE(d < 1);
CHECK_THAT(sp2.m_x, Catch::Matchers::WithinAbs(p.m_x, 0.5f));
CHECK_THAT(sp2.m_y, Catch::Matchers::WithinAbs(p.m_y, 0.5f));
CHECK_THAT(sp2.m_z, Catch::Matchers::WithinAbs(p.m_z, 0.5f));
}
TEST_CASE("symm_2bi3_1")
{
cif::file f(gTestDir / "2bi3.cif.gz");
@@ -538,6 +439,96 @@ TEST_CASE("symm_3bwh_1")
}
}
TEST_CASE("symm_476d")
{
cif::file f(gTestDir / "476d.cif.gz");
f.front().set_validator(cif::validator_factory::instance().get("mmcif_pdbx.dic"));
auto &db = f.front();
cif::mm::structure s(db);
cif::crystal c(db);
auto struct_conn = db["struct_conn"];
for (const auto &[asym1, seqid1, authseqid1, atomid1, symm1,
asym2, seqid2, authseqid2, atomid2, symm2,
dist] : struct_conn.find<std::string, int, std::string, std::string, std::string,
std::string, int, std::string, std::string, std::string,
float>(
cif::key("ptnr1_symmetry") != "1_555" or cif::key("ptnr2_symmetry") != "1_555",
"ptnr1_label_asym_id", "ptnr1_label_seq_id", "ptnr1_auth_seq_id", "ptnr1_label_atom_id", "ptnr1_symmetry",
"ptnr2_label_asym_id", "ptnr2_label_seq_id", "ptnr2_auth_seq_id", "ptnr2_label_atom_id", "ptnr2_symmetry",
"pdbx_dist_value"))
{
auto &r1 = s.get_residue(asym1, seqid1, authseqid1);
auto &r2 = s.get_residue(asym2, seqid2, authseqid2);
auto a1 = r1.get_atom_by_atom_id(atomid1);
auto a2 = r2.get_atom_by_atom_id(atomid2);
auto p1 = a1.get_location();
auto p2 = a2.get_location();
cif::sym_op so1(symm1);
cif::sym_op so2(symm2);
auto sa1 = c.symmetry_copy(p1, so1);
auto sa2 = c.symmetry_copy(p2, so2);
CHECK_THAT(cif::distance(sa1, sa2), Catch::Matchers::WithinAbs(dist, 0.01f));
}
}
TEST_CASE("symm-P_32_2_1_a")
{
cif::cell c{ 80, 80, 120, 90, 90, 120 };
cif::spacegroup sg{ 154 };
cif::crystal crystal{ c, sg };
cif::point a{ 1, 90, 1 };
cif::point p1{ 2, 2, 2 };
auto d = distance(a, p1);
auto [d2, p, so] = crystal.closest_symmetry_copy(a, p1);
std::cout << "d: " << d2 << " p: " << p << " so: " << so.string() << '\n';
auto p2 = crystal.symmetry_copy(p1, so);
auto d3 = distance(p2, a);
CHECK_THAT(cif::distance(p2, p), Catch::Matchers::WithinAbs(0.f, 0.01f));
CHECK_THAT(d3, Catch::Matchers::WithinAbs(d2, 0.01f));
CHECK(d2 <= d);
}
TEST_CASE("symm-P_32_2_1")
{
cif::cell c{ 82.162, 82.162, 135.202, 90, 90, 120 };
cif::spacegroup sg{ 154 };
cif::crystal crystal{ c, sg };
cif::point a{ 1.73727,89.1813,11.1388 };
cif::point p1{ -8.98574, 50.3861, -11.6447 };
auto d = distance(a, p1);
auto [d2, p, so] = crystal.closest_symmetry_copy(a, p1);
std::cout << "d: " << d2 << " p: " << p << " so: " << so.string() << '\n';
auto p2 = crystal.symmetry_copy(p1, so);
auto d3 = distance(p2, a);
CHECK_THAT(cif::distance(p2, p), Catch::Matchers::WithinAbs(0.f, 0.01f));
CHECK_THAT(d3, Catch::Matchers::WithinAbs(d2, 0.01f));
CHECK(d2 <= d);
}
TEST_CASE("volume_3bwh_1")
{
cif::file f(gTestDir / "1juh.cif.gz");