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