Actual source code: dscatter.c
1: /*
2: Contains the data structure for drawing scatter plots
3: graphs in a window with an axis. This is intended for scatter
4: plots that change dynamically.
5: */
7: #include <petscdraw.h>
8: #include <petsc/private/drawimpl.h>
10: PetscClassId PETSC_DRAWSP_CLASSID = 0;
12: /*@
13: PetscDrawSPCreate - Creates a scatter plot data structure.
15: Collective
17: Input Parameters:
18: + draw - the window where the graph will be made.
19: - dim - the number of sets of points which will be drawn
21: Output Parameter:
22: . drawsp - the scatter plot context
24: Level: intermediate
26: Notes:
27: Add points to the plot with `PetscDrawSPAddPoint()` or `PetscDrawSPAddPoints()`; the new points are not displayed until `PetscDrawSPDraw()` is called.
29: `PetscDrawSPReset()` removes all the points that have been added
31: `PetscDrawSPSetDimension()` determines how many point curves are being plotted.
33: The MPI communicator that owns the `PetscDraw` owns this `PetscDrawSP`, and each process can add points. All MPI ranks in the communicator must call `PetscDrawSPDraw()` to display the updated graph.
35: .seealso: `PetscDrawLGCreate()`, `PetscDrawLG`, `PetscDrawBarCreate()`, `PetscDrawBar`, `PetscDrawHGCreate()`, `PetscDrawHG`, `PetscDrawSPDestroy()`, `PetscDraw`, `PetscDrawSP`, `PetscDrawSPSetDimension()`, `PetscDrawSPReset()`,
36: `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawSPDraw()`, `PetscDrawSPSave()`, `PetscDrawSPSetLimits()`, `PetscDrawSPGetAxis()`, `PetscDrawAxis`, `PetscDrawSPGetDraw()`
37: @*/
38: PetscErrorCode PetscDrawSPCreate(PetscDraw draw, int dim, PetscDrawSP *drawsp)
39: {
40: PetscDrawSP sp;
42: PetscFunctionBegin;
44: PetscAssertPointer(drawsp, 3);
46: PetscCall(PetscHeaderCreate(sp, PETSC_DRAWSP_CLASSID, "DrawSP", "Scatter Plot", "Draw", PetscObjectComm((PetscObject)draw), PetscDrawSPDestroy, NULL));
47: PetscCall(PetscObjectReference((PetscObject)draw));
48: sp->win = draw;
49: sp->view = NULL;
50: sp->destroy = NULL;
51: sp->nopts = 0;
52: sp->dim = -1;
53: sp->xmin = (PetscReal)1.e20;
54: sp->ymin = (PetscReal)1.e20;
55: sp->zmin = (PetscReal)1.e20;
56: sp->xmax = (PetscReal)-1.e20;
57: sp->ymax = (PetscReal)-1.e20;
58: sp->zmax = (PetscReal)-1.e20;
59: sp->colorized = PETSC_FALSE;
60: sp->loc = 0;
62: PetscCall(PetscDrawSPSetDimension(sp, dim));
63: PetscCall(PetscDrawAxisCreate(draw, &sp->axis));
65: *drawsp = sp;
66: PetscFunctionReturn(PETSC_SUCCESS);
67: }
69: /*@
70: PetscDrawSPSetDimension - Change the number of points that are added at each `PetscDrawSPAddPoint()`
72: Not Collective
74: Input Parameters:
75: + sp - the scatter plot context.
76: - dim - the number of point curves on this process
78: Level: intermediate
80: .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`
81: @*/
82: PetscErrorCode PetscDrawSPSetDimension(PetscDrawSP sp, int dim)
83: {
84: PetscFunctionBegin;
86: if (sp->dim == dim) PetscFunctionReturn(PETSC_SUCCESS);
87: sp->dim = dim;
88: PetscCall(PetscFree3(sp->x, sp->y, sp->z));
89: PetscCall(PetscMalloc3(dim * PETSC_DRAW_SP_CHUNK_SIZE, &sp->x, dim * PETSC_DRAW_SP_CHUNK_SIZE, &sp->y, dim * PETSC_DRAW_SP_CHUNK_SIZE, &sp->z));
90: sp->len = dim * PETSC_DRAW_SP_CHUNK_SIZE;
91: PetscFunctionReturn(PETSC_SUCCESS);
92: }
94: /*@
95: PetscDrawSPGetDimension - Get the number of sets of points that are to be drawn at each `PetscDrawSPAddPoint()`
97: Not Collective
99: Input Parameter:
100: . sp - the scatter plot context.
102: Output Parameter:
103: . dim - the number of point curves on this process
105: Level: intermediate
107: .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`
108: @*/
109: PetscErrorCode PetscDrawSPGetDimension(PetscDrawSP sp, int *dim)
110: {
111: PetscFunctionBegin;
113: PetscAssertPointer(dim, 2);
114: *dim = sp->dim;
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: /*@
119: PetscDrawSPReset - Clears scatter plot to allow for reuse with new data.
121: Not Collective
123: Input Parameter:
124: . sp - the scatter plot context.
126: Level: intermediate
128: .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawSPDraw()`
129: @*/
130: PetscErrorCode PetscDrawSPReset(PetscDrawSP sp)
131: {
132: PetscFunctionBegin;
134: sp->xmin = (PetscReal)1.e20;
135: sp->ymin = (PetscReal)1.e20;
136: sp->zmin = (PetscReal)1.e20;
137: sp->xmax = (PetscReal)-1.e20;
138: sp->ymax = (PetscReal)-1.e20;
139: sp->zmax = (PetscReal)-1.e20;
140: sp->loc = 0;
141: sp->nopts = 0;
142: PetscFunctionReturn(PETSC_SUCCESS);
143: }
145: /*@
146: PetscDrawSPDestroy - Frees all space taken up by scatter plot data structure.
148: Collective
150: Input Parameter:
151: . sp - the scatter plot context
153: Level: intermediate
155: .seealso: `PetscDrawSPCreate()`, `PetscDrawSP`, `PetscDrawSPReset()`
156: @*/
157: PetscErrorCode PetscDrawSPDestroy(PetscDrawSP *sp)
158: {
159: PetscFunctionBegin;
160: if (!*sp) PetscFunctionReturn(PETSC_SUCCESS);
162: if (--((PetscObject)*sp)->refct > 0) {
163: *sp = NULL;
164: PetscFunctionReturn(PETSC_SUCCESS);
165: }
167: PetscCall(PetscFree3((*sp)->x, (*sp)->y, (*sp)->z));
168: PetscCall(PetscDrawAxisDestroy(&(*sp)->axis));
169: PetscCall(PetscDrawDestroy(&(*sp)->win));
170: PetscCall(PetscHeaderDestroy(sp));
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: /*@
175: PetscDrawSPAddPoint - Adds another point to each of the scatter plot point curves.
177: Not Collective
179: Input Parameters:
180: + sp - the scatter plot data structure
181: . x - the x coordinate values (of length dim) for the points of the curve
182: - y - the y coordinate values (of length dim) for the points of the curve
184: Level: intermediate
186: Note:
187: Here dim is the number of point curves passed to `PetscDrawSPCreate()`. The new points will
188: not be displayed until a call to `PetscDrawSPDraw()` is made.
190: .seealso: `PetscDrawSPAddPoints()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPointColorized()`
191: @*/
192: PetscErrorCode PetscDrawSPAddPoint(PetscDrawSP sp, PetscReal *x, PetscReal *y)
193: {
194: PetscInt i;
196: PetscFunctionBegin;
199: if (sp->loc + sp->dim >= sp->len) { /* allocate more space */
200: PetscReal *tmpx, *tmpy, *tmpz;
201: PetscCall(PetscMalloc3(sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpx, sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpy, sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpz));
202: PetscCall(PetscArraycpy(tmpx, sp->x, sp->len));
203: PetscCall(PetscArraycpy(tmpy, sp->y, sp->len));
204: PetscCall(PetscArraycpy(tmpz, sp->z, sp->len));
205: PetscCall(PetscFree3(sp->x, sp->y, sp->z));
206: sp->x = tmpx;
207: sp->y = tmpy;
208: sp->z = tmpz;
209: sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE;
210: }
211: for (i = 0; i < sp->dim; ++i) {
212: if (x[i] > sp->xmax) sp->xmax = x[i];
213: if (x[i] < sp->xmin) sp->xmin = x[i];
214: if (y[i] > sp->ymax) sp->ymax = y[i];
215: if (y[i] < sp->ymin) sp->ymin = y[i];
217: sp->x[sp->loc] = x[i];
218: sp->y[sp->loc++] = y[i];
219: }
220: ++sp->nopts;
221: PetscFunctionReturn(PETSC_SUCCESS);
222: }
224: /*@C
225: PetscDrawSPAddPoints - Adds several points to each of the scatter plot point curves.
227: Not Collective
229: Input Parameters:
230: + sp - the scatter plot context
231: . xx - array of pointers that point to arrays containing the new x coordinates for each curve.
232: . yy - array of pointers that point to arrays containing the new y points for each curve.
233: - n - number of points being added, each represents a subarray of length dim where dim is the value from `PetscDrawSPGetDimension()`
235: Level: intermediate
237: Note:
238: The new points will not be displayed until a call to `PetscDrawSPDraw()` is made
240: .seealso: `PetscDrawSPAddPoint()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPointColorized()`
241: @*/
242: PetscErrorCode PetscDrawSPAddPoints(PetscDrawSP sp, int n, PetscReal *xx[], PetscReal *yy[])
243: {
244: PetscInt i, j, k;
245: PetscReal *x, *y;
247: PetscFunctionBegin;
250: if (sp->loc + n * sp->dim >= sp->len) { /* allocate more space */
251: PetscReal *tmpx, *tmpy, *tmpz;
252: PetscInt chunk = PETSC_DRAW_SP_CHUNK_SIZE;
253: if (n > chunk) chunk = n;
254: PetscCall(PetscMalloc3(sp->len + sp->dim * chunk, &tmpx, sp->len + sp->dim * chunk, &tmpy, sp->len + sp->dim * chunk, &tmpz));
255: PetscCall(PetscArraycpy(tmpx, sp->x, sp->len));
256: PetscCall(PetscArraycpy(tmpy, sp->y, sp->len));
257: PetscCall(PetscArraycpy(tmpz, sp->z, sp->len));
258: PetscCall(PetscFree3(sp->x, sp->y, sp->z));
260: sp->x = tmpx;
261: sp->y = tmpy;
262: sp->z = tmpz;
263: sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE;
264: }
265: for (j = 0; j < sp->dim; ++j) {
266: x = xx[j];
267: y = yy[j];
268: k = sp->loc + j;
269: for (i = 0; i < n; ++i) {
270: if (x[i] > sp->xmax) sp->xmax = x[i];
271: if (x[i] < sp->xmin) sp->xmin = x[i];
272: if (y[i] > sp->ymax) sp->ymax = y[i];
273: if (y[i] < sp->ymin) sp->ymin = y[i];
275: sp->x[k] = x[i];
276: sp->y[k] = y[i];
277: k += sp->dim;
278: }
279: }
280: sp->loc += n * sp->dim;
281: sp->nopts += n;
282: PetscFunctionReturn(PETSC_SUCCESS);
283: }
285: /*@
286: PetscDrawSPAddPointColorized - Adds another point to each of the scatter plots as well as a numeric value to be used to colorize the scatter point.
288: Not Collective
290: Input Parameters:
291: + sp - the scatter plot data structure
292: . x - array of length dim containing the new x coordinate values for each of the point curves.
293: . y - array of length dim containing the new y coordinate values for each of the point curves.
294: - z - array of length dim containing the numeric values that will be mapped to [0,255] and used for scatter point colors.
296: Level: intermediate
298: Note:
299: The dimensions of the arrays is the number of point curves passed to `PetscDrawSPCreate()`.
300: The new points will not be displayed until a call to `PetscDrawSPDraw()` is made
302: .seealso: `PetscDrawSPAddPoints()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()`
303: @*/
304: PetscErrorCode PetscDrawSPAddPointColorized(PetscDrawSP sp, PetscReal *x, PetscReal *y, PetscReal *z)
305: {
306: PetscInt i;
308: PetscFunctionBegin;
310: sp->colorized = PETSC_TRUE;
311: if (sp->loc + sp->dim >= sp->len) { /* allocate more space */
312: PetscReal *tmpx, *tmpy, *tmpz;
313: PetscCall(PetscMalloc3(sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpx, sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpy, sp->len + sp->dim * PETSC_DRAW_SP_CHUNK_SIZE, &tmpz));
314: PetscCall(PetscArraycpy(tmpx, sp->x, sp->len));
315: PetscCall(PetscArraycpy(tmpy, sp->y, sp->len));
316: PetscCall(PetscArraycpy(tmpz, sp->z, sp->len));
317: PetscCall(PetscFree3(sp->x, sp->y, sp->z));
318: sp->x = tmpx;
319: sp->y = tmpy;
320: sp->z = tmpz;
321: sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE;
322: }
323: for (i = 0; i < sp->dim; ++i) {
324: if (x[i] > sp->xmax) sp->xmax = x[i];
325: if (x[i] < sp->xmin) sp->xmin = x[i];
326: if (y[i] > sp->ymax) sp->ymax = y[i];
327: if (y[i] < sp->ymin) sp->ymin = y[i];
328: if (z[i] < sp->zmin) sp->zmin = z[i];
329: if (z[i] > sp->zmax) sp->zmax = z[i];
330: // if (z[i] > sp->zmax && z[i] < 5.) sp->zmax = z[i];
332: sp->x[sp->loc] = x[i];
333: sp->y[sp->loc] = y[i];
334: sp->z[sp->loc++] = z[i];
335: }
336: ++sp->nopts;
337: PetscFunctionReturn(PETSC_SUCCESS);
338: }
340: /*@
341: PetscDrawSPDraw - Redraws a scatter plot.
343: Collective
345: Input Parameters:
346: + sp - the scatter plot context
347: - clear - clear the window before drawing the new plot
349: Level: intermediate
351: .seealso: `PetscDrawLGDraw()`, `PetscDrawLGSPDraw()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`
352: @*/
353: PetscErrorCode PetscDrawSPDraw(PetscDrawSP sp, PetscBool clear)
354: {
355: PetscDraw draw;
356: PetscBool isnull;
357: PetscMPIInt rank, size;
359: PetscFunctionBegin;
361: draw = sp->win;
362: PetscCall(PetscDrawIsNull(draw, &isnull));
363: if (isnull) PetscFunctionReturn(PETSC_SUCCESS);
364: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)sp), &rank));
365: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)sp), &size));
367: if (clear) {
368: PetscCall(PetscDrawCheckResizedWindow(draw));
369: PetscCall(PetscDrawClear(draw));
370: }
371: {
372: PetscReal lower[2] = {sp->xmin, sp->ymin}, glower[2];
373: PetscReal upper[2] = {sp->xmax, sp->ymax}, gupper[2];
374: PetscCallMPI(MPIU_Allreduce(lower, glower, 2, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)sp)));
375: PetscCallMPI(MPIU_Allreduce(upper, gupper, 2, MPIU_REAL, MPIU_MAX, PetscObjectComm((PetscObject)sp)));
376: PetscCall(PetscDrawAxisSetLimits(sp->axis, glower[0], gupper[0], glower[1], gupper[1]));
377: PetscCall(PetscDrawAxisDraw(sp->axis));
378: }
380: PetscDrawCollectiveBegin(draw);
381: {
382: const int dim = sp->dim, nopts = sp->nopts;
384: for (int i = 0; i < dim; ++i) {
385: for (int p = 0; p < nopts; ++p) {
386: int color = sp->colorized ? PetscDrawRealToColor(sp->z[p * dim], sp->zmin, sp->zmax) : (size > 1 ? PetscDrawRealToColor(rank, 0, size - 1) : PETSC_DRAW_RED);
388: PetscCall(PetscDrawPoint(draw, sp->x[p * dim + i], sp->y[p * dim + i], color));
389: }
390: }
391: }
392: PetscDrawCollectiveEnd(draw);
394: PetscCall(PetscDrawFlush(draw));
395: PetscCall(PetscDrawPause(draw));
396: PetscFunctionReturn(PETSC_SUCCESS);
397: }
399: /*@
400: PetscDrawSPSave - Saves a drawn image
402: Collective
404: Input Parameter:
405: . sp - the scatter plot context
407: Level: intermediate
409: .seealso: `PetscDrawSPCreate()`, `PetscDrawSPGetDraw()`, `PetscDrawSetSave()`, `PetscDrawSave()`
410: @*/
411: PetscErrorCode PetscDrawSPSave(PetscDrawSP sp)
412: {
413: PetscFunctionBegin;
415: PetscCall(PetscDrawSave(sp->win));
416: PetscFunctionReturn(PETSC_SUCCESS);
417: }
419: /*@
420: PetscDrawSPSetLimits - Sets the axis limits for a scatter plot. If more points are added after this call, the limits will be adjusted to include those additional points.
422: Not Collective
424: Input Parameters:
425: + sp - the line graph context
426: . x_min - the horizontal lower limit
427: . x_max - the horizontal upper limit
428: . y_min - the vertical lower limit
429: - y_max - the vertical upper limit
431: Level: intermediate
433: .seealso: `PetscDrawSP`, `PetscDrawAxis`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawSPGetAxis()`
434: @*/
435: PetscErrorCode PetscDrawSPSetLimits(PetscDrawSP sp, PetscReal x_min, PetscReal x_max, PetscReal y_min, PetscReal y_max)
436: {
437: PetscFunctionBegin;
439: sp->xmin = x_min;
440: sp->xmax = x_max;
441: sp->ymin = y_min;
442: sp->ymax = y_max;
443: PetscFunctionReturn(PETSC_SUCCESS);
444: }
446: /*@
447: PetscDrawSPGetAxis - Gets the axis context associated with a scatter plot
449: Not Collective
451: Input Parameter:
452: . sp - the scatter plot context
454: Output Parameter:
455: . axis - the axis context
457: Level: intermediate
459: Note:
460: This is useful if one wants to change some axis property, such as labels, color, etc. The axis context should not be destroyed by the application code.
462: .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawAxis`, `PetscDrawAxisCreate()`
463: @*/
464: PetscErrorCode PetscDrawSPGetAxis(PetscDrawSP sp, PetscDrawAxis *axis)
465: {
466: PetscFunctionBegin;
468: PetscAssertPointer(axis, 2);
469: *axis = sp->axis;
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@
474: PetscDrawSPGetDraw - Gets the draw context associated with a scatter plot
476: Not Collective
478: Input Parameter:
479: . sp - the scatter plot context
481: Output Parameter:
482: . draw - the draw context
484: Level: intermediate
486: .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDraw`
487: @*/
488: PetscErrorCode PetscDrawSPGetDraw(PetscDrawSP sp, PetscDraw *draw)
489: {
490: PetscFunctionBegin;
492: PetscAssertPointer(draw, 2);
493: *draw = sp->win;
494: PetscFunctionReturn(PETSC_SUCCESS);
495: }