Actual source code: hists.c
petsc-3.6.1 2015-08-06
2: /*
3: Contains the data structure for plotting a histogram in a window with an axis.
4: */
5: #include <petscdraw.h> /*I "petscdraw.h" I*/
6: #include <petsc/private/petscimpl.h> /*I "petscsys.h" I*/
7: #include <petscviewer.h> /*I "petscviewer.h" I*/
9: PetscClassId PETSC_DRAWHG_CLASSID = 0;
11: struct _p_PetscDrawHG {
12: PETSCHEADER(int);
13: PetscErrorCode (*destroy)(PetscDrawSP);
14: PetscErrorCode (*view)(PetscDrawSP,PetscViewer);
15: PetscDraw win;
16: PetscDrawAxis axis;
17: PetscReal xmin,xmax;
18: PetscReal ymin,ymax;
19: int numBins;
20: int maxBins;
21: PetscReal *bins;
22: int numValues;
23: int maxValues;
24: PetscReal *values;
25: int color;
26: PetscBool calcStats;
27: PetscBool integerBins;
28: };
30: #define CHUNKSIZE 100
34: /*@C
35: PetscDrawHGCreate - Creates a histogram data structure.
37: Collective over PetscDraw
39: Input Parameters:
40: + draw - The window where the graph will be made
41: - bins - The number of bins to use
43: Output Parameters:
44: . hist - The histogram context
46: Level: intermediate
48: Concepts: histogram^creating
50: .seealso: PetscDrawHGDestroy()
52: @*/
53: PetscErrorCode PetscDrawHGCreate(PetscDraw draw, int bins, PetscDrawHG *hist)
54: {
55: PetscDrawHG h;
56: PetscBool isnull;
62: PetscHeaderCreate(h, PETSC_DRAWHG_CLASSID, "PetscDrawHG", "Histogram", "Draw", PetscObjectComm((PetscObject)draw), PetscDrawHGDestroy, NULL);
64: h->view = NULL;
65: h->destroy = NULL;
66: h->win = draw;
68: PetscObjectReference((PetscObject) draw);
70: h->color = PETSC_DRAW_GREEN;
71: h->xmin = PETSC_MAX_REAL;
72: h->xmax = PETSC_MIN_REAL;
73: h->ymin = 0.;
74: h->ymax = 1.;
75: h->numBins = bins;
76: h->maxBins = bins;
78: PetscMalloc1(h->maxBins, &h->bins);
80: h->numValues = 0;
81: h->maxValues = CHUNKSIZE;
82: h->calcStats = PETSC_FALSE;
83: h->integerBins = PETSC_FALSE;
85: PetscMalloc1(h->maxValues, &h->values);
86: PetscLogObjectMemory((PetscObject)h, (h->maxBins + h->maxValues)*sizeof(PetscReal));
87: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
89: if (!isnull) {
90: PetscDrawAxisCreate(draw, &h->axis);
91: PetscLogObjectParent((PetscObject)h, (PetscObject)h->axis);
92: } else h->axis = NULL;
93: *hist = h;
94: return(0);
95: }
99: /*@
100: PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn.
102: Not Collective (ignored except on processor 0 of PetscDrawHG)
104: Input Parameter:
105: + hist - The histogram context.
106: - dim - The number of curves.
108: Level: intermediate
110: Concepts: histogram^setting number of bins
112: @*/
113: PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
114: {
119: if (hist->maxBins < bins) {
120: PetscFree(hist->bins);
121: PetscMalloc1(bins, &hist->bins);
122: PetscLogObjectMemory((PetscObject)hist, (bins - hist->maxBins) * sizeof(PetscReal));
124: hist->maxBins = bins;
125: }
126: hist->numBins = bins;
127: return(0);
128: }
132: /*@
133: PetscDrawHGReset - Clears histogram to allow for reuse with new data.
135: Not Collective (ignored except on processor 0 of PetscDrawHG)
137: Input Parameter:
138: . hist - The histogram context.
140: Level: intermediate
142: Concepts: histogram^resetting
143: @*/
144: PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
145: {
148: hist->xmin = PETSC_MAX_REAL;
149: hist->xmax = PETSC_MIN_REAL;
150: hist->ymin = 0.0;
151: hist->ymax = 0.0;
152: hist->numValues = 0;
153: return(0);
154: }
158: /*@C
159: PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
161: Collective over PetscDrawHG
163: Input Parameter:
164: . hist - The histogram context
166: Level: intermediate
168: .seealso: PetscDrawHGCreate()
169: @*/
170: PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
171: {
175: if (!*hist) return(0);
178: if (--((PetscObject)(*hist))->refct > 0) return(0);
179: PetscDrawAxisDestroy(&(*hist)->axis);
180: PetscDrawDestroy(&(*hist)->win);
181: PetscFree((*hist)->bins);
182: PetscFree((*hist)->values);
183: PetscHeaderDestroy(hist);
184: return(0);
185: }
189: /*@
190: PetscDrawHGAddValue - Adds another value to the histogram.
192: Not Collective (ignored except on processor 0 of PetscDrawHG)
194: Input Parameters:
195: + hist - The histogram
196: - value - The value
198: Level: intermediate
200: Concepts: histogram^adding values
202: .seealso: PetscDrawHGAddValues()
203: @*/
204: PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
205: {
208: /* Allocate more memory if necessary */
209: if (hist->numValues >= hist->maxValues) {
210: PetscReal *tmp;
213: PetscMalloc1(hist->maxValues+CHUNKSIZE, &tmp);
214: PetscLogObjectMemory((PetscObject)hist, CHUNKSIZE * sizeof(PetscReal));
215: PetscMemcpy(tmp, hist->values, hist->maxValues * sizeof(PetscReal));
216: PetscFree(hist->values);
218: hist->values = tmp;
219: hist->maxValues += CHUNKSIZE;
220: }
221: /* I disagree with the original Petsc implementation here. There should be no overshoot, but rather the
222: stated convention of using half-open intervals (always the way to go) */
223: if (!hist->numValues) {
224: hist->xmin = value;
225: hist->xmax = value;
226: #if 1
227: } else {
228: /* Update limits */
229: if (value > hist->xmax) hist->xmax = value;
230: if (value < hist->xmin) hist->xmin = value;
231: #else
232: } else if (hist->numValues == 1) {
233: /* Update limits -- We need to overshoot the largest value somewhat */
234: if (value > hist->xmax) hist->xmax = value + 0.001*(value - hist->xmin)/hist->numBins;
235: if (value < hist->xmin) {
236: hist->xmin = value;
237: hist->xmax = hist->xmax + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
238: }
239: } else {
240: /* Update limits -- We need to overshoot the largest value somewhat */
241: if (value > hist->xmax) hist->xmax = value + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
242: if (value < hist->xmin) hist->xmin = value;
243: #endif
244: }
246: hist->values[hist->numValues++] = value;
247: return(0);
248: }
252: /*@
253: PetscDrawHGDraw - Redraws a histogram.
255: Not Collective (ignored except on processor 0 of PetscDrawHG)
257: Input Parameter:
258: . hist - The histogram context
260: Level: intermediate
262: @*/
263: PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
264: {
265: PetscDraw draw = hist->win;
266: PetscBool isnull;
267: PetscReal xmin,xmax,ymin,ymax,*bins,*values,binSize,binLeft,binRight,maxHeight,mean,var;
268: char title[256];
269: char xlabel[256];
270: PetscInt numBins,numBinsOld,numValues,initSize,i,p,bcolor,color;
275: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
276: if (isnull) return(0);
277: if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
278: if (hist->numValues < 1) return(0);
280: #if 0
281: MPI_Comm_rank(PetscObjectComm((PetscObject)hist),&rank);
282: if (rank) return(0);
283: #endif
285: color = hist->color;
286: if (color == PETSC_DRAW_ROTATE) bcolor = 2;
287: else bcolor = color;
289: xmin = hist->xmin;
290: xmax = hist->xmax;
291: ymin = hist->ymin;
292: ymax = hist->ymax;
293: numValues = hist->numValues;
294: values = hist->values;
295: mean = 0.0;
296: var = 0.0;
298: PetscDrawClear(draw);
299: if (xmin == xmax) {
300: /* Calculate number of points in each bin */
301: bins = hist->bins;
302: bins[0] = 0.;
303: for (p = 0; p < numValues; p++) {
304: if (values[p] == xmin) bins[0]++;
305: mean += values[p];
306: var += values[p]*values[p];
307: }
308: maxHeight = bins[0];
309: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
310: xmax = xmin + 1;
311: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
312: if (hist->calcStats) {
313: mean /= numValues;
314: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
315: else var = 0.0;
316: PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var);
317: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
318: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
319: }
320: PetscDrawAxisDraw(hist->axis);
321: /* Draw bins */
322: binLeft = xmin;
323: binRight = xmax;
324: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[0],bcolor,bcolor,bcolor,bcolor);
326: if (color == PETSC_DRAW_ROTATE && bins[0] != 0.0) bcolor++;
327: if (bcolor > 31) bcolor = 2;
329: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[0],PETSC_DRAW_BLACK);
330: PetscDrawLine(draw,binRight,ymin,binRight,bins[0],PETSC_DRAW_BLACK);
331: PetscDrawLine(draw,binLeft,bins[0],binRight,bins[0],PETSC_DRAW_BLACK);
332: } else {
333: numBins = hist->numBins;
334: numBinsOld = hist->numBins;
335: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
336: initSize = (int) ((int) xmax - xmin)/numBins;
337: while (initSize*numBins != (int) xmax - xmin) {
338: initSize = PetscMax(initSize - 1, 1);
339: numBins = (int) ((int) xmax - xmin)/initSize;
340: PetscDrawHGSetNumberBins(hist, numBins);
341: }
342: }
343: binSize = (xmax - xmin)/numBins;
344: bins = hist->bins;
346: PetscMemzero(bins, numBins * sizeof(PetscReal));
348: maxHeight = 0.0;
349: for (i = 0; i < numBins; i++) {
350: binLeft = xmin + binSize*i;
351: binRight = xmin + binSize*(i+1);
352: for (p = 0; p < numValues; p++) {
353: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
354: /* Handle last bin separately */
355: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
356: if (!i) {
357: mean += values[p];
358: var += values[p]*values[p];
359: }
360: }
361: maxHeight = PetscMax(maxHeight, bins[i]);
362: }
363: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
365: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
366: if (hist->calcStats) {
367: mean /= numValues;
368: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
369: else var = 0.0;
371: PetscSNPrintf(title, 256,"Mean: %g Var: %g", (double)mean, (double)var);
372: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
373: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
374: }
375: PetscDrawAxisDraw(hist->axis);
376: /* Draw bins */
377: for (i = 0; i < numBins; i++) {
378: binLeft = xmin + binSize*i;
379: binRight = xmin + binSize*(i+1);
380: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[i],bcolor,bcolor,bcolor,bcolor);
381: if (color == PETSC_DRAW_ROTATE && bins[i]) bcolor++;
382: if (bcolor > 31) bcolor = 2;
383: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[i],PETSC_DRAW_BLACK);
384: PetscDrawLine(draw,binRight,ymin,binRight,bins[i],PETSC_DRAW_BLACK);
385: PetscDrawLine(draw,binLeft,bins[i],binRight,bins[i],PETSC_DRAW_BLACK);
386: }
387: PetscDrawHGSetNumberBins(hist, numBinsOld);
388: }
389: PetscDrawSynchronizedFlush(draw);
390: PetscDrawPause(draw);
391: return(0);
392: }
396: /*@
397: PetscDrawHGView - Prints the histogram information.
399: Not collective
401: Input Parameter:
402: . hist - The histogram context
404: Level: beginner
406: .keywords: draw, histogram
407: @*/
408: PetscErrorCode PetscDrawHGView(PetscDrawHG hist,PetscViewer viewer)
409: {
410: PetscReal xmax,xmin,*bins,*values,binSize,binLeft,binRight,mean,var;
412: PetscInt numBins,numBinsOld,numValues,initSize,i,p;
416: if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
417: if (hist->numValues < 1) return(0);
419: if (!viewer){
420: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)hist),&viewer);
421: }
422: PetscObjectPrintClassNamePrefixType((PetscObject)hist,viewer);
423: xmax = hist->xmax;
424: xmin = hist->xmin;
425: numValues = hist->numValues;
426: values = hist->values;
427: mean = 0.0;
428: var = 0.0;
429: if (xmax == xmin) {
430: /* Calculate number of points in the bin */
431: bins = hist->bins;
432: bins[0] = 0.;
433: for (p = 0; p < numValues; p++) {
434: if (values[p] == xmin) bins[0]++;
435: mean += values[p];
436: var += values[p]*values[p];
437: }
438: /* Draw bins */
439: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, (double)xmin, (double)xmax, (double)bins[0]);
440: } else {
441: numBins = hist->numBins;
442: numBinsOld = hist->numBins;
443: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
444: initSize = (int) ((int) xmax - xmin)/numBins;
445: while (initSize*numBins != (int) xmax - xmin) {
446: initSize = PetscMax(initSize - 1, 1);
447: numBins = (int) ((int) xmax - xmin)/initSize;
448: PetscDrawHGSetNumberBins(hist, numBins);
449: }
450: }
451: binSize = (xmax - xmin)/numBins;
452: bins = hist->bins;
454: /* Calculate number of points in each bin */
455: PetscMemzero(bins, numBins * sizeof(PetscReal));
456: for (i = 0; i < numBins; i++) {
457: binLeft = xmin + binSize*i;
458: binRight = xmin + binSize*(i+1);
459: for (p = 0; p < numValues; p++) {
460: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
461: /* Handle last bin separately */
462: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
463: if (!i) {
464: mean += values[p];
465: var += values[p]*values[p];
466: }
467: }
468: }
469: /* Draw bins */
470: for (i = 0; i < numBins; i++) {
471: binLeft = xmin + binSize*i;
472: binRight = xmin + binSize*(i+1);
473: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", (int)i, (double)binLeft, (double)binRight, (double)bins[i]);
474: }
475: PetscDrawHGSetNumberBins(hist, numBinsOld);
476: }
478: if (hist->calcStats) {
479: mean /= numValues;
480: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
481: else var = 0.0;
482: PetscViewerASCIIPrintf(viewer, "Mean: %g Var: %g\n", (double)mean, (double)var);
483: PetscViewerASCIIPrintf(viewer, "Total: %D\n", numValues);
484: }
485: return(0);
486: }
490: /*@
491: PetscDrawHGSetColor - Sets the color the bars will be drawn with.
493: Not Collective (ignored except on processor 0 of PetscDrawHG)
495: Input Parameters:
496: + hist - The histogram context
497: - color - one of the colors defined in petscdraw.h or PETSC_DRAW_ROTATE to make each bar a
498: different color
500: Level: intermediate
502: @*/
503: PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist, int color)
504: {
507: hist->color = color;
508: return(0);
509: }
513: /*@
514: PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
515: points are added after this call, the limits will be adjusted to
516: include those additional points.
518: Not Collective (ignored except on processor 0 of PetscDrawHG)
520: Input Parameters:
521: + hist - The histogram context
522: - x_min,x_max,y_min,y_max - The limits
524: Level: intermediate
526: Concepts: histogram^setting axis
527: @*/
528: PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
529: {
532: hist->xmin = x_min;
533: hist->xmax = x_max;
534: hist->ymin = y_min;
535: hist->ymax = y_max;
536: return(0);
537: }
541: /*@
542: PetscDrawHGCalcStats - Turns on calculation of descriptive statistics
544: Not collective
546: Input Parameters:
547: + hist - The histogram context
548: - calc - Flag for calculation
550: Level: intermediate
552: .keywords: draw, histogram, statistics
554: @*/
555: PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
556: {
559: hist->calcStats = calc;
560: return(0);
561: }
565: /*@
566: PetscDrawHGIntegerBins - Turns on integer width bins
568: Not collective
570: Input Parameters:
571: + hist - The histogram context
572: - ints - Flag for integer width bins
574: Level: intermediate
576: .keywords: draw, histogram, statistics
577: @*/
578: PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
579: {
582: hist->integerBins = ints;
583: return(0);
584: }
588: /*@C
589: PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
590: This is useful if one wants to change some axis property, such as
591: labels, color, etc. The axis context should not be destroyed by the
592: application code.
594: Not Collective (ignored except on processor 0 of PetscDrawHG)
596: Input Parameter:
597: . hist - The histogram context
599: Output Parameter:
600: . axis - The axis context
602: Level: intermediate
604: @*/
605: PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist, PetscDrawAxis *axis)
606: {
610: *axis = hist->axis;
611: return(0);
612: }
616: /*@C
617: PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
619: Not Collective, PetscDraw is parallel if PetscDrawHG is parallel
621: Input Parameter:
622: . hist - The histogram context
624: Output Parameter:
625: . win - The draw context
627: Level: intermediate
629: @*/
630: PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist, PetscDraw *win)
631: {
635: *win = hist->win;
636: return(0);
637: }