BRL-CAD
refract.c
Go to the documentation of this file.
1 /* R E F R A C T . C
2  * BRL-CAD
3  *
4  * Copyright (c) 1985-2014 United States Government as represented by
5  * the U.S. Army Research Laboratory.
6  *
7  * This library is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public License
9  * version 2.1 as published by the Free Software Foundation.
10  *
11  * This library is distributed in the hope that it will be useful, but
12  * WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14  * Lesser General Public License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with this file; see the file named COPYING for more
18  * information.
19  */
20 /** @file liboptical/refract.c
21  *
22  */
23 
24 #include "common.h"
25 
26 #include <stdio.h>
27 #include <string.h>
28 #include <math.h>
29 
30 #include "bu/sort.h"
31 #include "bu/parallel.h"
32 #include "vmath.h"
33 #include "mater.h"
34 #include "raytrace.h"
35 #include "optical.h"
36 #include "plot3.h"
37 
38 
39 extern int viewshade(struct application *ap,
40  register const struct partition *pp,
41  register struct shadework *swp);
42 
43 
44 int max_ireflect = 5; /* Maximum internal reflection level */
45 int max_bounces = 5; /* Maximum recursion level */
46 
47 #define MSG_PROLOGUE 20 /* # initial messages to see */
48 #define MSG_INTERVAL 4000 /* message interval thereafter */
49 
50 #define RI_AIR 1.0 /* Refractive index of air */
51 
52 #define AIR_GAP_TOL 0.01 /* Max permitted air gap for RI tracking */
53 
54 
55 #ifdef RT_MULTISPECTRAL
56 #include "spectrum.h"
57 extern struct bn_tabdata *background;
58 #else
59 extern vect_t background;
60 #endif
61 
62 HIDDEN int
63 rr_miss(struct application *ap)
64 {
65  RT_AP_CHECK(ap);
66  return 1; /* treat as escaping ray */
67 }
68 
69 
70 /*
71  * This routine is called when an internal reflection ray hits something
72  * (which is ordinarily the case).
73  *
74  * Generally, there will be one or two partitions on the hit list.
75  * The values for pt_outhit for the second partition should not be used,
76  * as a_onehit was set to 3, getting a maximum of 3 valid hit points.
77  *
78  * Explicit Returns -
79  * 0 dreadful internal error
80  * 1 treat as escaping ray & reshoot
81  * 2 Proper exit point determined, with Implicit Returns:
82  * a_uvec exit Point
83  * a_vvec exit Normal (inward pointing)
84  * a_refrac_index RI of *next* material
85  */
86 HIDDEN int
87 rr_hit(struct application *ap, struct partition *PartHeadp, struct seg *UNUSED(segp))
88 {
89  register struct partition *pp;
90  register struct hit *hitp;
91  register struct soltab *stp;
92  struct partition *psave = (struct partition *)NULL;
93  struct shadework sw;
94  struct application appl;
95  int ret;
96 
97  RT_AP_CHECK(ap);
98 
99  RT_APPLICATION_INIT(&appl);
100 
101  for (pp=PartHeadp->pt_forw; pp != PartHeadp; pp = pp->pt_forw)
102  if (pp->pt_outhit->hit_dist > 0.0) break;
103  if (pp == PartHeadp) {
105  bu_log("rr_hit: %d, %d no hit out front?\n",
106  ap->a_x, ap->a_y);
107  ret = 0; /* error */
108  goto out;
109  }
110  ret = 1; /* treat as escaping ray */
111  goto out;
112  }
113 
114  /*
115  * Ensure that the partition we are given is part of the same
116  * region that we started in. When the internal reflection
117  * is happening very near an edge or corner, this is not always
118  * the case, and either (a) a small sliver of some other region
119  * is found to be in the way, or (b) the ray completely misses the
120  * region that it started in, although not by much.
121  */
122  psave = pp;
123  if (R_DEBUG&RDEBUG_REFRACT) bu_log("rr_hit(%s)\n", psave->pt_regionp->reg_name);
124  for (; pp != PartHeadp; pp = pp->pt_forw)
125  if (pp->pt_regionp == (struct region *)(ap->a_uptr)) break;
126  if (pp == PartHeadp) {
127  if (R_DEBUG&(RDEBUG_SHOWERR|RDEBUG_REFRACT)) {
128  bu_log("rr_hit: %d, %d Ray internal to %s landed unexpectedly in %s\n",
129  ap->a_x, ap->a_y,
130  ((struct region *)(ap->a_uptr))->reg_name,
131  psave->pt_regionp->reg_name);
132  ret = 0; /* error */
133  goto out;
134  }
135  ret = 1; /* treat as escaping ray */
136  goto out;
137  }
138 
139  /*
140  * At one time, this was a check for pp->pt_inhit->hit_dist
141  * being NEAR zero. That was a mistake, because we may have
142  * been at the edge of a subtracted out center piece when
143  * internal reflection happened, except that floating point
144  * error (being right on the surface of the interior solid)
145  * prevented us from "seeing" that solid on the next ray,
146  * causing our ray endpoints to be quite far from the starting
147  * point, yet with the ray still validly inside the glass region.
148  *
149  * There is a major problem if the entry point
150  * is further ahead than the firing point, i.e., >0.
151  *
152  * Because this error has not yet been encountered, it is
153  * considered dreadful. Some recovery may be possible.
154  *
155  * For now, this seems to happen when a reflected ray starts outside
156  * the glass and doesn't even intersect the glass, so treat it as
157  * an escaping ray.
158  */
159 
160  if (pp->pt_inhit->hit_dist > 10) {
161  stp = pp->pt_inseg->seg_stp;
162  if (R_DEBUG&RDEBUG_REFRACT)
163  bu_log("rr_hit: %d, %d %s inhit %g > 10.0! (treating as escaping ray)\n",
164  ap->a_x, ap->a_y,
165  pp->pt_regionp->reg_name,
166  pp->pt_inhit->hit_dist);
167  ret = 1; /* treat as escaping ray */
168  goto out;
169  }
170 
171  /*
172  * If there is a very small crack in the glass, perhaps formed
173  * by a small error when taking the Union of two solids,
174  * attempt to find the real exit point.
175  * NOTE that this is usually taken care of inside librt
176  * in the bool_weave code, but it is inexpensive to check for it
177  * here. If this case is detected, push on, and log it.
178  * This code is not expected to be needed.
179  */
180  while (pp->pt_forw != PartHeadp) {
181  register fastf_t d;
182  d = pp->pt_forw->pt_inhit->hit_dist - pp->pt_outhit->hit_dist;
183  if (!NEAR_ZERO(d, AIR_GAP_TOL))
184  break;
185  if (pp->pt_forw->pt_regionp != pp->pt_regionp)
186  break;
187  if (R_DEBUG&(RDEBUG_SHOWERR|RDEBUG_REFRACT)) bu_log(
188  "rr_hit: %d, %d fusing small crack in glass %s\n",
189  ap->a_x, ap->a_y,
190  pp->pt_regionp->reg_name);
191  pp = pp->pt_forw;
192  }
193 
194  hitp = pp->pt_outhit;
195  stp = pp->pt_outseg->seg_stp;
196  if (hitp->hit_dist >= INFINITY) {
197  bu_log("rr_hit: %d, %d infinite glass (%g, %g) %s\n",
198  ap->a_x, ap->a_y,
199  pp->pt_inhit->hit_dist, hitp->hit_dist,
200  pp->pt_regionp->reg_name);
201  ret = 0; /* dreadful error */
202  goto out;
203  }
204  VJOIN1(hitp->hit_point, ap->a_ray.r_pt,
205  hitp->hit_dist, ap->a_ray.r_dir);
206  RT_HIT_NORMAL(ap->a_vvec, hitp, stp, &(ap->a_ray), pp->pt_outflip);
207 
208  /* For refraction, want exit normal to point inward. */
209  VREVERSE(ap->a_vvec, ap->a_vvec);
210  VMOVE(ap->a_uvec, hitp->hit_point);
211  ap->a_cumlen += (hitp->hit_dist - pp->pt_inhit->hit_dist);
212 
213  ap->a_refrac_index = RI_AIR; /* Default medium: air */
214 
215  /*
216  * Look ahead, and see if there is more glass to come.
217  * If so, obtain its refractive index, to enable correct
218  * calculation of the departing refraction angle.
219  */
220  if (pp->pt_forw != PartHeadp) {
221  register fastf_t d;
222  d = pp->pt_forw->pt_inhit->hit_dist - hitp->hit_dist;
223  if (NEAR_ZERO(d, AIR_GAP_TOL)) {
224  /*
225  * Make a private copy of the application struct,
226  * because viewshade() may change various fields.
227  */
228  appl = *ap; /* struct copy */
229 
230  memset((char *)&sw, 0, sizeof(sw));
231  sw.sw_transmit = sw.sw_reflect = 0.0;
232 
233  /* Set default in case shader doesn't fill this in. */
234  sw.sw_refrac_index = RI_AIR;
235 
236  /* Set special flag so that we get only shader
237  * parameters (refractive index, in this case).
238  * We don't even care about transmitted energy.
239  */
240  sw.sw_xmitonly = 2;
241  sw.sw_inputs = 0; /* no fields filled yet */
242 #ifdef RT_MULTISPECTRAL
243  sw.msw_color = bn_tabdata_get_constval(1.0, spectrum);
244  sw.msw_basecolor = bn_tabdata_get_constval(1.0, spectrum);
245 #else
246  VSETALL(sw.sw_color, 1);
247  VSETALL(sw.sw_basecolor, 1);
248 #endif
249 
250  if (R_DEBUG&(RDEBUG_SHADE|RDEBUG_REFRACT))
251  bu_log("rr_hit calling viewshade to discover refractive index\n");
252 
253  (void)viewshade(&appl, pp->pt_forw, &sw);
254 
255 #ifdef RT_MULTISPECTRAL
256  bu_free(sw.msw_color, "sw.msw_color");
257  bu_free(sw.msw_basecolor, "sw.msw_basecolor");
258 #endif
259 
260  if (R_DEBUG&(RDEBUG_SHADE|RDEBUG_REFRACT))
261  bu_log("rr_hit refractive index = %g\n", sw.sw_refrac_index);
262 
263  if (sw.sw_transmit > 0) {
264  ap->a_refrac_index = sw.sw_refrac_index;
265  if (R_DEBUG&RDEBUG_SHADE) {
266  bu_log("rr_hit a_refrac_index=%g (trans=%g)\n",
267  ap->a_refrac_index,
268  sw.sw_transmit);
269  }
270  ret= 3; /* OK -- more glass follows */
271  goto out;
272  }
273  }
274  }
275  ret = 2; /* OK -- no more glass */
276 out:
277  if (R_DEBUG&RDEBUG_REFRACT) bu_log("rr_hit(%s) return=%d\n",
278  psave ? psave->pt_regionp->reg_name : "",
279  ret);
280  return ret;
281 }
282 
283 
284 /*
285  * Compute the refracted ray 'v_2' from the incident ray 'v_1' with
286  * the refractive indices 'ri_2' and 'ri_1' respectively.
287  * Using Schnell's Law:
288  *
289  * theta_1 = angle of v_1 with surface normal
290  * theta_2 = angle of v_2 with reversed surface normal
291  * ri_1 * sin(theta_1) = ri_2 * sin(theta_2)
292  *
293  * sin(theta_2) = ri_1/ri_2 * sin(theta_1)
294  *
295  * The above condition is undefined for ri_1/ri_2 * sin(theta_1)
296  * being greater than 1, and this represents the condition for total
297  * reflection, the 'critical angle' is the angle theta_1 for which
298  * ri_1/ri_2 * sin(theta_1) equals 1.
299  *
300  * Returns TRUE if refracted, FALSE if reflected.
301  *
302  * Note: output (v_2) can be same storage as an input.
303  */
304 HIDDEN int
305 rr_refract(vect_t v_1, vect_t norml, double ri_1, double ri_2, vect_t v_2)
306 {
307  vect_t w, u;
308  fastf_t beta;
309 
310  if (NEAR_ZERO(ri_1, 0.0001) || NEAR_ZERO(ri_2, 0.0001)) {
311  bu_log("rr_refract:ri1=%g, ri2=%g\n", ri_1, ri_2);
312  beta = 1;
313  } else {
314  beta = ri_1/ri_2; /* temp */
315  if (beta > 10000) {
316  bu_log("rr_refract: beta=%g\n", beta);
317  beta = 1000;
318  }
319  }
320  VSCALE(w, v_1, beta);
321  VCROSS(u, w, norml);
322 
323  /*
324  * |w X norml| = |w||norml| * sin(theta_1)
325  * |u| = ri_1/ri_2 * sin(theta_1) = sin(theta_2)
326  */
327  if ((beta = VDOT(u, u)) > 1.0) {
328  /* Past critical angle, total reflection.
329  * Calculate reflected (bounced) incident ray.
330  */
331  if (R_DEBUG&RDEBUG_REFRACT) bu_log("rr_refract: reflected. ri1=%g ri2=%g beta=%g\n",
332  ri_1, ri_2, beta);
333  VREVERSE(u, v_1);
334  beta = 2 * VDOT(u, norml);
335  VSCALE(w, norml, beta);
336  VSUB2(v_2, w, u);
337  return 0; /* reflected */
338  } else {
339  /*
340  * 1 - beta = 1 - sin(theta_2)^^2
341  * = cos(theta_2)^^2.
342  * beta = -1.0 * cos(theta_2) - Dot(w, norml).
343  */
344  if (R_DEBUG&RDEBUG_REFRACT) bu_log("rr_refract: refracted. ri1=%g ri2=%g beta=%g\n",
345  ri_1, ri_2, beta);
346  beta = -sqrt(1.0 - beta) - VDOT(w, norml);
347  VSCALE(u, norml, beta);
348  VADD2(v_2, w, u);
349  return 1; /* refracted */
350  }
351  /* NOTREACHED */
352 }
353 
354 
355 int
356 rr_render(register struct application *ap,
357  const struct partition *pp,
358  struct shadework *swp)
359 {
360  struct application sub_ap;
361  vect_t work;
362  vect_t incident_dir;
363  fastf_t shader_fract;
364  fastf_t reflect;
365  fastf_t transmit;
366 
367 #ifdef RT_MULTISPECTRAL
368  struct bn_tabdata *ms_filter_color = BN_TABDATA_NULL;
369  struct bn_tabdata *ms_shader_color = BN_TABDATA_NULL;
370  struct bn_tabdata *ms_reflect_color = BN_TABDATA_NULL;
371  struct bn_tabdata *ms_transmit_color = BN_TABDATA_NULL;
372 #else
373  vect_t filter_color;
374  vect_t shader_color;
375  vect_t reflect_color;
376  vect_t transmit_color;
377 #endif
378 
379  fastf_t attenuation;
380  vect_t to_eye;
381  int code;
382 
383  RT_AP_CHECK(ap);
384 
385  RT_APPLICATION_INIT(&sub_ap);
386 
387 #ifdef RT_MULTISPECTRAL
388  sub_ap.a_spectrum = BN_TABDATA_NULL;
389  ms_reflect_color = bn_tabdata_get_constval(0.0, spectrum);
390 #endif
391 
392  /*
393  * sw_xmitonly is set primarily for light visibility rays.
394  * Need to compute (partial) transmission through to the light,
395  * or procedural shaders won't be able to cast shadows
396  * and light won't be able to get through glass
397  * (including "stained glass" and "filter glass").
398  *
399  * On the other hand, light visibility rays shouldn't be refracted,
400  * it is pointless to shoot at where the light isn't.
401  */
402  if (swp->sw_xmitonly) {
403  /* Caller wants transmission term only, don't fire reflected rays */
404  transmit = swp->sw_transmit + swp->sw_reflect; /* Don't loose energy */
405  reflect = 0;
406  } else {
407  reflect = swp->sw_reflect;
408  transmit = swp->sw_transmit;
409  }
410  if (R_DEBUG&RDEBUG_REFRACT) {
411  bu_log("rr_render(%s) START: lvl=%d reflect=%g, transmit=%g, xmitonly=%d\n",
412  pp->pt_regionp->reg_name,
413  ap->a_level,
414  reflect, transmit,
415  swp->sw_xmitonly);
416  }
417  if (reflect <= 0 && transmit <= 0)
418  goto out;
419 
420  if (ap->a_level > max_bounces) {
421  /* Nothing more to do for this ray */
422  static long count = 0; /* Not PARALLEL, should be OK */
423 
424  if ((R_DEBUG&(RDEBUG_SHOWERR|RDEBUG_REFRACT)) && (
425  count++ < MSG_PROLOGUE ||
426  (count%MSG_INTERVAL) == 3
427  )) {
428  bu_log("rr_render: %d, %d MAX BOUNCES=%d: %s\n",
429  ap->a_x, ap->a_y,
430  ap->a_level,
431  pp->pt_regionp->reg_name);
432  }
433 
434  /*
435  * Return the basic color of the object, ignoring the
436  * the fact that it is supposed to be
437  * filtering or reflecting light here.
438  * This is much better than returning just black,
439  * but something better might be done.
440  */
441 #ifdef RT_MULTISPECTRAL
442  BN_CK_TABDATA(swp->msw_color);
443  BN_CK_TABDATA(swp->msw_basecolor);
444  bn_tabdata_copy(swp->msw_color, swp->msw_basecolor);
445 #else
446  VMOVE(swp->sw_color, swp->sw_basecolor);
447 #endif
448  ap->a_cumlen += pp->pt_inhit->hit_dist;
449  goto out;
450  }
451 #ifdef RT_MULTISPECTRAL
452  BN_CK_TABDATA(swp->msw_basecolor);
453  ms_filter_color = bn_tabdata_dup(swp->msw_basecolor);
454 
455 #else
456  VMOVE(filter_color, swp->sw_basecolor);
457 #endif
458 
459  if ((swp->sw_inputs & (MFI_HIT|MFI_NORMAL)) != (MFI_HIT|MFI_NORMAL))
460  shade_inputs(ap, pp, swp, MFI_HIT|MFI_NORMAL);
461 
462  /*
463  * If this ray is being fired from the exit point of
464  * an object, and is directly entering another object,
465  * (i.e., there is no intervening air-gap), and
466  * the two refractive indices match, then do not fire a
467  * reflected ray -- just take the transmission contribution.
468  * This is important, e.g., for glass gun tubes projecting
469  * through a glass armor plate. :-)
470  */
472  && ZERO(ap->a_refrac_index - swp->sw_refrac_index))
473  {
474  transmit += reflect;
475  reflect = 0;
476  }
477 
478  /*
479  * Diminish base color appropriately, and add in
480  * contributions from mirror reflection & transparency
481  */
482  shader_fract = 1 - (reflect + transmit);
483  if (shader_fract < 0) {
484  shader_fract = 0;
485  } else if (shader_fract >= 1) {
486  goto out;
487  }
488  if (R_DEBUG&RDEBUG_REFRACT) {
489  bu_log("rr_render: lvl=%d start shader=%g, reflect=%g, transmit=%g %s\n",
490  ap->a_level,
491  shader_fract, reflect, transmit,
492  pp->pt_regionp->reg_name);
493  }
494 #ifdef RT_MULTISPECTRAL
495  BN_GET_TABDATA(ms_shader_color, swp->msw_color->table);
496  bn_tabdata_scale(ms_shader_color, swp->msw_color, shader_fract);
497 #else
498  VSCALE(shader_color, swp->sw_color, shader_fract);
499 #endif
500 
501  /*
502  * Compute transmission through an object.
503  * There may be a mirror reflection, which will be handled
504  * by the reflection code later
505  */
506  if (transmit > 0) {
507  if (R_DEBUG&RDEBUG_REFRACT) {
508  bu_log("rr_render: lvl=%d transmit=%g. Calculate refraction at entrance to %s.\n",
509  ap->a_level, transmit,
510  pp->pt_regionp->reg_name);
511  }
512  /*
513  * Calculate refraction at entrance.
514  */
515  sub_ap = *ap; /* struct copy */
516 #ifdef RT_MULTISPECTRAL
517  sub_ap.a_spectrum = bn_tabdata_dup((struct bn_tabdata *)ap->a_spectrum);
518 #endif
519  sub_ap.a_level = 0; /* # of internal reflections */
520  sub_ap.a_cumlen = 0; /* distance through the glass */
521  sub_ap.a_user = -1; /* sanity */
522  sub_ap.a_rbeam = ap->a_rbeam + swp->sw_hit.hit_dist * ap->a_diverge;
523  sub_ap.a_diverge = 0.0;
524  sub_ap.a_uptr = (void *)(pp->pt_regionp);
525  VMOVE(sub_ap.a_ray.r_pt, swp->sw_hit.hit_point);
526  VMOVE(incident_dir, ap->a_ray.r_dir);
527 
528  /* If there is an air gap, reset ray's RI to air */
529  if (pp->pt_inhit->hit_dist > AIR_GAP_TOL)
530  sub_ap.a_refrac_index = RI_AIR;
531 
532  if (!ZERO(sub_ap.a_refrac_index - swp->sw_refrac_index)
533  && !rr_refract(incident_dir, /* input direction */
534  swp->sw_hit.hit_normal, /* exit normal */
535  sub_ap.a_refrac_index, /* current RI */
536  swp->sw_refrac_index, /* next RI */
537  sub_ap.a_ray.r_dir /* output direction */
538  ))
539  {
540  /*
541  * Ray was mirror reflected back outside solid.
542  * Just add contribution to reflection,
543  * and quit.
544  */
545  reflect += transmit;
546  transmit = 0;
547 #ifdef RT_MULTISPECTRAL
548  ms_transmit_color = bn_tabdata_get_constval(0.0, spectrum);
549 #else
550  VSETALL(transmit_color, 0);
551 #endif
552  if (R_DEBUG&RDEBUG_REFRACT) {
553  bu_log("rr_render: lvl=%d change xmit into reflection %s\n",
554  ap->a_level,
555  pp->pt_regionp->reg_name);
556  }
557  goto do_reflection;
558  }
559  if (R_DEBUG&RDEBUG_REFRACT) {
560  bu_log("rr_render: lvl=%d begin transmission through %s.\n",
561  ap->a_level,
562  pp->pt_regionp->reg_name);
563  }
564 
565  /*
566  * Find new exit point from the inside.
567  * We will iterate, but not recurse, due to the special
568  * (non-recursing) hit and miss routines used here for
569  * internal reflection.
570  *
571  * a_onehit is set to 3, so that where possible,
572  * rr_hit() will be given three accurate hit points:
573  * the entry and exit points of this glass region,
574  * and the entry point into the next region.
575  * This permits calculation of the departing
576  * refraction angle based on the RI of the current and
577  * *next* regions along the ray.
578  */
579  sub_ap.a_purpose = "rr first glass transmission ray";
580  sub_ap.a_flag = 0;
581  do_inside:
582  sub_ap.a_hit = rr_hit;
583  sub_ap.a_miss = rr_miss;
584  sub_ap.a_logoverlap = ap->a_logoverlap;
585  sub_ap.a_onehit = 3;
586  sub_ap.a_rbeam = ap->a_rbeam + swp->sw_hit.hit_dist * ap->a_diverge;
587  sub_ap.a_diverge = 0.0;
588  switch (code = rt_shootray(&sub_ap)) {
589  case 3:
590  /* More glass to come.
591  * uvec=exit_pt, vvec=N, a_refrac_index = next RI.
592  */
593  break;
594  case 2:
595  /* No more glass to come.
596  * uvec=exit_pt, vvec=N, a_refrac_index = next RI.
597  */
598  break;
599  case 1:
600  /* Treat as escaping ray */
601  if (R_DEBUG&RDEBUG_REFRACT)
602  bu_log("rr_refract: Treating as escaping ray\n");
603  goto do_exit;
604  case 0:
605  default:
606  /* Dreadful error */
607 #ifdef RT_MULTISPECTRAL
608  bu_bomb("rr_refract: Stuck in glass. Very green pixel, unsupported in multi-spectral mode\n");
609 #else
610  VSET(swp->sw_color, 0, 99, 0); /* very green */
611 #endif
612  goto out; /* abandon hope */
613  }
614 
615  if (R_DEBUG&RDEBUG_REFRACT) {
616  bu_log("rr_render: calculating refraction @ exit from %s (green)\n", pp->pt_regionp->reg_name);
617  bu_log("Start point to exit point:\n\
618 vdraw open rr;vdraw params c 00ff00; vdraw write n 0 %g %g %g; vdraw wwrite n 1 %g %g %g; vdraw send\n",
619  V3ARGS(sub_ap.a_ray.r_pt),
620  V3ARGS(sub_ap.a_uvec));
621  }
622  /* NOTE: rr_hit returns EXIT Point in sub_ap.a_uvec,
623  * and returns EXIT Normal in sub_ap.a_vvec,
624  * and returns next RI in sub_ap.a_refrac_index
625  */
626  if (R_DEBUG&RDEBUG_RAYWRITE) {
627  wraypts(sub_ap.a_ray.r_pt,
628  sub_ap.a_ray.r_dir,
629  sub_ap.a_uvec,
630  2, ap, stdout); /* 2 = ?? */
631  }
632  if (R_DEBUG&RDEBUG_RAYPLOT) {
633  /* plotfp */
635  pl_color(stdout, 0, 255, 0);
636  pdv_3line(stdout,
637  sub_ap.a_ray.r_pt,
638  sub_ap.a_uvec);
640  }
641  /* Advance. Exit point becomes new start point */
642  VMOVE(sub_ap.a_ray.r_pt, sub_ap.a_uvec);
643  VMOVE(incident_dir, sub_ap.a_ray.r_dir);
644 
645  /*
646  * Calculate refraction at exit point.
647  * Use "look ahead" RI value from rr_hit.
648  */
649  if (!ZERO(sub_ap.a_refrac_index - swp->sw_refrac_index)
650  && !rr_refract(incident_dir, /* input direction */
651  sub_ap.a_vvec, /* exit normal */
652  swp->sw_refrac_index, /* current RI */
653  sub_ap.a_refrac_index, /* next RI */
654  sub_ap.a_ray.r_dir /* output direction */
655  ))
656  {
657  static long count = 0; /* not PARALLEL, should be OK */
658 
659  /* Reflected internally -- keep going */
660  if ((++sub_ap.a_level) <= max_ireflect) {
661  sub_ap.a_purpose = "rr reflected internal ray, probing for glass exit point";
662  sub_ap.a_flag = 0;
663  goto do_inside;
664  }
665 
666  /*
667  * Internal Reflection limit exceeded -- just let
668  * the ray escape, continuing on current course.
669  * This will cause some energy from somewhere in the
670  * scene to be received through this glass,
671  * which is much better than just returning
672  * grey or black, as before.
673  */
674  if ((R_DEBUG&(RDEBUG_SHOWERR|RDEBUG_REFRACT)) && (
675  count++ < MSG_PROLOGUE ||
676  (count%MSG_INTERVAL) == 3
677  )) {
678  bu_log("rr_render: %d, %d Int.reflect=%d: %s lvl=%d\n",
679  sub_ap.a_x, sub_ap.a_y,
680  sub_ap.a_level,
681  pp->pt_regionp->reg_name,
682  ap->a_level);
683  }
684  VMOVE(sub_ap.a_ray.r_dir, incident_dir);
685  goto do_exit;
686  }
687  do_exit:
688  /*
689  * Compute internal spectral transmittance.
690  * Bouger's law. pg 30 of "color science"
691  *
692  * Apply attenuation factor due to thickness of the glass.
693  * sw_extinction is in terms of fraction of light absorbed
694  * per linear meter of glass. a_cumlen is in mm.
695  */
696 /* XXX extinction should be a spectral curve, not scalar */
697  if (swp->sw_extinction > 0 && sub_ap.a_cumlen > 0) {
698  attenuation = pow(10.0, -1.0e-3 * sub_ap.a_cumlen *
699  swp->sw_extinction);
700  } else {
701  attenuation = 1;
702  }
703 
704  /*
705  * Process the escaping refracted ray.
706  * This is the only place we might recurse dangerously,
707  * so we are careful to use our caller's recursion level+1.
708  * NOTE: point & direction already filled in
709  */
710  sub_ap.a_hit = ap->a_hit;
711  sub_ap.a_miss = ap->a_miss;
712  sub_ap.a_logoverlap = ap->a_logoverlap;
713  sub_ap.a_onehit = ap->a_onehit;
714  sub_ap.a_level = ap->a_level+1;
715  sub_ap.a_uptr = ap->a_uptr;
716  sub_ap.a_rbeam = ap->a_rbeam + swp->sw_hit.hit_dist * ap->a_diverge;
717  sub_ap.a_diverge = 0.0;
718  if (code == 3) {
719  sub_ap.a_purpose = "rr recurse on next glass";
720  sub_ap.a_flag = 0;
721  } else {
722  sub_ap.a_purpose = "rr recurse on escaping internal ray";
723  sub_ap.a_flag = 1;
724  sub_ap.a_onehit = sub_ap.a_onehit > -3 ? -3 : sub_ap.a_onehit;
725  }
726  /* sub_ap.a_refrac_index was set to RI of next material by rr_hit().
727  */
728  sub_ap.a_cumlen = 0;
729  (void) rt_shootray(&sub_ap);
730 
731  /* a_user has hit/miss flag! */
732  if (sub_ap.a_user == 0) {
733 #ifdef RT_MULTISPECTRAL
734  ms_transmit_color = bn_tabdata_dup(background);
735 #else
736  VMOVE(transmit_color, background);
737 #endif
738  sub_ap.a_cumlen = 0;
739  } else {
740 #ifdef RT_MULTISPECTRAL
741  ms_transmit_color = bn_tabdata_dup(sub_ap.a_spectrum);
742 #else
743  VMOVE(transmit_color, sub_ap.a_color);
744 #endif
745  }
746  transmit *= attenuation;
747 #ifdef RT_MULTISPECTRAL
748  bn_tabdata_mul(ms_transmit_color, ms_filter_color, ms_transmit_color);
749 #else
750  VELMUL(transmit_color, filter_color, transmit_color);
751 #endif
752  if (R_DEBUG&RDEBUG_REFRACT) {
753  bu_log("rr_render: lvl=%d end of xmit through %s\n",
754  ap->a_level,
755  pp->pt_regionp->reg_name);
756  }
757  } else {
758 #ifdef RT_MULTISPECTRAL
759  ms_transmit_color = bn_tabdata_get_constval(0.0, spectrum);
760 #else
761  VSETALL(transmit_color, 0);
762 #endif
763  }
764 
765  /*
766  * Handle any reflection, including mirror reflections
767  * detected by the transmission code, above.
768  */
769 do_reflection:
770 #ifdef RT_MULTISPECTRAL
771  if (sub_ap.a_spectrum) {
772  bu_free(sub_ap.a_spectrum, "rr_render: sub_ap.a_spectrum bn_tabdata*");
773  sub_ap.a_spectrum = BN_TABDATA_NULL;
774  }
775 #endif
776  if (reflect > 0) {
777  register fastf_t f;
778 
779  /* Mirror reflection */
780  if (R_DEBUG&RDEBUG_REFRACT)
781  bu_log("rr_render: calculating mirror reflection off of %s\n", pp->pt_regionp->reg_name);
782  sub_ap = *ap; /* struct copy */
783 #ifdef RT_MULTISPECTRAL
784  sub_ap.a_spectrum = bn_tabdata_dup((struct bn_tabdata *)ap->a_spectrum);
785 #endif
786  sub_ap.a_rbeam = ap->a_rbeam + swp->sw_hit.hit_dist * ap->a_diverge;
787  sub_ap.a_diverge = 0.0;
788  sub_ap.a_level = ap->a_level+1;
789  sub_ap.a_onehit = -1; /* Require at least one non-air hit */
790  VMOVE(sub_ap.a_ray.r_pt, swp->sw_hit.hit_point);
791  VREVERSE(to_eye, ap->a_ray.r_dir);
792  f = 2 * VDOT(to_eye, swp->sw_hit.hit_normal);
793  VSCALE(work, swp->sw_hit.hit_normal, f);
794  /* I have been told this has unit length */
795  VSUB2(sub_ap.a_ray.r_dir, work, to_eye);
796  sub_ap.a_purpose = "rr reflected ray";
797  sub_ap.a_flag = 0;
798 
799  if (R_DEBUG&(RDEBUG_RAYPLOT|RDEBUG_REFRACT)) {
800  point_t endpt;
801  /* Plot the surface normal -- green/blue */
802  /* plotfp */
803  f = sub_ap.a_rt_i->rti_radius * 0.02;
804  VJOIN1(endpt, sub_ap.a_ray.r_pt,
805  f, swp->sw_hit.hit_normal);
806  if (R_DEBUG&RDEBUG_RAYPLOT) {
808  pl_color(stdout, 0, 255, 255);
809  pdv_3line(stdout, sub_ap.a_ray.r_pt, endpt);
811  }
812  bu_log("Surface normal for reflection:\n\
813 vdraw open rrnorm;vdraw params c 00ffff;vdraw write n 0 %g %g %g;vdraw write n 1 %g %g %g;vdraw send\n",
814  V3ARGS(sub_ap.a_ray.r_pt),
815  V3ARGS(endpt));
816 
817  }
818 
819  (void)rt_shootray(&sub_ap);
820 
821  /* a_user has hit/miss flag! */
822  if (sub_ap.a_user == 0) {
823  /* MISS */
824 #ifdef RT_MULTISPECTRAL
825  bn_tabdata_copy(ms_reflect_color, background);
826 #else
827  VMOVE(reflect_color, background);
828 #endif
829  } else {
830  ap->a_cumlen += sub_ap.a_cumlen;
831 #ifdef RT_MULTISPECTRAL
832  bn_tabdata_copy(ms_reflect_color, sub_ap.a_spectrum);
833 #else
834  VMOVE(reflect_color, sub_ap.a_color);
835 #endif
836  }
837  } else {
838 #ifdef RT_MULTISPECTRAL
839  bn_tabdata_constval(ms_reflect_color, 0.0);
840 #else
841  VSETALL(reflect_color, 0);
842 #endif
843  }
844 
845  /*
846  * Collect the contributions to the final color
847  */
848 #ifdef RT_MULTISPECTRAL
849  bn_tabdata_join2(swp->msw_color, ms_shader_color,
850  reflect, ms_reflect_color,
851  transmit, ms_transmit_color);
852 #else
853  VJOIN2(swp->sw_color, shader_color,
854  reflect, reflect_color,
855  transmit, transmit_color);
856 #endif
857  if (R_DEBUG&RDEBUG_REFRACT) {
858  bu_log("rr_render: lvl=%d end shader=%g reflect=%g, transmit=%g %s\n",
859  ap->a_level,
860  shader_fract, reflect, transmit,
861  pp->pt_regionp->reg_name);
862 #ifdef RT_MULTISPECTRAL
863  {
864  struct bu_vls str = BU_VLS_INIT_ZERO;
865 
866  bu_vls_strcat(&str, "ms_shader_color: ");
867  bn_tabdata_to_tcl(&str, ms_shader_color);
868  bu_vls_strcat(&str, "\nms_reflect_color: ");
869  bn_tabdata_to_tcl(&str, ms_reflect_color);
870  bu_vls_strcat(&str, "\nms_transmit_color: ");
871  bn_tabdata_to_tcl(&str, ms_transmit_color);
872  bu_log("rr_render: %s\n", bu_vls_addr(&str));
873  bu_vls_free(&str);
874  }
875 #else
876  VPRINT("shader ", shader_color);
877  VPRINT("reflect ", reflect_color);
878  VPRINT("transmit", transmit_color);
879 #endif
880  }
881 out:
882  if (R_DEBUG&RDEBUG_REFRACT) {
883 #ifdef RT_MULTISPECTRAL
884  {
885  struct bu_vls str = BU_VLS_INIT_ZERO;
886 
887  bu_vls_strcat(&str, "final swp->msw_color: ");
888  bn_tabdata_to_tcl(&str, swp->msw_color);
889  bu_log("rr_render: %s\n", bu_vls_addr(&str));
890  bu_vls_free(&str);
891  }
892 #else
893  VPRINT("final ", swp->sw_color);
894 #endif
895  }
896 
897  /* Release all the dynamic spectral curves */
898 #ifdef RT_MULTISPECTRAL
899  if (ms_filter_color) bu_free(ms_filter_color, "rr_render: ms_filter_color bn_tabdata*");
900  if (ms_shader_color) bu_free(ms_shader_color, "rr_render: ms_shader_color bn_tabdata*");
901  if (ms_reflect_color) bu_free(ms_reflect_color, "rr_render: ms_reflect_color bn_tabdata*");
902  if (sub_ap.a_spectrum) bu_free(sub_ap.a_spectrum, "rr_render: sub_ap.a_spectrum bn_tabdata*");
903 #endif
904 
905  return 1;
906 }
907 
908 
909 /*
910  * Local Variables:
911  * mode: C
912  * tab-width: 8
913  * indent-tabs-mode: t
914  * c-file-style: "stroustrup"
915  * End:
916  * ex: shiftwidth=4 tabstop=8
917  */
#define BN_TABDATA_NULL
Definition: tabdata.h:107
struct xray a_ray
Actual ray to be shot.
Definition: raytrace.h:1583
fastf_t a_cumlen
cumulative length of ray
Definition: raytrace.h:1625
void bu_log(const char *,...) _BU_ATTR_PRINTF12
Definition: log.c:176
struct region * pt_regionp
ptr to containing region
Definition: raytrace.h:580
struct hit * pt_outhit
OUT hit ptr.
Definition: raytrace.h:579
double rti_radius
radius of model bounding sphere
Definition: raytrace.h:1773
vect_t background
Definition: init.c:46
#define BN_GET_TABDATA(_data, _table)
Definition: tabdata.h:114
void wraypts(vect_t in, vect_t inorm, vect_t out, int id, struct application *app, FILE *fp)
Definition: wray.c:154
#define R_DEBUG
Definition: optical.h:115
struct soltab * seg_stp
pointer back to soltab
Definition: raytrace.h:372
int viewshade(struct application *ap, register const struct partition *pp, register struct shadework *swp)
int a_flag
application-specific flag
Definition: raytrace.h:1626
void bu_vls_strcat(struct bu_vls *vp, const char *s)
Definition: vls.c:368
#define VSET(a, b, c, d)
Definition: color.c:53
#define VSETALL(a, s)
Definition: color.c:54
void bu_semaphore_acquire(unsigned int i)
Definition: semaphore.c:180
Definition: raytrace.h:368
Definition: raytrace.h:248
Header file for the BRL-CAD common definitions.
int max_bounces
Definition: refract.c:45
struct seg * pt_outseg
OUT seg pointer.
Definition: raytrace.h:578
const char * reg_name
Identifying string.
Definition: raytrace.h:539
char pt_outflip
flip outhit->hit_normal
Definition: raytrace.h:582
HIDDEN int rr_miss(struct application *ap)
Definition: refract.c:63
int a_onehit
flag to stop on first hit
Definition: raytrace.h:1586
#define RDEBUG_RAYPLOT
Definition: optical.h:141
int(* a_hit)(struct application *, struct partition *, struct seg *)
called when shot hits model
Definition: raytrace.h:1584
struct hit * pt_inhit
IN hit pointer.
Definition: raytrace.h:577
#define HIDDEN
Definition: common.h:86
#define BN_CK_TABDATA(_p)
Definition: tabdata.h:106
void bn_tabdata_copy(struct bn_tabdata *out, const struct bn_tabdata *in)
Definition: tabdata.c:1025
#define RT_AP_CHECK(_ap)
Definition: raytrace.h:1685
HIDDEN int rr_refract(vect_t v_1, vect_t norml, double ri_1, double ri_2, vect_t v_2)
Definition: refract.c:305
void bn_tabdata_to_tcl(struct bu_vls *vp, const struct bn_tabdata *data)
Definition: tabdata.c:1087
void bu_vls_free(struct bu_vls *vp)
Definition: vls.c:248
fastf_t a_color[3]
application-specific color
Definition: raytrace.h:1620
struct rt_i * a_rt_i
this librt instance
Definition: raytrace.h:1588
void * memset(void *s, int c, size_t n)
void bn_tabdata_scale(struct bn_tabdata *out, const struct bn_tabdata *in1, double scale)
#define RDEBUG_SHOWERR
Definition: optical.h:128
fastf_t a_diverge
slope of beam divergence/mm
Definition: raytrace.h:1600
int a_user
application-specific value
Definition: raytrace.h:1617
#define V3ARGS(a)
Definition: color.c:56
int a_x
Screen X of ray, if applicable.
Definition: raytrace.h:1596
#define NEAR_ZERO(val, epsilon)
Definition: color.c:55
void bn_tabdata_mul(struct bn_tabdata *out, const struct bn_tabdata *in1, const struct bn_tabdata *in2)
Definition: tabdata.c:149
struct seg * pt_inseg
IN seg ptr (gives stp)
Definition: raytrace.h:576
point_t hit_point
DEPRECATED: Intersection point, use VJOIN1 hit_dist.
Definition: raytrace.h:251
struct bn_tabdata * a_spectrum
application-specific bn_tabdata pointer
Definition: raytrace.h:1619
int rr_render(register struct application *ap, const struct partition *pp, struct shadework *swp)
Definition: refract.c:356
int a_level
recursion level (for printing)
Definition: raytrace.h:1595
const char * a_purpose
Debug string: purpose of ray.
Definition: raytrace.h:1598
#define MSG_PROLOGUE
Definition: refract.c:47
#define AIR_GAP_TOL
Definition: refract.c:52
struct bn_tabdata * bn_tabdata_get_constval(double val, const struct bn_table *tabp)
Definition: tabdata.c:1054
#define UNUSED(parameter)
Definition: common.h:239
vect_t a_vvec
application-specific vector
Definition: raytrace.h:1623
#define RDEBUG_RAYWRITE
Definition: optical.h:140
goto out
Definition: nmg_mod.c:3846
char * bu_vls_addr(const struct bu_vls *vp)
Definition: vls.c:111
void shade_inputs(struct application *app, const struct partition *pp, struct shadework *swp, int want)
Definition: shade.c:118
void bn_tabdata_constval(struct bn_tabdata *data, double val)
Definition: tabdata.c:1072
void bu_semaphore_release(unsigned int i)
Definition: semaphore.c:218
void pl_color(register FILE *plotfp, int r, int g, int b)
Definition: plot3.c:325
vect_t r_dir
Direction of ray (UNIT Length)
Definition: raytrace.h:219
void pdv_3line(register FILE *plotfp, const fastf_t *a, const fastf_t *b)
Definition: plot3.c:642
int max_ireflect
Definition: refract.c:44
#define ZERO(val)
Definition: units.c:38
HIDDEN int code(fastf_t x, fastf_t y)
Definition: clip.c:43
void(* a_logoverlap)(struct application *, const struct partition *, const struct bu_ptbl *, const struct partition *)
called to log overlaps
Definition: raytrace.h:1594
point_t r_pt
Point at which ray starts.
Definition: raytrace.h:218
vect_t a_uvec
application-specific vector
Definition: raytrace.h:1622
#define RDEBUG_SHADE
Definition: optical.h:130
#define BU_SEM_SYSCALL
Definition: parallel.h:178
#define RI_AIR
Definition: refract.c:50
int a_y
Screen Y of ray, if applicable.
Definition: raytrace.h:1597
HIDDEN int rr_hit(struct application *ap, struct partition *PartHeadp, struct seg *segp)
Definition: refract.c:87
struct bn_tabdata * bn_tabdata_dup(const struct bn_tabdata *in)
Definition: tabdata.c:1041
#define MSG_INTERVAL
Definition: refract.c:48
int rt_shootray(struct application *ap)
void * a_uptr
application-specific pointer
Definition: raytrace.h:1618
void bu_free(void *ptr, const char *str)
Definition: malloc.c:328
fastf_t a_rbeam
initial beam radius (mm)
Definition: raytrace.h:1599
struct bn_table * spectrum
Definition: init.c:41
int(* a_miss)(struct application *)
called when shot misses
Definition: raytrace.h:1585
#define BU_VLS_INIT_ZERO
Definition: vls.h:84
#define RT_HIT_NORMAL(_normal, _hitp, _stp, _unused, _flipflag)
Definition: raytrace.h:273
struct partition * pt_forw
forwards link
Definition: raytrace.h:574
#define RDEBUG_REFRACT
Definition: optical.h:133
#define RT_APPLICATION_INIT(_p)
Definition: raytrace.h:1676
fastf_t hit_dist
dist from r_pt to hit_point
Definition: raytrace.h:250
fastf_t a_refrac_index
current index of refraction
Definition: raytrace.h:1624
Definition: vls.h:56
void bu_bomb(const char *str) _BU_ATTR_NORETURN
Definition: bomb.c:91
double fastf_t
Definition: defines.h:300
#define VPRINT(a, b)
Definition: raytrace.h:1881
void bn_tabdata_join2(struct bn_tabdata *out, const struct bn_tabdata *in1, double scale2, const struct bn_tabdata *in2, double scale3, const struct bn_tabdata *in3)
Header file for the BRL-CAD Optical Library, LIBOPTICAL.
const struct bn_table * table
Up pointer to definition of X axis.
Definition: tabdata.h:103