Actual source code: hists.c
petsc-3.3-p7 2013-05-11
2: /*
3: Contains the data structure for plotting a histogram in a window with an axis.
4: */
6: #include <petscsys.h> /*I "petscsys.h" I*/
8: PetscClassId PETSC_DRAWHG_CLASSID = 0;
10: struct _p_PetscDrawHG {
11: PETSCHEADER(int);
12: PetscErrorCode (*destroy)(PetscDrawSP);
13: PetscErrorCode (*view)(PetscDrawSP,PetscViewer);
14: PetscDraw win;
15: PetscDrawAxis axis;
16: PetscReal xmin,xmax;
17: PetscReal ymin,ymax;
18: int numBins;
19: int maxBins;
20: PetscReal *bins;
21: int numValues;
22: int maxValues;
23: PetscReal *values;
24: int color;
25: PetscBool calcStats;
26: PetscBool integerBins;
27: };
29: #define CHUNKSIZE 100
33: /*@C
34: PetscDrawHGCreate - Creates a histogram data structure.
36: Collective over PetscDraw
38: Input Parameters:
39: + draw - The window where the graph will be made
40: - bins - The number of bins to use
42: Output Parameters:
43: . hist - The histogram context
45: Level: intermediate
47: Concepts: histogram^creating
49: .seealso: PetscDrawHGDestroy()
51: @*/
52: PetscErrorCode PetscDrawHGCreate(PetscDraw draw, int bins, PetscDrawHG *hist) {
53: PetscDrawHG h;
54: MPI_Comm comm;
55: PetscBool isnull;
61: PetscObjectGetComm((PetscObject) draw, &comm);
62: PetscHeaderCreate(h, _p_PetscDrawHG, int, PETSC_DRAWHG_CLASSID, 0, "PetscDrawHG", "Histogram", "Draw", comm, PetscDrawHGDestroy, PETSC_NULL);
63: h->view = PETSC_NULL;
64: h->destroy = PETSC_NULL;
65: h->win = draw;
66: PetscObjectReference((PetscObject) draw);
67: h->color = PETSC_DRAW_GREEN;
68: h->xmin = PETSC_MAX_REAL;
69: h->xmax = PETSC_MIN_REAL;
70: h->ymin = 0.;
71: h->ymax = 1.;
72: h->numBins = bins;
73: h->maxBins = bins;
74: PetscMalloc(h->maxBins * sizeof(PetscReal), &h->bins);
75: h->numValues = 0;
76: h->maxValues = CHUNKSIZE;
77: h->calcStats = PETSC_FALSE;
78: h->integerBins = PETSC_FALSE;
79: PetscMalloc(h->maxValues * sizeof(PetscReal), &h->values);
80: PetscLogObjectMemory(h, (h->maxBins + h->maxValues)*sizeof(PetscReal));
81: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
82: if (!isnull) {
83: PetscDrawAxisCreate(draw, &h->axis);
84: PetscLogObjectParent(h, h->axis);
85: } else {
86: h->axis = PETSC_NULL;
87: }
88: *hist = h;
89: return(0);
90: }
94: /*@
95: PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn.
97: Not Collective (ignored except on processor 0 of PetscDrawHG)
99: Input Parameter:
100: + hist - The histogram context.
101: - dim - The number of curves.
103: Level: intermediate
105: Concepts: histogram^setting number of bins
107: @*/
108: PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
109: {
114: if (hist->maxBins < bins) {
115: PetscFree(hist->bins);
116: PetscMalloc(bins * sizeof(PetscReal), &hist->bins);
117: PetscLogObjectMemory(hist, (bins - hist->maxBins) * sizeof(PetscReal));
118: hist->maxBins = bins;
119: }
120: hist->numBins = bins;
121: return(0);
122: }
126: /*@
127: PetscDrawHGReset - Clears histogram to allow for reuse with new data.
129: Not Collective (ignored except on processor 0 of PetscDrawHG)
131: Input Parameter:
132: . hist - The histogram context.
134: Level: intermediate
136: Concepts: histogram^resetting
137: @*/
138: PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
139: {
142: hist->xmin = PETSC_MAX_REAL;
143: hist->xmax = PETSC_MIN_REAL;
144: hist->ymin = 0.0;
145: hist->ymax = 0.0;
146: hist->numValues = 0;
147: return(0);
148: }
152: /*@C
153: PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
155: Collective over PetscDrawHG
157: Input Parameter:
158: . hist - The histogram context
160: Level: intermediate
162: .seealso: PetscDrawHGCreate()
163: @*/
164: PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
165: {
169: if (!*hist) return(0);
172: if (--((PetscObject)(*hist))->refct > 0) return(0);
173: PetscDrawAxisDestroy(&(*hist)->axis);
174: PetscDrawDestroy(&(*hist)->win);
175: PetscFree((*hist)->bins);
176: PetscFree((*hist)->values);
177: PetscHeaderDestroy(hist);
178: return(0);
179: }
183: /*@
184: PetscDrawHGAddValue - Adds another value to the histogram.
186: Not Collective (ignored except on processor 0 of PetscDrawHG)
188: Input Parameters:
189: + hist - The histogram
190: - value - The value
192: Level: intermediate
194: Concepts: histogram^adding values
196: .seealso: PetscDrawHGAddValues()
197: @*/
198: PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
199: {
202: /* Allocate more memory if necessary */
203: if (hist->numValues >= hist->maxValues) {
204: PetscReal *tmp;
207: PetscMalloc((hist->maxValues+CHUNKSIZE) * sizeof(PetscReal), &tmp);
208: PetscLogObjectMemory(hist, CHUNKSIZE * sizeof(PetscReal));
209: PetscMemcpy(tmp, hist->values, hist->maxValues * sizeof(PetscReal));
210: PetscFree(hist->values);
211: hist->values = tmp;
212: hist->maxValues += CHUNKSIZE;
213: }
214: /* I disagree with the original Petsc implementation here. There should be no overshoot, but rather the
215: stated convention of using half-open intervals (always the way to go) */
216: if (!hist->numValues) {
217: hist->xmin = value;
218: hist->xmax = value;
219: #if 1
220: } else {
221: /* Update limits */
222: if (value > hist->xmax)
223: hist->xmax = value;
224: if (value < hist->xmin)
225: hist->xmin = value;
226: #else
227: } else if (hist->numValues == 1) {
228: /* Update limits -- We need to overshoot the largest value somewhat */
229: if (value > hist->xmax) {
230: hist->xmax = value + 0.001*(value - hist->xmin)/hist->numBins;
231: }
232: if (value < hist->xmin) {
233: hist->xmin = value;
234: hist->xmax = hist->xmax + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
235: }
236: } else {
237: /* Update limits -- We need to overshoot the largest value somewhat */
238: if (value > hist->xmax) {
239: hist->xmax = value + 0.001*(hist->xmax - hist->xmin)/hist->numBins;
240: }
241: if (value < hist->xmin) {
242: hist->xmin = value;
243: }
244: #endif
245: }
247: hist->values[hist->numValues++] = value;
248: return(0);
249: }
253: /*@
254: PetscDrawHGDraw - Redraws a histogram.
256: Not Collective (ignored except on processor 0 of PetscDrawHG)
258: Input Parameter:
259: . hist - The histogram context
261: Level: intermediate
263: @*/
264: PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
265: {
266: PetscDraw draw = hist->win;
267: PetscBool isnull;
268: PetscReal xmin,xmax,ymin,ymax,*bins,*values,binSize,binLeft,binRight,maxHeight,mean,var;
269: char title[256];
270: char xlabel[256];
271: PetscInt numBins,numBinsOld,numValues,initSize,i,p,bcolor,color;
276: PetscObjectTypeCompare((PetscObject) draw, PETSC_DRAW_NULL, &isnull);
277: if (isnull) return(0);
278: if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
279: if (hist->numValues < 1) return(0);
281: #if 0
282: MPI_Comm_rank(((PetscObject)hist)->comm,&rank);
283: if (rank) return(0);
284: #endif
286: color = hist->color;
287: if (color == PETSC_DRAW_ROTATE) {bcolor = 2;} else {bcolor = color;}
288: xmin = hist->xmin;
289: xmax = hist->xmax;
290: ymin = hist->ymin;
291: ymax = hist->ymax;
292: numValues = hist->numValues;
293: values = hist->values;
294: mean = 0.0;
295: var = 0.0;
296:
297: PetscDrawClear(draw);
298: if (xmin == xmax) {
299: /* Calculate number of points in each bin */
300: bins = hist->bins;
301: bins[0] = 0.;
302: for(p = 0; p < numValues; p++) {
303: if (values[p] == xmin) bins[0]++;
304: mean += values[p];
305: var += values[p]*values[p];
306: }
307: maxHeight = bins[0];
308: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
309: xmax = xmin + 1;
310: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
311: if (hist->calcStats) {
312: mean /= numValues;
313: if (numValues > 1) {
314: var = (var - numValues*mean*mean) / (numValues-1);
315: } else {
316: var = 0.0;
317: }
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, PETSC_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);
327: if (color == PETSC_DRAW_ROTATE && bins[0] != 0.0) bcolor++; if (bcolor > 31) bcolor = 2;
328: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[0],PETSC_DRAW_BLACK);
329: PetscDrawLine(draw,binRight,ymin,binRight,bins[0],PETSC_DRAW_BLACK);
330: PetscDrawLine(draw,binLeft,bins[0],binRight,bins[0],PETSC_DRAW_BLACK);
331: } else {
332: numBins = hist->numBins;
333: numBinsOld = hist->numBins;
334: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
335: initSize = (int) ((int) xmax - xmin)/numBins;
336: while (initSize*numBins != (int) xmax - xmin) {
337: initSize = PetscMax(initSize - 1, 1);
338: numBins = (int) ((int) xmax - xmin)/initSize;
339: PetscDrawHGSetNumberBins(hist, numBins);
340: }
341: }
342: binSize = (xmax - xmin)/numBins;
343: bins = hist->bins;
345: PetscMemzero(bins, numBins * sizeof(PetscReal));
346: maxHeight = 0.0;
347: for (i = 0; i < numBins; i++) {
348: binLeft = xmin + binSize*i;
349: binRight = xmin + binSize*(i+1);
350: for(p = 0; p < numValues; p++) {
351: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
352: /* Handle last bin separately */
353: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
354: if (!i) {
355: mean += values[p];
356: var += values[p]*values[p];
357: }
358: }
359: maxHeight = PetscMax(maxHeight, bins[i]);
360: }
361: if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
363: PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax);
364: if (hist->calcStats) {
365: mean /= numValues;
366: if (numValues > 1) {
367: var = (var - numValues*mean*mean) / (numValues-1);
368: } else {
369: var = 0.0;
370: }
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, PETSC_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++; if (bcolor > 31) bcolor = 2;
382: PetscDrawLine(draw,binLeft,ymin,binLeft,bins[i],PETSC_DRAW_BLACK);
383: PetscDrawLine(draw,binRight,ymin,binRight,bins[i],PETSC_DRAW_BLACK);
384: PetscDrawLine(draw,binLeft,bins[i],binRight,bins[i],PETSC_DRAW_BLACK);
385: }
386: PetscDrawHGSetNumberBins(hist, numBinsOld);
387: }
388: PetscDrawSynchronizedFlush(draw);
389: PetscDrawPause(draw);
390: return(0);
391: }
395: /*@
396: PetscDrawHGPrint - Prints the histogram information.
398: Not collective
400: Input Parameter:
401: . hist - The histogram context
403: Level: beginner
405: .keywords: draw, histogram
406: @*/
407: PetscErrorCode PetscDrawHGPrint(PetscDrawHG hist)
408: {
409: PetscReal xmax,xmin,*bins,*values,binSize,binLeft,binRight,mean,var;
411: PetscInt numBins,numBinsOld,numValues,initSize,i,p;
415: if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) return(0);
416: if (hist->numValues < 1) return(0);
418: xmax = hist->xmax;
419: xmin = hist->xmin;
420: numValues = hist->numValues;
421: values = hist->values;
422: mean = 0.0;
423: var = 0.0;
424: if (xmax == xmin) {
425: /* Calculate number of points in the bin */
426: bins = hist->bins;
427: bins[0] = 0.;
428: for(p = 0; p < numValues; p++) {
429: if (values[p] == xmin) bins[0]++;
430: mean += values[p];
431: var += values[p]*values[p];
432: }
433: /* Draw bins */
434: PetscPrintf(((PetscObject)hist)->comm, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, xmin, xmax, bins[0]);
435: } else {
436: numBins = hist->numBins;
437: numBinsOld = hist->numBins;
438: if (hist->integerBins && (((int) xmax - xmin) + 1.0e-05 > xmax - xmin)) {
439: initSize = (int) ((int) xmax - xmin)/numBins;
440: while (initSize*numBins != (int) xmax - xmin) {
441: initSize = PetscMax(initSize - 1, 1);
442: numBins = (int) ((int) xmax - xmin)/initSize;
443: PetscDrawHGSetNumberBins(hist, numBins);
444: }
445: }
446: binSize = (xmax - xmin)/numBins;
447: bins = hist->bins;
449: /* Calculate number of points in each bin */
450: PetscMemzero(bins, numBins * sizeof(PetscReal));
451: for (i = 0; i < numBins; i++) {
452: binLeft = xmin + binSize*i;
453: binRight = xmin + binSize*(i+1);
454: for(p = 0; p < numValues; p++) {
455: if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i]++;
456: /* Handle last bin separately */
457: if ((i == numBins-1) && (values[p] == binRight)) bins[i]++;
458: if (!i) {
459: mean += values[p];
460: var += values[p]*values[p];
461: }
462: }
463: }
464: /* Draw bins */
465: for (i = 0; i < numBins; i++) {
466: binLeft = xmin + binSize*i;
467: binRight = xmin + binSize*(i+1);
468: PetscPrintf(((PetscObject)hist)->comm, "Bin %2d (%6.2g - %6.2g): %.0g\n", i, binLeft, binRight, bins[i]);
469: }
470: PetscDrawHGSetNumberBins(hist, numBinsOld);
471: }
473: if (hist->calcStats) {
474: mean /= numValues;
475: if (numValues > 1) {
476: var = (var - numValues*mean*mean) / (numValues-1);
477: } else {
478: var = 0.0;
479: }
480: PetscPrintf(((PetscObject)hist)->comm, "Mean: %G Var: %G\n", mean, var);
481: PetscPrintf(((PetscObject)hist)->comm, "Total: %d\n", numValues);
482: }
483: return(0);
484: }
485:
488: /*@
489: PetscDrawHGSetColor - Sets the color the bars will be drawn with.
491: Not Collective (ignored except on processor 0 of PetscDrawHG)
493: Input Parameters:
494: + hist - The histogram context
495: - color - one of the colors defined in petscdraw.h or PETSC_DRAW_ROTATE to make each bar a
496: different color
498: Level: intermediate
500: @*/
501: PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist, int color)
502: {
505: hist->color = color;
506: return(0);
507: }
511: /*@
512: PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
513: points are added after this call, the limits will be adjusted to
514: include those additional points.
516: Not Collective (ignored except on processor 0 of PetscDrawHG)
518: Input Parameters:
519: + hist - The histogram context
520: - x_min,x_max,y_min,y_max - The limits
522: Level: intermediate
524: Concepts: histogram^setting axis
525: @*/
526: PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
527: {
530: hist->xmin = x_min;
531: hist->xmax = x_max;
532: hist->ymin = y_min;
533: hist->ymax = y_max;
534: return(0);
535: }
539: /*@
540: PetscDrawHGCalcStats - Turns on calculation of descriptive statistics
542: Not collective
544: Input Parameters:
545: + hist - The histogram context
546: - calc - Flag for calculation
548: Level: intermediate
550: .keywords: draw, histogram, statistics
552: @*/
553: PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
554: {
557: hist->calcStats = calc;
558: return(0);
559: }
563: /*@
564: PetscDrawHGIntegerBins - Turns on integer width bins
566: Not collective
568: Input Parameters:
569: + hist - The histogram context
570: - ints - Flag for integer width bins
572: Level: intermediate
574: .keywords: draw, histogram, statistics
575: @*/
576: PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
577: {
580: hist->integerBins = ints;
581: return(0);
582: }
586: /*@C
587: PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
588: This is useful if one wants to change some axis property, such as
589: labels, color, etc. The axis context should not be destroyed by the
590: application code.
592: Not Collective (ignored except on processor 0 of PetscDrawHG)
594: Input Parameter:
595: . hist - The histogram context
597: Output Parameter:
598: . axis - The axis context
600: Level: intermediate
602: @*/
603: PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist, PetscDrawAxis *axis)
604: {
608: *axis = hist->axis;
609: return(0);
610: }
614: /*@C
615: PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
617: Not Collective, PetscDraw is parallel if PetscDrawHG is parallel
619: Input Parameter:
620: . hist - The histogram context
622: Output Parameter:
623: . win - The draw context
625: Level: intermediate
627: @*/
628: PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist, PetscDraw *win)
629: {
633: *win = hist->win;
634: return(0);
635: }