CoordSystem2D.cpp

Go to the documentation of this file.

//  ************************************************************************************************
//
//  BornAgain: simulate and fit reflection and scattering
//
//! @file      Device/Coord/CoordSystem2D.cpp
//! @brief     Implements ICoordSystem classes.
//!
//! @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 "Device/Coord/CoordSystem2D.h"
#include "Base/Axis/FixedBinAxis.h"
#include "Base/Const/Units.h"
#include "Base/Math/Constants.h"
#include "Base/Pixel/RectangularPixel.h"
#include "Base/Util/Assert.h"
#include "Base/Vector/Direction.h"
#include "Device/Coord/AxisNames.h"
#include <algorithm>
#include <cmath>
 
namespace {
 
double getQ(double wavelength, double angle)
{
    return 4.0 * M_PI * std::sin(angle) / wavelength;
}
 
double axisAngle(size_t i_axis, R3 k_f)
{
    if (i_axis == 0)
        return R3Util::phi(k_f);
    if (i_axis == 1)
        return M_PI_2 - R3Util::theta(k_f);
    ASSERT(0);
}
 
} // namespace
 
//  ************************************************************************************************
//  class CoordSystem2D
//  ************************************************************************************************
 
CoordSystem2D::CoordSystem2D(const Direction& direction, double wavelength)
    : CoordSystem2D({}, direction, wavelength)
{
}
 
CoordSystem2D::CoordSystem2D(const OwningVector<IAxis>& axes, const Direction& direction,
                             double wavelength)
    : m_axes(axes)
    , m_wavelength(wavelength)
    , m_alpha_i(-direction.alpha())
    , m_phi_i(direction.phi())
{
}
 
CoordSystem2D::CoordSystem2D(const CoordSystem2D& other)
    : m_axes(other.m_axes)
    , m_wavelength(other.m_wavelength)
    , m_alpha_i(other.m_alpha_i)
    , m_phi_i(other.m_phi_i)
{
}
 
void CoordSystem2D::addAxisData(std::string name, double min, double max, size_t nbins)
{
    m_axes.emplace_back(new FixedBinAxis(name, min, max, nbins));
}
 
double CoordSystem2D::calculateMin(size_t i_axis, Coords units) const
{
    ASSERT(i_axis < rank());
    units = substituteDefaultUnits(units);
    if (units == Coords::NBINS)
        return 0.0;
    return calculateValue(i_axis, units, m_axes[i_axis]->min());
}
 
double CoordSystem2D::calculateMax(size_t i_axis, Coords units) const
{
    ASSERT(i_axis < rank());
    units = substituteDefaultUnits(units);
    if (units == Coords::NBINS)
        return m_axes[i_axis]->size();
    return calculateValue(i_axis, units, m_axes[i_axis]->max());
}
 
size_t CoordSystem2D::axisSize(size_t i_axis) const
{
    ASSERT(i_axis < rank());
    return m_axes[i_axis]->size();
}
 
std::vector<Coords> CoordSystem2D::availableUnits() const
{
    return {Coords::NBINS, Coords::RADIANS, Coords::DEGREES};
}
 
IAxis* CoordSystem2D::createConvertedAxis(size_t i_axis, Coords units) const
{
    const double min = calculateMin(i_axis, units);
    const double max = calculateMax(i_axis, units);
    const auto& axis_name = axisName(i_axis, units);
    const auto axis_size = axisSize(i_axis);
    return new FixedBinAxis(axis_name, axis_size, min, max);
}
 
//  ************************************************************************************************
//  class SphericalCoords
//  ************************************************************************************************
 
SphericalCoords::SphericalCoords(const OwningVector<IAxis>& axes, const Direction& direction,
                                 double wavelength)
    : CoordSystem2D(axes, direction, wavelength)
{
    ASSERT(axes.size() == 2);
}
 
SphericalCoords::SphericalCoords(const SphericalCoords& other)
    : CoordSystem2D(other)
{
}
 
SphericalCoords::~SphericalCoords() = default;
 
SphericalCoords* SphericalCoords::clone() const
{
    return new SphericalCoords(*this);
}
 
std::vector<Coords> SphericalCoords::availableUnits() const
{
    auto result = CoordSystem2D::availableUnits();
    result.push_back(Coords::QSPACE);
    return result;
}
 
double SphericalCoords::calculateValue(size_t i_axis, Coords units, double value) const
{
    switch (units) {
    case Coords::RADIANS:
        return value;
    case Coords::DEGREES:
        return Units::rad2deg(value);
    case Coords::QSPACE: {
        const R3 k_i = vecOfLambdaAlphaPhi(m_wavelength, m_alpha_i, m_phi_i);
        if (i_axis == 0) {
            const R3 k_f = vecOfLambdaAlphaPhi(m_wavelength, 0.0, value);
            return (k_i - k_f).y();
        }
        if (i_axis == 1) {
            const R3 k_f = vecOfLambdaAlphaPhi(m_wavelength, value, 0.0);
            return (k_f - k_i).z();
        }
        ASSERT(0);
    }
    case Coords::QXQY: {
        const R3 k_i = vecOfLambdaAlphaPhi(m_wavelength, m_alpha_i, m_phi_i);
        if (i_axis == 0) {
            const R3 k_f = vecOfLambdaAlphaPhi(m_wavelength, 0.0, value);
            return (k_i - k_f).y();
        }
        if (i_axis == 1) {
            const R3 k_f = vecOfLambdaAlphaPhi(m_wavelength, value, 0.0);
            return (k_f - k_i).x();
        }
        ASSERT(0);
    }
    default:
        ASSERT(0);
    }
}
 
std::vector<std::map<Coords, std::string>> SphericalCoords::createNameMaps() const
{
    return {DataUtils::AxisNames::sphericalAxis0, DataUtils::AxisNames::sphericalAxis1};
}
 
//  ************************************************************************************************
//  class ImageCoords
//  ************************************************************************************************
 
ImageCoords::ImageCoords(const OwningVector<IAxis>& axes, RectangularPixel* regionOfInterestPixel,
                         const Direction& direction, double wavelength)
    : CoordSystem2D(axes, direction, wavelength)
    , m_detector_pixel(regionOfInterestPixel)
{
    ASSERT(axes.size() == 2);
}
 
ImageCoords::ImageCoords(const ImageCoords& other)
    : CoordSystem2D(other)
    , m_detector_pixel(other.m_detector_pixel->clone())
{
}
 
ImageCoords::~ImageCoords() = default;
 
ImageCoords* ImageCoords::clone() const
{
    return new ImageCoords(*this);
}
 
std::vector<Coords> ImageCoords::availableUnits() const
{
    auto result = CoordSystem2D::availableUnits();
    result.push_back(Coords::QSPACE);
    result.push_back(Coords::MM);
    return result;
}
 
double ImageCoords::calculateValue(size_t i_axis, Coords units, double value) const
{
    ASSERT(m_wavelength > 0);
    if (units == Coords::MM)
        return value;
    const auto k00 = m_detector_pixel->getPosition(0.0, 0.0);
    const auto k01 = m_detector_pixel->getPosition(0.0, 1.0);
    const auto k10 = m_detector_pixel->getPosition(1.0, 0.0);
    const auto& max_pos = i_axis == 0 ? k10 : k01; // position of max along given axis
    const double shift = value - m_axes[i_axis]->min();
    const R3 out_dir = k00 + shift * (max_pos - k00).unit();
    const R3 k_f = out_dir.unit() * M_TWOPI / m_wavelength;
 
    switch (units) {
    case Coords::RADIANS:
        return axisAngle(i_axis, k_f);
    case Coords::DEGREES:
        return Units::rad2deg(axisAngle(i_axis, k_f));
    case Coords::QSPACE: {
        const auto k_i = vecOfLambdaAlphaPhi(m_wavelength, m_alpha_i, m_phi_i);
        if (i_axis == 0)
            return (k_i - k_f).y();
        if (i_axis == 1)
            return (k_f - k_i).z();
    } break;
    case Coords::QXQY: {
        const auto k_i = vecOfLambdaAlphaPhi(m_wavelength, m_alpha_i, m_phi_i);
        if (i_axis == 0)
            return (k_i - k_f).y();
        if (i_axis == 1)
            return (k_f - k_i).x();
    } break;
    default:
        break;
    }
    ASSERT(0);
}
 
std::vector<std::map<Coords, std::string>> ImageCoords::createNameMaps() const
{
    return {DataUtils::AxisNames::rectangularAxis0, DataUtils::AxisNames::rectangularAxis1};
}
 
//  ************************************************************************************************
//  class OffspecCoordinates
//  ************************************************************************************************
 
OffspecCoordinates::OffspecCoordinates(const OwningVector<IAxis>& axes, const Direction& direction)
    : CoordSystem2D(axes, direction)
{
}
 
OffspecCoordinates::OffspecCoordinates(const OffspecCoordinates& other)
    : CoordSystem2D(other)
{
}
 
OffspecCoordinates* OffspecCoordinates::clone() const
{
    return new OffspecCoordinates(*this);
}
 
double OffspecCoordinates::calculateValue(size_t, Coords units, double value) const
{
    switch (units) {
    case Coords::RADIANS:
        return value;
    case Coords::DEGREES:
        return Units::rad2deg(value);
    default:
        ASSERT(0);
    }
}
 
std::vector<std::map<Coords, std::string>> OffspecCoordinates::createNameMaps() const
{
    return {DataUtils::AxisNames::offspecAxis0, DataUtils::AxisNames::offspecAxis1};
}
 
//  ************************************************************************************************
//  class DepthProbeCoordinates
//  ************************************************************************************************
 
const std::string z_axis_name = "Position [nm]";
 
DepthProbeCoordinates::DepthProbeCoordinates(const OwningVector<IAxis>& axes,
                                             const Direction& direction, double wavelength)
    : CoordSystem2D(axes, direction, wavelength)
{
}
 
DepthProbeCoordinates::DepthProbeCoordinates(const DepthProbeCoordinates& other)
    : CoordSystem2D(other)
{
}
 
DepthProbeCoordinates::~DepthProbeCoordinates() = default;
 
DepthProbeCoordinates* DepthProbeCoordinates::clone() const
{
    return new DepthProbeCoordinates(*this);
}
 
std::vector<Coords> DepthProbeCoordinates::availableUnits() const
{
    auto result = CoordSystem2D::availableUnits();
    result.push_back(Coords::QSPACE);
    return result;
}
 
double DepthProbeCoordinates::calculateValue(size_t i_axis, Coords units, double value) const
{
    const auto& available_units = availableUnits();
    if (std::find(available_units.begin(), available_units.end(), units) == available_units.cend())
        ASSERT(0);
 
    if (i_axis == 1)
        return value; // unit conversions are not applied to sample position axis
    switch (units) {
    case Coords::DEGREES:
        return Units::rad2deg(value);
    case Coords::QSPACE:
        return getQ(m_wavelength, value);
    default:
        return value;
    }
}
 
std::vector<std::map<Coords, std::string>> DepthProbeCoordinates::createNameMaps() const
{
    return {DataUtils::AxisNames::specAxis, DataUtils::AxisNames::sampleDepthAxis};
}