ReParticle.cpp

//  ************************************************************************************************
//
//  BornAgain: simulate and fit reflection and scattering
//
//! @file      Resample/Particle/ReParticle.cpp
//! @brief     Implements class interface ReParticle.
//!
//! @homepage  http://www.bornagainproject.org
//! @license   GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors   Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
//  ************************************************************************************************

#include "Resample/Particle/ReParticle.h"
#include "Base/Types/Span.h"
#include "Base/Util/Assert.h"
#include "Base/Vector/WavevectorInfo.h" // debug
#include "Sample/Material/Material.h"
#include "Sample/Material/MaterialFactoryFuncs.h"
#include "Sample/Particle/IFormFactor.h"
#include "Sample/Scattering/Rotations.h"

ReParticle::ReParticle(const std::optional<size_t>& i_layer, const IFormFactor* ff,
                       const Material* material, const Material* ambient_material,
                       const R3* position, const RotMatrix* rotMatrix)
    : IReParticle(i_layer)
    , m_ff(ff)
    , m_material(material)
    , m_ambient_material(ambient_material)
    , m_position(position)
    , m_rotMatrix(rotMatrix)
{
}

ReParticle::ReParticle(const IFormFactor* ff, const R3* position, const RotMatrix* rot)
    : ReParticle({}, ff, nullptr, nullptr, position, rot)
{
}

ReParticle::~ReParticle() = default;

ReParticle* ReParticle::clone() const
{
    return new ReParticle(i_layer(), m_ff->clone(),
                          m_material ? new Material(*m_material) : nullptr,
                          m_ambient_material ? new Material(*m_ambient_material) : nullptr,
                          m_position ? new R3(*m_position) : nullptr,
                          m_rotMatrix ? new RotMatrix(*m_rotMatrix) : nullptr);
}

void ReParticle::setMaterial(const Material& material)
{
    m_material = std::make_unique<Material>(material);
}

void ReParticle::setAmbientMaterial(const Material& ambient_material)
{
    m_ambient_material = std::make_unique<Material>(ambient_material);
}

double ReParticle::volume() const
{
    return m_ff->volume();
}

double ReParticle::radialExtension() const
{
    return m_ff->radialExtension();
}

const IFormFactor* ReParticle::iformfactor() const
{
    return m_ff.get();
}

complex_t ReParticle::theFF(const WavevectorInfo& wavevectors) const
{
    WavevectorInfo wavevectors2 =
        m_rotMatrix ? wavevectors.transformed(m_rotMatrix->Inverse()) : wavevectors;
    complex_t result = m_ff->theFF(wavevectors2);
    if (m_material && m_ambient_material)
        result = (m_material->scalarSubtrSLD(wavevectors2)
                  - m_ambient_material->scalarSubtrSLD(wavevectors2))
                 * result;
    if (m_position)
        result *= exp_I(m_position->dot(wavevectors.getQ()));
    return result;
}

SpinMatrix ReParticle::thePolFF(const WavevectorInfo& wavevectors) const
{
    WavevectorInfo wavevectors2 =
        m_rotMatrix ? wavevectors.transformed(m_rotMatrix->Inverse()) : wavevectors;
    SpinMatrix result = m_ff->thePolFF(wavevectors2);
    if (m_material && m_ambient_material) {
        // the conjugated linear part of time reversal operator T
        // (T=UK with K complex conjugate operator and U is linear)
        SpinMatrix time_reverse_conj(0, 1, -1, 0);
        // the interaction and time reversal taken together:
        SpinMatrix V_eff = time_reverse_conj
                           * (m_material->polarizedSubtrSLD(wavevectors2)
                              - m_ambient_material->polarizedSubtrSLD(wavevectors2));
        result *= V_eff;
    }
    if (m_position)
        return result *= exp_I(m_position->dot(wavevectors.getQ()));
    return result;
}

Span ReParticle::zSpan() const
{
    RotMatrix transform = m_rotMatrix ? *m_rotMatrix : RotMatrix();
    std::unique_ptr<const IRotation> total_rotation(IRotation::createRotation(transform));
    Span span = m_ff->spanZ(total_rotation.get());
    if (m_position)
        return span + m_position->z();
    return span;
}