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