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