Actual source code: hists.c
petsc-3.7.7 2017-09-25
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 on 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;
63: PetscHeaderCreate(h,PETSC_DRAWHG_CLASSID,"DrawHG","Histogram","Draw",PetscObjectComm((PetscObject)draw),PetscDrawHGDestroy,NULL);
64: PetscLogObjectParent((PetscObject)draw,(PetscObject)h);
66: PetscObjectReference((PetscObject)draw);
67: h->win = draw;
69: h->view = NULL;
70: h->destroy = NULL;
71: h->color = PETSC_DRAW_GREEN;
72: h->xmin = PETSC_MAX_REAL;
73: h->xmax = PETSC_MIN_REAL;
74: h->ymin = 0.;
75: h->ymax = 1.;
76: h->numBins = bins;
77: h->maxBins = bins;
79: PetscMalloc1(h->maxBins,&h->bins);
81: h->numValues = 0;
82: h->maxValues = CHUNKSIZE;
83: h->calcStats = PETSC_FALSE;
84: h->integerBins = PETSC_FALSE;
86: PetscMalloc1(h->maxValues,&h->values);
87: PetscLogObjectMemory((PetscObject)h,(h->maxBins + h->maxValues)*sizeof(PetscReal));
89: PetscDrawAxisCreate(draw,&h->axis);
90: PetscLogObjectParent((PetscObject)h,(PetscObject)h->axis);
92: *hist = h;
93: return(0);
94: }
98: /*@
99: PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn.
101: Logically Collective on PetscDrawHG
103: Input Parameter:
104: + hist - The histogram context.
105: - bins - The number of bins.
107: Level: intermediate
109: Concepts: histogram^setting number of bins
111: @*/
112: PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
113: {
120: if (hist->maxBins < bins) {
121: PetscFree(hist->bins);
122: PetscMalloc1(bins, &hist->bins);
123: 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: Logically Collective on PetscDrawHG
137: Input Parameter:
138: . hist - The histogram context.
140: Level: intermediate
142: Concepts: histogram^resetting
143: @*/
144: PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
145: {
149: hist->xmin = PETSC_MAX_REAL;
150: hist->xmax = PETSC_MIN_REAL;
151: hist->ymin = 0.0;
152: hist->ymax = 0.0;
153: hist->numValues = 0;
154: return(0);
155: }
159: /*@C
160: PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
162: Collective on PetscDrawHG
164: Input Parameter:
165: . hist - The histogram context
167: Level: intermediate
169: .seealso: PetscDrawHGCreate()
170: @*/
171: PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
172: {
176: if (!*hist) return(0);
178: if (--((PetscObject)(*hist))->refct > 0) {*hist = NULL; return(0);}
180: PetscFree((*hist)->bins);
181: PetscFree((*hist)->values);
182: PetscDrawAxisDestroy(&(*hist)->axis);
183: PetscDrawDestroy(&(*hist)->win);
184: PetscHeaderDestroy(hist);
185: return(0);
186: }
190: /*@
191: PetscDrawHGAddValue - Adds another value to the histogram.
193: Logically Collective on PetscDrawHG
195: Input Parameters:
196: + hist - The histogram
197: - value - The value
199: Level: intermediate
201: Concepts: histogram^adding values
203: .seealso: PetscDrawHGAddValues()
204: @*/
205: PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
206: {
210: /* Allocate more memory if necessary */
211: if (hist->numValues >= hist->maxValues) {
212: PetscReal *tmp;
215: PetscMalloc1(hist->maxValues+CHUNKSIZE, &tmp);
216: PetscLogObjectMemory((PetscObject)hist, CHUNKSIZE * sizeof(PetscReal));
217: PetscMemcpy(tmp, hist->values, hist->maxValues * sizeof(PetscReal));
218: PetscFree(hist->values);
220: hist->values = tmp;
221: hist->maxValues += CHUNKSIZE;
222: }
223: /* I disagree with the original Petsc implementation here. There should be no overshoot, but rather the
224: stated convention of using half-open intervals (always the way to go) */
225: if (!hist->numValues) {
226: hist->xmin = value;
227: hist->xmax = value;
228: #if 1
229: } else {
230: /* Update limits */
231: if (value > hist->xmax) hist->xmax = value;
232: if (value < hist->xmin) hist->xmin = value;
233: #else
234: } else if (hist->numValues == 1) {
235: /* Update limits -- We need to overshoot the largest value somewhat */
236: if (value > hist->xmax) hist->xmax = value + 0.001*(value - hist->xmin)/hist->numBins;
237: if (value < hist->xmin) {
238: hist->xmin = value;
239: hist->xmax = hist->xmax + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
240: }
241: } else {
242: /* Update limits -- We need to overshoot the largest value somewhat */
243: if (value > hist->xmax) hist->xmax = value + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
244: if (value < hist->xmin) hist->xmin = value;
245: #endif
246: }
248: hist->values[hist->numValues++] = value;
249: return(0);
250: }
254: /*@
255: PetscDrawHGDraw - Redraws a histogram.
257: Collective on PetscDrawHG
259: Input Parameter:
260: . hist - The histogram context
262: Level: intermediate
264: @*/
265: PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
266: {
267: PetscDraw draw;
268: PetscBool isnull;
269: PetscReal xmin,xmax,ymin,ymax,*bins,*values,binSize,binLeft,binRight,maxHeight,mean,var;
270: char title[256];
271: char xlabel[256];
272: PetscInt numBins,numBinsOld,numValues,initSize,i,p,bcolor,color;
273: PetscMPIInt rank;
278: PetscDrawIsNull(hist->win,&isnull);
279: if (isnull) return(0);
280: MPI_Comm_rank(PetscObjectComm((PetscObject)hist),&rank);
282: if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
283: if (hist->numValues < 1) return(0);
285: color = hist->color;
286: if (color == PETSC_DRAW_ROTATE) bcolor = PETSC_DRAW_BLACK+1;
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: draw = hist->win;
299: PetscDrawCheckResizedWindow(draw);
300: PetscDrawClear(draw);
302: if (xmin == xmax) {
303: /* Calculate number of points in each bin */
304: bins = hist->bins;
305: bins[0] = 0.;
306: for (p = 0; p < numValues; p++) {
307: if (values[p] == xmin) bins[0]++;
308: mean += values[p];
309: var += values[p]*values[p];
310: }
311: maxHeight = bins[0];
312: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
313: xmax = xmin + 1;
314: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
315: if (hist->calcStats) {
316: mean /= numValues;
317: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
318: else var = 0.0;
319: PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var);
320: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
321: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
322: }
323: PetscDrawAxisDraw(hist->axis);
324: PetscDrawCollectiveBegin(draw);
325: if (!rank) { /* Draw bins */
326: binLeft = xmin;
327: binRight = xmax;
328: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[0],bcolor,bcolor,bcolor,bcolor);
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: }
333: PetscDrawCollectiveEnd(draw);
334: } else {
335: numBins = hist->numBins;
336: numBinsOld = hist->numBins;
337: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
338: initSize = (int) ((int) xmax - xmin)/numBins;
339: while (initSize*numBins != (int) xmax - xmin) {
340: initSize = PetscMax(initSize - 1, 1);
341: numBins = (int) ((int) xmax - xmin)/initSize;
342: PetscDrawHGSetNumberBins(hist, numBins);
343: }
344: }
345: binSize = (xmax - xmin)/numBins;
346: bins = hist->bins;
348: PetscMemzero(bins, numBins * sizeof(PetscReal));
350: maxHeight = 0.0;
351: for (i = 0; i < numBins; i++) {
352: binLeft = xmin + binSize*i;
353: binRight = xmin + binSize*(i+1);
354: for (p = 0; p < numValues; p++) {
355: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
356: /* Handle last bin separately */
357: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
358: if (!i) {
359: mean += values[p];
360: var += values[p]*values[p];
361: }
362: }
363: maxHeight = PetscMax(maxHeight, bins[i]);
364: }
365: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
367: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
368: if (hist->calcStats) {
369: mean /= numValues;
370: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
371: else var = 0.0;
372: PetscSNPrintf(title, 256,"Mean: %g Var: %g", (double)mean, (double)var);
373: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
374: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
375: }
376: PetscDrawAxisDraw(hist->axis);
377: PetscDrawCollectiveBegin(draw);
378: if (!rank) { /* Draw bins */
379: for (i = 0; i < numBins; i++) {
380: binLeft = xmin + binSize*i;
381: binRight = xmin + binSize*(i+1);
382: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[i],bcolor,bcolor,bcolor,bcolor);
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: if (color == PETSC_DRAW_ROTATE && bins[i]) bcolor++;
387: if (bcolor > PETSC_DRAW_BASIC_COLORS-1) bcolor = PETSC_DRAW_BLACK+1;
388: }
389: }
390: PetscDrawCollectiveEnd(draw);
391: PetscDrawHGSetNumberBins(hist,numBinsOld);
392: }
394: PetscDrawFlush(draw);
395: PetscDrawPause(draw);
396: return(0);
397: }
401: /*@
402: PetscDrawHGSave - Saves a drawn image
404: Collective on PetscDrawHG
406: Input Parameter:
407: . hist - The histogram context
409: Level: intermediate
411: Concepts: histogram^saving
413: .seealso: PetscDrawHGCreate(), PetscDrawHGGetDraw(), PetscDrawSetSave(), PetscDrawSave()
414: @*/
415: PetscErrorCode PetscDrawHGSave(PetscDrawHG hg)
416: {
421: PetscDrawSave(hg->win);
422: return(0);
423: }
427: /*@
428: PetscDrawHGView - Prints the histogram information.
430: Not collective
432: Input Parameter:
433: . hist - The histogram context
435: Level: beginner
437: .keywords: draw, histogram
438: @*/
439: PetscErrorCode PetscDrawHGView(PetscDrawHG hist,PetscViewer viewer)
440: {
441: PetscReal xmax,xmin,*bins,*values,binSize,binLeft,binRight,mean,var;
443: PetscInt numBins,numBinsOld,numValues,initSize,i,p;
448: if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
449: if (hist->numValues < 1) return(0);
451: if (!viewer){
452: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)hist),&viewer);
453: }
454: PetscObjectPrintClassNamePrefixType((PetscObject)hist,viewer);
455: xmax = hist->xmax;
456: xmin = hist->xmin;
457: numValues = hist->numValues;
458: values = hist->values;
459: mean = 0.0;
460: var = 0.0;
461: if (xmax == xmin) {
462: /* Calculate number of points in the bin */
463: bins = hist->bins;
464: bins[0] = 0.;
465: for (p = 0; p < numValues; p++) {
466: if (values[p] == xmin) bins[0]++;
467: mean += values[p];
468: var += values[p]*values[p];
469: }
470: /* Draw bins */
471: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, (double)xmin, (double)xmax, (double)bins[0]);
472: } else {
473: numBins = hist->numBins;
474: numBinsOld = hist->numBins;
475: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
476: initSize = (int) ((int) xmax - xmin)/numBins;
477: while (initSize*numBins != (int) xmax - xmin) {
478: initSize = PetscMax(initSize - 1, 1);
479: numBins = (int) ((int) xmax - xmin)/initSize;
480: PetscDrawHGSetNumberBins(hist, numBins);
481: }
482: }
483: binSize = (xmax - xmin)/numBins;
484: bins = hist->bins;
486: /* Calculate number of points in each bin */
487: PetscMemzero(bins, numBins * sizeof(PetscReal));
488: for (i = 0; i < numBins; i++) {
489: binLeft = xmin + binSize*i;
490: binRight = xmin + binSize*(i+1);
491: for (p = 0; p < numValues; p++) {
492: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
493: /* Handle last bin separately */
494: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
495: if (!i) {
496: mean += values[p];
497: var += values[p]*values[p];
498: }
499: }
500: }
501: /* Draw bins */
502: for (i = 0; i < numBins; i++) {
503: binLeft = xmin + binSize*i;
504: binRight = xmin + binSize*(i+1);
505: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", (int)i, (double)binLeft, (double)binRight, (double)bins[i]);
506: }
507: PetscDrawHGSetNumberBins(hist, numBinsOld);
508: }
510: if (hist->calcStats) {
511: mean /= numValues;
512: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
513: else var = 0.0;
514: PetscViewerASCIIPrintf(viewer, "Mean: %g Var: %g\n", (double)mean, (double)var);
515: PetscViewerASCIIPrintf(viewer, "Total: %D\n", numValues);
516: }
517: return(0);
518: }
522: /*@
523: PetscDrawHGSetColor - Sets the color the bars will be drawn with.
525: Logically Collective on PetscDrawHG
527: Input Parameters:
528: + hist - The histogram context
529: - color - one of the colors defined in petscdraw.h or PETSC_DRAW_ROTATE to make each bar a
530: different color
532: Level: intermediate
534: @*/
535: PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist,int color)
536: {
540: hist->color = color;
541: return(0);
542: }
546: /*@
547: PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
548: points are added after this call, the limits will be adjusted to
549: include those additional points.
551: Logically Collective on PetscDrawHG
553: Input Parameters:
554: + hist - The histogram context
555: - x_min,x_max,y_min,y_max - The limits
557: Level: intermediate
559: Concepts: histogram^setting axis
560: @*/
561: PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
562: {
566: hist->xmin = x_min;
567: hist->xmax = x_max;
568: hist->ymin = y_min;
569: hist->ymax = y_max;
570: return(0);
571: }
575: /*@
576: PetscDrawHGCalcStats - Turns on calculation of descriptive statistics
578: Not collective
580: Input Parameters:
581: + hist - The histogram context
582: - calc - Flag for calculation
584: Level: intermediate
586: .keywords: draw, histogram, statistics
588: @*/
589: PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
590: {
594: hist->calcStats = calc;
595: return(0);
596: }
600: /*@
601: PetscDrawHGIntegerBins - Turns on integer width bins
603: Not collective
605: Input Parameters:
606: + hist - The histogram context
607: - ints - Flag for integer width bins
609: Level: intermediate
611: .keywords: draw, histogram, statistics
612: @*/
613: PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
614: {
618: hist->integerBins = ints;
619: return(0);
620: }
624: /*@C
625: PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
626: This is useful if one wants to change some axis property, such as
627: labels, color, etc. The axis context should not be destroyed by the
628: application code.
630: Not Collective, PetscDrawAxis is parallel if PetscDrawHG is parallel
632: Input Parameter:
633: . hist - The histogram context
635: Output Parameter:
636: . axis - The axis context
638: Level: intermediate
640: @*/
641: PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist,PetscDrawAxis *axis)
642: {
646: *axis = hist->axis;
647: return(0);
648: }
652: /*@C
653: PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
655: Not Collective, PetscDraw is parallel if PetscDrawHG is parallel
657: Input Parameter:
658: . hist - The histogram context
660: Output Parameter:
661: . draw - The draw context
663: Level: intermediate
665: @*/
666: PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist,PetscDraw *draw)
667: {
671: *draw = hist->win;
672: return(0);
673: }