1 /*
2 Contains the data structure for plotting a histogram in a window with an axis.
3 */
4 #include <petscdraw.h> /*I "petscdraw.h" I*/
5 #include <petsc/private/petscimpl.h> /*I "petscsys.h" I*/
6 #include <petscviewer.h> /*I "petscviewer.h" I*/
7
8 PetscClassId PETSC_DRAWHG_CLASSID = 0;
9
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 PetscReal *weights;
25 int color;
26 PetscBool calcStats;
27 PetscBool integerBins;
28 };
29
30 #define CHUNKSIZE 100
31
32 /*@
33 PetscDrawHGCreate - Creates a histogram data structure.
34
35 Collective
36
37 Input Parameters:
38 + draw - The window where the graph will be made
39 - bins - The number of bins to use
40
41 Output Parameter:
42 . hist - The histogram context
43
44 Level: intermediate
45
46 Notes:
47 The difference between a bar chart, `PetscDrawBar`, and a histogram, `PetscDrawHG`, is explained here <https://stattrek.com/statistics/charts/histogram.aspx?Tutorial=AP>
48
49 The histogram is only displayed when `PetscDrawHGDraw()` is called.
50
51 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
52 zeroth MPI process in the communicator. All MPI processes in the communicator must call `PetscDrawHGDraw()` to display the updated graph.
53
54 .seealso: `PetscDrawHGDestroy()`, `PetscDrawHG`, `PetscDrawBarCreate()`, `PetscDrawBar`, `PetscDrawLGCreate()`, `PetscDrawLG`, `PetscDrawSPCreate()`, `PetscDrawSP`,
55 `PetscDrawHGSetNumberBins()`, `PetscDrawHGReset()`, `PetscDrawHGAddValue()`, `PetscDrawHGDraw()`, `PetscDrawHGSave()`, `PetscDrawHGView()`, `PetscDrawHGSetColor()`,
56 `PetscDrawHGSetLimits()`, `PetscDrawHGCalcStats()`, `PetscDrawHGIntegerBins()`, `PetscDrawHGGetAxis()`, `PetscDrawAxis`, `PetscDrawHGGetDraw()`
57 @*/
PetscDrawHGCreate(PetscDraw draw,int bins,PetscDrawHG * hist)58 PetscErrorCode PetscDrawHGCreate(PetscDraw draw, int bins, PetscDrawHG *hist)
59 {
60 PetscDrawHG h;
61
62 PetscFunctionBegin;
63 PetscValidHeaderSpecific(draw, PETSC_DRAW_CLASSID, 1);
64 PetscValidLogicalCollectiveInt(draw, bins, 2);
65 PetscAssertPointer(hist, 3);
66
67 PetscCall(PetscHeaderCreate(h, PETSC_DRAWHG_CLASSID, "DrawHG", "Histogram", "Draw", PetscObjectComm((PetscObject)draw), PetscDrawHGDestroy, NULL));
68 PetscCall(PetscObjectReference((PetscObject)draw));
69 h->win = draw;
70 h->view = NULL;
71 h->destroy = NULL;
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.e-6;
77 h->numBins = bins;
78 h->maxBins = bins;
79 PetscCall(PetscMalloc1(h->maxBins, &h->bins));
80 h->numValues = 0;
81 h->maxValues = CHUNKSIZE;
82 h->calcStats = PETSC_FALSE;
83 h->integerBins = PETSC_FALSE;
84 PetscCall(PetscMalloc2(h->maxValues, &h->values, h->maxValues, &h->weights));
85 PetscCall(PetscDrawAxisCreate(draw, &h->axis));
86 *hist = h;
87 PetscFunctionReturn(PETSC_SUCCESS);
88 }
89
90 /*@
91 PetscDrawHGSetNumberBins - Change the number of bins that are to be drawn in the histogram
92
93 Logically Collective
94
95 Input Parameters:
96 + hist - The histogram context.
97 - bins - The number of bins.
98
99 Level: intermediate
100
101 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawHGDraw()`, `PetscDrawHGIntegerBins()`
102 @*/
PetscDrawHGSetNumberBins(PetscDrawHG hist,int bins)103 PetscErrorCode PetscDrawHGSetNumberBins(PetscDrawHG hist, int bins)
104 {
105 PetscFunctionBegin;
106 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
107 PetscValidLogicalCollectiveInt(hist, bins, 2);
108
109 if (hist->maxBins < bins) {
110 PetscCall(PetscFree(hist->bins));
111 PetscCall(PetscMalloc1(bins, &hist->bins));
112 hist->maxBins = bins;
113 }
114 hist->numBins = bins;
115 PetscFunctionReturn(PETSC_SUCCESS);
116 }
117
118 /*@
119 PetscDrawHGReset - Clears histogram to allow for reuse with new data.
120
121 Logically Collective
122
123 Input Parameter:
124 . hist - The histogram context.
125
126 Level: intermediate
127
128 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawHGDraw()`, `PetscDrawHGAddValue()`
129 @*/
PetscDrawHGReset(PetscDrawHG hist)130 PetscErrorCode PetscDrawHGReset(PetscDrawHG hist)
131 {
132 PetscFunctionBegin;
133 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
134
135 hist->xmin = PETSC_MAX_REAL;
136 hist->xmax = PETSC_MIN_REAL;
137 hist->ymin = 0.0;
138 hist->ymax = 1.e-6;
139 hist->numValues = 0;
140 PetscFunctionReturn(PETSC_SUCCESS);
141 }
142
143 /*@
144 PetscDrawHGDestroy - Frees all space taken up by histogram data structure.
145
146 Collective
147
148 Input Parameter:
149 . hist - The histogram context
150
151 Level: intermediate
152
153 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`
154 @*/
PetscDrawHGDestroy(PetscDrawHG * hist)155 PetscErrorCode PetscDrawHGDestroy(PetscDrawHG *hist)
156 {
157 PetscFunctionBegin;
158 if (!*hist) PetscFunctionReturn(PETSC_SUCCESS);
159 PetscValidHeaderSpecific(*hist, PETSC_DRAWHG_CLASSID, 1);
160 if (--((PetscObject)*hist)->refct > 0) {
161 *hist = NULL;
162 PetscFunctionReturn(PETSC_SUCCESS);
163 }
164
165 PetscCall(PetscFree((*hist)->bins));
166 PetscCall(PetscFree2((*hist)->values, (*hist)->weights));
167 PetscCall(PetscDrawAxisDestroy(&(*hist)->axis));
168 PetscCall(PetscDrawDestroy(&(*hist)->win));
169 PetscCall(PetscHeaderDestroy(hist));
170 PetscFunctionReturn(PETSC_SUCCESS);
171 }
172
173 /*@
174 PetscDrawHGAddValue - Adds another value to the histogram.
175
176 Logically Collective
177
178 Input Parameters:
179 + hist - The histogram
180 - value - The value
181
182 Level: intermediate
183
184 Note:
185 Calls to this function are used to create a standard histogram with integer bin heights. Use calls to `PetscDrawHGAddWeightedValue()` to create a histogram with non-integer bin heights.
186
187 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawHGDraw()`, `PetscDrawHGReset()`, `PetscDrawHGAddWeightedValue()`
188 @*/
PetscDrawHGAddValue(PetscDrawHG hist,PetscReal value)189 PetscErrorCode PetscDrawHGAddValue(PetscDrawHG hist, PetscReal value)
190 {
191 PetscFunctionBegin;
192 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
193
194 /* Allocate more memory if necessary */
195 if (hist->numValues >= hist->maxValues) {
196 PetscReal *tmp, *tmpw;
197
198 PetscCall(PetscMalloc2(hist->maxValues + CHUNKSIZE, &tmp, hist->maxValues + CHUNKSIZE, &tmpw));
199 PetscCall(PetscArraycpy(tmp, hist->values, hist->maxValues));
200 PetscCall(PetscArraycpy(tmpw, hist->weights, hist->maxValues));
201 PetscCall(PetscFree2(hist->values, hist->weights));
202
203 hist->values = tmp;
204 hist->weights = tmpw;
205 hist->maxValues += CHUNKSIZE;
206 }
207 /* I disagree with the original PETSc implementation here. There should be no overshoot, but rather the
208 stated convention of using half-open intervals (always the way to go) */
209 if (!hist->numValues && (hist->xmin == PETSC_MAX_REAL) && (hist->xmax == PETSC_MIN_REAL)) {
210 hist->xmin = value;
211 hist->xmax = value;
212 #if 1
213 } else {
214 /* Update limits */
215 if (value > hist->xmax) hist->xmax = value;
216 if (value < hist->xmin) hist->xmin = value;
217 #else
218 } else if (hist->numValues == 1) {
219 /* Update limits -- We need to overshoot the largest value somewhat */
220 if (value > hist->xmax) hist->xmax = value + 0.001 * (value - hist->xmin) / (PetscReal)hist->numBins;
221 if (value < hist->xmin) {
222 hist->xmin = value;
223 hist->xmax = hist->xmax + 0.001 * (hist->xmax - hist->xmin) / (PetscReal)hist->numBins;
224 }
225 } else {
226 /* Update limits -- We need to overshoot the largest value somewhat */
227 if (value > hist->xmax) hist->xmax = value + 0.001 * (hist->xmax - hist->xmin) / (PetscReal)hist->numBins;
228 if (value < hist->xmin) hist->xmin = value;
229 #endif
230 }
231
232 hist->values[hist->numValues] = value;
233 hist->weights[hist->numValues] = 1.;
234 ++hist->numValues;
235 PetscFunctionReturn(PETSC_SUCCESS);
236 }
237
238 /*@
239 PetscDrawHGAddWeightedValue - Adds another value to the histogram with a weight.
240
241 Logically Collective
242
243 Input Parameters:
244 + hist - The histogram
245 . value - The value
246 - weight - The value weight
247
248 Level: intermediate
249
250 Notes:
251 Calls to this function are used to create a histogram with non-integer bin heights. Use calls to `PetscDrawHGAddValue()` to create a standard histogram with integer bin heights.
252
253 This allows us to histogram frequency and probability distributions (<https://learnche.org/pid/univariate-review/histograms-and-probability-distributions>). We can use this to look at particle weight distributions in Particle-in-Cell (PIC) methods, for example.
254
255 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawHGDraw()`, `PetscDrawHGReset()`, `PetscDrawHGAddValue()`
256 @*/
PetscDrawHGAddWeightedValue(PetscDrawHG hist,PetscReal value,PetscReal weight)257 PetscErrorCode PetscDrawHGAddWeightedValue(PetscDrawHG hist, PetscReal value, PetscReal weight)
258 {
259 PetscFunctionBegin;
260 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
261
262 /* Allocate more memory if necessary */
263 if (hist->numValues >= hist->maxValues) {
264 PetscReal *tmp, *tmpw;
265
266 PetscCall(PetscMalloc2(hist->maxValues + CHUNKSIZE, &tmp, hist->maxValues + CHUNKSIZE, &tmpw));
267 PetscCall(PetscArraycpy(tmp, hist->values, hist->maxValues));
268 PetscCall(PetscArraycpy(tmpw, hist->weights, hist->maxValues));
269 PetscCall(PetscFree2(hist->values, hist->weights));
270
271 hist->values = tmp;
272 hist->weights = tmpw;
273 hist->maxValues += CHUNKSIZE;
274 }
275 /* I disagree with the original PETSc implementation here. There should be no overshoot, but rather the
276 stated convention of using half-open intervals (always the way to go) */
277 if (!hist->numValues && (hist->xmin == PETSC_MAX_REAL) && (hist->xmax == PETSC_MIN_REAL)) {
278 hist->xmin = value;
279 hist->xmax = value;
280 } else {
281 /* Update limits */
282 if (value > hist->xmax) hist->xmax = value;
283 if (value < hist->xmin) hist->xmin = value;
284 }
285
286 hist->values[hist->numValues] = value;
287 hist->weights[hist->numValues] = weight;
288 ++hist->numValues;
289 PetscFunctionReturn(PETSC_SUCCESS);
290 }
291
292 /*@
293 PetscDrawHGDraw - Redraws a histogram.
294
295 Collective
296
297 Input Parameter:
298 . hist - The histogram context
299
300 Level: intermediate
301
302 .seealso: `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawHGAddValue()`, `PetscDrawHGReset()`
303 @*/
PetscDrawHGDraw(PetscDrawHG hist)304 PetscErrorCode PetscDrawHGDraw(PetscDrawHG hist)
305 {
306 PetscDraw draw;
307 PetscBool isnull, usewt = PETSC_FALSE;
308 PetscReal xmin, xmax, ymin, ymax, *bins, *values, *weights, binSize, binLeft, binRight, maxHeight, mean, var, totwt = 0.;
309 char title[256];
310 char xlabel[256];
311 PetscInt numValues, initSize, i, p;
312 int bcolor, color, numBins, numBinsOld;
313 PetscMPIInt rank;
314
315 PetscFunctionBegin;
316 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
317 PetscCall(PetscDrawIsNull(hist->win, &isnull));
318 if (isnull) PetscFunctionReturn(PETSC_SUCCESS);
319 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)hist), &rank));
320
321 if ((hist->xmin >= hist->xmax) || (hist->ymin >= hist->ymax)) PetscFunctionReturn(PETSC_SUCCESS);
322 if (hist->numValues < 1) PetscFunctionReturn(PETSC_SUCCESS);
323
324 color = hist->color;
325 if (color == PETSC_DRAW_ROTATE) bcolor = PETSC_DRAW_BLACK + 1;
326 else bcolor = color;
327
328 xmin = hist->xmin;
329 xmax = hist->xmax;
330 ymin = hist->ymin;
331 ymax = hist->ymax;
332 numValues = hist->numValues;
333 values = hist->values;
334 weights = hist->weights;
335 mean = 0.0;
336 var = 0.0;
337
338 draw = hist->win;
339 PetscCall(PetscDrawCheckResizedWindow(draw));
340 PetscCall(PetscDrawClear(draw));
341
342 if (xmin == xmax) {
343 /* Calculate number of points in each bin */
344 bins = hist->bins;
345 bins[0] = 0.;
346 for (p = 0; p < numValues; p++) {
347 if (values[p] == xmin) bins[0] += weights[0];
348 mean += values[p] * weights[p];
349 var += values[p] * values[p] * weights[p];
350 totwt += weights[p];
351 if (weights[p] != 1.) usewt = PETSC_TRUE;
352 }
353 maxHeight = bins[0];
354 if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
355 xmax = xmin + 1;
356 PetscCall(PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax));
357 if (hist->calcStats) {
358 if (usewt) {
359 mean /= totwt;
360 var = (var - totwt * mean * mean) / totwt;
361
362 PetscCall(PetscSNPrintf(xlabel, 256, "Total Weight: %g", (double)totwt));
363 } else {
364 mean /= numValues;
365 if (numValues > 1) var = (var - numValues * mean * mean) / (PetscReal)(numValues - 1);
366 else var = 0.0;
367
368 PetscCall(PetscSNPrintf(xlabel, 256, "Total: %" PetscInt_FMT, numValues));
369 }
370 PetscCall(PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var));
371 PetscCall(PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL));
372 }
373 PetscCall(PetscDrawAxisDraw(hist->axis));
374 PetscDrawCollectiveBegin(draw);
375 if (rank == 0) { /* Draw bins */
376 binLeft = xmin;
377 binRight = xmax;
378 PetscCall(PetscDrawRectangle(draw, binLeft, ymin, binRight, bins[0], bcolor, bcolor, bcolor, bcolor));
379 PetscCall(PetscDrawLine(draw, binLeft, ymin, binLeft, bins[0], PETSC_DRAW_BLACK));
380 PetscCall(PetscDrawLine(draw, binRight, ymin, binRight, bins[0], PETSC_DRAW_BLACK));
381 PetscCall(PetscDrawLine(draw, binLeft, bins[0], binRight, bins[0], PETSC_DRAW_BLACK));
382 }
383 PetscDrawCollectiveEnd(draw);
384 } else {
385 numBins = hist->numBins;
386 numBinsOld = hist->numBins;
387 if (hist->integerBins && (((int)xmax - xmin) + 1.0e-05 > xmax - xmin)) {
388 initSize = ((int)(xmax - xmin)) / numBins;
389 while (initSize * numBins != (int)xmax - xmin) {
390 initSize = PetscMax(initSize - 1, 1);
391 numBins = ((int)(xmax - xmin)) / initSize;
392 PetscCall(PetscDrawHGSetNumberBins(hist, numBins));
393 }
394 }
395 binSize = (xmax - xmin) / (PetscReal)numBins;
396 bins = hist->bins;
397
398 PetscCall(PetscArrayzero(bins, numBins));
399
400 maxHeight = 0.0;
401 for (i = 0; i < numBins; i++) {
402 binLeft = xmin + binSize * i;
403 binRight = xmin + binSize * (i + 1);
404 for (p = 0; p < numValues; p++) {
405 if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i] += weights[p];
406 /* Handle last bin separately */
407 if ((i == numBins - 1) && (values[p] == binRight)) bins[i] += weights[p];
408 if (!i) {
409 mean += values[p] * weights[p];
410 var += values[p] * values[p] * weights[p];
411 totwt += weights[p];
412 if (weights[p] != 1.) usewt = PETSC_TRUE;
413 }
414 }
415 maxHeight = PetscMax(maxHeight, bins[i]);
416 }
417 if (maxHeight > ymax) ymax = hist->ymax = maxHeight;
418
419 PetscCall(PetscDrawAxisSetLimits(hist->axis, xmin, xmax, ymin, ymax));
420 if (hist->calcStats) {
421 if (usewt) {
422 mean /= totwt;
423 var = (var - totwt * mean * mean) / totwt;
424
425 PetscCall(PetscSNPrintf(xlabel, 256, "Total Weight: %g", (double)totwt));
426 } else {
427 mean /= numValues;
428 if (numValues > 1) var = (var - numValues * mean * mean) / (PetscReal)(numValues - 1);
429 else var = 0.0;
430
431 PetscCall(PetscSNPrintf(xlabel, 256, "Total: %" PetscInt_FMT, numValues));
432 }
433 PetscCall(PetscSNPrintf(title, 256, "Mean: %g Var: %g", (double)mean, (double)var));
434 PetscCall(PetscDrawAxisSetLabels(hist->axis, title, xlabel, NULL));
435 }
436 PetscCall(PetscDrawAxisDraw(hist->axis));
437 PetscDrawCollectiveBegin(draw);
438 if (rank == 0) { /* Draw bins */
439 for (i = 0; i < numBins; i++) {
440 binLeft = xmin + binSize * i;
441 binRight = xmin + binSize * (i + 1);
442 PetscCall(PetscDrawRectangle(draw, binLeft, ymin, binRight, bins[i], bcolor, bcolor, bcolor, bcolor));
443 PetscCall(PetscDrawLine(draw, binLeft, ymin, binLeft, bins[i], PETSC_DRAW_BLACK));
444 PetscCall(PetscDrawLine(draw, binRight, ymin, binRight, bins[i], PETSC_DRAW_BLACK));
445 PetscCall(PetscDrawLine(draw, binLeft, bins[i], binRight, bins[i], PETSC_DRAW_BLACK));
446 if (color == PETSC_DRAW_ROTATE && bins[i]) bcolor++;
447 if (bcolor > PETSC_DRAW_BASIC_COLORS - 1) bcolor = PETSC_DRAW_BLACK + 1;
448 }
449 }
450 PetscDrawCollectiveEnd(draw);
451 PetscCall(PetscDrawHGSetNumberBins(hist, numBinsOld));
452 }
453
454 PetscCall(PetscDrawFlush(draw));
455 PetscCall(PetscDrawPause(draw));
456 PetscFunctionReturn(PETSC_SUCCESS);
457 }
458
459 /*@
460 PetscDrawHGSave - Saves a drawn image
461
462 Collective
463
464 Input Parameter:
465 . hg - The histogram context
466
467 Level: intermediate
468
469 .seealso: `PetscDrawSave()`, `PetscDrawHGCreate()`, `PetscDrawHGGetDraw()`, `PetscDrawSetSave()`, `PetscDrawHGDraw()`
470 @*/
PetscDrawHGSave(PetscDrawHG hg)471 PetscErrorCode PetscDrawHGSave(PetscDrawHG hg)
472 {
473 PetscFunctionBegin;
474 PetscValidHeaderSpecific(hg, PETSC_DRAWHG_CLASSID, 1);
475 PetscCall(PetscDrawSave(hg->win));
476 PetscFunctionReturn(PETSC_SUCCESS);
477 }
478
479 /*@
480 PetscDrawHGView - Prints the histogram information to a viewer
481
482 Not Collective
483
484 Input Parameters:
485 + hist - The histogram context
486 - viewer - The viewer to view it with
487
488 Level: beginner
489
490 .seealso: `PetscDrawHG`, `PetscViewer`, `PetscDrawHGCreate()`, `PetscDrawHGGetDraw()`, `PetscDrawSetSave()`, `PetscDrawSave()`, `PetscDrawHGDraw()`
491 @*/
PetscDrawHGView(PetscDrawHG hist,PetscViewer viewer)492 PetscErrorCode PetscDrawHGView(PetscDrawHG hist, PetscViewer viewer)
493 {
494 PetscReal xmax, xmin, *bins, *values, *weights, binSize, binLeft, binRight, mean, var, totwt = 0.;
495 PetscInt numBins, numBinsOld, numValues, initSize, i, p;
496 PetscBool usewt = PETSC_FALSE;
497 int inumBins;
498
499 PetscFunctionBegin;
500 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
501
502 if ((hist->xmin > hist->xmax) || (hist->ymin >= hist->ymax)) PetscFunctionReturn(PETSC_SUCCESS);
503 if (hist->numValues < 1) PetscFunctionReturn(PETSC_SUCCESS);
504
505 if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)hist), &viewer));
506 PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)hist, viewer));
507 xmax = hist->xmax;
508 xmin = hist->xmin;
509 numValues = hist->numValues;
510 values = hist->values;
511 weights = hist->weights;
512 mean = 0.0;
513 var = 0.0;
514 if (xmax == xmin) {
515 /* Calculate number of points in the bin */
516 bins = hist->bins;
517 bins[0] = 0.;
518 for (p = 0; p < numValues; p++) {
519 if (values[p] == xmin) bins[0] += weights[0];
520 mean += values[p] * weights[p];
521 var += values[p] * values[p] * weights[p];
522 totwt += weights[p];
523 if (weights[p] != 1.) usewt = PETSC_TRUE;
524 }
525 /* Draw bins */
526 PetscCall(PetscViewerASCIIPrintf(viewer, "Bin %2d (%6.2g - %6.2g): %.0g\n", 0, (double)xmin, (double)xmax, (double)bins[0]));
527 } else {
528 numBins = hist->numBins;
529 numBinsOld = hist->numBins;
530 if (hist->integerBins && (((int)xmax - xmin) + 1.0e-05 > xmax - xmin)) {
531 initSize = (int)((int)xmax - xmin) / numBins;
532 while (initSize * numBins != (int)xmax - xmin) {
533 initSize = PetscMax(initSize - 1, 1);
534 numBins = (int)((int)xmax - xmin) / initSize;
535 PetscCall(PetscCIntCast(numBins, &inumBins));
536 PetscCall(PetscDrawHGSetNumberBins(hist, inumBins));
537 }
538 }
539 binSize = (xmax - xmin) / numBins;
540 bins = hist->bins;
541
542 /* Calculate number of points in each bin */
543 PetscCall(PetscArrayzero(bins, numBins));
544 for (i = 0; i < numBins; i++) {
545 binLeft = xmin + binSize * i;
546 binRight = xmin + binSize * (i + 1);
547 for (p = 0; p < numValues; p++) {
548 if ((values[p] >= binLeft) && (values[p] < binRight)) bins[i] += weights[p];
549 /* Handle last bin separately */
550 if ((i == numBins - 1) && (values[p] == binRight)) bins[i] += weights[p];
551 if (!i) {
552 mean += values[p] * weights[p];
553 var += values[p] * values[p] * weights[p];
554 totwt += weights[p];
555 if (weights[p] != 1.) usewt = PETSC_TRUE;
556 }
557 }
558 }
559 /* Draw bins */
560 for (i = 0; i < numBins; i++) {
561 binLeft = xmin + binSize * i;
562 binRight = xmin + binSize * (i + 1);
563 PetscCall(PetscViewerASCIIPrintf(viewer, "Bin %2" PetscInt_FMT " (%6.2g - %6.2g): %.0g\n", i, (double)binLeft, (double)binRight, (double)bins[i]));
564 }
565 PetscCall(PetscCIntCast(numBinsOld, &inumBins));
566 PetscCall(PetscDrawHGSetNumberBins(hist, inumBins));
567 }
568
569 if (hist->calcStats) {
570 if (usewt) {
571 mean /= totwt;
572 var = (var - totwt * mean * mean) / totwt;
573
574 PetscCall(PetscViewerASCIIPrintf(viewer, "Total Weight: %g\n", (double)totwt));
575 } else {
576 mean /= numValues;
577 if (numValues > 1) var = (var - numValues * mean * mean) / (PetscReal)(numValues - 1);
578 else var = 0.0;
579
580 PetscCall(PetscViewerASCIIPrintf(viewer, "Total: %" PetscInt_FMT "\n", numValues));
581 }
582 PetscCall(PetscViewerASCIIPrintf(viewer, "Mean: %g Var: %g\n", (double)mean, (double)var));
583 }
584 PetscFunctionReturn(PETSC_SUCCESS);
585 }
586
587 /*@
588 PetscDrawHGSetColor - Sets the color the bars will be drawn with.
589
590 Logically Collective
591
592 Input Parameters:
593 + hist - The histogram context
594 - color - one of the colors defined in petscdraw.h or `PETSC_DRAW_ROTATE` to make each bar a different color
595
596 Level: intermediate
597
598 .seealso: `PetscDrawHG`, `PetscDrawHGCreate()`, `PetscDrawHGGetDraw()`, `PetscDrawSetSave()`, `PetscDrawSave()`, `PetscDrawHGDraw()`, `PetscDrawHGGetAxis()`
599 @*/
PetscDrawHGSetColor(PetscDrawHG hist,int color)600 PetscErrorCode PetscDrawHGSetColor(PetscDrawHG hist, int color)
601 {
602 PetscFunctionBegin;
603 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
604
605 hist->color = color;
606 PetscFunctionReturn(PETSC_SUCCESS);
607 }
608
609 /*@
610 PetscDrawHGSetLimits - Sets the axis limits for a histogram. If more
611 points are added after this call, the limits will be adjusted to
612 include those additional points.
613
614 Logically Collective
615
616 Input Parameters:
617 + hist - The histogram context
618 . x_min - the horizontal lower limit
619 . x_max - the horizontal upper limit
620 . y_min - the vertical lower limit
621 - y_max - the vertical upper limit
622
623 Level: intermediate
624
625 .seealso: `PetscDrawHG`, `PetscDrawHGCreate()`, `PetscDrawHGGetDraw()`, `PetscDrawSetSave()`, `PetscDrawSave()`, `PetscDrawHGDraw()`, `PetscDrawHGGetAxis()`
626 @*/
PetscDrawHGSetLimits(PetscDrawHG hist,PetscReal x_min,PetscReal x_max,int y_min,int y_max)627 PetscErrorCode PetscDrawHGSetLimits(PetscDrawHG hist, PetscReal x_min, PetscReal x_max, int y_min, int y_max)
628 {
629 PetscFunctionBegin;
630 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
631
632 hist->xmin = x_min;
633 hist->xmax = x_max;
634 hist->ymin = y_min;
635 hist->ymax = y_max;
636 PetscFunctionReturn(PETSC_SUCCESS);
637 }
638
639 /*@
640 PetscDrawHGCalcStats - Turns on calculation of descriptive statistics associated with the histogram
641
642 Not Collective
643
644 Input Parameters:
645 + hist - The histogram context
646 - calc - Flag for calculation
647
648 Level: intermediate
649
650 .seealso: `PetscDrawHG`, `PetscDrawHGCreate()`, `PetscDrawHGAddValue()`, `PetscDrawHGView()`, `PetscDrawHGDraw()`
651 @*/
PetscDrawHGCalcStats(PetscDrawHG hist,PetscBool calc)652 PetscErrorCode PetscDrawHGCalcStats(PetscDrawHG hist, PetscBool calc)
653 {
654 PetscFunctionBegin;
655 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
656
657 hist->calcStats = calc;
658 PetscFunctionReturn(PETSC_SUCCESS);
659 }
660
661 /*@
662 PetscDrawHGIntegerBins - Turns on integer width bins
663
664 Not Collective
665
666 Input Parameters:
667 + hist - The histogram context
668 - ints - Flag for integer width bins
669
670 Level: intermediate
671
672 .seealso: `PetscDrawHG`, `PetscDrawHGCreate()`, `PetscDrawHGAddValue()`, `PetscDrawHGView()`, `PetscDrawHGDraw()`, `PetscDrawHGSetColor()`
673 @*/
PetscDrawHGIntegerBins(PetscDrawHG hist,PetscBool ints)674 PetscErrorCode PetscDrawHGIntegerBins(PetscDrawHG hist, PetscBool ints)
675 {
676 PetscFunctionBegin;
677 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
678
679 hist->integerBins = ints;
680 PetscFunctionReturn(PETSC_SUCCESS);
681 }
682
683 /*@
684 PetscDrawHGGetAxis - Gets the axis context associated with a histogram.
685 This is useful if one wants to change some axis property, such as
686 labels, color, etc. The axis context should not be destroyed by the
687 application code.
688
689 Not Collective, axis is parallel if hist is parallel
690
691 Input Parameter:
692 . hist - The histogram context
693
694 Output Parameter:
695 . axis - The axis context
696
697 Level: intermediate
698
699 .seealso: `PetscDrawHG`, `PetscDrawAxis`, `PetscDrawHGCreate()`, `PetscDrawHGAddValue()`, `PetscDrawHGView()`, `PetscDrawHGDraw()`, `PetscDrawHGSetColor()`, `PetscDrawHGSetLimits()`
700 @*/
PetscDrawHGGetAxis(PetscDrawHG hist,PetscDrawAxis * axis)701 PetscErrorCode PetscDrawHGGetAxis(PetscDrawHG hist, PetscDrawAxis *axis)
702 {
703 PetscFunctionBegin;
704 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
705 PetscAssertPointer(axis, 2);
706 *axis = hist->axis;
707 PetscFunctionReturn(PETSC_SUCCESS);
708 }
709
710 /*@
711 PetscDrawHGGetDraw - Gets the draw context associated with a histogram.
712
713 Not Collective, draw is parallel if hist is parallel
714
715 Input Parameter:
716 . hist - The histogram context
717
718 Output Parameter:
719 . draw - The draw context
720
721 Level: intermediate
722
723 .seealso: `PetscDraw`, `PetscDrawHG`, `PetscDrawHGCreate()`, `PetscDrawHGAddValue()`, `PetscDrawHGView()`, `PetscDrawHGDraw()`, `PetscDrawHGSetColor()`, `PetscDrawAxis`, `PetscDrawHGSetLimits()`
724 @*/
PetscDrawHGGetDraw(PetscDrawHG hist,PetscDraw * draw)725 PetscErrorCode PetscDrawHGGetDraw(PetscDrawHG hist, PetscDraw *draw)
726 {
727 PetscFunctionBegin;
728 PetscValidHeaderSpecific(hist, PETSC_DRAWHG_CLASSID, 1);
729 PetscAssertPointer(draw, 2);
730 *draw = hist->win;
731 PetscFunctionReturn(PETSC_SUCCESS);
732 }
733