Logo ROOT  
Reference Guide
RVec.hxx
Go to the documentation of this file.
1// Author: Enrico Guiraud, Enric Tejedor, Danilo Piparo CERN 04/2021
2// Implementation adapted from from llvm::SmallVector.
3// See /math/vecops/ARCHITECTURE.md for more information.
4
5/*************************************************************************
6 * Copyright (C) 1995-2021, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13#ifndef ROOT_RVEC
14#define ROOT_RVEC
15
16#if __cplusplus > 201402L
17#define R__RVEC_NODISCARD [[nodiscard]]
18#else
19#define R__RVEC_NODISCARD
20#endif
21
22#ifdef _WIN32
23 #ifndef M_PI
24 #ifndef _USE_MATH_DEFINES
25 #define _USE_MATH_DEFINES
26 #endif
27 #include <math.h>
28 #undef _USE_MATH_DEFINES
29 #endif
30 #define _VECOPS_USE_EXTERN_TEMPLATES false
31#else
32 #define _VECOPS_USE_EXTERN_TEMPLATES true
33#endif
34
35#include <Rtypes.h> // R__CLING_PTRCHECK
36#include <TError.h> // R__ASSERT
37
38#include <algorithm>
39#include <cmath>
40#include <cstring>
41#include <limits> // for numeric_limits
42#include <new>
43#include <numeric> // for inner_product
44#include <sstream>
45#include <stdexcept>
46#include <string>
47#include <tuple>
48#include <type_traits>
49#include <utility>
50#include <vector>
51
52#ifdef R__HAS_VDT
53#include <vdt/vdtMath.h>
54#endif
55
56
57namespace ROOT {
58
59namespace VecOps {
60template<typename T>
61class RVec;
62}
63
64namespace Internal {
65namespace VecOps {
66
67template<typename T>
69
70// clang-format off
71template <typename>
72struct IsRVec : std::false_type {};
73
74template <typename T>
75struct IsRVec<ROOT::VecOps::RVec<T>> : std::true_type {};
76// clang-format on
77
78constexpr bool All(const bool *vals, std::size_t size)
79{
80 for (auto i = 0u; i < size; ++i)
81 if (!vals[i])
82 return false;
83 return true;
84}
85
86template <typename... T>
87std::size_t GetVectorsSize(const std::string &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... RVecs>
104auto MapImpl(F &&f, RVecs &&... 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/// Return the next power of two (in 64-bits) that is strictly greater than A.
124/// Return zero on overflow.
125inline uint64_t NextPowerOf2(uint64_t A)
126{
127 A |= (A >> 1);
128 A |= (A >> 2);
129 A |= (A >> 4);
130 A |= (A >> 8);
131 A |= (A >> 16);
132 A |= (A >> 32);
133 return A + 1;
134}
135
136/// This is all the stuff common to all SmallVectors.
137class R__CLING_PTRCHECK(off) SmallVectorBase {
138public:
139 // This limits the maximum size of an RVec<char> to ~4GB but we don't expect this to ever be a problem,
140 // and we prefer the smaller Size_T to reduce the size of each RVec object.
141 using Size_T = int32_t;
142
143protected:
144 void *fBeginX;
145 /// Always >= 0.
146 // Type is signed only for consistency with fCapacity.
148 /// Always >= -1. fCapacity == -1 indicates the RVec is in "memory adoption" mode.
150
151 /// The maximum value of the Size_T used.
152 static constexpr size_t SizeTypeMax() { return std::numeric_limits<Size_T>::max(); }
153
154 SmallVectorBase() = delete;
155 SmallVectorBase(void *FirstEl, size_t TotalCapacity) : fBeginX(FirstEl), fCapacity(TotalCapacity) {}
156
157 /// This is an implementation of the grow() method which only works
158 /// on POD-like data types and is out of line to reduce code duplication.
159 /// This function will report a fatal error if it cannot increase capacity.
160 void grow_pod(void *FirstEl, size_t MinSize, size_t TSize);
161
162 /// Report that MinSize doesn't fit into this vector's size type. Throws
163 /// std::length_error or calls report_fatal_error.
164 static void report_size_overflow(size_t MinSize);
165 /// Report that this vector is already at maximum capacity. Throws
166 /// std::length_error or calls report_fatal_error.
167 static void report_at_maximum_capacity();
168
169 /// If true, the RVec is in "memory adoption" mode, i.e. it is acting as a view on a memory buffer it does not own.
170 bool Owns() const { return fCapacity != -1; }
171
172public:
173 size_t size() const { return fSize; }
174 size_t capacity() const noexcept { return Owns() ? fCapacity : fSize; }
175
176 R__RVEC_NODISCARD bool empty() const { return !fSize; }
177
178 /// Set the array size to \p N, which the current array must have enough
179 /// capacity for.
180 ///
181 /// This does not construct or destroy any elements in the vector.
182 ///
183 /// Clients can use this in conjunction with capacity() to write past the end
184 /// of the buffer when they know that more elements are available, and only
185 /// update the size later. This avoids the cost of value initializing elements
186 /// which will only be overwritten.
187 void set_size(size_t N)
188 {
189 if (N > capacity()) {
190 throw std::runtime_error("Setting size to a value greater than capacity.");
191 }
192 fSize = N;
193 }
194};
195
196/// Used to figure out the offset of the first element of an RVec
197template <class T>
199 alignas(SmallVectorBase) char Base[sizeof(SmallVectorBase)];
200 alignas(T) char FirstEl[sizeof(T)];
201};
202
203/// This is the part of SmallVectorTemplateBase which does not depend on whether the type T is a POD.
204template <typename T>
205class R__CLING_PTRCHECK(off) SmallVectorTemplateCommon : public SmallVectorBase {
207
208 /// Find the address of the first element. For this pointer math to be valid
209 /// with small-size of 0 for T with lots of alignment, it's important that
210 /// SmallVectorStorage is properly-aligned even for small-size of 0.
211 void *getFirstEl() const
212 {
213 return const_cast<void *>(reinterpret_cast<const void *>(reinterpret_cast<const char *>(this) +
214 offsetof(SmallVectorAlignmentAndSize<T>, FirstEl)));
215 }
216 // Space after 'FirstEl' is clobbered, do not add any instance vars after it.
217
218protected:
219 SmallVectorTemplateCommon(size_t Size) : Base(getFirstEl(), Size) {}
220
221 void grow_pod(size_t MinSize, size_t TSize) { Base::grow_pod(getFirstEl(), MinSize, TSize); }
222
223 /// Return true if this is a smallvector which has not had dynamic
224 /// memory allocated for it.
225 bool isSmall() const { return this->fBeginX == getFirstEl(); }
226
227 /// Put this vector in a state of being small.
229 {
230 this->fBeginX = getFirstEl();
231 // from the original LLVM implementation:
232 // FIXME: Setting fCapacity to 0 is suspect.
233 this->fSize = this->fCapacity = 0;
234 }
235
236public:
237 // note that fSize is a _signed_ integer, but we expose it as an unsigned integer for consistency with STL containers
238 // as well as backward-compatibility
239 using size_type = size_t;
240 using difference_type = ptrdiff_t;
241 using value_type = T;
242 using iterator = T *;
243 using const_iterator = const T *;
244
245 using const_reverse_iterator = std::reverse_iterator<const_iterator>;
246 using reverse_iterator = std::reverse_iterator<iterator>;
247
248 using reference = T &;
249 using const_reference = const T &;
250 using pointer = T *;
251 using const_pointer = const T *;
252
253 using Base::capacity;
254 using Base::empty;
256
257 // forward iterator creation methods.
258 iterator begin() noexcept { return (iterator)this->fBeginX; }
259 const_iterator begin() const noexcept { return (const_iterator)this->fBeginX; }
260 const_iterator cbegin() const noexcept { return (const_iterator)this->fBeginX; }
261 iterator end() noexcept { return begin() + size(); }
262 const_iterator end() const noexcept { return begin() + size(); }
263 const_iterator cend() const noexcept { return begin() + size(); }
264
265 // reverse iterator creation methods.
266 reverse_iterator rbegin() noexcept { return reverse_iterator(end()); }
267 const_reverse_iterator rbegin() const noexcept { return const_reverse_iterator(end()); }
268 const_reverse_iterator crbegin() const noexcept { return const_reverse_iterator(end()); }
269 reverse_iterator rend() noexcept { return reverse_iterator(begin()); }
270 const_reverse_iterator rend() const noexcept { return const_reverse_iterator(begin()); }
271 const_reverse_iterator crend() const noexcept { return const_reverse_iterator(begin()); }
272
273 size_type size_in_bytes() const { return size() * sizeof(T); }
274 size_type max_size() const noexcept { return std::min(this->SizeTypeMax(), size_type(-1) / sizeof(T)); }
275
276 size_t capacity_in_bytes() const { return capacity() * sizeof(T); }
277
278 /// Return a pointer to the vector's buffer, even if empty().
279 pointer data() noexcept { return pointer(begin()); }
280 /// Return a pointer to the vector's buffer, even if empty().
281 const_pointer data() const noexcept { return const_pointer(begin()); }
282
284 {
285 if (empty()) {
286 throw std::runtime_error("`front` called on an empty RVec");
287 }
288 return begin()[0];
289 }
290
292 {
293 if (empty()) {
294 throw std::runtime_error("`front` called on an empty RVec");
295 }
296 return begin()[0];
297 }
298
300 {
301 if (empty()) {
302 throw std::runtime_error("`back` called on an empty RVec");
303 }
304 return end()[-1];
305 }
306
308 {
309 if (empty()) {
310 throw std::runtime_error("`back` called on an empty RVec");
311 }
312 return end()[-1];
313 }
314};
315
316/// SmallVectorTemplateBase<TriviallyCopyable = false> - This is where we put
317/// method implementations that are designed to work with non-trivial T's.
318///
319/// We approximate is_trivially_copyable with trivial move/copy construction and
320/// trivial destruction. While the standard doesn't specify that you're allowed
321/// copy these types with memcpy, there is no way for the type to observe this.
322/// This catches the important case of std::pair<POD, POD>, which is not
323/// trivially assignable.
324template <typename T, bool = (std::is_trivially_copy_constructible<T>::value) &&
325 (std::is_trivially_move_constructible<T>::value) &&
326 std::is_trivially_destructible<T>::value>
327class R__CLING_PTRCHECK(off) SmallVectorTemplateBase : public SmallVectorTemplateCommon<T> {
328protected:
330
331 static void destroy_range(T *S, T *E)
332 {
333 while (S != E) {
334 --E;
335 E->~T();
336 }
337 }
338
339 /// Move the range [I, E) into the uninitialized memory starting with "Dest",
340 /// constructing elements as needed.
341 template <typename It1, typename It2>
342 static void uninitialized_move(It1 I, It1 E, It2 Dest)
343 {
344 std::uninitialized_copy(std::make_move_iterator(I), std::make_move_iterator(E), Dest);
345 }
346
347 /// Copy the range [I, E) onto the uninitialized memory starting with "Dest",
348 /// constructing elements as needed.
349 template <typename It1, typename It2>
350 static void uninitialized_copy(It1 I, It1 E, It2 Dest)
351 {
352 std::uninitialized_copy(I, E, Dest);
353 }
354
355 /// Grow the allocated memory (without initializing new elements), doubling
356 /// the size of the allocated memory. Guarantees space for at least one more
357 /// element, or MinSize more elements if specified.
358 void grow(size_t MinSize = 0);
359
360public:
361 void push_back(const T &Elt)
362 {
363 if (R__unlikely(this->size() >= this->capacity()))
364 this->grow();
365 ::new ((void *)this->end()) T(Elt);
366 this->set_size(this->size() + 1);
367 }
368
369 void push_back(T &&Elt)
370 {
371 if (R__unlikely(this->size() >= this->capacity()))
372 this->grow();
373 ::new ((void *)this->end()) T(::std::move(Elt));
374 this->set_size(this->size() + 1);
375 }
376
377 void pop_back()
378 {
379 this->set_size(this->size() - 1);
380 this->end()->~T();
381 }
382};
383
384// Define this out-of-line to dissuade the C++ compiler from inlining it.
385template <typename T, bool TriviallyCopyable>
387{
388 // Ensure we can fit the new capacity.
389 // This is only going to be applicable when the capacity is 32 bit.
390 if (MinSize > this->SizeTypeMax())
391 this->report_size_overflow(MinSize);
392
393 // Ensure we can meet the guarantee of space for at least one more element.
394 // The above check alone will not catch the case where grow is called with a
395 // default MinSize of 0, but the current capacity cannot be increased.
396 // This is only going to be applicable when the capacity is 32 bit.
397 if (this->capacity() == this->SizeTypeMax())
398 this->report_at_maximum_capacity();
399
400 // Always grow, even from zero.
401 size_t NewCapacity = size_t(NextPowerOf2(this->capacity() + 2));
402 NewCapacity = std::min(std::max(NewCapacity, MinSize), this->SizeTypeMax());
403 T *NewElts = static_cast<T *>(malloc(NewCapacity * sizeof(T)));
404 R__ASSERT(NewElts != nullptr);
405
406 // Move the elements over.
407 this->uninitialized_move(this->begin(), this->end(), NewElts);
408
409 if (this->Owns()) {
410 // Destroy the original elements.
411 destroy_range(this->begin(), this->end());
412
413 // If this wasn't grown from the inline copy, deallocate the old space.
414 if (!this->isSmall())
415 free(this->begin());
416 }
417
418 this->fBeginX = NewElts;
419 this->fCapacity = NewCapacity;
420}
421
422/// SmallVectorTemplateBase<TriviallyCopyable = true> - This is where we put
423/// method implementations that are designed to work with trivially copyable
424/// T's. This allows using memcpy in place of copy/move construction and
425/// skipping destruction.
426template <typename T>
427class R__CLING_PTRCHECK(off) SmallVectorTemplateBase<T, true> : public SmallVectorTemplateCommon<T> {
429
430protected:
432
433 // No need to do a destroy loop for POD's.
434 static void destroy_range(T *, T *) {}
435
436 /// Move the range [I, E) onto the uninitialized memory
437 /// starting with "Dest", constructing elements into it as needed.
438 template <typename It1, typename It2>
439 static void uninitialized_move(It1 I, It1 E, It2 Dest)
440 {
441 // Just do a copy.
442 uninitialized_copy(I, E, Dest);
443 }
444
445 /// Copy the range [I, E) onto the uninitialized memory
446 /// starting with "Dest", constructing elements into it as needed.
447 template <typename It1, typename It2>
448 static void uninitialized_copy(It1 I, It1 E, It2 Dest)
449 {
450 // Arbitrary iterator types; just use the basic implementation.
451 std::uninitialized_copy(I, E, Dest);
452 }
453
454 /// Copy the range [I, E) onto the uninitialized memory
455 /// starting with "Dest", constructing elements into it as needed.
456 template <typename T1, typename T2>
458 T1 *I, T1 *E, T2 *Dest,
459 typename std::enable_if<std::is_same<typename std::remove_const<T1>::type, T2>::value>::type * = nullptr)
460 {
461 // Use memcpy for PODs iterated by pointers (which includes SmallVector
462 // iterators): std::uninitialized_copy optimizes to memmove, but we can
463 // use memcpy here. Note that I and E are iterators and thus might be
464 // invalid for memcpy if they are equal.
465 if (I != E)
466 memcpy(reinterpret_cast<void *>(Dest), I, (E - I) * sizeof(T));
467 }
468
469 /// Double the size of the allocated memory, guaranteeing space for at
470 /// least one more element or MinSize if specified.
471 void grow(size_t MinSize = 0)
472 {
473 this->grow_pod(MinSize, sizeof(T));
474 }
475
476public:
481
482 void push_back(const T &Elt)
483 {
484 if (R__unlikely(this->size() >= this->capacity()))
485 this->grow();
486 memcpy(reinterpret_cast<void *>(this->end()), &Elt, sizeof(T));
487 this->set_size(this->size() + 1);
488 }
489
490 void pop_back() { this->set_size(this->size() - 1); }
491};
492
493/// Storage for the SmallVector elements. This is specialized for the N=0 case
494/// to avoid allocating unnecessary storage.
495template <typename T, unsigned N>
497 alignas(T) char InlineElts[N * sizeof(T)]{};
498};
499
500/// We need the storage to be properly aligned even for small-size of 0 so that
501/// the pointer math in \a SmallVectorTemplateCommon::getFirstEl() is
502/// well-defined.
503template <typename T>
504struct alignas(T) SmallVectorStorage<T, 0> {
505};
506
507/// The size of the inline storage of an RVec.
508/// Our policy is to allocate at least 8 elements (or more if they all fit into one cacheline)
509/// unless the size of the buffer with 8 elements would be over a certain maximum size.
510template <typename T>
512private:
513#ifdef R__HAS_HARDWARE_INTERFERENCE_SIZE
514 constexpr std::size_t cacheLineSize = std::hardware_destructive_interference_size;
515#else
516 // safe bet: assume the typical 64 bytes
517 static constexpr std::size_t cacheLineSize = 64;
518#endif
519 static constexpr unsigned elementsPerCacheLine = (cacheLineSize - sizeof(SmallVectorBase)) / sizeof(T);
520 static constexpr unsigned maxInlineByteSize = 1024;
521
522public:
523 static constexpr unsigned value =
524 elementsPerCacheLine >= 8 ? elementsPerCacheLine : (sizeof(T) * 8 > maxInlineByteSize ? 0 : 8);
525};
526
527} // namespace VecOps
528} // namespace Internal
529
530namespace Detail {
531namespace VecOps {
532
533/// This class consists of common code factored out of the SmallVector class to
534/// reduce code duplication based on the SmallVector 'N' template parameter.
535template <typename T>
536class R__CLING_PTRCHECK(off) RVecImpl : public Internal::VecOps::SmallVectorTemplateBase<T> {
538
539public:
544
545protected:
546 // Default ctor - Initialize to empty.
547 explicit RVecImpl(unsigned N) : ROOT::Internal::VecOps::SmallVectorTemplateBase<T>(N) {}
548
549public:
550 RVecImpl(const RVecImpl &) = delete;
551
553 {
554 // Subclass has already destructed this vector's elements.
555 // If this wasn't grown from the inline copy, deallocate the old space.
556 if (!this->isSmall() && this->Owns())
557 free(this->begin());
558 }
559
560 // also give up adopted memory if applicable
561 void clear()
562 {
563 if (this->Owns()) {
564 this->destroy_range(this->begin(), this->end());
565 this->fSize = 0;
566 } else {
567 this->resetToSmall();
568 }
569 }
570
572 {
573 if (N < this->size()) {
574 if (this->Owns())
575 this->destroy_range(this->begin() + N, this->end());
576 this->set_size(N);
577 } else if (N > this->size()) {
578 if (this->capacity() < N)
579 this->grow(N);
580 for (auto I = this->end(), E = this->begin() + N; I != E; ++I)
581 new (&*I) T();
582 this->set_size(N);
583 }
584 }
585
586 void resize(size_type N, const T &NV)
587 {
588 if (N < this->size()) {
589 if (this->Owns())
590 this->destroy_range(this->begin() + N, this->end());
591 this->set_size(N);
592 } else if (N > this->size()) {
593 if (this->capacity() < N)
594 this->grow(N);
595 std::uninitialized_fill(this->end(), this->begin() + N, NV);
596 this->set_size(N);
597 }
598 }
599
601 {
602 if (this->capacity() < N)
603 this->grow(N);
604 }
605
606 void pop_back_n(size_type NumItems)
607 {
608 if (this->size() < NumItems) {
609 throw std::runtime_error("Popping back more elements than those available.");
610 }
611 if (this->Owns())
612 this->destroy_range(this->end() - NumItems, this->end());
613 this->set_size(this->size() - NumItems);
614 }
615
617 {
618 T Result = ::std::move(this->back());
619 this->pop_back();
620 return Result;
621 }
622
623 void swap(RVecImpl &RHS);
624
625 /// Add the specified range to the end of the SmallVector.
626 template <typename in_iter,
627 typename = typename std::enable_if<std::is_convertible<
628 typename std::iterator_traits<in_iter>::iterator_category, std::input_iterator_tag>::value>::type>
629 void append(in_iter in_start, in_iter in_end)
630 {
631 size_type NumInputs = std::distance(in_start, in_end);
632 if (NumInputs > this->capacity() - this->size())
633 this->grow(this->size() + NumInputs);
634
635 this->uninitialized_copy(in_start, in_end, this->end());
636 this->set_size(this->size() + NumInputs);
637 }
638
639 /// Append \p NumInputs copies of \p Elt to the end.
640 void append(size_type NumInputs, const T &Elt)
641 {
642 if (NumInputs > this->capacity() - this->size())
643 this->grow(this->size() + NumInputs);
644
645 std::uninitialized_fill_n(this->end(), NumInputs, Elt);
646 this->set_size(this->size() + NumInputs);
647 }
648
649 void append(std::initializer_list<T> IL) { append(IL.begin(), IL.end()); }
650
651 // from the original LLVM implementation:
652 // FIXME: Consider assigning over existing elements, rather than clearing &
653 // re-initializing them - for all assign(...) variants.
654
655 void assign(size_type NumElts, const T &Elt)
656 {
657 clear();
658 if (this->capacity() < NumElts)
659 this->grow(NumElts);
660 this->set_size(NumElts);
661 std::uninitialized_fill(this->begin(), this->end(), Elt);
662 }
663
664 template <typename in_iter,
665 typename = typename std::enable_if<std::is_convertible<
666 typename std::iterator_traits<in_iter>::iterator_category, std::input_iterator_tag>::value>::type>
667 void assign(in_iter in_start, in_iter in_end)
668 {
669 clear();
670 append(in_start, in_end);
671 }
672
673 void assign(std::initializer_list<T> IL)
674 {
675 clear();
676 append(IL);
677 }
678
680 {
681 // Just cast away constness because this is a non-const member function.
682 iterator I = const_cast<iterator>(CI);
683
684 if (I < this->begin() || I >= this->end()) {
685 throw std::runtime_error("The iterator passed to `erase` is out of bounds.");
686 }
687
688 iterator N = I;
689 // Shift all elts down one.
690 std::move(I + 1, this->end(), I);
691 // Drop the last elt.
692 this->pop_back();
693 return (N);
694 }
695
697 {
698 // Just cast away constness because this is a non-const member function.
699 iterator S = const_cast<iterator>(CS);
700 iterator E = const_cast<iterator>(CE);
701
702 if (S < this->begin() || E > this->end() || S > E) {
703 throw std::runtime_error("Invalid start/end pair passed to `erase` (out of bounds or start > end).");
704 }
705
706 iterator N = S;
707 // Shift all elts down.
708 iterator I = std::move(E, this->end(), S);
709 // Drop the last elts.
710 if (this->Owns())
711 this->destroy_range(I, this->end());
712 this->set_size(I - this->begin());
713 return (N);
714 }
715
717 {
718 if (I == this->end()) { // Important special case for empty vector.
719 this->push_back(::std::move(Elt));
720 return this->end() - 1;
721 }
722
723 if (I < this->begin() || I > this->end()) {
724 throw std::runtime_error("The iterator passed to `insert` is out of bounds.");
725 }
726
727 if (this->size() >= this->capacity()) {
728 size_t EltNo = I - this->begin();
729 this->grow();
730 I = this->begin() + EltNo;
731 }
732
733 ::new ((void *)this->end()) T(::std::move(this->back()));
734 // Push everything else over.
735 std::move_backward(I, this->end() - 1, this->end());
736 this->set_size(this->size() + 1);
737
738 // If we just moved the element we're inserting, be sure to update
739 // the reference.
740 T *EltPtr = &Elt;
741 if (I <= EltPtr && EltPtr < this->end())
742 ++EltPtr;
743
744 *I = ::std::move(*EltPtr);
745 return I;
746 }
747
749 {
750 if (I == this->end()) { // Important special case for empty vector.
751 this->push_back(Elt);
752 return this->end() - 1;
753 }
754
755 if (I < this->begin() || I > this->end()) {
756 throw std::runtime_error("The iterator passed to `insert` is out of bounds.");
757 }
758
759 if (this->size() >= this->capacity()) {
760 size_t EltNo = I - this->begin();
761 this->grow();
762 I = this->begin() + EltNo;
763 }
764 ::new ((void *)this->end()) T(std::move(this->back()));
765 // Push everything else over.
766 std::move_backward(I, this->end() - 1, this->end());
767 this->set_size(this->size() + 1);
768
769 // If we just moved the element we're inserting, be sure to update
770 // the reference.
771 const T *EltPtr = &Elt;
772 if (I <= EltPtr && EltPtr < this->end())
773 ++EltPtr;
774
775 *I = *EltPtr;
776 return I;
777 }
778
779 iterator insert(iterator I, size_type NumToInsert, const T &Elt)
780 {
781 // Convert iterator to elt# to avoid invalidating iterator when we reserve()
782 size_t InsertElt = I - this->begin();
783
784 if (I == this->end()) { // Important special case for empty vector.
785 append(NumToInsert, Elt);
786 return this->begin() + InsertElt;
787 }
788
789 if (I < this->begin() || I > this->end()) {
790 throw std::runtime_error("The iterator passed to `insert` is out of bounds.");
791 }
792
793 // Ensure there is enough space.
794 reserve(this->size() + NumToInsert);
795
796 // Uninvalidate the iterator.
797 I = this->begin() + InsertElt;
798
799 // If there are more elements between the insertion point and the end of the
800 // range than there are being inserted, we can use a simple approach to
801 // insertion. Since we already reserved space, we know that this won't
802 // reallocate the vector.
803 if (size_t(this->end() - I) >= NumToInsert) {
804 T *OldEnd = this->end();
805 append(std::move_iterator<iterator>(this->end() - NumToInsert), std::move_iterator<iterator>(this->end()));
806
807 // Copy the existing elements that get replaced.
808 std::move_backward(I, OldEnd - NumToInsert, OldEnd);
809
810 std::fill_n(I, NumToInsert, Elt);
811 return I;
812 }
813
814 // Otherwise, we're inserting more elements than exist already, and we're
815 // not inserting at the end.
816
817 // Move over the elements that we're about to overwrite.
818 T *OldEnd = this->end();
819 this->set_size(this->size() + NumToInsert);
820 size_t NumOverwritten = OldEnd - I;
821 this->uninitialized_move(I, OldEnd, this->end() - NumOverwritten);
822
823 // Replace the overwritten part.
824 std::fill_n(I, NumOverwritten, Elt);
825
826 // Insert the non-overwritten middle part.
827 std::uninitialized_fill_n(OldEnd, NumToInsert - NumOverwritten, Elt);
828 return I;
829 }
830
831 template <typename ItTy,
832 typename = typename std::enable_if<std::is_convertible<
833 typename std::iterator_traits<ItTy>::iterator_category, std::input_iterator_tag>::value>::type>
834 iterator insert(iterator I, ItTy From, ItTy To)
835 {
836 // Convert iterator to elt# to avoid invalidating iterator when we reserve()
837 size_t InsertElt = I - this->begin();
838
839 if (I == this->end()) { // Important special case for empty vector.
840 append(From, To);
841 return this->begin() + InsertElt;
842 }
843
844 if (I < this->begin() || I > this->end()) {
845 throw std::runtime_error("The iterator passed to `insert` is out of bounds.");
846 }
847
848 size_t NumToInsert = std::distance(From, To);
849
850 // Ensure there is enough space.
851 reserve(this->size() + NumToInsert);
852
853 // Uninvalidate the iterator.
854 I = this->begin() + InsertElt;
855
856 // If there are more elements between the insertion point and the end of the
857 // range than there are being inserted, we can use a simple approach to
858 // insertion. Since we already reserved space, we know that this won't
859 // reallocate the vector.
860 if (size_t(this->end() - I) >= NumToInsert) {
861 T *OldEnd = this->end();
862 append(std::move_iterator<iterator>(this->end() - NumToInsert), std::move_iterator<iterator>(this->end()));
863
864 // Copy the existing elements that get replaced.
865 std::move_backward(I, OldEnd - NumToInsert, OldEnd);
866
867 std::copy(From, To, I);
868 return I;
869 }
870
871 // Otherwise, we're inserting more elements than exist already, and we're
872 // not inserting at the end.
873
874 // Move over the elements that we're about to overwrite.
875 T *OldEnd = this->end();
876 this->set_size(this->size() + NumToInsert);
877 size_t NumOverwritten = OldEnd - I;
878 this->uninitialized_move(I, OldEnd, this->end() - NumOverwritten);
879
880 // Replace the overwritten part.
881 for (T *J = I; NumOverwritten > 0; --NumOverwritten) {
882 *J = *From;
883 ++J;
884 ++From;
885 }
886
887 // Insert the non-overwritten middle part.
888 this->uninitialized_copy(From, To, OldEnd);
889 return I;
890 }
891
892 void insert(iterator I, std::initializer_list<T> IL) { insert(I, IL.begin(), IL.end()); }
893
894 template <typename... ArgTypes>
895 reference emplace_back(ArgTypes &&...Args)
896 {
897 if (R__unlikely(this->size() >= this->capacity()))
898 this->grow();
899 ::new ((void *)this->end()) T(std::forward<ArgTypes>(Args)...);
900 this->set_size(this->size() + 1);
901 return this->back();
902 }
903
905
907};
908
909template <typename T>
911{
912 if (this == &RHS)
913 return;
914
915 // We can only avoid copying elements if neither vector is small.
916 if (!this->isSmall() && !RHS.isSmall()) {
917 std::swap(this->fBeginX, RHS.fBeginX);
918 std::swap(this->fSize, RHS.fSize);
919 std::swap(this->fCapacity, RHS.fCapacity);
920 return;
921 }
922
923 // This block handles the swap of a small and a non-owning vector
924 // It is more efficient to first move the non-owning vector, hence the 2 cases
925 if (this->isSmall() && !RHS.Owns()) { // the right vector is non-owning
926 RVecImpl<T> temp(0);
927 temp = std::move(RHS);
928 RHS = std::move(*this);
929 *this = std::move(temp);
930 return;
931 } else if (RHS.isSmall() && !this->Owns()) { // the left vector is non-owning
932 RVecImpl<T> temp(0);
933 temp = std::move(*this);
934 *this = std::move(RHS);
935 RHS = std::move(temp);
936 return;
937 }
938
939 if (RHS.size() > this->capacity())
940 this->grow(RHS.size());
941 if (this->size() > RHS.capacity())
942 RHS.grow(this->size());
943
944 // Swap the shared elements.
945 size_t NumShared = this->size();
946 if (NumShared > RHS.size())
947 NumShared = RHS.size();
948 for (size_type i = 0; i != NumShared; ++i)
949 std::iter_swap(this->begin() + i, RHS.begin() + i);
950
951 // Copy over the extra elts.
952 if (this->size() > RHS.size()) {
953 size_t EltDiff = this->size() - RHS.size();
954 this->uninitialized_copy(this->begin() + NumShared, this->end(), RHS.end());
955 RHS.set_size(RHS.size() + EltDiff);
956 if (this->Owns())
957 this->destroy_range(this->begin() + NumShared, this->end());
958 this->set_size(NumShared);
959 } else if (RHS.size() > this->size()) {
960 size_t EltDiff = RHS.size() - this->size();
961 this->uninitialized_copy(RHS.begin() + NumShared, RHS.end(), this->end());
962 this->set_size(this->size() + EltDiff);
963 if (RHS.Owns())
964 this->destroy_range(RHS.begin() + NumShared, RHS.end());
965 RHS.set_size(NumShared);
966 }
967}
968
969template <typename T>
971{
972 // Avoid self-assignment.
973 if (this == &RHS)
974 return *this;
975
976 // If we already have sufficient space, assign the common elements, then
977 // destroy any excess.
978 size_t RHSSize = RHS.size();
979 size_t CurSize = this->size();
980 if (CurSize >= RHSSize) {
981 // Assign common elements.
982 iterator NewEnd;
983 if (RHSSize)
984 NewEnd = std::copy(RHS.begin(), RHS.begin() + RHSSize, this->begin());
985 else
986 NewEnd = this->begin();
987
988 // Destroy excess elements.
989 if (this->Owns())
990 this->destroy_range(NewEnd, this->end());
991
992 // Trim.
993 this->set_size(RHSSize);
994 return *this;
995 }
996
997 // If we have to grow to have enough elements, destroy the current elements.
998 // This allows us to avoid copying them during the grow.
999 // From the original LLVM implementation:
1000 // FIXME: don't do this if they're efficiently moveable.
1001 if (this->capacity() < RHSSize) {
1002 if (this->Owns()) {
1003 // Destroy current elements.
1004 this->destroy_range(this->begin(), this->end());
1005 }
1006 this->set_size(0);
1007 CurSize = 0;
1008 this->grow(RHSSize);
1009 } else if (CurSize) {
1010 // Otherwise, use assignment for the already-constructed elements.
1011 std::copy(RHS.begin(), RHS.begin() + CurSize, this->begin());
1012 }
1013
1014 // Copy construct the new elements in place.
1015 this->uninitialized_copy(RHS.begin() + CurSize, RHS.end(), this->begin() + CurSize);
1016
1017 // Set end.
1018 this->set_size(RHSSize);
1019 return *this;
1020}
1021
1022template <typename T>
1024{
1025 // Avoid self-assignment.
1026 if (this == &RHS)
1027 return *this;
1028
1029 // If the RHS isn't small, clear this vector and then steal its buffer.
1030 if (!RHS.isSmall()) {
1031 if (this->Owns()) {
1032 this->destroy_range(this->begin(), this->end());
1033 if (!this->isSmall())
1034 free(this->begin());
1035 }
1036 this->fBeginX = RHS.fBeginX;
1037 this->fSize = RHS.fSize;
1038 this->fCapacity = RHS.fCapacity;
1039 RHS.resetToSmall();
1040 return *this;
1041 }
1042
1043 // If we already have sufficient space, assign the common elements, then
1044 // destroy any excess.
1045 size_t RHSSize = RHS.size();
1046 size_t CurSize = this->size();
1047 if (CurSize >= RHSSize) {
1048 // Assign common elements.
1049 iterator NewEnd = this->begin();
1050 if (RHSSize)
1051 NewEnd = std::move(RHS.begin(), RHS.end(), NewEnd);
1052
1053 // Destroy excess elements and trim the bounds.
1054 if (this->Owns())
1055 this->destroy_range(NewEnd, this->end());
1056 this->set_size(RHSSize);
1057
1058 // Clear the RHS.
1059 RHS.clear();
1060
1061 return *this;
1062 }
1063
1064 // If we have to grow to have enough elements, destroy the current elements.
1065 // This allows us to avoid copying them during the grow.
1066 // From the original LLVM implementation:
1067 // FIXME: this may not actually make any sense if we can efficiently move
1068 // elements.
1069 if (this->capacity() < RHSSize) {
1070 if (this->Owns()) {
1071 // Destroy current elements.
1072 this->destroy_range(this->begin(), this->end());
1073 }
1074 this->set_size(0);
1075 CurSize = 0;
1076 this->grow(RHSSize);
1077 } else if (CurSize) {
1078 // Otherwise, use assignment for the already-constructed elements.
1079 std::move(RHS.begin(), RHS.begin() + CurSize, this->begin());
1080 }
1081
1082 // Move-construct the new elements in place.
1083 this->uninitialized_move(RHS.begin() + CurSize, RHS.end(), this->begin() + CurSize);
1084
1085 // Set end.
1086 this->set_size(RHSSize);
1087
1088 RHS.clear();
1089 return *this;
1090}
1091
1092template <typename T>
1094{
1095 return v.isSmall();
1096}
1097
1098template <typename T>
1100{
1101 return !v.Owns();
1102}
1103
1104} // namespace VecOps
1105} // namespace Detail
1106
1107namespace VecOps {
1108// Note that we open here with @{ the Doxygen group vecops and it is
1109// closed again at the end of the C++ namespace VecOps
1110/**
1111 * \defgroup vecops VecOps
1112 * A "std::vector"-like collection of values implementing handy operation to analyse them
1113 * @{
1114*/
1115
1116// From the original SmallVector code:
1117// This is a 'vector' (really, a variable-sized array), optimized
1118// for the case when the array is small. It contains some number of elements
1119// in-place, which allows it to avoid heap allocation when the actual number of
1120// elements is below that threshold. This allows normal "small" cases to be
1121// fast without losing generality for large inputs.
1122//
1123// Note that this does not attempt to be exception safe.
1124
1125template <typename T, unsigned int N>
1126class R__CLING_PTRCHECK(off) RVecN : public Detail::VecOps::RVecImpl<T>, Internal::VecOps::SmallVectorStorage<T, N> {
1127public:
1128 RVecN() : Detail::VecOps::RVecImpl<T>(N) {}
1129
1131 {
1132 if (this->Owns()) {
1133 // Destroy the constructed elements in the vector.
1134 this->destroy_range(this->begin(), this->end());
1135 }
1136 }
1137
1138 explicit RVecN(size_t Size, const T &Value) : Detail::VecOps::RVecImpl<T>(N) { this->assign(Size, Value); }
1139
1140 explicit RVecN(size_t Size) : Detail::VecOps::RVecImpl<T>(N)
1141 {
1142 if (Size > N)
1143 this->grow(Size);
1144 this->fSize = Size;
1145 std::uninitialized_fill(this->begin(), this->end(), T{});
1146 }
1147
1148 template <typename ItTy,
1149 typename = typename std::enable_if<std::is_convertible<
1150 typename std::iterator_traits<ItTy>::iterator_category, std::input_iterator_tag>::value>::type>
1151 RVecN(ItTy S, ItTy E) : Detail::VecOps::RVecImpl<T>(N)
1152 {
1153 this->append(S, E);
1154 }
1155
1156 RVecN(std::initializer_list<T> IL) : Detail::VecOps::RVecImpl<T>(N) { this->assign(IL); }
1157
1158 RVecN(const RVecN &RHS) : Detail::VecOps::RVecImpl<T>(N)
1159 {
1160 if (!RHS.empty())
1162 }
1163
1164 RVecN &operator=(const RVecN &RHS)
1165 {
1167 return *this;
1168 }
1169
1170 RVecN(RVecN &&RHS) : Detail::VecOps::RVecImpl<T>(N)
1171 {
1172 if (!RHS.empty())
1174 }
1175
1176 RVecN(Detail::VecOps::RVecImpl<T> &&RHS) : Detail::VecOps::RVecImpl<T>(N)
1177 {
1178 if (!RHS.empty())
1180 }
1181
1182 RVecN(const std::vector<T> &RHS) : RVecN(RHS.begin(), RHS.end()) {}
1183
1185 {
1187 return *this;
1188 }
1189
1190 RVecN(T* p, size_t n) : Detail::VecOps::RVecImpl<T>(N)
1191 {
1192 this->fBeginX = p;
1193 this->fSize = n;
1194 this->fCapacity = -1;
1195 }
1196
1198 {
1200 return *this;
1201 }
1202
1203 RVecN &operator=(std::initializer_list<T> IL)
1204 {
1205 this->assign(IL);
1206 return *this;
1207 }
1208
1215
1217 {
1218 return begin()[idx];
1219 }
1220
1222 {
1223 return begin()[idx];
1224 }
1225
1227 RVecN operator[](const RVecN<V, M> &conds) const
1228 {
1229 const size_type n = conds.size();
1230
1231 if (n != this->size())
1232 throw std::runtime_error("Cannot index RVecN with condition vector of different size");
1233
1234 RVecN ret;
1235 ret.reserve(n);
1236 size_type j = 0u;
1237 for (size_type i = 0u; i < n; ++i) {
1238 if (conds[i]) {
1239 // the begin() is to go around the R__ASSERT in operator[]
1240 ret.begin()[j] = this->operator[](i);
1241 ++j;
1242 }
1243 }
1244 ret.set_size(j);
1245 return ret;
1246 }
1247
1248 // conversion
1250 operator RVecN<U, M>() const
1251 {
1252 return RVecN<U, M>(this->begin(), this->end());
1253 }
1254
1256 {
1257 if (pos >= size_type(this->fSize))
1258 throw std::out_of_range("RVecN");
1259 return this->operator[](pos);
1260 }
1261
1263 {
1264 if (pos >= size_type(this->fSize))
1265 throw std::out_of_range("RVecN");
1266 return this->operator[](pos);
1267 }
1268
1269 /// No exception thrown. The user specifies the desired value in case the RVecN is shorter than `pos`.
1271 {
1272 if (pos >= size_type(this->fSize))
1273 return fallback;
1274 return this->operator[](pos);
1275 }
1276
1277 /// No exception thrown. The user specifies the desired value in case the RVecN is shorter than `pos`.
1278 value_type at(size_type pos, value_type fallback) const
1279 {
1280 if (pos >= size_type(this->fSize))
1281 return fallback;
1282 return this->operator[](pos);
1283 }
1284};
1285
1286// clang-format off
1287/**
1288\class ROOT::VecOps::RVec
1289\brief A "std::vector"-like collection of values implementing handy operation to analyse them
1290\tparam T The type of the contained objects
1291
1292A RVec is a container designed to make analysis of values' collections fast and easy.
1293Its storage is contiguous in memory and its interface is designed such to resemble to the one
1294of the stl vector. In addition the interface features methods and
1295[external functions](https://root.cern/doc/master/namespaceROOT_1_1VecOps.html) to ease the manipulation and analysis
1296of the data in the RVec.
1297
1298\note ROOT::VecOps::RVec can also be spelled simply ROOT::RVec. Shorthand aliases such as ROOT::RVecI or ROOT::RVecD
1299are also available as template instantiations of RVec of fundamental types. The full list of available aliases:
1300- RVecB (`bool`)
1301- RVecC (`char`)
1302- RVecD (`double`)
1303- RVecF (`float`)
1304- RVecI (`int`)
1305- RVecL (`long`)
1306- RVecLL (`long long`)
1307- RVecU (`unsigned`)
1308- RVecUL (`unsigned long`)
1309- RVecULL (`unsigned long long`)
1310
1311\note RVec does not attempt to be exception safe. Exceptions thrown by element constructors during insertions, swaps or
1312other operations will be propagated potentially leaving the RVec object in an invalid state.
1313
1314\note RVec methods (e.g. `at` or `size`) follow the STL naming convention instead of the ROOT naming convention in order
1315to make RVec a drop-in replacement for `std::vector`.
1316
1317\htmlonly
1318<a href="https://doi.org/10.5281/zenodo.1253756"><img src="https://zenodo.org/badge/DOI/10.5281/zenodo.1253756.svg" alt="DOI"></a>
1319\endhtmlonly
1320
1321## Table of Contents
1322- [Example](\ref example)
1323- [Owning and adopting memory](\ref owningandadoptingmemory)
1324- [Sorting and manipulation of indices](\ref sorting)
1325- [Usage in combination with RDataFrame](\ref usagetdataframe)
1326- [Reference for the RVec class](\ref RVecdoxyref)
1327- [Reference for RVec helper functions](https://root.cern/doc/master/namespaceROOT_1_1VecOps.html)
1328
1329\anchor example
1330## Example
1331Suppose to have an event featuring a collection of muons with a certain pseudorapidity,
1332momentum and charge, e.g.:
1333~~~{.cpp}
1334std::vector<short> mu_charge {1, 1, -1, -1, -1, 1, 1, -1};
1335std::vector<float> mu_pt {56, 45, 32, 24, 12, 8, 7, 6.2};
1336std::vector<float> mu_eta {3.1, -.2, -1.1, 1, 4.1, 1.6, 2.4, -.5};
1337~~~
1338Suppose you want to extract the transverse momenta of the muons satisfying certain
1339criteria, for example consider only negatively charged muons with a pseudorapidity
1340smaller or equal to 2 and with a transverse momentum greater than 10 GeV.
1341Such a selection would require, among the other things, the management of an explicit
1342loop, for example:
1343~~~{.cpp}
1344std::vector<float> goodMuons_pt;
1345const auto size = mu_charge.size();
1346for (size_t i=0; i < size; ++i) {
1347 if (mu_pt[i] > 10 && abs(mu_eta[i]) <= 2. && mu_charge[i] == -1) {
1348 goodMuons_pt.emplace_back(mu_pt[i]);
1349 }
1350}
1351~~~
1352These operations become straightforward with RVec - we just need to *write what
1353we mean*:
1354~~~{.cpp}
1355auto goodMuons_pt = mu_pt[ (mu_pt > 10.f && abs(mu_eta) <= 2.f && mu_charge == -1) ]
1356~~~
1357Now the clean collection of transverse momenta can be used within the rest of the data analysis, for
1358example to fill a histogram.
1359
1360\anchor owningandadoptingmemory
1361## Owning and adopting memory
1362RVec has contiguous memory associated to it. It can own it or simply adopt it. In the latter case,
1363it can be constructed with the address of the memory associated to it and its length. For example:
1364~~~{.cpp}
1365std::vector<int> myStlVec {1,2,3};
1366RVec<int> myRVec(myStlVec.data(), myStlVec.size());
1367~~~
1368In this case, the memory associated to myStlVec and myRVec is the same, myRVec simply "adopted it".
1369If any method which implies a re-allocation is called, e.g. *emplace_back* or *resize*, the adopted
1370memory is released and new one is allocated. The previous content is copied in the new memory and
1371preserved.
1372
1373\anchor sorting
1374## Sorting and manipulation of indices
1375
1376### Sorting
1377RVec complies to the STL interfaces when it comes to iterations. As a result, standard algorithms
1378can be used, for example sorting:
1379~~~{.cpp}
1380RVec<double> v{6., 4., 5.};
1381std::sort(v.begin(), v.end());
1382~~~
1383
1384For convenience, helpers are provided too:
1385~~~{.cpp}
1386auto sorted_v = Sort(v);
1387auto reversed_v = Reverse(v);
1388~~~
1389
1390### Manipulation of indices
1391
1392It is also possible to manipulated the RVecs acting on their indices. For example,
1393the following syntax
1394~~~{.cpp}
1395RVecD v0 {9., 7., 8.};
1396auto v1 = Take(v0, {1, 2, 0});
1397~~~
1398will yield a new RVec<double> the content of which is the first, second and zeroth element of
1399v0, i.e. `{7., 8., 9.}`.
1400
1401The `Argsort` and `StableArgsort` helper extracts the indices which order the content of a `RVec`.
1402For example, this snippet accomplishes in a more expressive way what we just achieved:
1403~~~{.cpp}
1404auto v1_indices = Argsort(v0); // The content of v1_indices is {1, 2, 0}.
1405v1 = Take(v0, v1_indices);
1406~~~
1407
1408The `Take` utility allows to extract portions of the `RVec`. The content to be *taken*
1409can be specified with an `RVec` of indices or an integer. If the integer is negative,
1410elements will be picked starting from the end of the container:
1411~~~{.cpp}
1412RVecF vf {1.f, 2.f, 3.f, 4.f};
1413auto vf_1 = Take(vf, {1, 3}); // The content is {2.f, 4.f}
1414auto vf_2 = Take(vf, 2); // The content is {1.f, 2.f}
1415auto vf_3 = Take(vf, -3); // The content is {2.f, 3.f, 4.f}
1416~~~
1417
1418\anchor usagetdataframe
1419## Usage in combination with RDataFrame
1420RDataFrame leverages internally RVecs. Suppose to have a dataset stored in a
1421TTree which holds these columns (here we choose C arrays to represent the
1422collections, they could be as well std::vector instances):
1423~~~{.bash}
1424 nPart "nPart/I" An integer representing the number of particles
1425 px "px[nPart]/D" The C array of the particles' x component of the momentum
1426 py "py[nPart]/D" The C array of the particles' y component of the momentum
1427 E "E[nPart]/D" The C array of the particles' Energy
1428~~~
1429Suppose you'd like to plot in a histogram the transverse momenta of all particles
1430for which the energy is greater than 200 MeV.
1431The code required would just be:
1432~~~{.cpp}
1433RDataFrame d("mytree", "myfile.root");
1434auto cutPt = [](RVecD &pxs, RVecD &pys, RVecD &Es) {
1435 auto all_pts = sqrt(pxs * pxs + pys * pys);
1436 auto good_pts = all_pts[Es > 200.];
1437 return good_pts;
1438 };
1439
1440auto hpt = d.Define("pt", cutPt, {"px", "py", "E"})
1441 .Histo1D("pt");
1442hpt->Draw();
1443~~~
1444And if you'd like to express your selection as a string:
1445~~~{.cpp}
1446RDataFrame d("mytree", "myfile.root");
1447auto hpt = d.Define("pt", "sqrt(pxs * pxs + pys * pys)[E>200]")
1448 .Histo1D("pt");
1449hpt->Draw();
1450~~~
1451\anchor RVecdoxyref
1452**/
1453// clang-format on
1454
1455template <typename T>
1456class R__CLING_PTRCHECK(off) RVec : public RVecN<T, Internal::VecOps::RVecInlineStorageSize<T>::value> {
1458public:
1463 using SuperClass::begin;
1464 using SuperClass::size;
1465
1466 RVec() {}
1467
1468 explicit RVec(size_t Size, const T &Value) : SuperClass(Size, Value) {}
1469
1470 explicit RVec(size_t Size) : SuperClass(Size) {}
1471
1472 template <typename ItTy,
1473 typename = typename std::enable_if<std::is_convertible<
1474 typename std::iterator_traits<ItTy>::iterator_category, std::input_iterator_tag>::value>::type>
1475 RVec(ItTy S, ItTy E) : SuperClass(S, E)
1476 {
1477 }
1478
1479 RVec(std::initializer_list<T> IL) : SuperClass(IL) {}
1480
1481 RVec(const RVec &RHS) : SuperClass(RHS) {}
1482
1483 RVec &operator=(const RVec &RHS)
1484 {
1486 return *this;
1487 }
1488
1489 RVec(RVec &&RHS) : SuperClass(std::move(RHS)) {}
1490
1492 {
1493 SuperClass::operator=(std::move(RHS));
1494 return *this;
1495 }
1496
1497 RVec(Detail::VecOps::RVecImpl<T> &&RHS) : SuperClass(std::move(RHS)) {}
1498
1499 template <unsigned N>
1500 RVec(RVecN<T, N> &&RHS) : SuperClass(std::move(RHS)) {}
1501
1502 template <unsigned N>
1503 RVec(const RVecN<T, N> &RHS) : SuperClass(RHS) {}
1504
1505 RVec(const std::vector<T> &RHS) : SuperClass(RHS) {}
1506
1507 RVec(T* p, size_t n) : SuperClass(p, n) {}
1508
1509 // conversion
1511 operator RVec<U>() const
1512 {
1513 return RVec<U>(this->begin(), this->end());
1514 }
1515
1516 using SuperClass::operator[];
1517
1519 RVec operator[](const RVec<V> &conds) const
1520 {
1521 return RVec(SuperClass::operator[](conds));
1522 }
1523
1524 using SuperClass::at;
1525
1526 friend bool ROOT::Detail::VecOps::IsSmall<T>(const RVec<T> &v);
1527
1528 friend bool ROOT::Detail::VecOps::IsAdopting<T>(const RVec<T> &v);
1529};
1530
1531template <typename T, unsigned N>
1532inline size_t CapacityInBytes(const RVecN<T, N> &X)
1533{
1534 return X.capacity_in_bytes();
1535}
1536
1537///@name RVec Unary Arithmetic Operators
1538///@{
1539
1540#define RVEC_UNARY_OPERATOR(OP) \
1541template <typename T> \
1542RVec<T> operator OP(const RVec<T> &v) \
1543{ \
1544 RVec<T> ret(v); \
1545 for (auto &x : ret) \
1546 x = OP x; \
1547return ret; \
1548} \
1549
1554#undef RVEC_UNARY_OPERATOR
1555
1556///@}
1557///@name RVec Binary Arithmetic Operators
1558///@{
1559
1560#define ERROR_MESSAGE(OP) \
1561 "Cannot call operator " #OP " on vectors of different sizes."
1562
1563#define RVEC_BINARY_OPERATOR(OP) \
1564template <typename T0, typename T1> \
1565auto operator OP(const RVec<T0> &v, const T1 &y) \
1566 -> RVec<decltype(v[0] OP y)> \
1567{ \
1568 RVec<decltype(v[0] OP y)> ret(v.size()); \
1569 auto op = [&y](const T0 &x) { return x OP y; }; \
1570 std::transform(v.begin(), v.end(), ret.begin(), op); \
1571 return ret; \
1572} \
1573 \
1574template <typename T0, typename T1> \
1575auto operator OP(const T0 &x, const RVec<T1> &v) \
1576 -> RVec<decltype(x OP v[0])> \
1577{ \
1578 RVec<decltype(x OP v[0])> ret(v.size()); \
1579 auto op = [&x](const T1 &y) { return x OP y; }; \
1580 std::transform(v.begin(), v.end(), ret.begin(), op); \
1581 return ret; \
1582} \
1583 \
1584template <typename T0, typename T1> \
1585auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
1586 -> RVec<decltype(v0[0] OP v1[0])> \
1587{ \
1588 if (v0.size() != v1.size()) \
1589 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1590 \
1591 RVec<decltype(v0[0] OP v1[0])> ret(v0.size()); \
1592 auto op = [](const T0 &x, const T1 &y) { return x OP y; }; \
1593 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
1594 return ret; \
1595} \
1596
1605#undef RVEC_BINARY_OPERATOR
1606
1607///@}
1608///@name RVec Assignment Arithmetic Operators
1609///@{
1610
1611#define RVEC_ASSIGNMENT_OPERATOR(OP) \
1612template <typename T0, typename T1> \
1613RVec<T0>& operator OP(RVec<T0> &v, const T1 &y) \
1614{ \
1615 auto op = [&y](T0 &x) { return x OP y; }; \
1616 std::transform(v.begin(), v.end(), v.begin(), op); \
1617 return v; \
1618} \
1619 \
1620template <typename T0, typename T1> \
1621RVec<T0>& operator OP(RVec<T0> &v0, const RVec<T1> &v1) \
1622{ \
1623 if (v0.size() != v1.size()) \
1624 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1625 \
1626 auto op = [](T0 &x, const T1 &y) { return x OP y; }; \
1627 std::transform(v0.begin(), v0.end(), v1.begin(), v0.begin(), op); \
1628 return v0; \
1629} \
1630
1641#undef RVEC_ASSIGNMENT_OPERATOR
1642
1643///@}
1644///@name RVec Comparison and Logical Operators
1645///@{
1646
1647#define RVEC_LOGICAL_OPERATOR(OP) \
1648template <typename T0, typename T1> \
1649auto operator OP(const RVec<T0> &v, const T1 &y) \
1650 -> RVec<int> /* avoid std::vector<bool> */ \
1651{ \
1652 RVec<int> ret(v.size()); \
1653 auto op = [y](const T0 &x) -> int { return x OP y; }; \
1654 std::transform(v.begin(), v.end(), ret.begin(), op); \
1655 return ret; \
1656} \
1657 \
1658template <typename T0, typename T1> \
1659auto operator OP(const T0 &x, const RVec<T1> &v) \
1660 -> RVec<int> /* avoid std::vector<bool> */ \
1661{ \
1662 RVec<int> ret(v.size()); \
1663 auto op = [x](const T1 &y) -> int { return x OP y; }; \
1664 std::transform(v.begin(), v.end(), ret.begin(), op); \
1665 return ret; \
1666} \
1667 \
1668template <typename T0, typename T1> \
1669auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
1670 -> RVec<int> /* avoid std::vector<bool> */ \
1671{ \
1672 if (v0.size() != v1.size()) \
1673 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1674 \
1675 RVec<int> ret(v0.size()); \
1676 auto op = [](const T0 &x, const T1 &y) -> int { return x OP y; }; \
1677 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
1678 return ret; \
1679} \
1680
1689#undef RVEC_LOGICAL_OPERATOR
1690
1691///@}
1692///@name RVec Standard Mathematical Functions
1693///@{
1694
1695/// \cond
1696template <typename T> struct PromoteTypeImpl;
1697
1698template <> struct PromoteTypeImpl<float> { using Type = float; };
1699template <> struct PromoteTypeImpl<double> { using Type = double; };
1700template <> struct PromoteTypeImpl<long double> { using Type = long double; };
1701
1702template <typename T> struct PromoteTypeImpl { using Type = double; };
1703
1704template <typename T>
1705using PromoteType = typename PromoteTypeImpl<T>::Type;
1706
1707template <typename U, typename V>
1708using PromoteTypes = decltype(PromoteType<U>() + PromoteType<V>());
1709
1710/// \endcond
1711
1712#define RVEC_UNARY_FUNCTION(NAME, FUNC) \
1713 template <typename T> \
1714 RVec<PromoteType<T>> NAME(const RVec<T> &v) \
1715 { \
1716 RVec<PromoteType<T>> ret(v.size()); \
1717 auto f = [](const T &x) { return FUNC(x); }; \
1718 std::transform(v.begin(), v.end(), ret.begin(), f); \
1719 return ret; \
1720 }
1721
1722#define RVEC_BINARY_FUNCTION(NAME, FUNC) \
1723 template <typename T0, typename T1> \
1724 RVec<PromoteTypes<T0, T1>> NAME(const T0 &x, const RVec<T1> &v) \
1725 { \
1726 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
1727 auto f = [&x](const T1 &y) { return FUNC(x, y); }; \
1728 std::transform(v.begin(), v.end(), ret.begin(), f); \
1729 return ret; \
1730 } \
1731 \
1732 template <typename T0, typename T1> \
1733 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v, const T1 &y) \
1734 { \
1735 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
1736 auto f = [&y](const T1 &x) { return FUNC(x, y); }; \
1737 std::transform(v.begin(), v.end(), ret.begin(), f); \
1738 return ret; \
1739 } \
1740 \
1741 template <typename T0, typename T1> \
1742 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v0, const RVec<T1> &v1) \
1743 { \
1744 if (v0.size() != v1.size()) \
1745 throw std::runtime_error(ERROR_MESSAGE(NAME)); \
1746 \
1747 RVec<PromoteTypes<T0, T1>> ret(v0.size()); \
1748 auto f = [](const T0 &x, const T1 &y) { return FUNC(x, y); }; \
1749 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), f); \
1750 return ret; \
1751 } \
1752
1753#define RVEC_STD_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, std::F)
1754#define RVEC_STD_BINARY_FUNCTION(F) RVEC_BINARY_FUNCTION(F, std::F)
1755
1760
1764
1769
1774
1782
1789
1796
1801#undef RVEC_STD_UNARY_FUNCTION
1802
1803///@}
1804///@name RVec Fast Mathematical Functions with Vdt
1805///@{
1806
1807#ifdef R__HAS_VDT
1808#define RVEC_VDT_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, vdt::F)
1809
1810RVEC_VDT_UNARY_FUNCTION(fast_expf)
1811RVEC_VDT_UNARY_FUNCTION(fast_logf)
1812RVEC_VDT_UNARY_FUNCTION(fast_sinf)
1813RVEC_VDT_UNARY_FUNCTION(fast_cosf)
1814RVEC_VDT_UNARY_FUNCTION(fast_tanf)
1815RVEC_VDT_UNARY_FUNCTION(fast_asinf)
1816RVEC_VDT_UNARY_FUNCTION(fast_acosf)
1817RVEC_VDT_UNARY_FUNCTION(fast_atanf)
1818
1819RVEC_VDT_UNARY_FUNCTION(fast_exp)
1820RVEC_VDT_UNARY_FUNCTION(fast_log)
1821RVEC_VDT_UNARY_FUNCTION(fast_sin)
1822RVEC_VDT_UNARY_FUNCTION(fast_cos)
1823RVEC_VDT_UNARY_FUNCTION(fast_tan)
1824RVEC_VDT_UNARY_FUNCTION(fast_asin)
1825RVEC_VDT_UNARY_FUNCTION(fast_acos)
1826RVEC_VDT_UNARY_FUNCTION(fast_atan)
1827#undef RVEC_VDT_UNARY_FUNCTION
1828
1829#endif // R__HAS_VDT
1830
1831#undef RVEC_UNARY_FUNCTION
1832
1833///@}
1834
1835/// Inner product
1836///
1837/// Example code, at the ROOT prompt:
1838/// ~~~{.cpp}
1839/// using namespace ROOT::VecOps;
1840/// RVec<float> v1 {1., 2., 3.};
1841/// RVec<float> v2 {4., 5., 6.};
1842/// auto v1_dot_v2 = Dot(v1, v2);
1843/// v1_dot_v2
1844/// // (float) 32.0000f
1845/// ~~~
1846template <typename T, typename V>
1847auto Dot(const RVec<T> &v0, const RVec<V> &v1) -> decltype(v0[0] * v1[0])
1848{
1849 if (v0.size() != v1.size())
1850 throw std::runtime_error("Cannot compute inner product of vectors of different sizes");
1851 return std::inner_product(v0.begin(), v0.end(), v1.begin(), decltype(v0[0] * v1[0])(0));
1852}
1853
1854/// Sum elements of an RVec
1855///
1856/// Example code, at the ROOT prompt:
1857/// ~~~{.cpp}
1858/// using namespace ROOT::VecOps;
1859/// RVecF v {1.f, 2.f, 3.f};
1860/// auto v_sum = Sum(v);
1861/// v_sum
1862/// // (float) 6.f
1863/// auto v_sum_d = Sum(v, 0.);
1864/// v_sum_d
1865/// // (double) 6.0000000
1866/// ~~~
1867/// ~~~{.cpp}
1868/// using namespace ROOT::VecOps;
1869/// const ROOT::Math::PtEtaPhiMVector lv0 {15.5f, .3f, .1f, 105.65f},
1870/// lv1 {34.32f, 2.2f, 3.02f, 105.65f},
1871/// lv2 {12.95f, 1.32f, 2.2f, 105.65f};
1872/// RVec<ROOT::Math::PtEtaPhiMVector> v {lv0, lv1, lv2};
1873/// auto v_sum_lv = Sum(v, ROOT::Math::PtEtaPhiMVector());
1874/// v_sum_lv
1875/// // (ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double> > &) (30.8489,2.46534,2.58947,361.084)
1876/// ~~~
1877template <typename T>
1878T Sum(const RVec<T> &v, const T zero = T(0))
1879{
1880 return std::accumulate(v.begin(), v.end(), zero);
1881}
1882
1883/// Product method analogous to Sum:
1884/// Multiply the elements of a vector
1885template <typename T>
1886T Product(const RVec<T> &v, const T init = T(1)) // initialize with identity
1887{
1888 return std::accumulate(v.begin(), v.end(), init, std::multiplies<T>());
1889}
1890
1891/// Get the mean of the elements of an RVec
1892///
1893/// The return type is a double precision floating point number.
1894///
1895/// Example code, at the ROOT prompt:
1896/// ~~~{.cpp}
1897/// using namespace ROOT::VecOps;
1898/// RVecF v {1.f, 2.f, 4.f};
1899/// auto v_mean = Mean(v);
1900/// v_mean
1901/// // (double) 2.3333333
1902/// ~~~
1903template <typename T>
1904double Mean(const RVec<T> &v)
1905{
1906 if (v.empty()) return 0.;
1907 return double(Sum(v)) / v.size();
1908}
1909
1910/// Get the mean of the elements of an RVec with custom initial value
1911///
1912/// The return type will be deduced from the `zero` parameter
1913///
1914/// Example code, at the ROOT prompt:
1915/// ~~~{.cpp}
1916/// using namespace ROOT::VecOps;
1917/// RVecF v {1.f, 2.f, 4.f};
1918/// auto v_mean_f = Mean(v, 0.f);
1919/// v_mean_f
1920/// // (float) 2.33333f
1921/// auto v_mean_d = Mean(v, 0.);
1922/// v_mean_d
1923/// // (double) 2.3333333
1924/// ~~~
1925/// ~~~{.cpp}
1926/// using namespace ROOT::VecOps;
1927/// const ROOT::Math::PtEtaPhiMVector lv0 {15.5f, .3f, .1f, 105.65f},
1928/// lv1 {34.32f, 2.2f, 3.02f, 105.65f},
1929/// lv2 {12.95f, 1.32f, 2.2f, 105.65f};
1930/// RVec<ROOT::Math::PtEtaPhiMVector> v {lv0, lv1, lv2};
1931/// auto v_mean_lv = Mean(v, ROOT::Math::PtEtaPhiMVector());
1932/// v_mean_lv
1933/// // (ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double> > &) (10.283,2.46534,2.58947,120.361)
1934/// ~~~
1935template <typename T, typename R = T>
1936R Mean(const RVec<T> &v, const R zero)
1937{
1938 if (v.empty()) return zero;
1939 return Sum(v, zero) / v.size();
1940}
1941
1942/// Get the greatest element of an RVec
1943///
1944/// Example code, at the ROOT prompt:
1945/// ~~~~{.cpp}
1946/// using namespace ROOT::VecOps;
1947/// RVecF v {1.f, 2.f, 4.f};
1948/// auto v_max = Max(v);
1949/// v_max
1950/// (float) 4.00000f
1951/// ~~~~
1952template <typename T>
1954{
1955 return *std::max_element(v.begin(), v.end());
1956}
1957
1958/// Get the smallest element of an RVec
1959///
1960/// Example code, at the ROOT prompt:
1961/// ~~~~{.cpp}
1962/// using namespace ROOT::VecOps;
1963/// RVecF v {1.f, 2.f, 4.f};
1964/// auto v_min = Min(v);
1965/// v_min
1966/// (float) 1.00000f
1967/// ~~~~
1968template <typename T>
1970{
1971 return *std::min_element(v.begin(), v.end());
1972}
1973
1974/// Get the index of the greatest element of an RVec
1975/// In case of multiple occurrences of the maximum values,
1976/// the index corresponding to the first occurrence is returned.
1977///
1978/// Example code, at the ROOT prompt:
1979/// ~~~~{.cpp}
1980/// using namespace ROOT::VecOps;
1981/// RVecF v {1.f, 2.f, 4.f};
1982/// auto v_argmax = ArgMax(v);
1983/// v_argmax
1984/// // (unsigned long) 2
1985/// ~~~~
1986template <typename T>
1987std::size_t ArgMax(const RVec<T> &v)
1988{
1989 return std::distance(v.begin(), std::max_element(v.begin(), v.end()));
1990}
1991
1992/// Get the index of the smallest element of an RVec
1993/// In case of multiple occurrences of the minimum values,
1994/// the index corresponding to the first occurrence is returned.
1995///
1996/// Example code, at the ROOT prompt:
1997/// ~~~~{.cpp}
1998/// using namespace ROOT::VecOps;
1999/// RVecF v {1.f, 2.f, 4.f};
2000/// auto v_argmin = ArgMin(v);
2001/// v_argmin
2002/// // (unsigned long) 0
2003/// ~~~~
2004template <typename T>
2005std::size_t ArgMin(const RVec<T> &v)
2006{
2007 return std::distance(v.begin(), std::min_element(v.begin(), v.end()));
2008}
2009
2010/// Get the variance of the elements of an RVec
2011///
2012/// The return type is a double precision floating point number.
2013/// Example code, at the ROOT prompt:
2014/// ~~~{.cpp}
2015/// using namespace ROOT::VecOps;
2016/// RVecF v {1.f, 2.f, 4.f};
2017/// auto v_var = Var(v);
2018/// v_var
2019/// // (double) 2.3333333
2020/// ~~~
2021template <typename T>
2022double Var(const RVec<T> &v)
2023{
2024 const std::size_t size = v.size();
2025 if (size < std::size_t(2)) return 0.;
2026 T sum_squares(0), squared_sum(0);
2027 auto pred = [&sum_squares, &squared_sum](const T& x) {sum_squares+=x*x; squared_sum+=x;};
2028 std::for_each(v.begin(), v.end(), pred);
2029 squared_sum *= squared_sum;
2030 const auto dsize = (double) size;
2031 return 1. / (dsize - 1.) * (sum_squares - squared_sum / dsize );
2032}
2033
2034/// Get the standard deviation of the elements of an RVec
2035///
2036/// The return type is a double precision floating point number.
2037/// Example code, at the ROOT prompt:
2038/// ~~~{.cpp}
2039/// using namespace ROOT::VecOps;
2040/// RVecF v {1.f, 2.f, 4.f};
2041/// auto v_sd = StdDev(v);
2042/// v_sd
2043/// // (double) 1.5275252
2044/// ~~~
2045template <typename T>
2046double StdDev(const RVec<T> &v)
2047{
2048 return std::sqrt(Var(v));
2049}
2050
2051/// Create new collection applying a callable to the elements of the input collection
2052///
2053/// Example code, at the ROOT prompt:
2054/// ~~~{.cpp}
2055/// using namespace ROOT::VecOps;
2056/// RVecF v {1.f, 2.f, 4.f};
2057/// auto v_square = Map(v, [](float f){return f* 2.f;});
2058/// v_square
2059/// // (ROOT::VecOps::RVec<float> &) { 2.00000f, 4.00000f, 8.00000f }
2060///
2061/// RVecF x({1.f, 2.f, 3.f});
2062/// RVecF y({4.f, 5.f, 6.f});
2063/// RVecF z({7.f, 8.f, 9.f});
2064/// auto mod = [](float x, float y, float z) { return sqrt(x * x + y * y + z * z); };
2065/// auto v_mod = Map(x, y, z, mod);
2066/// v_mod
2067/// // (ROOT::VecOps::RVec<float> &) { 8.12404f, 9.64365f, 11.2250f }
2068/// ~~~
2069template <typename... Args>
2070auto Map(Args &&... args)
2071{
2072 /*
2073 Here the strategy in order to generalise the previous implementation of Map, i.e.
2074 `RVec Map(RVec, F)`, here we need to move the last parameter of the pack in first
2075 position in order to be able to invoke the Map function with automatic type deduction.
2076 This is achieved in two steps:
2077 1. Forward as tuple the pack to MapFromTuple
2078 2. Invoke the MapImpl helper which has the signature `template<...T, F> RVec MapImpl(F &&f, RVec<T>...)`
2079 */
2080
2081 // check the first N - 1 arguments are RVecs
2082 constexpr auto nArgs = sizeof...(Args);
2084 static_assert(ROOT::Internal::VecOps::All(isRVec, nArgs - 1),
2085 "Map: the first N-1 arguments must be RVecs or references to RVecs");
2086
2087 return ROOT::Internal::VecOps::MapFromTuple(std::forward_as_tuple(args...),
2088 std::make_index_sequence<sizeof...(args) - 1>());
2089}
2090
2091/// Create a new collection with the elements passing the filter expressed by the predicate
2092///
2093/// Example code, at the ROOT prompt:
2094/// ~~~{.cpp}
2095/// using namespace ROOT::VecOps;
2096/// RVecI v {1, 2, 4};
2097/// auto v_even = Filter(v, [](int i){return 0 == i%2;});
2098/// v_even
2099/// // (ROOT::VecOps::RVec<int> &) { 2, 4 }
2100/// ~~~
2101template <typename T, typename F>
2103{
2104 const auto thisSize = v.size();
2105 RVec<T> w;
2106 w.reserve(thisSize);
2107 for (auto &&val : v) {
2108 if (f(val))
2109 w.emplace_back(val);
2110 }
2111 return w;
2112}
2113
2114/// Return true if any of the elements equates to true, return false otherwise.
2115///
2116/// Example code, at the ROOT prompt:
2117/// ~~~{.cpp}
2118/// using namespace ROOT::VecOps;
2119/// RVecI v {0, 1, 0};
2120/// auto anyTrue = Any(v);
2121/// anyTrue
2122/// // (bool) true
2123/// ~~~
2124template <typename T>
2125auto Any(const RVec<T> &v) -> decltype(v[0] == true)
2126{
2127 for (auto &&e : v)
2128 if (static_cast<bool>(e) == true)
2129 return true;
2130 return false;
2131}
2132
2133/// Return true if all of the elements equate to true, return false otherwise.
2134///
2135/// Example code, at the ROOT prompt:
2136/// ~~~{.cpp}
2137/// using namespace ROOT::VecOps;
2138/// RVecI v {0, 1, 0};
2139/// auto allTrue = All(v);
2140/// allTrue
2141/// // (bool) false
2142/// ~~~
2143template <typename T>
2144auto All(const RVec<T> &v) -> decltype(v[0] == false)
2145{
2146 for (auto &&e : v)
2147 if (static_cast<bool>(e) == false)
2148 return false;
2149 return true;
2150}
2151
2152template <typename T>
2153void swap(RVec<T> &lhs, RVec<T> &rhs)
2154{
2155 lhs.swap(rhs);
2156}
2157
2158/// Return an RVec of indices that sort the input RVec
2159///
2160/// Example code, at the ROOT prompt:
2161/// ~~~{.cpp}
2162/// using namespace ROOT::VecOps;
2163/// RVecD v {2., 3., 1.};
2164/// auto sortIndices = Argsort(v)
2165/// // (ROOT::VecOps::RVec<unsigned long> &) { 2, 0, 1 }
2166/// auto values = Take(v, sortIndices)
2167/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 3.0000000 }
2168/// ~~~
2169template <typename T>
2171{
2172 using size_type = typename RVec<T>::size_type;
2173 RVec<size_type> i(v.size());
2174 std::iota(i.begin(), i.end(), 0);
2175 std::sort(i.begin(), i.end(), [&v](size_type i1, size_type i2) { return v[i1] < v[i2]; });
2176 return i;
2177}
2178
2179/// Return an RVec of indices that sort the input RVec based on a comparison function.
2180///
2181/// Example code, at the ROOT prompt:
2182/// ~~~{.cpp}
2183/// using namespace ROOT::VecOps;
2184/// RVecD v {2., 3., 1.};
2185/// auto sortIndices = Argsort(v, [](double x, double y) {return x > y;})
2186/// // (ROOT::VecOps::RVec<unsigned long> &) { 1, 0, 2 }
2187/// auto values = Take(v, sortIndices)
2188/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 1.0000000 }
2189/// ~~~
2190template <typename T, typename Compare>
2192{
2193 using size_type = typename RVec<T>::size_type;
2194 RVec<size_type> i(v.size());
2195 std::iota(i.begin(), i.end(), 0);
2196 std::sort(i.begin(), i.end(),
2197 [&v, &c](size_type i1, size_type i2) { return c(v[i1], v[i2]); });
2198 return i;
2199}
2200
2201/// Return an RVec of indices that sort the input RVec
2202/// while keeping the order of equal elements.
2203/// This is the stable variant of `Argsort`.
2204///
2205/// Example code, at the ROOT prompt:
2206/// ~~~{.cpp}
2207/// using namespace ROOT::VecOps;
2208/// RVecD v {2., 3., 2., 1.};
2209/// auto sortIndices = StableArgsort(v)
2210/// // (ROOT::VecOps::RVec<unsigned long> &) { 3, 0, 2, 1 }
2211/// auto values = Take(v, sortIndices)
2212/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 2.0000000, 3.0000000 }
2213/// ~~~
2214template <typename T>
2216{
2217 using size_type = typename RVec<T>::size_type;
2218 RVec<size_type> i(v.size());
2219 std::iota(i.begin(), i.end(), 0);
2220 std::stable_sort(i.begin(), i.end(), [&v](size_type i1, size_type i2) { return v[i1] < v[i2]; });
2221 return i;
2222}
2223
2224/// Return an RVec of indices that sort the input RVec based on a comparison function
2225/// while keeping the order of equal elements.
2226/// This is the stable variant of `Argsort`.
2227///
2228/// Example code, at the ROOT prompt:
2229/// ~~~{.cpp}
2230/// using namespace ROOT::VecOps;
2231/// RVecD v {2., 3., 2., 1.};
2232/// auto sortIndices = StableArgsort(v, [](double x, double y) {return x > y;})
2233/// // (ROOT::VecOps::RVec<unsigned long> &) { 1, 0, 2, 3 }
2234/// auto values = Take(v, sortIndices)
2235/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 2.0000000, 1.0000000 }
2236/// ~~~
2237template <typename T, typename Compare>
2239{
2240 using size_type = typename RVec<T>::size_type;
2241 RVec<size_type> i(v.size());
2242 std::iota(i.begin(), i.end(), 0);
2243 std::stable_sort(i.begin(), i.end(), [&v, &c](size_type i1, size_type i2) { return c(v[i1], v[i2]); });
2244 return i;
2245}
2246
2247/// Return elements of a vector at given indices
2248///
2249/// Example code, at the ROOT prompt:
2250/// ~~~{.cpp}
2251/// using namespace ROOT::VecOps;
2252/// RVecD v {2., 3., 1.};
2253/// auto vTaken = Take(v, {0,2});
2254/// vTaken
2255/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 1.0000000 }
2256/// ~~~
2257
2258template <typename T>
2259RVec<T> Take(const RVec<T> &v, const RVec<typename RVec<T>::size_type> &i)
2260{
2261 using size_type = typename RVec<T>::size_type;
2262 const size_type isize = i.size();
2263 RVec<T> r(isize);
2264 for (size_type k = 0; k < isize; k++)
2265 r[k] = v[i[k]];
2266 return r;
2267}
2268
2269/// Take version that defaults to (user-specified) output value if some index is out of range
2270template <typename T>
2271RVec<T> Take(const RVec<T> &v, const RVec<typename RVec<T>::size_type> &i, const T default_val)
2272{
2273 using size_type = typename RVec<T>::size_type;
2274 const size_type isize = i.size();
2275 RVec<T> r(isize);
2276 for (size_type k = 0; k < isize; k++)
2277 {
2278 if (k < v.size()){
2279 r[k] = v[i[k]];
2280 }
2281 else {
2282 r[k] = default_val;
2283 }
2284 }
2285 return r;
2286}
2287
2288
2289/// Return first or last `n` elements of an RVec
2290///
2291/// if `n > 0` and last elements if `n < 0`.
2292///
2293/// Example code, at the ROOT prompt:
2294/// ~~~{.cpp}
2295/// using namespace ROOT::VecOps;
2296/// RVecD v {2., 3., 1.};
2297/// auto firstTwo = Take(v, 2);
2298/// firstTwo
2299/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 3.0000000 }
2300/// auto lastOne = Take(v, -1);
2301/// lastOne
2302/// // (ROOT::VecOps::RVec<double>) { 1.0000000 }
2303/// ~~~
2304template <typename T>
2305RVec<T> Take(const RVec<T> &v, const int n)
2306{
2307 using size_type = typename RVec<T>::size_type;
2308 const size_type size = v.size();
2309 const size_type absn = std::abs(n);
2310 if (absn > size) {
2311 std::stringstream ss;
2312 ss << "Try to take " << absn << " elements but vector has only size " << size << ".";
2313 throw std::runtime_error(ss.str());
2314 }
2315 RVec<T> r(absn);
2316 if (n < 0) {
2317 for (size_type k = 0; k < absn; k++)
2318 r[k] = v[size - absn + k];
2319 } else {
2320 for (size_type k = 0; k < absn; k++)
2321 r[k] = v[k];
2322 }
2323 return r;
2324}
2325
2326/// Return a copy of the container without the elements at the specified indices.
2327///
2328/// Duplicated and out-of-range indices in idxs are ignored.
2329template <typename T>
2331{
2332 // clean up input indices
2333 std::sort(idxs.begin(), idxs.end());
2334 idxs.erase(std::unique(idxs.begin(), idxs.end()), idxs.end());
2335
2336 RVec<T> r;
2337 if (v.size() > idxs.size())
2338 r.reserve(v.size() - idxs.size());
2339
2340 auto discardIt = idxs.begin();
2341 using sz_t = typename RVec<T>::size_type;
2342 for (sz_t i = 0u; i < v.size(); ++i) {
2343 if (discardIt != idxs.end() && i == *discardIt)
2344 ++discardIt;
2345 else
2346 r.emplace_back(v[i]);
2347 }
2348
2349 return r;
2350}
2351
2352/// Return copy of reversed vector
2353///
2354/// Example code, at the ROOT prompt:
2355/// ~~~{.cpp}
2356/// using namespace ROOT::VecOps;
2357/// RVecD v {2., 3., 1.};
2358/// auto v_reverse = Reverse(v);
2359/// v_reverse
2360/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 3.0000000, 2.0000000 }
2361/// ~~~
2362template <typename T>
2364{
2365 RVec<T> r(v);
2366 std::reverse(r.begin(), r.end());
2367 return r;
2368}
2369
2370/// Return copy of RVec with elements sorted in ascending order
2371///
2372/// This helper is different from Argsort since it does not return an RVec of indices,
2373/// but an RVec of values.
2374///
2375/// Example code, at the ROOT prompt:
2376/// ~~~{.cpp}
2377/// using namespace ROOT::VecOps;
2378/// RVecD v {2., 3., 1.};
2379/// auto v_sorted = Sort(v);
2380/// v_sorted
2381/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 3.0000000 }
2382/// ~~~
2383template <typename T>
2385{
2386 RVec<T> r(v);
2387 std::sort(r.begin(), r.end());
2388 return r;
2389}
2390
2391/// Return copy of RVec with elements sorted based on a comparison operator
2392///
2393/// The comparison operator has to fullfill the same requirements of the
2394/// predicate of by std::sort.
2395///
2396///
2397/// This helper is different from Argsort since it does not return an RVec of indices,
2398/// but an RVec of values.
2399///
2400/// Example code, at the ROOT prompt:
2401/// ~~~{.cpp}
2402/// using namespace ROOT::VecOps;
2403/// RVecD v {2., 3., 1.};
2404/// auto v_sorted = Sort(v, [](double x, double y) {return 1/x < 1/y;});
2405/// v_sorted
2406/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 1.0000000 }
2407/// ~~~
2408template <typename T, typename Compare>
2410{
2411 RVec<T> r(v);
2412 std::sort(r.begin(), r.end(), std::forward<Compare>(c));
2413 return r;
2414}
2415
2416/// Return copy of RVec with elements sorted in ascending order
2417/// while keeping the order of equal elements.
2418///
2419/// This is the stable variant of `Sort`.
2420///
2421/// This helper is different from StableArgsort since it does not return an RVec of indices,
2422/// but an RVec of values.
2423///
2424/// Example code, at the ROOT prompt:
2425/// ~~~{.cpp}
2426/// using namespace ROOT::VecOps;
2427/// RVecD v {2., 3., 2, 1.};
2428/// auto v_sorted = StableSort(v);
2429/// v_sorted
2430/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 2.0000000, 3.0000000 }
2431/// ~~~
2432template <typename T>
2434{
2435 RVec<T> r(v);
2436 std::stable_sort(r.begin(), r.end());
2437 return r;
2438}
2439
2440// clang-format off
2441/// Return copy of RVec with elements sorted based on a comparison operator
2442/// while keeping the order of equal elements.
2443///
2444/// The comparison operator has to fullfill the same requirements of the
2445/// predicate of std::stable_sort.
2446///
2447/// This helper is different from StableArgsort since it does not return an RVec of indices,
2448/// but an RVec of values.
2449///
2450/// This is the stable variant of `Sort`.
2451///
2452/// Example code, at the ROOT prompt:
2453/// ~~~{.cpp}
2454/// using namespace ROOT::VecOps;
2455/// RVecD v {2., 3., 2., 1.};
2456/// auto v_sorted = StableSort(v, [](double x, double y) {return 1/x < 1/y;});
2457/// v_sorted
2458/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 2.0000000, 1.0000000 }
2459/// ~~~
2460/// ~~~{.cpp}
2461/// using namespace ROOT::VecOps;
2462/// RVec<RVecD> v {{2., 4.}, {3., 1.}, {2, 1.}, {1., 4.}};
2463/// auto v_sorted = StableSort(StableSort(v, [](const RVecD &x, const RVecD &y) {return x[1] < y[1];}), [](const RVecD &x, const RVecD &y) {return x[0] < y[0];});
2464/// v_sorted
2465/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<double> > &) { { 1.0000000, 4.0000000 }, { 2.0000000, 1.0000000 }, { 2.0000000, 4.0000000 }, { 3.0000000, 1.0000000 } }
2466/// ~~~
2467// clang-format off
2468template <typename T, typename Compare>
2470{
2471 RVec<T> r(v);
2472 std::stable_sort(r.begin(), r.end(), std::forward<Compare>(c));
2473 return r;
2474}
2475
2476/// Return the indices that represent all combinations of the elements of two
2477/// RVecs.
2478///
2479/// The type of the return value is an RVec of two RVecs containing indices.
2480///
2481/// Example code, at the ROOT prompt:
2482/// ~~~{.cpp}
2483/// using namespace ROOT::VecOps;
2484/// auto comb_idx = Combinations(3, 2);
2485/// comb_idx
2486/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 1, 1, 2, 2 }, { 0, 1, 0, 1, 0, 1 } }
2487/// ~~~
2488inline RVec<RVec<std::size_t>> Combinations(const std::size_t size1, const std::size_t size2)
2489{
2490 using size_type = std::size_t;
2492 r[0].resize(size1*size2);
2493 r[1].resize(size1*size2);
2494 size_type c = 0;
2495 for(size_type i=0; i<size1; i++) {
2496 for(size_type j=0; j<size2; j++) {
2497 r[0][c] = i;
2498 r[1][c] = j;
2499 c++;
2500 }
2501 }
2502 return r;
2503}
2504
2505/// Return the indices that represent all combinations of the elements of two
2506/// RVecs.
2507///
2508/// The type of the return value is an RVec of two RVecs containing indices.
2509///
2510/// Example code, at the ROOT prompt:
2511/// ~~~{.cpp}
2512/// using namespace ROOT::VecOps;
2513/// RVecD v1 {1., 2., 3.};
2514/// RVecD v2 {-4., -5.};
2515/// auto comb_idx = Combinations(v1, v2);
2516/// comb_idx
2517/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 1, 1, 2, 2 }, { 0, 1, 0, 1, 0, 1 } }
2518/// ~~~
2519template <typename T1, typename T2>
2521{
2522 return Combinations(v1.size(), v2.size());
2523}
2524
2525/// Return the indices that represent all unique combinations of the
2526/// elements of a given RVec.
2527///
2528/// ~~~{.cpp}
2529/// using namespace ROOT::VecOps;
2530/// RVecD v {1., 2., 3., 4.};
2531/// auto v_1 = Combinations(v, 1);
2532/// v_1
2533/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 1, 2, 3 } }
2534/// auto v_2 = Combinations(v, 2);
2535/// v_2
2536/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 0, 1, 1, 2 }, { 1, 2, 3, 2, 3, 3 } }
2537/// auto v_3 = Combinations(v, 3);
2538/// v_3
2539/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 0, 1 }, { 1, 1, 2, 2 }, { 2, 3, 3, 3 } }
2540/// auto v_4 = Combinations(v, 4);
2541/// v_4
2542/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0 }, { 1 }, { 2 }, { 3 } }
2543/// ~~~
2544template <typename T>
2546{
2547 using size_type = typename RVec<T>::size_type;
2548 const size_type s = v.size();
2549 if (n > s) {
2550 throw std::runtime_error("Cannot make unique combinations of size " + std::to_string(n) +
2551 " from vector of size " + std::to_string(s) + ".");
2552 }
2553
2555 for(size_type k=0; k<s; k++)
2556 indices[k] = k;
2557
2558 const auto innersize = [=] {
2559 size_type inners = s - n + 1;
2560 for (size_type m = s - n + 2; m <= s; ++m)
2561 inners *= m;
2562
2563 size_type factn = 1;
2564 for (size_type i = 2; i <= n; ++i)
2565 factn *= i;
2566 inners /= factn;
2567
2568 return inners;
2569 }();
2570
2572 size_type inneridx = 0;
2573 for (size_type k = 0; k < n; k++)
2574 c[k][inneridx] = indices[k];
2575 ++inneridx;
2576
2577 while (true) {
2578 bool run_through = true;
2579 long i = n - 1;
2580 for (; i>=0; i--) {
2581 if (indices[i] != i + s - n){
2582 run_through = false;
2583 break;
2584 }
2585 }
2586 if (run_through) {
2587 return c;
2588 }
2589 indices[i]++;
2590 for (long j=i+1; j<(long)n; j++)
2591 indices[j] = indices[j-1] + 1;
2592 for (size_type k = 0; k < n; k++)
2593 c[k][inneridx] = indices[k];
2594 ++inneridx;
2595 }
2596}
2597
2598/// Return the indices of the elements which are not zero
2599///
2600/// Example code, at the ROOT prompt:
2601/// ~~~{.cpp}
2602/// using namespace ROOT::VecOps;
2603/// RVecD v {2., 0., 3., 0., 1.};
2604/// auto nonzero_idx = Nonzero(v);
2605/// nonzero_idx
2606/// // (ROOT::VecOps::RVec<unsigned long> &) { 0, 2, 4 }
2607/// ~~~
2608template <typename T>
2610{
2611 using size_type = typename RVec<T>::size_type;
2613 const auto size = v.size();
2614 r.reserve(size);
2615 for(size_type i=0; i<size; i++) {
2616 if(v[i] != 0) {
2617 r.emplace_back(i);
2618 }
2619 }
2620 return r;
2621}
2622
2623/// Return the intersection of elements of two RVecs.
2624///
2625/// Each element of v1 is looked up in v2 and added to the returned vector if
2626/// found. Following, the order of v1 is preserved. If v2 is already sorted, the
2627/// optional argument v2_is_sorted can be used to toggle of the internal sorting
2628/// step, therewith optimising runtime.
2629///
2630/// Example code, at the ROOT prompt:
2631/// ~~~{.cpp}
2632/// using namespace ROOT::VecOps;
2633/// RVecD v1 {1., 2., 3.};
2634/// RVecD v2 {-4., -5., 2., 1.};
2635/// auto v1_intersect_v2 = Intersect(v1, v2);
2636/// v1_intersect_v2
2637/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000 }
2638/// ~~~
2639template <typename T>
2640RVec<T> Intersect(const RVec<T>& v1, const RVec<T>& v2, bool v2_is_sorted = false)
2641{
2642 RVec<T> v2_sorted;
2643 if (!v2_is_sorted) v2_sorted = Sort(v2);
2644 const auto v2_begin = v2_is_sorted ? v2.begin() : v2_sorted.begin();
2645 const auto v2_end = v2_is_sorted ? v2.end() : v2_sorted.end();
2646 RVec<T> r;
2647 const auto size = v1.size();
2648 r.reserve(size);
2649 using size_type = typename RVec<T>::size_type;
2650 for(size_type i=0; i<size; i++) {
2651 if (std::binary_search(v2_begin, v2_end, v1[i])) {
2652 r.emplace_back(v1[i]);
2653 }
2654 }
2655 return r;
2656}
2657
2658/// Return the elements of v1 if the condition c is true and v2 if the
2659/// condition c is false.
2660///
2661/// Example code, at the ROOT prompt:
2662/// ~~~{.cpp}
2663/// using namespace ROOT::VecOps;
2664/// RVecD v1 {1., 2., 3.};
2665/// RVecD v2 {-1., -2., -3.};
2666/// auto c = v1 > 1;
2667/// c
2668/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2669/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2670/// if_c_v1_else_v2
2671/// // (ROOT::VecOps::RVec<double> &) { -1.0000000, 2.0000000, 3.0000000 }
2672/// ~~~
2673template <typename T>
2674RVec<T> Where(const RVec<int>& c, const RVec<T>& v1, const RVec<T>& v2)
2675{
2676 using size_type = typename RVec<T>::size_type;
2677 const size_type size = c.size();
2678 RVec<T> r;
2679 r.reserve(size);
2680 for (size_type i=0; i<size; i++) {
2681 r.emplace_back(c[i] != 0 ? v1[i] : v2[i]);
2682 }
2683 return r;
2684}
2685
2686/// Return the elements of v1 if the condition c is true and sets the value v2
2687/// if the condition c is false.
2688///
2689/// Example code, at the ROOT prompt:
2690/// ~~~{.cpp}
2691/// using namespace ROOT::VecOps;
2692/// RVecD v1 {1., 2., 3.};
2693/// double v2 = 4.;
2694/// auto c = v1 > 1;
2695/// c
2696/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2697/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2698/// if_c_v1_else_v2
2699/// // (ROOT::VecOps::RVec<double>) { 4.0000000, 2.0000000, 3.0000000 }
2700/// ~~~
2701template <typename T>
2703{
2704 using size_type = typename RVec<T>::size_type;
2705 const size_type size = c.size();
2706 RVec<T> r;
2707 r.reserve(size);
2708 for (size_type i=0; i<size; i++) {
2709 r.emplace_back(c[i] != 0 ? v1[i] : v2);
2710 }
2711 return r;
2712}
2713
2714/// Return the elements of v2 if the condition c is false and sets the value v1
2715/// if the condition c is true.
2716///
2717/// Example code, at the ROOT prompt:
2718/// ~~~{.cpp}
2719/// using namespace ROOT::VecOps;
2720/// double v1 = 4.;
2721/// RVecD v2 {1., 2., 3.};
2722/// auto c = v2 > 1;
2723/// c
2724/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2725/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2726/// if_c_v1_else_v2
2727/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 4.0000000, 4.0000000 }
2728/// ~~~
2729template <typename T>
2731{
2732 using size_type = typename RVec<T>::size_type;
2733 const size_type size = c.size();
2734 RVec<T> r;
2735 r.reserve(size);
2736 for (size_type i=0; i<size; i++) {
2737 r.emplace_back(c[i] != 0 ? v1 : v2[i]);
2738 }
2739 return r;
2740}
2741
2742/// Return a vector with the value v2 if the condition c is false and sets the
2743/// value v1 if the condition c is true.
2744///
2745/// Example code, at the ROOT prompt:
2746/// ~~~{.cpp}
2747/// using namespace ROOT::VecOps;
2748/// double v1 = 4.;
2749/// double v2 = 2.;
2750/// RVecI c {0, 1, 1};
2751/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2752/// if_c_v1_else_v2
2753/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 4.0000000, 4.0000000 }
2754/// ~~~
2755template <typename T>
2757{
2758 using size_type = typename RVec<T>::size_type;
2759 const size_type size = c.size();
2760 RVec<T> r;
2761 r.reserve(size);
2762 for (size_type i=0; i<size; i++) {
2763 r.emplace_back(c[i] != 0 ? v1 : v2);
2764 }
2765 return r;
2766}
2767
2768/// Return the concatenation of two RVecs.
2769///
2770/// Example code, at the ROOT prompt:
2771/// ~~~{.cpp}
2772/// using namespace ROOT::VecOps;
2773/// RVecF rvf {0.f, 1.f, 2.f};
2774/// RVecI rvi {7, 8, 9};
2775/// Concatenate(rvf, rvi)
2776/// // (ROOT::VecOps::RVec<float>) { 0.00000f, 1.00000f, 2.00000f, 7.00000f, 8.00000f, 9.00000f }
2777/// ~~~
2780{
2781 RVec<Common_t> res;
2782 res.reserve(v0.size() + v1.size());
2783 std::copy(v0.begin(), v0.end(), std::back_inserter(res));
2784 std::copy(v1.begin(), v1.end(), std::back_inserter(res));
2785 return res;
2786}
2787
2788/// Return the angle difference \f$\Delta \phi\f$ of two scalars.
2789///
2790/// The function computes the closest angle from v1 to v2 with sign and is
2791/// therefore in the range \f$[-\pi, \pi]\f$.
2792/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2793/// to degrees \f$c = 180\f$.
2794template <typename T>
2795T DeltaPhi(T v1, T v2, const T c = M_PI)
2796{
2798 "DeltaPhi must be called with floating point values.");
2799 auto r = std::fmod(v2 - v1, 2.0 * c);
2800 if (r < -c) {
2801 r += 2.0 * c;
2802 }
2803 else if (r > c) {
2804 r -= 2.0 * c;
2805 }
2806 return r;
2807}
2808
2809/// Return the angle difference \f$\Delta \phi\f$ in radians of two vectors.
2810///
2811/// The function computes the closest angle from v1 to v2 with sign and is
2812/// therefore in the range \f$[-\pi, \pi]\f$.
2813/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2814/// to degrees \f$c = 180\f$.
2815template <typename T>
2816RVec<T> DeltaPhi(const RVec<T>& v1, const RVec<T>& v2, const T c = M_PI)
2817{
2818 using size_type = typename RVec<T>::size_type;
2819 const size_type size = v1.size();
2820 auto r = RVec<T>(size);
2821 for (size_type i = 0; i < size; i++) {
2822 r[i] = DeltaPhi(v1[i], v2[i], c);
2823 }
2824 return r;
2825}
2826
2827/// Return the angle difference \f$\Delta \phi\f$ in radians of a vector and a scalar.
2828///
2829/// The function computes the closest angle from v1 to v2 with sign and is
2830/// therefore in the range \f$[-\pi, \pi]\f$.
2831/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2832/// to degrees \f$c = 180\f$.
2833template <typename T>
2834RVec<T> DeltaPhi(const RVec<T>& v1, T v2, const T c = M_PI)
2835{
2836 using size_type = typename RVec<T>::size_type;
2837 const size_type size = v1.size();
2838 auto r = RVec<T>(size);
2839 for (size_type i = 0; i < size; i++) {
2840 r[i] = DeltaPhi(v1[i], v2, c);
2841 }
2842 return r;
2843}
2844
2845/// Return the angle difference \f$\Delta \phi\f$ in radians of a scalar and a vector.
2846///
2847/// The function computes the closest angle from v1 to v2 with sign and is
2848/// therefore in the range \f$[-\pi, \pi]\f$.
2849/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2850/// to degrees \f$c = 180\f$.
2851template <typename T>
2852RVec<T> DeltaPhi(T v1, const RVec<T>& v2, const T c = M_PI)
2853{
2854 using size_type = typename RVec<T>::size_type;
2855 const size_type size = v2.size();
2856 auto r = RVec<T>(size);
2857 for (size_type i = 0; i < size; i++) {
2858 r[i] = DeltaPhi(v1, v2[i], c);
2859 }
2860 return r;
2861}
2862
2863/// Return the square of the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2864/// the collections eta1, eta2, phi1 and phi2.
2865///
2866/// The function computes \f$\Delta R^2 = (\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2\f$
2867/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2868/// be set to radian or degrees using the optional argument c, see the documentation
2869/// of the DeltaPhi helper.
2870template <typename T>
2871RVec<T> DeltaR2(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
2872{
2873 const auto dphi = DeltaPhi(phi1, phi2, c);
2874 return (eta1 - eta2) * (eta1 - eta2) + dphi * dphi;
2875}
2876
2877/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2878/// the collections eta1, eta2, phi1 and phi2.
2879///
2880/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
2881/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2882/// be set to radian or degrees using the optional argument c, see the documentation
2883/// of the DeltaPhi helper.
2884template <typename T>
2885RVec<T> DeltaR(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
2886{
2887 return sqrt(DeltaR2(eta1, eta2, phi1, phi2, c));
2888}
2889
2890/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2891/// the scalars eta1, eta2, phi1 and phi2.
2892///
2893/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
2894/// of the given scalars eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2895/// be set to radian or degrees using the optional argument c, see the documentation
2896/// of the DeltaPhi helper.
2897template <typename T>
2898T DeltaR(T eta1, T eta2, T phi1, T phi2, const T c = M_PI)
2899{
2900 const auto dphi = DeltaPhi(phi1, phi2, c);
2901 return std::sqrt((eta1 - eta2) * (eta1 - eta2) + dphi * dphi);
2902}
2903
2904/// Return the invariant mass of two particles given the collections of the quantities
2905/// transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
2906///
2907/// The function computes the invariant mass of two particles with the four-vectors
2908/// (pt1, eta2, phi1, mass1) and (pt2, eta2, phi2, mass2).
2909template <typename T>
2911 const RVec<T>& pt1, const RVec<T>& eta1, const RVec<T>& phi1, const RVec<T>& mass1,
2912 const RVec<T>& pt2, const RVec<T>& eta2, const RVec<T>& phi2, const RVec<T>& mass2)
2913{
2914 std::size_t size = pt1.size();
2915
2916 R__ASSERT(eta1.size() == size && phi1.size() == size && mass1.size() == size);
2917 R__ASSERT(pt2.size() == size && phi2.size() == size && mass2.size() == size);
2918
2919 RVec<T> inv_masses(size);
2920
2921 for (std::size_t i = 0u; i < size; ++i) {
2922 // Conversion from (pt, eta, phi, mass) to (x, y, z, e) coordinate system
2923 const auto x1 = pt1[i] * std::cos(phi1[i]);
2924 const auto y1 = pt1[i] * std::sin(phi1[i]);
2925 const auto z1 = pt1[i] * std::sinh(eta1[i]);
2926 const auto e1 = std::sqrt(x1 * x1 + y1 * y1 + z1 * z1 + mass1[i] * mass1[i]);
2927
2928 const auto x2 = pt2[i] * std::cos(phi2[i]);
2929 const auto y2 = pt2[i] * std::sin(phi2[i]);
2930 const auto z2 = pt2[i] * std::sinh(eta2[i]);
2931 const auto e2 = std::sqrt(x2 * x2 + y2 * y2 + z2 * z2 + mass2[i] * mass2[i]);
2932
2933 // Addition of particle four-vector elements
2934 const auto e = e1 + e2;
2935 const auto x = x1 + x2;
2936 const auto y = y1 + y2;
2937 const auto z = z1 + z2;
2938
2939 inv_masses[i] = std::sqrt(e * e - x * x - y * y - z * z);
2940 }
2941
2942 // Return invariant mass with (+, -, -, -) metric
2943 return inv_masses;
2944}
2945
2946/// Return the invariant mass of multiple particles given the collections of the
2947/// quantities transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
2948///
2949/// The function computes the invariant mass of multiple particles with the
2950/// four-vectors (pt, eta, phi, mass).
2951template <typename T>
2952T InvariantMass(const RVec<T>& pt, const RVec<T>& eta, const RVec<T>& phi, const RVec<T>& mass)
2953{
2954 const std::size_t size = pt.size();
2955
2956 R__ASSERT(eta.size() == size && phi.size() == size && mass.size() == size);
2957
2958 T x_sum = 0.;
2959 T y_sum = 0.;
2960 T z_sum = 0.;
2961 T e_sum = 0.;
2962
2963 for (std::size_t i = 0u; i < size; ++ i) {
2964 // Convert to (e, x, y, z) coordinate system and update sums
2965 const auto x = pt[i] * std::cos(phi[i]);
2966 x_sum += x;
2967 const auto y = pt[i] * std::sin(phi[i]);
2968 y_sum += y;
2969 const auto z = pt[i] * std::sinh(eta[i]);
2970 z_sum += z;
2971 const auto e = std::sqrt(x * x + y * y + z * z + mass[i] * mass[i]);
2972 e_sum += e;
2973 }
2974
2975 // Return invariant mass with (+, -, -, -) metric
2976 return std::sqrt(e_sum * e_sum - x_sum * x_sum - y_sum * y_sum - z_sum * z_sum);
2977}
2978
2979////////////////////////////////////////////////////////////////////////////
2980/// \brief Build an RVec of objects starting from RVecs of input to their constructors.
2981/// \tparam T Type of the objects contained in the created RVec.
2982/// \tparam Args_t Pack of types templating the input RVecs.
2983/// \param[in] args The RVecs containing the values used to initialise the output objects.
2984/// \return The RVec of objects initialised with the input parameters.
2985///
2986/// Example code, at the ROOT prompt:
2987/// ~~~{.cpp}
2988/// using namespace ROOT::VecOps;
2989/// RVecF pts = {15.5, 34.32, 12.95};
2990/// RVecF etas = {0.3, 2.2, 1.32};
2991/// RVecF phis = {0.1, 3.02, 2.2};
2992/// RVecF masses = {105.65, 105.65, 105.65};
2993/// auto fourVecs = Construct<ROOT::Math::PtEtaPhiMVector>(pts, etas, phis, masses);
2994/// cout << fourVecs << endl;
2995/// // { (15.5,0.3,0.1,105.65), (34.32,2.2,3.02,105.65), (12.95,1.32,2.2,105.65) }
2996/// ~~~
2997template <typename T, typename... Args_t>
2999{
3000 const auto size = ::ROOT::Internal::VecOps::GetVectorsSize("Construct", args...);
3001 RVec<T> ret;
3002 ret.reserve(size);
3003 for (auto i = 0UL; i < size; ++i) {
3004 ret.emplace_back(args[i]...);
3005 }
3006 return ret;
3007}
3008
3009/// Enumerate method:
3010/// For any Rvec v produce another RVec with entries starting from 0, and incrementing by 1 until a N = v.size() is reached.
3011/// Example code, at the ROOT prompt:
3012/// ~~~{.cpp}
3013/// using namespace ROOT::VecOps;
3014/// RVecF v = {1., 2., 3.};
3015/// cout << Enumerate(v1) << "\n";
3016/// // { 0, 1, 2 }
3017template <typename T>
3019{
3020 const auto size = v.size();
3021 RVec<T> ret;
3022 ret.reserve(size);
3023 for (auto i = 0UL; i < size; ++i) {
3024 ret.emplace_back(i);
3025 }
3026 return ret;
3027}
3028
3029/// Range method:
3030/// Produce RVec with entries starting from 0, and incrementing by 1 until a user-specified N is reached.
3031/// Example code, at the ROOT prompt:
3032/// ~~~{.cpp}
3033/// using namespace ROOT::VecOps;
3034/// cout << Range(3) << "\n";
3035/// // { 0, 1, 2 }
3036inline RVec<std::size_t> Range(std::size_t length)
3037{
3039 ret.reserve(length);
3040 for (auto i = 0UL; i < length; ++i) {
3041 ret.emplace_back(i);
3042 }
3043 return ret;
3044}
3045
3046////////////////////////////////////////////////////////////////////////////////
3047/// Print a RVec at the prompt:
3048template <class T>
3049std::ostream &operator<<(std::ostream &os, const RVec<T> &v)
3050{
3051 // In order to print properly, convert to 64 bit int if this is a char
3052 constexpr bool mustConvert = std::is_same<char, T>::value || std::is_same<signed char, T>::value ||
3056 os << "{ ";
3057 auto size = v.size();
3058 if (size) {
3059 for (std::size_t i = 0; i < size - 1; ++i) {
3060 os << (Print_t)v[i] << ", ";
3061 }
3062 os << (Print_t)v[size - 1];
3063 }
3064 os << " }";
3065 return os;
3066}
3067
3068#if (_VECOPS_USE_EXTERN_TEMPLATES)
3069
3070#define RVEC_EXTERN_UNARY_OPERATOR(T, OP) \
3071 extern template RVec<T> operator OP<T>(const RVec<T> &);
3072
3073#define RVEC_EXTERN_BINARY_OPERATOR(T, OP) \
3074 extern template auto operator OP<T, T>(const T &x, const RVec<T> &v) \
3075 -> RVec<decltype(x OP v[0])>; \
3076 extern template auto operator OP<T, T>(const RVec<T> &v, const T &y) \
3077 -> RVec<decltype(v[0] OP y)>; \
3078 extern template auto operator OP<T, T>(const RVec<T> &v0, const RVec<T> &v1)\
3079 -> RVec<decltype(v0[0] OP v1[0])>;
3080
3081#define RVEC_EXTERN_ASSIGN_OPERATOR(T, OP) \
3082 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const T &); \
3083 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const RVec<T> &);
3084
3085#define RVEC_EXTERN_LOGICAL_OPERATOR(T, OP) \
3086 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const T &); \
3087 extern template RVec<int> operator OP<T, T>(const T &, const RVec<T> &); \
3088 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const RVec<T> &);
3089
3090#define RVEC_EXTERN_FLOAT_TEMPLATE(T) \
3091 extern template class RVec<T>; \
3092 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
3093 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
3094 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
3095 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
3096 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
3097 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
3098 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
3099 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
3100 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
3101 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
3102 RVEC_EXTERN_ASSIGN_OPERATOR(T, /=) \
3103 RVEC_EXTERN_LOGICAL_OPERATOR(T, <) \
3104 RVEC_EXTERN_LOGICAL_OPERATOR(T, >) \
3105 RVEC_EXTERN_LOGICAL_OPERATOR(T, ==) \
3106 RVEC_EXTERN_LOGICAL_OPERATOR(T, !=) \
3107 RVEC_EXTERN_LOGICAL_OPERATOR(T, <=) \
3108 RVEC_EXTERN_LOGICAL_OPERATOR(T, >=) \
3109 RVEC_EXTERN_LOGICAL_OPERATOR(T, &&) \
3110 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
3111
3112#define RVEC_EXTERN_INTEGER_TEMPLATE(T) \
3113 extern template class RVec<T>; \
3114 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
3115 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
3116 RVEC_EXTERN_UNARY_OPERATOR(T, ~) \
3117 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
3118 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
3119 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
3120 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
3121 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
3122 RVEC_EXTERN_BINARY_OPERATOR(T, %) \
3123 RVEC_EXTERN_BINARY_OPERATOR(T, &) \
3124 RVEC_EXTERN_BINARY_OPERATOR(T, |) \
3125 RVEC_EXTERN_BINARY_OPERATOR(T, ^) \
3126 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
3127 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
3128 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
3129 RVEC_EXTERN_ASSIGN_OPERATOR(T, /=) \
3130 RVEC_EXTERN_ASSIGN_OPERATOR(T, %=) \
3131 RVEC_EXTERN_ASSIGN_OPERATOR(T, &=) \
3132 RVEC_EXTERN_ASSIGN_OPERATOR(T, |=) \
3133 RVEC_EXTERN_ASSIGN_OPERATOR(T, ^=) \
3134 RVEC_EXTERN_ASSIGN_OPERATOR(T, >>=) \
3135 RVEC_EXTERN_ASSIGN_OPERATOR(T, <<=) \
3136 RVEC_EXTERN_LOGICAL_OPERATOR(T, <) \
3137 RVEC_EXTERN_LOGICAL_OPERATOR(T, >) \
3138 RVEC_EXTERN_LOGICAL_OPERATOR(T, ==) \
3139 RVEC_EXTERN_LOGICAL_OPERATOR(T, !=) \
3140 RVEC_EXTERN_LOGICAL_OPERATOR(T, <=) \
3141 RVEC_EXTERN_LOGICAL_OPERATOR(T, >=) \
3142 RVEC_EXTERN_LOGICAL_OPERATOR(T, &&) \
3143 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
3144
3145RVEC_EXTERN_INTEGER_TEMPLATE(char)
3146RVEC_EXTERN_INTEGER_TEMPLATE(short)
3147RVEC_EXTERN_INTEGER_TEMPLATE(int)
3148RVEC_EXTERN_INTEGER_TEMPLATE(long)
3149//RVEC_EXTERN_INTEGER_TEMPLATE(long long)
3150
3151RVEC_EXTERN_INTEGER_TEMPLATE(unsigned char)
3152RVEC_EXTERN_INTEGER_TEMPLATE(unsigned short)
3153RVEC_EXTERN_INTEGER_TEMPLATE(unsigned int)
3154RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long)
3155//RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long long)
3156
3157RVEC_EXTERN_FLOAT_TEMPLATE(float)
3158RVEC_EXTERN_FLOAT_TEMPLATE(double)
3159
3160#undef RVEC_EXTERN_UNARY_OPERATOR
3161#undef RVEC_EXTERN_BINARY_OPERATOR
3162#undef RVEC_EXTERN_ASSIGN_OPERATOR
3163#undef RVEC_EXTERN_LOGICAL_OPERATOR
3164#undef RVEC_EXTERN_INTEGER_TEMPLATE
3165#undef RVEC_EXTERN_FLOAT_TEMPLATE
3166
3167#define RVEC_EXTERN_UNARY_FUNCTION(T, NAME, FUNC) \
3168 extern template RVec<PromoteType<T>> NAME(const RVec<T> &);
3169
3170#define RVEC_EXTERN_STD_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, std::F)
3171
3172#define RVEC_EXTERN_BINARY_FUNCTION(T0, T1, NAME, FUNC) \
3173 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const T1 &); \
3174 extern template RVec<PromoteTypes<T0, T1>> NAME(const T0 &, const RVec<T1> &); \
3175 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const RVec<T1> &);
3176
3177#define RVEC_EXTERN_STD_BINARY_FUNCTION(T, F) RVEC_EXTERN_BINARY_FUNCTION(T, T, F, std::F)
3178
3179#define RVEC_EXTERN_STD_FUNCTIONS(T) \
3180 RVEC_EXTERN_STD_UNARY_FUNCTION(T, abs) \
3181 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fdim) \
3182 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fmod) \
3183 RVEC_EXTERN_STD_BINARY_FUNCTION(T, remainder) \
3184 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp) \
3185 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp2) \
3186 RVEC_EXTERN_STD_UNARY_FUNCTION(T, expm1) \
3187 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log) \
3188 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log10) \
3189 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log2) \
3190 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log1p) \
3191 RVEC_EXTERN_STD_BINARY_FUNCTION(T, pow) \
3192 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sqrt) \
3193 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cbrt) \
3194 RVEC_EXTERN_STD_BINARY_FUNCTION(T, hypot) \
3195 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sin) \
3196 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cos) \
3197 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tan) \
3198 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asin) \
3199 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acos) \
3200 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atan) \
3201 RVEC_EXTERN_STD_BINARY_FUNCTION(T, atan2) \
3202 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sinh) \
3203 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cosh) \
3204 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tanh) \
3205 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asinh) \
3206 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acosh) \
3207 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atanh) \
3208 RVEC_EXTERN_STD_UNARY_FUNCTION(T, floor) \
3209 RVEC_EXTERN_STD_UNARY_FUNCTION(T, ceil) \
3210 RVEC_EXTERN_STD_UNARY_FUNCTION(T, trunc) \
3211 RVEC_EXTERN_STD_UNARY_FUNCTION(T, round) \
3212 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erf) \
3213 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erfc) \
3214 RVEC_EXTERN_STD_UNARY_FUNCTION(T, lgamma) \
3215 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tgamma) \
3216
3217RVEC_EXTERN_STD_FUNCTIONS(float)
3218RVEC_EXTERN_STD_FUNCTIONS(double)
3219#undef RVEC_EXTERN_STD_UNARY_FUNCTION
3220#undef RVEC_EXTERN_STD_BINARY_FUNCTION
3221#undef RVEC_EXTERN_STD_UNARY_FUNCTIONS
3222
3223#ifdef R__HAS_VDT
3224
3225#define RVEC_EXTERN_VDT_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, vdt::F)
3226
3227RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_expf)
3228RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_logf)
3229RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_sinf)
3230RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_cosf)
3231RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_tanf)
3232RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_asinf)
3233RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_acosf)
3234RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_atanf)
3235
3236RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_exp)
3237RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_log)
3238RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_sin)
3239RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_cos)
3240RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_tan)
3241RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_asin)
3242RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_acos)
3243RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_atan)
3244
3245#endif // R__HAS_VDT
3246
3247#endif // _VECOPS_USE_EXTERN_TEMPLATES
3248
3249/** @} */ // end of Doxygen group vecops
3250
3251} // End of VecOps NS
3252
3253// Allow to use RVec as ROOT::RVec
3254using ROOT::VecOps::RVec;
3255
3266
3267} // End of ROOT NS
3268
3269#endif // ROOT_RVEC
size_t fSize
#define R__unlikely(expr)
Definition: RConfig.hxx:595
#define f(i)
Definition: RSha256.hxx:104
#define c(i)
Definition: RSha256.hxx:101
#define e(i)
Definition: RSha256.hxx:103
#define R__RVEC_NODISCARD
Definition: RVec.hxx:19
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#define M_PI
Definition: Rotated.cxx:105
#define R__CLING_PTRCHECK(ONOFF)
Definition: Rtypes.h:498
static Double_t Product(const Double_t *x, const Float_t *y)
Product.
Definition: TCTUB.cxx:101
#define R__ASSERT(e)
Definition: TError.h:118
#define N
Double_t Dot(const TGLVector3 &v1, const TGLVector3 &v2)
Definition: TGLUtil.h:317
Int_t Compare(const void *item1, const void *item2)
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h length
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char y2
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
Option_t Option_t TPoint TPoint const char y1
Binding & operator=(OUT(*fun)(void))
#define free
Definition: civetweb.c:1539
#define malloc
Definition: civetweb.c:1536
This class consists of common code factored out of the SmallVector class to reduce code duplication b...
Definition: RVec.hxx:536
void assign(size_type NumElts, const T &Elt)
Definition: RVec.hxx:655
void append(in_iter in_start, in_iter in_end)
Add the specified range to the end of the SmallVector.
Definition: RVec.hxx:629
iterator insert(iterator I, T &&Elt)
Definition: RVec.hxx:716
void resize(size_type N)
Definition: RVec.hxx:571
void assign(std::initializer_list< T > IL)
Definition: RVec.hxx:673
void resize(size_type N, const T &NV)
Definition: RVec.hxx:586
void reserve(size_type N)
Definition: RVec.hxx:600
iterator insert(iterator I, ItTy From, ItTy To)
Definition: RVec.hxx:834
reference emplace_back(ArgTypes &&...Args)
Definition: RVec.hxx:895
void assign(in_iter in_start, in_iter in_end)
Definition: RVec.hxx:667
iterator insert(iterator I, const T &Elt)
Definition: RVec.hxx:748
void swap(RVecImpl &RHS)
Definition: RVec.hxx:910
iterator insert(iterator I, size_type NumToInsert, const T &Elt)
Definition: RVec.hxx:779
RVecImpl & operator=(const RVecImpl &RHS)
Definition: RVec.hxx:970
iterator erase(const_iterator CS, const_iterator CE)
Definition: RVec.hxx:696
void append(size_type NumInputs, const T &Elt)
Append NumInputs copies of Elt to the end.
Definition: RVec.hxx:640
iterator erase(const_iterator CI)
Definition: RVec.hxx:679
RVecImpl & operator=(RVecImpl &&RHS)
Definition: RVec.hxx:1023
void pop_back_n(size_type NumItems)
Definition: RVec.hxx:606
RVecImpl(const RVecImpl &)=delete
void append(std::initializer_list< T > IL)
Definition: RVec.hxx:649
void insert(iterator I, std::initializer_list< T > IL)
Definition: RVec.hxx:892
This is all the stuff common to all SmallVectors.
Definition: RVec.hxx:137
SmallVectorBase(void *FirstEl, size_t TotalCapacity)
Definition: RVec.hxx:155
static constexpr size_t SizeTypeMax()
The maximum value of the Size_T used.
Definition: RVec.hxx:152
Size_T fCapacity
Always >= -1. fCapacity == -1 indicates the RVec is in "memory adoption" mode.
Definition: RVec.hxx:149
bool Owns() const
If true, the RVec is in "memory adoption" mode, i.e. it is acting as a view on a memory buffer it doe...
Definition: RVec.hxx:170
size_t capacity() const noexcept
Definition: RVec.hxx:174
void set_size(size_t N)
Set the array size to N, which the current array must have enough capacity for.
Definition: RVec.hxx:187
void grow(size_t MinSize=0)
Double the size of the allocated memory, guaranteeing space for at least one more element or MinSize ...
Definition: RVec.hxx:471
static void uninitialized_move(It1 I, It1 E, It2 Dest)
Move the range [I, E) onto the uninitialized memory starting with "Dest", constructing elements into ...
Definition: RVec.hxx:439
static void uninitialized_copy(T1 *I, T1 *E, T2 *Dest, typename std::enable_if< std::is_same< typename std::remove_const< T1 >::type, T2 >::value >::type *=nullptr)
Copy the range [I, E) onto the uninitialized memory starting with "Dest", constructing elements into ...
Definition: RVec.hxx:457
static void uninitialized_copy(It1 I, It1 E, It2 Dest)
Copy the range [I, E) onto the uninitialized memory starting with "Dest", constructing elements into ...
Definition: RVec.hxx:448
SmallVectorTemplateBase<TriviallyCopyable = false> - This is where we put method implementations that...
Definition: RVec.hxx:327
void grow(size_t MinSize=0)
Grow the allocated memory (without initializing new elements), doubling the size of the allocated mem...
static void uninitialized_move(It1 I, It1 E, It2 Dest)
Move the range [I, E) into the uninitialized memory starting with "Dest", constructing elements as ne...
Definition: RVec.hxx:342
static void uninitialized_copy(It1 I, It1 E, It2 Dest)
Copy the range [I, E) onto the uninitialized memory starting with "Dest", constructing elements as ne...
Definition: RVec.hxx:350
This is the part of SmallVectorTemplateBase which does not depend on whether the type T is a POD.
Definition: RVec.hxx:205
const_iterator cbegin() const noexcept
Definition: RVec.hxx:260
void grow_pod(size_t MinSize, size_t TSize)
Definition: RVec.hxx:221
reverse_iterator rbegin() noexcept
Definition: RVec.hxx:266
const_iterator cend() const noexcept
Definition: RVec.hxx:263
void resetToSmall()
Put this vector in a state of being small.
Definition: RVec.hxx:228
std::reverse_iterator< iterator > reverse_iterator
Definition: RVec.hxx:246
bool isSmall() const
Return true if this is a smallvector which has not had dynamic memory allocated for it.
Definition: RVec.hxx:225
const_reverse_iterator crend() const noexcept
Definition: RVec.hxx:271
const_iterator end() const noexcept
Definition: RVec.hxx:262
const_reverse_iterator crbegin() const noexcept
Definition: RVec.hxx:268
pointer data() noexcept
Return a pointer to the vector's buffer, even if empty().
Definition: RVec.hxx:279
const_reverse_iterator rbegin() const noexcept
Definition: RVec.hxx:267
std::reverse_iterator< const_iterator > const_reverse_iterator
Definition: RVec.hxx:245
const_iterator begin() const noexcept
Definition: RVec.hxx:259
const_pointer data() const noexcept
Return a pointer to the vector's buffer, even if empty().
Definition: RVec.hxx:281
void * getFirstEl() const
Find the address of the first element.
Definition: RVec.hxx:211
const_reverse_iterator rend() const noexcept
Definition: RVec.hxx:270
RVecN(size_t Size)
Definition: RVec.hxx:1140
RVecN(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1176
reference operator[](size_type idx)
Definition: RVec.hxx:1216
RVecN operator[](const RVecN< V, M > &conds) const
Definition: RVec.hxx:1227
RVecN(std::initializer_list< T > IL)
Definition: RVec.hxx:1156
const_reference at(size_type pos) const
Definition: RVec.hxx:1262
RVecN(const RVecN &RHS)
Definition: RVec.hxx:1158
RVecN & operator=(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1197
RVecN & operator=(RVecN &&RHS)
Definition: RVec.hxx:1184
value_type at(size_type pos, value_type fallback) const
No exception thrown. The user specifies the desired value in case the RVecN is shorter than pos.
Definition: RVec.hxx:1278
RVecN & operator=(std::initializer_list< T > IL)
Definition: RVec.hxx:1203
RVecN & operator=(const RVecN &RHS)
Definition: RVec.hxx:1164
RVecN(const std::vector< T > &RHS)
Definition: RVec.hxx:1182
RVecN(size_t Size, const T &Value)
Definition: RVec.hxx:1138
RVecN(RVecN &&RHS)
Definition: RVec.hxx:1170
RVecN(ItTy S, ItTy E)
Definition: RVec.hxx:1151
reference at(size_type pos)
Definition: RVec.hxx:1255
value_type at(size_type pos, value_type fallback)
No exception thrown. The user specifies the desired value in case the RVecN is shorter than pos.
Definition: RVec.hxx:1270
RVecN(T *p, size_t n)
Definition: RVec.hxx:1190
const_reference operator[](size_type idx) const
Definition: RVec.hxx:1221
A "std::vector"-like collection of values implementing handy operation to analyse them.
Definition: RVec.hxx:1456
RVec(RVecN< T, N > &&RHS)
Definition: RVec.hxx:1500
RVec(const RVecN< T, N > &RHS)
Definition: RVec.hxx:1503
RVec(size_t Size, const T &Value)
Definition: RVec.hxx:1468
RVec(const RVec &RHS)
Definition: RVec.hxx:1481
RVec & operator=(RVec &&RHS)
Definition: RVec.hxx:1491
RVec(T *p, size_t n)
Definition: RVec.hxx:1507
RVec operator[](const RVec< V > &conds) const
Definition: RVec.hxx:1519
RVec(std::initializer_list< T > IL)
Definition: RVec.hxx:1479
RVec(size_t Size)
Definition: RVec.hxx:1470
RVec(ItTy S, ItTy E)
Definition: RVec.hxx:1475
RVec(const std::vector< T > &RHS)
Definition: RVec.hxx:1505
typename SuperClass::size_type size_type
Definition: RVec.hxx:1461
RVec(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1497
RVec(RVec &&RHS)
Definition: RVec.hxx:1489
typename SuperClass::value_type value_type
Definition: RVec.hxx:1462
RVec & operator=(const RVec &RHS)
Definition: RVec.hxx:1483
TPaveText * pt
Type
enumeration specifying the integration types.
T Max(const RVec< T > &v)
Get the greatest element of an RVec.
Definition: RVec.hxx:1953
RVec< T > Reverse(const RVec< T > &v)
Return copy of reversed vector.
Definition: RVec.hxx:2363
RVec< T > DeltaPhi(T v1, const RVec< T > &v2, const T c=M_PI)
Return the angle difference in radians of a scalar and a vector.
Definition: RVec.hxx:2852
RVec< T > Where(const RVec< int > &c, T v1, T v2)
Return a vector with the value v2 if the condition c is false and sets the value v1 if the condition ...
Definition: RVec.hxx:2756
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:2640
RVec< std::size_t > Range(std::size_t length)
Enumerate method: For any Rvec v produce another RVec with entries starting from 0,...
Definition: RVec.hxx:3036
RVec< RVec< typename RVec< T >::size_type > > Combinations(const RVec< T > &v, const typename RVec< T >::size_type n)
Return the indices that represent all unique combinations of the elements of a given RVec.
Definition: RVec.hxx:2545
RVec< PromoteType< T > > asinh(const RVec< T > &v)
Definition: RVec.hxx:1786
RVec< typename RVec< T >::size_type > Argsort(const RVec< T > &v, Compare &&c)
Return an RVec of indices that sort the input RVec based on a comparison function.
Definition: RVec.hxx:2191
RVec< PromoteType< T > > expm1(const RVec< T > &v)
Definition: RVec.hxx:1763
RVec< PromoteType< T > > lgamma(const RVec< T > &v)
Definition: RVec.hxx:1799
RVec< typename RVec< T >::size_type > Nonzero(const RVec< T > &v)
Return the indices of the elements which are not zero.
Definition: RVec.hxx:2609
#define RVEC_UNARY_OPERATOR(OP)
Definition: RVec.hxx:1540
RVec< PromoteType< T > > cosh(const RVec< T > &v)
Definition: RVec.hxx:1784
R Mean(const RVec< T > &v, const R zero)
Get the mean of the elements of an RVec with custom initial value.
Definition: RVec.hxx:1936
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1756
#define RVEC_ASSIGNMENT_OPERATOR(OP)
Definition: RVec.hxx:1611
RVec< PromoteType< T > > round(const RVec< T > &v)
Definition: RVec.hxx:1793
T InvariantMass(const RVec< T > &pt, const RVec< T > &eta, const RVec< T > &phi, const RVec< T > &mass)
Return the invariant mass of multiple particles given the collections of the quantities transverse mo...
Definition: RVec.hxx:2952
RVec< PromoteType< T > > tgamma(const RVec< T > &v)
Definition: RVec.hxx:1800
RVec< PromoteType< T > > log2(const RVec< T > &v)
Definition: RVec.hxx:1767
RVec< Common_t > Concatenate(const RVec< T0 > &v0, const RVec< T1 > &v1)
Return the concatenation of two RVecs.
Definition: RVec.hxx:2779
RVec< PromoteType< T > > tan(const RVec< T > &v)
Definition: RVec.hxx:1777
RVec< PromoteType< T > > erf(const RVec< T > &v)
Definition: RVec.hxx:1797
RVec< T > StableSort(const RVec< T > &v, Compare &&c)
Return copy of RVec with elements sorted based on a comparison operator while keeping the order of eq...
Definition: RVec.hxx:2469
RVec< PromoteType< T > > cos(const RVec< T > &v)
Definition: RVec.hxx:1776
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:2910
RVec< PromoteType< T > > floor(const RVec< T > &v)
Definition: RVec.hxx:1790
RVec< PromoteType< T > > atanh(const RVec< T > &v)
Definition: RVec.hxx:1788
RVec< PromoteType< T > > acosh(const RVec< T > &v)
Definition: RVec.hxx:1787
RVec< T > DeltaR2(const RVec< T > &eta1, const RVec< T > &eta2, const RVec< T > &phi1, const RVec< T > &phi2, const T c=M_PI)
Return the square of the distance on the - plane ( ) from the collections eta1, eta2,...
Definition: RVec.hxx:2871
RVec< PromoteTypes< T0, T1 > > remainder(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1759
RVec< PromoteType< T > > ceil(const RVec< T > &v)
Definition: RVec.hxx:1791
RVec< PromoteType< T > > cbrt(const RVec< T > &v)
Definition: RVec.hxx:1772
T Sum(const RVec< T > &v, const T zero=T(0))
Sum elements of an RVec.
Definition: RVec.hxx:1878
RVec< PromoteType< T > > asin(const RVec< T > &v)
Definition: RVec.hxx:1778
void swap(RVec< T > &lhs, RVec< T > &rhs)
Definition: RVec.hxx:2153
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1765
RVec< T > Construct(const RVec< Args_t > &... args)
Build an RVec of objects starting from RVecs of input to their constructors.
Definition: RVec.hxx:2998
RVec< PromoteType< T > > trunc(const RVec< T > &v)
Definition: RVec.hxx:1792
RVec< PromoteType< T > > llround(const RVec< T > &v)
Definition: RVec.hxx:1795
RVec< PromoteTypes< T0, T1 > > pow(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1770
RVec< PromoteTypes< T0, T1 > > hypot(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1773
#define RVEC_STD_BINARY_FUNCTION(F)
Definition: RVec.hxx:1754
RVec< PromoteType< T > > lround(const RVec< T > &v)
Definition: RVec.hxx:1794
#define RVEC_BINARY_OPERATOR(OP)
Definition: RVec.hxx:1563
RVec< T > Drop(const RVec< T > &v, RVec< typename RVec< T >::size_type > idxs)
Return a copy of the container without the elements at the specified indices.
Definition: RVec.hxx:2330
RVec< typename RVec< T >::size_type > StableArgsort(const RVec< T > &v, Compare &&c)
Return an RVec of indices that sort the input RVec based on a comparison function while keeping the o...
Definition: RVec.hxx:2238
T Min(const RVec< T > &v)
Get the smallest element of an RVec.
Definition: RVec.hxx:1969
size_t CapacityInBytes(const RVecN< T, N > &X)
Definition: RVec.hxx:1532
#define RVEC_LOGICAL_OPERATOR(OP)
Definition: RVec.hxx:1647
RVec< T > Take(const RVec< T > &v, const int n)
Return first or last n elements of an RVec.
Definition: RVec.hxx:2305
RVec< PromoteType< T > > sinh(const RVec< T > &v)
Definition: RVec.hxx:1783
RVec< PromoteType< T > > sqrt(const RVec< T > &v)
Definition: RVec.hxx:1771
RVec< PromoteType< T > > erfc(const RVec< T > &v)
Definition: RVec.hxx:1798
RVec< T > Sort(const RVec< T > &v, Compare &&c)
Return copy of RVec with elements sorted based on a comparison operator.
Definition: RVec.hxx:2409
#define RVEC_STD_UNARY_FUNCTION(F)
Definition: RVec.hxx:1753
RVec< PromoteTypes< T0, T1 > > atan2(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1781
RVec< PromoteTypes< T0, T1 > > fmod(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1758
RVec< PromoteType< T > > log10(const RVec< T > &v)
Definition: RVec.hxx:1766
auto Map(Args &&... args)
Create new collection applying a callable to the elements of the input collection.
Definition: RVec.hxx:2070
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1761
RVec< PromoteTypes< T0, T1 > > fdim(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1757
double StdDev(const RVec< T > &v)
Get the standard deviation of the elements of an RVec.
Definition: RVec.hxx:2046
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:2125
RVec< PromoteType< T > > tanh(const RVec< T > &v)
Definition: RVec.hxx:1785
RVec< PromoteType< T > > acos(const RVec< T > &v)
Definition: RVec.hxx:1779
RVec< PromoteType< T > > log1p(const RVec< T > &v)
Definition: RVec.hxx:1768
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:2005
RVec< PromoteType< T > > atan(const RVec< T > &v)
Definition: RVec.hxx:1780
RVec< PromoteType< T > > sin(const RVec< T > &v)
Definition: RVec.hxx:1775
double Var(const RVec< T > &v)
Get the variance of the elements of an RVec.
Definition: RVec.hxx:2022
RVec< PromoteType< T > > exp2(const RVec< T > &v)
Definition: RVec.hxx:1762
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:2102
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:1987
T DeltaR(T eta1, T eta2, T phi1, T phi2, const T c=M_PI)
Return the distance on the - plane ( ) from the scalars eta1, eta2, phi1 and phi2.
Definition: RVec.hxx:2898
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 T2
Definition: md5.inl:147
#define F(x, y, z)
#define I(x, y, z)
#define T1
Definition: md5.inl:146
bool IsSmall(const ROOT::VecOps::RVec< T > &v)
Definition: RVec.hxx:1093
bool IsAdopting(const ROOT::VecOps::RVec< T > &v)
Definition: RVec.hxx:1099
auto MapImpl(F &&f, RVecs &&... vs) -> RVec< decltype(f(vs[0]...))>
Definition: RVec.hxx:104
uint64_t NextPowerOf2(uint64_t A)
Return the next power of two (in 64-bits) that is strictly greater than A.
Definition: RVec.hxx:125
constexpr bool All(const bool *vals, std::size_t size)
Definition: RVec.hxx:78
std::size_t GetVectorsSize(const std::string &id, const RVec< T > &... vs)
Definition: RVec.hxx:87
ROOT::VecOps::RVec< T > RVec
Definition: RVec.hxx:68
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
std::ostream & operator<<(std::ostream &os, const RConcurrentHashColl::HashValue &h)
static double A[]
double T(double x)
Definition: ChebyshevPol.h:34
double inner_product(const LAVector &, const LAVector &)
This file contains a specialised ROOT message handler to test for diagnostic in unit tests.
__roodevice__ double fast_exp(double x)
Definition: RooVDTHeaders.h:53
void init()
Inspect hardware capabilities, and load the optimal library for RooFit computations.
__roodevice__ double fast_log(double x)
Definition: RooVDTHeaders.h:57
RooArgSet S(Args_t &&... args)
Definition: RooArgSet.h:240
static constexpr double s
constexpr Double_t E()
Base of natural log: .
Definition: TMath.h:93
const char * Size
Definition: TXMLSetup.cxx:56
The size of the inline storage of an RVec.
Definition: RVec.hxx:511
Used to figure out the offset of the first element of an RVec.
Definition: RVec.hxx:198
Storage for the SmallVector elements.
Definition: RVec.hxx:496
TMarker m
Definition: textangle.C:8