From c3dceb92fbc77797c9d298c12c5469bf3a425c71 Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Fri, 8 Aug 2025 10:42:16 +0200 Subject: [PATCH] [Math] Fix bug and warning in `TMath::ModeHalfSample` MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The funciton de-duplicates the input values, but the logic was wrong. `std::lower_bound` returns a hit even if there is no exact match, it just looks for the first element that is larger. Therefore, we need an additional value check. Also, fixes the following compiler warning that relates to doing `data() + 1` which is not valid if the vector would be empty: ``` In member function ‘void std::__new_allocator<_Tp>::deallocate(_Tp*, size_type) [with _Tp = double]’, inlined from ‘static void std::allocator_traits >::deallocate(allocator_type&, pointer, size_type) [with _Tp = double]’ at /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/alloc_traits.h:550:23, inlined from ‘void std::_Vector_base<_Tp, _Alloc>::_M_deallocate(pointer, std::size_t) [with _Tp = double; _Alloc = std::allocator]’ at /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/stl_vector.h:392:19, inlined from ‘void std::_Vector_base<_Tp, _Alloc>::_M_deallocate(pointer, std::size_t) [with _Tp = double; _Alloc = std::allocator]’ at /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/stl_vector.h:388:7, inlined from ‘std::_Vector_base<_Tp, _Alloc>::~_Vector_base() [with _Tp = double; _Alloc = std::allocator]’ at /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/stl_vector.h:371:15, inlined from ‘std::vector<_Tp, _Alloc>::~vector() [with _Tp = double; _Alloc = std::allocator]’ at /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/stl_vector.h:751:7, inlined from ‘Double_t TMath::ModeHalfSample(Long64_t, const T*, const Double_t*) [with T = short int]’ at /home/rembserj/code/root/root_build/include/TMath.h:1540:1: /nix/store/eeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee-gcc-14.3.0/include/c++/14.3.0/bits/new_allocator.h:172:26: warning: ‘void operator delete(void*, std::size_t)’ called on a pointer to an unallocated object ‘18446744073709551608’ [-Wfree-nonheap-object] ``` --- math/mathcore/inc/TMath.h | 42 ++++++++++++++++++++++++--------------- 1 file changed, 26 insertions(+), 16 deletions(-) diff --git a/math/mathcore/inc/TMath.h b/math/mathcore/inc/TMath.h index d102d649f3420..2d37cc34d0808 100644 --- a/math/mathcore/inc/TMath.h +++ b/math/mathcore/inc/TMath.h @@ -1451,23 +1451,33 @@ template Double_t TMath::ModeHalfSample(Long64_t n, const T *a, con if (w && std::adjacent_find(w, w + n, std::not_equal_to<>()) == w + n) // If all weights are equal, ignore those w = nullptr; - // Sort the array and remove duplicates (if w != nullptr) + // Sort the array std::vector values; std::vector weights; - values.reserve(n); - weights.reserve(n); - for (Long64_t i = 0; i < n; ++i) { - const T value = a[i]; - const Double_t weight = w ? w[i] : 1; - const auto vstart = values.begin(); - const auto vstop = values.end(); - const auto viter = std::lower_bound(vstart, vstop, value); - const auto vidx = std::distance(vstart, viter); - if (viter != vstop && w) { - weights[vidx] += weight; - } else { - values.insert(viter, value); - weights.insert(weights.begin() + vidx, weight); + if (w == nullptr) { + values.assign(a, a + n); + weights.assign(n, 1.0); + std::sort(values.begin(), values.end()); + } else { + // Remove duplicates (if w != nullptr). + // The Bickel algorithm expects all values have the same weight of 1. + // The deduplication is relevant for varying weights, where we fall back + // to getting the Mode of the sample. + values.reserve(n); + weights.reserve(n); + for (Long64_t i = 0; i < n; ++i) { + const T value = a[i]; + const auto vstart = values.begin(); + const auto viter = std::lower_bound(vstart, values.end(), value); + const auto vidx = std::distance(vstart, viter); + if (viter != values.end() && *viter == value) { + // increase weight if value already added + weights[vidx] += w[i]; + } else { + // otherwise, insert at the right position + values.insert(viter, value); + weights.insert(weights.begin() + vidx, w[i]); + } } } @@ -1483,7 +1493,7 @@ template Double_t TMath::ModeHalfSample(Long64_t n, const T *a, con const double d1_0 = values[1] - values[0]; const double d2_1 = values[2] - values[1]; if (d2_1 < d1_0) - return TMath::Mean(2, values.data() + 1, weights.data() + 1); + return TMath::Mean(std::next(values.begin()), values.end(), std::next(weights.begin())); else if (d2_1 > d1_0) return TMath::Mean(2, values.data(), weights.data()); else