diff --git a/broker/OpenAPIConvert.cpp b/broker/OpenAPIConvert.cpp index ff273c16..cc296d52 100644 --- a/broker/OpenAPIConvert.cpp +++ b/broker/OpenAPIConvert.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../common/JFJochMath.h" #include "OpenAPIConvert.h" // From https://en.cppreference.com/w/cpp/string/byte/tolower @@ -460,8 +461,8 @@ ROIDefinition Convert(const org::openapitools::server::model::Roi_definitions& i output.circles.emplace_back(ROICircle(i.getName(), i.getCenterXPxl(), i.getCenterYPxl(), i.getRadiusPxl())); for (const auto &i: input.getAzim().getRois()) output.azimuthal.emplace_back(ROIAzimuthal(i.getName(), - (i.getQMaxRecipA() == 0.0) ? 0.0 : 2.0f * M_PI / i.getQMaxRecipA(), - (i.getQMinRecipA() == 0.0) ? 0.0 : 2.0f * M_PI / i.getQMinRecipA())); + (i.getQMaxRecipA() == 0.0) ? 0.0 : 2.0f * PI / i.getQMaxRecipA(), + (i.getQMinRecipA() == 0.0) ? 0.0 : 2.0f * PI / i.getQMinRecipA())); return output; } diff --git a/common/AzimuthalIntegrationMapping.cpp b/common/AzimuthalIntegrationMapping.cpp index 4c558f02..b72337ef 100644 --- a/common/AzimuthalIntegrationMapping.cpp +++ b/common/AzimuthalIntegrationMapping.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include #include #include @@ -145,8 +146,8 @@ void AzimuthalIntegrationMapping::SetupPixel(const DiffractionGeometry &geom, corrections[pxl] = corr; if (d > 0) { - float q = 2.0f * static_cast(M_PI) / d; - pixel_to_bin[pxl] = settings.GetBin(q, phi_rad * 180.0 / M_PI); + float q = 2.0f * static_cast(PI) / d; + pixel_to_bin[pxl] = settings.GetBin(q, phi_rad * 180.0 / PI); } } @@ -187,9 +188,9 @@ void AzimuthalIntegrationMapping::UpdateMaxBinNumber() { for (int j = 0; j < settings.GetAzimuthalBinCount(); j++) { for (int i = 0; i < settings.GetQBinCount(); i++) { bin_to_q[j * settings.GetQBinCount() + i] = static_cast(settings.GetQSpacing_recipA() * (i + 0.5) + settings.GetLowQ_recipA()); - bin_to_d[j * settings.GetQBinCount() + i] = 2.0f * static_cast(M_PI) / bin_to_q[j * settings.GetQBinCount() + i]; - bin_to_2theta[j * settings.GetQBinCount() + i] = 2.0f * asinf(bin_to_q[i] * wavelength / (4.0f * static_cast(M_PI))) * 180.0f / - static_cast(M_PI); + bin_to_d[j * settings.GetQBinCount() + i] = 2.0f * static_cast(PI) / bin_to_q[j * settings.GetQBinCount() + i]; + bin_to_2theta[j * settings.GetQBinCount() + i] = 2.0f * asinf(bin_to_q[i] * wavelength / (4.0f * static_cast(PI))) * 180.0f / + static_cast(PI); bin_to_phi[j * settings.GetQBinCount() + i] = static_cast(j) * 360.0f / static_cast(settings.GetAzimuthalBinCount()); } } diff --git a/common/AzimuthalIntegrationSettings.h b/common/AzimuthalIntegrationSettings.h index eafb97f5..7f7ee8e2 100644 --- a/common/AzimuthalIntegrationSettings.h +++ b/common/AzimuthalIntegrationSettings.h @@ -3,6 +3,7 @@ #pragma once +#include "JFJochMath.h" #include #include @@ -11,8 +12,8 @@ class AzimuthalIntegrationSettings { bool polarization_correction = true; float high_q_recipA = 5.0; float low_q_recipA = 0.1; - float bkg_estimate_high_q_recipA = 2.0f * M_PI / 3.0; - float bkg_estimate_low_q_recipA = 2.0f * M_PI / 5.0; + float bkg_estimate_high_q_recipA = 2.0f * PI / 3.0; + float bkg_estimate_low_q_recipA = 2.0f * PI / 5.0; float q_spacing = 0.05; int32_t azim_bins = 1; diff --git a/common/Coord.cpp b/common/Coord.cpp index 54809077..8595c40a 100644 --- a/common/Coord.cpp +++ b/common/Coord.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include #include "Coord.h" #include "JFJochException.h" @@ -143,7 +144,7 @@ float angle_deg(const Coord &c1, const Coord &c2) { else if (cos_ang < -1.0f && cos_ang > -1.0f - epsilon) cos_ang = -1.0f; - return acosf(cos_ang) * (180.0f / static_cast(M_PI)); + return acosf(cos_ang) * (180.0f / static_cast(PI)); } void Coord::swap(Coord &other) noexcept { diff --git a/common/CrystalLattice.cpp b/common/CrystalLattice.cpp index 17dbca65..ac3583c5 100644 --- a/common/CrystalLattice.cpp +++ b/common/CrystalLattice.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include #include "CrystalLattice.h" #include "JFJochException.h" @@ -9,7 +10,7 @@ #include "gemmi/unitcell.hpp" #include "gemmi/cellred.hpp" -#define DEG_TO_RAD static_cast(M_PI/180.0) +#define DEG_TO_RAD static_cast(PI/180.0) CrystalLattice::CrystalLattice(const UnitCell &cell) { vec[0] = {cell.a, 0, 0}; diff --git a/common/DiffractionGeometry.cpp b/common/DiffractionGeometry.cpp index 51163a59..1234ffa8 100644 --- a/common/DiffractionGeometry.cpp +++ b/common/DiffractionGeometry.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include #include "DiffractionGeometry.h" @@ -48,7 +49,7 @@ float DiffractionGeometry::Phi_rad(float x, float y) const { auto lab = LabCoord(x, y); auto v = atan2f(lab.y, lab.x); if (v < 0) - v += 2.0f * M_PI; + v += 2.0f * PI; return v; } @@ -58,7 +59,7 @@ float DiffractionGeometry::PxlToRes(float x, float y) const { } float DiffractionGeometry::PxlToQ(float x, float y) const { - return 2.0f * M_PI / PxlToRes(x,y); + return 2.0f * PI / PxlToRes(x,y); } float DiffractionGeometry::PxlToRes(float dist_pxl) const { @@ -87,7 +88,7 @@ float DiffractionGeometry::DistFromEwaldSphere(const Coord &recip) const { } float DiffractionGeometry::CalcAzIntSolidAngleCorr(float q) const { - float sin_theta = q * wavelength_A / (4 * static_cast(M_PI)); + float sin_theta = q * wavelength_A / (4 * static_cast(PI)); float cos_2theta = 1.0f - 2.0f * sin_theta * sin_theta; // cos(2*alpha) = 1 - 2 * sin(alpha)^2 float cos_2theta_3 = cos_2theta * cos_2theta * cos_2theta; return cos_2theta_3; diff --git a/common/GoniometerAxis.cpp b/common/GoniometerAxis.cpp index bfe9dedf..bd56d392 100644 --- a/common/GoniometerAxis.cpp +++ b/common/GoniometerAxis.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include #include @@ -138,6 +139,6 @@ std::vector GoniometerAxis::GetAngleContainerEnd(int64_t max_image_numbe } RotMatrix GoniometerAxis::GetTransformationAngle(float angle_deg) const { - auto angle_rad = angle_deg / 180.0f * static_cast(M_PI); + auto angle_rad = angle_deg / 180.0f * static_cast(PI); return {angle_rad, axis}; } \ No newline at end of file diff --git a/common/JFJochMath.h b/common/JFJochMath.h new file mode 100644 index 00000000..d6d9d579 --- /dev/null +++ b/common/JFJochMath.h @@ -0,0 +1,10 @@ +// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute +// SPDX-License-Identifier: GPL-3.0-only + +#pragma once + +// pi as a plain constexpr so it is usable everywhere: +// - host code on MSVC, which does not define M_PI without _USE_MATH_DEFINES; +// - CUDA device code, which is compiled as C++17 here, so PI +// (C++20) is not available. +constexpr double PI = 3.14159265358979323846; diff --git a/common/ROIAzimuthal.cpp b/common/ROIAzimuthal.cpp index 60dbc6a9..d9321b68 100644 --- a/common/ROIAzimuthal.cpp +++ b/common/ROIAzimuthal.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "JFJochMath.h" #include "ROIAzimuthal.h" #include "JFJochException.h" @@ -31,11 +32,11 @@ bool ROIAzimuthal::CheckROI(int64_t x, int64_t y, float resolution) const { } float ROIAzimuthal::GetQMax_recipA() const { - return 2.0f * M_PI / d_min_A; + return 2.0f * PI / d_min_A; } float ROIAzimuthal::GetQMin_recipA() const { - return 2.0f * M_PI / d_max_A; + return 2.0f * PI / d_max_A; } ROIConfig ROIAzimuthal::ExportMetadata() const { diff --git a/image_analysis/bragg_prediction/BraggPrediction.cpp b/image_analysis/bragg_prediction/BraggPrediction.cpp index 4c134e38..70ebea71 100644 --- a/image_analysis/bragg_prediction/BraggPrediction.cpp +++ b/image_analysis/bragg_prediction/BraggPrediction.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "BraggPrediction.h" #include "../bragg_integration/SystematicAbsence.h" @@ -46,7 +47,7 @@ int BraggPrediction::Calc(const DiffractionExperiment &experiment, const Crystal const float epsilon = 1e-5f; const float s0_sq = S0 * S0; - const float rad_to_deg = 180.0f / static_cast(M_PI); + const float rad_to_deg = 180.0f / static_cast(PI); int i = 0; diff --git a/image_analysis/bragg_prediction/BraggPredictionGPU.cu b/image_analysis/bragg_prediction/BraggPredictionGPU.cu index cf8f8cc6..04e4da64 100644 --- a/image_analysis/bragg_prediction/BraggPredictionGPU.cu +++ b/image_analysis/bragg_prediction/BraggPredictionGPU.cu @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "BraggPredictionGPU.h" #ifdef JFJOCH_USE_CUDA @@ -16,7 +17,7 @@ namespace { __device__ inline float angle_from_ewald_sphere_deg(const Coord &S0, float recip_x, float recip_y, float recip_z, float recip_sq) { const float epsilon = 1e-5f; - const float rad_to_deg = 180.0f / static_cast(M_PI); + const float rad_to_deg = 180.0f / static_cast(PI); const float s0_sq = S0.x * S0.x + S0.y * S0.y + S0.z * S0.z; const float s0_p0 = S0.x * recip_x + S0.y * recip_y + S0.z * recip_z; diff --git a/image_analysis/bragg_prediction/BraggPredictionRot.cpp b/image_analysis/bragg_prediction/BraggPredictionRot.cpp index c349fe3d..1743736e 100644 --- a/image_analysis/bragg_prediction/BraggPredictionRot.cpp +++ b/image_analysis/bragg_prediction/BraggPredictionRot.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "BraggPredictionRot.h" #include "../bragg_integration/SystematicAbsence.h" @@ -46,8 +47,8 @@ int BraggPredictionRot::Calc(const DiffractionExperiment &experiment, const Crys int i = 0; - const float mos_angle_rad = settings.mosaicity_deg * static_cast(M_PI) / 180.f; - const float half_wedge_angle_rad = settings.wedge_deg * static_cast(M_PI) / 180.f / 2.0f ; + const float mos_angle_rad = settings.mosaicity_deg * static_cast(PI) / 180.f; + const float half_wedge_angle_rad = settings.wedge_deg * static_cast(PI) / 180.f / 2.0f ; for (int h = -settings.max_hkl; h <= settings.max_hkl; h++) { // Precompute A* h contribution @@ -134,7 +135,7 @@ int BraggPredictionRot::Calc(const DiffractionExperiment &experiment, const Crys .h = h, .k = k, .l = l, - .delta_phi_deg = phi * 180.0f / static_cast(M_PI), + .delta_phi_deg = phi * 180.0f / static_cast(PI), .predicted_x = x, .predicted_y = y, .observed_x = NAN, diff --git a/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu b/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu index f86061f6..a06de212 100644 --- a/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu +++ b/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "BraggPredictionRotGPU.h" #ifdef JFJOCH_USE_CUDA @@ -150,7 +151,7 @@ namespace { out[count].h = h; out[count].k = k; out[count].l = l; - out[count].delta_phi_deg = phi * 180.0 / M_PI; + out[count].delta_phi_deg = phi * 180.0 / PI; out[count].predicted_x = x; out[count].predicted_y = y; out[count].observed_x = NAN; @@ -209,8 +210,8 @@ namespace { kc.one_over_wavelength = 1.0f / geom.GetWavelength_A(); // Store mosaicity and wedge in radians for partiality calculation - kc.mos_angle_rad = settings.mosaicity_deg * static_cast(M_PI) / 180.0f; - kc.wedge_angle_rad = settings.wedge_deg * static_cast(M_PI) / 180.0f; + kc.mos_angle_rad = settings.mosaicity_deg * static_cast(PI) / 180.0f; + kc.wedge_angle_rad = settings.wedge_deg * static_cast(PI) / 180.0f; kc.min_zeta = settings.min_zeta; kc.mosaicity_multiplier = settings.mosaicity_multiplier; diff --git a/image_analysis/geom_refinement/AssignSpotsToRings.cpp b/image_analysis/geom_refinement/AssignSpotsToRings.cpp index c50d0ef0..01691ddd 100644 --- a/image_analysis/geom_refinement/AssignSpotsToRings.cpp +++ b/image_analysis/geom_refinement/AssignSpotsToRings.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "AssignSpotsToRings.h" #include @@ -163,7 +164,7 @@ std::vector CalculateXtalRings(const UnitCell &cell, int hkl_max) { for (int l = 0; l <= hkl_max; l++) { if (h == 0 && k == 0 && l == 0) continue; auto p = Astar * h + Bstar * k + Cstar * l; - float Q = 2.0f * M_PI * p.Length(); + float Q = 2.0f * PI * p.Length(); u.push_back(Q); } } @@ -222,7 +223,7 @@ void GuessGeometry(DiffractionGeometry &geom, const std::vector &v, std::vector ring_q = CalculateXtalRings(calibrant_uc); if (ring_q.empty()) throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Unit cell problem"); - auto cluster_annot = GuessInitialGeometry(geom, v, 2 * M_PI / ring_q[0]); + auto cluster_annot = GuessInitialGeometry(geom, v, 2 * PI / ring_q[0]); std::vector optimizer_input; @@ -231,7 +232,7 @@ void GuessGeometry(DiffractionGeometry &geom, const std::vector &v, while (cluster_idx < cluster_annot.size() && ring_idx < ring_q.size()) { - float obs_q = 2 * M_PI / geom.PxlToRes(cluster_annot[cluster_idx].R_obs); + float obs_q = 2 * PI / geom.PxlToRes(cluster_annot[cluster_idx].R_obs); if (std::fabs(ring_q[ring_idx] - obs_q) < 0.1) { std::cout << "Found ring " << ring_idx << " with q_ring = " << ring_q[ring_idx] @@ -263,7 +264,7 @@ void OptimizeGeometry(DiffractionGeometry &geom, const std::vector & std::vector optimizer_input; for (const auto& s: v) { - float q_obs = 2 * M_PI / geom.PxlToRes(s.x, s.y); + float q_obs = 2 * PI / geom.PxlToRes(s.x, s.y); for (const auto &q : q_ring) { if (std::fabs(q - q_obs) < 0.1) { optimizer_input.push_back({s.x, s.y, q}); diff --git a/image_analysis/geom_refinement/LatticeReduction.cpp b/image_analysis/geom_refinement/LatticeReduction.cpp index 7da10261..d0cd4836 100644 --- a/image_analysis/geom_refinement/LatticeReduction.cpp +++ b/image_analysis/geom_refinement/LatticeReduction.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "LatticeReduction.h" #include @@ -251,7 +252,7 @@ CrystalLattice AngleAxisAndCellToLattice(const double rod[3], CrystalLattice AngleAxisAndLengthsToLattice(const double rod[3], const double lengths[3], bool hex) { return AngleAxisAndCellToLattice(rod, lengths, - /*alpha=*/M_PI / 2.0, - /*beta =*/M_PI / 2.0, - /*gamma=*/hex ? 2.0 * M_PI / 3.0 : M_PI / 2.0); + /*alpha=*/PI / 2.0, + /*beta =*/PI / 2.0, + /*gamma=*/hex ? 2.0 * PI / 3.0 : PI / 2.0); } diff --git a/image_analysis/geom_refinement/RingOptimizer.cpp b/image_analysis/geom_refinement/RingOptimizer.cpp index af214bf4..626f66f8 100644 --- a/image_analysis/geom_refinement/RingOptimizer.cpp +++ b/image_analysis/geom_refinement/RingOptimizer.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "RingOptimizer.h" #include "ceres/ceres.h" @@ -11,7 +12,7 @@ struct RingResidual { : obs_x(x), obs_y(y), lambda(lambda), pixel_size(pixel_size), - expected_len_recip_sq(expected_q * expected_q / (4.0 * M_PI * M_PI)) {} + expected_len_recip_sq(expected_q * expected_q / (4.0 * PI * PI)) {} template bool operator()(const T* const center_x, const T* const center_y, diff --git a/image_analysis/geom_refinement/XtalOptimizer.cpp b/image_analysis/geom_refinement/XtalOptimizer.cpp index 03d5ce08..a4676450 100644 --- a/image_analysis/geom_refinement/XtalOptimizer.cpp +++ b/image_analysis/geom_refinement/XtalOptimizer.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include #include "XtalOptimizer.h" @@ -293,9 +294,9 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, // Triclinic: initialize a,b,c and α,β,γ from current unit cell LatticeToRodriguesAndLengths_GS(data.latt, latt_vec0, latt_vec1); auto uc = data.latt.GetUnitCell(); - latt_vec2[0] = uc.alpha * M_PI / 180.0; - latt_vec2[1] = uc.beta * M_PI / 180.0; - latt_vec2[2] = uc.gamma * M_PI / 180.0; + latt_vec2[0] = uc.alpha * PI / 180.0; + latt_vec2[1] = uc.beta * PI / 180.0; + latt_vec2[2] = uc.gamma * PI / 180.0; break; } @@ -316,7 +317,7 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, if (data.axis) { const float angle_deg = data.axis->GetAngle_deg(i) + data.axis->GetWedge_deg() / 2.0; - angle_rad = angle_deg * M_PI / 180.0; + angle_rad = angle_deg * PI / 180.0; rot_matr = data.axis->GetTransformationAngle(angle_deg); } @@ -380,7 +381,7 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, if (!data.refine_detector_angles) { problem.SetParameterBlockConstant(detector_rot); } else { - const double rot_range = 3.0 / 180.0 * M_PI; + const double rot_range = 3.0 / 180.0 * PI; for (int i = 0; i < 2; ++i) { problem.SetParameterLowerBound(detector_rot, i, detector_rot[i] - rot_range); problem.SetParameterUpperBound(detector_rot, i, detector_rot[i] + rot_range); @@ -403,14 +404,14 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, } if (data.crystal_system == gemmi::CrystalSystem::Monoclinic) { - const double beta_lo = std::max(1e-6, M_PI * (data.min_angle_deg / 180.0)); - const double beta_hi = std::min(M_PI - 1e-6, M_PI * (data.max_angle_deg / 180.0)); + const double beta_lo = std::max(1e-6, PI * (data.min_angle_deg / 180.0)); + const double beta_hi = std::min(PI - 1e-6, PI * (data.max_angle_deg / 180.0)); problem.SetParameterLowerBound(latt_vec2, 0, beta_lo); problem.SetParameterUpperBound(latt_vec2, 0, beta_hi); } else if (data.crystal_system == gemmi::CrystalSystem::Triclinic) { // α, β, γ bounds (radians) - const double alo = M_PI * (data.min_angle_deg / 180.0); - const double ahi = M_PI * (data.max_angle_deg / 180.0); + const double alo = PI * (data.min_angle_deg / 180.0); + const double ahi = PI * (data.max_angle_deg / 180.0); for (int i = 0; i < 3; ++i) { problem.SetParameterLowerBound(latt_vec2, i, alo); problem.SetParameterUpperBound(latt_vec2, i, ahi); @@ -450,19 +451,19 @@ bool XtalOptimizerInternal(XtalOptimizerData &data, data.axis.value().Axis(Coord(rot_vec[0], rot_vec[1], rot_vec[2])); if (data.crystal_system == gemmi::CrystalSystem::Orthorhombic) - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, M_PI / 2.0, M_PI / 2.0, M_PI / 2.0); + data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); else if (data.crystal_system == gemmi::CrystalSystem::Tetragonal) { latt_vec1[1] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, M_PI / 2.0, M_PI / 2.0, M_PI / 2.0); + data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); } else if (data.crystal_system == gemmi::CrystalSystem::Cubic) { latt_vec1[1] = latt_vec1[0]; latt_vec1[2] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, M_PI / 2.0, M_PI / 2.0, M_PI / 2.0); + data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, PI / 2.0, PI / 2.0); } else if (data.crystal_system == gemmi::CrystalSystem::Hexagonal) { latt_vec1[1] = latt_vec1[0]; - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1,M_PI / 2.0, M_PI / 2.0, 2.0 * M_PI / 3.0); + data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1,PI / 2.0, PI / 2.0, 2.0 * PI / 3.0); } else if (data.crystal_system == gemmi::CrystalSystem::Monoclinic) { - data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, M_PI / 2.0, latt_vec2[0], M_PI / 2.0); + data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, PI / 2.0, latt_vec2[0], PI / 2.0); } else { // Triclinic via the same generic builder data.latt = AngleAxisAndCellToLattice(latt_vec0, latt_vec1, latt_vec2[0], latt_vec2[1], latt_vec2[2]); diff --git a/image_analysis/indexing/AnalyzeIndexing.cpp b/image_analysis/indexing/AnalyzeIndexing.cpp index 2a11ffd1..5da8f709 100644 --- a/image_analysis/indexing/AnalyzeIndexing.cpp +++ b/image_analysis/indexing/AnalyzeIndexing.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include #include #include @@ -20,11 +21,11 @@ namespace { } inline float deg_to_rad(float deg) { - return deg * (static_cast(M_PI) / 180.0f); + return deg * (static_cast(PI) / 180.0f); } inline float rad_to_deg(float rad) { - return rad * (180.0f / static_cast(M_PI)); + return rad * (180.0f / static_cast(PI)); } // Wrap to [-180, 180] (useful for residuals) @@ -139,7 +140,7 @@ namespace { float s1 = phi_ref + delta; float s2 = phi_ref - delta; - const float two_pi = 2.0f * static_cast(M_PI); + const float two_pi = 2.0f * static_cast(PI); auto shift_near = [&](float x) { if (!std::isfinite(x)) return std::numeric_limits::quiet_NaN(); @@ -153,9 +154,9 @@ namespace { return std::numeric_limits::quiet_NaN(); // fmod can return negative values; normalize to [-π, π] - if (shifted < -static_cast(M_PI)) + if (shifted < -static_cast(PI)) shifted += two_pi; - else if (shifted > static_cast(M_PI)) + else if (shifted > static_cast(PI)) shifted -= two_pi; return shifted + span_center; diff --git a/image_analysis/indexing/FFTIndexer.cpp b/image_analysis/indexing/FFTIndexer.cpp index 2578d589..c5fbb0f9 100644 --- a/image_analysis/indexing/FFTIndexer.cpp +++ b/image_analysis/indexing/FFTIndexer.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "FFTIndexer.h" #include #include "PostIndexingRefinement.h" @@ -23,7 +24,7 @@ FFTIndexer::FFTIndexer(const IndexingSettings &settings) // from when the extent was mistakenly written as 2*pi/d (the Q convention); it is kept // because the exact amount sets which marginal frames index - rounding it to a nearby // integer shifts the indexing rate ~0.5-1% (validated on the lyso datasets). - const float oversampling = 2.0f * static_cast(M_PI); + const float oversampling = 2.0f * static_cast(PI); const float one_over_d_max = 1.0f / settings.GetFFT_HighResolution_A(); histogram_spacing = 1.0f / (2.0f * max_length_A); @@ -52,7 +53,7 @@ void FFTIndexer::SetupDirectionVectors() { direction_vectors.reserve(static_cast(nDirections)); const double phi = (1.0 + std::sqrt(5.0)) / 2.0; // Golden ratio - const double golden_angle = 2.0 * M_PI / phi; + const double golden_angle = 2.0 * PI / phi; for (int i = 0; i < nDirections; i++) { // Half-sphere distribution (z in [0,1]) @@ -148,7 +149,7 @@ std::vector FFTIndexer::FilterFFTResults() const { std::vector ret; // Remove vectors less than 5 deg apart, as most likely these are colinear - const float COS_5_DEG = std::cos(5.0f * M_PI / 180.0f); + const float COS_5_DEG = std::cos(5.0f * PI / 180.0f); // Minimum relative amplitude to accept a shorter vector (fundamental) // over a longer one (harmonic). diff --git a/image_analysis/indexing/PostIndexingRefinement.cpp b/image_analysis/indexing/PostIndexingRefinement.cpp index 45960eca..75e21f58 100644 --- a/image_analysis/indexing/PostIndexingRefinement.cpp +++ b/image_analysis/indexing/PostIndexingRefinement.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "PostIndexingRefinement.h" #include @@ -175,9 +176,9 @@ std::vector Refine(const std::vector &in_spots, } // Filter for wrong angles - float alpha = std::acos(cell_rows.row(1).normalized().dot(cell_rows.row(2).normalized())) * 180.0f / M_PI; - float beta = std::acos(cell_rows.row(0).normalized().dot(cell_rows.row(2).normalized())) * 180.0f / M_PI; - float gamma = std::acos(cell_rows.row(0).normalized().dot(cell_rows.row(1).normalized())) * 180.0f / M_PI; + float alpha = std::acos(cell_rows.row(1).normalized().dot(cell_rows.row(2).normalized())) * 180.0f / PI; + float beta = std::acos(cell_rows.row(0).normalized().dot(cell_rows.row(2).normalized())) * 180.0f / PI; + float gamma = std::acos(cell_rows.row(0).normalized().dot(cell_rows.row(1).normalized())) * 180.0f / PI; if (alpha < p.min_angle_deg || alpha > p.max_angle_deg || beta < p.min_angle_deg || beta > p.max_angle_deg || diff --git a/image_analysis/lattice_search/LatticeSearch.cpp b/image_analysis/lattice_search/LatticeSearch.cpp index c2257a68..1648338a 100644 --- a/image_analysis/lattice_search/LatticeSearch.cpp +++ b/image_analysis/lattice_search/LatticeSearch.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "LatticeSearch.h" #include @@ -315,9 +316,9 @@ LatticeSearchResult LatticeSearch(const CrystalLattice &L, double dist_tolerance if (c.cond_BC && fabs((uc_reduced.b - uc_reduced.c) / (0.5 * (uc_reduced.b + uc_reduced.c))) > dist_tolerance) ok = false; - double expected_alpha = acos(c.cond_D / sqrt(B*C)) * 180 / M_PI; - double expected_beta = acos(c.cond_E / sqrt(A*C)) * 180 / M_PI; - double expected_gamma = acos(c.cond_F / sqrt(A*B)) * 180 / M_PI; + double expected_alpha = acos(c.cond_D / sqrt(B*C)) * 180 / PI; + double expected_beta = acos(c.cond_E / sqrt(A*C)) * 180 / PI; + double expected_gamma = acos(c.cond_F / sqrt(A*B)) * 180 / PI; if (fabs(expected_alpha - uc_reduced.alpha) > angle_tolerance) ok = false; diff --git a/image_analysis/pixel_refinement/PixelRefine.cpp b/image_analysis/pixel_refinement/PixelRefine.cpp index b9014d82..59d92de2 100644 --- a/image_analysis/pixel_refinement/PixelRefine.cpp +++ b/image_analysis/pixel_refinement/PixelRefine.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "PixelRefine.h" #include @@ -499,7 +500,7 @@ void PixelRefine::Run(const T *image, beam, dist_mm, detector_rot, latt_vec0, latt_vec1, latt_vec2, q_sq, eps_r, eps_t_sq)) continue; - const double P_t = std::exp(-eps_t_sq / (R1 * R1)) / (M_PI * R1 * R1); + const double P_t = std::exp(-eps_t_sq / (R1 * R1)) / (PI * R1 * R1); const double R0_eff_sq = R0 * R0 + g.R_bw_sq; const double P_rad = std::exp(-eps_r * eps_r / R0_eff_sq); const double v = std::max(g.Ibkg, 1.0); @@ -597,7 +598,7 @@ void PixelRefine::Run(const T *image, // Tangential profile shape (area-normalized) -> the fit template, using the // per-reflection R1_eff (Term 2). const double R1 = g.R1_eff; - const double P_t = std::exp(-eps_t_sq / (R1 * R1)) / (M_PI * R1 * R1); + const double P_t = std::exp(-eps_t_sq / (R1 * R1)) / (PI * R1 * R1); prof_full += P_t; // whole shoebox, on- or off-detector // Only real, unmasked detector pixels carry signal. diff --git a/image_analysis/rotation_indexer/RotationIndexer.cpp b/image_analysis/rotation_indexer/RotationIndexer.cpp index f004959a..66f0dc98 100644 --- a/image_analysis/rotation_indexer/RotationIndexer.cpp +++ b/image_analysis/rotation_indexer/RotationIndexer.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "RotationIndexer.h" #include "../geom_refinement/XtalOptimizer.h" @@ -84,7 +85,7 @@ void RotationIndexer::RunIndexing() { break; // Ignore lattices oriented by less than 3.0 degree - if (l.rotation_vector.Length() < 3.0 * M_PI / 180.0) + if (l.rotation_vector.Length() < 3.0 * PI / 180.0) continue; RotMatrix rot(l.rotation_vector.Length(), l.rotation_vector.Normalize()); diff --git a/image_analysis/spot_finding/SpotUtils.cpp b/image_analysis/spot_finding/SpotUtils.cpp index 1605b0f1..8dcf3254 100644 --- a/image_analysis/spot_finding/SpotUtils.cpp +++ b/image_analysis/spot_finding/SpotUtils.cpp @@ -1,6 +1,7 @@ // SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute // SPDX-License-Identifier: GPL-3.0-only +#include "../../common/JFJochMath.h" #include "SpotUtils.h" #include "../../common/ResolutionShells.h" @@ -25,10 +26,10 @@ void MarkIceRings(std::vector &spots, float tolerance_q_recipA) { std::vector ice_rings_q; for (const auto &i: ICE_RING_RES_A) - ice_rings_q.push_back(2 * M_PI / i); + ice_rings_q.push_back(2 * PI / i); for (auto &s: spots) { - auto spot_q = 2 * M_PI / s.d_A; + auto spot_q = 2 * PI / s.d_A; bool tmp = false; for (const auto &q: ice_rings_q) tmp |= (fabs(spot_q - q) < tolerance_q_recipA); @@ -65,8 +66,8 @@ void FilterSpuriousHighResolutionSpots(std::vector &spots, float thr // Avoid division by zero; d_A should be > 0 in valid data if (d1 <= 0.0f || d2 <= 0.0f) continue; - float q1 = 2 * M_PI / d1; - float q2 = 2 * M_PI / d2; + float q1 = 2 * PI / d1; + float q2 = 2 * PI / d2; float dq = q2 - q1; // should be >= 0 due to sorting if (dq > threshold) { cut_index = i + 1; // keep up to i inclusive diff --git a/viewer/JFJochImageReadingWorker.cpp b/viewer/JFJochImageReadingWorker.cpp index 67496478..f6aaedef 100644 --- a/viewer/JFJochImageReadingWorker.cpp +++ b/viewer/JFJochImageReadingWorker.cpp @@ -3,7 +3,7 @@ #include #include -#include +#include "../common/JFJochMath.h" #ifndef _WIN32 #include #include @@ -539,7 +539,7 @@ void JFJochImageReadingWorker::FindCenter(const UnitCell& calibrant, bool guess) QVector rings; for (int i = 0; i < 15 && i < ring_Q.size(); i++) { - rings.push_back(2 * std::numbers::pi / ring_Q[i]); + rings.push_back(2 * PI / ring_Q[i]); } emit setRings(rings); } diff --git a/viewer/image_viewer/JFJochAzIntImage.cpp b/viewer/image_viewer/JFJochAzIntImage.cpp index 74c86733..4588de53 100644 --- a/viewer/image_viewer/JFJochAzIntImage.cpp +++ b/viewer/image_viewer/JFJochAzIntImage.cpp @@ -5,7 +5,7 @@ #include #include "../../common/JFJochException.h" -#include +#include "../../common/JFJochMath.h" JFJochAzIntImage::JFJochAzIntImage(QWidget *parent) : JFJochImage(parent) { @@ -110,7 +110,7 @@ void JFJochAzIntImage::mouseDoubleClickEvent(QMouseEvent *event) { float phi = image->Dataset().az_int_bin_to_phi[idx]; auto geom = image->Dataset().experiment.GetDiffractionGeometry(); - auto coord = geom.ResPhiToPxl(2 * std::numbers::pi / q, phi / 180.0 * std::numbers::pi); + auto coord = geom.ResPhiToPxl(2 * PI / q, phi / 180.0 * PI); emit zoomOnBin(QPointF(coord.first, coord.second)); } } \ No newline at end of file diff --git a/viewer/image_viewer/JFJochDiffractionImage.cpp b/viewer/image_viewer/JFJochDiffractionImage.cpp index 51dde0f3..9b744348 100644 --- a/viewer/image_viewer/JFJochDiffractionImage.cpp +++ b/viewer/image_viewer/JFJochDiffractionImage.cpp @@ -3,7 +3,7 @@ #include "JFJochDiffractionImage.h" #include "../../common/DiffractionGeometry.h" -#include +#include "../../common/JFJochMath.h" #include #include @@ -240,9 +240,9 @@ void JFJochDiffractionImage::DrawResolutionRings() { continue; auto [x1,y1] = geom.ResPhiToPxl(d, 0); - auto [x2,y2] = geom.ResPhiToPxl(d, std::numbers::pi / 2); - auto [x3,y3] = geom.ResPhiToPxl(d, std::numbers::pi); - auto [x4,y4] = geom.ResPhiToPxl(d, 3.0 * std::numbers::pi / 2); + auto [x2,y2] = geom.ResPhiToPxl(d, PI / 2); + auto [x3,y3] = geom.ResPhiToPxl(d, PI); + auto [x4,y4] = geom.ResPhiToPxl(d, 3.0 * PI / 2); auto x_min = std::min({x1, x2, x3, x4}); auto x_max = std::max({x1, x2, x3, x4}); @@ -253,9 +253,9 @@ void JFJochDiffractionImage::DrawResolutionRings() { addOverlayItem(scene()->addEllipse(boundingRect, pen)); auto [x5,y5] = geom.ResPhiToPxl(d, phi_offset + 0); - auto [x6,y6] = geom.ResPhiToPxl(d, phi_offset + std::numbers::pi / 2); - auto [x7,y7] = geom.ResPhiToPxl(d, phi_offset + std::numbers::pi); - auto [x8,y8] = geom.ResPhiToPxl(d, phi_offset + 3.0 * std::numbers::pi / 2); + auto [x6,y6] = geom.ResPhiToPxl(d, phi_offset + PI / 2); + auto [x7,y7] = geom.ResPhiToPxl(d, phi_offset + PI); + auto [x8,y8] = geom.ResPhiToPxl(d, phi_offset + 3.0 * PI / 2); QPointF point_1(x5, y5); QPointF point_2(x6, y6); @@ -286,7 +286,7 @@ void JFJochDiffractionImage::DrawResolutionRings() { scene()->addItem(textItem); addOverlayItem(textItem); } - phi_offset += 4.0 / 180.0 * std::numbers::pi; + phi_offset += 4.0 / 180.0 * PI; } } diff --git a/viewer/widgets/PowderCalibrationWidget.cpp b/viewer/widgets/PowderCalibrationWidget.cpp index b7a78076..286f3c8a 100644 --- a/viewer/widgets/PowderCalibrationWidget.cpp +++ b/viewer/widgets/PowderCalibrationWidget.cpp @@ -5,7 +5,7 @@ #include #include #include -#include +#include "../../common/JFJochMath.h" #include "../image_analysis/geom_refinement/AssignSpotsToRings.h" @@ -26,7 +26,7 @@ PowderCalibrationWidget::PowderCalibrationWidget(QWidget *parent) : QWidget(pare std::vector rings = CalculateXtalRings(GetCalibrant(), 10); QVector q_rings; for (float ring : rings) { - q_rings.append(2 * std::numbers::pi / ring); + q_rings.append(2 * PI / ring); } emit ringsFromCalibration(q_rings); }); diff --git a/viewer/windows/JFJochViewerMetadataWindow.cpp b/viewer/windows/JFJochViewerMetadataWindow.cpp index 97db7e02..47f214ad 100644 --- a/viewer/windows/JFJochViewerMetadataWindow.cpp +++ b/viewer/windows/JFJochViewerMetadataWindow.cpp @@ -10,7 +10,7 @@ #include #include "JFJochViewerMetadataWindow.h" -#include +#include "../../common/JFJochMath.h" JFJochViewerMetadataWindow::JFJochViewerMetadataWindow(QWidget *parent) : JFJochHelperWindow(parent) { @@ -136,8 +136,8 @@ void JFJochViewerMetadataWindow::datasetUpdate() { tmp.BeamX_pxl(beam_center_x->value()); tmp.BeamY_pxl(beam_center_y->value()); tmp.IncidentEnergy_keV(WVL_1A_IN_KEV / wavelength_A->value()); - tmp.PoniRot1_rad(poni_rot1_deg->value() * std::numbers::pi / 180.0); - tmp.PoniRot2_rad(poni_rot2_deg->value() * std::numbers::pi / 180.0); + tmp.PoniRot1_rad(poni_rot1_deg->value() * PI / 180.0); + tmp.PoniRot2_rad(poni_rot2_deg->value() * PI / 180.0); tmp.DetectIceRings(detect_ice_rings->isChecked()); if (unit_cell_enabled->isChecked()) { @@ -177,8 +177,8 @@ void JFJochViewerMetadataWindow::datasetLoaded(std::shared_ptrsetValue(dataset->experiment.GetBeamX_pxl()); beam_center_y->setValue(dataset->experiment.GetBeamY_pxl()); wavelength_A->setValue(dataset->experiment.GetWavelength_A()); - poni_rot1_deg->setValue(dataset->experiment.GetPoniRot1_rad() * 180.0 / std::numbers::pi); - poni_rot2_deg->setValue(dataset->experiment.GetPoniRot2_rad() * 180.0 / std::numbers::pi); + poni_rot1_deg->setValue(dataset->experiment.GetPoniRot1_rad() * 180.0 / PI); + poni_rot2_deg->setValue(dataset->experiment.GetPoniRot2_rad() * 180.0 / PI); detect_ice_rings->setChecked(dataset->experiment.IsDetectIceRings()); auto unit_cell = dataset->experiment.GetUnitCell();