1bf3afd5fSEmilio G. Cota /*
2bf3afd5fSEmilio G. Cota * qdist.c - QEMU helpers for handling frequency distributions of data.
3bf3afd5fSEmilio G. Cota *
4bf3afd5fSEmilio G. Cota * Copyright (C) 2016, Emilio G. Cota <cota@braap.org>
5bf3afd5fSEmilio G. Cota *
6bf3afd5fSEmilio G. Cota * License: GNU GPL, version 2 or later.
7bf3afd5fSEmilio G. Cota * See the COPYING file in the top-level directory.
8bf3afd5fSEmilio G. Cota */
9e9abfcb5SPaolo Bonzini #include "qemu/osdep.h"
10bf3afd5fSEmilio G. Cota #include "qemu/qdist.h"
11bf3afd5fSEmilio G. Cota
12bf3afd5fSEmilio G. Cota #include <math.h>
13bf3afd5fSEmilio G. Cota #ifndef NAN
14bf3afd5fSEmilio G. Cota #define NAN (0.0 / 0.0)
15bf3afd5fSEmilio G. Cota #endif
16bf3afd5fSEmilio G. Cota
1711b7b07fSEmilio G. Cota #define QDIST_EMPTY_STR "(empty)"
1811b7b07fSEmilio G. Cota
qdist_init(struct qdist * dist)19bf3afd5fSEmilio G. Cota void qdist_init(struct qdist *dist)
20bf3afd5fSEmilio G. Cota {
21071d4054SEmilio G. Cota dist->entries = g_new(struct qdist_entry, 1);
22bf3afd5fSEmilio G. Cota dist->size = 1;
23bf3afd5fSEmilio G. Cota dist->n = 0;
24bf3afd5fSEmilio G. Cota }
25bf3afd5fSEmilio G. Cota
qdist_destroy(struct qdist * dist)26bf3afd5fSEmilio G. Cota void qdist_destroy(struct qdist *dist)
27bf3afd5fSEmilio G. Cota {
28bf3afd5fSEmilio G. Cota g_free(dist->entries);
29bf3afd5fSEmilio G. Cota }
30bf3afd5fSEmilio G. Cota
qdist_cmp_double(double a,double b)31bf3afd5fSEmilio G. Cota static inline int qdist_cmp_double(double a, double b)
32bf3afd5fSEmilio G. Cota {
33bf3afd5fSEmilio G. Cota if (a > b) {
34bf3afd5fSEmilio G. Cota return 1;
35bf3afd5fSEmilio G. Cota } else if (a < b) {
36bf3afd5fSEmilio G. Cota return -1;
37bf3afd5fSEmilio G. Cota }
38bf3afd5fSEmilio G. Cota return 0;
39bf3afd5fSEmilio G. Cota }
40bf3afd5fSEmilio G. Cota
qdist_cmp(const void * ap,const void * bp)41bf3afd5fSEmilio G. Cota static int qdist_cmp(const void *ap, const void *bp)
42bf3afd5fSEmilio G. Cota {
43bf3afd5fSEmilio G. Cota const struct qdist_entry *a = ap;
44bf3afd5fSEmilio G. Cota const struct qdist_entry *b = bp;
45bf3afd5fSEmilio G. Cota
46bf3afd5fSEmilio G. Cota return qdist_cmp_double(a->x, b->x);
47bf3afd5fSEmilio G. Cota }
48bf3afd5fSEmilio G. Cota
qdist_add(struct qdist * dist,double x,long count)49bf3afd5fSEmilio G. Cota void qdist_add(struct qdist *dist, double x, long count)
50bf3afd5fSEmilio G. Cota {
51bf3afd5fSEmilio G. Cota struct qdist_entry *entry = NULL;
52bf3afd5fSEmilio G. Cota
53bf3afd5fSEmilio G. Cota if (dist->n) {
54bf3afd5fSEmilio G. Cota struct qdist_entry e;
55bf3afd5fSEmilio G. Cota
56bf3afd5fSEmilio G. Cota e.x = x;
57bf3afd5fSEmilio G. Cota entry = bsearch(&e, dist->entries, dist->n, sizeof(e), qdist_cmp);
58bf3afd5fSEmilio G. Cota }
59bf3afd5fSEmilio G. Cota
60bf3afd5fSEmilio G. Cota if (entry) {
61bf3afd5fSEmilio G. Cota entry->count += count;
62bf3afd5fSEmilio G. Cota return;
63bf3afd5fSEmilio G. Cota }
64bf3afd5fSEmilio G. Cota
65bf3afd5fSEmilio G. Cota if (unlikely(dist->n == dist->size)) {
66bf3afd5fSEmilio G. Cota dist->size *= 2;
67071d4054SEmilio G. Cota dist->entries = g_renew(struct qdist_entry, dist->entries, dist->size);
68bf3afd5fSEmilio G. Cota }
69bf3afd5fSEmilio G. Cota dist->n++;
70bf3afd5fSEmilio G. Cota entry = &dist->entries[dist->n - 1];
71bf3afd5fSEmilio G. Cota entry->x = x;
72bf3afd5fSEmilio G. Cota entry->count = count;
73bf3afd5fSEmilio G. Cota qsort(dist->entries, dist->n, sizeof(*entry), qdist_cmp);
74bf3afd5fSEmilio G. Cota }
75bf3afd5fSEmilio G. Cota
qdist_inc(struct qdist * dist,double x)76bf3afd5fSEmilio G. Cota void qdist_inc(struct qdist *dist, double x)
77bf3afd5fSEmilio G. Cota {
78bf3afd5fSEmilio G. Cota qdist_add(dist, x, 1);
79bf3afd5fSEmilio G. Cota }
80bf3afd5fSEmilio G. Cota
81bf3afd5fSEmilio G. Cota /*
82bf3afd5fSEmilio G. Cota * Unicode for block elements. See:
83bf3afd5fSEmilio G. Cota * https://en.wikipedia.org/wiki/Block_Elements
84bf3afd5fSEmilio G. Cota */
85bf3afd5fSEmilio G. Cota static const gunichar qdist_blocks[] = {
86bf3afd5fSEmilio G. Cota 0x2581,
87bf3afd5fSEmilio G. Cota 0x2582,
88bf3afd5fSEmilio G. Cota 0x2583,
89bf3afd5fSEmilio G. Cota 0x2584,
90bf3afd5fSEmilio G. Cota 0x2585,
91bf3afd5fSEmilio G. Cota 0x2586,
92bf3afd5fSEmilio G. Cota 0x2587,
93bf3afd5fSEmilio G. Cota 0x2588
94bf3afd5fSEmilio G. Cota };
95bf3afd5fSEmilio G. Cota
96bf3afd5fSEmilio G. Cota #define QDIST_NR_BLOCK_CODES ARRAY_SIZE(qdist_blocks)
97bf3afd5fSEmilio G. Cota
98bf3afd5fSEmilio G. Cota /*
99bf3afd5fSEmilio G. Cota * Print a distribution into a string.
100bf3afd5fSEmilio G. Cota *
101bf3afd5fSEmilio G. Cota * This function assumes that appropriate binning has been done on the input;
102bf3afd5fSEmilio G. Cota * see qdist_bin__internal() and qdist_pr_plain().
103bf3afd5fSEmilio G. Cota *
104bf3afd5fSEmilio G. Cota * Callers must free the returned string with g_free().
105bf3afd5fSEmilio G. Cota */
qdist_pr_internal(const struct qdist * dist)106bf3afd5fSEmilio G. Cota static char *qdist_pr_internal(const struct qdist *dist)
107bf3afd5fSEmilio G. Cota {
108bf3afd5fSEmilio G. Cota double min, max;
109bf3afd5fSEmilio G. Cota GString *s = g_string_new("");
110bf3afd5fSEmilio G. Cota size_t i;
111bf3afd5fSEmilio G. Cota
112bf3afd5fSEmilio G. Cota /* if only one entry, its printout will be either full or empty */
113bf3afd5fSEmilio G. Cota if (dist->n == 1) {
114bf3afd5fSEmilio G. Cota if (dist->entries[0].count) {
115bf3afd5fSEmilio G. Cota g_string_append_unichar(s, qdist_blocks[QDIST_NR_BLOCK_CODES - 1]);
116bf3afd5fSEmilio G. Cota } else {
117bf3afd5fSEmilio G. Cota g_string_append_c(s, ' ');
118bf3afd5fSEmilio G. Cota }
119bf3afd5fSEmilio G. Cota goto out;
120bf3afd5fSEmilio G. Cota }
121bf3afd5fSEmilio G. Cota
122bf3afd5fSEmilio G. Cota /* get min and max counts */
123bf3afd5fSEmilio G. Cota min = dist->entries[0].count;
124bf3afd5fSEmilio G. Cota max = min;
125bf3afd5fSEmilio G. Cota for (i = 0; i < dist->n; i++) {
126bf3afd5fSEmilio G. Cota struct qdist_entry *e = &dist->entries[i];
127bf3afd5fSEmilio G. Cota
128bf3afd5fSEmilio G. Cota if (e->count < min) {
129bf3afd5fSEmilio G. Cota min = e->count;
130bf3afd5fSEmilio G. Cota }
131bf3afd5fSEmilio G. Cota if (e->count > max) {
132bf3afd5fSEmilio G. Cota max = e->count;
133bf3afd5fSEmilio G. Cota }
134bf3afd5fSEmilio G. Cota }
135bf3afd5fSEmilio G. Cota
136bf3afd5fSEmilio G. Cota for (i = 0; i < dist->n; i++) {
137bf3afd5fSEmilio G. Cota struct qdist_entry *e = &dist->entries[i];
138bf3afd5fSEmilio G. Cota int index;
139bf3afd5fSEmilio G. Cota
140bf3afd5fSEmilio G. Cota /* make an exception with 0; instead of using block[0], print a space */
141bf3afd5fSEmilio G. Cota if (e->count) {
142bf3afd5fSEmilio G. Cota /* divide first to avoid loss of precision when e->count == max */
143bf3afd5fSEmilio G. Cota index = (e->count - min) / (max - min) * (QDIST_NR_BLOCK_CODES - 1);
144bf3afd5fSEmilio G. Cota g_string_append_unichar(s, qdist_blocks[index]);
145bf3afd5fSEmilio G. Cota } else {
146bf3afd5fSEmilio G. Cota g_string_append_c(s, ' ');
147bf3afd5fSEmilio G. Cota }
148bf3afd5fSEmilio G. Cota }
149bf3afd5fSEmilio G. Cota out:
150bf3afd5fSEmilio G. Cota return g_string_free(s, FALSE);
151bf3afd5fSEmilio G. Cota }
152bf3afd5fSEmilio G. Cota
153bf3afd5fSEmilio G. Cota /*
154bf3afd5fSEmilio G. Cota * Bin the distribution in @from into @n bins of consecutive, non-overlapping
155bf3afd5fSEmilio G. Cota * intervals, copying the result to @to.
156bf3afd5fSEmilio G. Cota *
157bf3afd5fSEmilio G. Cota * This function is internal to qdist: only this file and test code should
158bf3afd5fSEmilio G. Cota * ever call it.
159bf3afd5fSEmilio G. Cota *
160bf3afd5fSEmilio G. Cota * Note: calling this function on an already-binned qdist is a bug.
161bf3afd5fSEmilio G. Cota *
162bf3afd5fSEmilio G. Cota * If @n == 0 or @from->n == 1, use @from->n.
163bf3afd5fSEmilio G. Cota */
qdist_bin__internal(struct qdist * to,const struct qdist * from,size_t n)164bf3afd5fSEmilio G. Cota void qdist_bin__internal(struct qdist *to, const struct qdist *from, size_t n)
165bf3afd5fSEmilio G. Cota {
166bf3afd5fSEmilio G. Cota double xmin, xmax;
167bf3afd5fSEmilio G. Cota double step;
168bf3afd5fSEmilio G. Cota size_t i, j;
169bf3afd5fSEmilio G. Cota
170bf3afd5fSEmilio G. Cota qdist_init(to);
171bf3afd5fSEmilio G. Cota
172bf3afd5fSEmilio G. Cota if (from->n == 0) {
173bf3afd5fSEmilio G. Cota return;
174bf3afd5fSEmilio G. Cota }
175bf3afd5fSEmilio G. Cota if (n == 0 || from->n == 1) {
176bf3afd5fSEmilio G. Cota n = from->n;
177bf3afd5fSEmilio G. Cota }
178bf3afd5fSEmilio G. Cota
179bf3afd5fSEmilio G. Cota /* set equally-sized bins between @from's left and right */
180bf3afd5fSEmilio G. Cota xmin = qdist_xmin(from);
181bf3afd5fSEmilio G. Cota xmax = qdist_xmax(from);
182bf3afd5fSEmilio G. Cota step = (xmax - xmin) / n;
183bf3afd5fSEmilio G. Cota
184bf3afd5fSEmilio G. Cota if (n == from->n) {
185bf3afd5fSEmilio G. Cota /* if @from's entries are equally spaced, no need to re-bin */
186bf3afd5fSEmilio G. Cota for (i = 0; i < from->n; i++) {
187bf3afd5fSEmilio G. Cota if (from->entries[i].x != xmin + i * step) {
188bf3afd5fSEmilio G. Cota goto rebin;
189bf3afd5fSEmilio G. Cota }
190bf3afd5fSEmilio G. Cota }
191bf3afd5fSEmilio G. Cota /* they're equally spaced, so copy the dist and bail out */
192071d4054SEmilio G. Cota to->entries = g_renew(struct qdist_entry, to->entries, n);
193bf3afd5fSEmilio G. Cota to->n = from->n;
194bf3afd5fSEmilio G. Cota memcpy(to->entries, from->entries, sizeof(*to->entries) * to->n);
195bf3afd5fSEmilio G. Cota return;
196bf3afd5fSEmilio G. Cota }
197bf3afd5fSEmilio G. Cota
198bf3afd5fSEmilio G. Cota rebin:
199bf3afd5fSEmilio G. Cota j = 0;
200bf3afd5fSEmilio G. Cota for (i = 0; i < n; i++) {
201bf3afd5fSEmilio G. Cota double x;
202bf3afd5fSEmilio G. Cota double left, right;
203bf3afd5fSEmilio G. Cota
204bf3afd5fSEmilio G. Cota left = xmin + i * step;
205bf3afd5fSEmilio G. Cota right = xmin + (i + 1) * step;
206bf3afd5fSEmilio G. Cota
207bf3afd5fSEmilio G. Cota /* Add x, even if it might not get any counts later */
208bf3afd5fSEmilio G. Cota x = left;
209bf3afd5fSEmilio G. Cota qdist_add(to, x, 0);
210bf3afd5fSEmilio G. Cota
211bf3afd5fSEmilio G. Cota /*
212bf3afd5fSEmilio G. Cota * To avoid double-counting we capture [left, right) ranges, except for
213*d02d06f8SMichael Tokarev * the rightmost bin, which captures a [left, right] range.
214bf3afd5fSEmilio G. Cota */
215bf3afd5fSEmilio G. Cota while (j < from->n && (from->entries[j].x < right || i == n - 1)) {
216bf3afd5fSEmilio G. Cota struct qdist_entry *o = &from->entries[j];
217bf3afd5fSEmilio G. Cota
218bf3afd5fSEmilio G. Cota qdist_add(to, x, o->count);
219bf3afd5fSEmilio G. Cota j++;
220bf3afd5fSEmilio G. Cota }
221bf3afd5fSEmilio G. Cota }
222bf3afd5fSEmilio G. Cota }
223bf3afd5fSEmilio G. Cota
224bf3afd5fSEmilio G. Cota /*
225bf3afd5fSEmilio G. Cota * Print @dist into a string, after re-binning it into @n bins of consecutive,
226bf3afd5fSEmilio G. Cota * non-overlapping intervals.
227bf3afd5fSEmilio G. Cota *
228bf3afd5fSEmilio G. Cota * If @n == 0, use @orig->n.
229bf3afd5fSEmilio G. Cota *
230bf3afd5fSEmilio G. Cota * Callers must free the returned string with g_free().
231bf3afd5fSEmilio G. Cota */
qdist_pr_plain(const struct qdist * dist,size_t n)232bf3afd5fSEmilio G. Cota char *qdist_pr_plain(const struct qdist *dist, size_t n)
233bf3afd5fSEmilio G. Cota {
234bf3afd5fSEmilio G. Cota struct qdist binned;
235bf3afd5fSEmilio G. Cota char *ret;
236bf3afd5fSEmilio G. Cota
237bf3afd5fSEmilio G. Cota if (dist->n == 0) {
23811b7b07fSEmilio G. Cota return g_strdup(QDIST_EMPTY_STR);
239bf3afd5fSEmilio G. Cota }
240bf3afd5fSEmilio G. Cota qdist_bin__internal(&binned, dist, n);
241bf3afd5fSEmilio G. Cota ret = qdist_pr_internal(&binned);
242bf3afd5fSEmilio G. Cota qdist_destroy(&binned);
243bf3afd5fSEmilio G. Cota return ret;
244bf3afd5fSEmilio G. Cota }
245bf3afd5fSEmilio G. Cota
qdist_pr_label(const struct qdist * dist,size_t n_bins,uint32_t opt,bool is_left)246bf3afd5fSEmilio G. Cota static char *qdist_pr_label(const struct qdist *dist, size_t n_bins,
247bf3afd5fSEmilio G. Cota uint32_t opt, bool is_left)
248bf3afd5fSEmilio G. Cota {
249bf3afd5fSEmilio G. Cota const char *percent;
250bf3afd5fSEmilio G. Cota const char *lparen;
251bf3afd5fSEmilio G. Cota const char *rparen;
252bf3afd5fSEmilio G. Cota GString *s;
253bf3afd5fSEmilio G. Cota double x1, x2, step;
254bf3afd5fSEmilio G. Cota double x;
255bf3afd5fSEmilio G. Cota double n;
256bf3afd5fSEmilio G. Cota int dec;
257bf3afd5fSEmilio G. Cota
258bf3afd5fSEmilio G. Cota s = g_string_new("");
259bf3afd5fSEmilio G. Cota if (!(opt & QDIST_PR_LABELS)) {
260bf3afd5fSEmilio G. Cota goto out;
261bf3afd5fSEmilio G. Cota }
262bf3afd5fSEmilio G. Cota
263bf3afd5fSEmilio G. Cota dec = opt & QDIST_PR_NODECIMAL ? 0 : 1;
264bf3afd5fSEmilio G. Cota percent = opt & QDIST_PR_PERCENT ? "%" : "";
265bf3afd5fSEmilio G. Cota
266bf3afd5fSEmilio G. Cota n = n_bins ? n_bins : dist->n;
267bf3afd5fSEmilio G. Cota x = is_left ? qdist_xmin(dist) : qdist_xmax(dist);
268bf3afd5fSEmilio G. Cota step = (qdist_xmax(dist) - qdist_xmin(dist)) / n;
269bf3afd5fSEmilio G. Cota
270bf3afd5fSEmilio G. Cota if (opt & QDIST_PR_100X) {
271bf3afd5fSEmilio G. Cota x *= 100.0;
272bf3afd5fSEmilio G. Cota step *= 100.0;
273bf3afd5fSEmilio G. Cota }
274bf3afd5fSEmilio G. Cota if (opt & QDIST_PR_NOBINRANGE) {
275bf3afd5fSEmilio G. Cota lparen = rparen = "";
276bf3afd5fSEmilio G. Cota x1 = x;
277bf3afd5fSEmilio G. Cota x2 = x; /* unnecessary, but a dumb compiler might not figure it out */
278bf3afd5fSEmilio G. Cota } else {
279bf3afd5fSEmilio G. Cota lparen = "[";
280bf3afd5fSEmilio G. Cota rparen = is_left ? ")" : "]";
281bf3afd5fSEmilio G. Cota if (is_left) {
282bf3afd5fSEmilio G. Cota x1 = x;
283bf3afd5fSEmilio G. Cota x2 = x + step;
284bf3afd5fSEmilio G. Cota } else {
285bf3afd5fSEmilio G. Cota x1 = x - step;
286bf3afd5fSEmilio G. Cota x2 = x;
287bf3afd5fSEmilio G. Cota }
288bf3afd5fSEmilio G. Cota }
289bf3afd5fSEmilio G. Cota g_string_append_printf(s, "%s%.*f", lparen, dec, x1);
290bf3afd5fSEmilio G. Cota if (!(opt & QDIST_PR_NOBINRANGE)) {
291bf3afd5fSEmilio G. Cota g_string_append_printf(s, ",%.*f%s", dec, x2, rparen);
292bf3afd5fSEmilio G. Cota }
293bf3afd5fSEmilio G. Cota g_string_append(s, percent);
294bf3afd5fSEmilio G. Cota out:
295bf3afd5fSEmilio G. Cota return g_string_free(s, FALSE);
296bf3afd5fSEmilio G. Cota }
297bf3afd5fSEmilio G. Cota
298bf3afd5fSEmilio G. Cota /*
299bf3afd5fSEmilio G. Cota * Print the distribution's histogram into a string.
300bf3afd5fSEmilio G. Cota *
301bf3afd5fSEmilio G. Cota * See also: qdist_pr_plain().
302bf3afd5fSEmilio G. Cota *
303bf3afd5fSEmilio G. Cota * Callers must free the returned string with g_free().
304bf3afd5fSEmilio G. Cota */
qdist_pr(const struct qdist * dist,size_t n_bins,uint32_t opt)305bf3afd5fSEmilio G. Cota char *qdist_pr(const struct qdist *dist, size_t n_bins, uint32_t opt)
306bf3afd5fSEmilio G. Cota {
307bf3afd5fSEmilio G. Cota const char *border = opt & QDIST_PR_BORDER ? "|" : "";
308bf3afd5fSEmilio G. Cota char *llabel, *rlabel;
309bf3afd5fSEmilio G. Cota char *hgram;
310bf3afd5fSEmilio G. Cota GString *s;
311bf3afd5fSEmilio G. Cota
312bf3afd5fSEmilio G. Cota if (dist->n == 0) {
31311b7b07fSEmilio G. Cota return g_strdup(QDIST_EMPTY_STR);
314bf3afd5fSEmilio G. Cota }
315bf3afd5fSEmilio G. Cota
316bf3afd5fSEmilio G. Cota s = g_string_new("");
317bf3afd5fSEmilio G. Cota
318bf3afd5fSEmilio G. Cota llabel = qdist_pr_label(dist, n_bins, opt, true);
319bf3afd5fSEmilio G. Cota rlabel = qdist_pr_label(dist, n_bins, opt, false);
320bf3afd5fSEmilio G. Cota hgram = qdist_pr_plain(dist, n_bins);
321bf3afd5fSEmilio G. Cota g_string_append_printf(s, "%s%s%s%s%s",
322bf3afd5fSEmilio G. Cota llabel, border, hgram, border, rlabel);
323bf3afd5fSEmilio G. Cota g_free(llabel);
324bf3afd5fSEmilio G. Cota g_free(rlabel);
325bf3afd5fSEmilio G. Cota g_free(hgram);
326bf3afd5fSEmilio G. Cota
327bf3afd5fSEmilio G. Cota return g_string_free(s, FALSE);
328bf3afd5fSEmilio G. Cota }
329bf3afd5fSEmilio G. Cota
qdist_x(const struct qdist * dist,int index)330bf3afd5fSEmilio G. Cota static inline double qdist_x(const struct qdist *dist, int index)
331bf3afd5fSEmilio G. Cota {
332bf3afd5fSEmilio G. Cota if (dist->n == 0) {
333bf3afd5fSEmilio G. Cota return NAN;
334bf3afd5fSEmilio G. Cota }
335bf3afd5fSEmilio G. Cota return dist->entries[index].x;
336bf3afd5fSEmilio G. Cota }
337bf3afd5fSEmilio G. Cota
qdist_xmin(const struct qdist * dist)338bf3afd5fSEmilio G. Cota double qdist_xmin(const struct qdist *dist)
339bf3afd5fSEmilio G. Cota {
340bf3afd5fSEmilio G. Cota return qdist_x(dist, 0);
341bf3afd5fSEmilio G. Cota }
342bf3afd5fSEmilio G. Cota
qdist_xmax(const struct qdist * dist)343bf3afd5fSEmilio G. Cota double qdist_xmax(const struct qdist *dist)
344bf3afd5fSEmilio G. Cota {
345bf3afd5fSEmilio G. Cota return qdist_x(dist, dist->n - 1);
346bf3afd5fSEmilio G. Cota }
347bf3afd5fSEmilio G. Cota
qdist_unique_entries(const struct qdist * dist)348bf3afd5fSEmilio G. Cota size_t qdist_unique_entries(const struct qdist *dist)
349bf3afd5fSEmilio G. Cota {
350bf3afd5fSEmilio G. Cota return dist->n;
351bf3afd5fSEmilio G. Cota }
352bf3afd5fSEmilio G. Cota
qdist_sample_count(const struct qdist * dist)353bf3afd5fSEmilio G. Cota unsigned long qdist_sample_count(const struct qdist *dist)
354bf3afd5fSEmilio G. Cota {
355bf3afd5fSEmilio G. Cota unsigned long count = 0;
356bf3afd5fSEmilio G. Cota size_t i;
357bf3afd5fSEmilio G. Cota
358bf3afd5fSEmilio G. Cota for (i = 0; i < dist->n; i++) {
359bf3afd5fSEmilio G. Cota struct qdist_entry *e = &dist->entries[i];
360bf3afd5fSEmilio G. Cota
361bf3afd5fSEmilio G. Cota count += e->count;
362bf3afd5fSEmilio G. Cota }
363bf3afd5fSEmilio G. Cota return count;
364bf3afd5fSEmilio G. Cota }
365bf3afd5fSEmilio G. Cota
qdist_pairwise_avg(const struct qdist * dist,size_t index,size_t n,unsigned long count)366bf3afd5fSEmilio G. Cota static double qdist_pairwise_avg(const struct qdist *dist, size_t index,
367bf3afd5fSEmilio G. Cota size_t n, unsigned long count)
368bf3afd5fSEmilio G. Cota {
369bf3afd5fSEmilio G. Cota /* amortize the recursion by using a base case > 2 */
370bf3afd5fSEmilio G. Cota if (n <= 8) {
371bf3afd5fSEmilio G. Cota size_t i;
372bf3afd5fSEmilio G. Cota double ret = 0;
373bf3afd5fSEmilio G. Cota
374bf3afd5fSEmilio G. Cota for (i = 0; i < n; i++) {
375bf3afd5fSEmilio G. Cota struct qdist_entry *e = &dist->entries[index + i];
376bf3afd5fSEmilio G. Cota
377bf3afd5fSEmilio G. Cota ret += e->x * e->count / count;
378bf3afd5fSEmilio G. Cota }
379bf3afd5fSEmilio G. Cota return ret;
380bf3afd5fSEmilio G. Cota } else {
381bf3afd5fSEmilio G. Cota size_t n2 = n / 2;
382bf3afd5fSEmilio G. Cota
383bf3afd5fSEmilio G. Cota return qdist_pairwise_avg(dist, index, n2, count) +
384bf3afd5fSEmilio G. Cota qdist_pairwise_avg(dist, index + n2, n - n2, count);
385bf3afd5fSEmilio G. Cota }
386bf3afd5fSEmilio G. Cota }
387bf3afd5fSEmilio G. Cota
qdist_avg(const struct qdist * dist)388bf3afd5fSEmilio G. Cota double qdist_avg(const struct qdist *dist)
389bf3afd5fSEmilio G. Cota {
390bf3afd5fSEmilio G. Cota unsigned long count;
391bf3afd5fSEmilio G. Cota
392bf3afd5fSEmilio G. Cota count = qdist_sample_count(dist);
393bf3afd5fSEmilio G. Cota if (!count) {
394bf3afd5fSEmilio G. Cota return NAN;
395bf3afd5fSEmilio G. Cota }
396bf3afd5fSEmilio G. Cota return qdist_pairwise_avg(dist, 0, dist->n, count);
397bf3afd5fSEmilio G. Cota }
398