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 */ 6 7 #include <petscdraw.h> /*I "petscdraw.h" I*/ 8 #include <petsc/private/drawimpl.h> /*I "petscsys.h" I*/ 9 10 PetscClassId PETSC_DRAWSP_CLASSID = 0; 11 12 /*@C 13 PetscDrawSPCreate - Creates a scatter plot data structure. 14 15 Collective on draw 16 17 Input Parameters: 18 + win - the window where the graph will be made. 19 - dim - the number of sets of points which will be drawn 20 21 Output Parameters: 22 . drawsp - the scatter plot context 23 24 Level: intermediate 25 26 Notes: 27 Add points to the plot with `PetscDrawSPAddPoint()` or `PetscDrawSPAddPoints()`; the new points are not displayed until `PetscDrawSPDraw()` is called. 28 29 `PetscDrawSPReset()` removes all the points that have been added 30 31 `PetscDrawSPSetDimension()` determines how many point curves are being plotted. 32 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. 34 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 PetscDrawSP sp; 40 41 PetscFunctionBegin; 42 PetscValidHeaderSpecific(draw, PETSC_DRAW_CLASSID, 1); 43 PetscValidPointer(drawsp, 3); 44 45 PetscCall(PetscHeaderCreate(sp, PETSC_DRAWSP_CLASSID, "DrawSP", "Scatter Plot", "Draw", PetscObjectComm((PetscObject)draw), PetscDrawSPDestroy, NULL)); 46 PetscCall(PetscObjectReference((PetscObject)draw)); 47 sp->win = draw; 48 sp->view = NULL; 49 sp->destroy = NULL; 50 sp->nopts = 0; 51 sp->dim = -1; 52 sp->xmin = 1.e20; 53 sp->ymin = 1.e20; 54 sp->zmin = 1.e20; 55 sp->xmax = -1.e20; 56 sp->ymax = -1.e20; 57 sp->zmax = -1.e20; 58 sp->colorized = PETSC_FALSE; 59 sp->loc = 0; 60 61 PetscCall(PetscDrawSPSetDimension(sp, dim)); 62 PetscCall(PetscDrawAxisCreate(draw, &sp->axis)); 63 64 *drawsp = sp; 65 PetscFunctionReturn(0); 66 } 67 68 /*@ 69 PetscDrawSPSetDimension - Change the number of points that are added at each `PetscDrawSPAddPoint()` 70 71 Not collective 72 73 Input Parameters: 74 + sp - the scatter plot context. 75 - dim - the number of point curves on this process 76 77 Level: intermediate 78 79 .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()` 80 @*/ 81 PetscErrorCode PetscDrawSPSetDimension(PetscDrawSP sp, int dim) { 82 PetscFunctionBegin; 83 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 84 if (sp->dim == dim) PetscFunctionReturn(0); 85 sp->dim = dim; 86 PetscCall(PetscFree3(sp->x, sp->y, sp->z)); 87 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)); 88 sp->len = dim * PETSC_DRAW_SP_CHUNK_SIZE; 89 PetscFunctionReturn(0); 90 } 91 92 /*@ 93 PetscDrawSPGetDimension - Get the number of sets of points that are to be drawn at each `PetscDrawSPAddPoint()` 94 95 Not collective 96 97 Input Parameters: 98 . sp - the scatter plot context. 99 100 Output Parameter: 101 . dim - the number of point curves on this process 102 103 Level: intermediate 104 105 .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()` 106 @*/ 107 PetscErrorCode PetscDrawSPGetDimension(PetscDrawSP sp, int *dim) { 108 PetscFunctionBegin; 109 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 110 PetscValidPointer(dim, 2); 111 *dim = sp->dim; 112 PetscFunctionReturn(0); 113 } 114 115 /*@ 116 PetscDrawSPReset - Clears scatter plot to allow for reuse with new data. 117 118 Not collective 119 120 Input Parameter: 121 . sp - the scatter plot context. 122 123 Level: intermediate 124 125 .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawSPDraw()` 126 @*/ 127 PetscErrorCode PetscDrawSPReset(PetscDrawSP sp) { 128 PetscFunctionBegin; 129 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 130 sp->xmin = 1.e20; 131 sp->ymin = 1.e20; 132 sp->zmin = 1.e20; 133 sp->xmax = -1.e20; 134 sp->ymax = -1.e20; 135 sp->zmax = -1.e20; 136 sp->loc = 0; 137 sp->nopts = 0; 138 PetscFunctionReturn(0); 139 } 140 141 /*@ 142 PetscDrawSPDestroy - Frees all space taken up by scatter plot data structure. 143 144 Collective on sp 145 146 Input Parameter: 147 . sp - the scatter plot context 148 149 Level: intermediate 150 151 .seealso: `PetscDrawSPCreate()`, `PetscDrawSP`, `PetscDrawSPReset()` 152 @*/ 153 PetscErrorCode PetscDrawSPDestroy(PetscDrawSP *sp) { 154 PetscFunctionBegin; 155 if (!*sp) PetscFunctionReturn(0); 156 PetscValidHeaderSpecific(*sp, PETSC_DRAWSP_CLASSID, 1); 157 if (--((PetscObject)(*sp))->refct > 0) { 158 *sp = NULL; 159 PetscFunctionReturn(0); 160 } 161 162 PetscCall(PetscFree3((*sp)->x, (*sp)->y, (*sp)->z)); 163 PetscCall(PetscDrawAxisDestroy(&(*sp)->axis)); 164 PetscCall(PetscDrawDestroy(&(*sp)->win)); 165 PetscCall(PetscHeaderDestroy(sp)); 166 PetscFunctionReturn(0); 167 } 168 169 /*@ 170 PetscDrawSPAddPoint - Adds another point to each of the scatter plot point curves. 171 172 Not collective 173 174 Input Parameters: 175 + sp - the scatter plot data structure 176 - x, y - two arrays of length dim containing the new x and y coordinate values for each of the point curves. Here dim is the number of point curves passed to PetscDrawSPCreate() 177 178 Level: intermediate 179 180 Note: 181 The new points will not be displayed until a call to `PetscDrawSPDraw()` is made 182 183 .seealso: `PetscDrawSPAddPoints()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPointColorized()` 184 @*/ 185 PetscErrorCode PetscDrawSPAddPoint(PetscDrawSP sp, PetscReal *x, PetscReal *y) { 186 PetscInt i; 187 188 PetscFunctionBegin; 189 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 190 191 if (sp->loc + sp->dim >= sp->len) { /* allocate more space */ 192 PetscReal *tmpx, *tmpy, *tmpz; 193 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)); 194 PetscCall(PetscArraycpy(tmpx, sp->x, sp->len)); 195 PetscCall(PetscArraycpy(tmpy, sp->y, sp->len)); 196 PetscCall(PetscArraycpy(tmpz, sp->z, sp->len)); 197 PetscCall(PetscFree3(sp->x, sp->y, sp->z)); 198 sp->x = tmpx; 199 sp->y = tmpy; 200 sp->z = tmpz; 201 sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE; 202 } 203 for (i = 0; i < sp->dim; ++i) { 204 if (x[i] > sp->xmax) sp->xmax = x[i]; 205 if (x[i] < sp->xmin) sp->xmin = x[i]; 206 if (y[i] > sp->ymax) sp->ymax = y[i]; 207 if (y[i] < sp->ymin) sp->ymin = y[i]; 208 209 sp->x[sp->loc] = x[i]; 210 sp->y[sp->loc++] = y[i]; 211 } 212 ++sp->nopts; 213 PetscFunctionReturn(0); 214 } 215 216 /*@C 217 PetscDrawSPAddPoints - Adds several points to each of the scatter plot point curves. 218 219 Not collective 220 221 Input Parameters: 222 + sp - the scatter plot context 223 . xx,yy - points to two arrays of pointers that point to arrays containing the new x and y points for each curve. 224 - n - number of points being added, each represents a subarray of length dim where dim is the value from `PetscDrawSPGetDimension()` 225 226 Level: intermediate 227 228 Note: 229 The new points will not be displayed until a call to `PetscDrawSPDraw()` is made 230 231 .seealso: `PetscDrawSPAddPoint()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPointColorized()` 232 @*/ 233 PetscErrorCode PetscDrawSPAddPoints(PetscDrawSP sp, int n, PetscReal **xx, PetscReal **yy) { 234 PetscInt i, j, k; 235 PetscReal *x, *y; 236 237 PetscFunctionBegin; 238 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 239 240 if (sp->loc + n * sp->dim >= sp->len) { /* allocate more space */ 241 PetscReal *tmpx, *tmpy, *tmpz; 242 PetscInt chunk = PETSC_DRAW_SP_CHUNK_SIZE; 243 if (n > chunk) chunk = n; 244 PetscCall(PetscMalloc3(sp->len + sp->dim * chunk, &tmpx, sp->len + sp->dim * chunk, &tmpy, sp->len + sp->dim * chunk, &tmpz)); 245 PetscCall(PetscArraycpy(tmpx, sp->x, sp->len)); 246 PetscCall(PetscArraycpy(tmpy, sp->y, sp->len)); 247 PetscCall(PetscArraycpy(tmpz, sp->z, sp->len)); 248 PetscCall(PetscFree3(sp->x, sp->y, sp->z)); 249 250 sp->x = tmpx; 251 sp->y = tmpy; 252 sp->z = tmpz; 253 sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE; 254 } 255 for (j = 0; j < sp->dim; ++j) { 256 x = xx[j]; 257 y = yy[j]; 258 k = sp->loc + j; 259 for (i = 0; i < n; ++i) { 260 if (x[i] > sp->xmax) sp->xmax = x[i]; 261 if (x[i] < sp->xmin) sp->xmin = x[i]; 262 if (y[i] > sp->ymax) sp->ymax = y[i]; 263 if (y[i] < sp->ymin) sp->ymin = y[i]; 264 265 sp->x[k] = x[i]; 266 sp->y[k] = y[i]; 267 k += sp->dim; 268 } 269 } 270 sp->loc += n * sp->dim; 271 sp->nopts += n; 272 PetscFunctionReturn(0); 273 } 274 275 /*@ 276 PetscDrawSPAddPointColorized - Adds another point to each of the scatter plots as well as a numeric value to be used to colorize the scatter point. 277 278 Not collective 279 280 Input Parameters: 281 + sp - the scatter plot data structure 282 . x, y - two arrays of length dim containing the new x and y coordinate values for each of the point curves. Here dim is the number of point curves passed to `PetscDrawSPCreate()` 283 - z - array of length dim containing the numeric values that will be mapped to [0,255] and used for scatter point colors. 284 285 Level: intermediate 286 287 Note: 288 The new points will not be displayed until a call to `PetscDrawSPDraw()` is made 289 290 .seealso: `PetscDrawSPAddPoints()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()` 291 @*/ 292 PetscErrorCode PetscDrawSPAddPointColorized(PetscDrawSP sp, PetscReal *x, PetscReal *y, PetscReal *z) { 293 PetscInt i; 294 295 PetscFunctionBegin; 296 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 297 sp->colorized = PETSC_TRUE; 298 if (sp->loc + sp->dim >= sp->len) { /* allocate more space */ 299 PetscReal *tmpx, *tmpy, *tmpz; 300 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)); 301 PetscCall(PetscArraycpy(tmpx, sp->x, sp->len)); 302 PetscCall(PetscArraycpy(tmpy, sp->y, sp->len)); 303 PetscCall(PetscArraycpy(tmpz, sp->z, sp->len)); 304 PetscCall(PetscFree3(sp->x, sp->y, sp->z)); 305 sp->x = tmpx; 306 sp->y = tmpy; 307 sp->z = tmpz; 308 sp->len += sp->dim * PETSC_DRAW_SP_CHUNK_SIZE; 309 } 310 for (i = 0; i < sp->dim; ++i) { 311 if (x[i] > sp->xmax) sp->xmax = x[i]; 312 if (x[i] < sp->xmin) sp->xmin = x[i]; 313 if (y[i] > sp->ymax) sp->ymax = y[i]; 314 if (y[i] < sp->ymin) sp->ymin = y[i]; 315 if (z[i] < sp->zmin) sp->zmin = z[i]; 316 if (z[i] > sp->zmax) sp->zmax = z[i]; 317 // if (z[i] > sp->zmax && z[i] < 5.) sp->zmax = z[i]; 318 319 sp->x[sp->loc] = x[i]; 320 sp->y[sp->loc] = y[i]; 321 sp->z[sp->loc++] = z[i]; 322 } 323 ++sp->nopts; 324 PetscFunctionReturn(0); 325 } 326 327 /*@ 328 PetscDrawSPDraw - Redraws a scatter plot. 329 330 Collective on sp 331 332 Input Parameters: 333 + sp - the scatter plot context 334 - clear - clear the window before drawing the new plot 335 336 Level: intermediate 337 338 .seealso: `PetscDrawLGDraw()`, `PetscDrawLGSPDraw()`, `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPReset()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()` 339 @*/ 340 PetscErrorCode PetscDrawSPDraw(PetscDrawSP sp, PetscBool clear) { 341 PetscDraw draw; 342 PetscBool isnull; 343 PetscMPIInt rank, size; 344 345 PetscFunctionBegin; 346 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 347 draw = sp->win; 348 PetscCall(PetscDrawIsNull(draw, &isnull)); 349 if (isnull) PetscFunctionReturn(0); 350 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)sp), &rank)); 351 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)sp), &size)); 352 353 if (clear) { 354 PetscCall(PetscDrawCheckResizedWindow(draw)); 355 PetscCall(PetscDrawClear(draw)); 356 } 357 { 358 PetscReal lower[2] = {sp->xmin, sp->ymin}, glower[2]; 359 PetscReal upper[2] = {sp->xmax, sp->ymax}, gupper[2]; 360 PetscCall(MPIU_Allreduce(lower, glower, 2, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)sp))); 361 PetscCall(MPIU_Allreduce(upper, gupper, 2, MPIU_REAL, MPIU_MAX, PetscObjectComm((PetscObject)sp))); 362 PetscCall(PetscDrawAxisSetLimits(sp->axis, glower[0], gupper[0], glower[1], gupper[1])); 363 PetscCall(PetscDrawAxisDraw(sp->axis)); 364 } 365 366 PetscDrawCollectiveBegin(draw); 367 { 368 const int dim = sp->dim, nopts = sp->nopts; 369 370 for (int i = 0; i < dim; ++i) { 371 for (int p = 0; p < nopts; ++p) { 372 PetscInt color = sp->colorized ? PetscDrawRealToColor(sp->z[p * dim], sp->zmin, sp->zmax) : (size > 1 ? PetscDrawRealToColor(rank, 0, size - 1) : PETSC_DRAW_RED); 373 374 PetscCall(PetscDrawPoint(draw, sp->x[p * dim + i], sp->y[p * dim + i], color)); 375 } 376 } 377 } 378 PetscDrawCollectiveEnd(draw); 379 380 PetscCall(PetscDrawFlush(draw)); 381 PetscCall(PetscDrawPause(draw)); 382 PetscFunctionReturn(0); 383 } 384 385 /*@ 386 PetscDrawSPSave - Saves a drawn image 387 388 Collective on sp 389 390 Input Parameter: 391 . sp - the scatter plot context 392 393 Level: intermediate 394 395 .seealso: `PetscDrawSPSave()`, `PetscDrawSPCreate()`, `PetscDrawSPGetDraw()`, `PetscDrawSetSave()`, `PetscDrawSave()` 396 @*/ 397 PetscErrorCode PetscDrawSPSave(PetscDrawSP sp) { 398 PetscFunctionBegin; 399 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 400 PetscCall(PetscDrawSave(sp->win)); 401 PetscFunctionReturn(0); 402 } 403 404 /*@ 405 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. 406 407 Not collective 408 409 Input Parameters: 410 + xsp - the line graph context 411 - x_min,x_max,y_min,y_max - the limits 412 413 Level: intermediate 414 415 .seealso: `PetscDrawSP`, `PetscDrawAxis`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawSPGetAxis()` 416 @*/ 417 PetscErrorCode PetscDrawSPSetLimits(PetscDrawSP sp, PetscReal x_min, PetscReal x_max, PetscReal y_min, PetscReal y_max) { 418 PetscFunctionBegin; 419 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 420 sp->xmin = x_min; 421 sp->xmax = x_max; 422 sp->ymin = y_min; 423 sp->ymax = y_max; 424 PetscFunctionReturn(0); 425 } 426 427 /*@ 428 PetscDrawSPGetAxis - Gets the axis context associated with a scatter plot 429 430 Not Collective 431 432 Input Parameter: 433 . sp - the scatter plot context 434 435 Output Parameter: 436 . axis - the axis context 437 438 Note: 439 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. 440 441 Level: intermediate 442 443 .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDrawSPAddPoint()`, `PetscDrawSPAddPoints()`, `PetscDrawAxis`, `PetscDrawAxisCreate()` 444 @*/ 445 PetscErrorCode PetscDrawSPGetAxis(PetscDrawSP sp, PetscDrawAxis *axis) { 446 PetscFunctionBegin; 447 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 448 PetscValidPointer(axis, 2); 449 *axis = sp->axis; 450 PetscFunctionReturn(0); 451 } 452 453 /*@ 454 PetscDrawSPGetDraw - Gets the draw context associated with a scatter plot 455 456 Not Collective 457 458 Input Parameter: 459 . sp - the scatter plot context 460 461 Output Parameter: 462 . draw - the draw context 463 464 Level: intermediate 465 466 .seealso: `PetscDrawSP`, `PetscDrawSPCreate()`, `PetscDrawSPDraw()`, `PetscDraw` 467 @*/ 468 PetscErrorCode PetscDrawSPGetDraw(PetscDrawSP sp, PetscDraw *draw) { 469 PetscFunctionBegin; 470 PetscValidHeaderSpecific(sp, PETSC_DRAWSP_CLASSID, 1); 471 PetscValidPointer(draw, 2); 472 *draw = sp->win; 473 PetscFunctionReturn(0); 474 } 475