Skip to content

L
libcifpp
Commit ea28ebdd by Maarten L. Hekkelman
Options

optimized caching of items in mmcif::Atom

parent 3ba46893
Showing with 267 additions and 208 deletions
include/cif++/Cif++.hpp
...@@ -36,6 +36,7 @@ ...@@ -36,6 +36,7 @@
#include <regex> #include <regex>
#include <set> #include <set>
#include <sstream> #include <sstream>
#include <shared_mutex>
#include "cif++/CifUtils.hpp" #include "cif++/CifUtils.hpp"
...@@ -142,13 +143,13 @@ class Item ...@@ -142,13 +143,13 @@ class Item
Item() {} Item() {}
template <typename T, std::enable_if_t<std::is_arithmetic_v<T>, int> = 0> template <typename T, std::enable_if_t<std::is_arithmetic_v<T>, int> = 0>
Item(const std::string &name, const T &value) Item(const std::string_view name, const T &value)
: mName(name) : mName(name)
, mValue(std::to_string(value)) , mValue(std::to_string(value))
{ {
} }
Item(const std::string &name, const std::string &value) Item(const std::string_view name, const std::string_view value)
: mName(name) : mName(name)
, mValue(value) , mValue(value)
{ {
...@@ -221,7 +222,7 @@ class Datablock ...@@ -221,7 +222,7 @@ class Datablock
using iterator = CategoryList::iterator; using iterator = CategoryList::iterator;
using const_iterator = CategoryList::const_iterator; using const_iterator = CategoryList::const_iterator;
Datablock(const std::string &name); Datablock(const std::string_view name);
~Datablock(); ~Datablock();
Datablock(const Datablock &) = delete; Datablock(const Datablock &) = delete;
...@@ -230,8 +231,6 @@ class Datablock ...@@ -230,8 +231,6 @@ class Datablock
std::string getName() const { return mName; } std::string getName() const { return mName; }
void setName(const std::string &n) { mName = n; } void setName(const std::string &n) { mName = n; }
std::string firstItem(const std::string &tag) const;
iterator begin() { return mCategories.begin(); } iterator begin() { return mCategories.begin(); }
iterator end() { return mCategories.end(); } iterator end() { return mCategories.end(); }
...@@ -256,7 +255,7 @@ class Datablock ...@@ -256,7 +255,7 @@ class Datablock
void write(std::ostream &os); void write(std::ostream &os);
// convenience function, add a line to the software category // convenience function, add a line to the software category
void add_software(const std::string &name, const std::string &classification, void add_software(const std::string_view name, const std::string &classification,
const std::string &versionNr, const std::string &versionDate); const std::string &versionNr, const std::string &versionDate);
friend bool operator==(const Datablock &lhs, const Datablock &rhs); friend bool operator==(const Datablock &lhs, const Datablock &rhs);
...@@ -264,7 +263,8 @@ class Datablock ...@@ -264,7 +263,8 @@ class Datablock
friend std::ostream& operator<<(std::ostream &os, const Datablock &data); friend std::ostream& operator<<(std::ostream &os, const Datablock &data);
private: private:
std::list<Category> mCategories; CategoryList mCategories; // LRU
mutable std::shared_mutex mLock;
std::string mName; std::string mName;
Validator *mValidator; Validator *mValidator;
Datablock *mNext; Datablock *mNext;
...@@ -1816,12 +1816,12 @@ class Category ...@@ -1816,12 +1816,12 @@ class Category
friend class Row; friend class Row;
friend class detail::ItemReference; friend class detail::ItemReference;
Category(Datablock &db, const std::string &name, Validator *Validator); Category(Datablock &db, const std::string_view name, Validator *Validator);
Category(const Category &) = delete; Category(const Category &) = delete;
Category &operator=(const Category &) = delete; Category &operator=(const Category &) = delete;
~Category(); ~Category();
const std::string name() const { return mName; } const std::string &name() const { return mName; }
using iterator = iterator_impl<Row>; using iterator = iterator_impl<Row>;
using const_iterator = iterator_impl<const Row>; using const_iterator = iterator_impl<const Row>;
......
include/cif++/Point.hpp
/*- /*-
* SPDX-License-Identifier: BSD-2-Clause * SPDX-License-Identifier: BSD-2-Clause
* *
* Copyright (c) 2020 NKI/AVL, Netherlands Cancer Institute * Copyright (c) 2020 NKI/AVL, Netherlands Cancer Institute
* *
* Redistribution and use in source and binary forms, with or without * Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are met: * modification, are permitted provided that the following conditions are met:
* *
* 1. Redistributions of source code must retain the above copyright notice, this * 1. Redistributions of source code must retain the above copyright notice, this
* list of conditions and the following disclaimer * list of conditions and the following disclaimer
* 2. Redistributions in binary form must reproduce the above copyright notice, * 2. Redistributions in binary form must reproduce the above copyright notice,
* this list of conditions and the following disclaimer in the documentation * this list of conditions and the following disclaimer in the documentation
* and/or other materials provided with the distribution. * and/or other materials provided with the distribution.
* *
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
* ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
* WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
...@@ -37,7 +37,7 @@ ...@@ -37,7 +37,7 @@
namespace mmcif namespace mmcif
{ {
typedef boost::math::quaternion<float> Quaternion; typedef boost::math::quaternion<float> Quaternion;
const double const double
kPI = 3.141592653589793238462643383279502884; kPI = 3.141592653589793238462643383279502884;
...@@ -51,26 +51,43 @@ const double ...@@ -51,26 +51,43 @@ const double
// float x, y, z; // float x, y, z;
// tie(x, y, z) = atom.loc(); // tie(x, y, z) = atom.loc();
template<typename F> template <typename F>
struct PointF struct PointF
{ {
typedef F FType; typedef F FType;
FType mX, mY, mZ; FType mX, mY, mZ;
PointF() : mX(0), mY(0), mZ(0) {}
PointF(FType x, FType y, FType z) : mX(x), mY(y), mZ(z) {}
template<typename PF> PointF()
PointF(const PointF<PF>& pt) : mX(0)
, mY(0)
, mZ(0)
{
}
PointF(FType x, FType y, FType z)
: mX(x)
, mY(y)
, mZ(z)
{
}
template <typename PF>
PointF(const PointF<PF> &pt)
: mX(static_cast<F>(pt.mX)) : mX(static_cast<F>(pt.mX))
, mY(static_cast<F>(pt.mY)) , mY(static_cast<F>(pt.mY))
, mZ(static_cast<F>(pt.mZ)) {} , mZ(static_cast<F>(pt.mZ))
{
}
#if HAVE_LIBCLIPPER #if HAVE_LIBCLIPPER
PointF(const clipper::Coord_orth& pt): mX(pt[0]), mY(pt[1]), mZ(pt[2]) {} PointF(const clipper::Coord_orth &pt)
: mX(pt[0])
, mY(pt[1])
, mZ(pt[2])
{
}
PointF& operator=(const clipper::Coord_orth& rhs) PointF &operator=(const clipper::Coord_orth &rhs)
{ {
mX = rhs[0]; mX = rhs[0];
mY = rhs[1]; mY = rhs[1];
...@@ -79,72 +96,72 @@ struct PointF ...@@ -79,72 +96,72 @@ struct PointF
} }
#endif #endif
template<typename PF> template <typename PF>
PointF& operator=(const PointF<PF>& rhs) PointF &operator=(const PointF<PF> &rhs)
{ {
mX = static_cast<F>(rhs.mX); mX = static_cast<F>(rhs.mX);
mY = static_cast<F>(rhs.mY); mY = static_cast<F>(rhs.mY);
mZ = static_cast<F>(rhs.mZ); mZ = static_cast<F>(rhs.mZ);
return *this; return *this;
} }
FType& getX() { return mX; } FType &getX() { return mX; }
FType getX() const { return mX; } FType getX() const { return mX; }
void setX(FType x) { mX = x; } void setX(FType x) { mX = x; }
FType& getY() { return mY; } FType &getY() { return mY; }
FType getY() const { return mY; } FType getY() const { return mY; }
void setY(FType y) { mY = y; } void setY(FType y) { mY = y; }
FType& getZ() { return mZ; } FType &getZ() { return mZ; }
FType getZ() const { return mZ; } FType getZ() const { return mZ; }
void setZ(FType z) { mZ = z; } void setZ(FType z) { mZ = z; }
PointF& operator+=(const PointF& rhs) PointF &operator+=(const PointF &rhs)
{ {
mX += rhs.mX; mX += rhs.mX;
mY += rhs.mY; mY += rhs.mY;
mZ += rhs.mZ; mZ += rhs.mZ;
return *this; return *this;
} }
PointF& operator+=(FType d) PointF &operator+=(FType d)
{ {
mX += d; mX += d;
mY += d; mY += d;
mZ += d; mZ += d;
return *this; return *this;
} }
PointF& operator-=(const PointF& rhs) PointF &operator-=(const PointF &rhs)
{ {
mX -= rhs.mX; mX -= rhs.mX;
mY -= rhs.mY; mY -= rhs.mY;
mZ -= rhs.mZ; mZ -= rhs.mZ;
return *this; return *this;
} }
PointF& operator-=(FType d) PointF &operator-=(FType d)
{ {
mX -= d; mX -= d;
mY -= d; mY -= d;
mZ -= d; mZ -= d;
return *this; return *this;
} }
PointF& operator*=(FType rhs) PointF &operator*=(FType rhs)
{ {
mX *= rhs; mX *= rhs;
mY *= rhs; mY *= rhs;
mZ *= rhs; mZ *= rhs;
return *this; return *this;
} }
PointF& operator/=(FType rhs) PointF &operator/=(FType rhs)
{ {
mX /= rhs; mX /= rhs;
mY /= rhs; mY /= rhs;
...@@ -162,18 +179,18 @@ struct PointF ...@@ -162,18 +179,18 @@ struct PointF
} }
return length; return length;
} }
void rotate(const boost::math::quaternion<FType>& q) void rotate(const boost::math::quaternion<FType> &q)
{ {
boost::math::quaternion<FType> p(0, mX, mY, mZ); boost::math::quaternion<FType> p(0, mX, mY, mZ);
p = q * p * boost::math::conj(q); p = q * p * boost::math::conj(q);
mX = p.R_component_2(); mX = p.R_component_2();
mY = p.R_component_3(); mY = p.R_component_3();
mZ = p.R_component_4(); mZ = p.R_component_4();
} }
#if HAVE_LIBCLIPPER #if HAVE_LIBCLIPPER
operator clipper::Coord_orth() const operator clipper::Coord_orth() const
{ {
...@@ -181,21 +198,21 @@ struct PointF ...@@ -181,21 +198,21 @@ struct PointF
} }
#endif #endif
operator std::tuple<const FType&, const FType&, const FType&>() const operator std::tuple<const FType &, const FType &, const FType &>() const
{ {
return std::make_tuple(std::ref(mX), std::ref(mY), std::ref(mZ)); return std::make_tuple(std::ref(mX), std::ref(mY), std::ref(mZ));
} }
operator std::tuple<FType&,FType&,FType&>() operator std::tuple<FType &, FType &, FType &>()
{ {
return std::make_tuple(std::ref(mX), std::ref(mY), std::ref(mZ)); return std::make_tuple(std::ref(mX), std::ref(mY), std::ref(mZ));
} }
bool operator==(const PointF& rhs) const bool operator==(const PointF &rhs) const
{ {
return mX == rhs.mX and mY == rhs.mY and mZ == rhs.mZ; return mX == rhs.mX and mY == rhs.mY and mZ == rhs.mZ;
} }
// consider point as a vector... perhaps I should rename Point? // consider point as a vector... perhaps I should rename Point?
FType lengthsq() const FType lengthsq() const
{ {
...@@ -211,45 +228,45 @@ struct PointF ...@@ -211,45 +228,45 @@ struct PointF
typedef PointF<float> Point; typedef PointF<float> Point;
typedef PointF<double> DPoint; typedef PointF<double> DPoint;
template<typename F> template <typename F>
inline std::ostream& operator<<(std::ostream& os, const PointF<F>& pt) inline std::ostream &operator<<(std::ostream &os, const PointF<F> &pt)
{ {
os << '(' << pt.mX << ',' << pt.mY << ',' << pt.mZ << ')'; os << '(' << pt.mX << ',' << pt.mY << ',' << pt.mZ << ')';
return os; return os;
} }
template<typename F> template <typename F>
inline PointF<F> operator+(const PointF<F>& lhs, const PointF<F>& rhs) inline PointF<F> operator+(const PointF<F> &lhs, const PointF<F> &rhs)
{ {
return PointF<F>(lhs.mX + rhs.mX, lhs.mY + rhs.mY, lhs.mZ + rhs.mZ); return PointF<F>(lhs.mX + rhs.mX, lhs.mY + rhs.mY, lhs.mZ + rhs.mZ);
} }
template<typename F> template <typename F>
inline PointF<F> operator-(const PointF<F>& lhs, const PointF<F>& rhs) inline PointF<F> operator-(const PointF<F> &lhs, const PointF<F> &rhs)
{ {
return PointF<F>(lhs.mX - rhs.mX, lhs.mY - rhs.mY, lhs.mZ - rhs.mZ); return PointF<F>(lhs.mX - rhs.mX, lhs.mY - rhs.mY, lhs.mZ - rhs.mZ);
} }
template<typename F> template <typename F>
inline PointF<F> operator-(const PointF<F>& pt) inline PointF<F> operator-(const PointF<F> &pt)
{ {
return PointF<F>(-pt.mX, -pt.mY, -pt.mZ); return PointF<F>(-pt.mX, -pt.mY, -pt.mZ);
} }
template<typename F> template <typename F>
inline PointF<F> operator*(const PointF<F>& pt, F f) inline PointF<F> operator*(const PointF<F> &pt, F f)
{ {
return PointF<F>(pt.mX * f, pt.mY * f, pt.mZ * f); return PointF<F>(pt.mX * f, pt.mY * f, pt.mZ * f);
} }
template<typename F> template <typename F>
inline PointF<F> operator*(F f, const PointF<F>& pt) inline PointF<F> operator*(F f, const PointF<F> &pt)
{ {
return PointF<F>(pt.mX * f, pt.mY * f, pt.mZ * f); return PointF<F>(pt.mX * f, pt.mY * f, pt.mZ * f);
} }
template<typename F> template <typename F>
inline PointF<F> operator/(const PointF<F>& pt, F f) inline PointF<F> operator/(const PointF<F> &pt, F f)
{ {
return PointF<F>(pt.mX / f, pt.mY / f, pt.mZ / f); return PointF<F>(pt.mX / f, pt.mY / f, pt.mZ / f);
} }
...@@ -257,17 +274,16 @@ inline PointF<F> operator/(const PointF<F>& pt, F f) ...@@ -257,17 +274,16 @@ inline PointF<F> operator/(const PointF<F>& pt, F f)
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// several standard 3d operations // several standard 3d operations
template<typename F> template <typename F>
inline double DistanceSquared(const PointF<F>& a, const PointF<F>& b) inline double DistanceSquared(const PointF<F> &a, const PointF<F> &b)
{ {
return return (a.mX - b.mX) * (a.mX - b.mX) +
(a.mX - b.mX) * (a.mX - b.mX) + (a.mY - b.mY) * (a.mY - b.mY) +
(a.mY - b.mY) * (a.mY - b.mY) + (a.mZ - b.mZ) * (a.mZ - b.mZ);
(a.mZ - b.mZ) * (a.mZ - b.mZ);
} }
template<typename F> template <typename F>
inline double Distance(const PointF<F>& a, const PointF<F>& b) inline double Distance(const PointF<F> &a, const PointF<F> &b)
{ {
return sqrt( return sqrt(
(a.mX - b.mX) * (a.mX - b.mX) + (a.mX - b.mX) * (a.mX - b.mX) +
...@@ -275,44 +291,44 @@ inline double Distance(const PointF<F>& a, const PointF<F>& b) ...@@ -275,44 +291,44 @@ inline double Distance(const PointF<F>& a, const PointF<F>& b)
(a.mZ - b.mZ) * (a.mZ - b.mZ)); (a.mZ - b.mZ) * (a.mZ - b.mZ));
} }
template<typename F> template <typename F>
inline F DotProduct(const PointF<F>& a, const PointF<F>& b) inline F DotProduct(const PointF<F> &a, const PointF<F> &b)
{ {
return a.mX * b.mX + a.mY * b.mY + a.mZ * b.mZ; return a.mX * b.mX + a.mY * b.mY + a.mZ * b.mZ;
} }
template<typename F> template <typename F>
inline PointF<F> CrossProduct(const PointF<F>& a, const PointF<F>& b) inline PointF<F> CrossProduct(const PointF<F> &a, const PointF<F> &b)
{ {
return PointF<F>(a.mY * b.mZ - b.mY * a.mZ, return PointF<F>(a.mY * b.mZ - b.mY * a.mZ,
a.mZ * b.mX - b.mZ * a.mX, a.mZ * b.mX - b.mZ * a.mX,
a.mX * b.mY - b.mX * a.mY); a.mX * b.mY - b.mX * a.mY);
} }
template<typename F> template <typename F>
double Angle(const PointF<F>& p1, const PointF<F>& p2, const PointF<F>& p3) double Angle(const PointF<F> &p1, const PointF<F> &p2, const PointF<F> &p3)
{ {
PointF<F> v1 = p1 - p2; PointF<F> v1 = p1 - p2;
PointF<F> v2 = p3 - p2; PointF<F> v2 = p3 - p2;
return std::acos(DotProduct(v1, v2) / (v1.length() * v2.length())) * 180 / kPI; return std::acos(DotProduct(v1, v2) / (v1.length() * v2.length())) * 180 / kPI;
} }
template<typename F> template <typename F>
double DihedralAngle(const PointF<F>& p1, const PointF<F>& p2, const PointF<F>& p3, const PointF<F>& p4) double DihedralAngle(const PointF<F> &p1, const PointF<F> &p2, const PointF<F> &p3, const PointF<F> &p4)
{ {
PointF<F> v12 = p1 - p2; // vector from p2 to p1 PointF<F> v12 = p1 - p2; // vector from p2 to p1
PointF<F> v43 = p4 - p3; // vector from p3 to p4 PointF<F> v43 = p4 - p3; // vector from p3 to p4
PointF<F> z = p2 - p3; // vector from p3 to p2 PointF<F> z = p2 - p3; // vector from p3 to p2
PointF<F> p = CrossProduct(z, v12); PointF<F> p = CrossProduct(z, v12);
PointF<F> x = CrossProduct(z, v43); PointF<F> x = CrossProduct(z, v43);
PointF<F> y = CrossProduct(z, x); PointF<F> y = CrossProduct(z, x);
double u = DotProduct(x, x); double u = DotProduct(x, x);
double v = DotProduct(y, y); double v = DotProduct(y, y);
double result = 360; double result = 360;
if (u > 0 and v > 0) if (u > 0 and v > 0)
{ {
...@@ -321,33 +337,33 @@ double DihedralAngle(const PointF<F>& p1, const PointF<F>& p2, const PointF<F>& ...@@ -321,33 +337,33 @@ double DihedralAngle(const PointF<F>& p1, const PointF<F>& p2, const PointF<F>&
if (u != 0 or v != 0) if (u != 0 or v != 0)
result = atan2(v, u) * 180 / kPI; result = atan2(v, u) * 180 / kPI;
} }
return result; return result;
} }
template<typename F> template <typename F>
double CosinusAngle(const PointF<F>& p1, const PointF<F>& p2, const PointF<F>& p3, const PointF<F>& p4) double CosinusAngle(const PointF<F> &p1, const PointF<F> &p2, const PointF<F> &p3, const PointF<F> &p4)
{ {
PointF<F> v12 = p1 - p2; PointF<F> v12 = p1 - p2;
PointF<F> v34 = p3 - p4; PointF<F> v34 = p3 - p4;
double result = 0; double result = 0;
double x = DotProduct(v12, v12) * DotProduct(v34, v34); double x = DotProduct(v12, v12) * DotProduct(v34, v34);
if (x > 0) if (x > 0)
result = DotProduct(v12, v34) / sqrt(x); result = DotProduct(v12, v34) / sqrt(x);
return result; return result;
} }
template<typename F> template <typename F>
auto DistancePointToLine(const PointF<F> &l1, const PointF<F> &l2, const PointF<F> &p) auto DistancePointToLine(const PointF<F> &l1, const PointF<F> &l2, const PointF<F> &p)
{ {
auto line = l2 - l1; auto line = l2 - l1;
auto p_to_l1 = p - l1; auto p_to_l1 = p - l1;
auto p_to_l2 = p - l2; auto p_to_l2 = p - l2;
auto cross = CrossProduct(p_to_l1, p_to_l2); auto cross = CrossProduct(p_to_l1, p_to_l2);
return cross.length() / line.length(); return cross.length() / line.length();
} }
// -------------------------------------------------------------------- // --------------------------------------------------------------------
...@@ -355,7 +371,7 @@ auto DistancePointToLine(const PointF<F> &l1, const PointF<F> &l2, const PointF< ...@@ -355,7 +371,7 @@ auto DistancePointToLine(const PointF<F> &l1, const PointF<F> &l2, const PointF<
// a random direction with a distance randomly chosen from a normal // a random direction with a distance randomly chosen from a normal
// distribution with a stddev of offset. // distribution with a stddev of offset.
template<typename F> template <typename F>
PointF<F> Nudge(PointF<F> p, F offset); PointF<F> Nudge(PointF<F> p, F offset);
// -------------------------------------------------------------------- // --------------------------------------------------------------------
...@@ -363,66 +379,68 @@ PointF<F> Nudge(PointF<F> p, F offset); ...@@ -363,66 +379,68 @@ PointF<F> Nudge(PointF<F> p, F offset);
Quaternion Normalize(Quaternion q); Quaternion Normalize(Quaternion q);
std::tuple<double,Point> QuaternionToAngleAxis(Quaternion q); std::tuple<double, Point> QuaternionToAngleAxis(Quaternion q);
Point Centroid(std::vector<Point>& Points); Point Centroid(const std::vector<Point> &Points);
Point CenterPoints(std::vector<Point>& Points); Point CenterPoints(std::vector<Point> &Points);
Quaternion AlignPoints(const std::vector<Point>& a, const std::vector<Point>& b); Quaternion AlignPoints(const std::vector<Point> &a, const std::vector<Point> &b);
double RMSd(const std::vector<Point>& a, const std::vector<Point>& b); double RMSd(const std::vector<Point> &a, const std::vector<Point> &b);
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// Helper class to generate evenly divided Points on a sphere // Helper class to generate evenly divided Points on a sphere
// we use a fibonacci sphere to calculate even distribution of the dots // we use a fibonacci sphere to calculate even distribution of the dots
template<int N> template <int N>
class SphericalDots class SphericalDots
{ {
public: public:
enum { P = 2 * N + 1 }; enum
typedef typename std::array<Point,P> array_type; {
typedef typename array_type::const_iterator iterator; P = 2 * N + 1
};
typedef typename std::array<Point, P> array_type;
typedef typename array_type::const_iterator iterator;
static SphericalDots& instance() static SphericalDots &instance()
{ {
static SphericalDots sInstance; static SphericalDots sInstance;
return sInstance; return sInstance;
} }
size_t size() const { return mPoints.size(); }
const Point operator[](uint32_t inIx) const { return mPoints[inIx]; }
iterator begin() const { return mPoints.begin(); }
iterator end() const { return mPoints.end(); }
double weight() const { return mWeight; } size_t size() const { return mPoints.size(); }
const Point operator[](uint32_t inIx) const { return mPoints[inIx]; }
iterator begin() const { return mPoints.begin(); }
iterator end() const { return mPoints.end(); }
double weight() const { return mWeight; }
SphericalDots() SphericalDots()
{ {
const double const double
kGoldenRatio = (1 + std::sqrt(5.0)) / 2; kGoldenRatio = (1 + std::sqrt(5.0)) / 2;
mWeight = (4 * kPI) / P; mWeight = (4 * kPI) / P;
auto p = mPoints.begin(); auto p = mPoints.begin();
for (int32_t i = -N; i <= N; ++i) for (int32_t i = -N; i <= N; ++i)
{ {
double lat = std::asin((2.0 * i) / P); double lat = std::asin((2.0 * i) / P);
double lon = std::fmod(i, kGoldenRatio) * 2 * kPI / kGoldenRatio; double lon = std::fmod(i, kGoldenRatio) * 2 * kPI / kGoldenRatio;
p->mX = sin(lon) * cos(lat); p->mX = sin(lon) * cos(lat);
p->mY = cos(lon) * cos(lat); p->mY = cos(lon) * cos(lat);
p->mZ = sin(lat); p->mZ = sin(lat);
++p; ++p;
} }
} }
private: private:
array_type mPoints;
array_type mPoints; double mWeight;
double mWeight;
}; };
typedef SphericalDots<50> SphericalDots_50; typedef SphericalDots<50> SphericalDots_50;
} } // namespace mmcif
include/cif++/Structure.hpp
...@@ -109,12 +109,12 @@ class Atom ...@@ -109,12 +109,12 @@ class Atom
float occupancy() const; float occupancy() const;
template <typename T> template <typename T>
T property(const std::string &name) const; T property(const std::string_view name) const;
void property(const std::string &name, const std::string &value); void property(const std::string_view name, const std::string &value);
template <typename T, std::enable_if_t<std::is_arithmetic_v<T>, int> = 0> template <typename T, std::enable_if_t<std::is_arithmetic_v<T>, int> = 0>
void property(const std::string &name, const T &value) void property(const std::string_view name, const T &value)
{ {
property(name, std::to_string(value)); property(name, std::to_string(value));
} }
...@@ -404,7 +404,7 @@ class File : public std::enable_shared_from_this<File> ...@@ -404,7 +404,7 @@ class File : public std::enable_shared_from_this<File>
File(const File &) = delete; File(const File &) = delete;
File &operator=(const File &) = delete; File &operator=(const File &) = delete;
cif::Datablock& createDatablock(const std::string &name); cif::Datablock& createDatablock(const std::string_view name);
void load(const std::filesystem::path &path); void load(const std::filesystem::path &path);
void save(const std::filesystem::path &path); void save(const std::filesystem::path &path);
......
src/Cif++.cpp
...@@ -33,6 +33,7 @@ ...@@ -33,6 +33,7 @@
#include <stack> #include <stack>
#include <tuple> #include <tuple>
#include <unordered_map> #include <unordered_map>
#include <shared_mutex>
#include <filesystem> #include <filesystem>
...@@ -351,7 +352,7 @@ namespace detail ...@@ -351,7 +352,7 @@ namespace detail
// -------------------------------------------------------------------- // --------------------------------------------------------------------
// Datablock implementation // Datablock implementation
Datablock::Datablock(const std::string &name) Datablock::Datablock(const std::string_view name)
: mName(name) : mName(name)
, mValidator(nullptr) , mValidator(nullptr)
, mNext(nullptr) , mNext(nullptr)
...@@ -363,43 +364,57 @@ Datablock::~Datablock() ...@@ -363,43 +364,57 @@ Datablock::~Datablock()
delete mNext; delete mNext;
} }
std::string Datablock::firstItem(const std::string &tag) const auto Datablock::emplace(std::string_view name) -> std::tuple<iterator, bool>
{ {
std::string result; // LRU code
std::string catName, itemName; std::shared_lock lock(mLock);
std::tie(catName, itemName) = splitTagName(tag);
for (auto &cat : mCategories) bool isNew = true;
auto i = begin();
while (i != end())
{ {
if (iequals(cat.name(), catName)) if (iequals(name, i->name()))
{ {
for (auto row : cat) isNew = false;
if (i != begin())
{ {
result = row[itemName].as<std::string>(); auto n = std::next(i);
break; mCategories.splice(begin(), mCategories, i, n);
} }
break; break;
} }
++i;
} }
return result; if (isNew)
} mCategories.emplace(begin(), *this, std::string(name), mValidator);
auto Datablock::emplace(std::string_view name) -> std::tuple<iterator, bool> return std::make_tuple(begin(), isNew);
{
bool isNew = false;
iterator i = find_if(begin(), end(), [name](const Category &cat) -> bool
{ return iequals(cat.name(), name); });
if (i == end())
{
isNew = true;
i = mCategories.emplace(end(), *this, std::string(name), mValidator);
}
return std::make_tuple(i, isNew); // bool isNew = false;
// iterator i = find_if(begin(), end(), [name](const Category &cat) -> bool
// { return iequals(cat.name(), name); });
// if (i == end())
// {
// isNew = true;
// i = mCategories.emplace(begin(), *this, std::string(name), mValidator);
// }
// else if (i != begin())
// {
// auto n = std::next(i);
// mCategories.splice(begin(), mCategories, i, n);
// i = mCategories.begin();
// }
// return std::make_tuple(i, isNew);
} }
Category &Datablock::operator[](std::string_view name) Category &Datablock::operator[](std::string_view name)
...@@ -411,14 +426,13 @@ Category &Datablock::operator[](std::string_view name) ...@@ -411,14 +426,13 @@ Category &Datablock::operator[](std::string_view name)
Category *Datablock::get(std::string_view name) Category *Datablock::get(std::string_view name)
{ {
auto i = find_if(begin(), end(), [name](const Category &cat) -> bool return &operator[](name);
{ return iequals(cat.name(), name); });
return i == end() ? nullptr : &*i;
} }
const Category *Datablock::get(std::string_view name) const const Category *Datablock::get(std::string_view name) const
{ {
std::shared_lock lock(mLock);
auto i = find_if(begin(), end(), [name](const Category &cat) -> bool auto i = find_if(begin(), end(), [name](const Category &cat) -> bool
{ return iequals(cat.name(), name); }); { return iequals(cat.name(), name); });
...@@ -427,6 +441,8 @@ const Category *Datablock::get(std::string_view name) const ...@@ -427,6 +441,8 @@ const Category *Datablock::get(std::string_view name) const
bool Datablock::isValid() bool Datablock::isValid()
{ {
std::shared_lock lock(mLock);
if (mValidator == nullptr) if (mValidator == nullptr)
throw std::runtime_error("Validator not specified"); throw std::runtime_error("Validator not specified");
...@@ -438,20 +454,26 @@ bool Datablock::isValid() ...@@ -438,20 +454,26 @@ bool Datablock::isValid()
void Datablock::validateLinks() const void Datablock::validateLinks() const
{ {
std::shared_lock lock(mLock);
for (auto &cat : *this) for (auto &cat : *this)
cat.validateLinks(); cat.validateLinks();
} }
void Datablock::setValidator(Validator *v) void Datablock::setValidator(Validator *v)
{ {
std::shared_lock lock(mLock);
mValidator = v; mValidator = v;
for (auto &cat : *this) for (auto &cat : *this)
cat.setValidator(v); cat.setValidator(v);
} }
void Datablock::add_software(const std::string &name, const std::string &classification, const std::string &versionNr, const std::string &versionDate) void Datablock::add_software(const std::string_view name, const std::string &classification, const std::string &versionNr, const std::string &versionDate)
{ {
std::shared_lock lock(mLock);
Category &cat = operator[]("software"); Category &cat = operator[]("software");
auto ordNr = cat.size() + 1; auto ordNr = cat.size() + 1;
// TODO: should we check this ordinal number??? // TODO: should we check this ordinal number???
...@@ -465,12 +487,16 @@ void Datablock::add_software(const std::string &name, const std::string &classif ...@@ -465,12 +487,16 @@ void Datablock::add_software(const std::string &name, const std::string &classif
void Datablock::getTagOrder(std::vector<std::string> &tags) const void Datablock::getTagOrder(std::vector<std::string> &tags) const
{ {
std::shared_lock lock(mLock);
for (auto &cat : *this) for (auto &cat : *this)
cat.getTagOrder(tags); cat.getTagOrder(tags);
} }
void Datablock::write(std::ostream &os) void Datablock::write(std::ostream &os)
{ {
std::shared_lock lock(mLock);
os << "data_" << mName << std::endl os << "data_" << mName << std::endl
<< "# " << std::endl; << "# " << std::endl;
...@@ -505,6 +531,8 @@ void Datablock::write(std::ostream &os) ...@@ -505,6 +531,8 @@ void Datablock::write(std::ostream &os)
void Datablock::write(std::ostream &os, const std::vector<std::string> &order) void Datablock::write(std::ostream &os, const std::vector<std::string> &order)
{ {
std::shared_lock lock(mLock);
os << "data_" << mName << std::endl os << "data_" << mName << std::endl
<< "# " << std::endl; << "# " << std::endl;
...@@ -580,6 +608,9 @@ void Datablock::write(std::ostream &os, const std::vector<std::string> &order) ...@@ -580,6 +608,9 @@ void Datablock::write(std::ostream &os, const std::vector<std::string> &order)
bool operator==(const cif::Datablock &dbA, const cif::Datablock &dbB) bool operator==(const cif::Datablock &dbA, const cif::Datablock &dbB)
{ {
std::shared_lock lockA(dbA.mLock);
std::shared_lock lockB(dbB.mLock);
std::vector<std::string> catA, catB; std::vector<std::string> catA, catB;
for (auto &cat : dbA) for (auto &cat : dbA)
...@@ -1311,7 +1342,7 @@ RowSet &RowSet::orderBy(std::initializer_list<std::string> items) ...@@ -1311,7 +1342,7 @@ RowSet &RowSet::orderBy(std::initializer_list<std::string> items)
// -------------------------------------------------------------------- // --------------------------------------------------------------------
Category::Category(Datablock &db, const std::string &name, Validator *Validator) Category::Category(Datablock &db, const std::string_view name, Validator *Validator)
: mDb(db) : mDb(db)
, mName(name) , mName(name)
, mValidator(Validator) , mValidator(Validator)
......
src/CifParser.cpp
...@@ -976,7 +976,8 @@ void DictParser::parseSaveFrame() ...@@ -976,7 +976,8 @@ void DictParser::parseSaveFrame()
if (isCategorySaveFrame) if (isCategorySaveFrame)
{ {
std::string category = dict.firstItem("_category.id"); std::string category;
cif::tie(category) = dict["category"].front().get("id");
std::vector<std::string> keys; std::vector<std::string> keys;
for (auto k: dict["category_key"]) for (auto k: dict["category_key"])
...@@ -991,7 +992,8 @@ void DictParser::parseSaveFrame() ...@@ -991,7 +992,8 @@ void DictParser::parseSaveFrame()
else else
{ {
// if the type code is missing, this must be a pointer, just skip it // if the type code is missing, this must be a pointer, just skip it
std::string typeCode = dict.firstItem("_item_type.code"); std::string typeCode;
cif::tie(typeCode) = dict["item_type"].front().get("code");
const ValidateType* tv = nullptr; const ValidateType* tv = nullptr;
if (not (typeCode.empty() or typeCode == "?")) if (not (typeCode.empty() or typeCode == "?"))
...@@ -1001,7 +1003,8 @@ void DictParser::parseSaveFrame() ...@@ -1001,7 +1003,8 @@ void DictParser::parseSaveFrame()
for (auto e: dict["item_enumeration"]) for (auto e: dict["item_enumeration"])
ess.insert(e["value"].as<std::string>()); ess.insert(e["value"].as<std::string>());
std::string defaultValue = dict.firstItem("_item_default.value"); std::string defaultValue;
cif::tie(defaultValue) = dict["item_default"].front().get("value");
bool defaultIsNull = false; bool defaultIsNull = false;
if (defaultValue.empty()) if (defaultValue.empty())
{ {
......
src/Point.cpp
...@@ -102,14 +102,14 @@ Point CenterPoints(std::vector<Point>& Points) ...@@ -102,14 +102,14 @@ Point CenterPoints(std::vector<Point>& Points)
return t; return t;
} }
Point Centroid(std::vector<Point>& Points) Point Centroid(const std::vector<Point>& pts)
{ {
Point result; Point result;
for (Point& pt : Points) for (auto &pt : pts)
result += pt; result += pt;
result /= static_cast<float>(Points.size()); result /= static_cast<float>(pts.size());
return result; return result;
} }
......
src/Structure.cpp
...@@ -216,9 +216,9 @@ struct AtomImpl ...@@ -216,9 +216,9 @@ struct AtomImpl
, mLocation(i.mLocation) , mLocation(i.mLocation)
, mRefcount(1) , mRefcount(1)
, mRow(i.mRow) , mRow(i.mRow)
, mCachedRefs(i.mCachedRefs)
, mCompound(i.mCompound) , mCompound(i.mCompound)
, mRadius(i.mRadius) , mRadius(i.mRadius)
, mCachedProperties(i.mCachedProperties)
, mSymmetryCopy(i.mSymmetryCopy) , mSymmetryCopy(i.mSymmetryCopy)
, mClone(true) , mClone(true)
// , mRTop(i.mRTop), mD(i.mD) // , mRTop(i.mRTop), mD(i.mD)
...@@ -270,9 +270,9 @@ struct AtomImpl ...@@ -270,9 +270,9 @@ struct AtomImpl
, mLocation(loc) , mLocation(loc)
, mRefcount(1) , mRefcount(1)
, mRow(impl.mRow) , mRow(impl.mRow)
, mCachedRefs(impl.mCachedRefs)
, mCompound(impl.mCompound) , mCompound(impl.mCompound)
, mRadius(impl.mRadius) , mRadius(impl.mRadius)
, mCachedProperties(impl.mCachedProperties)
, mSymmetryCopy(true) , mSymmetryCopy(true)
, mSymmetryOperator(sym_op) , mSymmetryOperator(sym_op)
{ {
...@@ -338,9 +338,9 @@ struct AtomImpl ...@@ -338,9 +338,9 @@ struct AtomImpl
if (not mClone) if (not mClone)
{ {
mRow["Cartn_x"] = p.getX(); property("Cartn_x", std::to_string(p.getX()));
mRow["Cartn_y"] = p.getY(); property("Cartn_y", std::to_string(p.getY()));
mRow["Cartn_z"] = p.getZ(); property("Cartn_z", std::to_string(p.getZ()));
} }
// boost::format kPosFmt("%.3f"); // boost::format kPosFmt("%.3f");
...@@ -382,26 +382,31 @@ struct AtomImpl ...@@ -382,26 +382,31 @@ struct AtomImpl
return mRadius; return mRadius;
} }
const std::string &property(const std::string &name) const const std::string property(const std::string_view name) const
{ {
static std::string kEmptyString; for (auto &&[tag, ref] : mCachedRefs)
auto i = mCachedProperties.find(name);
if (i == mCachedProperties.end())
{ {
auto v = mRow[name]; if (tag == name)
if (v.empty()) return ref.as<std::string>();
return kEmptyString;
return mCachedProperties[name] = v.as<std::string>();
} }
else
return i->second; mCachedRefs.emplace_back(name, mRow[name]);
return std::get<1>(mCachedRefs.back()).as<std::string>();
} }
void property(const std::string &name, const std::string &value) void property(const std::string_view name, const std::string &value)
{ {
mRow[name] = value; for (auto &&[tag, ref] : mCachedRefs)
{
if (tag != name)
continue;
ref = value;
return;
}
mCachedRefs.emplace_back(name, mRow[name]);
std::get<1>(mCachedRefs.back()) = value;
} }
int compare(const AtomImpl &b) const int compare(const AtomImpl &b) const
...@@ -432,9 +437,11 @@ struct AtomImpl ...@@ -432,9 +437,11 @@ struct AtomImpl
Point mLocation; Point mLocation;
int mRefcount; int mRefcount;
cif::Row mRow; cif::Row mRow;
mutable std::vector<std::tuple<std::string,cif::detail::ItemReference>> mCachedRefs;
mutable const Compound *mCompound = nullptr; mutable const Compound *mCompound = nullptr;
float mRadius = std::nanf("4"); float mRadius = std::nanf("4");
mutable std::map<std::string, std::string> mCachedProperties;
bool mSymmetryCopy = false; bool mSymmetryCopy = false;
bool mClone = false; bool mClone = false;
...@@ -533,25 +540,25 @@ const cif::Row Atom::getRowAniso() const ...@@ -533,25 +540,25 @@ const cif::Row Atom::getRowAniso() const
} }
template <> template <>
std::string Atom::property<std::string>(const std::string &name) const std::string Atom::property<std::string>(const std::string_view name) const
{ {
return impl()->property(name); return impl()->property(name);
} }
template <> template <>
int Atom::property<int>(const std::string &name) const int Atom::property<int>(const std::string_view name) const
{ {
auto v = impl()->property(name); auto v = impl()->property(name);
return v.empty() ? 0 : stoi(v); return v.empty() ? 0 : stoi(v);
} }
template <> template <>
float Atom::property<float>(const std::string &name) const float Atom::property<float>(const std::string_view name) const
{ {
return stof(impl()->property(name)); return stof(impl()->property(name));
} }
void Atom::property(const std::string &name, const std::string &value) void Atom::property(const std::string_view name, const std::string &value)
{ {
impl()->property(name, value); impl()->property(name, value);
} }
...@@ -1736,7 +1743,7 @@ File::~File() ...@@ -1736,7 +1743,7 @@ File::~File()
delete mImpl; delete mImpl;
} }
cif::Datablock& File::createDatablock(const std::string &name) cif::Datablock& File::createDatablock(const std::string_view name)
{ {
auto db = new cif::Datablock(name); auto db = new cif::Datablock(name);
...@@ -1807,9 +1814,9 @@ Structure::Structure(File &f, size_t modelNr, StructureOpenOptions options) ...@@ -1807,9 +1814,9 @@ Structure::Structure(File &f, size_t modelNr, StructureOpenOptions options)
} }
if (mAtoms.empty()) if (mAtoms.empty())
throw std::runtime_error("No atoms loaded, refuse to continue"); std::cerr << "Warning: no atoms loaded" << std::endl;
else
loadData(); loadData();
} }
void Structure::loadAtomsForModel(StructureOpenOptions options) void Structure::loadAtomsForModel(StructureOpenOptions options)
......
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