Actual source code: dscatter.c
petsc-3.10.5 2019-03-28
2: /*
3: Contains the data structure for drawing scatter plots
4: graphs in a window with an axis. This is intended for scatter
5: plots that change dynamically.
6: */
8: #include <petscdraw.h>
9: #include <petsc/private/petscimpl.h>
11: PetscClassId PETSC_DRAWSP_CLASSID = 0;
13: struct _p_PetscDrawSP {
14: PETSCHEADER(int);
15: PetscErrorCode (*destroy)(PetscDrawSP);
16: PetscErrorCode (*view)(PetscDrawSP,PetscViewer);
17: int len,loc;
18: PetscDraw win;
19: PetscDrawAxis axis;
20: PetscReal xmin,xmax,ymin,ymax,*x,*y;
21: int nopts,dim;
22: };
24: #define CHUNCKSIZE 100
26: /*@C
27: PetscDrawSPCreate - Creates a scatter plot data structure.
29: Collective on PetscDraw
31: Input Parameters:
32: + win - the window where the graph will be made.
33: - dim - the number of sets of points which will be drawn
35: Output Parameters:
36: . drawsp - the scatter plot context
38: Level: intermediate
40: Notes:
41: Add points to the plot with PetscDrawSPAddPoint() or PetscDrawSPAddPoints(); the new points are not displayed until PetscDrawSPDraw() is called.
43: PetscDrawSPReset() removes all the points that have been added
45: The MPI communicator that owns the PetscDraw owns this PetscDrawSP, but the calls to set options and add points are ignored on all processes except the
46: zeroth MPI process in the communicator. All MPI processes in the communicator must call PetscDrawSPDraw() to display the updated graph.
48: Concepts: scatter plot^creating
50: .seealso: PetscDrawLGCreate(), PetscDrawLG, PetscDrawBarCreate(), PetscDrawBar, PetscDrawHGCreate(), PetscDrawHG, PetscDrawSPDestroy(), PetscDraw, PetscDrawSP, PetscDrawSPSetDimension(), PetscDrawSPReset(),
51: PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPDraw(), PetscDrawSPSave(), PetscDrawSPSetLimits(), PetscDrawSPGetAxis(),PetscDrawAxis, PetscDrawSPGetDraw()
52: @*/
53: PetscErrorCode PetscDrawSPCreate(PetscDraw draw,int dim,PetscDrawSP *drawsp)
54: {
55: PetscDrawSP sp;
63: PetscHeaderCreate(sp,PETSC_DRAWSP_CLASSID,"DrawSP","Scatter Plot","Draw",PetscObjectComm((PetscObject)draw),PetscDrawSPDestroy,NULL);
64: PetscLogObjectParent((PetscObject)draw,(PetscObject)sp);
66: PetscObjectReference((PetscObject)draw);
67: sp->win = draw;
69: sp->view = NULL;
70: sp->destroy = NULL;
71: sp->nopts = 0;
72: sp->dim = dim;
73: sp->xmin = 1.e20;
74: sp->ymin = 1.e20;
75: sp->xmax = -1.e20;
76: sp->ymax = -1.e20;
78: PetscMalloc2(dim*CHUNCKSIZE,&sp->x,dim*CHUNCKSIZE,&sp->y);
79: PetscLogObjectMemory((PetscObject)sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
81: sp->len = dim*CHUNCKSIZE;
82: sp->loc = 0;
84: PetscDrawAxisCreate(draw,&sp->axis);
85: PetscLogObjectParent((PetscObject)sp,(PetscObject)sp->axis);
87: *drawsp = sp;
88: return(0);
89: }
91: /*@
92: PetscDrawSPSetDimension - Change the number of sets of points that are to be drawn.
94: Logically Collective on PetscDrawSP
96: Input Parameter:
97: + sp - the line graph context.
98: - dim - the number of curves.
100: Level: intermediate
102: Concepts: scatter plot^setting number of data types
104: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints()
106: @*/
107: PetscErrorCode PetscDrawSPSetDimension(PetscDrawSP sp,int dim)
108: {
114: if (sp->dim == dim) return(0);
116: PetscFree2(sp->x,sp->y);
117: sp->dim = dim;
118: PetscMalloc2(dim*CHUNCKSIZE,&sp->x,dim*CHUNCKSIZE,&sp->y);
119: PetscLogObjectMemory((PetscObject)sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
120: sp->len = dim*CHUNCKSIZE;
121: return(0);
122: }
124: /*@
125: PetscDrawSPReset - Clears line graph to allow for reuse with new data.
127: Logically Collective on PetscDrawSP
129: Input Parameter:
130: . sp - the line graph context.
132: Level: intermediate
134: Concepts: scatter plot^resetting
136: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPDraw()
137: @*/
138: PetscErrorCode PetscDrawSPReset(PetscDrawSP sp)
139: {
142: sp->xmin = 1.e20;
143: sp->ymin = 1.e20;
144: sp->xmax = -1.e20;
145: sp->ymax = -1.e20;
146: sp->loc = 0;
147: sp->nopts = 0;
148: return(0);
149: }
151: /*@C
152: PetscDrawSPDestroy - Frees all space taken up by scatter plot data structure.
154: Collective on PetscDrawSP
156: Input Parameter:
157: . sp - the line graph context
159: Level: intermediate
161: .seealso: PetscDrawSPCreate(), PetscDrawSP, PetscDrawSPReset()
163: @*/
164: PetscErrorCode PetscDrawSPDestroy(PetscDrawSP *sp)
165: {
169: if (!*sp) return(0);
171: if (--((PetscObject)(*sp))->refct > 0) {*sp = NULL; return(0);}
173: PetscFree2((*sp)->x,(*sp)->y);
174: PetscDrawAxisDestroy(&(*sp)->axis);
175: PetscDrawDestroy(&(*sp)->win);
176: PetscHeaderDestroy(sp);
177: return(0);
178: }
180: /*@
181: PetscDrawSPAddPoint - Adds another point to each of the scatter plots.
183: Logically Collective on PetscDrawSP
185: Input Parameters:
186: + sp - the scatter plot data structure
187: - x, y - the points to two vectors containing the new x and y
188: point for each curve.
190: Level: intermediate
192: Notes:
193: the new points will not be displayed until a call to PetscDrawSPDraw() is made
195: Concepts: scatter plot^adding points
197: .seealso: PetscDrawSPAddPoints(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPDraw()
199: @*/
200: PetscErrorCode PetscDrawSPAddPoint(PetscDrawSP sp,PetscReal *x,PetscReal *y)
201: {
203: PetscInt i;
208: if (sp->loc+sp->dim >= sp->len) { /* allocate more space */
209: PetscReal *tmpx,*tmpy;
210: PetscMalloc2(sp->len+sp->dim*CHUNCKSIZE,&tmpx,sp->len+sp->dim*CHUNCKSIZE,&tmpy);
211: PetscLogObjectMemory((PetscObject)sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
212: PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
213: PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
214: PetscFree2(sp->x,sp->y);
215: sp->x = tmpx;
216: sp->y = tmpy;
217: sp->len += sp->dim*CHUNCKSIZE;
218: }
219: for (i=0; i<sp->dim; i++) {
220: if (x[i] > sp->xmax) sp->xmax = x[i];
221: if (x[i] < sp->xmin) sp->xmin = x[i];
222: if (y[i] > sp->ymax) sp->ymax = y[i];
223: if (y[i] < sp->ymin) sp->ymin = y[i];
225: sp->x[sp->loc] = x[i];
226: sp->y[sp->loc++] = y[i];
227: }
228: sp->nopts++;
229: return(0);
230: }
233: /*@C
234: PetscDrawSPAddPoints - Adds several points to each of the scatter plots.
236: Logically Collective on PetscDrawSP
238: Input Parameters:
239: + sp - the LineGraph data structure
240: . xx,yy - points to two arrays of pointers that point to arrays
241: containing the new x and y points for each curve.
242: - n - number of points being added
244: Level: intermediate
246: Notes:
247: the new points will not be displayed until a call to PetscDrawSPDraw() is made
249: Concepts: scatter plot^adding points
251: .seealso: PetscDrawSPAddPoint(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPDraw()
252: @*/
253: PetscErrorCode PetscDrawSPAddPoints(PetscDrawSP sp,int n,PetscReal **xx,PetscReal **yy)
254: {
256: PetscInt i,j,k;
257: PetscReal *x,*y;
262: if (sp->loc+n*sp->dim >= sp->len) { /* allocate more space */
263: PetscReal *tmpx,*tmpy;
264: PetscInt chunk = CHUNCKSIZE;
265: if (n > chunk) chunk = n;
266: PetscMalloc2(sp->len+sp->dim*chunk,&tmpx,sp->len+sp->dim*chunk,&tmpy);
267: PetscLogObjectMemory((PetscObject)sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
268: PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
269: PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
270: PetscFree2(sp->x,sp->y);
272: sp->x = tmpx;
273: sp->y = tmpy;
274: sp->len += sp->dim*CHUNCKSIZE;
275: }
276: for (j=0; j<sp->dim; j++) {
277: x = xx[j]; y = yy[j];
278: k = sp->loc + j;
279: for (i=0; i<n; i++) {
280: if (x[i] > sp->xmax) sp->xmax = x[i];
281: if (x[i] < sp->xmin) sp->xmin = x[i];
282: if (y[i] > sp->ymax) sp->ymax = y[i];
283: if (y[i] < sp->ymin) sp->ymin = y[i];
285: sp->x[k] = x[i];
286: sp->y[k] = y[i];
287: k += sp->dim;
288: }
289: }
290: sp->loc += n*sp->dim;
291: sp->nopts += n;
292: return(0);
293: }
295: /*@
296: PetscDrawSPDraw - Redraws a scatter plot.
298: Collective on PetscDrawSP
300: Input Parameter:
301: + sp - the line graph context
302: - clear - clear the window before drawing the new plot
304: Level: intermediate
306: .seealso: PetscDrawLGDraw(), PetscDrawLGSPDraw(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints()
308: @*/
309: PetscErrorCode PetscDrawSPDraw(PetscDrawSP sp, PetscBool clear)
310: {
311: PetscReal xmin,xmax,ymin,ymax;
313: PetscMPIInt rank;
314: PetscBool isnull;
315: PetscDraw draw;
319: PetscDrawIsNull(sp->win,&isnull);
320: if (isnull) return(0);
321: MPI_Comm_rank(PetscObjectComm((PetscObject)sp),&rank);
323: if (sp->xmin > sp->xmax || sp->ymin > sp->ymax) return(0);
324: if (sp->nopts < 1) return(0);
326: draw = sp->win;
327: if (clear) {
328: PetscDrawCheckResizedWindow(draw);
329: PetscDrawClear(draw);
330: }
332: xmin = sp->xmin; xmax = sp->xmax; ymin = sp->ymin; ymax = sp->ymax;
333: PetscDrawAxisSetLimits(sp->axis,xmin,xmax,ymin,ymax);
334: PetscDrawAxisDraw(sp->axis);
336: PetscDrawCollectiveBegin(draw);
337: if (!rank) {
338: int i,j,dim=sp->dim,nopts=sp->nopts;
339: for (i=0; i<dim; i++) {
340: for (j=0; j<nopts; j++) {
341: PetscDrawPoint(draw,sp->x[j*dim+i],sp->y[j*dim+i],PETSC_DRAW_RED);
342: }
343: }
344: }
345: PetscDrawCollectiveEnd(draw);
347: PetscDrawFlush(draw);
348: PetscDrawPause(draw);
349: return(0);
350: }
352: /*@
353: PetscDrawSPSave - Saves a drawn image
355: Collective on PetscDrawSP
357: Input Parameter:
358: . sp - the scatter plot context
360: Level: intermediate
362: Concepts: scatter plot^saving
364: .seealso: PetscDrawSPCreate(), PetscDrawSPGetDraw(), PetscDrawSetSave(), PetscDrawSave()
365: @*/
366: PetscErrorCode PetscDrawSPSave(PetscDrawSP sp)
367: {
372: PetscDrawSave(sp->win);
373: return(0);
374: }
376: /*@
377: PetscDrawSPSetLimits - Sets the axis limits for a scatter plot If more
378: points are added after this call, the limits will be adjusted to
379: include those additional points.
381: Logically Collective on PetscDrawSP
383: Input Parameters:
384: + xsp - the line graph context
385: - x_min,x_max,y_min,y_max - the limits
387: Level: intermediate
389: Concepts: scatter plot^setting axis
391: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPGetAxis()
392: @*/
393: PetscErrorCode PetscDrawSPSetLimits(PetscDrawSP sp,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
394: {
397: sp->xmin = x_min;
398: sp->xmax = x_max;
399: sp->ymin = y_min;
400: sp->ymax = y_max;
401: return(0);
402: }
404: /*@C
405: PetscDrawSPGetAxis - Gets the axis context associated with a line graph.
406: This is useful if one wants to change some axis property, such as
407: labels, color, etc. The axis context should not be destroyed by the
408: application code.
410: Not Collective, if PetscDrawSP is parallel then PetscDrawAxis is parallel
412: Input Parameter:
413: . sp - the line graph context
415: Output Parameter:
416: . axis - the axis context
418: Level: intermediate
420: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawAxis, PetscDrawAxisCreate()
422: @*/
423: PetscErrorCode PetscDrawSPGetAxis(PetscDrawSP sp,PetscDrawAxis *axis)
424: {
428: *axis = sp->axis;
429: return(0);
430: }
432: /*@C
433: PetscDrawSPGetDraw - Gets the draw context associated with a line graph.
435: Not Collective, PetscDraw is parallel if PetscDrawSP is parallel
437: Input Parameter:
438: . sp - the line graph context
440: Output Parameter:
441: . draw - the draw context
443: Level: intermediate
445: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDraw
446: @*/
447: PetscErrorCode PetscDrawSPGetDraw(PetscDrawSP sp,PetscDraw *draw)
448: {
452: *draw = sp->win;
453: return(0);
454: }