Actual source code: dscatter.c

  1: #define PETSC_DLL
  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 petscsys.h

 10: PetscCookie DRAWSP_COOKIE = 0;

 12: struct _p_DrawSP {
 13:   PETSCHEADER(int);
 14:   PetscErrorCode (*destroy)(PetscDrawSP);
 15:   PetscErrorCode (*view)(PetscDrawSP,PetscViewer);
 16:   int           len,loc;
 17:   PetscDraw     win;
 18:   PetscDrawAxis axis;
 19:   PetscReal     xmin,xmax,ymin,ymax,*x,*y;
 20:   int           nopts,dim;
 21: };

 23: #define CHUNCKSIZE 100

 27: /*@C
 28:     PetscDrawSPCreate - Creates a scatter plot data structure.

 30:     Collective over PetscDraw

 32:     Input Parameters:
 33: +   win - the window where the graph will be made.
 34: -   dim - the number of sets of points which will be drawn

 36:     Output Parameters:
 37: .   drawsp - the scatter plot context

 39:    Level: intermediate

 41:    Concepts: scatter plot^creating

 43: .seealso:  PetscDrawSPDestroy()
 44: @*/
 45: PetscErrorCode  PetscDrawSPCreate(PetscDraw draw,int dim,PetscDrawSP *drawsp)
 46: {
 48:   PetscTruth     isnull;
 49:   PetscObject    obj = (PetscObject)draw;
 50:   PetscDrawSP    sp;

 55:   PetscTypeCompare(obj,PETSC_DRAW_NULL,&isnull);
 56:   if (isnull) {
 57:     PetscDrawOpenNull(((PetscObject)obj)->comm,(PetscDraw*)drawsp);
 58:     return(0);
 59:   }
 60:   PetscHeaderCreate(sp,_p_DrawSP,int,DRAWSP_COOKIE,0,"PetscDrawSP",((PetscObject)obj)->comm,PetscDrawSPDestroy,0);
 61:   sp->view    = 0;
 62:   sp->destroy = 0;
 63:   sp->nopts   = 0;
 64:   sp->win     = draw;
 65:   sp->dim     = dim;
 66:   sp->xmin    = 1.e20;
 67:   sp->ymin    = 1.e20;
 68:   sp->xmax    = -1.e20;
 69:   sp->ymax    = -1.e20;
 70:   PetscMalloc2(dim*CHUNCKSIZE,PetscReal,&sp->x,dim*CHUNCKSIZE,PetscReal,&sp->y);
 71:   PetscLogObjectMemory(sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
 72:   sp->len     = dim*CHUNCKSIZE;
 73:   sp->loc     = 0;
 74:   PetscDrawAxisCreate(draw,&sp->axis);
 75:   PetscLogObjectParent(sp,sp->axis);
 76:   *drawsp = sp;
 77:   return(0);
 78: }

 82: /*@
 83:    PetscDrawSPSetDimension - Change the number of sets of points  that are to be drawn.

 85:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

 87:    Input Parameter:
 88: +  sp - the line graph context.
 89: -  dim - the number of curves.

 91:    Level: intermediate

 93:    Concepts: scatter plot^setting number of data types

 95: @*/
 96: PetscErrorCode  PetscDrawSPSetDimension(PetscDrawSP sp,int dim)
 97: {

101:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);
103:   if (sp->dim == dim) return(0);

105:   PetscFree2(sp->x,sp->y);
106:   sp->dim     = dim;
107:   PetscMalloc2(dim*CHUNCKSIZE,PetscReal,&sp->x,dim*CHUNCKSIZE,PetscReal,&sp->y);
108:   PetscLogObjectMemory(sp,2*dim*CHUNCKSIZE*sizeof(PetscReal));
109:   sp->len     = dim*CHUNCKSIZE;
110:   return(0);
111: }

115: /*@
116:    PetscDrawSPReset - Clears line graph to allow for reuse with new data.

118:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

120:    Input Parameter:
121: .  sp - the line graph context.

123:    Level: intermediate

125:   Concepts: scatter plot^resetting

127: @*/
128: PetscErrorCode  PetscDrawSPReset(PetscDrawSP sp)
129: {
131:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);
133:   sp->xmin  = 1.e20;
134:   sp->ymin  = 1.e20;
135:   sp->xmax  = -1.e20;
136:   sp->ymax  = -1.e20;
137:   sp->loc   = 0;
138:   sp->nopts = 0;
139:   return(0);
140: }

144: /*@C
145:    PetscDrawSPDestroy - Frees all space taken up by scatter plot data structure.

147:    Collective over PetscDrawSP

149:    Input Parameter:
150: .  sp - the line graph context

152:    Level: intermediate

154: .seealso:  PetscDrawSPCreate()
155: @*/
156: PetscErrorCode  PetscDrawSPDestroy(PetscDrawSP sp)
157: {


163:   if (--((PetscObject)sp)->refct > 0) return(0);
164:   if (((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE){
165:     PetscDrawDestroy((PetscDraw) sp);
166:     return(0);
167:   }
168:   PetscDrawAxisDestroy(sp->axis);
169:   PetscFree2(sp->x,sp->y);
170:   PetscHeaderDestroy(sp);
171:   return(0);
172: }

176: /*@
177:    PetscDrawSPAddPoint - Adds another point to each of the scatter plots.

179:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

181:    Input Parameters:
182: +  sp - the scatter plot data structure
183: -  x, y - the points to two vectors containing the new x and y 
184:           point for each curve.

186:    Level: intermediate

188:    Concepts: scatter plot^adding points

190: .seealso: PetscDrawSPAddPoints()
191: @*/
192: PetscErrorCode  PetscDrawSPAddPoint(PetscDrawSP sp,PetscReal *x,PetscReal *y)
193: {
195:   PetscInt       i;

198:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);

201:   if (sp->loc+sp->dim >= sp->len) { /* allocate more space */
202:     PetscReal *tmpx,*tmpy;
203:     PetscMalloc2(sp->len+sp->dim*CHUNCKSIZE,PetscReal,&tmpx,sp->len+sp->dim*CHUNCKSIZE,PetscReal,&tmpy);
204:     PetscLogObjectMemory(sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
205:     PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
206:     PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
207:     PetscFree2(sp->x,sp->y);
208:     sp->x = tmpx;
209:     sp->y = tmpy;
210:     sp->len += sp->dim*CHUNCKSIZE;
211:   }
212:   for (i=0; i<sp->dim; i++) {
213:     if (x[i] > sp->xmax) sp->xmax = x[i];
214:     if (x[i] < sp->xmin) sp->xmin = x[i];
215:     if (y[i] > sp->ymax) sp->ymax = y[i];
216:     if (y[i] < sp->ymin) sp->ymin = y[i];

218:     sp->x[sp->loc]   = x[i];
219:     sp->y[sp->loc++] = y[i];
220:   }
221:   sp->nopts++;
222:   return(0);
223: }


228: /*@C
229:    PetscDrawSPAddPoints - Adds several points to each of the scatter plots.

231:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

233:    Input Parameters:
234: +  sp - the LineGraph data structure
235: .  xx,yy - points to two arrays of pointers that point to arrays 
236:            containing the new x and y points for each curve.
237: -  n - number of points being added

239:    Level: intermediate

241:    Concepts: scatter plot^adding points

243: .seealso: PetscDrawSPAddPoint()
244: @*/
245: PetscErrorCode  PetscDrawSPAddPoints(PetscDrawSP sp,int n,PetscReal **xx,PetscReal **yy)
246: {
248:   PetscInt       i,j,k;
249:   PetscReal      *x,*y;

252:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);

255:   if (sp->loc+n*sp->dim >= sp->len) { /* allocate more space */
256:     PetscReal *tmpx,*tmpy;
257:     PetscInt  chunk = CHUNCKSIZE;
258:     if (n > chunk) chunk = n;
259:     PetscMalloc2(sp->len+sp->dim*chunk,PetscReal,&tmpx,sp->len+sp->dim*chunk,PetscReal,&tmpy);
260:     PetscLogObjectMemory(sp,2*sp->dim*CHUNCKSIZE*sizeof(PetscReal));
261:     PetscMemcpy(tmpx,sp->x,sp->len*sizeof(PetscReal));
262:     PetscMemcpy(tmpy,sp->y,sp->len*sizeof(PetscReal));
263:     PetscFree2(sp->x,sp->y);
264:     sp->x    = tmpx;
265:     sp->y    = tmpy;
266:     sp->len += sp->dim*CHUNCKSIZE;
267:   }
268:   for (j=0; j<sp->dim; j++) {
269:     x = xx[j]; y = yy[j];
270:     k = sp->loc + j;
271:     for (i=0; i<n; i++) {
272:       if (x[i] > sp->xmax) sp->xmax = x[i];
273:       if (x[i] < sp->xmin) sp->xmin = x[i];
274:       if (y[i] > sp->ymax) sp->ymax = y[i];
275:       if (y[i] < sp->ymin) sp->ymin = y[i];

277:       sp->x[k] = x[i];
278:       sp->y[k] = y[i];
279:       k += sp->dim;
280:     }
281:   }
282:   sp->loc   += n*sp->dim;
283:   sp->nopts += n;
284:   return(0);
285: }

289: /*@
290:    PetscDrawSPDraw - Redraws a scatter plot.

292:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

294:    Input Parameter:
295: .  sp - the line graph context

297:    Level: intermediate

299: .seealso: PetscDrawLGDraw(), PetscDrawLGSPDraw()

301: @*/
302: PetscErrorCode  PetscDrawSPDraw(PetscDrawSP sp)
303: {
304:   PetscReal      xmin=sp->xmin,xmax=sp->xmax,ymin=sp->ymin,ymax=sp->ymax;
306:   PetscInt       i,j,dim = sp->dim,nopts = sp->nopts;
307:   PetscMPIInt    rank;
308:   PetscDraw      draw = sp->win;

311:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);

314:   if (nopts < 1) return(0);
315:   if (xmin > xmax || ymin > ymax) return(0);
316:   PetscDrawClear(draw);
317:   PetscDrawAxisSetLimits(sp->axis,xmin,xmax,ymin,ymax);
318:   PetscDrawAxisDraw(sp->axis);
319: 
320:   MPI_Comm_rank(((PetscObject)sp)->comm,&rank);
321:   if (!rank) {
322:     for (i=0; i<dim; i++) {
323:       for (j=0; j<nopts; j++) {
324:         PetscDrawString(draw,sp->x[j*dim+i],sp->y[j*dim+i],PETSC_DRAW_RED,"x");
325:       }
326:     }
327:   }
328:   PetscDrawFlush(sp->win);
329:   PetscDrawPause(sp->win);
330:   return(0);
331: }
332: 
335: /*@
336:    PetscDrawSPSetLimits - Sets the axis limits for a line graph. If more
337:    points are added after this call, the limits will be adjusted to
338:    include those additional points.

340:    Not Collective (ignored on all processors except processor 0 of PetscDrawSP)

342:    Input Parameters:
343: +  xsp - the line graph context
344: -  x_min,x_max,y_min,y_max - the limits

346:    Level: intermediate

348:    Concepts: scatter plot^setting axis

350: @*/
351: PetscErrorCode  PetscDrawSPSetLimits(PetscDrawSP sp,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
352: {
354:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) return(0);
356:   sp->xmin = x_min;
357:   sp->xmax = x_max;
358:   sp->ymin = y_min;
359:   sp->ymax = y_max;
360:   return(0);
361: }
362: 
365: /*@C
366:    PetscDrawSPGetAxis - Gets the axis context associated with a line graph.
367:    This is useful if one wants to change some axis property, such as
368:    labels, color, etc. The axis context should not be destroyed by the
369:    application code.

371:    Not Collective (except PetscDrawAxis can only be used on processor 0 of PetscDrawSP)

373:    Input Parameter:
374: .  sp - the line graph context

376:    Output Parameter:
377: .  axis - the axis context

379:    Level: intermediate

381: @*/
382: PetscErrorCode  PetscDrawSPGetAxis(PetscDrawSP sp,PetscDrawAxis *axis)
383: {
385:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) {
386:     *axis = 0;
387:     return(0);
388:   }
390:   *axis = sp->axis;
391:   return(0);
392: }

396: /*@C
397:    PetscDrawSPGetDraw - Gets the draw context associated with a line graph.

399:    Not Collective, PetscDraw is parallel if PetscDrawSP is parallel

401:    Input Parameter:
402: .  sp - the line graph context

404:    Output Parameter:
405: .  draw - the draw context

407:    Level: intermediate

409: @*/
410: PetscErrorCode  PetscDrawSPGetDraw(PetscDrawSP sp,PetscDraw *draw)
411: {
415:   if (sp && ((PetscObject)sp)->cookie == PETSC_DRAW_COOKIE) {
416:     *draw = (PetscDraw)sp;
417:   } else {
418:     *draw = sp->win;
419:   }
420:   return(0);
421: }