scatteringobservables2d.cpp

Go to the documentation of this file.

/**
 * @file scatteringobservables2d.cpp
 * @brief Implements the scatteringobservables2d physics/numerical component used by QPE experiments.
 * @ingroup qpe_quantum
 *
 * Source: @c core/quantum/scatteringobservables2d.cpp. This file is part of Quantum Phenomena Explorer (QPE),
 * maintained by Officina Turini. See LICENSE_IT.txt / LICENSE.txt for the
 * QPE freeware source-available terms and THIRD_PARTY_NOTICES.txt for
 * third-party components.
 */
 
#include "scatteringobservables2d.h"
 
#include "core/quantum/gaussianpacket2d.h"
 
#include <cmath>
#include <utility>
 
namespace qpe {
 
ScatteringRegions2D ScatteringObservables2D::measureX(const ComplexField2D &field,
                                                       double leftBoundary,
                                                       double rightBoundary) noexcept
{
    ScatteringRegions2D out;
    if (field.empty()) {
        return out;
    }
 
    if (leftBoundary > rightBoundary) {
        std::swap(leftBoundary, rightBoundary);
    }
 
    const double cellArea = static_cast<double>(field.dx()) *
                            static_cast<double>(field.dy());
 
    for (std::size_t y = 0; y < field.height(); ++y) {
        for (std::size_t x = 0; x < field.width(); ++x) {
            const auto &v = field(x, y);
            const double re = static_cast<double>(v.real());
            const double im = static_cast<double>(v.imag());
            const double p = (re * re + im * im) * cellArea;
            const double px = GaussianPacket2D::xCoordinate(field, x);
 
            out.total += p;
            if (px < leftBoundary) {
                out.left += p;
            } else if (px > rightBoundary) {
                out.right += p;
            } else {
                out.interaction += p;
            }
        }
    }
 
    if (out.total > 0.0 && std::isfinite(out.total)) {
        const double inv = 1.0 / out.total;
        out.leftFraction = out.left * inv;
        out.interactionFraction = out.interaction * inv;
        out.rightFraction = out.right * inv;
    }
 
    return out;
}
 
} // namespace qpe