Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ref: re-template algebra #893

Merged
merged 1 commit into from
Dec 18, 2024
Merged

ref: re-template algebra #893

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 5 additions & 2 deletions core/include/detray/builders/bin_fillers.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -122,10 +122,13 @@ struct bin_associator {
const transform_container &transforms,
const mask_container &masks,
const context_t ctx, Args &&...) const -> void {

using scalar_t = dscalar<typename grid_t::algebra_type>;

// Fill the surfaces into the grid by matching their contour onto the
// grid bins
bin_association(ctx, surfaces, transforms, masks, grid, {0.1f, 0.1f},
false);
bin_association(ctx, surfaces, transforms, masks, grid,
std::array<scalar_t, 2>{0.1f, 0.1f}, false);
niermann999 marked this conversation as resolved.
Show resolved Hide resolved
}
};

Expand Down
13 changes: 7 additions & 6 deletions core/include/detray/builders/cuboid_portal_generator.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,13 +30,14 @@ template <typename detector_t>
class cuboid_portal_generator final
: public surface_factory_interface<detector_t> {

using scalar_type = typename detector_t::scalar_type;
using transform3_type = typename detector_t::transform3_type;
using algebra_type = typename detector_t::algebra_type;
using scalar_type = dscalar<algebra_type>;
using transform3_type = dtransform3D<algebra_type>;
niermann999 marked this conversation as resolved.
Show resolved Hide resolved

/// A functor to construct global bounding boxes around masks
struct bounding_box_creator {

using aabb_t = axis_aligned_bounding_volume<cuboid3D>;
using aabb_t = axis_aligned_bounding_volume<cuboid3D, algebra_type>;

template <typename mask_group_t, typename index_t>
DETRAY_HOST_DEVICE inline void operator()(
Expand Down Expand Up @@ -86,15 +87,15 @@ class cuboid_portal_generator final
typename detector_t::geometry_context ctx = {})
-> dindex_range override {

using point3_t = typename detector_t::point3_type;
using vector3_t = typename detector_t::vector3_type;
using point3_t = dpoint3D<algebra_type>;
using vector3_t = dvector3D<algebra_type>;
niermann999 marked this conversation as resolved.
Show resolved Hide resolved

using surface_t = typename detector_t::surface_type;
using nav_link_t = typename surface_t::navigation_link;
using mask_link_t = typename surface_t::mask_link;
using material_link_t = typename surface_t::material_link;

using aabb_t = axis_aligned_bounding_volume<cuboid3D>;
using aabb_t = axis_aligned_bounding_volume<cuboid3D, algebra_type>;

constexpr auto invalid_src_link{detail::invalid_value<std::uint64_t>()};

Expand Down
19 changes: 10 additions & 9 deletions core/include/detray/builders/cylinder_portal_generator.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -110,15 +110,16 @@ template <typename detector_t>
class cylinder_portal_generator final
: public surface_factory_interface<detector_t> {

using scalar_t = typename detector_t::scalar_type;
using point3_t = typename detector_t::point3_type;
using vector3_t = typename detector_t::vector3_type;
using transform3_t = typename detector_t::transform3_type;
using algebra_t = typename detector_t::algebra_type;
using scalar_t = dscalar<algebra_t>;
using point3_t = dpoint3D<algebra_t>;
using vector3_t = dvector3D<algebra_t>;
using transform3_t = dtransform3D<algebra_t>;
niermann999 marked this conversation as resolved.
Show resolved Hide resolved

/// A functor to construct global bounding boxes around masks
struct bounding_box_creator {

using aabb_t = axis_aligned_bounding_volume<cuboid3D>;
using aabb_t = axis_aligned_bounding_volume<cuboid3D, algebra_t>;

template <typename mask_group_t, typename index_t>
DETRAY_HOST_DEVICE inline void operator()(
Expand Down Expand Up @@ -182,7 +183,7 @@ class cylinder_portal_generator final
typename detector_t::geometry_context ctx = {})
-> dindex_range override {

using aabb_t = axis_aligned_bounding_volume<cuboid3D>;
using aabb_t = axis_aligned_bounding_volume<cuboid3D, algebra_t>;

// Only build portals for cylinder volumes
assert(volume.id() == volume_id::e_cylinder);
Expand Down Expand Up @@ -227,9 +228,9 @@ class cylinder_portal_generator final

// Get the half lengths for the cylinder height and disc translation
const point3_t h_lengths = 0.5f * (box_max - box_min);
const scalar h_x{math::fabs(h_lengths[0])};
const scalar h_y{math::fabs(h_lengths[1])};
const scalar h_z{math::fabs(h_lengths[2])};
const scalar_t h_x{math::fabs(h_lengths[0])};
const scalar_t h_y{math::fabs(h_lengths[1])};
const scalar_t h_z{math::fabs(h_lengths[2])};

const scalar_t outer_r_min{math::max(h_x, h_y)};
const scalar_t mean_radius{get_mean_radius(surfaces, transforms)};
Expand Down
100 changes: 54 additions & 46 deletions core/include/detray/builders/detail/associator.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,31 +18,34 @@

namespace detray::detail {

/** Struct that assigns the center of gravity to a rectangular bin */
/// Struct that assigns the center of gravity to a rectangular bin
template <typename algebra_t>
struct center_of_gravity_rectangle {
/** Call operator to the struct, allows to chain several chain operators
* together
*
* @param bin_contour The contour description of the bin -> target
* @param surface_contour The contour description of the surface -> test
*
* @note the bin_contour is asummed to be a rectangle
*
* @return whether this should be associated
*/
template <typename point2_t>
bool operator()(const std::vector<point2_t> &bin_contour,
const std::vector<point2_t> &surface_contour) {
/// Call operator to the struct, allows to chain several chain operators
/// together
///
/// @param bin_contour The contour description of the bin -> target
/// @param surface_contour The contour description of the surface -> test
///
/// @note the bin_contour is asummed to be a rectangle
///
/// @return whether this should be associated
bool operator()(const std::vector<dpoint2D<algebra_t>> &bin_contour,
const std::vector<dpoint2D<algebra_t>> &surface_contour) {

using scalar_t = dscalar<algebra_t>;
using point2_t = dpoint2D<algebra_t>;

// Check if centre of gravity is inside bin
point2_t cgs = {0.f, 0.f};
for (const auto &svtx : surface_contour) {
cgs = cgs + svtx;
}
cgs = 1.f / static_cast<scalar>(surface_contour.size()) * cgs;
scalar min_l0 = std::numeric_limits<scalar>::max();
scalar max_l0 = -std::numeric_limits<scalar>::max();
scalar min_l1 = std::numeric_limits<scalar>::max();
scalar max_l1 = -std::numeric_limits<scalar>::max();
cgs = 1.f / static_cast<scalar_t>(surface_contour.size()) * cgs;
scalar_t min_l0 = std::numeric_limits<scalar_t>::max();
scalar_t max_l0 = -std::numeric_limits<scalar_t>::max();
scalar_t min_l1 = std::numeric_limits<scalar_t>::max();
scalar_t max_l1 = -std::numeric_limits<scalar_t>::max();
for (const auto &b : bin_contour) {
min_l0 = math::min(b[0], min_l0);
max_l0 = math::max(b[0], max_l0);
Expand All @@ -59,25 +62,28 @@ struct center_of_gravity_rectangle {
}
};

/** Check if center of mass is inside a generic polygon bin */
/// Check if center of mass is inside a generic polygon bin
template <typename algebra_t>
struct center_of_gravity_generic {
/** Call operator to the struct, allows to chain several chain operators
* together
*
* @param bin_contour The contour description of the bin -> target
* @param surface_contour The contour description of the surface -> test
*
* @return whether this should be associated
*/
template <typename point2_t>
bool operator()(const std::vector<point2_t> &bin_contour,
const std::vector<point2_t> &surface_contour) {
/// Call operator to the struct, allows to chain several chain operators
/// together
///
/// @param bin_contour The contour description of the bin -> target
/// @param surface_contour The contour description of the surface -> test
///
/// @return whether this should be associated
bool operator()(const std::vector<dpoint2D<algebra_t>> &bin_contour,
const std::vector<dpoint2D<algebra_t>> &surface_contour) {

using scalar_t = dscalar<algebra_t>;
using point2_t = dpoint2D<algebra_t>;

// Check if centre of gravity is inside bin
point2_t cgs = {0.f, 0.f};
for (const auto &svtx : surface_contour) {
cgs = cgs + svtx;
}
cgs = 1.f / static_cast<scalar>(surface_contour.size()) * cgs;
cgs = 1.f / static_cast<scalar_t>(surface_contour.size()) * cgs;

std::size_t i = 0u;
std::size_t j = 0u;
Expand All @@ -97,25 +103,27 @@ struct center_of_gravity_generic {
}
};

/** Check if the egdes of the bin and surface contour overlap */
/// Check if the egdes of the bin and surface contour overlap
template <typename algebra_t>
struct edges_intersect_generic {

/** Call operator to the struct, allows to chain several chain operators
* together
*
* @param bin_contour The contour description of the bin -> target
* @param surface_contour The contour description of the surface -> test
*
* @return whether this should be associated
*/
template <typename point2_t>
bool operator()(const std::vector<point2_t> &bin_contour,
const std::vector<point2_t> &surface_contour) {
/// Call operator to the struct, allows to chain several chain operators
/// together
///
/// @param bin_contour The contour description of the bin -> target
/// @param surface_contour The contour description of the surface -> test
///
/// @return whether this should be associated
bool operator()(const std::vector<dpoint2D<algebra_t>> &bin_contour,
const std::vector<dpoint2D<algebra_t>> &surface_contour) {

using scalar_t = dscalar<algebra_t>;
using point2_t = dpoint2D<algebra_t>;

auto intersect = [](const point2_t &pi, const point2_t &pj,
const point2_t &pk, const point2_t &pl) {
scalar d = (pj[0] - pi[0]) * (pl[1] - pk[1]) -
(pj[1] - pi[1]) * (pl[0] - pk[0]);
scalar_t d = (pj[0] - pi[0]) * (pl[1] - pk[1]) -
(pj[1] - pi[1]) * (pl[0] - pk[0]);

if (d != 0.f) {
double r = ((pi[1] - pk[1]) * (pl[0] - pk[0]) -
Expand Down
65 changes: 33 additions & 32 deletions core/include/detray/builders/detail/bin_association.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,13 +37,13 @@ namespace detray::detail {
/// taken absolute or relative
template <typename context_t, typename surface_container_t,
typename transform_container_t, typename mask_container_t,
concepts::surface_grid grid_t>
concepts::surface_grid grid_t, typename scalar_t>
static inline void bin_association(const context_t & /*context*/,
const surface_container_t &surfaces,
const transform_container_t &transforms,
const mask_container_t &surface_masks,
grid_t &grid,
const std::array<scalar, 2> &bin_tolerance,
const std::array<scalar_t, 2> &bin_tolerance,
bool absolute_tolerance = true) {

using algebra_t = typename grid_t::local_frame_type::algebra_type;
Expand All @@ -58,8 +58,8 @@ static inline void bin_association(const context_t & /*context*/,
polar2D<algebra_t>>) {
// Run with two different associators: center of gravity and edge
// intersection
center_of_gravity_generic cgs_assoc;
edges_intersect_generic edges_assoc;
center_of_gravity_generic<algebra_t> cgs_assoc;
edges_intersect_generic<algebra_t> edges_assoc;

// Loop over all bins and associate the surfaces
for (unsigned int bin_0 = 0; bin_0 < axis_0.nbins(); ++bin_0) {
Expand All @@ -68,18 +68,18 @@ static inline void bin_association(const context_t & /*context*/,
auto r_borders = axis_0.bin_edges(bin_0);
auto phi_borders = axis_1.bin_edges(bin_1);

scalar r_add =
scalar_t r_add =
absolute_tolerance
? bin_tolerance[0]
: bin_tolerance[0] * (r_borders[1] - r_borders[0]);
scalar phi_add =
scalar_t phi_add =
absolute_tolerance
? bin_tolerance[1]
: bin_tolerance[1] * (phi_borders[1] - phi_borders[0]);

// Create a contour for the bin
std::vector<point2_t> bin_contour =
detail::r_phi_polygon<scalar, point2_t>(
detail::r_phi_polygon<scalar_t, point2_t>(
r_borders[0] - r_add, r_borders[1] + r_add,
phi_borders[0] - phi_add, phi_borders[1] + phi_add);

Expand Down Expand Up @@ -125,8 +125,8 @@ static inline void bin_association(const context_t & /*context*/,
std::is_same_v<typename grid_t::local_frame_type,
concentric_cylindrical2D<algebra_t>>) {

center_of_gravity_rectangle cgs_assoc;
edges_intersect_generic edges_assoc;
center_of_gravity_rectangle<algebra_t> cgs_assoc;
edges_intersect_generic<algebra_t> edges_assoc;

// Loop over all bins and associate the surfaces
for (unsigned int bin_0 = 0; bin_0 < axis_0.nbins(); ++bin_0) {
Expand All @@ -135,19 +135,19 @@ static inline void bin_association(const context_t & /*context*/,
auto z_borders = axis_0.bin_edges(bin_0);
auto phi_borders = axis_1.bin_edges(bin_1);

scalar z_add =
scalar_t z_add =
absolute_tolerance
? bin_tolerance[0]
: bin_tolerance[0] * (z_borders[1] - z_borders[0]);
scalar phi_add =
scalar_t phi_add =
absolute_tolerance
? bin_tolerance[1]
: bin_tolerance[1] * (phi_borders[1] - phi_borders[0]);

scalar z_min = z_borders[0];
scalar z_max = z_borders[1];
scalar phi_min_rep = phi_borders[0];
scalar phi_max_rep = phi_borders[1];
scalar_t z_min = z_borders[0];
scalar_t z_max = z_borders[1];
scalar_t phi_min_rep = phi_borders[0];
scalar_t phi_max_rep = phi_borders[1];

point2_t p0_bin = {z_min - z_add, phi_min_rep - phi_add};
point2_t p1_bin = {z_min - z_add, phi_max_rep + phi_add};
Expand Down Expand Up @@ -177,25 +177,26 @@ static inline void bin_association(const context_t & /*context*/,
// Create a surface contour
std::vector<point2_t> surface_contour;
surface_contour.reserve(vertices.size());
scalar phi_min = std::numeric_limits<scalar>::max();
scalar phi_max =
-std::numeric_limits<scalar>::max();
scalar_t phi_min =
std::numeric_limits<scalar_t>::max();
scalar_t phi_max =
-std::numeric_limits<scalar_t>::max();
// We poentially need the split vertices
std::vector<point2_t> s_c_neg;
std::vector<point2_t> s_c_pos;
scalar z_min_neg =
std::numeric_limits<scalar>::max();
scalar z_max_neg =
-std::numeric_limits<scalar>::max();
scalar z_min_pos =
std::numeric_limits<scalar>::max();
scalar z_max_pos =
-std::numeric_limits<scalar>::max();
scalar_t z_min_neg =
std::numeric_limits<scalar_t>::max();
scalar_t z_max_neg =
-std::numeric_limits<scalar_t>::max();
scalar_t z_min_pos =
std::numeric_limits<scalar_t>::max();
scalar_t z_max_pos =
-std::numeric_limits<scalar_t>::max();

for (const auto &v : vertices) {
const point3_t vg =
transform.point_to_global(v);
scalar phi = math::atan2(vg[1], vg[0]);
scalar_t phi = math::atan2(vg[1], vg[0]);
phi_min = math::min(phi, phi_min);
phi_max = math::max(phi, phi_max);
surface_contour.push_back({vg[2], phi});
Expand All @@ -211,16 +212,16 @@ static inline void bin_association(const context_t & /*context*/,
}
// Check for phi wrapping
std::vector<std::vector<point2_t>> surface_contours;
if (phi_max - phi_min > constant<scalar>::pi &&
if (phi_max - phi_min > constant<scalar_t>::pi &&
phi_max * phi_min < 0.) {
s_c_neg.push_back(
{z_max_neg, -constant<scalar>::pi});
{z_max_neg, -constant<scalar_t>::pi});
s_c_neg.push_back(
{z_min_neg, -constant<scalar>::pi});
{z_min_neg, -constant<scalar_t>::pi});
s_c_pos.push_back(
{z_max_pos, constant<scalar>::pi});
{z_max_pos, constant<scalar_t>::pi});
s_c_pos.push_back(
{z_min_pos, constant<scalar>::pi});
{z_min_pos, constant<scalar_t>::pi});
surface_contours = {s_c_neg, s_c_pos};
} else {
surface_contours = {surface_contour};
Expand Down
2 changes: 1 addition & 1 deletion core/include/detray/builders/detail/portal_accessor.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ namespace detray::detail {
template <typename detector_t>
auto get_cylinder_portals(const tracking_volume<detector_t> &vol) {

using scalar_t = typename detector_t::scalar_type;
using scalar_t = dscalar<typename detector_t::algebra_type>;

std::vector<const typename detector_t::volume_type &> inner_pt{},
outer_pt{}, lower_pt{}, upper_pr{};
Expand Down
2 changes: 1 addition & 1 deletion core/include/detray/builders/detail/volume_connector.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ template <typename detector_t,
void connect_cylindrical_volumes(
detector_t &d, const typename detector_t::volume_finder &volume_grid) {

using scalar_t = typename detector_t::scalar_type;
using scalar_t = dscalar<typename detector_t::algebra_type>;
typename detector_t::context default_context = {};

// The grid is populated, now create portal surfaces
Expand Down
Loading
Loading