Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RVec.hxx
Go to the documentation of this file.
1// Author: Enrico Guiraud, Enric Tejedor, Danilo Piparo CERN 01/2018
2
3/*************************************************************************
4 * Copyright (C) 1995-2018, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11/**
12 \defgroup vecops VecOps
13*/
14
15#ifndef ROOT_RVEC
16#define ROOT_RVEC
17
18#ifdef _WIN32
19 #ifndef M_PI
20 #ifndef _USE_MATH_DEFINES
21 #define _USE_MATH_DEFINES
22 #endif
23 #include <math.h>
24 #undef _USE_MATH_DEFINES
25 #endif
26 #define _VECOPS_USE_EXTERN_TEMPLATES false
27#else
28 #define _VECOPS_USE_EXTERN_TEMPLATES true
29#endif
30
33#include <ROOT/RStringView.hxx>
34#include <TError.h> // R__ASSERT
35#include <ROOT/TypeTraits.hxx>
36
37#include <algorithm>
38#include <cmath>
39#include <numeric> // for inner_product
40#include <sstream>
41#include <stdexcept>
42#include <type_traits>
43#include <vector>
44#include <utility>
45#include <tuple>
46
47// TODO remove when fData is removed (6.26)
48#include <map>
49#include <string>
50
51#ifdef R__HAS_VDT
52#include <vdt/vdtMath.h>
53#endif
54
55
56// TODO remove when fData is removed (6.26)
57class TTree;
58class TBranch;
59
60namespace ROOT {
61
62namespace VecOps {
63template<typename T>
64class RVec;
65}
66
67// TODO remove when fData is removed (6.26)
68namespace Internal {
69namespace RDF {
70class BoolArray;
71class RBranchSet;
72using BoolArrayMap = std::map<std::string, BoolArray>;
73template <typename T>
74void SetBranchesHelper(BoolArrayMap &boolArrays, TTree *inputTree, TTree &outputTree, const std::string &inName,
75 const std::string &outName, TBranch *&branch, void *&branchAddress, ROOT::VecOps::RVec<T> *ab,
76 RBranchSet &outputBranches);
77}
78}
79
80namespace Detail {
81namespace VecOps {
82
83template<typename T>
85
86template <typename... T>
87std::size_t GetVectorsSize(std::string_view id, const RVec<T> &... vs)
88{
89 constexpr const auto nArgs = sizeof...(T);
90 const std::size_t sizes[] = {vs.size()...};
91 if (nArgs > 1) {
92 for (auto i = 1UL; i < nArgs; i++) {
93 if (sizes[0] == sizes[i])
94 continue;
95 std::string msg(id);
96 msg += ": input RVec instances have different lengths!";
97 throw std::runtime_error(msg);
98 }
99 }
100 return sizes[0];
101}
102
103template <typename F, typename... T>
104auto MapImpl(F &&f, const RVec<T> &... vs) -> RVec<decltype(f(vs[0]...))>
105{
106 const auto size = GetVectorsSize("Map", vs...);
107 RVec<decltype(f(vs[0]...))> ret(size);
108
109 for (auto i = 0UL; i < size; i++)
110 ret[i] = f(vs[i]...);
111
112 return ret;
113}
114
115template <typename Tuple_t, std::size_t... Is>
116auto MapFromTuple(Tuple_t &&t, std::index_sequence<Is...>)
117 -> decltype(MapImpl(std::get<std::tuple_size<Tuple_t>::value - 1>(t), std::get<Is>(t)...))
118{
119 constexpr const auto tupleSizeM1 = std::tuple_size<Tuple_t>::value - 1;
120 return MapImpl(std::get<tupleSizeM1>(t), std::get<Is>(t)...);
121}
122
123}
124}
125
126namespace Internal {
127namespace VecOps {
128
129// We use this helper to workaround a limitation of compilers such as
130// gcc 4.8 amd clang on osx 10.14 for which std::vector<bool>::emplace_back
131// is not defined.
132template <typename T, typename... Args>
133void EmplaceBack(T &v, Args &&... args)
134{
135 v.emplace_back(std::forward<Args>(args)...);
136}
137
138template <typename... Args>
139void EmplaceBack(std::vector<bool> &v, Args &&... args)
140{
141 v.push_back(std::forward<Args>(args)...);
142}
143
144} // End of VecOps NS
145} // End of Internal NS
146
147namespace VecOps {
148// clang-format off
149/**
150\class ROOT::VecOps::RVec
151\ingroup vecops
152\brief A "std::vector"-like collection of values implementing handy operation to analyse them
153\tparam T The type of the contained objects
154
155A RVec is a container designed to make analysis of values' collections fast and easy.
156Its storage is contiguous in memory and its interface is designed such to resemble to the one
157of the stl vector. In addition the interface features methods and external functions to ease
158the manipulation and analysis of the data in the RVec.
159
160\htmlonly
161<a href="https://doi.org/10.5281/zenodo.1253756"><img src="https://zenodo.org/badge/DOI/10.5281/zenodo.1253756.svg" alt="DOI"></a>
162\endhtmlonly
163
164## Table of Contents
165- [Example](#example)
166- [Owning and adopting memory](#owningandadoptingmemory)
167- [Sorting and manipulation of indices](#sorting)
168- [Usage in combination with RDataFrame](#usagetdataframe)
169- [Reference for the RVec class](#RVecdoxyref)
170
171Also see the [reference for RVec helper functions](https://root.cern/doc/master/namespaceROOT_1_1VecOps.html).
172
173## <a name="example"></a>Example
174Suppose to have an event featuring a collection of muons with a certain pseudorapidity,
175momentum and charge, e.g.:
176~~~{.cpp}
177std::vector<short> mu_charge {1, 1, -1, -1, -1, 1, 1, -1};
178std::vector<float> mu_pt {56, 45, 32, 24, 12, 8, 7, 6.2};
179std::vector<float> mu_eta {3.1, -.2, -1.1, 1, 4.1, 1.6, 2.4, -.5};
180~~~
181Suppose you want to extract the transverse momenta of the muons satisfying certain
182criteria, for example consider only negatively charged muons with a pseudorapidity
183smaller or equal to 2 and with a transverse momentum greater than 10 GeV.
184Such a selection would require, among the other things, the management of an explicit
185loop, for example:
186~~~{.cpp}
187std::vector<float> goodMuons_pt;
188const auto size = mu_charge.size();
189for (size_t i=0; i < size; ++i) {
190 if (mu_pt[i] > 10 && abs(mu_eta[i]) <= 2. && mu_charge[i] == -1) {
191 goodMuons_pt.emplace_back(mu_pt[i]);
192 }
193}
194~~~
195These operations become straightforward with RVec - we just need to *write what
196we mean*:
197~~~{.cpp}
198auto goodMuons_pt = mu_pt[ (mu_pt > 10.f && abs(mu_eta) <= 2.f && mu_charge == -1) ]
199~~~
200Now the clean collection of transverse momenta can be used within the rest of the data analysis, for
201example to fill a histogram.
202
203## <a name="owningandadoptingmemory"></a>Owning and adopting memory
204RVec has contiguous memory associated to it. It can own it or simply adopt it. In the latter case,
205it can be constructed with the address of the memory associated to it and its length. For example:
206~~~{.cpp}
207std::vector<int> myStlVec {1,2,3};
208RVec<int> myRVec(myStlVec.data(), myStlVec.size());
209~~~
210In this case, the memory associated to myStlVec and myRVec is the same, myRVec simply "adopted it".
211If any method which implies a re-allocation is called, e.g. *emplace_back* or *resize*, the adopted
212memory is released and new one is allocated. The previous content is copied in the new memory and
213preserved.
214
215## <a name="#sorting"></a>Sorting and manipulation of indices
216
217### Sorting
218RVec complies to the STL interfaces when it comes to iterations. As a result, standard algorithms
219can be used, for example sorting:
220~~~{.cpp}
221RVec<double> v{6., 4., 5.};
222std::sort(v.begin(), v.end());
223~~~
224
225For convinience, helpers are provided too:
226~~~{.cpp}
227auto sorted_v = Sort(v);
228auto reversed_v = Reverse(v);
229~~~
230
231### Manipulation of indices
232
233It is also possible to manipulated the RVecs acting on their indices. For example,
234the following syntax
235~~~{.cpp}
236RVec<double> v0 {9., 7., 8.};
237auto v1 = Take(v0, {1, 2, 0});
238~~~
239will yield a new RVec<double> the content of which is the first, second and zeroth element of
240v0, i.e. `{7., 8., 9.}`.
241
242The `Argsort` helper extracts the indices which order the content of a `RVec`. For example,
243this snippet accomplish in a more expressive way what we just achieved:
244~~~{.cpp}
245auto v1_indices = Argsort(v0); // The content of v1_indices is {1, 2, 0}.
246v1 = Take(v0, v1_indices);
247~~~
248
249The `Take` utility allows to extract portions of the `RVec`. The content to be *taken*
250can be specified with an `RVec` of indices or an integer. If the integer is negative,
251elements will be picked starting from the end of the container:
252~~~{.cpp}
253RVec<float> vf {1.f, 2.f, 3.f, 4.f};
254auto vf_1 = Take(vf, {1, 3}); // The content is {2.f, 4.f}
255auto vf_2 = Take(vf, 2); // The content is {1.f, 2.f}
256auto vf_3 = Take(vf, -3); // The content is {2.f, 3.f, 4.f}
257~~~
258
259## <a name="usagetdataframe"></a>Usage in combination with RDataFrame
260RDataFrame leverages internally RVecs. Suppose to have a dataset stored in a
261TTree which holds these columns (here we choose C arrays to represent the
262collections, they could be as well std::vector instances):
263~~~{.bash}
264 nPart "nPart/I" An integer representing the number of particles
265 px "px[nPart]/D" The C array of the particles' x component of the momentum
266 py "py[nPart]/D" The C array of the particles' y component of the momentum
267 E "E[nPart]/D" The C array of the particles' Energy
268~~~
269Suppose you'd like to plot in a histogram the transverse momenta of all particles
270for which the energy is greater than 200 MeV.
271The code required would just be:
272~~~{.cpp}
273RDataFrame d("mytree", "myfile.root");
274using doubles = RVec<double>;
275auto cutPt = [](doubles &pxs, doubles &pys, doubles &Es) {
276 auto all_pts = sqrt(pxs * pxs + pys * pys);
277 auto good_pts = all_pts[Es > 200.];
278 return good_pts;
279 };
280
281auto hpt = d.Define("pt", cutPt, {"px", "py", "E"})
282 .Histo1D("pt");
283hpt->Draw();
284~~~
285And if you'd like to express your selection as a string:
286~~~{.cpp}
287RDataFrame d("mytree", "myfile.root");
288auto hpt = d.Define("pt", "sqrt(pxs * pxs + pys * pys)[E>200]")
289 .Histo1D("pt");
290hpt->Draw();
291~~~
292<a name="RVecdoxyref"></a>
293**/
294// clang-format on
295template <typename T>
296class RVec {
297 // TODO remove this friendship in 6.26, when fData is removed
298 friend void ::ROOT::Internal::RDF::SetBranchesHelper<T>(ROOT::Internal::RDF::BoolArrayMap &boolArrays,
299 TTree *inputTree, TTree &outputTree,
300 const std::string &inName, const std::string &outName,
301 TBranch *&branch, void *&branchAddress, RVec<T> *ab,
302 ROOT::Internal::RDF::RBranchSet &outputBranches);
303
304 // Here we check if T is a bool. This is done in order to decide what type
305 // to use as a storage. If T is anything but bool, we use a vector<T, RAdoptAllocator<T>>.
306 // If T is a bool, we opt for a plain vector<bool> otherwise we'll not be able
307 // to write the data type given the shortcomings of TCollectionProxy design.
308 static constexpr const auto IsVecBool = std::is_same<bool, T>::value;
309public:
310 using Impl_t = typename std::conditional<IsVecBool, std::vector<bool>, std::vector<T, ::ROOT::Detail::VecOps::RAdoptAllocator<T>>>::type;
311 using value_type = typename Impl_t::value_type;
312 using size_type = typename Impl_t::size_type;
313 using difference_type = typename Impl_t::difference_type;
314 using reference = typename Impl_t::reference;
315 using const_reference = typename Impl_t::const_reference;
316 using pointer = typename Impl_t::pointer;
317 using const_pointer = typename Impl_t::const_pointer;
318 // The data_t and const_data_t types are chosen to be void in case T is a bool.
319 // This way we can as elegantly as in the STL return void upon calling the data() method.
320 using data_t = typename std::conditional<IsVecBool, void, typename Impl_t::pointer>::type;
321 using const_data_t = typename std::conditional<IsVecBool, void, typename Impl_t::const_pointer>::type;
322 using iterator = typename Impl_t::iterator;
323 using const_iterator = typename Impl_t::const_iterator;
324 using reverse_iterator = typename Impl_t::reverse_iterator;
325 using const_reverse_iterator = typename Impl_t::const_reverse_iterator;
326
327private:
329
330public:
331 // constructors
332 RVec() {}
333
334 explicit RVec(size_type count) : fData(count) {}
335
336 RVec(size_type count, const T &value) : fData(count, value) {}
337
338 RVec(const RVec<T> &v) : fData(v.fData) {}
339
340 RVec(RVec<T> &&v) : fData(std::move(v.fData)) {}
341
342 RVec(const std::vector<T> &v) : fData(v.cbegin(), v.cend()) {}
343
344 RVec(pointer p, size_type n) : fData(n, T(), ROOT::Detail::VecOps::RAdoptAllocator<T>(p)) {}
345
346 template <class InputIt>
347 RVec(InputIt first, InputIt last) : fData(first, last) {}
348
349 RVec(std::initializer_list<T> init) : fData(init) {}
350
351 // assignment
353 {
354 fData = v.fData;
355 return *this;
356 }
357
359 {
360 std::swap(fData, v.fData);
361 return *this;
362 }
363
364 RVec<T> &operator=(std::initializer_list<T> ilist)
365 {
366 fData = ilist;
367 return *this;
368 }
369
370 // conversion
371 template <typename U, typename = std::enable_if<std::is_convertible<T, U>::value>>
372 operator RVec<U>() const
373 {
374 RVec<U> ret(size());
375 std::copy(begin(), end(), ret.begin());
376 return ret;
377 }
378
379 R__DEPRECATED(6, 26, "Please use `std::vector<T>(rvec.begin(), rvec.end())` instead.")
380 const Impl_t &AsVector() const { return fData; }
381 R__DEPRECATED(6, 26, "Please use `std::vector<T>(rvec.begin(), rvec.end())` instead.") Impl_t &AsVector()
382 {
383 return fData;
384 }
385
386 // accessors
387 reference at(size_type pos) { return fData.at(pos); }
388 const_reference at(size_type pos) const { return fData.at(pos); }
389 /// No exception thrown. The user specifies the desired value in case the RVec is shorter than `pos`.
390 value_type at(size_type pos, value_type fallback) { return pos < fData.size() ? fData[pos] : fallback; }
391 /// No exception thrown. The user specifies the desired value in case the RVec is shorter than `pos`.
392 value_type at(size_type pos, value_type fallback) const { return pos < fData.size() ? fData[pos] : fallback; }
393 reference operator[](size_type pos) { return fData[pos]; }
394 const_reference operator[](size_type pos) const { return fData[pos]; }
395
396 template <typename V, typename = std::enable_if<std::is_convertible<V, bool>::value>>
397 RVec operator[](const RVec<V> &conds) const
398 {
399 const size_type n = conds.size();
400
401 if (n != size())
402 throw std::runtime_error("Cannot index RVec with condition vector of different size");
403
404 RVec<T> ret;
405 ret.reserve(n);
406 for (size_type i = 0; i < n; ++i)
407 if (conds[i])
408 ret.emplace_back(fData[i]);
409 return ret;
410 }
411
412 reference front() { return fData.front(); }
413 const_reference front() const { return fData.front(); }
414 reference back() { return fData.back(); }
415 const_reference back() const { return fData.back(); }
416 data_t data() noexcept { return fData.data(); }
417 const_data_t data() const noexcept { return fData.data(); }
418 // iterators
419 iterator begin() noexcept { return fData.begin(); }
420 const_iterator begin() const noexcept { return fData.begin(); }
421 const_iterator cbegin() const noexcept { return fData.cbegin(); }
422 iterator end() noexcept { return fData.end(); }
423 const_iterator end() const noexcept { return fData.end(); }
424 const_iterator cend() const noexcept { return fData.cend(); }
425 reverse_iterator rbegin() noexcept { return fData.rbegin(); }
426 const_reverse_iterator rbegin() const noexcept { return fData.rbegin(); }
427 const_reverse_iterator crbegin() const noexcept { return fData.crbegin(); }
428 reverse_iterator rend() noexcept { return fData.rend(); }
429 const_reverse_iterator rend() const noexcept { return fData.rend(); }
430 const_reverse_iterator crend() const noexcept { return fData.crend(); }
431 // capacity
432 bool empty() const noexcept { return fData.empty(); }
433 size_type size() const noexcept { return fData.size(); }
434 size_type max_size() const noexcept { return fData.size(); }
435 void reserve(size_type new_cap) { fData.reserve(new_cap); }
436 size_type capacity() const noexcept { return fData.capacity(); }
437 void shrink_to_fit() { fData.shrink_to_fit(); };
438 // modifiers
439 void clear() noexcept { fData.clear(); }
440 iterator erase(iterator pos) { return fData.erase(pos); }
441 iterator erase(iterator first, iterator last) { return fData.erase(first, last); }
442 void push_back(T &&value) { fData.push_back(std::forward<T>(value)); }
443 void push_back(const value_type& value) { fData.push_back(value); };
444 template <class... Args>
445 reference emplace_back(Args &&... args)
446 {
447 ROOT::Internal::VecOps::EmplaceBack(fData, std::forward<Args>(args)...);
448 return fData.back();
449 }
450 /// This method is intended only for arithmetic types unlike the std::vector
451 /// corresponding one which is generic.
452 template<typename U = T, typename std::enable_if<std::is_arithmetic<U>::value, int>* = nullptr>
454 {
455 return fData.emplace(pos, value);
456 }
457 void pop_back() { fData.pop_back(); }
458 void resize(size_type count) { fData.resize(count); }
459 void resize(size_type count, const value_type &value) { fData.resize(count, value); }
460 void swap(RVec<T> &other) { std::swap(fData, other.fData); }
461};
462
463///@name RVec Unary Arithmetic Operators
464///@{
465
466#define RVEC_UNARY_OPERATOR(OP) \
467template <typename T> \
468RVec<T> operator OP(const RVec<T> &v) \
469{ \
470 RVec<T> ret(v); \
471 for (auto &x : ret) \
472 x = OP x; \
473return ret; \
474} \
475
480#undef RVEC_UNARY_OPERATOR
481
482///@}
483///@name RVec Binary Arithmetic Operators
484///@{
485
486#define ERROR_MESSAGE(OP) \
487 "Cannot call operator " #OP " on vectors of different sizes."
488
489#define RVEC_BINARY_OPERATOR(OP) \
490template <typename T0, typename T1> \
491auto operator OP(const RVec<T0> &v, const T1 &y) \
492 -> RVec<decltype(v[0] OP y)> \
493{ \
494 RVec<decltype(v[0] OP y)> ret(v.size()); \
495 auto op = [&y](const T0 &x) { return x OP y; }; \
496 std::transform(v.begin(), v.end(), ret.begin(), op); \
497 return ret; \
498} \
499 \
500template <typename T0, typename T1> \
501auto operator OP(const T0 &x, const RVec<T1> &v) \
502 -> RVec<decltype(x OP v[0])> \
503{ \
504 RVec<decltype(x OP v[0])> ret(v.size()); \
505 auto op = [&x](const T1 &y) { return x OP y; }; \
506 std::transform(v.begin(), v.end(), ret.begin(), op); \
507 return ret; \
508} \
509 \
510template <typename T0, typename T1> \
511auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
512 -> RVec<decltype(v0[0] OP v1[0])> \
513{ \
514 if (v0.size() != v1.size()) \
515 throw std::runtime_error(ERROR_MESSAGE(OP)); \
516 \
517 RVec<decltype(v0[0] OP v1[0])> ret(v0.size()); \
518 auto op = [](const T0 &x, const T1 &y) { return x OP y; }; \
519 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
520 return ret; \
521} \
522
531#undef RVEC_BINARY_OPERATOR
532
533///@}
534///@name RVec Assignment Arithmetic Operators
535///@{
536
537#define RVEC_ASSIGNMENT_OPERATOR(OP) \
538template <typename T0, typename T1> \
539RVec<T0>& operator OP(RVec<T0> &v, const T1 &y) \
540{ \
541 auto op = [&y](T0 &x) { return x OP y; }; \
542 std::transform(v.begin(), v.end(), v.begin(), op); \
543 return v; \
544} \
545 \
546template <typename T0, typename T1> \
547RVec<T0>& operator OP(RVec<T0> &v0, const RVec<T1> &v1) \
548{ \
549 if (v0.size() != v1.size()) \
550 throw std::runtime_error(ERROR_MESSAGE(OP)); \
551 \
552 auto op = [](T0 &x, const T1 &y) { return x OP y; }; \
553 std::transform(v0.begin(), v0.end(), v1.begin(), v0.begin(), op); \
554 return v0; \
555} \
556
567#undef RVEC_ASSIGNMENT_OPERATOR
568
569///@}
570///@name RVec Comparison and Logical Operators
571///@{
572
573#define RVEC_LOGICAL_OPERATOR(OP) \
574template <typename T0, typename T1> \
575auto operator OP(const RVec<T0> &v, const T1 &y) \
576 -> RVec<int> /* avoid std::vector<bool> */ \
577{ \
578 RVec<int> ret(v.size()); \
579 auto op = [y](const T0 &x) -> int { return x OP y; }; \
580 std::transform(v.begin(), v.end(), ret.begin(), op); \
581 return ret; \
582} \
583 \
584template <typename T0, typename T1> \
585auto operator OP(const T0 &x, const RVec<T1> &v) \
586 -> RVec<int> /* avoid std::vector<bool> */ \
587{ \
588 RVec<int> ret(v.size()); \
589 auto op = [x](const T1 &y) -> int { return x OP y; }; \
590 std::transform(v.begin(), v.end(), ret.begin(), op); \
591 return ret; \
592} \
593 \
594template <typename T0, typename T1> \
595auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
596 -> RVec<int> /* avoid std::vector<bool> */ \
597{ \
598 if (v0.size() != v1.size()) \
599 throw std::runtime_error(ERROR_MESSAGE(OP)); \
600 \
601 RVec<int> ret(v0.size()); \
602 auto op = [](const T0 &x, const T1 &y) -> int { return x OP y; }; \
603 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
604 return ret; \
605} \
606
615#undef RVEC_LOGICAL_OPERATOR
616
617///@}
618///@name RVec Standard Mathematical Functions
619///@{
620
621/// \cond
622template <typename T> struct PromoteTypeImpl;
623
624template <> struct PromoteTypeImpl<float> { using Type = float; };
625template <> struct PromoteTypeImpl<double> { using Type = double; };
626template <> struct PromoteTypeImpl<long double> { using Type = long double; };
627
628template <typename T> struct PromoteTypeImpl { using Type = double; };
629
630template <typename T>
631using PromoteType = typename PromoteTypeImpl<T>::Type;
632
633template <typename U, typename V>
634using PromoteTypes = decltype(PromoteType<U>() + PromoteType<V>());
635
636/// \endcond
637
638#define RVEC_UNARY_FUNCTION(NAME, FUNC) \
639 template <typename T> \
640 RVec<PromoteType<T>> NAME(const RVec<T> &v) \
641 { \
642 RVec<PromoteType<T>> ret(v.size()); \
643 auto f = [](const T &x) { return FUNC(x); }; \
644 std::transform(v.begin(), v.end(), ret.begin(), f); \
645 return ret; \
646 }
647
648#define RVEC_BINARY_FUNCTION(NAME, FUNC) \
649 template <typename T0, typename T1> \
650 RVec<PromoteTypes<T0, T1>> NAME(const T0 &x, const RVec<T1> &v) \
651 { \
652 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
653 auto f = [&x](const T1 &y) { return FUNC(x, y); }; \
654 std::transform(v.begin(), v.end(), ret.begin(), f); \
655 return ret; \
656 } \
657 \
658 template <typename T0, typename T1> \
659 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v, const T1 &y) \
660 { \
661 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
662 auto f = [&y](const T1 &x) { return FUNC(x, y); }; \
663 std::transform(v.begin(), v.end(), ret.begin(), f); \
664 return ret; \
665 } \
666 \
667 template <typename T0, typename T1> \
668 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v0, const RVec<T1> &v1) \
669 { \
670 if (v0.size() != v1.size()) \
671 throw std::runtime_error(ERROR_MESSAGE(NAME)); \
672 \
673 RVec<PromoteTypes<T0, T1>> ret(v0.size()); \
674 auto f = [](const T0 &x, const T1 &y) { return FUNC(x, y); }; \
675 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), f); \
676 return ret; \
677 } \
678
679#define RVEC_STD_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, std::F)
680#define RVEC_STD_BINARY_FUNCTION(F) RVEC_BINARY_FUNCTION(F, std::F)
681
686
690
695
700
708
715
722
727#undef RVEC_STD_UNARY_FUNCTION
728
729///@}
730///@name RVec Fast Mathematical Functions with Vdt
731///@{
732
733#ifdef R__HAS_VDT
734#define RVEC_VDT_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, vdt::F)
735
736RVEC_VDT_UNARY_FUNCTION(fast_expf)
737RVEC_VDT_UNARY_FUNCTION(fast_logf)
738RVEC_VDT_UNARY_FUNCTION(fast_sinf)
739RVEC_VDT_UNARY_FUNCTION(fast_cosf)
740RVEC_VDT_UNARY_FUNCTION(fast_tanf)
741RVEC_VDT_UNARY_FUNCTION(fast_asinf)
742RVEC_VDT_UNARY_FUNCTION(fast_acosf)
743RVEC_VDT_UNARY_FUNCTION(fast_atanf)
744
745RVEC_VDT_UNARY_FUNCTION(fast_exp)
746RVEC_VDT_UNARY_FUNCTION(fast_log)
747RVEC_VDT_UNARY_FUNCTION(fast_sin)
748RVEC_VDT_UNARY_FUNCTION(fast_cos)
749RVEC_VDT_UNARY_FUNCTION(fast_tan)
750RVEC_VDT_UNARY_FUNCTION(fast_asin)
751RVEC_VDT_UNARY_FUNCTION(fast_acos)
752RVEC_VDT_UNARY_FUNCTION(fast_atan)
753#undef RVEC_VDT_UNARY_FUNCTION
754
755#endif // R__HAS_VDT
756
757#undef RVEC_UNARY_FUNCTION
758
759///@}
760
761/// Inner product
762///
763/// Example code, at the ROOT prompt:
764/// ~~~{.cpp}
765/// using namespace ROOT::VecOps;
766/// RVec<float> v1 {1., 2., 3.};
767/// RVec<float> v2 {4., 5., 6.};
768/// auto v1_dot_v2 = Dot(v1, v2);
769/// v1_dot_v2
770/// // (float) 32.f
771/// ~~~
772template <typename T, typename V>
773auto Dot(const RVec<T> &v0, const RVec<V> &v1) -> decltype(v0[0] * v1[0])
774{
775 if (v0.size() != v1.size())
776 throw std::runtime_error("Cannot compute inner product of vectors of different sizes");
777 return std::inner_product(v0.begin(), v0.end(), v1.begin(), decltype(v0[0] * v1[0])(0));
778}
779
780/// Sum elements of an RVec
781///
782/// Example code, at the ROOT prompt:
783/// ~~~{.cpp}
784/// using namespace ROOT::VecOps;
785/// RVec<float> v {1.f, 2.f, 3.f};
786/// auto v_sum = Sum(v);
787/// v_sum
788/// // (float) 6.f
789/// ~~~
790template <typename T>
791T Sum(const RVec<T> &v)
792{
793 return std::accumulate(v.begin(), v.end(), T(0));
794}
795
796/// Get the mean of the elements of an RVec
797///
798/// The return type is a double precision floating point number.
799/// Example code, at the ROOT prompt:
800/// ~~~{.cpp}
801/// using namespace ROOT::VecOps;
802/// RVec<float> v {1.f, 2.f, 4.f};
803/// auto v_mean = Mean(v);
804/// v_mean
805/// // (double) 2.3333333
806/// ~~~
807template <typename T>
808double Mean(const RVec<T> &v)
809{
810 if (v.empty()) return 0.;
811 return double(Sum(v)) / v.size();
812}
813
814/// Get the greatest element of an RVec
815///
816/// Example code, at the ROOT prompt:
817/// ~~~~{.cpp}
818/// using namespace ROOT::VecOps;
819/// RVec<float> v {1.f, 2.f, 4.f};
820/// auto v_max = Max(v)
821/// v_max
822/// (float) 4.f
823/// ~~~~
824template <typename T>
825T Max(const RVec<T> &v)
826{
827 return *std::max_element(v.begin(), v.end());
828}
829
830/// Get the smallest element of an RVec
831///
832/// Example code, at the ROOT prompt:
833/// ~~~~{.cpp}
834/// using namespace ROOT::VecOps;
835/// RVec<float> v {1.f, 2.f, 4.f};
836/// auto v_min = Min(v)
837/// v_min
838/// (float) 1.f
839/// ~~~~
840template <typename T>
841T Min(const RVec<T> &v)
842{
843 return *std::min_element(v.begin(), v.end());
844}
845
846/// Get the index of the greatest element of an RVec
847/// In case of multiple occurrences of the maximum values,
848/// the index corresponding to the first occurrence is returned.
849///
850/// Example code, at the ROOT prompt:
851/// ~~~~{.cpp}
852/// using namespace ROOT::VecOps;
853/// RVec<float> v {1.f, 2.f, 4.f};
854/// auto v_argmax = ArgMax(v);
855/// v_argmax
856/// // (int) 2
857/// ~~~~
858template <typename T>
859std::size_t ArgMax(const RVec<T> &v)
860{
861 return std::distance(v.begin(), std::max_element(v.begin(), v.end()));
862}
863
864/// Get the index of the smallest element of an RVec
865/// In case of multiple occurrences of the minimum values,
866/// the index corresponding to the first occurrence is returned.
867///
868/// Example code, at the ROOT prompt:
869/// ~~~~{.cpp}
870/// using namespace ROOT::VecOps;
871/// RVec<float> v {1.f, 2.f, 4.f};
872/// auto v_argmin = ArgMin(v);
873/// v_argmin
874/// // (int) 0
875/// ~~~~
876template <typename T>
877std::size_t ArgMin(const RVec<T> &v)
878{
879 return std::distance(v.begin(), std::min_element(v.begin(), v.end()));
880}
881
882/// Get the variance of the elements of an RVec
883///
884/// The return type is a double precision floating point number.
885/// Example code, at the ROOT prompt:
886/// ~~~{.cpp}
887/// using namespace ROOT::VecOps;
888/// RVec<float> v {1.f, 2.f, 4.f};
889/// auto v_var = Var(v);
890/// v_var
891/// // (double) 2.3333333
892/// ~~~
893template <typename T>
894double Var(const RVec<T> &v)
895{
896 const std::size_t size = v.size();
897 if (size < std::size_t(2)) return 0.;
898 T sum_squares(0), squared_sum(0);
899 auto pred = [&sum_squares, &squared_sum](const T& x) {sum_squares+=x*x; squared_sum+=x;};
900 std::for_each(v.begin(), v.end(), pred);
901 squared_sum *= squared_sum;
902 const auto dsize = (double) size;
903 return 1. / (dsize - 1.) * (sum_squares - squared_sum / dsize );
904}
905
906/// Get the standard deviation of the elements of an RVec
907///
908/// The return type is a double precision floating point number.
909/// Example code, at the ROOT prompt:
910/// ~~~{.cpp}
911/// using namespace ROOT::VecOps;
912/// RVec<float> v {1.f, 2.f, 4.f};
913/// auto v_sd = StdDev(v);
914/// v_sd
915/// // (double) 1.5275252
916/// ~~~
917template <typename T>
918double StdDev(const RVec<T> &v)
919{
920 return std::sqrt(Var(v));
921}
922
923/// Create new collection applying a callable to the elements of the input collection
924///
925/// Example code, at the ROOT prompt:
926/// ~~~{.cpp}
927/// using namespace ROOT::VecOps;
928/// RVec<float> v {1.f, 2.f, 4.f};
929/// auto v_square = Map(v, [](float f){return f* 2.f;});
930/// v_square
931/// // (ROOT::VecOps::RVec<float> &) { 2.00000f, 4.00000f, 8.00000f }
932///
933/// RVec<float> x({1.f, 2.f, 3.f});
934/// RVec<float> y({4.f, 5.f, 6.f});
935/// RVec<float> z({7.f, 8.f, 9.f});
936/// auto mod = [](float x, float y, float z) { return sqrt(x * x + y * y + z * z); };
937/// auto v_mod = Map(x, y, z, mod);
938/// v_mod
939/// // (ROOT::VecOps::RVec<float> &) { 8.12404f, 9.64365f, 11.2250f }
940/// ~~~
941template <typename... Args>
942auto Map(Args &&... args)
943 -> decltype(ROOT::Detail::VecOps::MapFromTuple(std::forward_as_tuple(args...),
944 std::make_index_sequence<sizeof...(args) - 1>()))
945{
946 /*
947 Here the strategy in order to generalise the previous implementation of Map, i.e.
948 `RVec Map(RVec, F)`, here we need to move the last parameter of the pack in first
949 position in order to be able to invoke the Map function with automatic type deduction.
950 This is achieved in two steps:
951 1. Forward as tuple the pack to MapFromTuple
952 2. Invoke the MapImpl helper which has the signature `template<...T, F> RVec MapImpl(F &&f, RVec<T>...)`
953 NOTA BENE: the signature is very heavy but it is one of the lightest ways to manage in C++11
954 to build the return type based on the template args.
955 */
956 return ROOT::Detail::VecOps::MapFromTuple(std::forward_as_tuple(args...),
957 std::make_index_sequence<sizeof...(args) - 1>());
958}
959
960/// Create a new collection with the elements passing the filter expressed by the predicate
961///
962/// Example code, at the ROOT prompt:
963/// ~~~{.cpp}
964/// using namespace ROOT::VecOps;
965/// RVec<int> v {1, 2, 4};
966/// auto v_even = Filter(v, [](int i){return 0 == i%2;});
967/// v_even
968/// // (ROOT::VecOps::RVec<int> &) { 2, 4 }
969/// ~~~
970template <typename T, typename F>
972{
973 const auto thisSize = v.size();
974 RVec<T> w;
975 w.reserve(thisSize);
976 for (auto &&val : v) {
977 if (f(val))
978 w.emplace_back(val);
979 }
980 return w;
981}
982
983/// Return true if any of the elements equates to true, return false otherwise.
984///
985/// Example code, at the ROOT prompt:
986/// ~~~{.cpp}
987/// using namespace ROOT::VecOps;
988/// RVec<int> v {0, 1, 0};
989/// auto anyTrue = Any(v);
990/// anyTrue
991/// // (bool) true
992/// ~~~
993template <typename T>
994auto Any(const RVec<T> &v) -> decltype(v[0] == true)
995{
996 for (auto &&e : v)
997 if (e == true)
998 return true;
999 return false;
1000}
1001
1002/// Return true if all of the elements equate to true, return false otherwise.
1003///
1004/// Example code, at the ROOT prompt:
1005/// ~~~{.cpp}
1006/// using namespace ROOT::VecOps;
1007/// RVec<int> v {0, 1, 0};
1008/// auto allTrue = All(v);
1009/// allTrue
1010/// // (bool) false
1011/// ~~~
1012template <typename T>
1013auto All(const RVec<T> &v) -> decltype(v[0] == false)
1014{
1015 for (auto &&e : v)
1016 if (e == false)
1017 return false;
1018 return true;
1019}
1020
1021template <typename T>
1022void swap(RVec<T> &lhs, RVec<T> &rhs)
1023{
1024 lhs.swap(rhs);
1025}
1026
1027/// Return an RVec of indices that sort the input RVec
1028///
1029/// Example code, at the ROOT prompt:
1030/// ~~~{.cpp}
1031/// using namespace ROOT::VecOps;
1032/// RVec<double> v {2., 3., 1.};
1033/// auto sortIndices = Argsort(v);
1034/// sortIndices
1035/// // (ROOT::VecOps::RVec<unsigned long> &) { 2, 0, 1 }
1036/// ~~~
1037template <typename T>
1039{
1040 using size_type = typename RVec<T>::size_type;
1041 RVec<size_type> i(v.size());
1042 std::iota(i.begin(), i.end(), 0);
1043 std::sort(i.begin(), i.end(), [&v](size_type i1, size_type i2) { return v[i1] < v[i2]; });
1044 return i;
1045}
1046
1047/// Return elements of a vector at given indices
1048///
1049/// Example code, at the ROOT prompt:
1050/// ~~~{.cpp}
1051/// using namespace ROOT::VecOps;
1052/// RVec<double> v {2., 3., 1.};
1053/// auto vTaken = Take(v, {0,2});
1054/// vTaken
1055/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 1.0000000 }
1056/// ~~~
1057template <typename T>
1058RVec<T> Take(const RVec<T> &v, const RVec<typename RVec<T>::size_type> &i)
1059{
1060 using size_type = typename RVec<T>::size_type;
1061 const size_type isize = i.size();
1062 RVec<T> r(isize);
1063 for (size_type k = 0; k < isize; k++)
1064 r[k] = v[i[k]];
1065 return r;
1066}
1067
1068/// Return first or last `n` elements of an RVec
1069///
1070/// if `n > 0` and last elements if `n < 0`.
1071///
1072/// Example code, at the ROOT prompt:
1073/// ~~~{.cpp}
1074/// using namespace ROOT::VecOps;
1075/// RVec<double> v {2., 3., 1.};
1076/// auto firstTwo = Take(v, 2);
1077/// firstTwo
1078/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 3.0000000 }
1079/// auto lastOne = Take(v, -1);
1080/// lastOne
1081/// // (ROOT::VecOps::RVec<double>) { 1.0000000 }
1082/// ~~~
1083template <typename T>
1084RVec<T> Take(const RVec<T> &v, const int n)
1085{
1086 using size_type = typename RVec<T>::size_type;
1087 const size_type size = v.size();
1088 const size_type absn = std::abs(n);
1089 if (absn > size) {
1090 std::stringstream ss;
1091 ss << "Try to take " << absn << " elements but vector has only size " << size << ".";
1092 throw std::runtime_error(ss.str());
1093 }
1094 RVec<T> r(absn);
1095 if (n < 0) {
1096 for (size_type k = 0; k < absn; k++)
1097 r[k] = v[size - absn + k];
1098 } else {
1099 for (size_type k = 0; k < absn; k++)
1100 r[k] = v[k];
1101 }
1102 return r;
1103}
1104
1105/// Return copy of reversed vector
1106///
1107/// Example code, at the ROOT prompt:
1108/// ~~~{.cpp}
1109/// using namespace ROOT::VecOps;
1110/// RVec<double> v {2., 3., 1.};
1111/// auto v_reverse = Reverse(v);
1112/// v_reverse
1113/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 3.0000000, 2.0000000 }
1114/// ~~~
1115template <typename T>
1117{
1118 RVec<T> r(v);
1119 std::reverse(r.begin(), r.end());
1120 return r;
1121}
1122
1123/// Return copy of RVec with elements sorted in ascending order
1124///
1125/// This helper is different from ArgSort since it does not return an RVec of indices,
1126/// but an RVec of values.
1127///
1128/// Example code, at the ROOT prompt:
1129/// ~~~{.cpp}
1130/// using namespace ROOT::VecOps;
1131/// RVec<double> v {2., 3., 1.};
1132/// auto v_sorted = Sort(v);
1133/// v_sorted
1134/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 2.0000000, 3.0000000 }
1135/// ~~~
1136template <typename T>
1137RVec<T> Sort(const RVec<T> &v)
1138{
1139 RVec<T> r(v);
1140 std::sort(r.begin(), r.end());
1141 return r;
1142}
1143
1144/// Return copy of RVec with elements sorted based on a comparison operator
1145///
1146/// The comparison operator has to fullfill the same requirements of the
1147/// predicate of by std::sort.
1148///
1149///
1150/// This helper is different from ArgSort since it does not return an RVec of indices,
1151/// but an RVec of values.
1152///
1153/// Example code, at the ROOT prompt:
1154/// ~~~{.cpp}
1155/// using namespace ROOT::VecOps;
1156/// RVec<double> v {2., 3., 1.};
1157/// auto v_sorted = Sort(v, [](double x, double y) {return 1/x < 1/y;});
1158/// v_sorted
1159/// // (ROOT::VecOps::RVec<double>) { 3.0000000, 2.0000000, 1.0000000 }
1160/// ~~~
1161template <typename T, typename Compare>
1162RVec<T> Sort(const RVec<T> &v, Compare &&c)
1163{
1164 RVec<T> r(v);
1165 std::sort(r.begin(), r.end(), std::forward<Compare>(c));
1166 return r;
1167}
1168
1169/// Return the indices that represent all combinations of the elements of two
1170/// RVecs.
1171///
1172/// The type of the return value is an RVec of two RVecs containing indices.
1173///
1174/// Example code, at the ROOT prompt:
1175/// ~~~{.cpp}
1176/// using namespace ROOT::VecOps;
1177/// auto comb_idx = Combinations(3, 2);
1178/// comb_idx
1179/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0, 0, 1, 1, 2, 2 }, { 0, 1, 0, 1, 0, 1 } }
1180/// ~~~
1181inline RVec<RVec<std::size_t>> Combinations(const std::size_t size1, const std::size_t size2)
1182{
1183 using size_type = std::size_t;
1185 r[0].resize(size1*size2);
1186 r[1].resize(size1*size2);
1187 size_type c = 0;
1188 for(size_type i=0; i<size1; i++) {
1189 for(size_type j=0; j<size2; j++) {
1190 r[0][c] = i;
1191 r[1][c] = j;
1192 c++;
1193 }
1194 }
1195 return r;
1196}
1197
1198/// Return the indices that represent all combinations of the elements of two
1199/// RVecs.
1200///
1201/// The type of the return value is an RVec of two RVecs containing indices.
1202///
1203/// Example code, at the ROOT prompt:
1204/// ~~~{.cpp}
1205/// using namespace ROOT::VecOps;
1206/// RVec<double> v1 {1., 2., 3.};
1207/// RVec<double> v2 {-4., -5.};
1208/// auto comb_idx = Combinations(v1, v2);
1209/// comb_idx
1210/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0, 0, 1, 1, 2, 2 }, { 0, 1,
1211/// 0, 1, 0, 1 } }
1212/// ~~~
1213template <typename T1, typename T2>
1215{
1216 return Combinations(v1.size(), v2.size());
1217}
1218
1219/// Return the indices that represent all unique combinations of the
1220/// elements of a given RVec.
1221///
1222/// ~~~{.cpp}
1223/// using namespace ROOT::VecOps;
1224/// RVec<double> v {1., 2., 3., 4.};
1225/// auto v_1 = Combinations(v, 1);
1226/// v_1
1227/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0, 1, 2, 3 } }
1228/// auto v_2 = Combinations(v, 2);
1229/// auto v_2
1230/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0, 0, 0, 1, 1, 2 }, { 1, 2, 3, 2, 3, 3 } }
1231/// auto v_3 = Combinations(v, 3);
1232/// v_3
1233/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0, 0, 0, 1 }, { 1, 1, 2, 2 }, { 2, 3, 3, 3 } }
1234/// auto v_4 = Combinations(v, 4);
1235/// v_4
1236/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type> >) { { 0 }, { 1 }, { 2 }, { 3 } }
1237/// ~~~
1238template <typename T>
1240{
1241 using size_type = typename RVec<T>::size_type;
1242 const size_type s = v.size();
1243 if (n > s) {
1244 std::stringstream ss;
1245 ss << "Cannot make unique combinations of size " << n << " from vector of size " << s << ".";
1246 throw std::runtime_error(ss.str());
1247 }
1249 for(size_type k=0; k<s; k++)
1250 indices[k] = k;
1252 for(size_type k=0; k<n; k++)
1253 c[k].emplace_back(indices[k]);
1254 while (true) {
1255 bool run_through = true;
1256 long i = n - 1;
1257 for (; i>=0; i--) {
1258 if (indices[i] != i + s - n){
1259 run_through = false;
1260 break;
1261 }
1262 }
1263 if (run_through) {
1264 return c;
1265 }
1266 indices[i]++;
1267 for (long j=i+1; j<(long)n; j++)
1268 indices[j] = indices[j-1] + 1;
1269 for(size_type k=0; k<n; k++)
1270 c[k].emplace_back(indices[k]);
1271 }
1272}
1273
1274/// Return the indices of the elements which are not zero
1275///
1276/// Example code, at the ROOT prompt:
1277/// ~~~{.cpp}
1278/// using namespace ROOT::VecOps;
1279/// RVec<double> v {2., 0., 3., 0., 1.};
1280/// auto nonzero_idx = Nonzero(v);
1281/// nonzero_idx
1282/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<double>::size_type>) { 0, 2, 4 }
1283/// ~~~
1284template <typename T>
1286{
1287 using size_type = typename RVec<T>::size_type;
1289 const auto size = v.size();
1290 r.reserve(size);
1291 for(size_type i=0; i<size; i++) {
1292 if(v[i] != 0) {
1293 r.emplace_back(i);
1294 }
1295 }
1296 return r;
1297}
1298
1299/// Return the intersection of elements of two RVecs.
1300///
1301/// Each element of v1 is looked up in v2 and added to the returned vector if
1302/// found. Following, the order of v1 is preserved. If v2 is already sorted, the
1303/// optional argument v2_is_sorted can be used to toggle of the internal sorting
1304/// step, therewith optimising runtime.
1305///
1306/// Example code, at the ROOT prompt:
1307/// ~~~{.cpp}
1308/// using namespace ROOT::VecOps;
1309/// RVec<double> v1 {1., 2., 3.};
1310/// RVec<double> v2 {-4., -5., 2., 1.};
1311/// auto v1_intersect_v2 = Intersect(v1, v2);
1312/// v1_intersect_v2
1313/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 2.0000000 }
1314/// ~~~
1315template <typename T>
1316RVec<T> Intersect(const RVec<T>& v1, const RVec<T>& v2, bool v2_is_sorted = false)
1317{
1318 RVec<T> v2_sorted;
1319 if (!v2_is_sorted) v2_sorted = Sort(v2);
1320 const auto v2_begin = v2_is_sorted ? v2.begin() : v2_sorted.begin();
1321 const auto v2_end = v2_is_sorted ? v2.end() : v2_sorted.end();
1322 RVec<T> r;
1323 const auto size = v1.size();
1324 r.reserve(size);
1325 using size_type = typename RVec<T>::size_type;
1326 for(size_type i=0; i<size; i++) {
1327 if (std::binary_search(v2_begin, v2_end, v1[i])) {
1328 r.emplace_back(v1[i]);
1329 }
1330 }
1331 return r;
1332}
1333
1334/// Return the elements of v1 if the condition c is true and v2 if the
1335/// condition c is false.
1336///
1337/// Example code, at the ROOT prompt:
1338/// ~~~{.cpp}
1339/// using namespace ROOT::VecOps;
1340/// RVec<double> v1 {1., 2., 3.};
1341/// RVec<double> v2 {-1., -2., -3.};
1342/// auto c = v1 > 1;
1343/// c
1344/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
1345/// auto if_c_v1_else_v2 = Where(c, v1, v2);
1346/// if_c_v1_else_v2
1347/// // (ROOT::VecOps::RVec<double> &) { -1.0000000, 2.0000000, 3.0000000 }
1348/// ~~~
1349template <typename T>
1350RVec<T> Where(const RVec<int>& c, const RVec<T>& v1, const RVec<T>& v2)
1351{
1352 using size_type = typename RVec<T>::size_type;
1353 const size_type size = c.size();
1354 RVec<T> r;
1355 r.reserve(size);
1356 for (size_type i=0; i<size; i++) {
1357 r.emplace_back(c[i] != 0 ? v1[i] : v2[i]);
1358 }
1359 return r;
1360}
1361
1362/// Return the elements of v1 if the condition c is true and sets the value v2
1363/// if the condition c is false.
1364///
1365/// Example code, at the ROOT prompt:
1366/// ~~~{.cpp}
1367/// using namespace ROOT::VecOps;
1368/// RVec<double> v1 {1., 2., 3.};
1369/// double v2 = 4.;
1370/// auto c = v1 > 1;
1371/// c
1372/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
1373/// auto if_c_v1_else_v2 = Where(c, v1, v2);
1374/// if_c_v1_else_v2
1375/// // (ROOT::VecOps::RVec<double>) { 4.0000000, 2.0000000, 3.0000000 }
1376/// ~~~
1377template <typename T>
1379{
1380 using size_type = typename RVec<T>::size_type;
1381 const size_type size = c.size();
1382 RVec<T> r;
1383 r.reserve(size);
1384 for (size_type i=0; i<size; i++) {
1385 r.emplace_back(c[i] != 0 ? v1[i] : v2);
1386 }
1387 return r;
1388}
1389
1390/// Return the elements of v2 if the condition c is false and sets the value v1
1391/// if the condition c is true.
1392///
1393/// Example code, at the ROOT prompt:
1394/// ~~~{.cpp}
1395/// using namespace ROOT::VecOps;
1396/// double v1 = 4.;
1397/// RVec<double> v2 {1., 2., 3.};
1398/// auto c = v2 > 1;
1399/// c
1400/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
1401/// auto if_c_v1_else_v2 = Where(c, v1, v2);
1402/// if_c_v1_else_v2
1403/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 4.0000000, 4.0000000 }
1404/// ~~~
1405template <typename T>
1407{
1408 using size_type = typename RVec<T>::size_type;
1409 const size_type size = c.size();
1410 RVec<T> r;
1411 r.reserve(size);
1412 for (size_type i=0; i<size; i++) {
1413 r.emplace_back(c[i] != 0 ? v1 : v2[i]);
1414 }
1415 return r;
1416}
1417
1418/// Return a vector with the value v2 if the condition c is false and sets the
1419/// value v1 if the condition c is true.
1420///
1421/// Example code, at the ROOT prompt:
1422/// ~~~{.cpp}
1423/// using namespace ROOT::VecOps;
1424/// double v1 = 4.;
1425/// double v2 = 2.;
1426/// RVec<int> c {0, 1, 1};
1427/// auto if_c_v1_else_v2 = Where(c, v1, v2);
1428/// if_c_v1_else_v2
1429/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 4.0000000, 4.0000000 }
1430/// ~~~
1431template <typename T>
1433{
1434 using size_type = typename RVec<T>::size_type;
1435 const size_type size = c.size();
1436 RVec<T> r;
1437 r.reserve(size);
1438 for (size_type i=0; i<size; i++) {
1439 r.emplace_back(c[i] != 0 ? v1 : v2);
1440 }
1441 return r;
1442}
1443
1444/// Return the concatenation of two RVecs.
1445///
1446/// Example code, at the ROOT prompt:
1447/// ~~~{.cpp}
1448/// using namespace ROOT::VecOps;
1449/// RVec<float> rvf {0.f, 1.f, 2.f};
1450/// RVec<int> rvi {7, 8, 9};
1451/// Concatenate(rvf, rvi);
1452/// // (ROOT::VecOps::RVec<float>) { 2.0000000, 4.0000000, 4.0000000 }
1453/// ~~~
1454template <typename T0, typename T1, typename Common_t = typename std::common_type<T0, T1>::type>
1456{
1457 RVec<Common_t> res;
1458 res.reserve(v0.size() + v1.size());
1459 std::copy(v0.begin(), v0.end(), std::back_inserter(res));
1460 std::copy(v1.begin(), v1.end(), std::back_inserter(res));
1461 return res;
1462}
1463
1464/// Return the angle difference \f$\Delta \phi\f$ of two scalars.
1465///
1466/// The function computes the closest angle from v1 to v2 with sign and is
1467/// therefore in the range \f$[-\pi, \pi]\f$.
1468/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
1469/// to degrees \f$c = 180\f$.
1470template <typename T>
1471T DeltaPhi(T v1, T v2, const T c = M_PI)
1472{
1473 static_assert(std::is_floating_point<T>::value,
1474 "DeltaPhi must be called with floating point values.");
1475 auto r = std::fmod(v2 - v1, 2.0 * c);
1476 if (r < -c) {
1477 r += 2.0 * c;
1478 }
1479 else if (r > c) {
1480 r -= 2.0 * c;
1481 }
1482 return r;
1483}
1484
1485/// Return the angle difference \f$\Delta \phi\f$ in radians of two vectors.
1486///
1487/// The function computes the closest angle from v1 to v2 with sign and is
1488/// therefore in the range \f$[-\pi, \pi]\f$.
1489/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
1490/// to degrees \f$c = 180\f$.
1491template <typename T>
1492RVec<T> DeltaPhi(const RVec<T>& v1, const RVec<T>& v2, const T c = M_PI)
1493{
1494 using size_type = typename RVec<T>::size_type;
1495 const size_type size = v1.size();
1496 auto r = RVec<T>(size);
1497 for (size_type i = 0; i < size; i++) {
1498 r[i] = DeltaPhi(v1[i], v2[i], c);
1499 }
1500 return r;
1501}
1502
1503/// Return the angle difference \f$\Delta \phi\f$ in radians of a vector and a scalar.
1504///
1505/// The function computes the closest angle from v1 to v2 with sign and is
1506/// therefore in the range \f$[-\pi, \pi]\f$.
1507/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
1508/// to degrees \f$c = 180\f$.
1509template <typename T>
1510RVec<T> DeltaPhi(const RVec<T>& v1, T v2, const T c = M_PI)
1511{
1512 using size_type = typename RVec<T>::size_type;
1513 const size_type size = v1.size();
1514 auto r = RVec<T>(size);
1515 for (size_type i = 0; i < size; i++) {
1516 r[i] = DeltaPhi(v1[i], v2, c);
1517 }
1518 return r;
1519}
1520
1521/// Return the angle difference \f$\Delta \phi\f$ in radians of a scalar and a vector.
1522///
1523/// The function computes the closest angle from v1 to v2 with sign and is
1524/// therefore in the range \f$[-\pi, \pi]\f$.
1525/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
1526/// to degrees \f$c = 180\f$.
1527template <typename T>
1528RVec<T> DeltaPhi(T v1, const RVec<T>& v2, const T c = M_PI)
1529{
1530 using size_type = typename RVec<T>::size_type;
1531 const size_type size = v2.size();
1532 auto r = RVec<T>(size);
1533 for (size_type i = 0; i < size; i++) {
1534 r[i] = DeltaPhi(v1, v2[i], c);
1535 }
1536 return r;
1537}
1538
1539/// Return the square of the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
1540/// the collections eta1, eta2, phi1 and phi2.
1541///
1542/// The function computes \f$\Delta R^2 = (\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2\f$
1543/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
1544/// be set to radian or degrees using the optional argument c, see the documentation
1545/// of the DeltaPhi helper.
1546template <typename T>
1547RVec<T> DeltaR2(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
1548{
1549 const auto dphi = DeltaPhi(phi1, phi2, c);
1550 return (eta1 - eta2) * (eta1 - eta2) + dphi * dphi;
1551}
1552
1553/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
1554/// the collections eta1, eta2, phi1 and phi2.
1555///
1556/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
1557/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
1558/// be set to radian or degrees using the optional argument c, see the documentation
1559/// of the DeltaPhi helper.
1560template <typename T>
1561RVec<T> DeltaR(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
1562{
1563 return sqrt(DeltaR2(eta1, eta2, phi1, phi2, c));
1564}
1565
1566/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
1567/// the scalars eta1, eta2, phi1 and phi2.
1568///
1569/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
1570/// of the given scalars eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
1571/// be set to radian or degrees using the optional argument c, see the documentation
1572/// of the DeltaPhi helper.
1573template <typename T>
1574T DeltaR(T eta1, T eta2, T phi1, T phi2, const T c = M_PI)
1575{
1576 const auto dphi = DeltaPhi(phi1, phi2, c);
1577 return std::sqrt((eta1 - eta2) * (eta1 - eta2) + dphi * dphi);
1578}
1579
1580/// Return the invariant mass of two particles given the collections of the quantities
1581/// transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
1582///
1583/// The function computes the invariant mass of two particles with the four-vectors
1584/// (pt1, eta2, phi1, mass1) and (pt2, eta2, phi2, mass2).
1585template <typename T>
1587 const RVec<T>& pt1, const RVec<T>& eta1, const RVec<T>& phi1, const RVec<T>& mass1,
1588 const RVec<T>& pt2, const RVec<T>& eta2, const RVec<T>& phi2, const RVec<T>& mass2)
1589{
1590 std::size_t size = pt1.size();
1591
1592 R__ASSERT(eta1.size() == size && phi1.size() == size && mass1.size() == size);
1593 R__ASSERT(pt2.size() == size && phi2.size() == size && mass2.size() == size);
1594
1595 RVec<T> inv_masses(size);
1596
1597 for (std::size_t i = 0u; i < size; ++i) {
1598 // Conversion from (pt, eta, phi, mass) to (x, y, z, e) coordinate system
1599 const auto x1 = pt1[i] * std::cos(phi1[i]);
1600 const auto y1 = pt1[i] * std::sin(phi1[i]);
1601 const auto z1 = pt1[i] * std::sinh(eta1[i]);
1602 const auto e1 = std::sqrt(x1 * x1 + y1 * y1 + z1 * z1 + mass1[i] * mass1[i]);
1603
1604 const auto x2 = pt2[i] * std::cos(phi2[i]);
1605 const auto y2 = pt2[i] * std::sin(phi2[i]);
1606 const auto z2 = pt2[i] * std::sinh(eta2[i]);
1607 const auto e2 = std::sqrt(x2 * x2 + y2 * y2 + z2 * z2 + mass2[i] * mass2[i]);
1608
1609 // Addition of particle four-vector elements
1610 const auto e = e1 + e2;
1611 const auto x = x1 + x2;
1612 const auto y = y1 + y2;
1613 const auto z = z1 + z2;
1614
1615 inv_masses[i] = std::sqrt(e * e - x * x - y * y - z * z);
1616 }
1617
1618 // Return invariant mass with (+, -, -, -) metric
1619 return inv_masses;
1620}
1621
1622/// Return the invariant mass of multiple particles given the collections of the
1623/// quantities transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
1624///
1625/// The function computes the invariant mass of multiple particles with the
1626/// four-vectors (pt, eta, phi, mass).
1627template <typename T>
1628T InvariantMass(const RVec<T>& pt, const RVec<T>& eta, const RVec<T>& phi, const RVec<T>& mass)
1629{
1630 const std::size_t size = pt.size();
1631
1632 R__ASSERT(eta.size() == size && phi.size() == size && mass.size() == size);
1633
1634 T x_sum = 0.;
1635 T y_sum = 0.;
1636 T z_sum = 0.;
1637 T e_sum = 0.;
1638
1639 for (std::size_t i = 0u; i < size; ++ i) {
1640 // Convert to (e, x, y, z) coordinate system and update sums
1641 const auto x = pt[i] * std::cos(phi[i]);
1642 x_sum += x;
1643 const auto y = pt[i] * std::sin(phi[i]);
1644 y_sum += y;
1645 const auto z = pt[i] * std::sinh(eta[i]);
1646 z_sum += z;
1647 const auto e = std::sqrt(x * x + y * y + z * z + mass[i] * mass[i]);
1648 e_sum += e;
1649 }
1650
1651 // Return invariant mass with (+, -, -, -) metric
1652 return std::sqrt(e_sum * e_sum - x_sum * x_sum - y_sum * y_sum - z_sum * z_sum);
1653}
1654
1655////////////////////////////////////////////////////////////////////////////
1656/// \brief Build an RVec of objects starting from RVecs of input to their constructors.
1657/// \tparam T Type of the objects contained in the created RVec.
1658/// \tparam Args_t Pack of types templating the input RVecs.
1659/// \param[in] args The RVecs containing the values used to initialise the output objects.
1660/// \return The RVec of objects initialised with the input parameters.
1661///
1662/// Example code, at the ROOT prompt:
1663/// ~~~{.cpp}
1664/// using namespace ROOT::VecOps;
1665/// RVec<float> pts = {15.5, 34.32, 12.95};
1666/// RVec<float> etas = {0.3, 2.2, 1.32};
1667/// RVec<float> phis = {0.1, 3.02, 2.2};
1668/// RVec<float> masses = {105.65, 105.65, 105.65};
1669/// auto fourVecs = Construct<ROOT::Math::PtEtaPhiMVector>(pts, etas, phis, masses);
1670/// cout << fourVecs << endl;
1671/// // { (15.5,0.3,0.1,105.65), (34.32,2.2,3.02,105.65), (12.95,1.32,2.2,105.65) }
1672/// ~~~
1673template <typename T, typename... Args_t>
1675{
1676 const auto size = ::ROOT::Detail::VecOps::GetVectorsSize("Construct", args...);
1677 RVec<T> ret;
1678 ret.reserve(size);
1679 for (auto i = 0UL; i < size; ++i) {
1680 ret.emplace_back(args[i]...);
1681 }
1682 return ret;
1683}
1684
1685////////////////////////////////////////////////////////////////////////////////
1686/// Print a RVec at the prompt:
1687template <class T>
1688std::ostream &operator<<(std::ostream &os, const RVec<T> &v)
1689{
1690 // In order to print properly, convert to 64 bit int if this is a char
1691 constexpr bool mustConvert = std::is_same<char, T>::value || std::is_same<signed char, T>::value ||
1692 std::is_same<unsigned char, T>::value || std::is_same<wchar_t, T>::value ||
1693 std::is_same<char16_t, T>::value || std::is_same<char32_t, T>::value;
1694 using Print_t = typename std::conditional<mustConvert, long long int, T>::type;
1695 os << "{ ";
1696 auto size = v.size();
1697 if (size) {
1698 for (std::size_t i = 0; i < size - 1; ++i) {
1699 os << (Print_t)v[i] << ", ";
1700 }
1701 os << (Print_t)v[size - 1];
1702 }
1703 os << " }";
1704 return os;
1705}
1706
1707#if (_VECOPS_USE_EXTERN_TEMPLATES)
1708
1709#define RVEC_EXTERN_UNARY_OPERATOR(T, OP) \
1710 extern template RVec<T> operator OP<T>(const RVec<T> &);
1711
1712#define RVEC_EXTERN_BINARY_OPERATOR(T, OP) \
1713 extern template auto operator OP<T, T>(const T &x, const RVec<T> &v) \
1714 -> RVec<decltype(x OP v[0])>; \
1715 extern template auto operator OP<T, T>(const RVec<T> &v, const T &y) \
1716 -> RVec<decltype(v[0] OP y)>; \
1717 extern template auto operator OP<T, T>(const RVec<T> &v0, const RVec<T> &v1)\
1718 -> RVec<decltype(v0[0] OP v1[0])>;
1719
1720#define RVEC_EXTERN_ASSIGN_OPERATOR(T, OP) \
1721 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const T &); \
1722 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const RVec<T> &);
1723
1724#define RVEC_EXTERN_LOGICAL_OPERATOR(T, OP) \
1725 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const T &); \
1726 extern template RVec<int> operator OP<T, T>(const T &, const RVec<T> &); \
1727 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const RVec<T> &);
1728
1729#define RVEC_EXTERN_FLOAT_TEMPLATE(T) \
1730 extern template class RVec<T>; \
1731 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
1732 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
1733 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
1734 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
1735 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
1736 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
1737 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
1738 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
1739 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
1740 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
1741 RVEC_EXTERN_ASSIGN_OPERATOR(T, /=) \
1742 RVEC_EXTERN_LOGICAL_OPERATOR(T, <) \
1743 RVEC_EXTERN_LOGICAL_OPERATOR(T, >) \
1744 RVEC_EXTERN_LOGICAL_OPERATOR(T, ==) \
1745 RVEC_EXTERN_LOGICAL_OPERATOR(T, !=) \
1746 RVEC_EXTERN_LOGICAL_OPERATOR(T, <=) \
1747 RVEC_EXTERN_LOGICAL_OPERATOR(T, >=) \
1748 RVEC_EXTERN_LOGICAL_OPERATOR(T, &&) \
1749 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
1750
1751#define RVEC_EXTERN_INTEGER_TEMPLATE(T) \
1752 extern template class RVec<T>; \
1753 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
1754 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
1755 RVEC_EXTERN_UNARY_OPERATOR(T, ~) \
1756 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
1757 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
1758 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
1759 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
1760 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
1761 RVEC_EXTERN_BINARY_OPERATOR(T, %) \
1762 RVEC_EXTERN_BINARY_OPERATOR(T, &) \
1763 RVEC_EXTERN_BINARY_OPERATOR(T, |) \
1764 RVEC_EXTERN_BINARY_OPERATOR(T, ^) \
1765 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
1766 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
1767 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
1768 RVEC_EXTERN_ASSIGN_OPERATOR(T, /=) \
1769 RVEC_EXTERN_ASSIGN_OPERATOR(T, %=) \
1770 RVEC_EXTERN_ASSIGN_OPERATOR(T, &=) \
1771 RVEC_EXTERN_ASSIGN_OPERATOR(T, |=) \
1772 RVEC_EXTERN_ASSIGN_OPERATOR(T, ^=) \
1773 RVEC_EXTERN_ASSIGN_OPERATOR(T, >>=) \
1774 RVEC_EXTERN_ASSIGN_OPERATOR(T, <<=) \
1775 RVEC_EXTERN_LOGICAL_OPERATOR(T, <) \
1776 RVEC_EXTERN_LOGICAL_OPERATOR(T, >) \
1777 RVEC_EXTERN_LOGICAL_OPERATOR(T, ==) \
1778 RVEC_EXTERN_LOGICAL_OPERATOR(T, !=) \
1779 RVEC_EXTERN_LOGICAL_OPERATOR(T, <=) \
1780 RVEC_EXTERN_LOGICAL_OPERATOR(T, >=) \
1781 RVEC_EXTERN_LOGICAL_OPERATOR(T, &&) \
1782 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
1783
1784RVEC_EXTERN_INTEGER_TEMPLATE(char)
1785RVEC_EXTERN_INTEGER_TEMPLATE(short)
1786RVEC_EXTERN_INTEGER_TEMPLATE(int)
1787RVEC_EXTERN_INTEGER_TEMPLATE(long)
1788//RVEC_EXTERN_INTEGER_TEMPLATE(long long)
1789
1790RVEC_EXTERN_INTEGER_TEMPLATE(unsigned char)
1791RVEC_EXTERN_INTEGER_TEMPLATE(unsigned short)
1792RVEC_EXTERN_INTEGER_TEMPLATE(unsigned int)
1793RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long)
1794//RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long long)
1795
1796RVEC_EXTERN_FLOAT_TEMPLATE(float)
1797RVEC_EXTERN_FLOAT_TEMPLATE(double)
1798
1799#undef RVEC_EXTERN_UNARY_OPERATOR
1800#undef RVEC_EXTERN_BINARY_OPERATOR
1801#undef RVEC_EXTERN_ASSIGN_OPERATOR
1802#undef RVEC_EXTERN_LOGICAL_OPERATOR
1803#undef RVEC_EXTERN_INTEGER_TEMPLATE
1804#undef RVEC_EXTERN_FLOAT_TEMPLATE
1805
1806#define RVEC_EXTERN_UNARY_FUNCTION(T, NAME, FUNC) \
1807 extern template RVec<PromoteType<T>> NAME(const RVec<T> &);
1808
1809#define RVEC_EXTERN_STD_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, std::F)
1810
1811#define RVEC_EXTERN_BINARY_FUNCTION(T0, T1, NAME, FUNC) \
1812 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const T1 &); \
1813 extern template RVec<PromoteTypes<T0, T1>> NAME(const T0 &, const RVec<T1> &); \
1814 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const RVec<T1> &);
1815
1816#define RVEC_EXTERN_STD_BINARY_FUNCTION(T, F) RVEC_EXTERN_BINARY_FUNCTION(T, T, F, std::F)
1817
1818#define RVEC_EXTERN_STD_FUNCTIONS(T) \
1819 RVEC_EXTERN_STD_UNARY_FUNCTION(T, abs) \
1820 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fdim) \
1821 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fmod) \
1822 RVEC_EXTERN_STD_BINARY_FUNCTION(T, remainder) \
1823 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp) \
1824 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp2) \
1825 RVEC_EXTERN_STD_UNARY_FUNCTION(T, expm1) \
1826 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log) \
1827 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log10) \
1828 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log2) \
1829 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log1p) \
1830 RVEC_EXTERN_STD_BINARY_FUNCTION(T, pow) \
1831 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sqrt) \
1832 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cbrt) \
1833 RVEC_EXTERN_STD_BINARY_FUNCTION(T, hypot) \
1834 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sin) \
1835 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cos) \
1836 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tan) \
1837 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asin) \
1838 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acos) \
1839 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atan) \
1840 RVEC_EXTERN_STD_BINARY_FUNCTION(T, atan2) \
1841 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sinh) \
1842 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cosh) \
1843 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tanh) \
1844 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asinh) \
1845 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acosh) \
1846 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atanh) \
1847 RVEC_EXTERN_STD_UNARY_FUNCTION(T, floor) \
1848 RVEC_EXTERN_STD_UNARY_FUNCTION(T, ceil) \
1849 RVEC_EXTERN_STD_UNARY_FUNCTION(T, trunc) \
1850 RVEC_EXTERN_STD_UNARY_FUNCTION(T, round) \
1851 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erf) \
1852 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erfc) \
1853 RVEC_EXTERN_STD_UNARY_FUNCTION(T, lgamma) \
1854 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tgamma) \
1855
1856RVEC_EXTERN_STD_FUNCTIONS(float)
1857RVEC_EXTERN_STD_FUNCTIONS(double)
1858#undef RVEC_EXTERN_STD_UNARY_FUNCTION
1859#undef RVEC_EXTERN_STD_BINARY_FUNCTION
1860#undef RVEC_EXTERN_STD_UNARY_FUNCTIONS
1861
1862#ifdef R__HAS_VDT
1863
1864#define RVEC_EXTERN_VDT_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, vdt::F)
1865
1866RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_expf)
1867RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_logf)
1868RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_sinf)
1869RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_cosf)
1870RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_tanf)
1871RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_asinf)
1872RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_acosf)
1873RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_atanf)
1874
1875RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_exp)
1876RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_log)
1877RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_sin)
1878RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_cos)
1879RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_tan)
1880RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_asin)
1881RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_acos)
1882RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_atan)
1883
1884#endif // R__HAS_VDT
1885
1886#endif // _VECOPS_USE_EXTERN_TEMPLATES
1887
1888} // End of VecOps NS
1889
1890// Allow to use RVec as ROOT::RVec
1891using ROOT::VecOps::RVec;
1892
1893} // End of ROOT NS
1894
1895#endif // ROOT_RVEC
double
long
ROOT::R::TRInterface & r
Definition Object.C:4
#define f(i)
Definition RSha256.hxx:104
#define c(i)
Definition RSha256.hxx:101
#define e(i)
Definition RSha256.hxx:103
#define RVEC_UNARY_OPERATOR(OP)
Definition RVec.hxx:466
#define RVEC_ASSIGNMENT_OPERATOR(OP)
Definition RVec.hxx:537
#define RVEC_STD_BINARY_FUNCTION(F)
Definition RVec.hxx:680
#define RVEC_BINARY_OPERATOR(OP)
Definition RVec.hxx:489
#define RVEC_LOGICAL_OPERATOR(OP)
Definition RVec.hxx:573
#define RVEC_STD_UNARY_FUNCTION(F)
Definition RVec.hxx:679
static const double x2[5]
static const double x1[5]
#define M_PI
Definition Rotated.cxx:105
TBuffer & operator<<(TBuffer &buf, const Tmpl *obj)
Definition TBuffer.h:399
#define R__ASSERT(e)
Definition TError.h:120
Int_t Compare(const void *item1, const void *item2)
int type
Definition TGX11.cxx:121
double tanh(double)
double ceil(double)
double acos(double)
double log10(double)
double atan(double)
double sqrt(double)
double asin(double)
A "std::vector"-like collection of values implementing handy operation to analyse them.
Definition RVec.hxx:296
typename std::conditional< IsVecBool, std::vector< bool >, std::vector< T, ::ROOT::Detail::VecOps::RAdoptAllocator< T > > >::type Impl_t
Definition RVec.hxx:310
typename Impl_t::iterator iterator
Definition RVec.hxx:322
static constexpr const auto IsVecBool
Definition RVec.hxx:308
const_iterator begin() const noexcept
Definition RVec.hxx:420
RVec(const RVec< T > &v)
Definition RVec.hxx:338
reverse_iterator rbegin() noexcept
Definition RVec.hxx:425
R__DEPRECATED(6, 26, "Please use `std::vector<T>(rvec.begin(), rvec.end())` instead.") Impl_t &AsVector()
Definition RVec.hxx:381
const_iterator cbegin() const noexcept
Definition RVec.hxx:421
typename Impl_t::value_type value_type
Definition RVec.hxx:311
reference operator[](size_type pos)
Definition RVec.hxx:393
void clear() noexcept
Definition RVec.hxx:439
iterator erase(iterator pos)
Definition RVec.hxx:440
RVec< T > & operator=(std::initializer_list< T > ilist)
Definition RVec.hxx:364
typename Impl_t::const_pointer const_pointer
Definition RVec.hxx:317
RVec(std::initializer_list< T > init)
Definition RVec.hxx:349
void resize(size_type count)
Definition RVec.hxx:458
void push_back(T &&value)
Definition RVec.hxx:442
const_reference back() const
Definition RVec.hxx:415
typename Impl_t::difference_type difference_type
Definition RVec.hxx:313
const_reverse_iterator crbegin() const noexcept
Definition RVec.hxx:427
RVec operator[](const RVec< V > &conds) const
Definition RVec.hxx:397
iterator erase(iterator first, iterator last)
Definition RVec.hxx:441
RVec< T > & operator=(RVec< T > &&v)
Definition RVec.hxx:358
const_reverse_iterator crend() const noexcept
Definition RVec.hxx:430
const_reference at(size_type pos) const
Definition RVec.hxx:388
typename Impl_t::const_iterator const_iterator
Definition RVec.hxx:323
bool empty() const noexcept
Definition RVec.hxx:432
reference back()
Definition RVec.hxx:414
typename Impl_t::reference reference
Definition RVec.hxx:314
void push_back(const value_type &value)
Definition RVec.hxx:443
iterator end() noexcept
Definition RVec.hxx:422
size_type size() const noexcept
Definition RVec.hxx:433
RVec(size_type count, const T &value)
Definition RVec.hxx:336
const_iterator end() const noexcept
Definition RVec.hxx:423
iterator emplace(const_iterator pos, U value)
This method is intended only for arithmetic types unlike the std::vector corresponding one which is g...
Definition RVec.hxx:453
typename Impl_t::const_reference const_reference
Definition RVec.hxx:315
const_iterator cend() const noexcept
Definition RVec.hxx:424
RVec(RVec< T > &&v)
Definition RVec.hxx:340
RVec(size_type count)
Definition RVec.hxx:334
typename Impl_t::const_reverse_iterator const_reverse_iterator
Definition RVec.hxx:325
const_reference front() const
Definition RVec.hxx:413
R__DEPRECATED(6, 26, "Please use `std::vector<T>(rvec.begin(), rvec.end())` instead.") const Impl_t &AsVector() const
Definition RVec.hxx:379
void reserve(size_type new_cap)
Definition RVec.hxx:435
reference at(size_type pos)
Definition RVec.hxx:387
const_reverse_iterator rend() const noexcept
Definition RVec.hxx:429
typename std::conditional< IsVecBool, void, typename Impl_t::pointer >::type data_t
Definition RVec.hxx:320
reference front()
Definition RVec.hxx:412
typename Impl_t::reverse_iterator reverse_iterator
Definition RVec.hxx:324
const_reference operator[](size_type pos) const
Definition RVec.hxx:394
typename Impl_t::pointer pointer
Definition RVec.hxx:316
size_type capacity() const noexcept
Definition RVec.hxx:436
reverse_iterator rend() noexcept
Definition RVec.hxx:428
RVec< T > & operator=(const RVec< T > &v)
Definition RVec.hxx:352
typename std::conditional< IsVecBool, void, typename Impl_t::const_pointer >::type const_data_t
Definition RVec.hxx:321
void shrink_to_fit()
Definition RVec.hxx:437
value_type at(size_type pos, value_type fallback)
No exception thrown. The user specifies the desired value in case the RVec is shorter than pos.
Definition RVec.hxx:390
iterator begin() noexcept
Definition RVec.hxx:419
const_reverse_iterator rbegin() const noexcept
Definition RVec.hxx:426
value_type at(size_type pos, value_type fallback) const
No exception thrown. The user specifies the desired value in case the RVec is shorter than pos.
Definition RVec.hxx:392
typename Impl_t::size_type size_type
Definition RVec.hxx:312
void resize(size_type count, const value_type &value)
Definition RVec.hxx:459
reference emplace_back(Args &&... args)
Definition RVec.hxx:445
const_data_t data() const noexcept
Definition RVec.hxx:417
RVec(const std::vector< T > &v)
Definition RVec.hxx:342
data_t data() noexcept
Definition RVec.hxx:416
void swap(RVec< T > &other)
Definition RVec.hxx:460
RVec(InputIt first, InputIt last)
Definition RVec.hxx:347
RVec(pointer p, size_type n)
Definition RVec.hxx:344
size_type max_size() const noexcept
Definition RVec.hxx:434
A TTree is a list of TBranches.
Definition TBranch.h:89
A TTree represents a columnar dataset.
Definition TTree.h:79
TPaveText * pt
Type
enumeration specifying the integration types.
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
#define F(x, y, z)
std::size_t GetVectorsSize(std::string_view id, const RVec< T > &... vs)
Definition RVec.hxx:87
auto MapFromTuple(Tuple_t &&t, std::index_sequence< Is... >) -> decltype(MapImpl(std::get< std::tuple_size< Tuple_t >::value - 1 >(t), std::get< Is >(t)...))
Definition RVec.hxx:116
auto MapImpl(F &&f, const RVec< T > &... vs) -> RVec< decltype(f(vs[0]...))>
Definition RVec.hxx:104
void SetBranchesHelper(BoolArrayMap &boolArrays, TTree *inputTree, TTree &outputTree, const std::string &inName, const std::string &outName, TBranch *&branch, void *&branchAddress, ROOT::VecOps::RVec< T > *ab, RBranchSet &outputBranches)
std::map< std::string, BoolArray > BoolArrayMap
Definition RVec.hxx:72
void EmplaceBack(T &v, Args &&... args)
Definition RVec.hxx:133
RVec< T > Reverse(const RVec< T > &v)
Return copy of reversed vector.
Definition RVec.hxx:1116
RVec< T > Intersect(const RVec< T > &v1, const RVec< T > &v2, bool v2_is_sorted=false)
Return the intersection of elements of two RVecs.
Definition RVec.hxx:1316
auto All(const RVec< T > &v) -> decltype(v[0]==false)
Return true if all of the elements equate to true, return false otherwise.
Definition RVec.hxx:1013
T Sum(const RVec< T > &v)
Sum elements of an RVec.
Definition RVec.hxx:791
RVec< typename RVec< T >::size_type > Nonzero(const RVec< T > &v)
Return the indices of the elements which are not zero.
Definition RVec.hxx:1285
auto Map(Args &&... args) -> decltype(ROOT::Detail::VecOps::MapFromTuple(std::forward_as_tuple(args...), std::make_index_sequence< sizeof...(args) - 1 >()))
Create new collection applying a callable to the elements of the input collection.
Definition RVec.hxx:942
RVec< Common_t > Concatenate(const RVec< T0 > &v0, const RVec< T1 > &v1)
Return the concatenation of two RVecs.
Definition RVec.hxx:1455
RVec< T > InvariantMasses(const RVec< T > &pt1, const RVec< T > &eta1, const RVec< T > &phi1, const RVec< T > &mass1, const RVec< T > &pt2, const RVec< T > &eta2, const RVec< T > &phi2, const RVec< T > &mass2)
Return the invariant mass of two particles given the collections of the quantities transverse momentu...
Definition RVec.hxx:1586
RVec< T > Take(const RVec< T > &v, const RVec< typename RVec< T >::size_type > &i)
Return elements of a vector at given indices.
Definition RVec.hxx:1058
RVec< T > Construct(const RVec< Args_t > &... args)
Build an RVec of objects starting from RVecs of input to their constructors.
Definition RVec.hxx:1674
RVec< RVec< std::size_t > > Combinations(const std::size_t size1, const std::size_t size2)
Return the indices that represent all combinations of the elements of two RVecs.
Definition RVec.hxx:1181
RVec< T > Where(const RVec< int > &c, const RVec< T > &v1, const RVec< T > &v2)
Return the elements of v1 if the condition c is true and v2 if the condition c is false.
Definition RVec.hxx:1350
auto Any(const RVec< T > &v) -> decltype(v[0]==true)
Return true if any of the elements equates to true, return false otherwise.
Definition RVec.hxx:994
RVec< typename RVec< T >::size_type > Argsort(const RVec< T > &v)
Return an RVec of indices that sort the input RVec.
Definition RVec.hxx:1038
std::size_t ArgMin(const RVec< T > &v)
Get the index of the smallest element of an RVec In case of multiple occurrences of the minimum value...
Definition RVec.hxx:877
double Var(const RVec< T > &v)
Get the variance of the elements of an RVec.
Definition RVec.hxx:894
RVec< T > Filter(const RVec< T > &v, F &&f)
Create a new collection with the elements passing the filter expressed by the predicate.
Definition RVec.hxx:971
std::size_t ArgMax(const RVec< T > &v)
Get the index of the greatest element of an RVec In case of multiple occurrences of the maximum value...
Definition RVec.hxx:859
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
Definition first.py:1
#define Dot(u, v)
Definition normal.c:49