Fix calculating an elliptical arc.

develop
Roman Telezhynskyi 2022-08-27 16:46:25 +03:00
parent ab75b783b6
commit 1e344d6df0
12 changed files with 1509 additions and 86 deletions

View File

@ -20,6 +20,7 @@
- Fix segmenting a simple curve.
- Fix export tiled pdf on Mac OS.
- Improve visualization for tools.
- Fix calculating an elliptical arc.
# Valentina 0.7.51 April 18, 2022
- Z value change for a layout piece.

View File

@ -77,7 +77,6 @@ public:
* @brief Clear Clears the carrousel (removes everything)
*/
void Clear();
auto Layout() const -> VPLayoutWeakPtr;
public slots:

View File

@ -37,9 +37,171 @@
#include "../ifc/ifcdef.h"
#include "../ifc/exception/vexception.h"
#include "../vmisc/vabstractapplication.h"
#include "../vmisc/fpm/fixed.hpp"
#include "../vmisc/fpm/math.hpp"
#include "../vmisc/compatibility.h"
#include "vabstractcurve.h"
#include "vellipticalarc_p.h"
#include "../vmisc/vmath.h"
namespace
{
constexpr qreal tolerance = accuracyPointOnLine/32;
//---------------------------------------------------------------------------------------------------------------------
auto VLen(fpm::fixed_16_16 x, fpm::fixed_16_16 y) -> fpm::fixed_16_16
{
x = fpm::abs(x);
y = fpm::abs(y);
if (x > y)
{
return x + qMax(y/8, y/2 - x/8);
}
return y + qMax(x/8, x/2 - y/8);
}
//---------------------------------------------------------------------------------------------------------------------
auto AuxRadius(fpm::fixed_16_16 xP, fpm::fixed_16_16 yP, fpm::fixed_16_16 xQ, fpm::fixed_16_16 yQ) -> fpm::fixed_16_16
{
fpm::fixed_16_16 dP = VLen(xP, yP);
fpm::fixed_16_16 dQ = VLen(xQ, yQ);
fpm::fixed_16_16 dJ = VLen(xP + xQ, yP + yQ);
fpm::fixed_16_16 dK = VLen(xP - xQ, yP - yQ);
fpm::fixed_16_16 r1 = qMax(dP, dQ);
fpm::fixed_16_16 r2 = qMax(dJ, dK);
return qMax(r1 + r1/16, r2 - r2/4);
}
//---------------------------------------------------------------------------------------------------------------------
auto AngularInc(fpm::fixed_16_16 xP, fpm::fixed_16_16 yP, fpm::fixed_16_16 xQ, fpm::fixed_16_16 yQ,
fpm::fixed_16_16 flatness) -> int
{
fpm::fixed_16_16 r = AuxRadius(xP, yP, xQ, yQ);
fpm::fixed_16_16 err2{r >> 3};
// 2nd-order term
fpm::fixed_16_16 err4{r >> 7};
// 4th-order term
const int kmax = qRound(0.5 * std::log2(maxSceneSize / (8. * tolerance)));
for (int k = 0; k < kmax; ++k)
{
if (flatness >= err2 + err4)
{
return k;
}
err2 >>= 2;
err4 >>= 4;
}
return kmax;
}
//---------------------------------------------------------------------------------------------------------------------
inline void CircleGen(fpm::fixed_16_16& u, fpm::fixed_16_16& v, uint k)
{
u -= v >> k;
v += u >> k;
}
//---------------------------------------------------------------------------------------------------------------------
auto InitialValue(fpm::fixed_16_16 u0, fpm::fixed_16_16 v0, uint k) -> fpm::fixed_16_16
{
uint shift = 2*k + 3;
fpm::fixed_16_16 w {u0 >> shift};
fpm::fixed_16_16 U0 = u0 - w + (v0 >> (k + 1));
w >>= (shift + 1);
U0 -= w;
w >>= shift;
U0 -= w;
return U0;
}
//---------------------------------------------------------------------------------------------------------------------
auto EllipseCore(fpm::fixed_16_16 xC, fpm::fixed_16_16 yC, fpm::fixed_16_16 xP, fpm::fixed_16_16 yP,
fpm::fixed_16_16 xQ, fpm::fixed_16_16 yQ, fpm::fixed_16_16 sweep,
fpm::fixed_16_16 flatness) -> QVector<QPointF>
{
uint k = qMin(static_cast<uint>(AngularInc(xP, yP, xQ, yQ, flatness)), 16U);
const uint count = static_cast<std::uint32_t>(sweep.raw_value()) >> (16 - k);
QVector<QPointF> arc;
arc.reserve(static_cast<int>(count) + 1);
// Arc start point
arc.append({static_cast<qreal>(xP + xC), static_cast<qreal>(yP + yC)});
xQ = InitialValue(xQ, xP, k);
yQ = InitialValue(yQ, yP, k);
for (uint i = 0; i < count; ++i)
{
CircleGen(xQ, xP, k);
CircleGen(yQ, yP, k);
arc.append({static_cast<qreal>(xP + xC), static_cast<qreal>(yP + yC)});
}
return arc;
}
//---------------------------------------------------------------------------------------------------------------------
auto EllipticArcPoints(QPointF c, qreal radius1, qreal radius2, qreal astart, qreal asweep,
qreal approximationScale) -> QVector<QPointF>
{
fpm::fixed_16_16 xC{c.x()};
fpm::fixed_16_16 yC{c.y()};
fpm::fixed_16_16 xP{c.x() + radius1};
fpm::fixed_16_16 yP{c.y()};
fpm::fixed_16_16 xQ{c.x()};
fpm::fixed_16_16 yQ{c.y() - radius2};
xP -= xC;
yP -= yC;
xQ -= xC;
yQ -= yC;
if (not qFuzzyIsNull(astart))
{
// Set new conjugate diameter end points P and Q
fpm::fixed_16_16 cosa {cos(astart)};
fpm::fixed_16_16 sina {sin(astart)};
fpm::fixed_16_16 x {xP * cosa + xQ * sina};
fpm::fixed_16_16 y {yP * cosa + yQ * sina};
xQ = xQ * cosa - xP * sina;
yQ = yQ * cosa - yP * sina;
xP = x;
yP = y;
}
// If sweep angle is negative, switch direction
if (asweep < 0)
{
xQ = -xQ;
yQ = -yQ;
asweep = -asweep;
}
if(approximationScale < minCurveApproximationScale || approximationScale > maxCurveApproximationScale)
{
approximationScale = VAbstractApplication::VApp()->Settings()->GetCurveApproximationScale();
}
fpm::fixed_16_16 flatness {maxCurveApproximationScale / approximationScale * tolerance};
fpm::fixed_16_16 swangle{asweep};
QVector<QPointF> arc = EllipseCore(xC, yC, xP, yP, xQ, yQ, swangle, flatness);
// Arc end point
fpm::fixed_16_16 cosb {qCos(asweep)};
fpm::fixed_16_16 sinb {qSin(asweep)};
xP = xP*cosb + xQ*sinb;
yP = yP*cosb + yQ*sinb;
arc.append({static_cast<qreal>(xP+xC), static_cast<qreal>(yP+yC)});
return arc;
}
} // namespace
//---------------------------------------------------------------------------------------------------------------------
/**
@ -59,11 +221,11 @@ VEllipticalArc::VEllipticalArc()
* @param f2 end angle (degree).
*/
VEllipticalArc::VEllipticalArc (const VPointF &center, qreal radius1, qreal radius2, const QString &formulaRadius1,
const QString &formulaRadius2, qreal f1, const QString &formulaF1, qreal f2,
const QString &formulaF2, qreal rotationAngle, const QString &formulaRotationAngle,
quint32 idObject, Draw mode)
const QString &formulaRadius2, qreal f1, const QString &formulaF1, qreal f2,
const QString &formulaF2, qreal rotationAngle, const QString &formulaRotationAngle,
quint32 idObject, Draw mode)
: VAbstractArc(GOType::EllipticalArc, center, f1, formulaF1, f2, formulaF2, idObject, mode),
d (new VEllipticalArcData(radius1, radius2, formulaRadius1, formulaRadius2, rotationAngle, formulaRotationAngle))
d (new VEllipticalArcData(radius1, radius2, formulaRadius1, formulaRadius2, rotationAngle, formulaRotationAngle))
{
CreateName();
}
@ -72,7 +234,7 @@ VEllipticalArc::VEllipticalArc (const VPointF &center, qreal radius1, qreal radi
VEllipticalArc::VEllipticalArc(const VPointF &center, qreal radius1, qreal radius2, qreal f1, qreal f2,
qreal rotationAngle)
: VAbstractArc(GOType::EllipticalArc, center, f1, f2, NULL_ID, Draw::Calculation),
d (new VEllipticalArcData(radius1, radius2, rotationAngle))
d (new VEllipticalArcData(radius1, radius2, rotationAngle))
{
CreateName();
}
@ -83,7 +245,7 @@ VEllipticalArc::VEllipticalArc(qreal length, const QString &formulaLength, const
const QString &formulaF1, qreal rotationAngle, const QString &formulaRotationAngle,
quint32 idObject, Draw mode)
: VAbstractArc(GOType::EllipticalArc, formulaLength, center, f1, formulaF1, idObject, mode),
d (new VEllipticalArcData(radius1, radius2, formulaRadius1, formulaRadius2, rotationAngle, formulaRotationAngle))
d (new VEllipticalArcData(radius1, radius2, formulaRadius1, formulaRadius2, rotationAngle, formulaRotationAngle))
{
CreateName();
FindF2(length);
@ -93,7 +255,7 @@ VEllipticalArc::VEllipticalArc(qreal length, const QString &formulaLength, const
VEllipticalArc::VEllipticalArc(qreal length, const VPointF &center, qreal radius1, qreal radius2, qreal f1,
qreal rotationAngle)
: VAbstractArc(GOType::EllipticalArc, center, f1, NULL_ID, Draw::Calculation),
d (new VEllipticalArcData(radius1, radius2, rotationAngle))
d (new VEllipticalArcData(radius1, radius2, rotationAngle))
{
CreateName();
FindF2(length);
@ -193,7 +355,7 @@ auto VEllipticalArc::Move(qreal length, qreal angle, const QString &prefix) cons
const VPointF center = oldCenter.Move(length, angle);
const QPointF position = d->m_transform.inverted().map(center.toQPointF()) -
d->m_transform.inverted().map(oldCenter.toQPointF());
d->m_transform.inverted().map(oldCenter.toQPointF());
QTransform t = d->m_transform;
t.translate(position.x(), position.y());
@ -285,45 +447,19 @@ auto VEllipticalArc::GetCenter() const -> VPointF
auto VEllipticalArc::GetPoints() const -> QVector<QPointF>
{
const QPointF center = VAbstractArc::GetCenter().toQPointF();
QRectF box(center.x() - d->radius1, center.y() - d->radius2, d->radius1*2, d->radius2*2);
QLineF startLine(center.x(), center.y(), center.x() + d->radius1, center.y());
QLineF endLine = startLine;
startLine.setAngle(VAbstractArc::GetStartAngle());
endLine.setAngle(RealEndAngle());
qreal sweepAngle = startLine.angleTo(endLine);
if (qFuzzyIsNull(sweepAngle))
{
sweepAngle = 360;
}
QPainterPath path;
path.moveTo(GetP1());
path.arcTo(box, VAbstractArc::GetStartAngle(), sweepAngle);
path.moveTo(GetP2());
// Generate complete ellipse because angles are not correct and have to be fixed manually
QVector<QPointF> points = EllipticArcPoints(center, d->radius1, d->radius2, 0.0, M_2PI, GetApproximationScale());
points = ArcPoints(points);
QTransform t = d->m_transform;
t.translate(center.x(), center.y());
t.rotate(-GetRotationAngle());
t.translate(-center.x(), -center.y());
path = t.map(path);
std::transform(points.begin(), points.end(), points.begin(), [t](const QPointF &point) { return t.map(point); });
QPolygonF polygon;
const QList<QPolygonF> sub = path.toSubpathPolygons();
if (not sub.isEmpty())
{
polygon = ConstFirst (path.toSubpathPolygons());
if (not polygon.isEmpty() && not VFuzzyComparePoints(GetP1(), ConstFirst<QPointF> (polygon)))
{
polygon.removeFirst(); // remove point (0;0)
}
}
return static_cast<QVector<QPointF>>(polygon);
return points;
}
//---------------------------------------------------------------------------------------------------------------------
@ -360,7 +496,7 @@ auto VEllipticalArc::CutArc(const qreal &length, VEllipticalArc &arc1, VElliptic
const QString errorMsg = QObject::tr("Unable to cut curve '%1'. The curve is too short.").arg(name());
VAbstractApplication::VApp()->IsPedantic() ? throw VException(errorMsg) :
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
return {};
}
@ -380,10 +516,10 @@ auto VEllipticalArc::CutArc(const qreal &length, VEllipticalArc &arc1, VElliptic
else
{
errorMsg = QObject::tr("Curve '%1'. Length of a cut segment is too small. Optimize it to minimal value.")
.arg(name());
.arg(name());
}
VAbstractApplication::VApp()->IsPedantic() ? throw VException(errorMsg) :
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
}
else if (length > maxLength)
{
@ -393,15 +529,15 @@ auto VEllipticalArc::CutArc(const qreal &length, VEllipticalArc &arc1, VElliptic
if (not pointName.isEmpty())
{
errorMsg = QObject::tr("Curve '%1'. Length of a cut segment (%2) is too big. Optimize it to maximal value.")
.arg(name(), pointName);
.arg(name(), pointName);
}
else
{
errorMsg = QObject::tr("Curve '%1'. Length of a cut segment is too big. Optimize it to maximal value.")
.arg(name());
.arg(name());
}
VAbstractApplication::VApp()->IsPedantic() ? throw VException(errorMsg) :
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
qWarning() << VAbstractApplication::warningMessageSignature + errorMsg;
}
else
{
@ -410,13 +546,13 @@ auto VEllipticalArc::CutArc(const qreal &length, VEllipticalArc &arc1, VElliptic
// the first arc has given length and startAngle just like in the origin arc
arc1 = VEllipticalArc (len, QString().setNum(length), GetCenter(), d->radius1, d->radius2,
d->formulaRadius1, d->formulaRadius2, GetStartAngle(), GetFormulaF1(), d->rotationAngle,
GetFormulaRotationAngle(), getIdObject(), getMode());
d->formulaRadius1, d->formulaRadius2, GetStartAngle(), GetFormulaF1(), d->rotationAngle,
GetFormulaRotationAngle(), getIdObject(), getMode());
// the second arc has startAngle just like endAngle of the first arc
// and it has endAngle just like endAngle of the origin arc
arc2 = VEllipticalArc (GetCenter(), d->radius1, d->radius2, d->formulaRadius1, d->formulaRadius2,
arc1.GetEndAngle(), arc1.GetFormulaF2(), GetEndAngle(), GetFormulaF2(), d->rotationAngle,
GetFormulaRotationAngle(), getIdObject(), getMode());
arc1.GetEndAngle(), arc1.GetFormulaF2(), GetEndAngle(), GetFormulaF2(), d->rotationAngle,
GetFormulaRotationAngle(), getIdObject(), getMode());
return arc1.GetP1();
}
@ -432,7 +568,7 @@ auto VEllipticalArc::CutArc(const qreal &length, const QString &pointName) const
//---------------------------------------------------------------------------------------------------------------------
void VEllipticalArc::CreateName()
{
QString name = ELARC_ + QString("%1").arg(this->GetCenter().name());
QString name = ELARC_ + QStringLiteral("%1").arg(this->GetCenter().name());
const QString nameStr = QStringLiteral("_%1");
if (getMode() == Draw::Modeling && getIdObject() != NULL_ID)
@ -555,23 +691,101 @@ auto VEllipticalArc::GetP(qreal angle) const -> QPointF
}
//---------------------------------------------------------------------------------------------------------------------
auto VEllipticalArc::RealEndAngle() const -> qreal
auto VEllipticalArc::ArcPoints(QVector<QPointF> points) const -> QVector<QPointF>
{
qreal endAngle = VEllipticalArc::OptimizeAngle(VAbstractArc::GetEndAngle());
if (qFuzzyIsNull(endAngle) ||
VFuzzyComparePossibleNulls(endAngle, 90) ||
VFuzzyComparePossibleNulls(endAngle, 180) ||
VFuzzyComparePossibleNulls(endAngle, 270) ||
VFuzzyComparePossibleNulls(endAngle, 360))
if (points.size() < 2 || VFuzzyComparePossibleNulls(VAbstractArc::GetStartAngle(), VAbstractArc::GetEndAngle()))
{
return endAngle;
return points;
}
endAngle = qRadiansToDegrees(qAtan2(d->radius1 * qSin(qDegreesToRadians(endAngle)),
d->radius2 * qCos(qDegreesToRadians(endAngle))));
QPointF center = static_cast<QPointF>(GetCenter());
qreal radius = qMax(d->radius1, d->radius2) * 2;
return endAngle;
QLineF start(center.x(), center.y(), center.x() + radius, center.y());
start.setAngle(VAbstractArc::GetStartAngle());
QLineF end(center.x(), center.y(), center.x() + radius, center.y());
end.setAngle(VAbstractArc::GetEndAngle());
if (VFuzzyComparePossibleNulls(start.angle(), end.angle()))
{
return points;
}
auto IsBoundedIntersection = [](QLineF::IntersectType type, QPointF p, const QLineF &segment1,
const QLineF &segment2)
{
return type == QLineF::BoundedIntersection ||
(VGObject::IsPointOnLineSegment (p, segment1.p1(), segment2.p1()) &&
VGObject::IsPointOnLineSegment (p, segment2.p1(), segment2.p2()));
};
bool begin = true;
if (start.angle() > end.angle())
{
for (int i=0; i < points.size()-1; ++i)
{
QLineF edge(points.at(i), points.at(i+1));
QPointF p;
QLineF::IntersectType type = Intersects(start, edge, &p);
// QLineF::intersects not always accurate on edge cases
if (IsBoundedIntersection(type, p, edge, start))
{
QVector<QPointF> head = points.mid(0, i+1);
QVector<QPointF> tail = points.mid(i+1, -1);
tail.prepend(p);
points = tail + head;
begin = false;
break;
}
}
}
QVector<QPointF> arc;
arc.reserve(points.size());
for (int i=0; i < points.size()-1; ++i)
{
QLineF edge(points.at(i), points.at(i+1));
if (begin)
{
QPointF p;
QLineF::IntersectType type = Intersects(start, edge, &p);
// QLineF::intersects not always accurate on edge cases
if (IsBoundedIntersection(type, p, edge, start))
{
arc.append(p);
begin = false;
}
}
else
{
QPointF p;
QLineF::IntersectType type = Intersects(end, edge, &p);
// QLineF::intersects not always accurate on edge cases
if (IsBoundedIntersection(type, p, edge, end))
{
arc.append(points.at(i));
arc.append(p);
break;
}
arc.append(points.at(i));
}
}
if (arc.isEmpty())
{
return points;
}
return arc;
}
//---------------------------------------------------------------------------------------------------------------------

View File

@ -114,10 +114,8 @@ private:
QSharedDataPointer<VEllipticalArcData> d;
auto MaxLength() const -> qreal;
auto GetP(qreal angle) const -> QPointF;
auto RealEndAngle() const -> qreal;
auto ArcPoints(QVector<QPointF> points) const -> QVector<QPointF>;
};
Q_DECLARE_METATYPE(VEllipticalArc) // NOLINT

View File

@ -64,13 +64,6 @@
#include "../ifc/exception/vexception.h"
#include "literals.h"
const qreal defCurveApproximationScale = 0.5;
const qreal minCurveApproximationScale = 0.2;
const qreal maxCurveApproximationScale = 10.0;
const int minLabelFontSize = 5;
const int maxLabelFontSize = 100;
//---------------------------------------------------------------------------------------------------------------------
QPixmap QPixmapFromCache(const QString &pixmapPath)
{

View File

@ -74,13 +74,14 @@ class QMarginsF;
class VTranslateMeasurements;
class QGraphicsItem;
#define SceneSize 50000
extern const qreal defCurveApproximationScale;
extern const qreal minCurveApproximationScale;
extern const qreal maxCurveApproximationScale;
constexpr qreal maxSceneSize = ((20.0 * 1000.0) / 25.4) * PrintDPI; // 20 meters in pixels
extern const int minLabelFontSize;
extern const int maxLabelFontSize;
constexpr qreal defCurveApproximationScale = 0.5;
constexpr qreal minCurveApproximationScale = 0.2;
constexpr qreal maxCurveApproximationScale = 10.0;
constexpr int minLabelFontSize = 5;
constexpr int maxLabelFontSize = 100;
enum class NodeDetail : qint8 { Contour, Modeling };
enum class SceneObject : qint8 { Point, Line, Spline, Arc, ElArc, SplinePath, Detail, Unknown };

View File

@ -0,0 +1,21 @@
MIT License
Copyright (c) 2019 Mike Lankamp
Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.

View File

@ -0,0 +1,497 @@
#ifndef FPM_FIXED_HPP
#define FPM_FIXED_HPP
#include <cassert>
#include <cmath>
#include <cstdint>
#include <functional>
#include <limits>
#include <type_traits>
#include <QtGlobal>
#if QT_VERSION < QT_VERSION_CHECK(5, 5, 0)
#include "../diagnostic.h"
#endif // QT_VERSION < QT_VERSION_CHECK(5, 5, 0)
namespace fpm
{
//! Fixed-point number type
//! \tparam BaseType the base integer type used to store the fixed-point number. This can be a signed or
//! unsigned type.
//! \tparam IntermediateType the integer type used to store intermediate results during calculations.
//! \tparam FractionBits the number of bits of the BaseType used to store the fraction
template <typename BaseType, typename IntermediateType, unsigned int FractionBits>
class fixed
{
static_assert(std::is_integral<BaseType>::value, "BaseType must be an integral type");
static_assert(FractionBits > 0, "FractionBits must be greater than zero");
static_assert(FractionBits <= sizeof(BaseType) * 8 - 1, "BaseType must at least be able to contain entire "
"fraction, with space for at least one integral bit");
static_assert(sizeof(IntermediateType) > sizeof(BaseType), "IntermediateType must be larger than BaseType");
static_assert(std::is_signed<IntermediateType>::value == std::is_signed<BaseType>::value,
"IntermediateType must have same signedness as BaseType");
// Although this value fits in the BaseType in terms of bits, if there's only one integral bit, this value
// is incorrect (flips from positive to negative), so we must extend the size to IntermediateType.
static constexpr IntermediateType FRACTION_MULT = IntermediateType(1) << FractionBits;
struct raw_construct_tag {};
constexpr inline fixed(BaseType val, raw_construct_tag /*unused*/) noexcept : m_value(val) {}
public:
inline fixed() noexcept = default;
// Converts an integral number to the fixed-point type.
// Like static_cast, this truncates bits that don't fit.
template <typename T, typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline explicit fixed(T val) noexcept
: m_value(static_cast<BaseType>(val * FRACTION_MULT))
{}
// Converts an floating-point number to the fixed-point type.
// Like static_cast, this truncates bits that don't fit.
template <typename T, typename std::enable_if<std::is_floating_point<T>::value>::type* = nullptr>
constexpr inline explicit fixed(T val) noexcept
: m_value(static_cast<BaseType>((val >= 0.0) ? (val * FRACTION_MULT + T{0.5}) : (val * FRACTION_MULT - T{0.5})))
{}
// Constructs from another fixed-point type with possibly different underlying representation.
// Like static_cast, this truncates bits that don't fit.
template <typename B, typename I, unsigned int F>
constexpr inline explicit fixed(fixed<B,I,F> val) noexcept
: m_value(from_fixed_point<F>(val.raw_value()).raw_value())
{}
// Explicit conversion to a floating-point type
template <typename T, typename std::enable_if<std::is_floating_point<T>::value>::type* = nullptr>
constexpr inline explicit operator T() const noexcept
{
return static_cast<T>(m_value) / FRACTION_MULT;
}
// Explicit conversion to an integral type
template <typename T, typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline explicit operator T() const noexcept
{
return static_cast<T>(m_value / FRACTION_MULT);
}
// Returns the raw underlying value of this type.
// Do not use this unless you know what you're doing.
constexpr inline auto raw_value() const noexcept -> BaseType
{
return m_value;
}
//! Constructs a fixed-point number from another fixed-point number.
//! \tparam NumFractionBits the number of bits used by the fraction in \a value.
//! \param value the integer fixed-point number
template <unsigned int NumFractionBits, typename T,
typename std::enable_if<(NumFractionBits > FractionBits)>::type* = nullptr>
static constexpr inline auto from_fixed_point(T value) noexcept -> fixed
{
// To correctly round the last bit in the result, we need one more bit of information.
// We do this by multiplying by two before dividing and adding the LSB to the real result.
return fixed(static_cast<BaseType>(
value / (T(1) << (NumFractionBits - FractionBits)) +
(value / (T(1) << (NumFractionBits - FractionBits - 1)) % 2)),
raw_construct_tag{});
}
template <unsigned int NumFractionBits, typename T,
typename std::enable_if<(NumFractionBits <= FractionBits)>::type* = nullptr>
static constexpr inline auto from_fixed_point(T value) noexcept -> fixed
{
return fixed(static_cast<BaseType>(
value * (T(1) << (FractionBits - NumFractionBits))),
raw_construct_tag{});
}
// Constructs a fixed-point number from its raw underlying value.
// Do not use this unless you know what you're doing.
static constexpr inline auto from_raw_value(BaseType value) noexcept -> fixed
{
return fixed(value, raw_construct_tag{});
}
//
// Constants
//
static constexpr auto e() -> fixed { return from_fixed_point<61>(6267931151224907085LL); }
static constexpr auto pi() -> fixed { return from_fixed_point<61>(7244019458077122842LL); }
static constexpr auto half_pi() -> fixed { return from_fixed_point<62>(7244019458077122842LL); }
static constexpr auto two_pi() -> fixed { return from_fixed_point<60>(7244019458077122842LL); }
//
// Arithmetic member operators
//
constexpr inline auto operator-() const noexcept -> fixed
{
return fixed::from_raw_value(-m_value);
}
inline auto operator+=(const fixed& y) noexcept -> fixed&
{
m_value += y.m_value;
return *this;
}
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator+=(I y) noexcept -> fixed&
{
m_value += y * FRACTION_MULT;
return *this;
}
inline auto operator-=(const fixed& y) noexcept -> fixed&
{
m_value -= y.m_value;
return *this;
}
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator-=(I y) noexcept -> fixed&
{
m_value -= y * FRACTION_MULT;
return *this;
}
inline auto operator*=(const fixed& y) noexcept -> fixed&
{
// Normal fixed-point multiplication is: x * y / 2**FractionBits.
// To correctly round the last bit in the result, we need one more bit of information.
// We do this by multiplying by two before dividing and adding the LSB to the real result.
auto value = static_cast<IntermediateType>(
static_cast<double>(static_cast<IntermediateType>(m_value) * y.m_value) / (FRACTION_MULT / 2));
m_value = static_cast<BaseType>((value / 2) + (value % 2));
return *this;
}
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator*=(I y) noexcept -> fixed&
{
m_value *= y;
return *this;
}
inline auto operator/=(const fixed& y) noexcept -> fixed&
{
assert(y.m_value != 0);
// Normal fixed-point division is: x * 2**FractionBits / y.
// To correctly round the last bit in the result, we need one more bit of information.
// We do this by multiplying by two before dividing and adding the LSB to the real result.
auto value = (static_cast<IntermediateType>(m_value) * FRACTION_MULT * 2) / y.m_value;
m_value = static_cast<BaseType>((value / 2) + (value % 2));
return *this;
}
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator/=(I y) noexcept -> fixed&
{
m_value /= y;
return *this;
}
//
// Shift operators
//
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator>>=(I y) noexcept -> fixed&
{
m_value >>= y;
return *this;
}
template <typename I, typename std::enable_if<std::is_integral<I>::value>::type* = nullptr>
inline auto operator<<=(I y) noexcept -> fixed&
{
m_value <<= y;
return *this;
}
private:
BaseType m_value;
};
//
// Convenience typedefs
//
using fixed_16_16 = fixed<std::int32_t, std::int64_t, 16>;
using fixed_24_8 = fixed<std::int32_t, std::int64_t, 8>;
using fixed_8_24 = fixed<std::int32_t, std::int64_t, 24>;
//
// Addition
//
template <typename B, typename I, unsigned int F>
constexpr inline auto operator+(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) += y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator+(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) += y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator+(T x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(y) += x;
}
//
// Subtraction
//
template <typename B, typename I, unsigned int F>
constexpr inline auto operator-(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) -= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator-(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) -= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator-(T x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) -= y;
}
//
// Multiplication
//
template <typename B, typename I, unsigned int F>
constexpr inline auto operator*(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) *= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator*(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) *= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator*(T x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(y) *= x;
}
//
// Division
//
template <typename B, typename I, unsigned int F>
constexpr inline auto operator/(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) /= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator/(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) /= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator/(T x, const fixed<B, I, F>& y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) /= y;
}
//
// Shift operators
//
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator>>(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) >>= y;
}
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
constexpr inline auto operator<<(const fixed<B, I, F>& x, T y) noexcept -> fixed<B, I, F>
{
return fixed<B, I, F>(x) <<= y;
}
//
// Comparison operators
//
template <typename B, typename I, unsigned int F>
constexpr inline auto operator==(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() == y.raw_value();
}
template <typename B, typename I, unsigned int F>
constexpr inline auto operator!=(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() != y.raw_value();
}
template <typename B, typename I, unsigned int F>
constexpr inline auto operator<(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() < y.raw_value();
}
template <typename B, typename I, unsigned int F>
constexpr inline auto operator>(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() > y.raw_value();
}
template <typename B, typename I, unsigned int F>
constexpr inline auto operator<=(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() <= y.raw_value();
}
template <typename B, typename I, unsigned int F>
constexpr inline auto operator>=(const fixed<B, I, F>& x, const fixed<B, I, F>& y) noexcept -> bool
{
return x.raw_value() >= y.raw_value();
}
namespace detail
{
QT_WARNING_PUSH
QT_WARNING_DISABLE_CLANG("-Wunneeded-internal-declaration")
// Number of base-10 digits required to fully represent a number of bits
static constexpr auto max_digits10(int bits) -> int
{
// 8.24 fixed-point equivalent of (int)ceil(bits * std::log10(2));
using T = long long; // NOLINT(google-runtime-int)
return static_cast<int>((T{bits} * 5050445 + (T{1} << 24) - 1) >> 24); // NOLINT(hicpp-signed-bitwise)
}
// Number of base-10 digits that can be fully represented by a number of bits
static constexpr auto digits10(int bits) -> int
{
// 8.24 fixed-point equivalent of (int)(bits * std::log10(2));
using T = long long; // NOLINT(google-runtime-int)
return static_cast<int>((T{bits} * 5050445) >> 24); //NOLINT(hicpp-signed-bitwise)
}
QT_WARNING_POP
} // namespace detail
} // namespace fpm
// Specializations for customization points
namespace std // NOLINT(cert-dcl58-cpp)
{
template <typename B, typename I, unsigned int F>
struct hash<fpm::fixed<B,I,F>>
{
using argument_type = fpm::fixed<B, I, F>;
using result_type = std::size_t;
auto operator()(const argument_type &arg) const noexcept(noexcept(std::declval<std::hash<B>>()(arg.raw_value())))
-> result_type
{
return m_hash(arg.raw_value());
}
private:
std::hash<B> m_hash;
};
template <typename B, typename I, unsigned int F>
struct numeric_limits<fpm::fixed<B,I,F>>
{
static constexpr bool is_specialized = true;
static constexpr bool is_signed = std::numeric_limits<B>::is_signed;
static constexpr bool is_integer = false;
static constexpr bool is_exact = true;
static constexpr bool has_infinity = false;
static constexpr bool has_quiet_NaN = false;
static constexpr bool has_signaling_NaN = false;
static constexpr bool has_denorm = std::denorm_absent;
static constexpr bool has_denorm_loss = false;
static constexpr std::float_round_style round_style = std::round_to_nearest;
static constexpr bool is_iec_559 = false;
static constexpr bool is_bounded = true;
static constexpr bool is_modulo = std::numeric_limits<B>::is_modulo;
static constexpr int digits = std::numeric_limits<B>::digits;
// Any number with `digits10` significant base-10 digits (that fits in
// the range of the type) is guaranteed to be convertible from text and
// back without change. Worst case, this is 0.000...001, so we can only
// guarantee this case. Nothing more.
static constexpr int digits10 = 1;
// This is equal to max_digits10 for the integer and fractional part together.
static constexpr int max_digits10 =
fpm::detail::max_digits10(std::numeric_limits<B>::digits - F) + fpm::detail::max_digits10(F);
static constexpr int radix = 2;
static constexpr int min_exponent = 1 - F;
static constexpr int min_exponent10 = -fpm::detail::digits10(F);
static constexpr int max_exponent = std::numeric_limits<B>::digits - F;
static constexpr int max_exponent10 = fpm::detail::digits10(std::numeric_limits<B>::digits - F);
static constexpr bool traps = true;
static constexpr bool tinyness_before = false;
static constexpr auto lowest() noexcept -> fpm::fixed<B,I,F>
{
return fpm::fixed<B,I,F>::from_raw_value(std::numeric_limits<B>::lowest());
}
static constexpr auto min() noexcept -> fpm::fixed<B,I,F>
{
return lowest();
}
static constexpr auto max() noexcept -> fpm::fixed<B,I,F>
{
return fpm::fixed<B,I,F>::from_raw_value(std::numeric_limits<B>::max());
}
static constexpr auto epsilon() noexcept -> fpm::fixed<B,I,F>
{
return fpm::fixed<B,I,F>::from_raw_value(1);
}
static constexpr auto round_error() noexcept -> fpm::fixed<B,I,F>
{
return fpm::fixed<B,I,F>(1) / 2;
}
static constexpr auto denorm_min() noexcept -> fpm::fixed<B,I,F>
{
return min();
}
};
} // namespace std
#endif

700
src/libs/vmisc/fpm/math.hpp Normal file
View File

@ -0,0 +1,700 @@
#ifndef FPM_MATH_HPP
#define FPM_MATH_HPP
#include "fixed.hpp"
#include <cmath>
#ifdef _MSC_VER
#include <intrin.h>
#endif
#include <QtGlobal>
#if QT_VERSION < QT_VERSION_CHECK(5, 5, 0)
#include "../diagnostic.h"
#endif // QT_VERSION < QT_VERSION_CHECK(5, 5, 0)
namespace fpm
{
//
// Helper functions
//
namespace detail
{
QT_WARNING_PUSH
QT_WARNING_DISABLE_CLANG("-Wsign-conversion")
// Returns the index of the most-signifcant set bit
inline auto find_highest_bit(unsigned long long value) noexcept -> long
{
assert(value != 0);
#if defined(_MSC_VER)
unsigned long index;
#if defined(_WIN64)
_BitScanReverse64(&index, value);
#else
if (_BitScanReverse(&index, static_cast<unsigned long>(value >> 32)) != 0) {
index += 32;
} else {
_BitScanReverse(&index, static_cast<unsigned long>(value & 0xfffffffflu));
}
#endif
return index;
#elif defined(__GNUC__) || defined(__clang__)
return sizeof(value) * 8 - 1 - __builtin_clzll(value); // NOLINT
#else
# error "your platform does not support find_highest_bit()"
#endif
}
QT_WARNING_POP
} // namespace detail
//
// Classification methods
//
template <typename B, typename I, unsigned int F>
constexpr inline auto fpclassify(fixed<B, I, F> x) noexcept -> int
{
return (x.raw_value() == 0) ? FP_ZERO : FP_NORMAL;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isfinite(fixed<B, I, F> /*unused*/) noexcept -> bool
{
return true;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isinf(fixed<B, I, F> /*unused*/) noexcept -> bool
{
return false;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isnan(fixed<B, I, F> /*unused*/) noexcept -> bool
{
return false;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isnormal(fixed<B, I, F> x) noexcept -> bool
{
return x.raw_value() != 0;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto signbit(fixed<B, I, F> x) noexcept -> bool
{
return x.raw_value() < 0;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isgreater(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> bool
{
return x > y;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isgreaterequal(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> bool
{
return x >= y;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isless(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> bool
{
return x < y;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto islessequal(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> bool
{
return x <= y;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto islessgreater(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> bool
{
return x != y;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto isunordered(fixed<B, I, F> /*x*/, fixed<B, I, F> /*y*/) noexcept -> bool
{
return false;
}
//
// Nearest integer operations
//
template <typename B, typename I, unsigned int F>
inline auto ceil(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
constexpr auto FRAC = B(1) << F;
auto value = x.raw_value();
if (value > 0)
{
value += FRAC - 1;
}
return fixed<B, I, F>::from_raw_value(value / FRAC * FRAC);
}
template <typename B, typename I, unsigned int F>
inline auto floor(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
constexpr auto FRAC = B(1) << F;
auto value = x.raw_value();
if (value < 0)
{
value -= FRAC - 1;
}
return fixed<B, I, F>::from_raw_value(value / FRAC * FRAC);
}
template <typename B, typename I, unsigned int F>
inline auto trunc(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
constexpr auto FRAC = B(1) << F;
return fixed<B, I, F>::from_raw_value(x.raw_value() / FRAC * FRAC);
}
template <typename B, typename I, unsigned int F>
inline auto round(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
constexpr auto FRAC = B(1) << F;
auto value = x.raw_value() / (FRAC / 2);
return fixed<B, I, F>::from_raw_value(((value / 2) + (value % 2)) * FRAC);
}
template <typename B, typename I, unsigned int F>
auto nearbyint(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
// Rounding mode is assumed to be FE_TONEAREST
constexpr auto FRAC = B(1) << F;
auto value = x.raw_value();
const bool is_half = std::abs(value % FRAC) == FRAC / 2;
value /= FRAC / 2;
value = (value / 2) + (value % 2);
value -= (value % 2) * is_half;
return fixed<B, I, F>::from_raw_value(value * FRAC);
}
template <typename B, typename I, unsigned int F>
constexpr inline auto rint(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
// Rounding mode is assumed to be FE_TONEAREST
return nearbyint(x);
}
//
// Mathematical functions
//
template <typename B, typename I, unsigned int F>
constexpr inline auto abs(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
return (x >= fixed<B, I, F>{0}) ? x : -x;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto fmod(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> fixed<B, I, F>
{
return
assert(y.raw_value() != 0),
fixed<B, I, F>::from_raw_value(x.raw_value() % y.raw_value());
}
template <typename B, typename I, unsigned int F>
constexpr inline auto remainder(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> fixed<B, I, F>
{
return
assert(y.raw_value() != 0),
x - nearbyint(x / y) * y;
}
template <typename B, typename I, unsigned int F>
inline auto remquo(fixed<B, I, F> x, fixed<B, I, F> y, int* quo) noexcept -> fixed<B, I, F>
{
assert(y.raw_value() != 0);
assert(quo != nullptr);
*quo = x.raw_value() / y.raw_value();
return fixed<B, I, F>::from_raw_value(x.raw_value() % y.raw_value());
}
//
// Manipulation functions
//
template <typename B, typename I, unsigned int F, typename C, typename J, unsigned int G>
constexpr inline auto copysign(fixed<B, I, F> x, fixed<C, J, G> y) noexcept -> fixed<B, I, F>
{
return
x = abs(x),
(y >= fixed<C, J, G>{0}) ? x : -x;
}
template <typename B, typename I, unsigned int F>
constexpr inline auto nextafter(fixed<B, I, F> from, fixed<B, I, F> to) noexcept -> fixed<B, I, F>
{
return from == to ? to :
to > from ? fixed<B, I, F>::from_raw_value(from.raw_value() + 1)
: fixed<B, I, F>::from_raw_value(from.raw_value() - 1);
}
template <typename B, typename I, unsigned int F>
constexpr inline auto nexttoward(fixed<B, I, F> from, fixed<B, I, F> to) noexcept -> fixed<B, I, F>
{
return nextafter(from, to);
}
template <typename B, typename I, unsigned int F>
inline auto modf(fixed<B, I, F> x, fixed<B, I, F>* iptr) noexcept -> fixed<B, I, F>
{
const auto raw = x.raw_value();
constexpr auto FRAC = B{1} << F;
*iptr = fixed<B, I, F>::from_raw_value(raw / FRAC * FRAC);
return fixed<B, I, F>::from_raw_value(raw % FRAC);
}
//
// Power functions
//
template <typename B, typename I, unsigned int F, typename T,
typename std::enable_if<std::is_integral<T>::value>::type* = nullptr>
auto pow(fixed<B, I, F> base, T exp) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (base == Fixed(0)) {
assert(exp > 0);
return Fixed(0);
}
Fixed result {1};
if (exp < 0)
{
for (Fixed intermediate = base; exp != 0; exp /= 2, intermediate *= intermediate)
{
if ((exp % 2) != 0)
{
result /= intermediate;
}
}
}
else
{
for (Fixed intermediate = base; exp != 0; exp /= 2, intermediate *= intermediate)
{
if ((exp % 2) != 0)
{
result *= intermediate;
}
}
}
return result;
}
template <typename B, typename I, unsigned int F>
auto pow(fixed<B, I, F> base, fixed<B, I, F> exp) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (base == Fixed(0)) {
assert(exp > Fixed(0));
return Fixed(0);
}
if (exp < Fixed(0))
{
return 1 / pow(base, -exp);
}
constexpr auto FRAC = B(1) << F;
if (exp.raw_value() % FRAC == 0)
{
// Non-fractional exponents are easier to calculate
return pow(base, exp.raw_value() / FRAC);
}
// For negative bases we do not support fractional exponents.
// Technically fractions with odd denominators could work,
// but that's too much work to figure out.
assert(base > Fixed(0));
return exp2(log2(base) * exp);
}
template <typename B, typename I, unsigned int F>
auto exp(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (x < Fixed(0)) {
return 1 / exp(-x);
}
constexpr auto FRAC = B(1) << F;
const B x_int = x.raw_value() / FRAC;
x -= x_int;
assert(x >= Fixed(0) && x < Fixed(1));
constexpr auto fA = Fixed::template from_fixed_point<63>( 128239257017632854LL); // 1.3903728105644451e-2
constexpr auto fB = Fixed::template from_fixed_point<63>( 320978614890280666LL); // 3.4800571158543038e-2
constexpr auto fC = Fixed::template from_fixed_point<63>(1571680799599592947LL); // 1.7040197373796334e-1
constexpr auto fD = Fixed::template from_fixed_point<63>(4603349000587966862LL); // 4.9909609871464493e-1
constexpr auto fE = Fixed::template from_fixed_point<62>(4612052447974689712LL); // 1.0000794567422495
constexpr auto fF = Fixed::template from_fixed_point<63>(9223361618412247875LL); // 9.9999887043019773e-1
return pow(Fixed::e(), x_int) * (((((fA * x + fB) * x + fC) * x + fD) * x + fE) * x + fF);
}
template <typename B, typename I, unsigned int F>
auto exp2(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (x < Fixed(0)) {
return 1 / exp2(-x);
}
constexpr auto FRAC = B(1) << F;
const B x_int = x.raw_value() / FRAC;
x -= x_int;
assert(x >= Fixed(0) && x < Fixed(1));
constexpr auto fA = Fixed::template from_fixed_point<63>( 17491766697771214LL); // 1.8964611454333148e-3
constexpr auto fB = Fixed::template from_fixed_point<63>( 82483038782406547LL); // 8.9428289841091295e-3
constexpr auto fC = Fixed::template from_fixed_point<63>( 515275173969157690LL); // 5.5866246304520701e-2
constexpr auto fD = Fixed::template from_fixed_point<63>(2214897896212987987LL); // 2.4013971109076949e-1
constexpr auto fE = Fixed::template from_fixed_point<63>(6393224161192452326LL); // 6.9315475247516736e-1
constexpr auto fF = Fixed::template from_fixed_point<63>(9223371050976163566LL); // 9.9999989311082668e-1
return Fixed(1 << x_int) * (((((fA * x + fB) * x + fC) * x + fD) * x + fE) * x + fF);
}
template <typename B, typename I, unsigned int F>
auto expm1(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
return exp(x) - 1; // cppcheck-suppress unpreciseMathCall
}
template <typename B, typename I, unsigned int F>
auto log2(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x > Fixed(0));
// Normalize input to the [1:2] domain
B value = x.raw_value();
const long highest = detail::find_highest_bit(value);
if (highest >= F)
{
value >>= (highest - F);
}
else
{
value <<= (F - highest);
}
x = Fixed::from_raw_value(value);
assert(x >= Fixed(1) && x < Fixed(2));
constexpr auto fA = Fixed::template from_fixed_point<63>( 413886001457275979LL); // 4.4873610194131727e-2
constexpr auto fB = Fixed::template from_fixed_point<63>(-3842121857793256941LL); // -4.1656368651734915e-1
constexpr auto fC = Fixed::template from_fixed_point<62>( 7522345947206307744LL); // 1.6311487636297217
constexpr auto fD = Fixed::template from_fixed_point<61>(-8187571043052183818LL); // -3.5507929249026341
constexpr auto fE = Fixed::template from_fixed_point<60>( 5870342889289496598LL); // 5.0917108110420042
constexpr auto fF = Fixed::template from_fixed_point<61>(-6457199832668582866LL); // -2.8003640347009253
return Fixed(highest - F) + (((((fA * x + fB) * x + fC) * x + fD) * x + fE) * x + fF);
}
template <typename B, typename I, unsigned int F>
auto log(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
return log2(x) / log2(Fixed::e());
}
template <typename B, typename I, unsigned int F>
auto log10(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
return log2(x) / log2(Fixed(10));
}
template <typename B, typename I, unsigned int F>
auto log1p(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
return log(1 + x); // cppcheck-suppress unpreciseMathCall
}
template <typename B, typename I, unsigned int F>
auto cbrt(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (x == Fixed(0))
{
return x;
}
if (x < Fixed(0))
{
return -cbrt(-x);
}
assert(x >= Fixed(0));
// Finding the cube root of an integer, taken from Hacker's Delight,
// based on the square root algorithm.
// We start at the greatest power of eight that's less than the argument.
int ofs = ((detail::find_highest_bit(x.raw_value()) + 2*F) / 3 * 3);
I num = I{x.raw_value()};
I res = 0;
const auto do_round = [&]
{
for (; ofs >= 0; ofs -= 3)
{
res += res;
const I val = (3*res*(res + 1) + 1) << ofs;
if (num >= val)
{
num -= val;
res++;
}
}
};
// We should shift by 2*F (since there are two multiplications), but that
// could overflow even the intermediate type, so we have to split the
// algorithm up in two rounds of F bits each. Each round will deplete
// 'num' digit by digit, so after a round we can shift it again.
num <<= F;
ofs -= F;
do_round();
num <<= F;
ofs += F;
do_round();
return Fixed::from_raw_value(static_cast<B>(res));
}
template <typename B, typename I, unsigned int F>
auto sqrt(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x >= Fixed(0));
if (x == Fixed(0))
{
return x;
}
// Finding the square root of an integer in base-2, from:
// https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Binary_numeral_system_.28base_2.29
// Shift by F first because it's fixed-point.
I num = I{x.raw_value()} << F;
I res = 0;
// "bit" starts at the greatest power of four that's less than the argument.
for (I bit = I{1} << ((detail::find_highest_bit(x.raw_value()) + F) / 2 * 2); bit != 0; bit >>= 2)
{
const I val = res + bit;
res >>= 1;
if (num >= val)
{
num -= val;
res += bit;
}
}
// Round the last digit up if necessary
if (num > res)
{
res++;
}
return Fixed::from_raw_value(static_cast<B>(res));
}
template <typename B, typename I, unsigned int F>
auto hypot(fixed<B, I, F> x, fixed<B, I, F> y) noexcept -> fixed<B, I, F>
{
assert(x != 0 || y != 0);
return sqrt(x*x + y*y);
}
//
// Trigonometry functions
//
template <typename B, typename I, unsigned int F>
auto sin(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
// This sine uses a fifth-order curve-fitting approximation originally
// described by Jasper Vijn on coranac.com which has a worst-case
// relative error of 0.07% (over [-pi:pi]).
using Fixed = fixed<B, I, F>;
// Turn x from [0..2*PI] domain into [0..4] domain
x = fmod(x, Fixed::two_pi());
x = x / Fixed::half_pi();
// Take x modulo one rotation, so [-4..+4].
if (x < Fixed(0)) {
x += Fixed(4);
}
int sign = +1;
if (x > Fixed(2)) {
// Reduce domain to [0..2].
sign = -1;
x -= Fixed(2);
}
if (x > Fixed(1)) {
// Reduce domain to [0..1].
x = Fixed(2) - x;
}
const Fixed x2 = x*x;
return sign * x * (Fixed::pi() - x2*(Fixed::two_pi() - 5 - x2*(Fixed::pi() - 3)))/2;
}
template <typename B, typename I, unsigned int F>
inline auto cos(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
return sin(fixed<B, I, F>::half_pi() + x);
}
template <typename B, typename I, unsigned int F>
inline auto tan(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
auto cx = cos(x);
// Tangent goes to infinity at 90 and -90 degrees.
// We can't represent that with fixed-point maths.
assert(abs(cx).raw_value() > 1);
return sin(x) / cx;
}
namespace detail {
// Calculates atan(x) assuming that x is in the range [0,1]
template <typename B, typename I, unsigned int F>
auto atan_sanitized(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x >= Fixed(0) && x <= Fixed(1));
constexpr auto fA = Fixed::template from_fixed_point<63>( 716203666280654660LL); // 0.0776509570923569
constexpr auto fB = Fixed::template from_fixed_point<63>(-2651115102768076601LL); // -0.287434475393028
constexpr auto fC = Fixed::template from_fixed_point<63>( 9178930894564541004LL); // 0.995181681698119 (PI/4 - A - B)
const auto xx = x * x;
return ((fA*xx + fB)*xx + fC)*x;
}
// Calculate atan(y / x), assuming x != 0.
//
// If x is very, very small, y/x can easily overflow the fixed-point range.
// If q = y/x and q > 1, atan(q) would calculate atan(1/q) as intermediate step
// anyway. We can shortcut that here and avoid the loss of information, thus
// improving the accuracy of atan(y/x) for very small x.
template <typename B, typename I, unsigned int F>
auto atan_div(fixed<B, I, F> y, fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x != Fixed(0));
// Make sure y and x are positive.
// If y / x is negative (when y or x, but not both, are negative), negate the result to
// keep the correct outcome.
if (y < Fixed(0)) {
if (x < Fixed(0)) {
return atan_div(-y, -x);
}
return -atan_div(-y, x);
}
if (x < Fixed(0)) {
return -atan_div(y, -x);
}
assert(y >= Fixed(0));
assert(x > Fixed(0));
if (y > x) {
return Fixed::half_pi() - detail::atan_sanitized(x / y);
}
return detail::atan_sanitized(y / x);
}
} // namespace detail
template <typename B, typename I, unsigned int F>
auto atan(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (x < Fixed(0))
{
return -atan(-x);
}
if (x > Fixed(1))
{
return Fixed::half_pi() - detail::atan_sanitized(Fixed(1) / x);
}
return detail::atan_sanitized(x);
}
template <typename B, typename I, unsigned int F>
auto asin(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x >= Fixed(-1) && x <= Fixed(+1));
const auto yy = Fixed(1) - x * x;
if (yy == Fixed(0))
{
return copysign(Fixed::half_pi(), x);
}
return detail::atan_div(x, sqrt(yy));
}
template <typename B, typename I, unsigned int F>
auto acos(fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
assert(x >= Fixed(-1) && x <= Fixed(+1));
if (x == Fixed(-1))
{
return Fixed::pi();
}
const auto yy = Fixed(1) - x * x;
return Fixed(2)*detail::atan_div(sqrt(yy), Fixed(1) + x);
}
template <typename B, typename I, unsigned int F>
auto atan2(fixed<B, I, F> y, fixed<B, I, F> x) noexcept -> fixed<B, I, F>
{
using Fixed = fixed<B, I, F>;
if (x == Fixed(0))
{
assert(y != Fixed(0));
return (y > Fixed(0)) ? Fixed::half_pi() : -Fixed::half_pi();
}
auto ret = detail::atan_div(y, x);
if (x < Fixed(0))
{
return (y >= Fixed(0)) ? ret + Fixed::pi() : ret - Fixed::pi();
}
return ret;
}
} // namespace fpm
#endif

View File

@ -55,7 +55,9 @@ HEADERS += \
$$PWD/vmodifierkey.h \
$$PWD/typedef.h \
$$PWD/backport/qscopeguard.h \
$$PWD/dialogs/dialogselectlanguage.h
$$PWD/dialogs/dialogselectlanguage.h \
$$PWD/fpm/fixed.hpp \
$$PWD/fpm/math.hpp
contains(DEFINES, APPIMAGE) {
SOURCES += \

View File

@ -53,8 +53,8 @@ DialogSinglePoint::DialogSinglePoint(const VContainer *data, quint32 toolId, QWi
ui->lineEditName->setClearButtonEnabled(true);
ui->doubleSpinBoxX->setRange(0, VAbstractValApplication::VApp()->fromPixel(SceneSize));
ui->doubleSpinBoxY->setRange(0, VAbstractValApplication::VApp()->fromPixel(SceneSize));
ui->doubleSpinBoxX->setRange(0, VAbstractValApplication::VApp()->fromPixel(maxSceneSize));
ui->doubleSpinBoxY->setRange(0, VAbstractValApplication::VApp()->fromPixel(maxSceneSize));
InitOkCancel(ui);
connect(ui->lineEditName, &QLineEdit::textChanged, this, [this]()

View File

@ -51,21 +51,18 @@
#include <QOpenGLWidget>
#include <QMimeData>
#include <QMimeDatabase>
#include <QtMath>
#include "../vmisc/def.h"
#include "../vmisc/vmath.h"
#include "vmaingraphicsscene.h"
#include "vsimplecurve.h"
#include "vcontrolpointspline.h"
#include "../vmisc/vabstractapplication.h"
#include "../vmisc/vcommonsettings.h"
#include "../vmisc/literals.h"
#include "vabstractmainwindow.h"
#include "global.h"
#include "../ifc/xml/utils.h"
const qreal maxSceneSize = ((20.0 * 1000.0) / 25.4) * PrintDPI; // 20 meters in pixels
namespace
{
qreal ScrollingSteps(QWheelEvent* wheel_event)