Actual source code: hists.c
petsc-3.11.4 2019-09-28
2: /*
3: Contains the data structure for plotting a histogram in a window with an axis.
4: */
5: #include <petscdraw.h>
6: #include <petsc/private/petscimpl.h>
7: #include <petscviewer.h>
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
32: /*@C
33: PetscDrawHGCreate - Creates a histogram data structure.
35: Collective on PetscDraw
37: Input Parameters:
38: + draw - The window where the graph will be made
39: - bins - The number of bins to use
41: Output Parameters:
42: . hist - The histogram context
44: Notes:
45: The difference between a bar chart, PetscDrawBar, and a histogram, PetscDrawHG, is explained here http://stattrek.com/statistics/charts/histogram.aspx?Tutorial=AP
47: The histogram is only displayed when PetscDrawHGDraw() is called.
49: The MPI communicator that owns the PetscDraw owns this PetscDrawHG, but the calls to set options and add data are ignored on all processes except the
50: zeroth MPI process in the communicator. All MPI processes in the communicator must call PetscDrawHGDraw() to display the updated graph.
52: Level: intermediate
54: Concepts: histogram^creating
56: .seealso: PetscDrawHGDestroy(), PetscDrawHG, PetscDrawBarCreate(), PetscDrawBar, PetscDrawLGCreate(), PetscDrawLG, PetscDrawSPCreate(), PetscDrawSP,
57: PetscDrawHGSetNumberBins(), PetscDrawHGReset(), PetscDrawHGAddValue(), PetscDrawHGDraw(), PetscDrawHGSave(), PetscDrawHGView(), PetscDrawHGSetColor(),
58: PetscDrawHGSetLimits(), PetscDrawHGCalcStats(), PetscDrawHGIntegerBins(), PetscDrawHGGetAxis(), PetscDrawAxis, PetscDrawHGGetDraw()
60: @*/
61: PetscErrorCode PetscDrawHGCreate(PetscDraw draw,int bins,PetscDrawHG *hist)
62: {
63: PetscDrawHG h;
71: PetscHeaderCreate(h,PETSC_DRAWHG_CLASSID,"DrawHG","Histogram","Draw",PetscObjectComm((PetscObject)draw),PetscDrawHGDestroy,NULL);
72: PetscLogObjectParent((PetscObject)draw,(PetscObject)h);
74: PetscObjectReference((PetscObject)draw);
75: h->win = draw;
77: h->view = NULL;
78: h->destroy = NULL;
79: h->color = PETSC_DRAW_GREEN;
80: h->xmin = PETSC_MAX_REAL;
81: h->xmax = PETSC_MIN_REAL;
82: h->ymin = 0.;
83: h->ymax = 1.;
84: h->numBins = bins;
85: h->maxBins = bins;
87: PetscMalloc1(h->maxBins,&h->bins);
89: h->numValues = 0;
90: h->maxValues = CHUNKSIZE;
91: h->calcStats = PETSC_FALSE;
92: h->integerBins = PETSC_FALSE;
94: PetscMalloc1(h->maxValues,&h->values);
95: PetscLogObjectMemory((PetscObject)h,(h->maxBins + h->maxValues)*sizeof(PetscReal));
97: PetscDrawAxisCreate(draw,&h->axis);
98: PetscLogObjectParent((PetscObject)h,(PetscObject)h->axis);
100: *hist = h;
101: return(0);
102: }
104: /*@
105: PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn.
107: Logically Collective on PetscDrawHG
109: Input Parameter:
110: + hist - The histogram context.
111: - bins - The number of bins.
113: Level: intermediate
115: Concepts: histogram^setting number of bins
117: .seealso: PetscDrawHGCreate(), PetscDrawHG, PetscDrawHGDraw(), PetscDrawHGIntegerBins()
119: @*/
120: PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
121: {
128: if (hist->maxBins < bins) {
129: PetscFree(hist->bins);
130: PetscMalloc1(bins, &hist->bins);
131: PetscLogObjectMemory((PetscObject)hist, (bins - hist->maxBins) * sizeof(PetscReal));
132: hist->maxBins = bins;
133: }
134: hist->numBins = bins;
135: return(0);
136: }
138: /*@
139: PetscDrawHGReset - Clears histogram to allow for reuse with new data.
141: Logically Collective on PetscDrawHG
143: Input Parameter:
144: . hist - The histogram context.
146: Level: intermediate
148: Concepts: histogram^resetting
150: .seealso: PetscDrawHGCreate(), PetscDrawHG, PetscDrawHGDraw(), PetscDrawHGAddValue()
152: @*/
153: PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
154: {
158: hist->xmin = PETSC_MAX_REAL;
159: hist->xmax = PETSC_MIN_REAL;
160: hist->ymin = 0.0;
161: hist->ymax = 0.0;
162: hist->numValues = 0;
163: return(0);
164: }
166: /*@C
167: PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
169: Collective on PetscDrawHG
171: Input Parameter:
172: . hist - The histogram context
174: Level: intermediate
176: .seealso: PetscDrawHGCreate(), PetscDrawHG
177: @*/
178: PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
179: {
183: if (!*hist) return(0);
185: if (--((PetscObject)(*hist))->refct > 0) {*hist = NULL; return(0);}
187: PetscFree((*hist)->bins);
188: PetscFree((*hist)->values);
189: PetscDrawAxisDestroy(&(*hist)->axis);
190: PetscDrawDestroy(&(*hist)->win);
191: PetscHeaderDestroy(hist);
192: return(0);
193: }
195: /*@
196: PetscDrawHGAddValue - Adds another value to the histogram.
198: Logically Collective on PetscDrawHG
200: Input Parameters:
201: + hist - The histogram
202: - value - The value
204: Level: intermediate
206: Concepts: histogram^adding values
208: .seealso: PetscDrawHGCreate(), PetscDrawHG, PetscDrawHGDraw(), PetscDrawHGAddValue(), PetscDrawHGReset()
209: @*/
210: PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
211: {
215: /* Allocate more memory if necessary */
216: if (hist->numValues >= hist->maxValues) {
217: PetscReal *tmp;
220: PetscMalloc1(hist->maxValues+CHUNKSIZE, &tmp);
221: PetscLogObjectMemory((PetscObject)hist, CHUNKSIZE * sizeof(PetscReal));
222: PetscMemcpy(tmp, hist->values, hist->maxValues * sizeof(PetscReal));
223: PetscFree(hist->values);
225: hist->values = tmp;
226: hist->maxValues += CHUNKSIZE;
227: }
228: /* I disagree with the original Petsc implementation here. There should be no overshoot, but rather the
229: stated convention of using half-open intervals (always the way to go) */
230: if (!hist->numValues) {
231: hist->xmin = value;
232: hist->xmax = value;
233: #if 1
234: } else {
235: /* Update limits */
236: if (value > hist->xmax) hist->xmax = value;
237: if (value < hist->xmin) hist->xmin = value;
238: #else
239: } else if (hist->numValues == 1) {
240: /* Update limits -- We need to overshoot the largest value somewhat */
241: if (value > hist->xmax) hist->xmax = value + 0.001*(value - hist->xmin)/hist->numBins;
242: if (value < hist->xmin) {
243: hist->xmin = value;
244: hist->xmax = hist->xmax + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
245: }
246: } else {
247: /* Update limits -- We need to overshoot the largest value somewhat */
248: if (value > hist->xmax) hist->xmax = value + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
249: if (value < hist->xmin) hist->xmin = value;
250: #endif
251: }
253: hist->values[hist->numValues++] = value;
254: return(0);
255: }
257: /*@
258: PetscDrawHGDraw - Redraws a histogram.
260: Collective on PetscDrawHG
262: Input Parameter:
263: . hist - The histogram context
265: Level: intermediate
267: .seealso: PetscDrawHGCreate(), PetscDrawHG, PetscDrawHGDraw(), PetscDrawHGAddValue(), PetscDrawHGReset()
269: @*/
270: PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
271: {
272: PetscDraw draw;
273: PetscBool isnull;
274: PetscReal xmin,xmax,ymin,ymax,*bins,*values,binSize,binLeft,binRight,maxHeight,mean,var;
275: char title[256];
276: char xlabel[256];
277: PetscInt numBins,numBinsOld,numValues,initSize,i,p,bcolor,color;
278: PetscMPIInt rank;
283: PetscDrawIsNull(hist->win,&isnull);
284: if (isnull) return(0);
285: MPI_Comm_rank(PetscObjectComm((PetscObject)hist),&rank);
287: if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
288: if (hist->numValues < 1) return(0);
290: color = hist->color;
291: if (color == PETSC_DRAW_ROTATE) bcolor = PETSC_DRAW_BLACK+1;
292: else bcolor = color;
294: xmin = hist->xmin;
295: xmax = hist->xmax;
296: ymin = hist->ymin;
297: ymax = hist->ymax;
298: numValues = hist->numValues;
299: values = hist->values;
300: mean = 0.0;
301: var = 0.0;
303: draw = hist->win;
304: PetscDrawCheckResizedWindow(draw);
305: PetscDrawClear(draw);
307: if (xmin == xmax) {
308: /* Calculate number of points in each bin */
309: bins = hist->bins;
310: bins[0] = 0.;
311: for (p = 0; p < numValues; p++) {
312: if (values[p] == xmin) bins[0]++;
313: mean += values[p];
314: var += values[p]*values[p];
315: }
316: maxHeight = bins[0];
317: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
318: xmax = xmin + 1;
319: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
320: if (hist->calcStats) {
321: mean /= numValues;
322: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
323: else var = 0.0;
324: PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var);
325: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
326: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
327: }
328: PetscDrawAxisDraw(hist->axis);
329: PetscDrawCollectiveBegin(draw);
330: if (!rank) { /* Draw bins */
331: binLeft = xmin;
332: binRight = xmax;
333: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[0],bcolor,bcolor,bcolor,bcolor);
334: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[0],PETSC_DRAW_BLACK);
335: PetscDrawLine(draw,binRight,ymin,binRight,bins[0],PETSC_DRAW_BLACK);
336: PetscDrawLine(draw,binLeft,bins[0],binRight,bins[0],PETSC_DRAW_BLACK);
337: }
338: PetscDrawCollectiveEnd(draw);
339: } else {
340: numBins = hist->numBins;
341: numBinsOld = hist->numBins;
342: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
343: initSize = (int) ((int) xmax - xmin)/numBins;
344: while (initSize*numBins != (int) xmax - xmin) {
345: initSize = PetscMax(initSize - 1, 1);
346: numBins = (int) ((int) xmax - xmin)/initSize;
347: PetscDrawHGSetNumberBins(hist, numBins);
348: }
349: }
350: binSize = (xmax - xmin)/numBins;
351: bins = hist->bins;
353: PetscMemzero(bins, numBins * sizeof(PetscReal));
355: maxHeight = 0.0;
356: for (i = 0; i < numBins; i++) {
357: binLeft = xmin + binSize*i;
358: binRight = xmin + binSize*(i+1);
359: for (p = 0; p < numValues; p++) {
360: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
361: /* Handle last bin separately */
362: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
363: if (!i) {
364: mean += values[p];
365: var += values[p]*values[p];
366: }
367: }
368: maxHeight = PetscMax(maxHeight, bins[i]);
369: }
370: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
372: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
373: if (hist->calcStats) {
374: mean /= numValues;
375: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
376: else var = 0.0;
377: PetscSNPrintf(title, 256,"Mean: %g Var: %g", (double)mean, (double)var);
378: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
379: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
380: }
381: PetscDrawAxisDraw(hist->axis);
382: PetscDrawCollectiveBegin(draw);
383: if (!rank) { /* Draw bins */
384: for (i = 0; i < numBins; i++) {
385: binLeft = xmin + binSize*i;
386: binRight = xmin + binSize*(i+1);
387: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[i],bcolor,bcolor,bcolor,bcolor);
388: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[i],PETSC_DRAW_BLACK);
389: PetscDrawLine(draw,binRight,ymin,binRight,bins[i],PETSC_DRAW_BLACK);
390: PetscDrawLine(draw,binLeft,bins[i],binRight,bins[i],PETSC_DRAW_BLACK);
391: if (color == PETSC_DRAW_ROTATE && bins[i]) bcolor++;
392: if (bcolor > PETSC_DRAW_BASIC_COLORS-1) bcolor = PETSC_DRAW_BLACK+1;
393: }
394: }
395: PetscDrawCollectiveEnd(draw);
396: PetscDrawHGSetNumberBins(hist,numBinsOld);
397: }
399: PetscDrawFlush(draw);
400: PetscDrawPause(draw);
401: return(0);
402: }
404: /*@
405: PetscDrawHGSave - Saves a drawn image
407: Collective on PetscDrawHG
409: Input Parameter:
410: . hist - The histogram context
412: Level: intermediate
414: Concepts: histogram^saving
416: .seealso: PetscDrawHGCreate(), PetscDrawHGGetDraw(), PetscDrawSetSave(), PetscDrawSave(), PetscDrawHGDraw()
417: @*/
418: PetscErrorCode PetscDrawHGSave(PetscDrawHG hg)
419: {
424: PetscDrawSave(hg->win);
425: return(0);
426: }
428: /*@
429: PetscDrawHGView - Prints the histogram information.
431: Not collective
433: Input Parameter:
434: . hist - The histogram context
436: Level: beginner
438: .seealso: PetscDrawHGCreate(), PetscDrawHGGetDraw(), PetscDrawSetSave(), PetscDrawSave(), PetscDrawHGDraw()
440: .keywords: draw, histogram
441: @*/
442: PetscErrorCode PetscDrawHGView(PetscDrawHG hist,PetscViewer viewer)
443: {
444: PetscReal xmax,xmin,*bins,*values,binSize,binLeft,binRight,mean,var;
446: PetscInt numBins,numBinsOld,numValues,initSize,i,p;
451: if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
452: if (hist->numValues < 1) return(0);
454: if (!viewer){
455: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)hist),&viewer);
456: }
457: PetscObjectPrintClassNamePrefixType((PetscObject)hist,viewer);
458: xmax = hist->xmax;
459: xmin = hist->xmin;
460: numValues = hist->numValues;
461: values = hist->values;
462: mean = 0.0;
463: var = 0.0;
464: if (xmax == xmin) {
465: /* Calculate number of points in the bin */
466: bins = hist->bins;
467: bins[0] = 0.;
468: for (p = 0; p < numValues; p++) {
469: if (values[p] == xmin) bins[0]++;
470: mean += values[p];
471: var += values[p]*values[p];
472: }
473: /* Draw bins */
474: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, (double)xmin, (double)xmax, (double)bins[0]);
475: } else {
476: numBins = hist->numBins;
477: numBinsOld = hist->numBins;
478: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
479: initSize = (int) ((int) xmax - xmin)/numBins;
480: while (initSize*numBins != (int) xmax - xmin) {
481: initSize = PetscMax(initSize - 1, 1);
482: numBins = (int) ((int) xmax - xmin)/initSize;
483: PetscDrawHGSetNumberBins(hist, numBins);
484: }
485: }
486: binSize = (xmax - xmin)/numBins;
487: bins = hist->bins;
489: /* Calculate number of points in each bin */
490: PetscMemzero(bins, numBins * sizeof(PetscReal));
491: for (i = 0; i < numBins; i++) {
492: binLeft = xmin + binSize*i;
493: binRight = xmin + binSize*(i+1);
494: for (p = 0; p < numValues; p++) {
495: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
496: /* Handle last bin separately */
497: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
498: if (!i) {
499: mean += values[p];
500: var += values[p]*values[p];
501: }
502: }
503: }
504: /* Draw bins */
505: for (i = 0; i < numBins; i++) {
506: binLeft = xmin + binSize*i;
507: binRight = xmin + binSize*(i+1);
508: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", (int)i, (double)binLeft, (double)binRight, (double)bins[i]);
509: }
510: PetscDrawHGSetNumberBins(hist, numBinsOld);
511: }
513: if (hist->calcStats) {
514: mean /= numValues;
515: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
516: else var = 0.0;
517: PetscViewerASCIIPrintf(viewer, "Mean: %g Var: %g\n", (double)mean, (double)var);
518: PetscViewerASCIIPrintf(viewer, "Total: %D\n", numValues);
519: }
520: return(0);
521: }
523: /*@
524: PetscDrawHGSetColor - Sets the color the bars will be drawn with.
526: Logically Collective on PetscDrawHG
528: Input Parameters:
529: + hist - The histogram context
530: - color - one of the colors defined in petscdraw.h or PETSC_DRAW_ROTATE to make each bar a
531: different color
533: Level: intermediate
535: .seealso: PetscDrawHGCreate(), PetscDrawHGGetDraw(), PetscDrawSetSave(), PetscDrawSave(), PetscDrawHGDraw(), PetscDrawHGGetAxis()
537: @*/
538: PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist,int color)
539: {
543: hist->color = color;
544: return(0);
545: }
547: /*@
548: PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
549: points are added after this call, the limits will be adjusted to
550: include those additional points.
552: Logically Collective on PetscDrawHG
554: Input Parameters:
555: + hist - The histogram context
556: - x_min,x_max,y_min,y_max - The limits
558: Level: intermediate
560: Concepts: histogram^setting axis
562: .seealso: PetscDrawHGCreate(), PetscDrawHGGetDraw(), PetscDrawSetSave(), PetscDrawSave(), PetscDrawHGDraw(), PetscDrawHGGetAxis()
564: @*/
565: PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
566: {
570: hist->xmin = x_min;
571: hist->xmax = x_max;
572: hist->ymin = y_min;
573: hist->ymax = y_max;
574: return(0);
575: }
577: /*@
578: PetscDrawHGCalcStats - Turns on calculation of descriptive statistics
580: Not collective
582: Input Parameters:
583: + hist - The histogram context
584: - calc - Flag for calculation
586: Level: intermediate
588: .keywords: draw, histogram, statistics
590: .seealso: PetscDrawHGCreate(), PetscDrawHGAddValue(), PetscDrawHGView(), PetscDrawHGDraw()
592: @*/
593: PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
594: {
598: hist->calcStats = calc;
599: return(0);
600: }
602: /*@
603: PetscDrawHGIntegerBins - Turns on integer width bins
605: Not collective
607: Input Parameters:
608: + hist - The histogram context
609: - ints - Flag for integer width bins
611: Level: intermediate
613: .keywords: draw, histogram, statistics
615: .seealso: PetscDrawHGCreate(), PetscDrawHGAddValue(), PetscDrawHGView(), PetscDrawHGDraw(), PetscDrawHGSetColor()
617: @*/
618: PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
619: {
623: hist->integerBins = ints;
624: return(0);
625: }
627: /*@C
628: PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
629: This is useful if one wants to change some axis property, such as
630: labels, color, etc. The axis context should not be destroyed by the
631: application code.
633: Not Collective, PetscDrawAxis is parallel if PetscDrawHG is parallel
635: Input Parameter:
636: . hist - The histogram context
638: Output Parameter:
639: . axis - The axis context
641: Level: intermediate
643: .seealso: PetscDrawHGCreate(), PetscDrawHGAddValue(), PetscDrawHGView(), PetscDrawHGDraw(), PetscDrawHGSetColor(), PetscDrawAxis, PetscDrawHGSetLimits()
645: @*/
646: PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist,PetscDrawAxis *axis)
647: {
651: *axis = hist->axis;
652: return(0);
653: }
655: /*@C
656: PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
658: Not Collective, PetscDraw is parallel if PetscDrawHG is parallel
660: Input Parameter:
661: . hist - The histogram context
663: Output Parameter:
664: . draw - The draw context
666: Level: intermediate
668: .seealso: PetscDrawHGCreate(), PetscDrawHGAddValue(), PetscDrawHGView(), PetscDrawHGDraw(), PetscDrawHGSetColor(), PetscDrawAxis, PetscDrawHGSetLimits()
670: @*/
671: PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist,PetscDraw *draw)
672: {
676: *draw = hist->win;
677: return(0);
678: }