1 | /* |
2 | * Copyright 2018-present Facebook, Inc. |
3 | * |
4 | * Licensed under the Apache License, Version 2.0 (the "License"); |
5 | * you may not use this file except in compliance with the License. |
6 | * You may obtain a copy of the License at |
7 | * |
8 | * http://www.apache.org/licenses/LICENSE-2.0 |
9 | * |
10 | * Unless required by applicable law or agreed to in writing, software |
11 | * distributed under the License is distributed on an "AS IS" BASIS, |
12 | * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. |
13 | * See the License for the specific language governing permissions and |
14 | * limitations under the License. |
15 | */ |
16 | |
17 | #include <folly/stats/detail/DoubleRadixSort.h> |
18 | |
19 | #include <algorithm> |
20 | #include <cstring> |
21 | |
22 | namespace folly { |
23 | namespace detail { |
24 | |
25 | // Convert floats to something comparable via radix sort. |
26 | // See http://stereopsis.com/radix.html for details. |
27 | static uint8_t getRadixBucket(double* f, uint8_t shift) { |
28 | uint64_t val; |
29 | memcpy(&val, f, sizeof(double)); |
30 | uint64_t mask = -int64_t(val >> 63) | 0x8000000000000000; |
31 | auto adjusted = val ^ mask; |
32 | return (adjusted >> (64 - 8 - shift)) & 0xFF; |
33 | } |
34 | |
35 | // MSB radix sort for doubles. |
36 | static void double_radix_sort_rec( |
37 | uint64_t n, |
38 | uint64_t* buckets, |
39 | uint8_t shift, |
40 | bool inout, |
41 | double* in, |
42 | double* out) { |
43 | // First pass: calculate bucket counts. |
44 | memset(buckets, 0, 256 * sizeof(uint64_t)); |
45 | for (uint64_t i = 0; i < n; i++) { |
46 | buckets[getRadixBucket(&in[i], shift)]++; |
47 | } |
48 | |
49 | // Second pass: calculate bucket start positions. |
50 | uint64_t tot = 0; |
51 | for (uint64_t i = 0; i < 256; i++) { |
52 | auto prev = tot; |
53 | tot += buckets[i]; |
54 | buckets[i + 256] = prev; |
55 | } |
56 | |
57 | // Third pass: Move based on radix counts. |
58 | for (uint64_t i = 0; i < n; i++) { |
59 | auto pos = buckets[getRadixBucket(&in[i], shift) + 256]++; |
60 | out[pos] = in[i]; |
61 | } |
62 | |
63 | // If we haven't used up all input bytes, recurse and sort. if the |
64 | // bucket is too small, use std::sort instead, and copy output to |
65 | // correct array. |
66 | if (shift < 56) { |
67 | tot = 0; |
68 | for (int i = 0; i < 256; i++) { |
69 | if (buckets[i] < 256) { |
70 | std::sort(out + tot, out + tot + buckets[i]); |
71 | if (!inout) { |
72 | memcpy(in + tot, out + tot, buckets[i] * sizeof(double)); |
73 | } |
74 | } else { |
75 | double_radix_sort_rec( |
76 | buckets[i], buckets + 256, shift + 8, !inout, out + tot, in + tot); |
77 | } |
78 | tot += buckets[i]; |
79 | } |
80 | } |
81 | } |
82 | |
83 | void double_radix_sort(uint64_t n, uint64_t* buckets, double* in, double* tmp) { |
84 | // If array is too small, use std::sort directly. |
85 | if (n < 700) { |
86 | std::sort(in, in + n); |
87 | } else { |
88 | detail::double_radix_sort_rec(n, buckets, 0, false, in, tmp); |
89 | } |
90 | } |
91 | |
92 | } // namespace detail |
93 | } // namespace folly |
94 | |