Actual source code: dscatter.c
petsc-3.9.4 2018-09-11
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: Add points to the plot with PetscDrawSPAddPoint() or PetscDrawSPAddPoints(); the new points are not displayed until PetscDrawSPDraw() is called.
42: PetscDrawSPReset() removes all the points that have been added
44: 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
45: zeroth MPI process in the communicator. All MPI processes in the communicator must call PetscDrawSPDraw() to display the updated graph.
47: Concepts: scatter plot^creating
49: .seealso: PetscDrawLGCreate(), PetscDrawLG, PetscDrawBarCreate(), PetscDrawBar, PetscDrawHGCreate(), PetscDrawHG, PetscDrawSPDestroy(), PetscDraw, PetscDrawSP, PetscDrawSPSetDimension(), PetscDrawSPReset(),
50: PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPDraw(), PetscDrawSPSave(), PetscDrawSPSetLimits(), PetscDrawSPGetAxis(),PetscDrawAxis, PetscDrawSPGetDraw()
51: @*/
52: PetscErrorCode PetscDrawSPCreate(PetscDraw draw,int dim,PetscDrawSP *drawsp)
53: {
54: PetscDrawSP sp;
62: PetscHeaderCreate(sp,PETSC_DRAWSP_CLASSID,"DrawSP","Scatter Plot","Draw",PetscObjectComm((PetscObject)draw),PetscDrawSPDestroy,NULL);
63: PetscLogObjectParent((PetscObject)draw,(PetscObject)sp);
65: PetscObjectReference((PetscObject)draw);
66: sp->win = draw;
68: sp->view = NULL;
69: sp->destroy = NULL;
70: sp->nopts = 0;
71: sp->dim = dim;
72: sp->xmin = 1.e20;
73: sp->ymin = 1.e20;
74: sp->xmax = -1.e20;
75: sp->ymax = -1.e20;
77: PetscMalloc2(dim*CHUNCKSIZE,&sp->x,dim*CHUNCKSIZE,&sp->y);
78: PetscLogObjectMemory((PetscObject)sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
80: sp->len = dim*CHUNCKSIZE;
81: sp->loc = 0;
83: PetscDrawAxisCreate(draw,&sp->axis);
84: PetscLogObjectParent((PetscObject)sp,(PetscObject)sp->axis);
86: *drawsp = sp;
87: return(0);
88: }
90: /*@
91: PetscDrawSPSetDimension - Change the number of sets of points that are to be drawn.
93: Logically Collective on PetscDrawSP
95: Input Parameter:
96: + sp - the line graph context.
97: - dim - the number of curves.
99: Level: intermediate
101: Concepts: scatter plot^setting number of data types
103: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints()
105: @*/
106: PetscErrorCode PetscDrawSPSetDimension(PetscDrawSP sp,int dim)
107: {
113: if (sp->dim == dim) return(0);
115: PetscFree2(sp->x,sp->y);
116: sp->dim = dim;
117: PetscMalloc2(dim*CHUNCKSIZE,&sp->x,dim*CHUNCKSIZE,&sp->y);
118: PetscLogObjectMemory((PetscObject)sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
119: sp->len = dim*CHUNCKSIZE;
120: return(0);
121: }
123: /*@
124: PetscDrawSPReset - Clears line graph to allow for reuse with new data.
126: Logically Collective on PetscDrawSP
128: Input Parameter:
129: . sp - the line graph context.
131: Level: intermediate
133: Concepts: scatter plot^resetting
135: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPDraw()
136: @*/
137: PetscErrorCode PetscDrawSPReset(PetscDrawSP sp)
138: {
141: sp->xmin = 1.e20;
142: sp->ymin = 1.e20;
143: sp->xmax = -1.e20;
144: sp->ymax = -1.e20;
145: sp->loc = 0;
146: sp->nopts = 0;
147: return(0);
148: }
150: /*@C
151: PetscDrawSPDestroy - Frees all space taken up by scatter plot data structure.
153: Collective on PetscDrawSP
155: Input Parameter:
156: . sp - the line graph context
158: Level: intermediate
160: .seealso: PetscDrawSPCreate(), PetscDrawSP, PetscDrawSPReset()
162: @*/
163: PetscErrorCode PetscDrawSPDestroy(PetscDrawSP *sp)
164: {
168: if (!*sp) return(0);
170: if (--((PetscObject)(*sp))->refct > 0) {*sp = NULL; return(0);}
172: PetscFree2((*sp)->x,(*sp)->y);
173: PetscDrawAxisDestroy(&(*sp)->axis);
174: PetscDrawDestroy(&(*sp)->win);
175: PetscHeaderDestroy(sp);
176: return(0);
177: }
179: /*@
180: PetscDrawSPAddPoint - Adds another point to each of the scatter plots.
182: Logically Collective on PetscDrawSP
184: Input Parameters:
185: + sp - the scatter plot data structure
186: - x, y - the points to two vectors containing the new x and y
187: point for each curve.
189: Level: intermediate
191: Notes: the new points will not be displayed until a call to PetscDrawSPDraw() is made
193: Concepts: scatter plot^adding points
195: .seealso: PetscDrawSPAddPoints(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPDraw()
197: @*/
198: PetscErrorCode PetscDrawSPAddPoint(PetscDrawSP sp,PetscReal *x,PetscReal *y)
199: {
201: PetscInt i;
206: if (sp->loc+sp->dim >= sp->len) { /* allocate more space */
207: PetscReal *tmpx,*tmpy;
208: PetscMalloc2(sp->len+sp->dim*CHUNCKSIZE,&tmpx,sp->len+sp->dim*CHUNCKSIZE,&tmpy);
209: PetscLogObjectMemory((PetscObject)sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
210: PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
211: PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
212: PetscFree2(sp->x,sp->y);
213: sp->x = tmpx;
214: sp->y = tmpy;
215: sp->len += sp->dim*CHUNCKSIZE;
216: }
217: for (i=0; i<sp->dim; i++) {
218: if (x[i] > sp->xmax) sp->xmax = x[i];
219: if (x[i] < sp->xmin) sp->xmin = x[i];
220: if (y[i] > sp->ymax) sp->ymax = y[i];
221: if (y[i] < sp->ymin) sp->ymin = y[i];
223: sp->x[sp->loc] = x[i];
224: sp->y[sp->loc++] = y[i];
225: }
226: sp->nopts++;
227: return(0);
228: }
231: /*@C
232: PetscDrawSPAddPoints - Adds several points to each of the scatter plots.
234: Logically Collective on PetscDrawSP
236: Input Parameters:
237: + sp - the LineGraph data structure
238: . xx,yy - points to two arrays of pointers that point to arrays
239: containing the new x and y points for each curve.
240: - n - number of points being added
242: Level: intermediate
244: Notes: the new points will not be displayed until a call to PetscDrawSPDraw() is made
246: Concepts: scatter plot^adding points
248: .seealso: PetscDrawSPAddPoint(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPDraw()
249: @*/
250: PetscErrorCode PetscDrawSPAddPoints(PetscDrawSP sp,int n,PetscReal **xx,PetscReal **yy)
251: {
253: PetscInt i,j,k;
254: PetscReal *x,*y;
259: if (sp->loc+n*sp->dim >= sp->len) { /* allocate more space */
260: PetscReal *tmpx,*tmpy;
261: PetscInt chunk = CHUNCKSIZE;
262: if (n > chunk) chunk = n;
263: PetscMalloc2(sp->len+sp->dim*chunk,&tmpx,sp->len+sp->dim*chunk,&tmpy);
264: PetscLogObjectMemory((PetscObject)sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
265: PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
266: PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
267: PetscFree2(sp->x,sp->y);
269: sp->x = tmpx;
270: sp->y = tmpy;
271: sp->len += sp->dim*CHUNCKSIZE;
272: }
273: for (j=0; j<sp->dim; j++) {
274: x = xx[j]; y = yy[j];
275: k = sp->loc + j;
276: for (i=0; i<n; i++) {
277: if (x[i] > sp->xmax) sp->xmax = x[i];
278: if (x[i] < sp->xmin) sp->xmin = x[i];
279: if (y[i] > sp->ymax) sp->ymax = y[i];
280: if (y[i] < sp->ymin) sp->ymin = y[i];
282: sp->x[k] = x[i];
283: sp->y[k] = y[i];
284: k += sp->dim;
285: }
286: }
287: sp->loc += n*sp->dim;
288: sp->nopts += n;
289: return(0);
290: }
292: /*@
293: PetscDrawSPDraw - Redraws a scatter plot.
295: Collective on PetscDrawSP
297: Input Parameter:
298: + sp - the line graph context
299: - clear - clear the window before drawing the new plot
301: Level: intermediate
303: .seealso: PetscDrawLGDraw(), PetscDrawLGSPDraw(), PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPReset(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints()
305: @*/
306: PetscErrorCode PetscDrawSPDraw(PetscDrawSP sp, PetscBool clear)
307: {
308: PetscReal xmin,xmax,ymin,ymax;
310: PetscMPIInt rank;
311: PetscBool isnull;
312: PetscDraw draw;
316: PetscDrawIsNull(sp->win,&isnull);
317: if (isnull) return(0);
318: MPI_Comm_rank(PetscObjectComm((PetscObject)sp),&rank);
320: if (sp->xmin > sp->xmax || sp->ymin > sp->ymax) return(0);
321: if (sp->nopts < 1) return(0);
323: draw = sp->win;
324: if (clear) {
325: PetscDrawCheckResizedWindow(draw);
326: PetscDrawClear(draw);
327: }
329: xmin = sp->xmin; xmax = sp->xmax; ymin = sp->ymin; ymax = sp->ymax;
330: PetscDrawAxisSetLimits(sp->axis,xmin,xmax,ymin,ymax);
331: PetscDrawAxisDraw(sp->axis);
333: PetscDrawCollectiveBegin(draw);
334: if (!rank) {
335: int i,j,dim=sp->dim,nopts=sp->nopts;
336: for (i=0; i<dim; i++) {
337: for (j=0; j<nopts; j++) {
338: PetscDrawPoint(draw,sp->x[j*dim+i],sp->y[j*dim+i],PETSC_DRAW_RED);
339: }
340: }
341: }
342: PetscDrawCollectiveEnd(draw);
344: PetscDrawFlush(draw);
345: PetscDrawPause(draw);
346: return(0);
347: }
349: /*@
350: PetscDrawSPSave - Saves a drawn image
352: Collective on PetscDrawSP
354: Input Parameter:
355: . sp - the scatter plot context
357: Level: intermediate
359: Concepts: scatter plot^saving
361: .seealso: PetscDrawSPCreate(), PetscDrawSPGetDraw(), PetscDrawSetSave(), PetscDrawSave()
362: @*/
363: PetscErrorCode PetscDrawSPSave(PetscDrawSP sp)
364: {
369: PetscDrawSave(sp->win);
370: return(0);
371: }
373: /*@
374: PetscDrawSPSetLimits - Sets the axis limits for a scatter plot If more
375: points are added after this call, the limits will be adjusted to
376: include those additional points.
378: Logically Collective on PetscDrawSP
380: Input Parameters:
381: + xsp - the line graph context
382: - x_min,x_max,y_min,y_max - the limits
384: Level: intermediate
386: Concepts: scatter plot^setting axis
388: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawSPGetAxis()
389: @*/
390: PetscErrorCode PetscDrawSPSetLimits(PetscDrawSP sp,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
391: {
394: sp->xmin = x_min;
395: sp->xmax = x_max;
396: sp->ymin = y_min;
397: sp->ymax = y_max;
398: return(0);
399: }
401: /*@C
402: PetscDrawSPGetAxis - Gets the axis context associated with a line graph.
403: This is useful if one wants to change some axis property, such as
404: labels, color, etc. The axis context should not be destroyed by the
405: application code.
407: Not Collective, if PetscDrawSP is parallel then PetscDrawAxis is parallel
409: Input Parameter:
410: . sp - the line graph context
412: Output Parameter:
413: . axis - the axis context
415: Level: intermediate
417: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDrawSPAddPoint(), PetscDrawSPAddPoints(), PetscDrawAxis, PetscDrawAxisCreate()
419: @*/
420: PetscErrorCode PetscDrawSPGetAxis(PetscDrawSP sp,PetscDrawAxis *axis)
421: {
425: *axis = sp->axis;
426: return(0);
427: }
429: /*@C
430: PetscDrawSPGetDraw - Gets the draw context associated with a line graph.
432: Not Collective, PetscDraw is parallel if PetscDrawSP is parallel
434: Input Parameter:
435: . sp - the line graph context
437: Output Parameter:
438: . draw - the draw context
440: Level: intermediate
442: .seealso: PetscDrawSP, PetscDrawSPCreate(), PetscDrawSPDraw(), PetscDraw
443: @*/
444: PetscErrorCode PetscDrawSPGetDraw(PetscDrawSP sp,PetscDraw *draw)
445: {
449: *draw = sp->win;
450: return(0);
451: }