1 | /***************************************
2 | $Header: /cvsroot/petscgraphics/render.c,v 1.15 2004/11/30 14:35:54 hazelsct Exp $
3 |
4 | This file contains the rendering code for Illuminator, which renders 2-D or
5 | 3-D data into an RGB(A) unsigned char array (using perspective in 3-D).
6 | ***************************************/
7 |
8 | #include "illuminator.h"
9 |
10 | /* Build with -DDEBUG for debugging output */
11 | #undef DPRINTF
12 | #ifdef DEBUG
13 | #define DPRINTF(fmt, args...) PetscPrintf (PETSC_COMM_WORLD, "%s: " fmt, __FUNCT__, args)
14 | #else
15 | #define DPRINTF(fmt, args...)
16 | #endif
17 |
18 |
19 | #undef __FUNCT__
20 | #define __FUNCT__ "pseudocolor"
21 |
22 | /*++++++++++++++++++++++++++++++++++++++
23 | This little function converts a scalar value into an rgb color from red to
24 | blue.
25 |
26 | PetscScalar val Value to convert.
27 |
28 | PetscScalar* scale Array with minimum and maximum values in which to scale
29 | val.
30 |
31 | guchar *pixel Address in rgb buffer where this pixel should be painted.
32 | ++++++++++++++++++++++++++++++++++++++*/
33 |
34 | static inline void pseudocolor (PetscScalar val, PetscScalar* scale,
35 | guchar *pixel)
36 | {
37 | PetscScalar shade = (val - scale[0]) / (scale[1] - scale[0]);
38 | /* Old stuff *
39 | if (shade < 0.2) { /* red -> yellow *
40 | *pixel++ = 255;
41 | *pixel++ = 1275*shade;
42 | *pixel++ = 0; }
43 | else if (shade < 0.4) { /* yellow -> green *
44 | *pixel++ = 510-1275*shade;
45 | *pixel++ = 255;
46 | *pixel++ = 0; }
47 | else if (shade < 0.6) { /* green -> cyan *
48 | *pixel++ = 0;
49 | *pixel++ = 255;
50 | *pixel++ = 1275*shade-510; }
51 | else if (shade < 0.8) { /* cyan -> blue *
52 | *pixel++ = 0;
53 | *pixel++ = 1020-1275*shade;
54 | *pixel++ = 255; }
55 | else { /* blue -> magenta *
56 | *pixel++ = 1275*shade-1020;
57 | *pixel++ = 0;
58 | *pixel++ = 255; }
59 | /* New stuff */
60 | if (shade < 0.25) { /* red -> yellow */
61 | *pixel++ = 255;
62 | *pixel++ = 1020*shade;
63 | *pixel++ = 0; }
64 | else if (shade < 0.5) { /* yellow -> green */
65 | *pixel++ = 510-1020*shade;
66 | *pixel++ = 255;
67 | *pixel++ = 0; }
68 | else if (shade < 0.75) { /* green -> cyan */
69 | *pixel++ = 0;
70 | *pixel++ = 255;
71 | *pixel++ = 1020*shade-510; }
72 | else { /* cyan -> blue */
73 | *pixel++ = 0;
74 | *pixel++ = 1020-1020*shade;
75 | *pixel++ = 255; }
76 | }
77 |
78 |
79 | #undef __FUNCT__
80 | #define __FUNCT__ "pseudoternarycolor"
81 |
82 | /*++++++++++++++++++++++++++++++++++++++
83 | This little function converts two ternary fractions into an rgb color, with
84 | yellow, cyan and magenta indicating the corners.
85 |
86 | PetscScalar A First ternary fraction.
87 |
88 | PetscScalar B Second ternary fraction.
89 |
90 | PetscScalar *scale Array first and second ternary fractions of each of the
91 | three corner values for scaling.
92 |
93 | guchar *pixel Address in rgb buffer where this pixel should be painted.
94 | ++++++++++++++++++++++++++++++++++++++*/
95 |
96 | static inline void pseudoternarycolor
97 | (PetscScalar A, PetscScalar B, PetscScalar *scale, guchar *pixel)
98 | {
99 | PetscScalar x1, x2, inverse_det;
100 |
101 | /* Transform A,B into x1,x2 based on corners */
102 | inverse_det = 1./((scale[2]-scale[0])*(scale[5]-scale[1]) -
103 | (scale[4]-scale[0])*(scale[3]-scale[1]));
104 | x1 = ((A-scale[0])*(scale[5]-scale[1]) -
105 | (B-scale[1])*(scale[4]-scale[0])) * inverse_det;
106 | x2 = ((B-scale[1])*(scale[2]-scale[0]) -
107 | (A-scale[0])*(scale[3]-scale[1])) * inverse_det;
108 |
109 | /* Now colorize */
110 | *pixel++ = 255*(1.-x1);
111 | *pixel++ = 255*(x1+x2);
112 | *pixel++ = 255*(1.-x2);
113 | }
114 |
115 |
116 | #undef __FUNCT__
117 | #define __FUNCT__ "pseudohueintcolor"
118 |
119 | /*++++++++++++++++++++++++++++++++++++++
120 | This little function converts a vector into an rgb color with hue indicating
121 | direction (green, yellow, red, blue at 0, 90, 180, 270 degrees) and intensity
122 | indicating magnitude relative to reference magnitude in scale[1].
123 |
124 | PetscScalar vx Vector's
125 | +latex+$x$-component.
126 | +html+ <i>x</i>-component.
127 |
128 | PetscScalar vy Vector's
129 | +latex+$y$-component.
130 | +html+ <i>y</i>-component.
131 |
132 | PetscScalar *scale Array whose second entry has the reference magnitude.
133 |
134 | guchar *pixel Address in rgb buffer where this pixel should be painted.
135 | ++++++++++++++++++++++++++++++++++++++*/
136 |
137 | static inline void pseudohueintcolor
138 | (PetscScalar vx, PetscScalar vy, PetscScalar *scale, guchar *pixel)
139 | {
140 | PetscScalar mag=sqrt(vx*vx+vy*vy), theta=atan2(vy,vx), red, green, blue;
141 | if (scale[1] <= 0.)
142 | {
143 | *pixel = *(pixel+1) = *(pixel+2) = 0;
144 | return;
145 | }
146 | mag = (mag > scale[1]) ? 1.0 : mag/scale[1];
147 |
148 | red = 2./M_PI * ((theta<-M_PI/2.) ? -M_PI/2.-theta :
149 | ((theta<0.) ? 0. : ((theta<M_PI/2.) ? theta : M_PI/2.)));
150 | green = 2./M_PI * ((theta<-M_PI/2.) ? 0. :
151 | ((theta<0.) ? theta+M_PI/2. :
152 | ((theta<M_PI/2.) ? M_PI/2. : M_PI-theta)));
153 | blue = 2./M_PI * ((theta<-M_PI/2.) ? theta+M_PI :
154 | ((theta<0.) ? -theta : 0.));
155 |
156 | *pixel++ = 255*mag*red;
157 | *pixel++ = 255*mag*green;
158 | *pixel++ = 255*mag*blue;
159 | }
160 |
161 |
162 | #undef __FUNCT__
163 | #define __FUNCT__ "pseudoshearcolor"
164 |
165 | /*++++++++++++++++++++++++++++++++++++++
166 | This little function converts a shear tensor (symmetric, diagonals sum to
167 | zero) into an rgb color with hue indicating direction of the tensile stress
168 | (red, yellow, green, cyan, blue, magenta at 0, 30, 60, 90, 120, 150 degrees
169 | respectively; 180 is equivalent to zero for stress) and intensity
170 | indicating magnitude relative to reference magnitude in scale[0].
171 |
172 | PetscScalar gxx Tensor's
173 | +latex+$xx$-component.
174 | +html+ <i>xx</i>-component.
175 |
176 | PetscScalar gxy Tensor's
177 | +latex+$xy$-component.
178 | +html+ <i>xy</i>-component.
179 |
180 | PetscScalar *scale Array whose first entry has the reference magnitude.
181 |
182 | guchar *pixel Address in rgb buffer where this pixel should be painted.
183 | ++++++++++++++++++++++++++++++++++++++*/
184 |
185 | static inline void pseudoshearcolor
186 | (PetscScalar gxx, PetscScalar gxy, PetscScalar *scale, guchar *pixel)
187 | {
188 | PetscScalar mag=sqrt(gxx*gxx+gxy*gxy), theta=atan2(gxy,gxx),
189 | red,green,blue;
190 | if (scale[0] <= 0.)
191 | {
192 | *pixel = *(pixel+1) = *(pixel+2) = 0;
193 | return;
194 | }
195 | mag = (mag > scale[0]) ? 1.0 : mag/scale[0];
196 |
197 | red = 3./M_PI * ((theta < -2.*M_PI/3.) ? 0. :
198 | ((theta < -M_PI/3.) ? theta + 2.*M_PI/3. :
199 | ((theta < M_PI/3.) ? M_PI/3. :
200 | ((theta < 2.*M_PI/3.) ? 2.*M_PI/3. - theta: 0.))));
201 | green = 3./M_PI * ((theta < -M_PI/3.) ? M_PI/3. :
202 | ((theta < 0.) ? -theta :
203 | ((theta < 2.*M_PI/3.) ? 0. : theta - 2.*M_PI/3.)));
204 | blue = 3./M_PI * ((theta < -2.*M_PI/3) ? -2.*M_PI/3. - theta :
205 | ((theta < 0.) ? 0. :
206 | ((theta < M_PI/3.) ? theta : M_PI/3.)));
207 |
208 | *pixel++ = 255*mag*red;
209 | *pixel++ = 255*mag*green;
210 | *pixel++ = 255*mag*blue;
211 | }
212 |
213 |
214 | #undef __FUNCT__
215 | #define __FUNCT__ "render_scale_2d"
216 |
217 | /*++++++++++++++++++++++++++++++++++++++
218 | This draws a little rwidth x rheight image depicting the scale of a field
219 | variable.
220 |
221 | int render_scale_2d It returns zero or an error code.
222 |
223 | guchar *rgb Image to draw into.
224 |
225 | int rwidth Intended width of the image.
226 |
227 | int rheight Intended height of the image.
228 |
229 | int bytes_per_pixel Number of bytes per pixel in the image.
230 |
231 | field_plot_type fieldtype Type of plot.
232 |
233 | int symmetry Symmetry order for vector scale image.
234 | ++++++++++++++++++++++++++++++++++++++*/
235 |
236 | int render_scale_2d (guchar *rgb, int rwidth, int rheight, int bytes_per_pixel,
237 | field_plot_type fieldtype, int symmetry)
238 | {
239 | PetscScalar scale[6];
240 | int i, j, myheight;
241 |
242 | switch (fieldtype)
243 | {
244 | case FIELD_SCALAR:
245 | scale[0]=0.;
246 | scale[1]=1.;
247 | for (j=0; j<rheight; j++)
248 | for (i=0; i<rwidth; i++)
249 | pseudocolor ((PetscScalar) i/(rwidth-1), scale,
250 | rgb+(j*rwidth + i) * bytes_per_pixel);
251 | return 0;
252 |
253 | case FIELD_TERNARY:
254 | scale[0] = scale[1] = scale[3] = scale[4] = 0.;
255 | scale[2] = scale[5] = 1.;
256 | myheight = ((PetscScalar) rheight > 0.5*sqrt(3.) * rwidth) ?
257 | (int) (0.5*sqrt(3.)* rwidth) : rheight;
258 | for (j=0; j<myheight; j++)
259 | for (i=0; i<(int) (2./sqrt(3.)*(myheight-j)); i++)
260 | pseudoternarycolor
261 | ((PetscScalar) i/(myheight*2./sqrt(3.)),
262 | (PetscScalar) j/myheight, scale,
263 | rgb + (((rheight+myheight)/2 -j)*rwidth +i +(int)(sqrt(3.)/3.*j))*
264 | bytes_per_pixel);
265 | return 0;
266 |
267 | case FIELD_VECTOR:
268 | scale[0] = 0.;
269 | scale[1] = 1.;
270 | myheight = ((rheight > rwidth) ? rwidth : rheight)/2;
271 | for (j=-myheight; j<myheight; j++)
272 | for (i=-(int)(myheight*sqrt(1.-((PetscScalar)j*j/myheight/myheight)));
273 | i<(int)(myheight*sqrt(1.-((PetscScalar)j*j/myheight/myheight)));
274 | i++)
275 | {
276 | PetscScalar mag = (PetscScalar) sqrt (i*i+j*j)/myheight,
277 | theta = atan2 (j, i);
278 | pseudohueintcolor
279 | (mag * cos (symmetry*theta), mag * sin (symmetry*theta), scale,
280 | rgb + ((myheight-j)*rwidth + rwidth/2+i)*bytes_per_pixel);
281 | }
282 | return 0;
283 |
284 | case FIELD_TENSOR_SHEAR:
285 | scale[0] = 1.;
286 | myheight = ((rheight > rwidth) ? rwidth : rheight)/2;
287 | for (j=-myheight; j<myheight; j++)
288 | for (i=-(int)(myheight*sqrt(1.-((PetscScalar)j*j/myheight/myheight)));
289 | i<(int)(myheight*sqrt(1.-((PetscScalar)j*j/myheight/myheight)));
290 | i++)
291 | {
292 | PetscScalar mag = (PetscScalar) sqrt (i*i+j*j)/myheight,
293 | theta = atan2 (j, i);
294 | pseudoshearcolor
295 | (mag * cos (2*theta), mag * sin (2*theta), scale,
296 | rgb + ((myheight-j)*rwidth + rwidth/2+i)*bytes_per_pixel);
297 | }
298 | return 0;
299 | }
300 | SETERRQ (PETSC_ERR_ARG_OUTOFRANGE, "Field type not yet supported");
301 | }
302 |
303 |
304 | #undef __FUNCT__
305 | #define __FUNCT__ "render_rgb_local_2d"
306 |
307 | /*++++++++++++++++++++++++++++++++++++++
308 | Render data from global_array into local part of an RGB buffer. When running
309 | in parallel, these local buffers should be collected and layered to produce
310 | the full image.
311 |
312 | int render_rgb_local_2d Returns zero or an error code.
313 |
314 | guchar *rgb RGB buffer, or location in the RGB buffer, in which to render.
315 |
316 | int rwidth Total width of the section of the RGB buffer to draw.
317 |
318 | int rheight Total height of the section of the RGB buffer to draw.
319 |
320 | int rowskip Total width of the RGB buffer.
321 |
322 | int bytes_per_pixel Number of bytes per pixel in this RGB buffer (typically 3
323 | or 4).
324 |
325 | PetscScalar *global_array Local array of global vector data to render.
326 |
327 | int num_fields Number of field variables in the array.
328 |
329 | int display_field The (first) field we are rendering now.
330 |
331 | field_plot_type fieldtype The type of this field.
332 |
333 | PetscScalar *scale Array of minimum and maximum values to pass to the
334 | various pseudocolor functions; if NULL, call auto_scale to determine those
335 | values.
336 |
337 | int nx Width of the array.
338 |
339 | int ny Height of the array.
340 |
341 | int xs Starting
342 | +latex+$x$-coordinate
343 | +html+ <i>x</i>-coordinate
344 | of the local part of the global vector.
345 |
346 | int ys Starting
347 | +latex+$y$-coordinate
348 | +html+ <i>y</i>-coordinate
349 | of the local part of the global vector.
350 |
351 | int xm Width of the local part of the global vector.
352 |
353 | int ym Height of the local part of the global vector.
354 | ++++++++++++++++++++++++++++++++++++++*/
355 |
356 | int render_rgb_local_2d
357 | (guchar *rgb, int rwidth, int rheight, int rowskip, int bytes_per_pixel,
358 | PetscScalar *global_array, int num_fields, int display_field,
359 | field_plot_type fieldtype, PetscScalar *scale, int nx,int ny, int xs,int ys,
360 | int xm,int ym)
361 | {
362 | int ix,iy, ierr;
363 | PetscScalar local_scale[6];
364 |
365 | /* Determine default min and max if none are provided */
366 | if (scale == NULL)
367 | {
368 | scale = local_scale;
369 | ierr = auto_scale (global_array, xm*ym, num_fields, display_field,
370 | fieldtype, 2, scale); CHKERRQ (ierr);
371 | }
372 | DPRINTF ("Final scale: %g %g %g %g %g %g\n", scale[0],
373 | scale[1], scale[2], scale[3], scale[4], scale[5]);
374 |
375 | /* Do the rendering (note switch in inner loop, gotta get rid of that) */
376 | for (iy=rheight*ys/ny; iy<rheight*(ys+ym)/ny; iy++)
377 | for (ix=rwidth*xs/nx; ix<rwidth*(xs+xm)/nx; ix++)
378 | {
379 | int vecindex = (((rheight-iy-1)*ny/rheight)*nx +
380 | ix*nx/rwidth)*num_fields + display_field;
381 | guchar *pixel = rgb + bytes_per_pixel*(iy*rowskip + ix);
382 |
383 | switch (fieldtype)
384 | {
385 | case FIELD_SCALAR:
386 | case FIELD_SCALAR+1:
387 | {
388 | pseudocolor (global_array [vecindex], scale, pixel);
389 | break;
390 | }
391 | case FIELD_TERNARY:
392 | {
393 | pseudoternarycolor
394 | (global_array [vecindex], global_array [vecindex+1], scale,
395 | pixel);
396 | break;
397 | }
398 | case FIELD_VECTOR:
399 | case FIELD_VECTOR+1:
400 | {
401 | pseudohueintcolor
402 | (global_array [vecindex], global_array [vecindex+1], scale,
403 | pixel);
404 | break;
405 | }
406 | case FIELD_TENSOR_SHEAR:
407 | {
408 | pseudoshearcolor
409 | (global_array [vecindex], global_array [vecindex+1], scale,
410 | pixel);
411 | break;
412 | }
413 | default:
414 | SETERRQ (PETSC_ERR_ARG_OUTOFRANGE, "Field type not yet supported");
415 | }
416 | }
417 | return 0;
418 | }
419 |
420 |
421 | #undef __FUNCT__
422 | #define __FUNCT__ "render_rgb_local_3d"
423 |
424 | /*++++++++++++++++++++++++++++++++++++++
425 | Render triangle data into an RGB buffer. When called in parallel, the
426 | resulting images should be layered to give the complete picture. Zooming is
427 | done by adjusting the ratio of the dir vector to the right vector.
428 | +latex+\par
429 | +html+ <p>
430 | The coordinate transformation is pretty simple.
431 | +latex+A point $\vec{p}$ in space is transformed to $x,y$ on the screen by
432 | +latex+representing it as:
433 | +latex+\begin{equation}
434 | +latex+ \label{eq:transform}
435 | +latex+ \vec{p} = \vec{\imath} + a \vec{d} + ax \vec{r} + ay \vec{u},
436 | +latex+\end{equation}
437 | +latex+where $\vec{\imath}$ is the observer's point (passed in as {\tt eye}),
438 | +latex+$\vec{d}$ is the direction of observation (passed in as {\tt dir}),
439 | +latex+$\vec{r}$ is the rightward direction to the observer (passed in as
440 | +latex+{\tt right}), and $\vec{u}$ is the upward direction to the observer
441 | +latex+which is given the direction of $\vec{r}\times\vec{d}$ and the
442 | +latex+magnitude of $\vec{r}$.
443 | +html+ A point <u><i>p</i></u> in space is transformed to <i>x, y</i> on the
444 | +html+ screen by representing it as:
445 | +html+ <center><u><i>p</i></u> = <u><i>i</i></u> + <i>a <u>d</u></i> +
446 | +html+ <i>ax <u>r</u></i> + <i>ay <u>u</u></i>,</center>
447 | +html+ where <u><i>i</i></u> is the observer's point (passed in as
448 | +html+ <tt>eye</tt>), <u><i>d</i></u> is the direction of observation (passed
449 | +html+ in as <tt>dir</tt>), <u><i>r</i></u> is the rightward direction to the
450 | +html+ observer (passed in as <tt>right</tt>), and <u><i>u</i></u> is the
451 | +html+ upward direction to the observer which is given the direction of
452 | +html+ the cross product of <u><i>d</i></u> and <u><i>r</i></u> and the
453 | +html+ magnitude of <u><i>r</i></u>.
454 | +latex+\par
455 | +html+ <p>
456 | +latex+This system in equation \ref{eq:transform} can easily be solved for
457 | +latex+$x$ and $y$ by first making it into a matrix equation:
458 | +latex+\begin{equation}
459 | +latex+ \label{eq:matrix}
460 | +latex+ \left(\begin{array}{ccc} d_x&r_x&u_x \\ d_y&r_y&u_y \\ d_z&r_z&u_z
461 | +latex+ \end{array}\right)
462 | +latex+ \left(\begin{array}{c} a \\ ax \\ ay \end{array}\right) =
463 | +latex+ \left(\begin{array}{c} p_x-i_x \\ p_y-i_y \\ p_z-i_z
464 | +latex+ \end{array}\right).
465 | +latex+\end{equation}
466 | +latex+Calling the matrix $M$ the inverse of the matrix on the left side of
467 | +latex+equation \ref{eq:matrix} gives the result:
468 | +latex+\begin{eqnarray}
469 | +latex+ a&=& M_{00} (p_x-i_x) + M_{01} (p_y-i_y) + M_{02} (p_z-i_z), \\
470 | +latex+ x&=& \frac{1}{a}[M_{10}(p_x-i_x)+M_{11}(p_y-i_y)+M_{12}(p_z-i_z)],\\
471 | +latex+ y&=& \frac{1}{a}[M_{20}(p_x-i_x)+M_{21}(p_y-i_y)+M_{22}(p_z-i_z)].
472 | +latex+\end{eqnarray}
473 | +html+ This system can easily be solved for <i>x</i> and <i>y</i> by first
474 | +html+ making it into a matrix equation:
475 | +html+ <center>(<i><u>d</u> <u>r</u> <u>u</u></i>) (<i>a, ax, ay</i>) =
476 | +html+ <u><i>p</i></u> - <u><i>i</i></u>.</center>
477 | +html+ Calling the matrix <u><i>M</i></u> the inverse of the matrix on the
478 | +html+ left side gives the result:
479 | +html+ <center><table>
480 | +html+ <tr><td><i>a</i></td><td>=</td><td>
481 | +html+ <i>M</i><sub>00</sub> (<i>p<sub>x</sub></i>-<i>i<sub>x</sub></i>) +
482 | +html+ <i>M</i><sub>01</sub> (<i>p<sub>y</sub></i>-<i>i<sub>y</sub></i>) +
483 | +html+ <i>M</i><sub>02</sub> (<i>p<sub>z</sub></i>-<i>i<sub>z</sub></i>),
484 | +html+ </td></tr><tr><td><i>x</i></td><td>=</td><td>
485 | +html+ [<i>M</i><sub>10</sub> (<i>p<sub>x</sub></i>-<i>i<sub>x</sub></i>) +
486 | +html+ <i>M</i><sub>11</sub> (<i>p<sub>y</sub></i>-<i>i<sub>y</sub></i>) +
487 | +html+ <i>M</i><sub>12</sub> (<i>p<sub>z</sub></i>-<i>i<sub>z</sub></i>)]
488 | +html+ / <i>a</i>,</td></tr><tr><td><i>y</i></td><td>=</td><td>
489 | +html+ [<i>M</i><sub>00</sub> (<i>p<sub>x</sub></i>-<i>i<sub>x</sub></i>) +
490 | +html+ <i>M</i><sub>01</sub> (<i>p<sub>y</sub></i>-<i>i<sub>y</sub></i>) +
491 | +html+ <i>M</i><sub>02</sub> (<i>p<sub>z</sub></i>-<i>i<sub>z</sub></i>)]
492 | +html+ / <i>a</i>.</td></tr></table></center>
493 | +latex+\par
494 | +html+ <p>
495 | The triple product of the vector from the light source to the triangle
496 | centroid and the two vectors making up the triangle edges determines the
497 | cosine of the angle made by the triangle normal and the incident light, whose
498 | absolute value is used here to shade the triangle. At this point, the light
499 | source is taken as the observer location at
500 | +latex+``{\tt eye}'',
501 | +html+ "<tt>eye</tt>",
502 | but that can easily be modified to use one or more independent light sources.
503 |
504 | int render_rgb_local_3d Returns zero or an error code.
505 |
506 | guchar *rgb RGB buffer in which to render.
507 |
508 | int rwidth Total width of the RGB buffer.
509 |
510 | int rheight Total height of the RGB buffer.
511 |
512 | int bytes_per_pixel Number of bytes per pixel in this RGB buffer (typically 3
513 | or 4).
514 |
515 | int num_triangles Number of triangles to render.
516 |
517 | PetscScalar *vertices Table of coordinates (x1,y1,z1, x2,y2,z2, x3,y3,z3) and
518 | colors (RGBA 0-1) making thirteen values per triangle.
519 |
520 | PetscScalar *eye Point from where we're looking (x,y,z).
521 |
522 | PetscScalar *dir Direction we're looking (x,y,z).
523 |
524 | PetscScalar *right Rightward direction in physical space (x,y,z).
525 | ++++++++++++++++++++++++++++++++++++++*/
526 |
527 | int render_rgb_local_3d
528 | (guchar *rgb, int rwidth, int rheight, int bytes_per_pixel,
529 | int num_triangles, PetscScalar *vertices,
530 | PetscScalar *eye, PetscScalar *dir, PetscScalar *right)
531 | {
532 | #ifdef IMLIB2_EXISTS
533 | int *screen_coords, i;
534 | PetscScalar *shades, up[3], m00,m01,m02,m10,m11,m12,m20,m21,m22, a;
535 |
536 | if (bytes_per_pixel != 4)
537 | SETERRQ (PETSC_ERR_ARG_SIZ, "Imlib2 rendering requires 4 bytes per pixel");
538 |
539 | /* Allocate memory for the shades and integer coordinates */
540 | i = PetscMalloc (6*num_triangles * sizeof(int), &screen_coords); CHKERRQ (i);
541 | i = PetscMalloc (num_triangles * sizeof(PetscScalar), &shades); CHKERRQ (i);
542 |
543 | /* Calculate the "up" vector as described above, storing the magnitude of dir
544 | in a */
545 | a = sqrt (dir[0]*dir[0] + dir[1]*dir[1] + dir[2]*dir[2]);
546 | up[0] = (right[1]*dir[2] - right[2]*dir[1]) / a;
547 | up[1] = (right[2]*dir[0] - right[0]*dir[2]) / a;
548 | up[2] = (right[0]*dir[1] - right[1]*dir[0]) / a;
549 |
550 | /* Calculate the inverse of the dir, right, up columnar matrix, storing the
551 | determinant of the matrix in a */
552 | a = (dir[0] * (right[1]*up[2] - right[2]*up[1]) +
553 | dir[1] * (right[2]*up[0] - right[0]*up[2]) +
554 | dir[2] * (right[0]*up[1] - right[1]*up[0]));
555 | m00 = (right[1]*up[2] - right[2]*up[1]) / a;
556 | m01 = (right[2]*up[0] - right[0]*up[2]) / a;
557 | m02 = (right[0]*up[1] - right[1]*up[0]) / a;
558 | m10 = (up[1] *dir[2] - up[2] *dir[1]) / a;
559 | m11 = (up[2] *dir[0] - up[0] *dir[2]) / a;
560 | m12 = (up[0] *dir[1] - up[1] *dir[0]) / a;
561 | m20 = (dir[1] *right[2] - dir[2] *right[1]) / a;
562 | m21 = (dir[2] *right[0] - dir[0] *right[2]) / a;
563 | m22 = (dir[0] *right[1] - dir[1] *right[0]) / a;
564 |
565 | /* Loop over triangles, transforming their coordinates and setting shading */
566 | for (i=0; i<num_triangles; i++)
567 | {
568 | PetscScalar pmi00,pmi01,pmi02,pmi10,pmi11,pmi12,pmi20,pmi21,pmi22, c[3];
569 |
570 | /* Calculate coordinates relative to observer p-i (hence pmi) */
571 | pmi00 = vertices[13*i] - eye[0];
572 | pmi01 = vertices[13*i+1] - eye[1];
573 | pmi02 = vertices[13*i+2] - eye[2];
574 | pmi10 = vertices[13*i+3] - eye[0];
575 | pmi11 = vertices[13*i+4] - eye[1];
576 | pmi12 = vertices[13*i+5] - eye[2];
577 | pmi20 = vertices[13*i+6] - eye[0];
578 | pmi21 = vertices[13*i+7] - eye[1];
579 | pmi22 = vertices[13*i+8] - eye[2];
580 |
581 | /* Calculate x,y for point 0 */
582 | a = m00*pmi00 + m01*pmi01 + m02*pmi02;
583 | if (a < 0)
584 | screen_coords [6*i+0] = screen_coords [6*i+1] = -1;
585 | else
586 | {
587 | screen_coords [6*i+0] = rwidth/2 *
588 | (1. + (m10*pmi00 + m11*pmi01 + m12*pmi02) / a);
589 | screen_coords [6*i+1] = rheight/2 - rwidth/2 *
590 | (m20*pmi00 + m21*pmi01 + m22*pmi02) / a;
591 | }
592 |
593 | /* Calculate x,y for point 1 */
594 | a = m00*pmi10 + m01*pmi11 + m02*pmi12;
595 | if (a <= 0)
596 | screen_coords [6*i+2] = screen_coords [6*i+3] = -1;
597 | else
598 | {
599 | screen_coords [6*i+2] = rwidth/2 *
600 | (1. + (m10*pmi10 + m11*pmi11 + m12*pmi12) / a);
601 | screen_coords [6*i+3] = rheight/2 - rwidth/2 *
602 | (m20*pmi10 + m21*pmi11 + m22*pmi12) / a;
603 | }
604 |
605 | /* Calculate x,y for point 2 */
606 | a = m00*pmi20 + m01*pmi21 + m02*pmi22;
607 | if (a <= 0)
608 | screen_coords [6*i+4] = screen_coords [6*i+5] = -1;
609 | else
610 | {
611 | screen_coords [6*i+4] = rwidth/2 *
612 | (1. + (m10*pmi20 + m11*pmi21 + m12*pmi22) / a);
613 | screen_coords [6*i+5] = rheight/2 - rwidth/2 *
614 | (m20*pmi20 + m21*pmi21 + m22*pmi22) / a;
615 | }
616 |
617 | #ifdef DEBUG
618 | if (i==0)
619 | DPRINTF ("First triangle: (%g,%g,%g), (%g,%g,%g), (%g,%g,%g);\n"
620 | "2-D coordinates: (%d,%d), (%d,%d), (%d,%d)\n",
621 | vertices[0],vertices[1],vertices[2],vertices[3],vertices[4],
622 | vertices[5],vertices[6],vertices[7],vertices[8],
623 | screen_coords[0],screen_coords[1],screen_coords[3],
624 | screen_coords[3],screen_coords[4],screen_coords[5]);
625 | #endif
626 |
627 | /* Calculate relative triangle centroid, triangle arms, and shade */
628 | c[0] = (pmi00 + pmi10 + pmi20) / 3.;
629 | c[1] = (pmi01 + pmi11 + pmi21) / 3.;
630 | c[2] = (pmi02 + pmi12 + pmi22) / 3.;
631 | pmi10 -= pmi00; pmi11 -= pmi01; pmi12 -= pmi02;
632 | pmi20 -= pmi00; pmi21 -= pmi01; pmi22 -= pmi02;
633 | shades [i] = PetscAbsScalar (c[0] * (pmi11*pmi22 - pmi12*pmi21) +
634 | c[1] * (pmi12*pmi20 - pmi10*pmi22) +
635 | c[2] * (pmi10*pmi21 - pmi11*pmi20)) /
636 | sqrt (c[0]*c[0] + c[1]*c[1] + c[2]*c[2]) /
637 | sqrt (pmi10*pmi10 + pmi11*pmi11 + pmi12*pmi12) /
638 | sqrt (pmi20*pmi20 + pmi21*pmi21 + pmi22*pmi22);
639 | }
640 |
641 | /* Do the Imlib2 rendering */
642 | imlib2_render_triangles ((DATA32 *) rgb, rwidth, rheight, num_triangles,
643 | screen_coords, vertices+9, 13, shades, 1);
644 |
645 | i = PetscFree (shades); CHKERRQ (i);
646 | i = PetscFree (screen_coords); CHKERRQ (i);
647 | #else
648 | SETERRQ (PETSC_ERR_SUP, "3-D rendering requires Imlib2");
649 | #endif
650 | }