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