1/******************************************************************************
2 *
3 *   Copyright © International Business Machines  Corp., 2006, 2008
4 *
5 *   This program is free software;  you can redistribute it and/or modify
6 *   it under the terms of the GNU General Public License as published by
7 *   the Free Software Foundation; either version 2 of the License, or
8 *   (at your option) any later version.
9 *
10 *   This program is distributed in the hope that it will be useful,
11 *   but WITHOUT ANY WARRANTY;  without even the implied warranty of
12 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See
13 *   the GNU General Public License for more details.
14 *
15 *   You should have received a copy of the GNU General Public License
16 *   along with this program;  if not, write to the Free Software
17 *   Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 *
19 * NAME
20 *	   libstats.c
21 *
22 * DESCRIPTION
23 *	  Some basic statistical analysis convenience tools.
24 *
25 *
26 * USAGE:
27 *	  To be included in test cases
28 *
29 * AUTHOR
30 *		Darren Hart <dvhltc@us.ibm.com>
31 *
32 * HISTORY
33 *	  2006-Oct-17: Initial version by Darren Hart
34 *	  2009-Jul-22: Addition of stats_container_append function by Kiran Prakash
35 *
36 * TODO: the save routine for gnuplot plotting should be more modular...
37 *
38 *****************************************************************************/
39
40#include <stdio.h>
41#include <stdlib.h>
42#include <string.h>
43#include <errno.h>
44#include <unistd.h>
45#include <math.h>
46#include <libstats.h>
47#include <librttest.h>
48
49#include "../include/realtime_config.h"
50
51#ifndef HAVE_EXP10
52# define exp10(x) (exp((x) * log(10)))
53#endif
54
55int save_stats = 0;
56
57/* static helper functions */
58static int stats_record_compare(const void *a, const void *b)
59{
60	int ret = 0;
61	stats_record_t *rec_a = (stats_record_t *) a;
62	stats_record_t *rec_b = (stats_record_t *) b;
63	if (rec_a->y < rec_b->y)
64		ret = -1;
65	else if (rec_a->y > rec_b->y)
66		ret = 1;
67	return ret;
68}
69
70/* function implementations */
71int stats_container_init(stats_container_t * data, long size)
72{
73	data->size = size;
74	data->index = -1;
75	data->records = calloc(size, sizeof(stats_record_t));
76	if (!data->records)
77		return -1;
78	return 0;
79}
80
81int stats_container_append(stats_container_t * data, stats_record_t rec)
82{
83	int myindex = ++data->index;
84	if (myindex >= data->size) {
85		debug(DBG_ERR, "Number of elements cannot be more than %ld\n",
86		      data->size);
87		data->index--;
88		return -1;
89	}
90	data->records[myindex] = rec;
91	return myindex;
92}
93
94int stats_container_resize(stats_container_t * data, long size)
95{
96	stats_record_t *newrecords =
97	    realloc(data->records, size * sizeof(stats_record_t));
98	if (!newrecords)
99		return -1;
100	data->records = newrecords;
101	if (data->size < size)
102		memset(data->records + data->size, 0, size - data->size);
103	data->size = size;
104	return 0;
105}
106
107int stats_container_free(stats_container_t * data)
108{
109	free(data->records);
110	return 0;
111}
112
113int stats_sort(stats_container_t * data, enum stats_sort_method method)
114{
115	/* method not implemented, always ascending on y atm */
116	qsort(data->records, data->index + 1, sizeof(stats_record_t),
117	      stats_record_compare);
118	return 0;
119}
120
121float stats_stddev(stats_container_t * data)
122{
123	long i;
124	float sd, avg, sum, delta;
125	long n;
126
127	sd = 0.0;
128	n = data->index + 1;
129	sum = 0.0;
130
131	/* calculate the mean */
132	for (i = 0; i < n; i++) {
133		sum += data->records[i].y;
134	}
135	avg = sum / (float)n;
136
137	/* calculate the standard deviation */
138	sum = 0.0;
139	for (i = 0; i < n; i++) {
140		delta = (data->records[i].y - avg);
141		sum += delta * delta;
142	}
143	sum /= n;
144
145	sd = sqrt(sum);
146
147	return sd;
148}
149
150float stats_avg(stats_container_t * data)
151{
152	long i;
153	float avg, sum;
154	long n;
155
156	n = data->index + 1;
157	sum = 0.0;
158
159	/* calculate the mean */
160	for (i = 0; i < n; i++) {
161		sum += data->records[i].y;
162	}
163	avg = sum / (float)n;
164
165	return avg;
166}
167
168long stats_min(stats_container_t * data)
169{
170	long i;
171	long min;
172	long n;
173
174	n = data->index + 1;
175
176	/* calculate the mean */
177	min = data->records[0].y;
178	for (i = 1; i < n; i++) {
179		if (data->records[i].y < min) {
180			min = data->records[i].y;
181		}
182	}
183
184	return min;
185}
186
187long stats_max(stats_container_t * data)
188{
189	long i;
190	long max;
191	long n;
192
193	n = data->index + 1;
194
195	/* calculate the mean */
196	max = data->records[0].y;
197	for (i = 1; i < n; i++) {
198		if (data->records[i].y > max) {
199			max = data->records[i].y;
200		}
201	}
202
203	return max;
204}
205
206int stats_quantiles_init(stats_quantiles_t * quantiles, int nines)
207{
208	if (nines < 2) {
209		return -1;
210	}
211	quantiles->nines = nines;
212	/* allocate space for quantiles, starting with 0.99 (two nines) */
213	quantiles->quantiles = calloc(sizeof(long), (nines - 1));
214	if (!quantiles->quantiles) {
215		return -1;
216	}
217	return 0;
218}
219
220int stats_quantiles_free(stats_quantiles_t * quantiles)
221{
222	free(quantiles->quantiles);
223	return 0;
224}
225
226int stats_quantiles_calc(stats_container_t * data,
227			 stats_quantiles_t * quantiles)
228{
229	int i;
230	int size;
231	int index;
232
233	// check for sufficient data size of accurate calculation
234	if (data->index < 0 ||
235	    (data->index + 1) < (long)exp10(quantiles->nines)) {
236		return -1;
237	}
238
239	size = data->index + 1;
240	stats_sort(data, ASCENDING_ON_Y);
241
242	for (i = 2; i <= quantiles->nines; i++) {
243		index = size - size / exp10(i);
244		quantiles->quantiles[i - 2] = data->records[index].y;
245	}
246	return 0;
247}
248
249void stats_quantiles_print(stats_quantiles_t * quantiles)
250{
251	int i;
252	int fraction = 0;
253	for (i = 0; i <= quantiles->nines - 2; i++) {
254		if (i > 0)
255			fraction += 9 * exp10(i - 1);
256		printf("99.%d%% < %ld\n", fraction, quantiles->quantiles[i]);
257	}
258}
259
260int stats_hist(stats_container_t * hist, stats_container_t * data)
261{
262	int i;
263	int ret;
264	long min, max, width;
265	long y, b;
266
267	ret = 0;
268
269	if (hist->size <= 0 || data->index < 0) {
270		return -1;
271	}
272
273	/* calculate the range of dataset */
274	min = max = data->records[0].y;
275	for (i = 0; i <= data->index; i++) {
276		y = data->records[i].y;
277		if (y > max)
278			max = y;
279		if (y < min)
280			min = y;
281	}
282
283	/* define the bucket ranges */
284	width = MAX((max - min) / hist->size, 1);
285	hist->records[0].x = min;
286	for (i = 1; i < (hist->size); i++) {
287		hist->records[i].x = min + i * width;
288	}
289
290	/* fill in the counts */
291	for (i = 0; i <= data->index; i++) {
292		y = data->records[i].y;
293		b = MIN((y - min) / width, hist->size - 1);
294		hist->records[b].y++;
295	}
296
297	return 0;
298}
299
300void stats_hist_print(stats_container_t * hist)
301{
302	long i, x;
303	for (i = 0; i < hist->size; i++) {
304		x = hist->records[i].x;
305		if (i < hist->size - 1)
306			printf("[%ld,%ld) = %ld\n", x,
307			       hist->records[i + 1].x, hist->records[i].y);
308		else
309			printf("[%ld,-) = %ld\n", x, hist->records[i].y);
310	}
311}
312
313int stats_container_save(char *filename, char *title, char *xlabel,
314			 char *ylabel, stats_container_t * data, char *mode)
315{
316	int i;
317	int minx = 0, maxx = 0, miny = 0, maxy = 0;
318	FILE *dat_fd;
319	FILE *plt_fd;
320	char *datfile;
321	char *pltfile;
322	stats_record_t *rec;
323
324	if (!save_stats)
325		return 0;
326
327	/* generate the filenames */
328	if (asprintf(&datfile, "%s.dat", filename) == -1) {
329		fprintf(stderr,
330			"Failed to allocate string for data filename\n");
331		return -1;
332	}
333	if (asprintf(&pltfile, "%s.plt", filename) == -1) {
334		fprintf(stderr,
335			"Failed to allocate string for plot filename\n");
336		return -1;
337	}
338
339	/* generate the data file */
340	if ((dat_fd = fopen(datfile, "w")) == NULL) {
341		perror("Failed to open dat file");
342		return -1;
343	} else {
344		minx = maxx = data->records[0].x;
345		miny = maxy = data->records[0].y;
346		for (i = 0; i <= data->index; i++) {
347			rec = &data->records[i];
348			minx = MIN(minx, rec->x);
349			maxx = MAX(maxx, rec->x);
350			miny = MIN(miny, rec->y);
351			maxy = MAX(maxy, rec->y);
352			fprintf(dat_fd, "%ld %ld\n", rec->x, rec->y);
353		}
354		fclose(dat_fd);
355	}
356
357	/* generate the plt file */
358	if (!(plt_fd = fopen(pltfile, "w"))) {
359		perror("Failed to open plt file");
360		return -1;
361	} else {
362		fprintf(plt_fd, "set terminal png\n");
363		fprintf(plt_fd, "set output \"%s.png\"\n", pltfile);
364		fprintf(plt_fd, "set title \"%s\"\n", title);
365		fprintf(plt_fd, "set xlabel \"%s\"\n", xlabel);
366		fprintf(plt_fd, "set ylabel \"%s\"\n", ylabel);
367		fprintf(plt_fd, "plot [0:%d] [0:%d] \"%s\" with %s\n",
368			maxx + 1, maxy + 1, datfile, mode);
369		fclose(plt_fd);
370	}
371
372	return 0;
373}
374