actsvg/core/utils.hpp

0001 // This file is part of the actsvg package.
0002 //
0003 // Copyright (C) 2022 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include <cmath>
0012 #include <iomanip>
0013 #include <iostream>
0014 #include <sstream>
0015 #include <tuple>
0016
0017 namespace actsvg {
0018
0019 namespace utils {
0020
0021 /** Helper method to run enumerate(...) with structured binding
0022  * over STL type containers.
0023  *
0024  * @param iterable is a std-like iterable container type
0025  *
0026  **/
0027 template <typename container_type,
0028           typename container_type_iter =
0029               decltype(std::begin(std::declval<container_type>())),
0030           typename = decltype(std::end(std::declval<container_type>()))>
0031 constexpr auto enumerate(container_type &&iterable) {
0032     struct iterator {
0033         size_t i;
0034         container_type_iter iter;
0035
0036         bool operator!=(const iterator &rhs) const { return iter != rhs.iter; }
0037
0038         /** Increase index and iterator at once */
0039         void operator++() {
0040             ++i;
0041             ++iter;
0042         }
0043
0044         /** Tie them together for returning */
0045         auto operator*() const { return std::tie(i, *iter); }
0046     };
0047     struct iterable_wrapper {
0048         container_type iterable;
0049         auto begin() { return iterator{0, std::begin(iterable)}; }
0050         auto end() { return iterator{0, std::end(iterable)}; }
0051     };
0052     return iterable_wrapper{std::forward<container_type>(iterable)};
0053 }
0054
0055 /** Helper for perp
0056  *
0057  * @param p_ the point
0058  *
0059  * @return a scalar of the norm
0060  **/
0061 template <typename point2_type>
0062 scalar perp(const point2_type &p_) {
0063     return std::sqrt(p_[0] * p_[0] + p_[1] * p_[1]);
0064 }
0065
0066 /** Helper method for calculating the distance between two points
0067  *
0068  * @param x0_ the starting point
0069  * @param x1_ the ending point
0070  *
0071  * @return the distance
0072  **/
0073 template <typename point2_type>
0074 scalar distance(const point2_type &x0_, const point2_type &x1_) {
0075     return std::sqrt(std::pow(x1_[0] - x0_[0], 2) +
0076                      std::pow(x1_[1] - x0_[1], 2));
0077 }
0078
0079 /** Helper method to intersect two lines in 2-d
0080  * @param x0_ the starting point of line 0
0081  * @param d0_ the direction of line 0
0082  * @param x1_ the starting point of line 1
0083  * @param d1_ the direction of line 1
0084  *
0085  * @return the intersection point
0086  **/
0087 template <typename point2_type>
0088 point2_type intersect(const point2_type &x0_, const point2_type &d0_,
0089                       const point2_type &x1_, const point2_type &d1_) {
0090     point2_type p_int;
0091
0092     struct line {
0093         scalar a;
0094         scalar d;
0095
0096         line(const point2_type &x_, const point2_type &d_) {
0097             a = d_[1u] / d_[0u];
0098             d = x_[1u] - a * x_[0u];
0099         }
0100     };
0101
0102     if (d0_[0u] == 0) {
0103
0104         line l1(x1_, d1_);
0105         p_int[0u] = x0_[0u];
0106         p_int[1u] = l1.a * p_int[0u] + l1.d;
0107
0108     } else if (d1_[0u] == 0) {
0109
0110         line l0(x0_, d0_);
0111         p_int[0u] = x1_[0u];
0112         p_int[1u] = l0.a * p_int[0u] + l0.d;
0113
0114     } else {
0115
0116         line l0(x0_, d0_);
0117         line l1(x1_, d1_);
0118
0119         p_int[0u] = (l1.d - l0.d) / (l0.a - l1.a);
0120         p_int[1u] = l0.a * p_int[0u] + l0.d;
0121     }
0122
0123     return p_int;
0124 }
0125
0126 /** Helper from id to url
0127  * @param id_ the idnetification to be transformed
0128  **/
0129 std::string id_to_url(const std::string &id_);
0130
0131 /** Helper to format point2
0132  *
0133  * @param s_ the scalar
0134  * @param pr_ the precision
0135  *
0136  * @return a string
0137  */
0138 std::string to_string(const scalar &s_, size_t pr_ = 4);
0139
0140 /** Helper to format point2
0141  *
0142  * @param p_ the point
0143  * @param pr_ the precision
0144  *
0145  * @return a string
0146  */
0147 template <typename point2_type>
0148 std::string to_string(const point2_type &p_, size_t pr_ = 4) {
0149     std::stringstream sstream;
0150     sstream << std::setprecision(pr_) << "(" << p_[0] << "," << p_[1] << ")";
0151     return sstream.str();
0152 }
0153
0154 /** Check if a parameter is within range
0155  *
0156  * @param p_ parameter to be checked
0157  * @param range_ the range
0158  *
0159  * @note borders are included in this check
0160  **/
0161 template <typename range_type>
0162 bool inside_range(scalar p_, const range_type &range_) {
0163     return (range_[0] <= p_ and p_ <= range_[1]);
0164 }
0165
0166 /** Helper method to calculate the barycenter
0167  *
0168  * @param vs_ vertices
0169  *
0170  * @return a new barycenter
0171  **/
0172 template <typename point2_type>
0173 point2_type barycenter(const std::vector<point2_type> &vs_) {
0174     point2_type rv = {0., 0.};
0175     for (const auto &v : vs_) {
0176         rv[0] += v[0];
0177         rv[1] += v[1];
0178     }
0179     rv[0] /= vs_.size();
0180     rv[1] /= vs_.size();
0181     return rv;
0182 }
0183
0184 /** Helper method to rotate a 2-d point
0185  * @param p_ the point to be rotated
0186  * @param a_ the angle alpha
0187  *
0188  * @return the rotated point
0189  **/
0190 template <typename point2_type>
0191 point2_type rotate(const point2_type &p_, scalar a_) {
0192     point2_type p_rot;
0193     p_rot[0] = std::cos(a_) * p_[0] - std::sin(a_) * p_[1];
0194     p_rot[1] = std::sin(a_) * p_[0] + std::cos(a_) * p_[1];
0195     return p_rot;
0196 }
0197
0198 /** Helper method to scale object
0199  *
0200  * @param p0_ the point object and @param s_ the scale
0201  *
0202  * @return a new point_type object
0203  **/
0204 template <typename point_type, std::size_t kDIM = 3u>
0205 point_type scale(const point_type &p0_, scalar s_) {
0206     point_type scaled;
0207     scaled[0] = s_ * p0_[0];
0208     scaled[1] = s_ * p0_[1];
0209     if constexpr (kDIM == 3u) {
0210         scaled[2] = s_ * p0_[2];
0211     }
0212     return scaled;
0213 }
0214
0215 /** Helper method for making a unit vector
0216  *
0217  * @param x0_ the starting point
0218  * @param x1_ the ending point
0219  *
0220  * @return the distance
0221  **/
0222 template <typename point_type, std::size_t kDIM = 2u>
0223 point_type unit_direction(const point_type &x0_, const point_type &x1_) {
0224     // 2 dimensional case
0225     if constexpr (kDIM == 2u) {
0226         scalar length = std::sqrt(std::pow(x1_[0] - x0_[0], 2) +
0227                                   std::pow(x1_[1] - x0_[1], 2));
0228         scalar norml = 1. / length;
0229         return point_type{norml * (x1_[0] - x0_[0]), norml * (x1_[1] - x0_[1])};
0230     }
0231     // 3 dimensional case
0232     if constexpr (kDIM == 3u) {
0233         scalar length = std::sqrt(std::pow(x1_[0] - x0_[0], 2) +
0234                                   std::pow(x1_[1] - x0_[1], 2) +
0235                                   std::pow(x1_[2] - x0_[2], 2));
0236         scalar norml = 1. / length;
0237         return point_type{norml * (x1_[0] - x0_[0]), norml * (x1_[1] - x0_[1]),
0238                           norml * (x1_[2] - x0_[2])};
0239     }
0240     return point_type{};
0241 }
0242
0243 /** Helper method to add two point3 objects
0244  *
0245  * @param p0_ and @param p1_ the  points
0246  *
0247  * @return a new point_type object
0248  **/
0249 template <typename point_type, std::size_t kDIM = 3u>
0250 point_type add(const point_type &p0_, const point_type &p1_) {
0251     point_type added;
0252     added[0] = p0_[0] + p1_[0];
0253     added[1] = p0_[1] + p1_[1];
0254     if constexpr (kDIM == 3u) {
0255         added[2] = p0_[2] + p1_[2];
0256     }
0257     return added;
0258 }
0259
0260 /** Helper method to rotate a 3D point @param p_ under
0261  * rotation @param rt_
0262  *
0263  * @return new rotated point3
0264  **/
0265 template <typename point3_type>
0266 point3_type rotate(const std::array<point3_type, 3> &rt_,
0267                    const point3_type &p_) {
0268
0269     point3_type rotated;
0270     rotated[0] = rt_[0][0] * p_[0] + rt_[0][1] * p_[1] + rt_[0][2] * p_[2];
0271     rotated[1] = rt_[1][0] * p_[0] + rt_[1][1] * p_[1] + rt_[1][2] * p_[2];
0272     rotated[2] = rt_[2][0] * p_[0] + rt_[2][1] * p_[1] + rt_[2][2] * p_[2];
0273     return rotated;
0274 }
0275
0276 /** Helper method to place a 3D point @param p_ under
0277  * @param tr_ transform and @param rt_ rotation
0278  *
0279  * @return new transform
0280  **/
0281 template <typename point3_type>
0282 point3_type place(const point3_type &p_, const point3_type &tr_,
0283                   const std::array<point3_type, 3> &rt_) {
0284     point3_type placed = rotate(rt_, p_);
0285     placed[0] = p_[0] + tr_[0];
0286     placed[1] = p_[1] + tr_[1];
0287     placed[2] = p_[2] + tr_[2];
0288     return placed;
0289 }
0290
0291 /** Helper method to place a 3D points @param pc_ under
0292  * @param tr_ transform and @param rt_ rotation
0293  *
0294  * @return new transformed point3 collection
0295  **/
0296 template <typename point3_container, typename point3_type>
0297 point3_container place_vertices(const point3_container &pc_,
0298                                 const point3_type &tr_,
0299                                 const std::array<point3_type, 3> &rt_) {
0300     point3_container placed;
0301     for (const auto &p : pc_) {
0302         placed.push_back(place(p, tr_, rt_));
0303     }
0304     return placed;
0305 }
0306
0307 }  // namespace utils
0308
0309 }  // namespace actsvg