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
906 // This block handles the swap of a small and a non-owning vector
907 // It is more efficient to first move the non-owning vector, hence the 2 cases
908 if (this->isSmall() && !RHS.Owns()) { // the right vector is non-owning
909 RVecImpl<T> temp(0);
910 temp = std::move(RHS);
911 RHS = std::move(*this);
912 *this = std::move(temp);
913 return;
914 } else if (RHS.isSmall() && !this->Owns()) { // the left vector is non-owning
915 RVecImpl<T> temp(0);
916 temp = std::move(*this);
917 *this = std::move(RHS);
918 RHS = std::move(temp);
919 return;
920 }
921
922 if (RHS.size() > this->capacity())
923 this->grow(RHS.size());
924 if (this->size() > RHS.capacity())
925 RHS.grow(this->size());
926
927 // Swap the shared elements.
928 size_t NumShared = this->size();
929 if (NumShared > RHS.size())
930 NumShared = RHS.size();
931 for (size_type i = 0; i != NumShared; ++i)
932 std::iter_swap(this->begin() + i, RHS.begin() + i);
933
934 // Copy over the extra elts.
935 if (this->size() > RHS.size()) {
936 size_t EltDiff = this->size() - RHS.size();
937 this->uninitialized_copy(this->begin() + NumShared, this->end(), RHS.end());
938 RHS.set_size(RHS.size() + EltDiff);
939 if (this->Owns())
940 this->destroy_range(this->begin() + NumShared, this->end());
941 this->set_size(NumShared);
942 } else if (RHS.size() > this->size()) {
943 size_t EltDiff = RHS.size() - this->size();
944 this->uninitialized_copy(RHS.begin() + NumShared, RHS.end(), this->end());
945 this->set_size(this->size() + EltDiff);
946 if (RHS.Owns())
947 this->destroy_range(RHS.begin() + NumShared, RHS.end());
948 RHS.set_size(NumShared);
949 }
950}
951
952template <typename T>
954{
955 // Avoid self-assignment.
956 if (this == &RHS)
957 return *this;
958
959 // If we already have sufficient space, assign the common elements, then
960 // destroy any excess.
961 size_t RHSSize = RHS.size();
962 size_t CurSize = this->size();
963 if (CurSize >= RHSSize) {
964 // Assign common elements.
965 iterator NewEnd;
966 if (RHSSize)
967 NewEnd = std::copy(RHS.begin(), RHS.begin() + RHSSize, this->begin());
968 else
969 NewEnd = this->begin();
970
971 // Destroy excess elements.
972 if (this->Owns())
973 this->destroy_range(NewEnd, this->end());
974
975 // Trim.
976 this->set_size(RHSSize);
977 return *this;
978 }
979
980 // If we have to grow to have enough elements, destroy the current elements.
981 // This allows us to avoid copying them during the grow.
982 // From the original LLVM implementation:
983 // FIXME: don't do this if they're efficiently moveable.
984 if (this->capacity() < RHSSize) {
985 if (this->Owns()) {
986 // Destroy current elements.
987 this->destroy_range(this->begin(), this->end());
988 }
989 this->set_size(0);
990 CurSize = 0;
991 this->grow(RHSSize);
992 } else if (CurSize) {
993 // Otherwise, use assignment for the already-constructed elements.
994 std::copy(RHS.begin(), RHS.begin() + CurSize, this->begin());
995 }
996
997 // Copy construct the new elements in place.
998 this->uninitialized_copy(RHS.begin() + CurSize, RHS.end(), this->begin() + CurSize);
999
1000 // Set end.
1001 this->set_size(RHSSize);
1002 return *this;
1003}
1004
1005template <typename T>
1007{
1008 // Avoid self-assignment.
1009 if (this == &RHS)
1010 return *this;
1011
1012 // If the RHS isn't small, clear this vector and then steal its buffer.
1013 if (!RHS.isSmall()) {
1014 if (this->Owns()) {
1015 this->destroy_range(this->begin(), this->end());
1016 if (!this->isSmall())
1017 free(this->begin());
1018 }
1019 this->fBeginX = RHS.fBeginX;
1020 this->fSize = RHS.fSize;
1021 this->fCapacity = RHS.fCapacity;
1022 RHS.resetToSmall();
1023 return *this;
1024 }
1025
1026 // If we already have sufficient space, assign the common elements, then
1027 // destroy any excess.
1028 size_t RHSSize = RHS.size();
1029 size_t CurSize = this->size();
1030 if (CurSize >= RHSSize) {
1031 // Assign common elements.
1032 iterator NewEnd = this->begin();
1033 if (RHSSize)
1034 NewEnd = std::move(RHS.begin(), RHS.end(), NewEnd);
1035
1036 // Destroy excess elements and trim the bounds.
1037 if (this->Owns())
1038 this->destroy_range(NewEnd, this->end());
1039 this->set_size(RHSSize);
1040
1041 // Clear the RHS.
1042 RHS.clear();
1043
1044 return *this;
1045 }
1046
1047 // If we have to grow to have enough elements, destroy the current elements.
1048 // This allows us to avoid copying them during the grow.
1049 // From the original LLVM implementation:
1050 // FIXME: this may not actually make any sense if we can efficiently move
1051 // elements.
1052 if (this->capacity() < RHSSize) {
1053 if (this->Owns()) {
1054 // Destroy current elements.
1055 this->destroy_range(this->begin(), this->end());
1056 }
1057 this->set_size(0);
1058 CurSize = 0;
1059 this->grow(RHSSize);
1060 } else if (CurSize) {
1061 // Otherwise, use assignment for the already-constructed elements.
1062 std::move(RHS.begin(), RHS.begin() + CurSize, this->begin());
1063 }
1064
1065 // Move-construct the new elements in place.
1066 this->uninitialized_move(RHS.begin() + CurSize, RHS.end(), this->begin() + CurSize);
1067
1068 // Set end.
1069 this->set_size(RHSSize);
1070
1071 RHS.clear();
1072 return *this;
1073}
1074
1075template <typename T>
1077{
1078 return v.isSmall();
1079}
1080
1081template <typename T>
1083{
1084 return !v.Owns();
1085}
1086
1087} // namespace VecOps
1088} // namespace Detail
1089
1090namespace VecOps {
1091// Note that we open here with @{ the Doxygen group vecops and it is
1092// closed again at the end of the C++ namespace VecOps
1093/**
1094 * \defgroup vecops VecOps
1095 * A "std::vector"-like collection of values implementing handy operation to analyse them
1096 * @{
1097*/
1098
1099// From the original SmallVector code:
1100// This is a 'vector' (really, a variable-sized array), optimized
1101// for the case when the array is small. It contains some number of elements
1102// in-place, which allows it to avoid heap allocation when the actual number of
1103// elements is below that threshold. This allows normal "small" cases to be
1104// fast without losing generality for large inputs.
1105//
1106// Note that this does not attempt to be exception safe.
1107
1108template <typename T, unsigned int N>
1109class R__CLING_PTRCHECK(off) RVecN : public Detail::VecOps::RVecImpl<T>, Internal::VecOps::SmallVectorStorage<T, N> {
1110public:
1111 RVecN() : Detail::VecOps::RVecImpl<T>(N) {}
1112
1114 {
1115 if (this->Owns()) {
1116 // Destroy the constructed elements in the vector.
1117 this->destroy_range(this->begin(), this->end());
1118 }
1119 }
1120
1121 explicit RVecN(size_t Size, const T &Value) : Detail::VecOps::RVecImpl<T>(N) { this->assign(Size, Value); }
1122
1123 explicit RVecN(size_t Size) : Detail::VecOps::RVecImpl<T>(N)
1124 {
1125 if (Size > N)
1126 this->grow(Size);
1127 this->fSize = Size;
1128 std::uninitialized_fill(this->begin(), this->end(), T{});
1129 }
1130
1131 template <typename ItTy,
1132 typename = typename std::enable_if<std::is_convertible<
1133 typename std::iterator_traits<ItTy>::iterator_category, std::input_iterator_tag>::value>::type>
1134 RVecN(ItTy S, ItTy E) : Detail::VecOps::RVecImpl<T>(N)
1135 {
1136 this->append(S, E);
1137 }
1138
1139 RVecN(std::initializer_list<T> IL) : Detail::VecOps::RVecImpl<T>(N) { this->assign(IL); }
1140
1141 RVecN(const RVecN &RHS) : Detail::VecOps::RVecImpl<T>(N)
1142 {
1143 if (!RHS.empty())
1145 }
1146
1147 RVecN &operator=(const RVecN &RHS)
1148 {
1150 return *this;
1151 }
1152
1153 RVecN(RVecN &&RHS) : Detail::VecOps::RVecImpl<T>(N)
1154 {
1155 if (!RHS.empty())
1157 }
1158
1159 RVecN(Detail::VecOps::RVecImpl<T> &&RHS) : Detail::VecOps::RVecImpl<T>(N)
1160 {
1161 if (!RHS.empty())
1163 }
1164
1165 RVecN(const std::vector<T> &RHS) : RVecN(RHS.begin(), RHS.end()) {}
1166
1168 {
1170 return *this;
1171 }
1172
1173 RVecN(T* p, size_t n) : Detail::VecOps::RVecImpl<T>(N)
1174 {
1175 this->fBeginX = p;
1176 this->fSize = n;
1177 this->fCapacity = -1;
1178 }
1179
1181 {
1183 return *this;
1184 }
1185
1186 RVecN &operator=(std::initializer_list<T> IL)
1187 {
1188 this->assign(IL);
1189 return *this;
1190 }
1191
1198
1200 {
1201 return begin()[idx];
1202 }
1203
1205 {
1206 return begin()[idx];
1207 }
1208
1210 RVecN operator[](const RVecN<V, M> &conds) const
1211 {
1212 const size_type n = conds.size();
1213
1214 if (n != this->size())
1215 throw std::runtime_error("Cannot index RVecN with condition vector of different size");
1216
1217 RVecN ret;
1218 ret.reserve(n);
1219 size_type j = 0u;
1220 for (size_type i = 0u; i < n; ++i) {
1221 if (conds[i]) {
1222 // the begin() is to go around the R__ASSERT in operator[]
1223 ret.begin()[j] = this->operator[](i);
1224 ++j;
1225 }
1226 }
1227 ret.set_size(j);
1228 return ret;
1229 }
1230
1231 // conversion
1233 operator RVecN<U, M>() const
1234 {
1235 return RVecN<U, M>(this->begin(), this->end());
1236 }
1237
1239 {
1240 if (pos >= size_type(this->fSize))
1241 throw std::out_of_range("RVecN");
1242 return this->operator[](pos);
1243 }
1244
1246 {
1247 if (pos >= size_type(this->fSize))
1248 throw std::out_of_range("RVecN");
1249 return this->operator[](pos);
1250 }
1251
1252 /// No exception thrown. The user specifies the desired value in case the RVecN is shorter than `pos`.
1254 {
1255 if (pos >= size_type(this->fSize))
1256 return fallback;
1257 return this->operator[](pos);
1258 }
1259
1260 /// No exception thrown. The user specifies the desired value in case the RVecN is shorter than `pos`.
1261 value_type at(size_type pos, value_type fallback) const
1262 {
1263 if (pos >= size_type(this->fSize))
1264 return fallback;
1265 return this->operator[](pos);
1266 }
1267};
1268
1269// clang-format off
1270/**
1271\class ROOT::VecOps::RVec
1272\brief A "std::vector"-like collection of values implementing handy operation to analyse them
1273\tparam T The type of the contained objects
1274
1275A RVec is a container designed to make analysis of values' collections fast and easy.
1276Its storage is contiguous in memory and its interface is designed such to resemble to the one
1277of the stl vector. In addition the interface features methods and external functions to ease
1278the manipulation and analysis of the data in the RVec.
1279
1280\note ROOT::VecOps::RVec can also be spelled simply ROOT::RVec. Shorthand aliases such as ROOT::RVecI or ROOT::RVecD
1281are also available as template instantiations of RVec of fundamental types. The full list of available aliases:
1282- RVecB (`bool`)
1283- RVecC (`char`)
1284- RVecD (`double`)
1285- RVecF (`float`)
1286- RVecI (`int`)
1287- RVecL (`long`)
1288- RVecLL (`long long`)
1289- RVecU (`unsigned`)
1290- RVecUL (`unsigned long`)
1291- RVecULL (`unsigned long long`)
1292
1293\note RVec does not attempt to be exception safe. Exceptions thrown by element constructors during insertions, swaps or
1294other operations will be propagated potentially leaving the RVec object in an invalid state.
1295
1296\note RVec methods (e.g. `at` or `size`) follow the STL naming convention instead of the ROOT naming convention in order
1297to make RVec a drop-in replacement for `std::vector`.
1298
1299\htmlonly
1300<a href="https://doi.org/10.5281/zenodo.1253756"><img src="https://zenodo.org/badge/DOI/10.5281/zenodo.1253756.svg" alt="DOI"></a>
1301\endhtmlonly
1302
1303## Table of Contents
1304- [Example](\ref example)
1305- [Owning and adopting memory](\ref owningandadoptingmemory)
1306- [Sorting and manipulation of indices](\ref sorting)
1307- [Usage in combination with RDataFrame](\ref usagetdataframe)
1308- [Reference for the RVec class](\ref RVecdoxyref)
1309
1310Also see the [reference for RVec helper functions](https://root.cern/doc/master/namespaceROOT_1_1VecOps.html).
1311
1312\anchor example
1313## Example
1314Suppose to have an event featuring a collection of muons with a certain pseudorapidity,
1315momentum and charge, e.g.:
1316~~~{.cpp}
1317std::vector<short> mu_charge {1, 1, -1, -1, -1, 1, 1, -1};
1318std::vector<float> mu_pt {56, 45, 32, 24, 12, 8, 7, 6.2};
1319std::vector<float> mu_eta {3.1, -.2, -1.1, 1, 4.1, 1.6, 2.4, -.5};
1320~~~
1321Suppose you want to extract the transverse momenta of the muons satisfying certain
1322criteria, for example consider only negatively charged muons with a pseudorapidity
1323smaller or equal to 2 and with a transverse momentum greater than 10 GeV.
1324Such a selection would require, among the other things, the management of an explicit
1325loop, for example:
1326~~~{.cpp}
1327std::vector<float> goodMuons_pt;
1328const auto size = mu_charge.size();
1329for (size_t i=0; i < size; ++i) {
1330 if (mu_pt[i] > 10 && abs(mu_eta[i]) <= 2. && mu_charge[i] == -1) {
1331 goodMuons_pt.emplace_back(mu_pt[i]);
1332 }
1333}
1334~~~
1335These operations become straightforward with RVec - we just need to *write what
1336we mean*:
1337~~~{.cpp}
1338auto goodMuons_pt = mu_pt[ (mu_pt > 10.f && abs(mu_eta) <= 2.f && mu_charge == -1) ]
1339~~~
1340Now the clean collection of transverse momenta can be used within the rest of the data analysis, for
1341example to fill a histogram.
1342
1343\anchor owningandadoptingmemory
1344## Owning and adopting memory
1345RVec has contiguous memory associated to it. It can own it or simply adopt it. In the latter case,
1346it can be constructed with the address of the memory associated to it and its length. For example:
1347~~~{.cpp}
1348std::vector<int> myStlVec {1,2,3};
1349RVec<int> myRVec(myStlVec.data(), myStlVec.size());
1350~~~
1351In this case, the memory associated to myStlVec and myRVec is the same, myRVec simply "adopted it".
1352If any method which implies a re-allocation is called, e.g. *emplace_back* or *resize*, the adopted
1353memory is released and new one is allocated. The previous content is copied in the new memory and
1354preserved.
1355
1356\anchor sorting
1357## Sorting and manipulation of indices
1358
1359### Sorting
1360RVec complies to the STL interfaces when it comes to iterations. As a result, standard algorithms
1361can be used, for example sorting:
1362~~~{.cpp}
1363RVec<double> v{6., 4., 5.};
1364std::sort(v.begin(), v.end());
1365~~~
1366
1367For convenience, helpers are provided too:
1368~~~{.cpp}
1369auto sorted_v = Sort(v);
1370auto reversed_v = Reverse(v);
1371~~~
1372
1373### Manipulation of indices
1374
1375It is also possible to manipulated the RVecs acting on their indices. For example,
1376the following syntax
1377~~~{.cpp}
1378RVecD v0 {9., 7., 8.};
1379auto v1 = Take(v0, {1, 2, 0});
1380~~~
1381will yield a new RVec<double> the content of which is the first, second and zeroth element of
1382v0, i.e. `{7., 8., 9.}`.
1383
1384The `Argsort` and `StableArgsort` helper extracts the indices which order the content of a `RVec`.
1385For example, this snippet accomplishes in a more expressive way what we just achieved:
1386~~~{.cpp}
1387auto v1_indices = Argsort(v0); // The content of v1_indices is {1, 2, 0}.
1388v1 = Take(v0, v1_indices);
1389~~~
1390
1391The `Take` utility allows to extract portions of the `RVec`. The content to be *taken*
1392can be specified with an `RVec` of indices or an integer. If the integer is negative,
1393elements will be picked starting from the end of the container:
1394~~~{.cpp}
1395RVecF vf {1.f, 2.f, 3.f, 4.f};
1396auto vf_1 = Take(vf, {1, 3}); // The content is {2.f, 4.f}
1397auto vf_2 = Take(vf, 2); // The content is {1.f, 2.f}
1398auto vf_3 = Take(vf, -3); // The content is {2.f, 3.f, 4.f}
1399~~~
1400
1401\anchor usagetdataframe
1402## Usage in combination with RDataFrame
1403RDataFrame leverages internally RVecs. Suppose to have a dataset stored in a
1404TTree which holds these columns (here we choose C arrays to represent the
1405collections, they could be as well std::vector instances):
1406~~~{.bash}
1407 nPart "nPart/I" An integer representing the number of particles
1408 px "px[nPart]/D" The C array of the particles' x component of the momentum
1409 py "py[nPart]/D" The C array of the particles' y component of the momentum
1410 E "E[nPart]/D" The C array of the particles' Energy
1411~~~
1412Suppose you'd like to plot in a histogram the transverse momenta of all particles
1413for which the energy is greater than 200 MeV.
1414The code required would just be:
1415~~~{.cpp}
1416RDataFrame d("mytree", "myfile.root");
1417auto cutPt = [](RVecD &pxs, RVecD &pys, RVecD &Es) {
1418 auto all_pts = sqrt(pxs * pxs + pys * pys);
1419 auto good_pts = all_pts[Es > 200.];
1420 return good_pts;
1421 };
1422
1423auto hpt = d.Define("pt", cutPt, {"px", "py", "E"})
1424 .Histo1D("pt");
1425hpt->Draw();
1426~~~
1427And if you'd like to express your selection as a string:
1428~~~{.cpp}
1429RDataFrame d("mytree", "myfile.root");
1430auto hpt = d.Define("pt", "sqrt(pxs * pxs + pys * pys)[E>200]")
1431 .Histo1D("pt");
1432hpt->Draw();
1433~~~
1434\anchor RVecdoxyref
1435**/
1436// clang-format on
1437
1438template <typename T>
1439class R__CLING_PTRCHECK(off) RVec : public RVecN<T, Internal::VecOps::RVecInlineStorageSize<T>::value> {
1441public:
1446 using SuperClass::begin;
1447 using SuperClass::size;
1448
1449 RVec() {}
1450
1451 explicit RVec(size_t Size, const T &Value) : SuperClass(Size, Value) {}
1452
1453 explicit RVec(size_t Size) : SuperClass(Size) {}
1454
1455 template <typename ItTy,
1456 typename = typename std::enable_if<std::is_convertible<
1457 typename std::iterator_traits<ItTy>::iterator_category, std::input_iterator_tag>::value>::type>
1458 RVec(ItTy S, ItTy E) : SuperClass(S, E)
1459 {
1460 }
1461
1462 RVec(std::initializer_list<T> IL) : SuperClass(IL) {}
1463
1464 RVec(const RVec &RHS) : SuperClass(RHS) {}
1465
1466 RVec &operator=(const RVec &RHS)
1467 {
1469 return *this;
1470 }
1471
1472 RVec(RVec &&RHS) : SuperClass(std::move(RHS)) {}
1473
1475 {
1476 SuperClass::operator=(std::move(RHS));
1477 return *this;
1478 }
1479
1480 RVec(Detail::VecOps::RVecImpl<T> &&RHS) : SuperClass(std::move(RHS)) {}
1481
1482 template <unsigned N>
1483 RVec(RVecN<T, N> &&RHS) : SuperClass(std::move(RHS)) {}
1484
1485 template <unsigned N>
1486 RVec(const RVecN<T, N> &RHS) : SuperClass(RHS) {}
1487
1488 RVec(const std::vector<T> &RHS) : SuperClass(RHS) {}
1489
1490 RVec(T* p, size_t n) : SuperClass(p, n) {}
1491
1492 // conversion
1494 operator RVec<U>() const
1495 {
1496 return RVec<U>(this->begin(), this->end());
1497 }
1498
1499 using SuperClass::operator[];
1500
1502 RVec operator[](const RVec<V> &conds) const
1503 {
1504 return RVec(SuperClass::operator[](conds));
1505 }
1506
1507 using SuperClass::at;
1508
1509 friend bool ROOT::Detail::VecOps::IsSmall<T>(const RVec<T> &v);
1510
1511 friend bool ROOT::Detail::VecOps::IsAdopting<T>(const RVec<T> &v);
1512};
1513
1514template <typename T, unsigned N>
1515inline size_t CapacityInBytes(const RVecN<T, N> &X)
1516{
1517 return X.capacity_in_bytes();
1518}
1519
1520///@name RVec Unary Arithmetic Operators
1521///@{
1522
1523#define RVEC_UNARY_OPERATOR(OP) \
1524template <typename T> \
1525RVec<T> operator OP(const RVec<T> &v) \
1526{ \
1527 RVec<T> ret(v); \
1528 for (auto &x : ret) \
1529 x = OP x; \
1530return ret; \
1531} \
1532
1537#undef RVEC_UNARY_OPERATOR
1538
1539///@}
1540///@name RVec Binary Arithmetic Operators
1541///@{
1542
1543#define ERROR_MESSAGE(OP) \
1544 "Cannot call operator " #OP " on vectors of different sizes."
1545
1546#define RVEC_BINARY_OPERATOR(OP) \
1547template <typename T0, typename T1> \
1548auto operator OP(const RVec<T0> &v, const T1 &y) \
1549 -> RVec<decltype(v[0] OP y)> \
1550{ \
1551 RVec<decltype(v[0] OP y)> ret(v.size()); \
1552 auto op = [&y](const T0 &x) { return x OP y; }; \
1553 std::transform(v.begin(), v.end(), ret.begin(), op); \
1554 return ret; \
1555} \
1556 \
1557template <typename T0, typename T1> \
1558auto operator OP(const T0 &x, const RVec<T1> &v) \
1559 -> RVec<decltype(x OP v[0])> \
1560{ \
1561 RVec<decltype(x OP v[0])> ret(v.size()); \
1562 auto op = [&x](const T1 &y) { return x OP y; }; \
1563 std::transform(v.begin(), v.end(), ret.begin(), op); \
1564 return ret; \
1565} \
1566 \
1567template <typename T0, typename T1> \
1568auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
1569 -> RVec<decltype(v0[0] OP v1[0])> \
1570{ \
1571 if (v0.size() != v1.size()) \
1572 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1573 \
1574 RVec<decltype(v0[0] OP v1[0])> ret(v0.size()); \
1575 auto op = [](const T0 &x, const T1 &y) { return x OP y; }; \
1576 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
1577 return ret; \
1578} \
1579
1588#undef RVEC_BINARY_OPERATOR
1589
1590///@}
1591///@name RVec Assignment Arithmetic Operators
1592///@{
1593
1594#define RVEC_ASSIGNMENT_OPERATOR(OP) \
1595template <typename T0, typename T1> \
1596RVec<T0>& operator OP(RVec<T0> &v, const T1 &y) \
1597{ \
1598 auto op = [&y](T0 &x) { return x OP y; }; \
1599 std::transform(v.begin(), v.end(), v.begin(), op); \
1600 return v; \
1601} \
1602 \
1603template <typename T0, typename T1> \
1604RVec<T0>& operator OP(RVec<T0> &v0, const RVec<T1> &v1) \
1605{ \
1606 if (v0.size() != v1.size()) \
1607 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1608 \
1609 auto op = [](T0 &x, const T1 &y) { return x OP y; }; \
1610 std::transform(v0.begin(), v0.end(), v1.begin(), v0.begin(), op); \
1611 return v0; \
1612} \
1613
1624#undef RVEC_ASSIGNMENT_OPERATOR
1625
1626///@}
1627///@name RVec Comparison and Logical Operators
1628///@{
1629
1630#define RVEC_LOGICAL_OPERATOR(OP) \
1631template <typename T0, typename T1> \
1632auto operator OP(const RVec<T0> &v, const T1 &y) \
1633 -> RVec<int> /* avoid std::vector<bool> */ \
1634{ \
1635 RVec<int> ret(v.size()); \
1636 auto op = [y](const T0 &x) -> int { return x OP y; }; \
1637 std::transform(v.begin(), v.end(), ret.begin(), op); \
1638 return ret; \
1639} \
1640 \
1641template <typename T0, typename T1> \
1642auto operator OP(const T0 &x, const RVec<T1> &v) \
1643 -> RVec<int> /* avoid std::vector<bool> */ \
1644{ \
1645 RVec<int> ret(v.size()); \
1646 auto op = [x](const T1 &y) -> int { return x OP y; }; \
1647 std::transform(v.begin(), v.end(), ret.begin(), op); \
1648 return ret; \
1649} \
1650 \
1651template <typename T0, typename T1> \
1652auto operator OP(const RVec<T0> &v0, const RVec<T1> &v1) \
1653 -> RVec<int> /* avoid std::vector<bool> */ \
1654{ \
1655 if (v0.size() != v1.size()) \
1656 throw std::runtime_error(ERROR_MESSAGE(OP)); \
1657 \
1658 RVec<int> ret(v0.size()); \
1659 auto op = [](const T0 &x, const T1 &y) -> int { return x OP y; }; \
1660 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), op); \
1661 return ret; \
1662} \
1663
1672#undef RVEC_LOGICAL_OPERATOR
1673
1674///@}
1675///@name RVec Standard Mathematical Functions
1676///@{
1677
1678/// \cond
1679template <typename T> struct PromoteTypeImpl;
1680
1681template <> struct PromoteTypeImpl<float> { using Type = float; };
1682template <> struct PromoteTypeImpl<double> { using Type = double; };
1683template <> struct PromoteTypeImpl<long double> { using Type = long double; };
1684
1685template <typename T> struct PromoteTypeImpl { using Type = double; };
1686
1687template <typename T>
1688using PromoteType = typename PromoteTypeImpl<T>::Type;
1689
1690template <typename U, typename V>
1691using PromoteTypes = decltype(PromoteType<U>() + PromoteType<V>());
1692
1693/// \endcond
1694
1695#define RVEC_UNARY_FUNCTION(NAME, FUNC) \
1696 template <typename T> \
1697 RVec<PromoteType<T>> NAME(const RVec<T> &v) \
1698 { \
1699 RVec<PromoteType<T>> ret(v.size()); \
1700 auto f = [](const T &x) { return FUNC(x); }; \
1701 std::transform(v.begin(), v.end(), ret.begin(), f); \
1702 return ret; \
1703 }
1704
1705#define RVEC_BINARY_FUNCTION(NAME, FUNC) \
1706 template <typename T0, typename T1> \
1707 RVec<PromoteTypes<T0, T1>> NAME(const T0 &x, const RVec<T1> &v) \
1708 { \
1709 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
1710 auto f = [&x](const T1 &y) { return FUNC(x, y); }; \
1711 std::transform(v.begin(), v.end(), ret.begin(), f); \
1712 return ret; \
1713 } \
1714 \
1715 template <typename T0, typename T1> \
1716 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v, const T1 &y) \
1717 { \
1718 RVec<PromoteTypes<T0, T1>> ret(v.size()); \
1719 auto f = [&y](const T1 &x) { return FUNC(x, y); }; \
1720 std::transform(v.begin(), v.end(), ret.begin(), f); \
1721 return ret; \
1722 } \
1723 \
1724 template <typename T0, typename T1> \
1725 RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &v0, const RVec<T1> &v1) \
1726 { \
1727 if (v0.size() != v1.size()) \
1728 throw std::runtime_error(ERROR_MESSAGE(NAME)); \
1729 \
1730 RVec<PromoteTypes<T0, T1>> ret(v0.size()); \
1731 auto f = [](const T0 &x, const T1 &y) { return FUNC(x, y); }; \
1732 std::transform(v0.begin(), v0.end(), v1.begin(), ret.begin(), f); \
1733 return ret; \
1734 } \
1735
1736#define RVEC_STD_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, std::F)
1737#define RVEC_STD_BINARY_FUNCTION(F) RVEC_BINARY_FUNCTION(F, std::F)
1738
1743
1747
1752
1757
1765
1772
1779
1784#undef RVEC_STD_UNARY_FUNCTION
1785
1786///@}
1787///@name RVec Fast Mathematical Functions with Vdt
1788///@{
1789
1790#ifdef R__HAS_VDT
1791#define RVEC_VDT_UNARY_FUNCTION(F) RVEC_UNARY_FUNCTION(F, vdt::F)
1792
1793RVEC_VDT_UNARY_FUNCTION(fast_expf)
1794RVEC_VDT_UNARY_FUNCTION(fast_logf)
1795RVEC_VDT_UNARY_FUNCTION(fast_sinf)
1796RVEC_VDT_UNARY_FUNCTION(fast_cosf)
1797RVEC_VDT_UNARY_FUNCTION(fast_tanf)
1798RVEC_VDT_UNARY_FUNCTION(fast_asinf)
1799RVEC_VDT_UNARY_FUNCTION(fast_acosf)
1800RVEC_VDT_UNARY_FUNCTION(fast_atanf)
1801
1802RVEC_VDT_UNARY_FUNCTION(fast_exp)
1803RVEC_VDT_UNARY_FUNCTION(fast_log)
1804RVEC_VDT_UNARY_FUNCTION(fast_sin)
1805RVEC_VDT_UNARY_FUNCTION(fast_cos)
1806RVEC_VDT_UNARY_FUNCTION(fast_tan)
1807RVEC_VDT_UNARY_FUNCTION(fast_asin)
1808RVEC_VDT_UNARY_FUNCTION(fast_acos)
1809RVEC_VDT_UNARY_FUNCTION(fast_atan)
1810#undef RVEC_VDT_UNARY_FUNCTION
1811
1812#endif // R__HAS_VDT
1813
1814#undef RVEC_UNARY_FUNCTION
1815
1816///@}
1817
1818/// Inner product
1819///
1820/// Example code, at the ROOT prompt:
1821/// ~~~{.cpp}
1822/// using namespace ROOT::VecOps;
1823/// RVec<float> v1 {1., 2., 3.};
1824/// RVec<float> v2 {4., 5., 6.};
1825/// auto v1_dot_v2 = Dot(v1, v2);
1826/// v1_dot_v2
1827/// // (float) 32.0000f
1828/// ~~~
1829template <typename T, typename V>
1830auto Dot(const RVec<T> &v0, const RVec<V> &v1) -> decltype(v0[0] * v1[0])
1831{
1832 if (v0.size() != v1.size())
1833 throw std::runtime_error("Cannot compute inner product of vectors of different sizes");
1834 return std::inner_product(v0.begin(), v0.end(), v1.begin(), decltype(v0[0] * v1[0])(0));
1835}
1836
1837/// Sum elements of an RVec
1838///
1839/// Example code, at the ROOT prompt:
1840/// ~~~{.cpp}
1841/// using namespace ROOT::VecOps;
1842/// RVecF v {1.f, 2.f, 3.f};
1843/// auto v_sum = Sum(v);
1844/// v_sum
1845/// // (float) 6.f
1846/// auto v_sum_d = Sum(v, 0.);
1847/// v_sum_d
1848/// // (double) 6.0000000
1849/// ~~~
1850/// ~~~{.cpp}
1851/// using namespace ROOT::VecOps;
1852/// const ROOT::Math::PtEtaPhiMVector lv0 {15.5f, .3f, .1f, 105.65f},
1853/// lv1 {34.32f, 2.2f, 3.02f, 105.65f},
1854/// lv2 {12.95f, 1.32f, 2.2f, 105.65f};
1855/// RVec<ROOT::Math::PtEtaPhiMVector> v {lv0, lv1, lv2};
1856/// auto v_sum_lv = Sum(v, ROOT::Math::PtEtaPhiMVector());
1857/// v_sum_lv
1858/// // (ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double> > &) (30.8489,2.46534,2.58947,361.084)
1859/// ~~~
1860template <typename T, typename R = T>
1861R Sum(const RVec<T> &v, const R zero = R(0))
1862{
1863 return std::accumulate(v.begin(), v.end(), zero);
1864}
1865
1866/// Get the mean of the elements of an RVec
1867///
1868/// The return type is a double precision floating point number.
1869///
1870/// Example code, at the ROOT prompt:
1871/// ~~~{.cpp}
1872/// using namespace ROOT::VecOps;
1873/// RVecF v {1.f, 2.f, 4.f};
1874/// auto v_mean = Mean(v);
1875/// v_mean
1876/// // (double) 2.3333333
1877/// ~~~
1878template <typename T>
1879double Mean(const RVec<T> &v)
1880{
1881 if (v.empty()) return 0.;
1882 return double(Sum(v)) / v.size();
1883}
1884
1885/// Get the mean of the elements of an RVec with custom initial value
1886///
1887/// The return type will be deduced from the `zero` parameter
1888///
1889/// Example code, at the ROOT prompt:
1890/// ~~~{.cpp}
1891/// using namespace ROOT::VecOps;
1892/// RVecF v {1.f, 2.f, 4.f};
1893/// auto v_mean_f = Mean(v, 0.f);
1894/// v_mean_f
1895/// // (float) 2.33333f
1896/// auto v_mean_d = Mean(v, 0.);
1897/// v_mean_d
1898/// // (double) 2.3333333
1899/// ~~~
1900/// ~~~{.cpp}
1901/// using namespace ROOT::VecOps;
1902/// const ROOT::Math::PtEtaPhiMVector lv0 {15.5f, .3f, .1f, 105.65f},
1903/// lv1 {34.32f, 2.2f, 3.02f, 105.65f},
1904/// lv2 {12.95f, 1.32f, 2.2f, 105.65f};
1905/// RVec<ROOT::Math::PtEtaPhiMVector> v {lv0, lv1, lv2};
1906/// auto v_mean_lv = Mean(v, ROOT::Math::PtEtaPhiMVector());
1907/// v_mean_lv
1908/// // (ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double> > &) (10.283,2.46534,2.58947,120.361)
1909/// ~~~
1910template <typename T, typename R = T>
1911R Mean(const RVec<T> &v, const R zero)
1912{
1913 if (v.empty()) return zero;
1914 return Sum(v, zero) / v.size();
1915}
1916
1917/// Get the greatest element of an RVec
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_max = Max(v);
1924/// v_max
1925/// (float) 4.00000f
1926/// ~~~~
1927template <typename T>
1929{
1930 return *std::max_element(v.begin(), v.end());
1931}
1932
1933/// Get the smallest element of an RVec
1934///
1935/// Example code, at the ROOT prompt:
1936/// ~~~~{.cpp}
1937/// using namespace ROOT::VecOps;
1938/// RVecF v {1.f, 2.f, 4.f};
1939/// auto v_min = Min(v);
1940/// v_min
1941/// (float) 1.00000f
1942/// ~~~~
1943template <typename T>
1945{
1946 return *std::min_element(v.begin(), v.end());
1947}
1948
1949/// Get the index of the greatest element of an RVec
1950/// In case of multiple occurrences of the maximum values,
1951/// the index corresponding to the first occurrence is returned.
1952///
1953/// Example code, at the ROOT prompt:
1954/// ~~~~{.cpp}
1955/// using namespace ROOT::VecOps;
1956/// RVecF v {1.f, 2.f, 4.f};
1957/// auto v_argmax = ArgMax(v);
1958/// v_argmax
1959/// // (unsigned long) 2
1960/// ~~~~
1961template <typename T>
1962std::size_t ArgMax(const RVec<T> &v)
1963{
1964 return std::distance(v.begin(), std::max_element(v.begin(), v.end()));
1965}
1966
1967/// Get the index of the smallest element of an RVec
1968/// In case of multiple occurrences of the minimum values,
1969/// the index corresponding to the first occurrence is returned.
1970///
1971/// Example code, at the ROOT prompt:
1972/// ~~~~{.cpp}
1973/// using namespace ROOT::VecOps;
1974/// RVecF v {1.f, 2.f, 4.f};
1975/// auto v_argmin = ArgMin(v);
1976/// v_argmin
1977/// // (unsigned long) 0
1978/// ~~~~
1979template <typename T>
1980std::size_t ArgMin(const RVec<T> &v)
1981{
1982 return std::distance(v.begin(), std::min_element(v.begin(), v.end()));
1983}
1984
1985/// Get the variance of the elements of an RVec
1986///
1987/// The return type is a double precision floating point number.
1988/// Example code, at the ROOT prompt:
1989/// ~~~{.cpp}
1990/// using namespace ROOT::VecOps;
1991/// RVecF v {1.f, 2.f, 4.f};
1992/// auto v_var = Var(v);
1993/// v_var
1994/// // (double) 2.3333333
1995/// ~~~
1996template <typename T>
1997double Var(const RVec<T> &v)
1998{
1999 const std::size_t size = v.size();
2000 if (size < std::size_t(2)) return 0.;
2001 T sum_squares(0), squared_sum(0);
2002 auto pred = [&sum_squares, &squared_sum](const T& x) {sum_squares+=x*x; squared_sum+=x;};
2003 std::for_each(v.begin(), v.end(), pred);
2004 squared_sum *= squared_sum;
2005 const auto dsize = (double) size;
2006 return 1. / (dsize - 1.) * (sum_squares - squared_sum / dsize );
2007}
2008
2009/// Get the standard deviation of the elements of an RVec
2010///
2011/// The return type is a double precision floating point number.
2012/// Example code, at the ROOT prompt:
2013/// ~~~{.cpp}
2014/// using namespace ROOT::VecOps;
2015/// RVecF v {1.f, 2.f, 4.f};
2016/// auto v_sd = StdDev(v);
2017/// v_sd
2018/// // (double) 1.5275252
2019/// ~~~
2020template <typename T>
2021double StdDev(const RVec<T> &v)
2022{
2023 return std::sqrt(Var(v));
2024}
2025
2026/// Create new collection applying a callable to the elements of the input collection
2027///
2028/// Example code, at the ROOT prompt:
2029/// ~~~{.cpp}
2030/// using namespace ROOT::VecOps;
2031/// RVecF v {1.f, 2.f, 4.f};
2032/// auto v_square = Map(v, [](float f){return f* 2.f;});
2033/// v_square
2034/// // (ROOT::VecOps::RVec<float> &) { 2.00000f, 4.00000f, 8.00000f }
2035///
2036/// RVecF x({1.f, 2.f, 3.f});
2037/// RVecF y({4.f, 5.f, 6.f});
2038/// RVecF z({7.f, 8.f, 9.f});
2039/// auto mod = [](float x, float y, float z) { return sqrt(x * x + y * y + z * z); };
2040/// auto v_mod = Map(x, y, z, mod);
2041/// v_mod
2042/// // (ROOT::VecOps::RVec<float> &) { 8.12404f, 9.64365f, 11.2250f }
2043/// ~~~
2044template <typename... Args>
2045auto Map(Args &&... args)
2046{
2047 /*
2048 Here the strategy in order to generalise the previous implementation of Map, i.e.
2049 `RVec Map(RVec, F)`, here we need to move the last parameter of the pack in first
2050 position in order to be able to invoke the Map function with automatic type deduction.
2051 This is achieved in two steps:
2052 1. Forward as tuple the pack to MapFromTuple
2053 2. Invoke the MapImpl helper which has the signature `template<...T, F> RVec MapImpl(F &&f, RVec<T>...)`
2054 */
2055 return ROOT::Internal::VecOps::MapFromTuple(std::forward_as_tuple(args...),
2056 std::make_index_sequence<sizeof...(args) - 1>());
2057}
2058
2059/// Create a new collection with the elements passing the filter expressed by the predicate
2060///
2061/// Example code, at the ROOT prompt:
2062/// ~~~{.cpp}
2063/// using namespace ROOT::VecOps;
2064/// RVecI v {1, 2, 4};
2065/// auto v_even = Filter(v, [](int i){return 0 == i%2;});
2066/// v_even
2067/// // (ROOT::VecOps::RVec<int> &) { 2, 4 }
2068/// ~~~
2069template <typename T, typename F>
2071{
2072 const auto thisSize = v.size();
2073 RVec<T> w;
2074 w.reserve(thisSize);
2075 for (auto &&val : v) {
2076 if (f(val))
2077 w.emplace_back(val);
2078 }
2079 return w;
2080}
2081
2082/// Return true if any of the elements equates to true, return false otherwise.
2083///
2084/// Example code, at the ROOT prompt:
2085/// ~~~{.cpp}
2086/// using namespace ROOT::VecOps;
2087/// RVecI v {0, 1, 0};
2088/// auto anyTrue = Any(v);
2089/// anyTrue
2090/// // (bool) true
2091/// ~~~
2092template <typename T>
2093auto Any(const RVec<T> &v) -> decltype(v[0] == true)
2094{
2095 for (auto &&e : v)
2096 if (static_cast<bool>(e) == true)
2097 return true;
2098 return false;
2099}
2100
2101/// Return true if all of the elements equate to true, return false otherwise.
2102///
2103/// Example code, at the ROOT prompt:
2104/// ~~~{.cpp}
2105/// using namespace ROOT::VecOps;
2106/// RVecI v {0, 1, 0};
2107/// auto allTrue = All(v);
2108/// allTrue
2109/// // (bool) false
2110/// ~~~
2111template <typename T>
2112auto All(const RVec<T> &v) -> decltype(v[0] == false)
2113{
2114 for (auto &&e : v)
2115 if (static_cast<bool>(e) == false)
2116 return false;
2117 return true;
2118}
2119
2120template <typename T>
2121void swap(RVec<T> &lhs, RVec<T> &rhs)
2122{
2123 lhs.swap(rhs);
2124}
2125
2126/// Return an RVec of indices that sort the input RVec
2127///
2128/// Example code, at the ROOT prompt:
2129/// ~~~{.cpp}
2130/// using namespace ROOT::VecOps;
2131/// RVecD v {2., 3., 1.};
2132/// auto sortIndices = Argsort(v)
2133/// // (ROOT::VecOps::RVec<unsigned long> &) { 2, 0, 1 }
2134/// auto values = Take(v, sortIndices)
2135/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 3.0000000 }
2136/// ~~~
2137template <typename T>
2139{
2140 using size_type = typename RVec<T>::size_type;
2141 RVec<size_type> i(v.size());
2142 std::iota(i.begin(), i.end(), 0);
2143 std::sort(i.begin(), i.end(), [&v](size_type i1, size_type i2) { return v[i1] < v[i2]; });
2144 return i;
2145}
2146
2147/// Return an RVec of indices that sort the input RVec based on a comparison function.
2148///
2149/// Example code, at the ROOT prompt:
2150/// ~~~{.cpp}
2151/// using namespace ROOT::VecOps;
2152/// RVecD v {2., 3., 1.};
2153/// auto sortIndices = Argsort(v, [](double x, double y) {return x > y;})
2154/// // (ROOT::VecOps::RVec<unsigned long> &) { 1, 0, 2 }
2155/// auto values = Take(v, sortIndices)
2156/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 1.0000000 }
2157/// ~~~
2158template <typename T, typename Compare>
2160{
2161 using size_type = typename RVec<T>::size_type;
2162 RVec<size_type> i(v.size());
2163 std::iota(i.begin(), i.end(), 0);
2164 std::sort(i.begin(), i.end(),
2165 [&v, &c](size_type i1, size_type i2) { return c(v[i1], v[i2]); });
2166 return i;
2167}
2168
2169/// Return an RVec of indices that sort the input RVec
2170/// while keeping the order of equal elements.
2171/// This is the stable variant of `Argsort`.
2172///
2173/// Example code, at the ROOT prompt:
2174/// ~~~{.cpp}
2175/// using namespace ROOT::VecOps;
2176/// RVecD v {2., 3., 2., 1.};
2177/// auto sortIndices = StableArgsort(v)
2178/// // (ROOT::VecOps::RVec<unsigned long> &) { 3, 0, 2, 1 }
2179/// auto values = Take(v, sortIndices)
2180/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 2.0000000, 3.0000000 }
2181/// ~~~
2182template <typename T>
2184{
2185 using size_type = typename RVec<T>::size_type;
2186 RVec<size_type> i(v.size());
2187 std::iota(i.begin(), i.end(), 0);
2188 std::stable_sort(i.begin(), i.end(), [&v](size_type i1, size_type i2) { return v[i1] < v[i2]; });
2189 return i;
2190}
2191
2192/// Return an RVec of indices that sort the input RVec based on a comparison function
2193/// while keeping the order of equal elements.
2194/// This is the stable variant of `Argsort`.
2195///
2196/// Example code, at the ROOT prompt:
2197/// ~~~{.cpp}
2198/// using namespace ROOT::VecOps;
2199/// RVecD v {2., 3., 2., 1.};
2200/// auto sortIndices = StableArgsort(v, [](double x, double y) {return x > y;})
2201/// // (ROOT::VecOps::RVec<unsigned long> &) { 1, 0, 2, 3 }
2202/// auto values = Take(v, sortIndices)
2203/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 2.0000000, 1.0000000 }
2204/// ~~~
2205template <typename T, typename Compare>
2207{
2208 using size_type = typename RVec<T>::size_type;
2209 RVec<size_type> i(v.size());
2210 std::iota(i.begin(), i.end(), 0);
2211 std::stable_sort(i.begin(), i.end(), [&v, &c](size_type i1, size_type i2) { return c(v[i1], v[i2]); });
2212 return i;
2213}
2214
2215/// Return elements of a vector at given indices
2216///
2217/// Example code, at the ROOT prompt:
2218/// ~~~{.cpp}
2219/// using namespace ROOT::VecOps;
2220/// RVecD v {2., 3., 1.};
2221/// auto vTaken = Take(v, {0,2});
2222/// vTaken
2223/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 1.0000000 }
2224/// ~~~
2225template <typename T>
2226RVec<T> Take(const RVec<T> &v, const RVec<typename RVec<T>::size_type> &i)
2227{
2228 using size_type = typename RVec<T>::size_type;
2229 const size_type isize = i.size();
2230 RVec<T> r(isize);
2231 for (size_type k = 0; k < isize; k++)
2232 r[k] = v[i[k]];
2233 return r;
2234}
2235
2236/// Return first or last `n` elements of an RVec
2237///
2238/// if `n > 0` and last elements if `n < 0`.
2239///
2240/// Example code, at the ROOT prompt:
2241/// ~~~{.cpp}
2242/// using namespace ROOT::VecOps;
2243/// RVecD v {2., 3., 1.};
2244/// auto firstTwo = Take(v, 2);
2245/// firstTwo
2246/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 3.0000000 }
2247/// auto lastOne = Take(v, -1);
2248/// lastOne
2249/// // (ROOT::VecOps::RVec<double>) { 1.0000000 }
2250/// ~~~
2251template <typename T>
2252RVec<T> Take(const RVec<T> &v, const int n)
2253{
2254 using size_type = typename RVec<T>::size_type;
2255 const size_type size = v.size();
2256 const size_type absn = std::abs(n);
2257 if (absn > size) {
2258 std::stringstream ss;
2259 ss << "Try to take " << absn << " elements but vector has only size " << size << ".";
2260 throw std::runtime_error(ss.str());
2261 }
2262 RVec<T> r(absn);
2263 if (n < 0) {
2264 for (size_type k = 0; k < absn; k++)
2265 r[k] = v[size - absn + k];
2266 } else {
2267 for (size_type k = 0; k < absn; k++)
2268 r[k] = v[k];
2269 }
2270 return r;
2271}
2272
2273/// Return a copy of the container without the elements at the specified indices.
2274///
2275/// Duplicated and out-of-range indices in idxs are ignored.
2276template <typename T>
2278{
2279 // clean up input indices
2280 std::sort(idxs.begin(), idxs.end());
2281 idxs.erase(std::unique(idxs.begin(), idxs.end()), idxs.end());
2282
2283 RVec<T> r;
2284 if (v.size() > idxs.size())
2285 r.reserve(v.size() - idxs.size());
2286
2287 auto discardIt = idxs.begin();
2288 using sz_t = typename RVec<T>::size_type;
2289 for (sz_t i = 0u; i < v.size(); ++i) {
2290 if (discardIt != idxs.end() && i == *discardIt)
2291 ++discardIt;
2292 else
2293 r.emplace_back(v[i]);
2294 }
2295
2296 return r;
2297}
2298
2299/// Return copy of reversed vector
2300///
2301/// Example code, at the ROOT prompt:
2302/// ~~~{.cpp}
2303/// using namespace ROOT::VecOps;
2304/// RVecD v {2., 3., 1.};
2305/// auto v_reverse = Reverse(v);
2306/// v_reverse
2307/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 3.0000000, 2.0000000 }
2308/// ~~~
2309template <typename T>
2311{
2312 RVec<T> r(v);
2313 std::reverse(r.begin(), r.end());
2314 return r;
2315}
2316
2317/// Return copy of RVec with elements sorted in ascending order
2318///
2319/// This helper is different from Argsort since it does not return an RVec of indices,
2320/// but an RVec of values.
2321///
2322/// Example code, at the ROOT prompt:
2323/// ~~~{.cpp}
2324/// using namespace ROOT::VecOps;
2325/// RVecD v {2., 3., 1.};
2326/// auto v_sorted = Sort(v);
2327/// v_sorted
2328/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 3.0000000 }
2329/// ~~~
2330template <typename T>
2332{
2333 RVec<T> r(v);
2334 std::sort(r.begin(), r.end());
2335 return r;
2336}
2337
2338/// Return copy of RVec with elements sorted based on a comparison operator
2339///
2340/// The comparison operator has to fullfill the same requirements of the
2341/// predicate of by std::sort.
2342///
2343///
2344/// This helper is different from Argsort since it does not return an RVec of indices,
2345/// but an RVec of values.
2346///
2347/// Example code, at the ROOT prompt:
2348/// ~~~{.cpp}
2349/// using namespace ROOT::VecOps;
2350/// RVecD v {2., 3., 1.};
2351/// auto v_sorted = Sort(v, [](double x, double y) {return 1/x < 1/y;});
2352/// v_sorted
2353/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 1.0000000 }
2354/// ~~~
2355template <typename T, typename Compare>
2357{
2358 RVec<T> r(v);
2359 std::sort(r.begin(), r.end(), std::forward<Compare>(c));
2360 return r;
2361}
2362
2363/// Return copy of RVec with elements sorted in ascending order
2364/// while keeping the order of equal elements.
2365///
2366/// This is the stable variant of `Sort`.
2367///
2368/// This helper is different from StableArgsort since it does not return an RVec of indices,
2369/// but an RVec of values.
2370///
2371/// Example code, at the ROOT prompt:
2372/// ~~~{.cpp}
2373/// using namespace ROOT::VecOps;
2374/// RVecD v {2., 3., 2, 1.};
2375/// auto v_sorted = StableSort(v);
2376/// v_sorted
2377/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000, 2.0000000, 3.0000000 }
2378/// ~~~
2379template <typename T>
2381{
2382 RVec<T> r(v);
2383 std::stable_sort(r.begin(), r.end());
2384 return r;
2385}
2386
2387// clang-format off
2388/// Return copy of RVec with elements sorted based on a comparison operator
2389/// while keeping the order of equal elements.
2390///
2391/// The comparison operator has to fullfill the same requirements of the
2392/// predicate of std::stable_sort.
2393///
2394/// This helper is different from StableArgsort since it does not return an RVec of indices,
2395/// but an RVec of values.
2396///
2397/// This is the stable variant of `Sort`.
2398///
2399/// Example code, at the ROOT prompt:
2400/// ~~~{.cpp}
2401/// using namespace ROOT::VecOps;
2402/// RVecD v {2., 3., 2., 1.};
2403/// auto v_sorted = StableSort(v, [](double x, double y) {return 1/x < 1/y;});
2404/// v_sorted
2405/// // (ROOT::VecOps::RVec<double> &) { 3.0000000, 2.0000000, 2.0000000, 1.0000000 }
2406/// ~~~
2407/// ~~~{.cpp}
2408/// using namespace ROOT::VecOps;
2409/// RVec<RVecD> v {{2., 4.}, {3., 1.}, {2, 1.}, {1., 4.}};
2410/// 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];});
2411/// v_sorted
2412/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<double> > &) { { 1.0000000, 4.0000000 }, { 2.0000000, 1.0000000 }, { 2.0000000, 4.0000000 }, { 3.0000000, 1.0000000 } }
2413/// ~~~
2414// clang-format off
2415template <typename T, typename Compare>
2417{
2418 RVec<T> r(v);
2419 std::stable_sort(r.begin(), r.end(), std::forward<Compare>(c));
2420 return r;
2421}
2422
2423/// Return the indices that represent all combinations of the elements of two
2424/// RVecs.
2425///
2426/// The type of the return value is an RVec of two RVecs containing indices.
2427///
2428/// Example code, at the ROOT prompt:
2429/// ~~~{.cpp}
2430/// using namespace ROOT::VecOps;
2431/// auto comb_idx = Combinations(3, 2);
2432/// comb_idx
2433/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 1, 1, 2, 2 }, { 0, 1, 0, 1, 0, 1 } }
2434/// ~~~
2435inline RVec<RVec<std::size_t>> Combinations(const std::size_t size1, const std::size_t size2)
2436{
2437 using size_type = std::size_t;
2439 r[0].resize(size1*size2);
2440 r[1].resize(size1*size2);
2441 size_type c = 0;
2442 for(size_type i=0; i<size1; i++) {
2443 for(size_type j=0; j<size2; j++) {
2444 r[0][c] = i;
2445 r[1][c] = j;
2446 c++;
2447 }
2448 }
2449 return r;
2450}
2451
2452/// Return the indices that represent all combinations of the elements of two
2453/// RVecs.
2454///
2455/// The type of the return value is an RVec of two RVecs containing indices.
2456///
2457/// Example code, at the ROOT prompt:
2458/// ~~~{.cpp}
2459/// using namespace ROOT::VecOps;
2460/// RVecD v1 {1., 2., 3.};
2461/// RVecD v2 {-4., -5.};
2462/// auto comb_idx = Combinations(v1, v2);
2463/// comb_idx
2464/// // (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 1, 1, 2, 2 }, { 0, 1, 0, 1, 0, 1 } }
2465/// ~~~
2466template <typename T1, typename T2>
2468{
2469 return Combinations(v1.size(), v2.size());
2470}
2471
2472/// Return the indices that represent all unique combinations of the
2473/// elements of a given RVec.
2474///
2475/// ~~~{.cpp}
2476/// using namespace ROOT::VecOps;
2477/// RVecD v {1., 2., 3., 4.};
2478/// auto v_1 = Combinations(v, 1);
2479/// v_1
2480/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 1, 2, 3 } }
2481/// auto v_2 = Combinations(v, 2);
2482/// v_2
2483/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 0, 1, 1, 2 }, { 1, 2, 3, 2, 3, 3 } }
2484/// auto v_3 = Combinations(v, 3);
2485/// v_3
2486/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0, 0, 0, 1 }, { 1, 1, 2, 2 }, { 2, 3, 3, 3 } }
2487/// auto v_4 = Combinations(v, 4);
2488/// v_4
2489/// (ROOT::VecOps::RVec<ROOT::VecOps::RVec<unsigned long> > &) { { 0 }, { 1 }, { 2 }, { 3 } }
2490/// ~~~
2491template <typename T>
2493{
2494 using size_type = typename RVec<T>::size_type;
2495 const size_type s = v.size();
2496 if (n > s) {
2497 throw std::runtime_error("Cannot make unique combinations of size " + std::to_string(n) +
2498 " from vector of size " + std::to_string(s) + ".");
2499 }
2500
2502 for(size_type k=0; k<s; k++)
2503 indices[k] = k;
2504
2505 const auto innersize = [=] {
2506 size_type inners = s - n + 1;
2507 for (size_type m = s - n + 2; m <= s; ++m)
2508 inners *= m;
2509
2510 size_type factn = 1;
2511 for (size_type i = 2; i <= n; ++i)
2512 factn *= i;
2513 inners /= factn;
2514
2515 return inners;
2516 }();
2517
2519 size_type inneridx = 0;
2520 for (size_type k = 0; k < n; k++)
2521 c[k][inneridx] = indices[k];
2522 ++inneridx;
2523
2524 while (true) {
2525 bool run_through = true;
2526 long i = n - 1;
2527 for (; i>=0; i--) {
2528 if (indices[i] != i + s - n){
2529 run_through = false;
2530 break;
2531 }
2532 }
2533 if (run_through) {
2534 return c;
2535 }
2536 indices[i]++;
2537 for (long j=i+1; j<(long)n; j++)
2538 indices[j] = indices[j-1] + 1;
2539 for (size_type k = 0; k < n; k++)
2540 c[k][inneridx] = indices[k];
2541 ++inneridx;
2542 }
2543}
2544
2545/// Return the indices of the elements which are not zero
2546///
2547/// Example code, at the ROOT prompt:
2548/// ~~~{.cpp}
2549/// using namespace ROOT::VecOps;
2550/// RVecD v {2., 0., 3., 0., 1.};
2551/// auto nonzero_idx = Nonzero(v);
2552/// nonzero_idx
2553/// // (ROOT::VecOps::RVec<unsigned long> &) { 0, 2, 4 }
2554/// ~~~
2555template <typename T>
2557{
2558 using size_type = typename RVec<T>::size_type;
2560 const auto size = v.size();
2561 r.reserve(size);
2562 for(size_type i=0; i<size; i++) {
2563 if(v[i] != 0) {
2564 r.emplace_back(i);
2565 }
2566 }
2567 return r;
2568}
2569
2570/// Return the intersection of elements of two RVecs.
2571///
2572/// Each element of v1 is looked up in v2 and added to the returned vector if
2573/// found. Following, the order of v1 is preserved. If v2 is already sorted, the
2574/// optional argument v2_is_sorted can be used to toggle of the internal sorting
2575/// step, therewith optimising runtime.
2576///
2577/// Example code, at the ROOT prompt:
2578/// ~~~{.cpp}
2579/// using namespace ROOT::VecOps;
2580/// RVecD v1 {1., 2., 3.};
2581/// RVecD v2 {-4., -5., 2., 1.};
2582/// auto v1_intersect_v2 = Intersect(v1, v2);
2583/// v1_intersect_v2
2584/// // (ROOT::VecOps::RVec<double> &) { 1.0000000, 2.0000000 }
2585/// ~~~
2586template <typename T>
2587RVec<T> Intersect(const RVec<T>& v1, const RVec<T>& v2, bool v2_is_sorted = false)
2588{
2589 RVec<T> v2_sorted;
2590 if (!v2_is_sorted) v2_sorted = Sort(v2);
2591 const auto v2_begin = v2_is_sorted ? v2.begin() : v2_sorted.begin();
2592 const auto v2_end = v2_is_sorted ? v2.end() : v2_sorted.end();
2593 RVec<T> r;
2594 const auto size = v1.size();
2595 r.reserve(size);
2596 using size_type = typename RVec<T>::size_type;
2597 for(size_type i=0; i<size; i++) {
2598 if (std::binary_search(v2_begin, v2_end, v1[i])) {
2599 r.emplace_back(v1[i]);
2600 }
2601 }
2602 return r;
2603}
2604
2605/// Return the elements of v1 if the condition c is true and v2 if the
2606/// condition c is false.
2607///
2608/// Example code, at the ROOT prompt:
2609/// ~~~{.cpp}
2610/// using namespace ROOT::VecOps;
2611/// RVecD v1 {1., 2., 3.};
2612/// RVecD v2 {-1., -2., -3.};
2613/// auto c = v1 > 1;
2614/// c
2615/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2616/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2617/// if_c_v1_else_v2
2618/// // (ROOT::VecOps::RVec<double> &) { -1.0000000, 2.0000000, 3.0000000 }
2619/// ~~~
2620template <typename T>
2621RVec<T> Where(const RVec<int>& c, const RVec<T>& v1, const RVec<T>& v2)
2622{
2623 using size_type = typename RVec<T>::size_type;
2624 const size_type size = c.size();
2625 RVec<T> r;
2626 r.reserve(size);
2627 for (size_type i=0; i<size; i++) {
2628 r.emplace_back(c[i] != 0 ? v1[i] : v2[i]);
2629 }
2630 return r;
2631}
2632
2633/// Return the elements of v1 if the condition c is true and sets the value v2
2634/// if the condition c is false.
2635///
2636/// Example code, at the ROOT prompt:
2637/// ~~~{.cpp}
2638/// using namespace ROOT::VecOps;
2639/// RVecD v1 {1., 2., 3.};
2640/// double v2 = 4.;
2641/// auto c = v1 > 1;
2642/// c
2643/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2644/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2645/// if_c_v1_else_v2
2646/// // (ROOT::VecOps::RVec<double>) { 4.0000000, 2.0000000, 3.0000000 }
2647/// ~~~
2648template <typename T>
2650{
2651 using size_type = typename RVec<T>::size_type;
2652 const size_type size = c.size();
2653 RVec<T> r;
2654 r.reserve(size);
2655 for (size_type i=0; i<size; i++) {
2656 r.emplace_back(c[i] != 0 ? v1[i] : v2);
2657 }
2658 return r;
2659}
2660
2661/// Return the elements of v2 if the condition c is false and sets the value v1
2662/// if the condition c is true.
2663///
2664/// Example code, at the ROOT prompt:
2665/// ~~~{.cpp}
2666/// using namespace ROOT::VecOps;
2667/// double v1 = 4.;
2668/// RVecD v2 {1., 2., 3.};
2669/// auto c = v2 > 1;
2670/// c
2671/// // (ROOT::VecOps::RVec<int> &) { 0, 1, 1 }
2672/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2673/// if_c_v1_else_v2
2674/// // (ROOT::VecOps::RVec<double>) { 1.0000000, 4.0000000, 4.0000000 }
2675/// ~~~
2676template <typename T>
2678{
2679 using size_type = typename RVec<T>::size_type;
2680 const size_type size = c.size();
2681 RVec<T> r;
2682 r.reserve(size);
2683 for (size_type i=0; i<size; i++) {
2684 r.emplace_back(c[i] != 0 ? v1 : v2[i]);
2685 }
2686 return r;
2687}
2688
2689/// Return a vector with the value v2 if the condition c is false and sets the
2690/// value v1 if the condition c is true.
2691///
2692/// Example code, at the ROOT prompt:
2693/// ~~~{.cpp}
2694/// using namespace ROOT::VecOps;
2695/// double v1 = 4.;
2696/// double v2 = 2.;
2697/// RVecI c {0, 1, 1};
2698/// auto if_c_v1_else_v2 = Where(c, v1, v2);
2699/// if_c_v1_else_v2
2700/// // (ROOT::VecOps::RVec<double>) { 2.0000000, 4.0000000, 4.0000000 }
2701/// ~~~
2702template <typename T>
2704{
2705 using size_type = typename RVec<T>::size_type;
2706 const size_type size = c.size();
2707 RVec<T> r;
2708 r.reserve(size);
2709 for (size_type i=0; i<size; i++) {
2710 r.emplace_back(c[i] != 0 ? v1 : v2);
2711 }
2712 return r;
2713}
2714
2715/// Return the concatenation of two RVecs.
2716///
2717/// Example code, at the ROOT prompt:
2718/// ~~~{.cpp}
2719/// using namespace ROOT::VecOps;
2720/// RVecF rvf {0.f, 1.f, 2.f};
2721/// RVecI rvi {7, 8, 9};
2722/// Concatenate(rvf, rvi)
2723/// // (ROOT::VecOps::RVec<float>) { 0.00000f, 1.00000f, 2.00000f, 7.00000f, 8.00000f, 9.00000f }
2724/// ~~~
2727{
2728 RVec<Common_t> res;
2729 res.reserve(v0.size() + v1.size());
2730 std::copy(v0.begin(), v0.end(), std::back_inserter(res));
2731 std::copy(v1.begin(), v1.end(), std::back_inserter(res));
2732 return res;
2733}
2734
2735/// Return the angle difference \f$\Delta \phi\f$ of two scalars.
2736///
2737/// The function computes the closest angle from v1 to v2 with sign and is
2738/// therefore in the range \f$[-\pi, \pi]\f$.
2739/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2740/// to degrees \f$c = 180\f$.
2741template <typename T>
2742T DeltaPhi(T v1, T v2, const T c = M_PI)
2743{
2745 "DeltaPhi must be called with floating point values.");
2746 auto r = std::fmod(v2 - v1, 2.0 * c);
2747 if (r < -c) {
2748 r += 2.0 * c;
2749 }
2750 else if (r > c) {
2751 r -= 2.0 * c;
2752 }
2753 return r;
2754}
2755
2756/// Return the angle difference \f$\Delta \phi\f$ in radians of two vectors.
2757///
2758/// The function computes the closest angle from v1 to v2 with sign and is
2759/// therefore in the range \f$[-\pi, \pi]\f$.
2760/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2761/// to degrees \f$c = 180\f$.
2762template <typename T>
2763RVec<T> DeltaPhi(const RVec<T>& v1, const RVec<T>& v2, const T c = M_PI)
2764{
2765 using size_type = typename RVec<T>::size_type;
2766 const size_type size = v1.size();
2767 auto r = RVec<T>(size);
2768 for (size_type i = 0; i < size; i++) {
2769 r[i] = DeltaPhi(v1[i], v2[i], c);
2770 }
2771 return r;
2772}
2773
2774/// Return the angle difference \f$\Delta \phi\f$ in radians of a vector and a scalar.
2775///
2776/// The function computes the closest angle from v1 to v2 with sign and is
2777/// therefore in the range \f$[-\pi, \pi]\f$.
2778/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2779/// to degrees \f$c = 180\f$.
2780template <typename T>
2781RVec<T> DeltaPhi(const RVec<T>& v1, T v2, const T c = M_PI)
2782{
2783 using size_type = typename RVec<T>::size_type;
2784 const size_type size = v1.size();
2785 auto r = RVec<T>(size);
2786 for (size_type i = 0; i < size; i++) {
2787 r[i] = DeltaPhi(v1[i], v2, c);
2788 }
2789 return r;
2790}
2791
2792/// Return the angle difference \f$\Delta \phi\f$ in radians of a scalar and a vector.
2793///
2794/// The function computes the closest angle from v1 to v2 with sign and is
2795/// therefore in the range \f$[-\pi, \pi]\f$.
2796/// The computation is done per default in radians \f$c = \pi\f$ but can be switched
2797/// to degrees \f$c = 180\f$.
2798template <typename T>
2799RVec<T> DeltaPhi(T v1, const RVec<T>& v2, const T c = M_PI)
2800{
2801 using size_type = typename RVec<T>::size_type;
2802 const size_type size = v2.size();
2803 auto r = RVec<T>(size);
2804 for (size_type i = 0; i < size; i++) {
2805 r[i] = DeltaPhi(v1, v2[i], c);
2806 }
2807 return r;
2808}
2809
2810/// Return the square of the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2811/// the collections eta1, eta2, phi1 and phi2.
2812///
2813/// The function computes \f$\Delta R^2 = (\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2\f$
2814/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2815/// be set to radian or degrees using the optional argument c, see the documentation
2816/// of the DeltaPhi helper.
2817template <typename T>
2818RVec<T> DeltaR2(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
2819{
2820 const auto dphi = DeltaPhi(phi1, phi2, c);
2821 return (eta1 - eta2) * (eta1 - eta2) + dphi * dphi;
2822}
2823
2824/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2825/// the collections eta1, eta2, phi1 and phi2.
2826///
2827/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
2828/// of the given collections eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2829/// be set to radian or degrees using the optional argument c, see the documentation
2830/// of the DeltaPhi helper.
2831template <typename T>
2832RVec<T> DeltaR(const RVec<T>& eta1, const RVec<T>& eta2, const RVec<T>& phi1, const RVec<T>& phi2, const T c = M_PI)
2833{
2834 return sqrt(DeltaR2(eta1, eta2, phi1, phi2, c));
2835}
2836
2837/// Return the distance on the \f$\eta\f$-\f$\phi\f$ plane (\f$\Delta R\f$) from
2838/// the scalars eta1, eta2, phi1 and phi2.
2839///
2840/// The function computes \f$\Delta R = \sqrt{(\eta_1 - \eta_2)^2 + (\phi_1 - \phi_2)^2}\f$
2841/// of the given scalars eta1, eta2, phi1 and phi2. The angle \f$\phi\f$ can
2842/// be set to radian or degrees using the optional argument c, see the documentation
2843/// of the DeltaPhi helper.
2844template <typename T>
2845T DeltaR(T eta1, T eta2, T phi1, T phi2, const T c = M_PI)
2846{
2847 const auto dphi = DeltaPhi(phi1, phi2, c);
2848 return std::sqrt((eta1 - eta2) * (eta1 - eta2) + dphi * dphi);
2849}
2850
2851/// Return the invariant mass of two particles given the collections of the quantities
2852/// transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
2853///
2854/// The function computes the invariant mass of two particles with the four-vectors
2855/// (pt1, eta2, phi1, mass1) and (pt2, eta2, phi2, mass2).
2856template <typename T>
2858 const RVec<T>& pt1, const RVec<T>& eta1, const RVec<T>& phi1, const RVec<T>& mass1,
2859 const RVec<T>& pt2, const RVec<T>& eta2, const RVec<T>& phi2, const RVec<T>& mass2)
2860{
2861 std::size_t size = pt1.size();
2862
2863 R__ASSERT(eta1.size() == size && phi1.size() == size && mass1.size() == size);
2864 R__ASSERT(pt2.size() == size && phi2.size() == size && mass2.size() == size);
2865
2866 RVec<T> inv_masses(size);
2867
2868 for (std::size_t i = 0u; i < size; ++i) {
2869 // Conversion from (pt, eta, phi, mass) to (x, y, z, e) coordinate system
2870 const auto x1 = pt1[i] * std::cos(phi1[i]);
2871 const auto y1 = pt1[i] * std::sin(phi1[i]);
2872 const auto z1 = pt1[i] * std::sinh(eta1[i]);
2873 const auto e1 = std::sqrt(x1 * x1 + y1 * y1 + z1 * z1 + mass1[i] * mass1[i]);
2874
2875 const auto x2 = pt2[i] * std::cos(phi2[i]);
2876 const auto y2 = pt2[i] * std::sin(phi2[i]);
2877 const auto z2 = pt2[i] * std::sinh(eta2[i]);
2878 const auto e2 = std::sqrt(x2 * x2 + y2 * y2 + z2 * z2 + mass2[i] * mass2[i]);
2879
2880 // Addition of particle four-vector elements
2881 const auto e = e1 + e2;
2882 const auto x = x1 + x2;
2883 const auto y = y1 + y2;
2884 const auto z = z1 + z2;
2885
2886 inv_masses[i] = std::sqrt(e * e - x * x - y * y - z * z);
2887 }
2888
2889 // Return invariant mass with (+, -, -, -) metric
2890 return inv_masses;
2891}
2892
2893/// Return the invariant mass of multiple particles given the collections of the
2894/// quantities transverse momentum (pt), rapidity (eta), azimuth (phi) and mass.
2895///
2896/// The function computes the invariant mass of multiple particles with the
2897/// four-vectors (pt, eta, phi, mass).
2898template <typename T>
2899T InvariantMass(const RVec<T>& pt, const RVec<T>& eta, const RVec<T>& phi, const RVec<T>& mass)
2900{
2901 const std::size_t size = pt.size();
2902
2903 R__ASSERT(eta.size() == size && phi.size() == size && mass.size() == size);
2904
2905 T x_sum = 0.;
2906 T y_sum = 0.;
2907 T z_sum = 0.;
2908 T e_sum = 0.;
2909
2910 for (std::size_t i = 0u; i < size; ++ i) {
2911 // Convert to (e, x, y, z) coordinate system and update sums
2912 const auto x = pt[i] * std::cos(phi[i]);
2913 x_sum += x;
2914 const auto y = pt[i] * std::sin(phi[i]);
2915 y_sum += y;
2916 const auto z = pt[i] * std::sinh(eta[i]);
2917 z_sum += z;
2918 const auto e = std::sqrt(x * x + y * y + z * z + mass[i] * mass[i]);
2919 e_sum += e;
2920 }
2921
2922 // Return invariant mass with (+, -, -, -) metric
2923 return std::sqrt(e_sum * e_sum - x_sum * x_sum - y_sum * y_sum - z_sum * z_sum);
2924}
2925
2926////////////////////////////////////////////////////////////////////////////
2927/// \brief Build an RVec of objects starting from RVecs of input to their constructors.
2928/// \tparam T Type of the objects contained in the created RVec.
2929/// \tparam Args_t Pack of types templating the input RVecs.
2930/// \param[in] args The RVecs containing the values used to initialise the output objects.
2931/// \return The RVec of objects initialised with the input parameters.
2932///
2933/// Example code, at the ROOT prompt:
2934/// ~~~{.cpp}
2935/// using namespace ROOT::VecOps;
2936/// RVecF pts = {15.5, 34.32, 12.95};
2937/// RVecF etas = {0.3, 2.2, 1.32};
2938/// RVecF phis = {0.1, 3.02, 2.2};
2939/// RVecF masses = {105.65, 105.65, 105.65};
2940/// auto fourVecs = Construct<ROOT::Math::PtEtaPhiMVector>(pts, etas, phis, masses);
2941/// cout << fourVecs << endl;
2942/// // { (15.5,0.3,0.1,105.65), (34.32,2.2,3.02,105.65), (12.95,1.32,2.2,105.65) }
2943/// ~~~
2944template <typename T, typename... Args_t>
2946{
2947 const auto size = ::ROOT::Internal::VecOps::GetVectorsSize("Construct", args...);
2948 RVec<T> ret;
2949 ret.reserve(size);
2950 for (auto i = 0UL; i < size; ++i) {
2951 ret.emplace_back(args[i]...);
2952 }
2953 return ret;
2954}
2955
2956////////////////////////////////////////////////////////////////////////////////
2957/// Print a RVec at the prompt:
2958template <class T>
2959std::ostream &operator<<(std::ostream &os, const RVec<T> &v)
2960{
2961 // In order to print properly, convert to 64 bit int if this is a char
2962 constexpr bool mustConvert = std::is_same<char, T>::value || std::is_same<signed char, T>::value ||
2966 os << "{ ";
2967 auto size = v.size();
2968 if (size) {
2969 for (std::size_t i = 0; i < size - 1; ++i) {
2970 os << (Print_t)v[i] << ", ";
2971 }
2972 os << (Print_t)v[size - 1];
2973 }
2974 os << " }";
2975 return os;
2976}
2977
2978#if (_VECOPS_USE_EXTERN_TEMPLATES)
2979
2980#define RVEC_EXTERN_UNARY_OPERATOR(T, OP) \
2981 extern template RVec<T> operator OP<T>(const RVec<T> &);
2982
2983#define RVEC_EXTERN_BINARY_OPERATOR(T, OP) \
2984 extern template auto operator OP<T, T>(const T &x, const RVec<T> &v) \
2985 -> RVec<decltype(x OP v[0])>; \
2986 extern template auto operator OP<T, T>(const RVec<T> &v, const T &y) \
2987 -> RVec<decltype(v[0] OP y)>; \
2988 extern template auto operator OP<T, T>(const RVec<T> &v0, const RVec<T> &v1)\
2989 -> RVec<decltype(v0[0] OP v1[0])>;
2990
2991#define RVEC_EXTERN_ASSIGN_OPERATOR(T, OP) \
2992 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const T &); \
2993 extern template RVec<T> &operator OP<T, T>(RVec<T> &, const RVec<T> &);
2994
2995#define RVEC_EXTERN_LOGICAL_OPERATOR(T, OP) \
2996 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const T &); \
2997 extern template RVec<int> operator OP<T, T>(const T &, const RVec<T> &); \
2998 extern template RVec<int> operator OP<T, T>(const RVec<T> &, const RVec<T> &);
2999
3000#define RVEC_EXTERN_FLOAT_TEMPLATE(T) \
3001 extern template class RVec<T>; \
3002 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
3003 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
3004 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
3005 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
3006 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
3007 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
3008 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
3009 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
3010 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
3011 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
3012 RVEC_EXTERN_ASSIGN_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 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
3021
3022#define RVEC_EXTERN_INTEGER_TEMPLATE(T) \
3023 extern template class RVec<T>; \
3024 RVEC_EXTERN_UNARY_OPERATOR(T, +) \
3025 RVEC_EXTERN_UNARY_OPERATOR(T, -) \
3026 RVEC_EXTERN_UNARY_OPERATOR(T, ~) \
3027 RVEC_EXTERN_UNARY_OPERATOR(T, !) \
3028 RVEC_EXTERN_BINARY_OPERATOR(T, +) \
3029 RVEC_EXTERN_BINARY_OPERATOR(T, -) \
3030 RVEC_EXTERN_BINARY_OPERATOR(T, *) \
3031 RVEC_EXTERN_BINARY_OPERATOR(T, /) \
3032 RVEC_EXTERN_BINARY_OPERATOR(T, %) \
3033 RVEC_EXTERN_BINARY_OPERATOR(T, &) \
3034 RVEC_EXTERN_BINARY_OPERATOR(T, |) \
3035 RVEC_EXTERN_BINARY_OPERATOR(T, ^) \
3036 RVEC_EXTERN_ASSIGN_OPERATOR(T, +=) \
3037 RVEC_EXTERN_ASSIGN_OPERATOR(T, -=) \
3038 RVEC_EXTERN_ASSIGN_OPERATOR(T, *=) \
3039 RVEC_EXTERN_ASSIGN_OPERATOR(T, /=) \
3040 RVEC_EXTERN_ASSIGN_OPERATOR(T, %=) \
3041 RVEC_EXTERN_ASSIGN_OPERATOR(T, &=) \
3042 RVEC_EXTERN_ASSIGN_OPERATOR(T, |=) \
3043 RVEC_EXTERN_ASSIGN_OPERATOR(T, ^=) \
3044 RVEC_EXTERN_ASSIGN_OPERATOR(T, >>=) \
3045 RVEC_EXTERN_ASSIGN_OPERATOR(T, <<=) \
3046 RVEC_EXTERN_LOGICAL_OPERATOR(T, <) \
3047 RVEC_EXTERN_LOGICAL_OPERATOR(T, >) \
3048 RVEC_EXTERN_LOGICAL_OPERATOR(T, ==) \
3049 RVEC_EXTERN_LOGICAL_OPERATOR(T, !=) \
3050 RVEC_EXTERN_LOGICAL_OPERATOR(T, <=) \
3051 RVEC_EXTERN_LOGICAL_OPERATOR(T, >=) \
3052 RVEC_EXTERN_LOGICAL_OPERATOR(T, &&) \
3053 RVEC_EXTERN_LOGICAL_OPERATOR(T, ||)
3054
3055RVEC_EXTERN_INTEGER_TEMPLATE(char)
3056RVEC_EXTERN_INTEGER_TEMPLATE(short)
3057RVEC_EXTERN_INTEGER_TEMPLATE(int)
3058RVEC_EXTERN_INTEGER_TEMPLATE(long)
3059//RVEC_EXTERN_INTEGER_TEMPLATE(long long)
3060
3061RVEC_EXTERN_INTEGER_TEMPLATE(unsigned char)
3062RVEC_EXTERN_INTEGER_TEMPLATE(unsigned short)
3063RVEC_EXTERN_INTEGER_TEMPLATE(unsigned int)
3064RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long)
3065//RVEC_EXTERN_INTEGER_TEMPLATE(unsigned long long)
3066
3067RVEC_EXTERN_FLOAT_TEMPLATE(float)
3068RVEC_EXTERN_FLOAT_TEMPLATE(double)
3069
3070#undef RVEC_EXTERN_UNARY_OPERATOR
3071#undef RVEC_EXTERN_BINARY_OPERATOR
3072#undef RVEC_EXTERN_ASSIGN_OPERATOR
3073#undef RVEC_EXTERN_LOGICAL_OPERATOR
3074#undef RVEC_EXTERN_INTEGER_TEMPLATE
3075#undef RVEC_EXTERN_FLOAT_TEMPLATE
3076
3077#define RVEC_EXTERN_UNARY_FUNCTION(T, NAME, FUNC) \
3078 extern template RVec<PromoteType<T>> NAME(const RVec<T> &);
3079
3080#define RVEC_EXTERN_STD_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, std::F)
3081
3082#define RVEC_EXTERN_BINARY_FUNCTION(T0, T1, NAME, FUNC) \
3083 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const T1 &); \
3084 extern template RVec<PromoteTypes<T0, T1>> NAME(const T0 &, const RVec<T1> &); \
3085 extern template RVec<PromoteTypes<T0, T1>> NAME(const RVec<T0> &, const RVec<T1> &);
3086
3087#define RVEC_EXTERN_STD_BINARY_FUNCTION(T, F) RVEC_EXTERN_BINARY_FUNCTION(T, T, F, std::F)
3088
3089#define RVEC_EXTERN_STD_FUNCTIONS(T) \
3090 RVEC_EXTERN_STD_UNARY_FUNCTION(T, abs) \
3091 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fdim) \
3092 RVEC_EXTERN_STD_BINARY_FUNCTION(T, fmod) \
3093 RVEC_EXTERN_STD_BINARY_FUNCTION(T, remainder) \
3094 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp) \
3095 RVEC_EXTERN_STD_UNARY_FUNCTION(T, exp2) \
3096 RVEC_EXTERN_STD_UNARY_FUNCTION(T, expm1) \
3097 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log) \
3098 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log10) \
3099 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log2) \
3100 RVEC_EXTERN_STD_UNARY_FUNCTION(T, log1p) \
3101 RVEC_EXTERN_STD_BINARY_FUNCTION(T, pow) \
3102 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sqrt) \
3103 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cbrt) \
3104 RVEC_EXTERN_STD_BINARY_FUNCTION(T, hypot) \
3105 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sin) \
3106 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cos) \
3107 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tan) \
3108 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asin) \
3109 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acos) \
3110 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atan) \
3111 RVEC_EXTERN_STD_BINARY_FUNCTION(T, atan2) \
3112 RVEC_EXTERN_STD_UNARY_FUNCTION(T, sinh) \
3113 RVEC_EXTERN_STD_UNARY_FUNCTION(T, cosh) \
3114 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tanh) \
3115 RVEC_EXTERN_STD_UNARY_FUNCTION(T, asinh) \
3116 RVEC_EXTERN_STD_UNARY_FUNCTION(T, acosh) \
3117 RVEC_EXTERN_STD_UNARY_FUNCTION(T, atanh) \
3118 RVEC_EXTERN_STD_UNARY_FUNCTION(T, floor) \
3119 RVEC_EXTERN_STD_UNARY_FUNCTION(T, ceil) \
3120 RVEC_EXTERN_STD_UNARY_FUNCTION(T, trunc) \
3121 RVEC_EXTERN_STD_UNARY_FUNCTION(T, round) \
3122 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erf) \
3123 RVEC_EXTERN_STD_UNARY_FUNCTION(T, erfc) \
3124 RVEC_EXTERN_STD_UNARY_FUNCTION(T, lgamma) \
3125 RVEC_EXTERN_STD_UNARY_FUNCTION(T, tgamma) \
3126
3127RVEC_EXTERN_STD_FUNCTIONS(float)
3128RVEC_EXTERN_STD_FUNCTIONS(double)
3129#undef RVEC_EXTERN_STD_UNARY_FUNCTION
3130#undef RVEC_EXTERN_STD_BINARY_FUNCTION
3131#undef RVEC_EXTERN_STD_UNARY_FUNCTIONS
3132
3133#ifdef R__HAS_VDT
3134
3135#define RVEC_EXTERN_VDT_UNARY_FUNCTION(T, F) RVEC_EXTERN_UNARY_FUNCTION(T, F, vdt::F)
3136
3137RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_expf)
3138RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_logf)
3139RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_sinf)
3140RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_cosf)
3141RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_tanf)
3142RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_asinf)
3143RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_acosf)
3144RVEC_EXTERN_VDT_UNARY_FUNCTION(float, fast_atanf)
3145
3146RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_exp)
3147RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_log)
3148RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_sin)
3149RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_cos)
3150RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_tan)
3151RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_asin)
3152RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_acos)
3153RVEC_EXTERN_VDT_UNARY_FUNCTION(double, fast_atan)
3154
3155#endif // R__HAS_VDT
3156
3157#endif // _VECOPS_USE_EXTERN_TEMPLATES
3158
3159/** @} */ // end of Doxygen group vecops
3160
3161} // End of VecOps NS
3162
3163// Allow to use RVec as ROOT::RVec
3164using ROOT::VecOps::RVec;
3165
3176
3177} // End of ROOT NS
3178
3179#endif // ROOT_RVEC
size_t fSize
#define R__unlikely(expr)
Definition: RConfig.hxx:595
#define f(i)
Definition: RSha256.hxx:104
#define c(i)
Definition: RSha256.hxx:101
#define 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
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#define M_PI
Definition: Rotated.cxx:105
#define R__CLING_PTRCHECK(ONOFF)
Definition: Rtypes.h:498
#define R__ASSERT(e)
Definition: TError.h:118
#define N
Double_t Dot(const TGLVector3 &v1, const TGLVector3 &v2)
Definition: TGLUtil.h:317
Int_t Compare(const void *item1, const void *item2)
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char y2
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
Option_t Option_t TPoint TPoint const char y1
Binding & operator=(OUT(*fun)(void))
#define free
Definition: civetweb.c:1539
#define malloc
Definition: civetweb.c:1536
This class consists of common code factored out of the SmallVector class to reduce code duplication b...
Definition: RVec.hxx: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:953
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:1006
void pop_back_n(size_type NumItems)
Definition: RVec.hxx:589
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
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
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:1123
RVecN(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1159
reference operator[](size_type idx)
Definition: RVec.hxx:1199
RVecN operator[](const RVecN< V, M > &conds) const
Definition: RVec.hxx:1210
RVecN(std::initializer_list< T > IL)
Definition: RVec.hxx:1139
const_reference at(size_type pos) const
Definition: RVec.hxx:1245
RVecN(const RVecN &RHS)
Definition: RVec.hxx:1141
RVecN & operator=(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1180
RVecN & operator=(RVecN &&RHS)
Definition: RVec.hxx:1167
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:1261
RVecN & operator=(std::initializer_list< T > IL)
Definition: RVec.hxx:1186
RVecN & operator=(const RVecN &RHS)
Definition: RVec.hxx:1147
RVecN(const std::vector< T > &RHS)
Definition: RVec.hxx:1165
RVecN(size_t Size, const T &Value)
Definition: RVec.hxx:1121
RVecN(RVecN &&RHS)
Definition: RVec.hxx:1153
RVecN(ItTy S, ItTy E)
Definition: RVec.hxx:1134
reference at(size_type pos)
Definition: RVec.hxx:1238
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:1253
RVecN(T *p, size_t n)
Definition: RVec.hxx:1173
const_reference operator[](size_type idx) const
Definition: RVec.hxx:1204
A "std::vector"-like collection of values implementing handy operation to analyse them.
Definition: RVec.hxx:1439
RVec(RVecN< T, N > &&RHS)
Definition: RVec.hxx:1483
RVec(const RVecN< T, N > &RHS)
Definition: RVec.hxx:1486
RVec(size_t Size, const T &Value)
Definition: RVec.hxx:1451
RVec(const RVec &RHS)
Definition: RVec.hxx:1464
RVec & operator=(RVec &&RHS)
Definition: RVec.hxx:1474
RVec(T *p, size_t n)
Definition: RVec.hxx:1490
RVec operator[](const RVec< V > &conds) const
Definition: RVec.hxx:1502
RVec(std::initializer_list< T > IL)
Definition: RVec.hxx:1462
RVec(size_t Size)
Definition: RVec.hxx:1453
RVec(ItTy S, ItTy E)
Definition: RVec.hxx:1458
RVec(const std::vector< T > &RHS)
Definition: RVec.hxx:1488
typename SuperClass::size_type size_type
Definition: RVec.hxx:1444
RVec(Detail::VecOps::RVecImpl< T > &&RHS)
Definition: RVec.hxx:1480
RVec(RVec &&RHS)
Definition: RVec.hxx:1472
typename SuperClass::value_type value_type
Definition: RVec.hxx:1445
RVec & operator=(const RVec &RHS)
Definition: RVec.hxx:1466
TPaveText * pt
Type
enumeration specifying the integration types.
T Max(const RVec< T > &v)
Get the greatest element of an RVec.
Definition: RVec.hxx:1928
RVec< T > Reverse(const RVec< T > &v)
Return copy of reversed vector.
Definition: RVec.hxx:2310
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:2799
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:2703
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:2587
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:2112
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:2492
RVec< PromoteType< T > > asinh(const RVec< T > &v)
Definition: RVec.hxx:1769
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:2159
RVec< PromoteType< T > > expm1(const RVec< T > &v)
Definition: RVec.hxx:1746
RVec< PromoteType< T > > lgamma(const RVec< T > &v)
Definition: RVec.hxx:1782
RVec< typename RVec< T >::size_type > Nonzero(const RVec< T > &v)
Return the indices of the elements which are not zero.
Definition: RVec.hxx:2556
#define RVEC_UNARY_OPERATOR(OP)
Definition: RVec.hxx:1523
RVec< PromoteType< T > > cosh(const RVec< T > &v)
Definition: RVec.hxx:1767
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:1911
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1739
#define RVEC_ASSIGNMENT_OPERATOR(OP)
Definition: RVec.hxx:1594
RVec< PromoteType< T > > round(const RVec< T > &v)
Definition: RVec.hxx:1776
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:2899
RVec< PromoteType< T > > tgamma(const RVec< T > &v)
Definition: RVec.hxx:1783
RVec< PromoteType< T > > log2(const RVec< T > &v)
Definition: RVec.hxx:1750
RVec< Common_t > Concatenate(const RVec< T0 > &v0, const RVec< T1 > &v1)
Return the concatenation of two RVecs.
Definition: RVec.hxx:2726
R Sum(const RVec< T > &v, const R zero=R(0))
Sum elements of an RVec.
Definition: RVec.hxx:1861
RVec< PromoteType< T > > tan(const RVec< T > &v)
Definition: RVec.hxx:1760
RVec< PromoteType< T > > erf(const RVec< T > &v)
Definition: RVec.hxx:1780
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:2416
RVec< PromoteType< T > > cos(const RVec< T > &v)
Definition: RVec.hxx:1759
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:2857
RVec< PromoteType< T > > floor(const RVec< T > &v)
Definition: RVec.hxx:1773
RVec< PromoteType< T > > atanh(const RVec< T > &v)
Definition: RVec.hxx:1771
RVec< PromoteType< T > > acosh(const RVec< T > &v)
Definition: RVec.hxx:1770
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:2818
RVec< PromoteTypes< T0, T1 > > remainder(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1742
RVec< PromoteType< T > > ceil(const RVec< T > &v)
Definition: RVec.hxx:1774
RVec< PromoteType< T > > cbrt(const RVec< T > &v)
Definition: RVec.hxx:1755
RVec< PromoteType< T > > asin(const RVec< T > &v)
Definition: RVec.hxx:1761
void swap(RVec< T > &lhs, RVec< T > &rhs)
Definition: RVec.hxx:2121
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1748
RVec< T > Construct(const RVec< Args_t > &... args)
Build an RVec of objects starting from RVecs of input to their constructors.
Definition: RVec.hxx:2945
RVec< PromoteType< T > > trunc(const RVec< T > &v)
Definition: RVec.hxx:1775
RVec< PromoteType< T > > llround(const RVec< T > &v)
Definition: RVec.hxx:1778
RVec< PromoteTypes< T0, T1 > > pow(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1753
RVec< PromoteTypes< T0, T1 > > hypot(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1756
#define RVEC_STD_BINARY_FUNCTION(F)
Definition: RVec.hxx:1737
RVec< PromoteType< T > > lround(const RVec< T > &v)
Definition: RVec.hxx:1777
#define RVEC_BINARY_OPERATOR(OP)
Definition: RVec.hxx:1546
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:2277
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:2206
T Min(const RVec< T > &v)
Get the smallest element of an RVec.
Definition: RVec.hxx:1944
size_t CapacityInBytes(const RVecN< T, N > &X)
Definition: RVec.hxx:1515
#define RVEC_LOGICAL_OPERATOR(OP)
Definition: RVec.hxx:1630
RVec< T > Take(const RVec< T > &v, const int n)
Return first or last n elements of an RVec.
Definition: RVec.hxx:2252
RVec< PromoteType< T > > sinh(const RVec< T > &v)
Definition: RVec.hxx:1766
RVec< PromoteType< T > > sqrt(const RVec< T > &v)
Definition: RVec.hxx:1754
RVec< PromoteType< T > > erfc(const RVec< T > &v)
Definition: RVec.hxx:1781
RVec< T > Sort(const RVec< T > &v, Compare &&c)
Return copy of RVec with elements sorted based on a comparison operator.
Definition: RVec.hxx:2356
#define RVEC_STD_UNARY_FUNCTION(F)
Definition: RVec.hxx:1736
RVec< PromoteTypes< T0, T1 > > atan2(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1764
RVec< PromoteTypes< T0, T1 > > fmod(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1741
RVec< PromoteType< T > > log10(const RVec< T > &v)
Definition: RVec.hxx:1749
auto Map(Args &&... args)
Create new collection applying a callable to the elements of the input collection.
Definition: RVec.hxx:2045
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1744
RVec< PromoteTypes< T0, T1 > > fdim(const RVec< T0 > &v0, const RVec< T1 > &v1)
Definition: RVec.hxx:1740
double StdDev(const RVec< T > &v)
Get the standard deviation of the elements of an RVec.
Definition: RVec.hxx:2021
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:2093
RVec< PromoteType< T > > tanh(const RVec< T > &v)
Definition: RVec.hxx:1768
RVec< PromoteType< T > > acos(const RVec< T > &v)
Definition: RVec.hxx:1762
RVec< PromoteType< T > > log1p(const RVec< T > &v)
Definition: RVec.hxx:1751
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:1980
RVec< PromoteType< T > > atan(const RVec< T > &v)
Definition: RVec.hxx:1763
RVec< PromoteType< T > > sin(const RVec< T > &v)
Definition: RVec.hxx:1758
double Var(const RVec< T > &v)
Get the variance of the elements of an RVec.
Definition: RVec.hxx:1997
RVec< PromoteType< T > > exp2(const RVec< T > &v)
Definition: RVec.hxx:1745
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:2070
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:1962
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:2845
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
#define T2
Definition: md5.inl:147
#define F(x, y, z)
#define I(x, y, z)
#define T1
Definition: md5.inl:146
bool IsSmall(const ROOT::VecOps::RVec< T > &v)
Definition: RVec.hxx:1076
bool IsAdopting(const ROOT::VecOps::RVec< T > &v)
Definition: RVec.hxx:1082
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 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:241
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