xref: /petsc/src/sys/classes/draw/impls/x/ximage.c (revision 8067a7d5d1c29607366ed2aff819d529007817c0)
1*8067a7d5SLisandro Dalcin /*
2*8067a7d5SLisandro Dalcin     Code for getting raster images out of a X image or pixmap
3*8067a7d5SLisandro Dalcin */
4*8067a7d5SLisandro Dalcin 
5*8067a7d5SLisandro Dalcin #include <../src/sys/classes/draw/impls/x/ximpl.h>
6*8067a7d5SLisandro Dalcin 
7*8067a7d5SLisandro Dalcin PETSC_INTERN PetscErrorCode PetscDrawGetImage_X(PetscDraw,unsigned char[][3],unsigned int*,unsigned int*,unsigned char*[]);
8*8067a7d5SLisandro Dalcin 
9*8067a7d5SLisandro Dalcin /*
10*8067a7d5SLisandro Dalcin    Get RGB color entries out of the X colormap
11*8067a7d5SLisandro Dalcin */
12*8067a7d5SLisandro Dalcin #undef __FUNCT__
13*8067a7d5SLisandro Dalcin #define __FUNCT__ "PetscDrawXiGetColorsRGB"
14*8067a7d5SLisandro Dalcin static PetscErrorCode PetscDrawXiGetColorsRGB(PetscDraw_X *Xwin,unsigned char rgb[256][3])
15*8067a7d5SLisandro Dalcin {
16*8067a7d5SLisandro Dalcin   int    k;
17*8067a7d5SLisandro Dalcin   XColor colordef[256];
18*8067a7d5SLisandro Dalcin 
19*8067a7d5SLisandro Dalcin   PetscFunctionBegin;
20*8067a7d5SLisandro Dalcin   for (k=0; k<256; k++) {
21*8067a7d5SLisandro Dalcin     colordef[k].pixel = Xwin->cmapping[k];
22*8067a7d5SLisandro Dalcin     colordef[k].flags = DoRed|DoGreen|DoBlue;
23*8067a7d5SLisandro Dalcin   }
24*8067a7d5SLisandro Dalcin   XQueryColors(Xwin->disp,Xwin->cmap,colordef,256);
25*8067a7d5SLisandro Dalcin   for (k=0; k<256; k++) {
26*8067a7d5SLisandro Dalcin     rgb[k][0] = (unsigned char)(colordef[k].red   >> 8);
27*8067a7d5SLisandro Dalcin     rgb[k][1] = (unsigned char)(colordef[k].green >> 8);
28*8067a7d5SLisandro Dalcin     rgb[k][2] = (unsigned char)(colordef[k].blue  >> 8);
29*8067a7d5SLisandro Dalcin   }
30*8067a7d5SLisandro Dalcin   PetscFunctionReturn(0);
31*8067a7d5SLisandro Dalcin }
32*8067a7d5SLisandro Dalcin 
33*8067a7d5SLisandro Dalcin /*
34*8067a7d5SLisandro Dalcin    Map a pixel value to PETSc color value (index in the colormap)
35*8067a7d5SLisandro Dalcin */
36*8067a7d5SLisandro Dalcin #undef __FUNCT__
37*8067a7d5SLisandro Dalcin #define __FUNCT__ "PetscDrawXiPixelToColor"
38*8067a7d5SLisandro Dalcin PETSC_STATIC_INLINE int PetscDrawXiPixelToColor(PetscDraw_X *Xwin,PetscDrawXiPixVal pixval)
39*8067a7d5SLisandro Dalcin {
40*8067a7d5SLisandro Dalcin   int               color;
41*8067a7d5SLisandro Dalcin   PetscDrawXiPixVal *cmap = Xwin->cmapping;
42*8067a7d5SLisandro Dalcin 
43*8067a7d5SLisandro Dalcin   PetscFunctionBegin;
44*8067a7d5SLisandro Dalcin   for (color=0; color<256; color++)   /* slow linear search */
45*8067a7d5SLisandro Dalcin     if (cmap[color] == pixval) break; /* found color */
46*8067a7d5SLisandro Dalcin   if (PetscUnlikely(color == 256))    /* should not happen */
47*8067a7d5SLisandro Dalcin     color = PETSC_DRAW_BLACK;
48*8067a7d5SLisandro Dalcin   PetscFunctionReturn(color);
49*8067a7d5SLisandro Dalcin }
50*8067a7d5SLisandro Dalcin 
51*8067a7d5SLisandro Dalcin #undef __FUNCT__
52*8067a7d5SLisandro Dalcin #define __FUNCT__ "PetscDrawGetImage_X"
53*8067a7d5SLisandro Dalcin PetscErrorCode PetscDrawGetImage_X(PetscDraw draw,unsigned char palette[256][3],unsigned int *out_w,unsigned int *out_h,unsigned char *out_pixels[])
54*8067a7d5SLisandro Dalcin {
55*8067a7d5SLisandro Dalcin   PetscDraw_X      *Xwin = (PetscDraw_X*)draw->data;
56*8067a7d5SLisandro Dalcin   PetscMPIInt      rank;
57*8067a7d5SLisandro Dalcin   PetscErrorCode   ierr;
58*8067a7d5SLisandro Dalcin 
59*8067a7d5SLisandro Dalcin   PetscFunctionBegin;
60*8067a7d5SLisandro Dalcin   if (out_w)      *out_w      = 0;
61*8067a7d5SLisandro Dalcin   if (out_h)      *out_h      = 0;
62*8067a7d5SLisandro Dalcin   if (out_pixels) *out_pixels = NULL;
63*8067a7d5SLisandro Dalcin   ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)draw),&rank);CHKERRQ(ierr);
64*8067a7d5SLisandro Dalcin 
65*8067a7d5SLisandro Dalcin   /* make sure the X server processed requests from all processes */
66*8067a7d5SLisandro Dalcin   ierr = PetscDrawCollectiveBegin(draw);CHKERRQ(ierr);
67*8067a7d5SLisandro Dalcin   XSync(Xwin->disp,True);
68*8067a7d5SLisandro Dalcin   ierr = PetscDrawCollectiveEnd(draw);CHKERRQ(ierr);
69*8067a7d5SLisandro Dalcin   ierr = MPI_Barrier(PetscObjectComm((PetscObject)draw));CHKERRQ(ierr);
70*8067a7d5SLisandro Dalcin 
71*8067a7d5SLisandro Dalcin   /* only the first process handles the saving business */
72*8067a7d5SLisandro Dalcin   ierr = PetscDrawCollectiveBegin(draw);CHKERRQ(ierr);
73*8067a7d5SLisandro Dalcin   if (!rank) {
74*8067a7d5SLisandro Dalcin     Window        root;
75*8067a7d5SLisandro Dalcin     XImage        *ximage;
76*8067a7d5SLisandro Dalcin     unsigned char *pixels = NULL;
77*8067a7d5SLisandro Dalcin     unsigned int  w,h,dummy;
78*8067a7d5SLisandro Dalcin     int           x,y,p;
79*8067a7d5SLisandro Dalcin     /* get RGB colors out of the colormap */
80*8067a7d5SLisandro Dalcin     ierr = PetscDrawXiGetColorsRGB(Xwin,palette);CHKERRQ(ierr);
81*8067a7d5SLisandro Dalcin     /* get image out of the drawable */
82*8067a7d5SLisandro Dalcin     XGetGeometry(Xwin->disp,PetscDrawXiDrawable(Xwin),&root,&x,&y,&w,&h,&dummy,&dummy);
83*8067a7d5SLisandro Dalcin     ximage = XGetImage(Xwin->disp,PetscDrawXiDrawable(Xwin),0,0,w,h,AllPlanes,ZPixmap);
84*8067a7d5SLisandro Dalcin     if (!ximage) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Cannot XGetImage()");
85*8067a7d5SLisandro Dalcin     /* extract pixel values out of the image  */
86*8067a7d5SLisandro Dalcin     ierr = PetscMalloc1(w*h,&pixels);CHKERRQ(ierr);
87*8067a7d5SLisandro Dalcin     for (p=0,y=0; y<(int)h; y++)
88*8067a7d5SLisandro Dalcin       for (x=0; x<(int)w; x++)
89*8067a7d5SLisandro Dalcin         pixels[p++] = PetscDrawXiPixelToColor(Xwin,XGetPixel(ximage,x,y));
90*8067a7d5SLisandro Dalcin     XDestroyImage(ximage);
91*8067a7d5SLisandro Dalcin     *out_w      = w;
92*8067a7d5SLisandro Dalcin     *out_h      = h;
93*8067a7d5SLisandro Dalcin     *out_pixels = pixels;
94*8067a7d5SLisandro Dalcin   }
95*8067a7d5SLisandro Dalcin   ierr = PetscDrawCollectiveEnd(draw);CHKERRQ(ierr);
96*8067a7d5SLisandro Dalcin   PetscFunctionReturn(0);
97*8067a7d5SLisandro Dalcin }
98