Actual source code: hists.c
petsc-3.4.5 2014-06-29
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: MPI_Comm comm;
57: PetscBool isnull;
63: PetscObjectGetComm((PetscObject) draw, &comm);
64: PetscHeaderCreate(h, _p_PetscDrawHG, int, PETSC_DRAWHG_CLASSID, "PetscDrawHG", "Histogram", "Draw", comm, PetscDrawHGDestroy, NULL);
66: h->view = NULL;
67: h->destroy = NULL;
68: h->win = draw;
70: PetscObjectReference((PetscObject) draw);
72: h->color = PETSC_DRAW_GREEN;
73: h->xmin = PETSC_MAX_REAL;
74: h->xmax = PETSC_MIN_REAL;
75: h->ymin = 0.;
76: h->ymax = 1.;
77: h->numBins = bins;
78: h->maxBins = bins;
80: PetscMalloc(h->maxBins * sizeof(PetscReal), &h->bins);
82: h->numValues = 0;
83: h->maxValues = CHUNKSIZE;
84: h->calcStats = PETSC_FALSE;
85: h->integerBins = PETSC_FALSE;
87: PetscMalloc(h->maxValues * sizeof(PetscReal), &h->values);
88: PetscLogObjectMemory(h, (h->maxBins + h->maxValues)*sizeof(PetscReal));
89: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
91: if (!isnull) {
92: PetscDrawAxisCreate(draw, &h->axis);
93: PetscLogObjectParent(h, h->axis);
94: } else h->axis = NULL;
95: *hist = h;
96: return(0);
97: }
101: /*@
102: PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn.
104: Not Collective (ignored except on processor 0 of PetscDrawHG)
106: Input Parameter:
107: + hist - The histogram context.
108: - dim - The number of curves.
110: Level: intermediate
112: Concepts: histogram^setting number of bins
114: @*/
115: PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
116: {
121: if (hist->maxBins < bins) {
122: PetscFree(hist->bins);
123: PetscMalloc(bins * sizeof(PetscReal), &hist->bins);
124: PetscLogObjectMemory(hist, (bins - hist->maxBins) * sizeof(PetscReal));
126: hist->maxBins = bins;
127: }
128: hist->numBins = bins;
129: return(0);
130: }
134: /*@
135: PetscDrawHGReset - Clears histogram to allow for reuse with new data.
137: Not Collective (ignored except on processor 0 of PetscDrawHG)
139: Input Parameter:
140: . hist - The histogram context.
142: Level: intermediate
144: Concepts: histogram^resetting
145: @*/
146: PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
147: {
150: hist->xmin = PETSC_MAX_REAL;
151: hist->xmax = PETSC_MIN_REAL;
152: hist->ymin = 0.0;
153: hist->ymax = 0.0;
154: hist->numValues = 0;
155: return(0);
156: }
160: /*@C
161: PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
163: Collective over PetscDrawHG
165: Input Parameter:
166: . hist - The histogram context
168: Level: intermediate
170: .seealso: PetscDrawHGCreate()
171: @*/
172: PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
173: {
177: if (!*hist) return(0);
180: if (--((PetscObject)(*hist))->refct > 0) return(0);
181: PetscDrawAxisDestroy(&(*hist)->axis);
182: PetscDrawDestroy(&(*hist)->win);
183: PetscFree((*hist)->bins);
184: PetscFree((*hist)->values);
185: PetscHeaderDestroy(hist);
186: return(0);
187: }
191: /*@
192: PetscDrawHGAddValue - Adds another value to the histogram.
194: Not Collective (ignored except on processor 0 of PetscDrawHG)
196: Input Parameters:
197: + hist - The histogram
198: - value - The value
200: Level: intermediate
202: Concepts: histogram^adding values
204: .seealso: PetscDrawHGAddValues()
205: @*/
206: PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
207: {
210: /* Allocate more memory if necessary */
211: if (hist->numValues >= hist->maxValues) {
212: PetscReal *tmp;
215: PetscMalloc((hist->maxValues+CHUNKSIZE) * sizeof(PetscReal), &tmp);
216: PetscLogObjectMemory(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: Not Collective (ignored except on processor 0 of PetscDrawHG)
259: Input Parameter:
260: . hist - The histogram context
262: Level: intermediate
264: @*/
265: PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
266: {
267: PetscDraw draw = hist->win;
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;
277: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
278: if (isnull) return(0);
279: if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
280: if (hist->numValues < 1) return(0);
282: #if 0
283: MPI_Comm_rank(PetscObjectComm((PetscObject)hist),&rank);
284: if (rank) return(0);
285: #endif
287: color = hist->color;
288: if (color == PETSC_DRAW_ROTATE) bcolor = 2;
289: else bcolor = color;
291: xmin = hist->xmin;
292: xmax = hist->xmax;
293: ymin = hist->ymin;
294: ymax = hist->ymax;
295: numValues = hist->numValues;
296: values = hist->values;
297: mean = 0.0;
298: var = 0.0;
300: PetscDrawClear(draw);
301: if (xmin == xmax) {
302: /* Calculate number of points in each bin */
303: bins = hist->bins;
304: bins[0] = 0.;
305: for (p = 0; p < numValues; p++) {
306: if (values[p] == xmin) bins[0]++;
307: mean += values[p];
308: var += values[p]*values[p];
309: }
310: maxHeight = bins[0];
311: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
312: xmax = xmin + 1;
313: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
314: if (hist->calcStats) {
315: mean /= numValues;
316: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
317: else var = 0.0;
318: PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var);
319: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
320: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
321: }
322: PetscDrawAxisDraw(hist->axis);
323: /* Draw bins */
324: binLeft = xmin;
325: binRight = xmax;
326: PetscDrawRectangle(draw,binLeft,ymin,binRight,bins[0],bcolor,bcolor,bcolor,bcolor);
328: if (color == PETSC_DRAW_ROTATE && bins[0] != 0.0) bcolor++;
329: if (bcolor > 31) bcolor = 2;
331: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[0],PETSC_DRAW_BLACK);
332: PetscDrawLine(draw,binRight,ymin,binRight,bins[0],PETSC_DRAW_BLACK);
333: PetscDrawLine(draw,binLeft,bins[0],binRight,bins[0],PETSC_DRAW_BLACK);
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;
373: PetscSNPrintf(title, 256,"Mean: %g Var: %g", (double)mean, (double)var);
374: PetscSNPrintf(xlabel,256, "Total: %D", numValues);
375: PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL);
376: }
377: PetscDrawAxisDraw(hist->axis);
378: /* 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: if (color == PETSC_DRAW_ROTATE && bins[i]) bcolor++;
384: if (bcolor > 31) bcolor = 2;
385: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[i],PETSC_DRAW_BLACK);
386: PetscDrawLine(draw,binRight,ymin,binRight,bins[i],PETSC_DRAW_BLACK);
387: PetscDrawLine(draw,binLeft,bins[i],binRight,bins[i],PETSC_DRAW_BLACK);
388: }
389: PetscDrawHGSetNumberBins(hist, numBinsOld);
390: }
391: PetscDrawSynchronizedFlush(draw);
392: PetscDrawPause(draw);
393: return(0);
394: }
398: /*@
399: PetscDrawHGView - Prints the histogram information.
401: Not collective
403: Input Parameter:
404: . hist - The histogram context
406: Level: beginner
408: .keywords: draw, histogram
409: @*/
410: PetscErrorCode PetscDrawHGView(PetscDrawHG hist,PetscViewer viewer)
411: {
412: PetscReal xmax,xmin,*bins,*values,binSize,binLeft,binRight,mean,var;
414: PetscInt numBins,numBinsOld,numValues,initSize,i,p;
418: if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
419: if (hist->numValues < 1) return(0);
421: if (!viewer){
422: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)hist),&viewer);
423: }
424: xmax = hist->xmax;
425: xmin = hist->xmin;
426: numValues = hist->numValues;
427: values = hist->values;
428: mean = 0.0;
429: var = 0.0;
430: if (xmax == xmin) {
431: /* Calculate number of points in the bin */
432: bins = hist->bins;
433: bins[0] = 0.;
434: for (p = 0; p < numValues; p++) {
435: if (values[p] == xmin) bins[0]++;
436: mean += values[p];
437: var += values[p]*values[p];
438: }
439: /* Draw bins */
440: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, (double)xmin, (double)xmax, (double)bins[0]);
441: } else {
442: numBins = hist->numBins;
443: numBinsOld = hist->numBins;
444: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
445: initSize = (int) ((int) xmax - xmin)/numBins;
446: while (initSize*numBins != (int) xmax - xmin) {
447: initSize = PetscMax(initSize - 1, 1);
448: numBins = (int) ((int) xmax - xmin)/initSize;
449: PetscDrawHGSetNumberBins(hist, numBins);
450: }
451: }
452: binSize = (xmax - xmin)/numBins;
453: bins = hist->bins;
455: /* Calculate number of points in each bin */
456: PetscMemzero(bins, numBins * sizeof(PetscReal));
457: for (i = 0; i < numBins; i++) {
458: binLeft = xmin + binSize*i;
459: binRight = xmin + binSize*(i+1);
460: for (p = 0; p < numValues; p++) {
461: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
462: /* Handle last bin separately */
463: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
464: if (!i) {
465: mean += values[p];
466: var += values[p]*values[p];
467: }
468: }
469: }
470: /* Draw bins */
471: for (i = 0; i < numBins; i++) {
472: binLeft = xmin + binSize*i;
473: binRight = xmin + binSize*(i+1);
474: PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", (int)i, (double)binLeft, (double)binRight, (double)bins[i]);
475: }
476: PetscDrawHGSetNumberBins(hist, numBinsOld);
477: }
479: if (hist->calcStats) {
480: mean /= numValues;
481: if (numValues > 1) var = (var - numValues*mean*mean) / (numValues-1);
482: else var = 0.0;
483: PetscViewerASCIIPrintf(viewer, "Mean: %G Var: %G\n", (double)mean, (double)var);
484: PetscViewerASCIIPrintf(viewer, "Total: %D\n", numValues);
485: }
486: return(0);
487: }
491: /*@
492: PetscDrawHGSetColor - Sets the color the bars will be drawn with.
494: Not Collective (ignored except on processor 0 of PetscDrawHG)
496: Input Parameters:
497: + hist - The histogram context
498: - color - one of the colors defined in petscdraw.h or PETSC_DRAW_ROTATE to make each bar a
499: different color
501: Level: intermediate
503: @*/
504: PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist, int color)
505: {
508: hist->color = color;
509: return(0);
510: }
514: /*@
515: PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
516: points are added after this call, the limits will be adjusted to
517: include those additional points.
519: Not Collective (ignored except on processor 0 of PetscDrawHG)
521: Input Parameters:
522: + hist - The histogram context
523: - x_min,x_max,y_min,y_max - The limits
525: Level: intermediate
527: Concepts: histogram^setting axis
528: @*/
529: PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
530: {
533: hist->xmin = x_min;
534: hist->xmax = x_max;
535: hist->ymin = y_min;
536: hist->ymax = y_max;
537: return(0);
538: }
542: /*@
543: PetscDrawHGCalcStats - Turns on calculation of descriptive statistics
545: Not collective
547: Input Parameters:
548: + hist - The histogram context
549: - calc - Flag for calculation
551: Level: intermediate
553: .keywords: draw, histogram, statistics
555: @*/
556: PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
557: {
560: hist->calcStats = calc;
561: return(0);
562: }
566: /*@
567: PetscDrawHGIntegerBins - Turns on integer width bins
569: Not collective
571: Input Parameters:
572: + hist - The histogram context
573: - ints - Flag for integer width bins
575: Level: intermediate
577: .keywords: draw, histogram, statistics
578: @*/
579: PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
580: {
583: hist->integerBins = ints;
584: return(0);
585: }
589: /*@C
590: PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
591: This is useful if one wants to change some axis property, such as
592: labels, color, etc. The axis context should not be destroyed by the
593: application code.
595: Not Collective (ignored except on processor 0 of PetscDrawHG)
597: Input Parameter:
598: . hist - The histogram context
600: Output Parameter:
601: . axis - The axis context
603: Level: intermediate
605: @*/
606: PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist, PetscDrawAxis *axis)
607: {
611: *axis = hist->axis;
612: return(0);
613: }
617: /*@C
618: PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
620: Not Collective, PetscDraw is parallel if PetscDrawHG is parallel
622: Input Parameter:
623: . hist - The histogram context
625: Output Parameter:
626: . win - The draw context
628: Level: intermediate
630: @*/
631: PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist, PetscDraw *win)
632: {
636: *win = hist->win;
637: return(0);
638: }