Line data Source code
1 : // SPDX-License-Identifier: BSD-3-Clause
2 : // Copyright (c) 2023 Scipp contributors (https://github.com/scipp)
3 : /// @file
4 : /// @author Simon Heybrock
5 : #pragma once
6 :
7 : #include "scipp/common/numeric.h"
8 : #include "scipp/common/overloaded.h"
9 : #include "scipp/core/eigen.h"
10 : #include "scipp/core/element/arg_list.h"
11 : #include "scipp/core/spatial_transforms.h"
12 : #include "scipp/core/subbin_sizes.h"
13 : #include "scipp/core/transform_common.h"
14 : #include "scipp/units/except.h"
15 :
16 : namespace scipp::core::element {
17 :
18 : template <class... Extra>
19 : constexpr auto add_inplace_types =
20 : arg_list<double, float, int64_t, int32_t, Eigen::Vector3d,
21 : std::tuple<scipp::core::time_point, int64_t>,
22 : std::tuple<scipp::core::time_point, int32_t>,
23 : std::tuple<double, float>, std::tuple<float, double>,
24 : std::tuple<int64_t, int32_t>, std::tuple<int32_t, int64_t>,
25 : std::tuple<double, int64_t>, std::tuple<double, int32_t>,
26 : std::tuple<float, int64_t>, std::tuple<float, int32_t>,
27 : std::tuple<double, bool>, std::tuple<int64_t, bool>, Extra...>;
28 :
29 : constexpr auto add_equals = overloaded{add_inplace_types<SubbinSizes>,
30 51803682 : [](auto &&a, const auto &b) { a += b; }};
31 :
32 : constexpr auto nan_add_equals =
33 2249 : overloaded{add_inplace_types<>, [](auto &&a, const auto &b) {
34 : using numeric::isnan;
35 2249 : if (isnan(a))
36 0 : a = std::decay_t<decltype(a)>{0}; // Force zero
37 2249 : if (!isnan(b))
38 2055 : a += b;
39 2249 : }};
40 :
41 : constexpr auto subtract_equals =
42 2412 : overloaded{add_inplace_types<>, [](auto &&a, const auto &b) { a -= b; }};
43 :
44 : constexpr auto mul_inplace_types = arg_list<
45 : double, float, int64_t, int32_t, Eigen::Matrix3d, std::tuple<double, float>,
46 : std::tuple<float, double>, std::tuple<int64_t, int32_t>,
47 : std::tuple<int64_t, bool>, std::tuple<int32_t, int64_t>,
48 : std::tuple<double, int64_t>, std::tuple<double, int32_t>,
49 : std::tuple<float, int64_t>, std::tuple<float, int32_t>,
50 : std::tuple<Eigen::Vector3d, double>, std::tuple<Eigen::Vector3d, float>,
51 : std::tuple<Eigen::Vector3d, int64_t>, std::tuple<Eigen::Vector3d, int32_t>>;
52 :
53 : // Note that we do *not* support any integer type as left-hand-side, to match
54 : // Python 3 and numpy "truediv" behavior.
55 : constexpr auto div_inplace_types = arg_list<
56 : double, float, std::tuple<double, float>, std::tuple<float, double>,
57 : std::tuple<double, int64_t>, std::tuple<double, int32_t>,
58 : std::tuple<float, int64_t>, std::tuple<float, int32_t>,
59 : std::tuple<Eigen::Vector3d, double>, std::tuple<Eigen::Vector3d, float>,
60 : std::tuple<Eigen::Vector3d, int64_t>, std::tuple<Eigen::Vector3d, int32_t>>;
61 :
62 : constexpr auto floor_div_inplace_types =
63 : arg_list<double, float, int64_t, int32_t, std::tuple<double, float>,
64 : std::tuple<float, double>, std::tuple<double, int64_t>,
65 : std::tuple<double, int32_t>, std::tuple<float, int64_t>,
66 : std::tuple<float, int32_t>>;
67 :
68 : constexpr auto multiply_equals =
69 3150 : overloaded{mul_inplace_types, [](auto &&a, const auto &b) { a *= b; }};
70 : constexpr auto divide_equals =
71 2268 : overloaded{div_inplace_types, [](auto &&a, const auto &b) { a /= b; }};
72 : constexpr auto floor_divide_equals = overloaded{
73 : floor_div_inplace_types, transform_flags::expect_no_variance_arg<0>,
74 : transform_flags::expect_no_variance_arg<1>,
75 5 : [](auto &&a, const auto &b) { a = numeric::floor_divide(a, b); }};
76 :
77 : using arithmetic_and_matrix_type_pairs = decltype(std::tuple_cat(
78 : std::declval<arithmetic_type_pairs>(),
79 : std::tuple<std::tuple<Eigen::Vector3d, Eigen::Vector3d>>()));
80 :
81 : struct add_types_t {
82 : constexpr void operator()() const noexcept;
83 : using types = decltype(std::tuple_cat(
84 : std::declval<arithmetic_and_matrix_type_pairs>(),
85 : std::tuple<
86 : std::tuple<time_point, int64_t>, std::tuple<time_point, int32_t>,
87 : std::tuple<int64_t, time_point>, std::tuple<int32_t, time_point>>{}));
88 : };
89 :
90 : struct subtract_types_t {
91 : constexpr void operator()() const noexcept;
92 : using types = decltype(std::tuple_cat(
93 : std::declval<arithmetic_and_matrix_type_pairs>(),
94 : std::tuple<std::tuple<time_point, int64_t>,
95 : std::tuple<time_point, int32_t>,
96 : std::tuple<time_point, time_point>>{}));
97 : };
98 :
99 : struct multiplies_types_t {
100 : constexpr void operator()() const noexcept;
101 : using types = decltype(std::tuple_cat(
102 : std::declval<arithmetic_type_pairs_with_bool>(),
103 : std::tuple<std::tuple<double, Eigen::Vector3d>>(),
104 : std::tuple<std::tuple<float, Eigen::Vector3d>>(),
105 : std::tuple<std::tuple<int64_t, Eigen::Vector3d>>(),
106 : std::tuple<std::tuple<int32_t, Eigen::Vector3d>>(),
107 : std::tuple<std::tuple<Eigen::Vector3d, double>>(),
108 : std::tuple<std::tuple<Eigen::Vector3d, float>>(),
109 : std::tuple<std::tuple<Eigen::Vector3d, int64_t>>(),
110 : std::tuple<std::tuple<Eigen::Vector3d, int32_t>>(),
111 : std::declval<
112 : pair_product_t<Eigen::Matrix3d, Eigen::Affine3d,
113 : scipp::core::Quaternion, scipp::core::Translation>>(),
114 : std::tuple<std::tuple<scipp::core::Quaternion, Eigen::Vector3d>>(),
115 : std::tuple<std::tuple<Eigen::Matrix3d, Eigen::Vector3d>>()));
116 : };
117 :
118 : struct apply_spatial_transformation_t {
119 : constexpr void operator()() const noexcept;
120 : using types = decltype(std::tuple_cat(
121 : std::tuple<std::tuple<Eigen::Affine3d, Eigen::Vector3d>>(),
122 : std::tuple<std::tuple<scipp::core::Translation, Eigen::Vector3d>>(),
123 :
124 : std::tuple<
125 : std::tuple<scipp::core::Translation, scipp::core::Translation>>(),
126 : std::tuple<std::tuple<scipp::core::Translation, Eigen::Affine3d>>(),
127 :
128 : std::tuple<std::tuple<Eigen::Affine3d, Eigen::Affine3d>>(),
129 : std::tuple<std::tuple<Eigen::Affine3d, scipp::core::Translation>>()));
130 : };
131 :
132 : struct true_divide_types_t {
133 : constexpr void operator()() const noexcept;
134 : using types = decltype(std::tuple_cat(
135 : std::declval<arithmetic_type_pairs>(),
136 : std::tuple<std::tuple<Eigen::Vector3d, double>>(),
137 : std::tuple<std::tuple<Eigen::Vector3d, float>>(),
138 : std::tuple<std::tuple<Eigen::Vector3d, int64_t>>(),
139 : std::tuple<std::tuple<Eigen::Vector3d, int32_t>>()));
140 : };
141 :
142 : struct floor_divide_types_t {
143 : constexpr void operator()() const noexcept;
144 : using types = arithmetic_type_pairs;
145 : };
146 :
147 : struct remainder_types_t {
148 : constexpr void operator()() const noexcept;
149 : using types = arithmetic_type_pairs;
150 : };
151 :
152 : constexpr auto add =
153 13331665 : overloaded{add_types_t{}, [](const auto a, const auto b) { return a + b; }};
154 : constexpr auto subtract = overloaded{
155 22137155 : subtract_types_t{}, [](const auto a, const auto b) { return a - b; }};
156 : constexpr auto multiply = overloaded{
157 : multiplies_types_t{},
158 : transform_flags::expect_no_in_variance_if_out_cannot_have_variance,
159 3931043 : [](const auto a, const auto b) { return a * b; }};
160 :
161 : constexpr auto apply_spatial_transformation = overloaded{
162 : apply_spatial_transformation_t{},
163 : transform_flags::expect_no_in_variance_if_out_cannot_have_variance,
164 20 : [](const auto a, const auto b) { return a * b; },
165 20 : [](const units::Unit &a, const units::Unit &b) {
166 20 : if (a != b)
167 2 : throw except::UnitError(
168 : "Cannot apply spatial transform as the units of the transformation "
169 4 : "are not the same as the units of transformation or vector.");
170 : else
171 18 : return a;
172 : }};
173 :
174 : // truediv defined as in Python.
175 : constexpr auto divide = overloaded{
176 : true_divide_types_t{},
177 : transform_flags::expect_no_in_variance_if_out_cannot_have_variance,
178 656426 : [](const auto &a, const auto &b) { return numeric::true_divide(a, b); },
179 6247 : [](const units::Unit &a, const units::Unit &b) { return a / b; }};
180 :
181 : // floordiv defined as in Python. Complementary to mod.
182 : constexpr auto floor_divide = overloaded{
183 : floor_divide_types_t{}, transform_flags::expect_no_variance_arg<0>,
184 : transform_flags::expect_no_variance_arg<1>,
185 24103515 : [](const auto a, const auto b) { return numeric::floor_divide(a, b); },
186 110 : [](const units::Unit &a, const units::Unit &b) { return a / b; }};
187 :
188 : // remainder defined as in Python
189 : constexpr auto mod = overloaded{
190 : remainder_types_t{}, transform_flags::expect_no_variance_arg<0>,
191 : transform_flags::expect_no_variance_arg<1>,
192 24099583 : [](const auto a, const auto b) { return numeric::remainder(a, b); },
193 69 : [](const units::Unit &a, const units::Unit &b) { return a % b; }};
194 :
195 : constexpr auto remainder_inplace_types =
196 : arg_list<double, float, int64_t, int32_t, std::tuple<double, float>,
197 : std::tuple<float, double>, std::tuple<double, int64_t>,
198 : std::tuple<double, int32_t>, std::tuple<float, int64_t>,
199 : std::tuple<float, int32_t>, std::tuple<int32_t, int64_t>,
200 : std::tuple<int64_t, int32_t>>;
201 :
202 : constexpr auto mod_equals = overloaded{
203 : remainder_inplace_types, transform_flags::expect_no_variance_arg<0>,
204 : transform_flags::expect_no_variance_arg<1>,
205 18 : [](units::Unit &a, const units::Unit &b) { a %= b; },
206 23826824 : [](auto &&a, const auto &b) { a = mod(a, b); }};
207 :
208 : constexpr auto negative =
209 : overloaded{arg_list<double, float, int64_t, int32_t, Eigen::Vector3d>,
210 2370 : [](const auto x) { return -x; }};
211 :
212 : } // namespace scipp::core::element
|