diff --git a/Core/FormFactors/FormFactorPyramid.cpp b/Core/FormFactors/FormFactorPyramid.cpp
index c1e334d357a01e917370bd6ae89e49e0f2c67637..83d1a2bc328b1a28e4d06ab90ab9b7dd7203f152 100644
--- a/Core/FormFactors/FormFactorPyramid.cpp
+++ b/Core/FormFactors/FormFactorPyramid.cpp
@@ -7,7 +7,7 @@
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
-//! @copyright Forschungszentrum Jülich GmbH 2015
+//! @copyright Forschungszentrum Jülich GmbH 2015-16
//! @authors Scientific Computing Group at MLZ Garching
//! @authors C. Durniak, M. Ganeva, G. Pospelov, W. Van Herck, J. Wuttke
//
@@ -15,14 +15,11 @@
#include "FormFactorPyramid.h"
#include "BornAgainNamespace.h"
-#include "MathFunctions.h"
-using namespace BornAgain;
-
-FormFactorPyramid::FormFactorPyramid(
- double length, double height, double alpha)
+FormFactorPyramid::FormFactorPyramid(double length, double height, double alpha)
+ : FormFactorPolyhedron( polyhedral_faces( length, height, alpha ), 0. )
{
- setName(FFPyramidType);
+ setName(BornAgain::FFPyramidType);
m_length = length;
m_height = height;
m_alpha = alpha;
@@ -30,34 +27,46 @@ FormFactorPyramid::FormFactorPyramid(
init_parameters();
}
+FormFactorPyramid::~FormFactorPyramid() {}
+
+std::vector<PolyhedralFace> FormFactorPyramid::polyhedral_faces(
+ double length, double height, double alpha)
+{
+ double a = length/2;
+ double b = length/2 - height/std::tan(alpha);
+
+ kvector_t V[8] = {
+ // base:
+ { -a, -a, 0. },
+ { a, -a, 0. },
+ { a, a, 0. },
+ { -a, a, 0. },
+ // top:
+ { -b, -b, height },
+ { b, -b, height },
+ { b, b, height },
+ { -b, b, height } };
+ std::vector<PolyhedralFace> faces;
+ faces.push_back( PolyhedralFace( { V[3], V[2], V[1], V[0] }, true ) );
+ faces.push_back( PolyhedralFace( { V[0], V[1], V[5], V[4] } ) );
+ faces.push_back( PolyhedralFace( { V[1], V[2], V[6], V[5] } ) );
+ faces.push_back( PolyhedralFace( { V[2], V[3], V[7], V[6] } ) );
+ faces.push_back( PolyhedralFace( { V[3], V[0], V[4], V[7] } ) );
+ faces.push_back( PolyhedralFace( { V[4], V[5], V[6], V[7] }, true ) );
+
+ return faces;
+}
+
bool FormFactorPyramid::check_initialization() const
{
bool result(true);
- if(m_alpha<0.0 || m_alpha>Units::PID2) {
- std::ostringstream ostr;
- ostr << "FormFactorPyramid() -> Error in class initialization with parameters ";
- ostr << " length:" << m_length;
- ostr << " height:" << m_height;
- ostr << " alpha:" << m_alpha << "\n\n";
- ostr << "Check for '0 <= alpha <= pi/2' failed.";
- throw Exceptions::ClassInitializationException(ostr.str());
- }
- if(m_length<0.0 || m_height<0.0) {
- std::ostringstream ostr;
- ostr << "FormFactorPyramid() -> Error in class initialization with parameters ";
- ostr << " length:" << m_length;
- ostr << " height:" << m_height;
- ostr << " alpha:" << m_alpha << "\n\n";
- ostr << "Check for '0 <= length and 0 <= height' failed.";
- throw Exceptions::ClassInitializationException(ostr.str());
- }
- if(2.*m_height > m_length*std::tan(m_alpha)) {
+ if(m_height > m_length*std::tan(m_alpha)) {
std::ostringstream ostr;
- ostr << "FormFactorPyramid() -> Error in class initialization with parameters ";
+ ostr << "FormFactorPyramid() -> Error in class initialization with parameters";
ostr << " length:" << m_length;
ostr << " height:" << m_height;
- ostr << " alpha:" << m_alpha << "\n\n";
- ostr << "Check for '2.*height <= length*tan(alpha)' failed.";
+ ostr << " alpha[rad]:" << m_alpha << "\n\n";
+ ostr << "Check for 'height <= length*tan(alpha)' failed.";
throw Exceptions::ClassInitializationException(ostr.str());
}
return result;
@@ -66,102 +75,17 @@ bool FormFactorPyramid::check_initialization() const
void FormFactorPyramid::init_parameters()
{
clearParameterPool();
- registerParameter(Length, &m_length, AttLimits::n_positive());
- registerParameter(Height, &m_height, AttLimits::n_positive());
- registerParameter(Alpha, &m_alpha, AttLimits::n_positive());
+ registerParameter(BornAgain::Length, &m_length, AttLimits::n_positive());
+ registerParameter(BornAgain::Height, &m_height, AttLimits::n_positive());
+ registerParameter(BornAgain::Alpha, &m_alpha, AttLimits::n_positive());
}
FormFactorPyramid* FormFactorPyramid::clone() const
{
- return new FormFactorPyramid(m_length, m_height, m_alpha);
+ return new FormFactorPyramid(m_length, m_height, m_alpha);
}
void FormFactorPyramid::accept(ISampleVisitor *visitor) const
{
visitor->visit(this);
}
-
-double FormFactorPyramid::getRadius() const
-{
- return m_length/2.0;
-}
-
-complex_t FormFactorPyramid::evaluate_for_q(const cvector_t& q) const
-{
-
- double H = m_height;
- double R = m_length/2.;
- double tga = std::tan(m_alpha);
- if (m_height == 0.0) return 0.0;
- if (m_alpha == 0.0) return 0.0;
- if (m_length == 0.0) return 0.0;
-
- complex_t qx = q.x();
- complex_t qy = q.y();
- complex_t qz = q.z();
-
- const complex_t im(0, 1);
- complex_t full = fullPyramidPrimitive(qx/tga, qy/tga, qz, -R*tga);
- complex_t top = fullPyramidPrimitive(qx/tga, qy/tga, qz, H-R*tga);
- return std::exp(im*qz*R*tga) * (full-top) / (tga*tga);
-}
-
-complex_t FormFactorPyramid::fullPyramidPrimitive(complex_t a, complex_t b, complex_t c,
- double z) const
-{
- const complex_t im(0, 1);
- if (std::norm(a * z) > Numeric::double_epsilon && std::norm(b * z) > Numeric::double_epsilon) {
- if (std::abs((a - b) * (a - b) - c * c) * z * z > Numeric::double_epsilon
- && std::abs((a + b) * (a + b) - c * c) * z * z > Numeric::double_epsilon) {
- complex_t phase = std::exp(im * c * z);
- complex_t numerator
- = std::sin(a * z) * (b * (a * a - b * b + c * c) * std::cos(b * z)
- + im * c * (a * a + b * b - c * c) * std::sin(b * z))
- + a * std::cos(a * z) * ((-a * a + b * b + c * c) * std::sin(b * z)
- + 2.0 * im * b * c * std::cos(b * z));
- complex_t denominator = a * b * (a - b - c) * (a + b - c) * (a - b + c) * (a + b + c);
- return -4.0 * phase * numerator / denominator;
- } else {
- return 2.0 * (g(a - b, c, z) - g(a + b, c, z)) / (a * b);
- }
- } else if (std::norm(a * z) <= Numeric::double_epsilon
- && std::norm(b * z) <= Numeric::double_epsilon) {
- if (std::norm(c * z) <= Numeric::double_epsilon) {
- return -4.0 * std::pow(z, 3) / 3.0;
- } else
- return 4.0 * im
- * (2.0 + std::exp(im * c * z) * (c * c * z * z + 2.0 * im * c * z - 2.0))
- / std::pow(c, 3);
- } else {
- complex_t abmax;
- if (std::norm(b * z) <= Numeric::double_epsilon) {
- abmax = a;
- } else {
- abmax = b;
- }
- return 2.0 * (h(c - abmax, z) - h(c + abmax, z)) / abmax;
- }
-}
-
-complex_t FormFactorPyramid::g(complex_t x, complex_t c, double z) const
-{
- const complex_t im(0, 1);
- if (std::abs((x*x-c*c)*z*z) > Numeric::double_epsilon) {
- return (im*c - std::exp(im*c*z)*(x*std::sin(x*z) + im*c*std::cos(x*z))) / (x*x-c*c);
- } else {
- if (std::norm(c*z) > Numeric::double_epsilon) {
- return im * (std::exp(2.0*im*c*z) + 2.0*im*c*z - 1.0) / (4.0*c);
- } else {
- return -z;
- }
- }
-}
-
-complex_t FormFactorPyramid::h(complex_t x, double z) const
-{
- const complex_t im(0, 1);
- if (std::norm(x*z) > Numeric::double_epsilon) {
- return (im - (im+x*z)*std::exp(im*x*z))/(x*x);
- }
- return -im*z*z/2.0;
-}
diff --git a/Core/FormFactors/FormFactorPyramid.h b/Core/FormFactors/FormFactorPyramid.h
index ce0c04ca67440025a0703b10509b788d9ca4d8cc..95458cf92f0ffb98527fca498d3d58e7f98df5db 100644
--- a/Core/FormFactors/FormFactorPyramid.h
+++ b/Core/FormFactors/FormFactorPyramid.h
@@ -7,7 +7,7 @@
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
-//! @copyright Forschungszentrum Jülich GmbH 2015
+//! @copyright Forschungszentrum Jülich GmbH 2015-16
//! @authors Scientific Computing Group at MLZ Garching
//! @authors C. Durniak, M. Ganeva, G. Pospelov, W. Van Herck, J. Wuttke
//
@@ -16,13 +16,12 @@
#ifndef FORMFACTORPYRAMID_H
#define FORMFACTORPYRAMID_H
-#include "IFormFactorBorn.h"
+#include "FormFactorPolyhedron.h"
//! @class FormFactorPyramid
//! @ingroup formfactors
-//! @brief The form factor of pyramid.
-
-class BA_CORE_API_ FormFactorPyramid : public IFormFactorBorn
+//! @brief The formfactor of a cone6.
+class BA_CORE_API_ FormFactorPyramid : public FormFactorPolyhedron
{
public:
//! @brief Pyramid constructor
@@ -30,48 +29,33 @@ public:
//! @param height of Pyramid
//! @param angle in radians between base and facet
FormFactorPyramid(double length, double height, double alpha);
+ virtual ~FormFactorPyramid();
- virtual ~FormFactorPyramid() {}
- virtual FormFactorPyramid *clone() const;
+ static std::vector<PolyhedralFace> polyhedral_faces(
+ double length, double height, double alpha);
- virtual void accept(ISampleVisitor *visitor) const;
+ virtual FormFactorPyramid* clone() const;
- virtual double getRadius() const;
+ virtual void accept(ISampleVisitor *visitor) const;
+ virtual double getRadius() const final;
double getHeight() const;
-
double getLength() const;
-
double getAlpha() const;
- virtual complex_t evaluate_for_q(const cvector_t& q) const;
-
protected:
virtual bool check_initialization() const;
virtual void init_parameters();
private:
- complex_t fullPyramidPrimitive(complex_t a, complex_t b, complex_t c, double z) const;
- complex_t g(complex_t x, complex_t c, double z) const; // helper function
- complex_t h(complex_t x, double z) const; // helper function
double m_length;
double m_height;
double m_alpha;
};
-inline double FormFactorPyramid::getHeight() const
-{
- return m_height;
-}
-
-inline double FormFactorPyramid::getLength() const
-{
- return m_length;
-}
-
-inline double FormFactorPyramid::getAlpha() const
-{
- return m_alpha;
-}
+inline double FormFactorPyramid::getHeight() const { return m_height; }
+inline double FormFactorPyramid::getLength() const { return m_length; }
+inline double FormFactorPyramid::getAlpha() const { return m_alpha; }
+inline double FormFactorPyramid::getRadius() const { return m_length/2.0; }
#endif // FORMFACTORPYRAMID_H
diff --git a/Doc/UserManual/FormFactors.tex b/Doc/UserManual/FormFactors.tex
index 75129c7b2909218444737ff23c63fbab1ffca0b5..81b44c78438718706f25c920014909e3a1f08847 100644
--- a/Doc/UserManual/FormFactors.tex
+++ b/Doc/UserManual/FormFactors.tex
@@ -998,7 +998,7 @@ for four different angles~$\omega$ of rotation around the $z$ axis.}
\label{fig:FFprism3Ex}
\end{figure}
-\paragraph{References}\strut\\
+\paragraph{References and History}\strut\\
Has been validated against the \E{Prism3} form factor of \IsGISAXS\
\cite[Eq.~2.29]{Laz08} \cite[Eq.~219]{ReLL09}.
Note the different parameterization $L= 2 R_{\rm{\Code{IsGISAXS}}}$.
@@ -1114,23 +1114,9 @@ They must fulfill
\end{displaymath}
\paragraph{Form factor, volume, horizontal section}\strut\\
-Notation:
-\begin{displaymath}
- a\coloneqq q_x/\tan\alpha,\quad
- b\coloneqq q_y/\tan\alpha,\quad
- c\coloneqq q_z.
-\end{displaymath}
-Results:
-\begin{equation*}
- \DS F= \frac{\exp(iq_zL\tan\alpha/2)}{\tan^2\alpha} \Big\{ I_p(H-L\tan\alpha/2) - I_p(-L\tan\alpha/2) \Big\}
-\end{equation*}
-Where $I_p(z)$ is an antiderivative for the final z--integration in the Fourier transform of the non--truncated pyramid shape function:
\begin{equation*}
-\begin{array}{@{}l@{}l@{}}
- \DS I_p(z) = & 4\exp(icz)\Big\{ \sin(az)\big[b(a^2-b^2+c^2)\cos(bz) + ic(a^2+b^2-c^2)\sin(bz)\big] \\
- & +a\cos(az)\big[(-a^2+b^2+c^2)\sin(bz) + 2ibc\cos(bz) \big] \Big\} \\
- & / \Big\{ab\big[(a+b)^2-c^2\big]\big[(a-b)^2-c^2\big]\Big\}
-\end{array}
+ F \text{~: computed algebraically,
+ using the generic polyhedron form factor~\cite{ba:ffp},}
\end{equation*}
\begin{equation*}
V = \dfrac{1}{6} L^3 \tan\alpha\left[ 1
@@ -1152,12 +1138,15 @@ computed with $L=10$~nm, $H=4.2$~nm and $\alpha=60^{\circ}$,
for four different angles~$\omega$ of rotation around the $z$ axis.}
\end{figure}
-\paragraph{References}\strut\\
+\paragraph{References and History}\strut\\
Corresponds to \E{Pyramid} form factor of \IsGISAXS\
\cite[Eq.~2.31]{Laz08} \cite[Eq.~221]{ReLL09},
-with different parametrization $L=2R_{\rm{\Code{IsGISXAXS}}}$,
-with correction of a sign error,
-and with a more compact form of $F(\q)$.
+except for different parametrization $L=2R_{\rm{\Code{IsGISXAXS}}}$
+and a corrected sign.
+
+Reimplemented in \BornAgain-1.6 using the generic form factor
+of a polygonal prism \cite{ba:ffp},
+to achieve numerical stability near the removable singularity at $q\to0$.
%===============================================================================