50b9569c09f0ebd56427bd506d7989e78147cd37
[blender-staging.git] / source / blender / radiosity / intern / source / radrender.c
1 /* ***************************************
2  *
3  * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
4  *
5  * This program is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU General Public License
7  * as published by the Free Software Foundation; either version 2
8  * of the License, or (at your option) any later version. The Blender
9  * Foundation also sells licenses for use in proprietary software under
10  * the Blender License.  See http://www.blender.org/BL/ for information
11  * about this.
12  *
13  * This program is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with this program; if not, write to the Free Software Foundation,
20  * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
21  *
22  * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
23  * All rights reserved.
24  *
25  * The Original Code is: all of this file.
26  *
27  * Contributor(s): none yet.
28  *
29  * ***** END GPL/BL DUAL LICENSE BLOCK *****
30  */
31  
32 /* radrender.c, aug 2003
33  *
34  * Most of the code here is copied from radiosity code, to optimize for renderfaces.
35  * Shared function calls mostly reside in radfactors.c
36  * No adaptive subdivision takes place
37  *
38  * - do_radio_render();  main call, extern
39  *   - initradfaces(); add radface structs in render faces, init radio globals
40  *   - 
41  *   - initradiosity(); LUTs
42  *   - inithemiwindows();
43  *   - progressiverad(); main itteration loop
44  *     - hemi zbuffers
45  *     - calc rad factors
46  * 
47  *   - closehemiwindows();
48  *   - freeAllRad();
49  *   - make vertex colors
50  *
51  * - during render, materials use totrad as ambient replacement
52  * - free radfaces
53  */
54
55 #include <stdlib.h>
56 #include <string.h>
57 #include <math.h>
58
59 #include "MEM_guardedalloc.h"
60
61 #include "BLI_blenlib.h"
62 #include "BLI_arithb.h"
63 #include "BLI_rand.h"
64
65 #include "BKE_utildefines.h"
66 #include "BKE_global.h"
67 #include "BKE_main.h"
68
69 #include "BIF_screen.h"
70
71 #include "radio.h"
72 #include "render.h" 
73
74 #ifdef HAVE_CONFIG_H
75 #include <config.h>
76 #endif
77
78 /* only needed now for a print, if its useful move to RG */
79 static float maxenergy; 
80
81 /* find the face with maximum energy to become shooter */
82 /* nb: _rr means rad-render version of existing radio call */
83 static VlakRen *findshoot_rr(void)
84 {
85         RadFace *rf;
86         VlakRen *vlr=NULL, *shoot;
87         float energy;
88         int a;
89         
90         shoot= NULL;
91         maxenergy= 0.0;
92         
93         for(a=0; a<R.totvlak; a++) {
94                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
95                 if(vlr->radface) {
96                         rf= vlr->radface;
97                         rf->flag &= ~RAD_SHOOT;
98                         
99                         energy= rf->unshot[0]*rf->area;
100                         energy+= rf->unshot[1]*rf->area;
101                         energy+= rf->unshot[2]*rf->area;
102
103                         if(energy>maxenergy) {
104                                 shoot= vlr;
105                                 maxenergy= energy;
106                         }
107                 }
108         }
109
110         if(shoot) {
111                 maxenergy/= RG.totenergy;
112                 if(maxenergy<RG.convergence) return NULL;
113                 shoot->radface->flag |= RAD_SHOOT;
114         }
115
116         return shoot;
117 }
118
119 static void backface_test_rr(VlakRen *shoot)
120 {
121         VlakRen *vlr=NULL;
122         RadFace *rf;
123         float tvec[3];
124         int a;
125         
126         /* backface testing */
127         for(a=0; a<R.totvlak; a++) {
128                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
129                 if(vlr->radface) {
130                         if(vlr!=shoot) {
131                                 rf= vlr->radface;
132                                 VecSubf(tvec, shoot->radface->cent, rf->cent);
133                                 
134                                 if( tvec[0]*rf->norm[0]+ tvec[1]*rf->norm[1]+ tvec[2]*rf->norm[2] < 0.0) {              
135                                         rf->flag |= RAD_BACKFACE;
136                                 }
137                         }
138                 }
139         }
140 }
141
142 static void clear_backface_test_rr()
143 {
144         VlakRen *vlr=NULL;
145         RadFace *rf;
146         int a;
147         
148         /* backface flag clear */
149         for(a=0; a<R.totvlak; a++) {
150                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
151                 
152                 if(vlr->radface) {
153                         rf= vlr->radface;
154                         rf->flag &= ~RAD_BACKFACE;
155                 }
156         }
157 }
158
159 extern RadView hemitop, hemiside; // radfactors.c
160
161 /* hemi-zbuffering, delivers formfactors array */
162 static void makeformfactors_rr(VlakRen *shoot)
163 {
164         VlakRen *vlr=NULL;
165         RadFace *rf;
166         float len, vec[3], up[3], side[3], tar[5][3], *fp;
167         int a;
168
169         memset(RG.formfactors, 0, sizeof(float)*RG.totelem);
170
171         /* set up hemiview */
172         /* first: upvector for hemitop, we use diagonal hemicubes to prevent aliasing */
173         
174         VecSubf(vec, shoot->v1->co, shoot->radface->cent);
175         Crossf(up, shoot->radface->norm, vec);
176         len= Normalise(up);
177         
178         VECCOPY(hemitop.up, up);
179         VECCOPY(hemiside.up, shoot->radface->norm);
180
181         Crossf(side, shoot->radface->norm, up);
182
183         /* five targets */
184         VecAddf(tar[0], shoot->radface->cent, shoot->radface->norm);
185         VecAddf(tar[1], shoot->radface->cent, up);
186         VecSubf(tar[2], shoot->radface->cent, up);
187         VecAddf(tar[3], shoot->radface->cent, side);
188         VecSubf(tar[4], shoot->radface->cent, side);
189
190         /* camera */
191         VECCOPY(hemiside.cam, shoot->radface->cent);
192         VECCOPY(hemitop.cam, shoot->radface->cent);
193
194         /* do it! */
195         VECCOPY(hemitop.tar, tar[0]);
196         hemizbuf(&hemitop);
197
198         for(a=1; a<5; a++) {
199                 VECCOPY(hemiside.tar, tar[a]);
200                 hemizbuf(&hemiside);
201         }
202
203         /* convert factors to real radiosity */
204         fp= RG.formfactors;
205
206         for(a=0; a<R.totvlak; a++) {
207                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
208                 
209                 if(vlr->radface) {
210                         rf= vlr->radface;
211                         if(*fp!=0.0 && rf->area!=0.0) {
212                                 *fp *= shoot->radface->area/rf->area;
213                                 if(*fp>1.0) *fp= 1.0001;
214                         }
215                         fp++;
216                 }
217         }
218 }
219
220 /* based at RG.formfactors array, distribute shoot energy over other faces */
221 static void applyformfactors_rr(VlakRen *shoot)
222 {
223         VlakRen *vlr=NULL;
224         RadFace *rf;
225         float *fp, *ref, unr, ung, unb, r, g, b;
226         int a;
227
228         unr= shoot->radface->unshot[0];
229         ung= shoot->radface->unshot[1];
230         unb= shoot->radface->unshot[2];
231
232         fp= RG.formfactors;
233         
234         for(a=0; a<R.totvlak; a++) {
235                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
236                 
237                 if(vlr->radface) {
238                         rf= vlr->radface;
239                         if(*fp!= 0.0) {
240                                 
241                                 ref= &(vlr->mat->r);
242                                 
243                                 r= (*fp)*unr*ref[0];
244                                 g= (*fp)*ung*ref[1];
245                                 b= (*fp)*unb*ref[2];
246                                 
247                                 // if(rf->flag & RAD_BACKFACE) {
248                                 
249                                 rf->totrad[0]+= r;
250                                 rf->totrad[1]+= g;
251                                 rf->totrad[2]+= b;
252                                 
253                                 rf->unshot[0]+= r;
254                                 rf->unshot[1]+= g;
255                                 rf->unshot[2]+= b;
256                         }
257                         fp++;
258                 }
259         }
260         /* shoot energy has been shot */
261         shoot->radface->unshot[0]= shoot->radface->unshot[1]= shoot->radface->unshot[2]= 0.0;
262 }
263
264
265 /* main loop for itterations */
266 static void progressiverad_rr()
267 {
268         extern void RE_local_timecursor(int);   // RE_callbacks.c
269         VlakRen *shoot;
270         float unshot[3];
271         int it= 0;
272         
273         shoot= findshoot_rr();
274         while( shoot ) {
275                 
276                 /* backfaces receive no energy, but are zbuffered... */
277                 backface_test_rr(shoot);
278                 
279                 /* ...unless it's two sided */
280                 if(shoot->radface->flag & RAD_TWOSIDED) {
281                         VECCOPY(unshot, shoot->radface->unshot);
282                         VecMulf(shoot->radface->norm, -1.0);
283                         makeformfactors_rr(shoot);
284                         applyformfactors_rr(shoot);
285                         VecMulf(shoot->radface->norm, -1.0);
286                         VECCOPY(shoot->radface->unshot, unshot);
287                 }
288
289                 /* hemi-zbuffers */
290                 makeformfactors_rr(shoot);
291                 /* based at RG.formfactors array, distribute shoot energy over other faces */
292                 applyformfactors_rr(shoot);
293                 
294                 it++;
295                 RE_local_timecursor(it);
296                 
297                 clear_backface_test_rr();
298                 
299                 if(blender_test_break()) break;
300                 if(RG.maxiter && RG.maxiter<=it) break;
301                 
302                 shoot= findshoot_rr();
303         }
304         printf(" Unshot energy:%f\n", 1000.0*maxenergy);
305         
306         RE_local_timecursor((G.scene->r.cfra));
307 }
308
309 static RadFace *radfaces=NULL;
310
311 static void initradfaces(void)  
312 {
313         VlakRen *vlr= NULL;
314         RadFace *rf;
315         int a, b;
316         
317         /* globals */
318         RG.totenergy= 0.0;
319         RG.totpatch= 0;         // we count initial emittors here
320         RG.totelem= 0;          // total # faces are put here (so we can use radfactors.c calls)
321         /* size is needed for hemicube clipping */
322         RG.min[0]= RG.min[1]= RG.min[2]= 1.0e20;
323         RG.max[0]= RG.max[1]= RG.max[2]= -1.0e20;
324         
325         /* count first for fast malloc */
326         for(a=0; a<R.totvlak; a++) {
327                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
328                 
329                 if(vlr->mat->mode & MA_RADIO) {
330                         if(vlr->mat->emit > 0.0) {
331                                 RG.totpatch++;
332                         }
333                         RG.totelem++;
334                 }
335         }
336         
337 printf(" Rad elems: %d emittors %d\n", RG.totelem, RG.totpatch);        
338         if(RG.totelem==0 || RG.totpatch==0) return;
339
340         /* make/init radfaces */
341         rf=radfaces= MEM_callocN(RG.totelem*sizeof(RadFace), "radfaces");
342         for(a=0; a<R.totvlak; a++) {
343                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
344                 
345                 if(vlr->mat->mode & MA_RADIO) {
346                         
347                         /* during render, vlr->n gets flipped/corrected, we cannot have that */
348                         if(vlr->v4) CalcNormFloat4(vlr->v1->co, vlr->v2->co, vlr->v3->co, vlr->v4->co, rf->norm);
349                         else CalcNormFloat(vlr->v1->co, vlr->v2->co, vlr->v3->co, rf->norm);
350                         
351                         rf->totrad[0]= vlr->mat->emit*vlr->mat->r;
352                         rf->totrad[1]= vlr->mat->emit*vlr->mat->g;
353                         rf->totrad[2]= vlr->mat->emit*vlr->mat->b;
354                         VECCOPY(rf->unshot, rf->totrad);
355                         
356                         if(vlr->v4) {
357                                 rf->area= AreaQ3Dfl(vlr->v1->co, vlr->v2->co, vlr->v3->co, vlr->v4->co);
358                                 CalcCent4f(rf->cent, vlr->v1->co, vlr->v2->co, vlr->v3->co, vlr->v4->co);
359                         }
360                         else {
361                                 rf->area= AreaT3Dfl(vlr->v1->co, vlr->v2->co, vlr->v3->co);
362                                 CalcCent3f(rf->cent, vlr->v1->co, vlr->v2->co, vlr->v3->co);
363                         }
364                         
365                         RG.totenergy+= rf->unshot[0]*rf->area;
366                         RG.totenergy+= rf->unshot[1]*rf->area;
367                         RG.totenergy+= rf->unshot[2]*rf->area;
368                         
369                         for(b=0; b<3; b++) {
370                                 RG.min[b]= MIN2(RG.min[b], rf->cent[b]);
371                                 RG.max[b]= MAX2(RG.max[b], rf->cent[b]);
372                         }
373
374 // uncommented; this isnt satisfying, but i leave it in the code for now (ton)                  
375 //                      if(vlr->mat->translucency!=0.0) rf->flag |= RAD_TWOSIDED;
376                         
377                         vlr->radface= rf++;
378                 }
379         }
380         RG.size[0]= (RG.max[0]- RG.min[0]);
381         RG.size[1]= (RG.max[1]- RG.min[1]);
382         RG.size[2]= (RG.max[2]- RG.min[2]);
383         RG.maxsize= MAX3(RG.size[0],RG.size[1],RG.size[2]);
384
385         /* formfactor array */
386         if(RG.formfactors) MEM_freeN(RG.formfactors);
387         if(RG.totelem)
388                 RG.formfactors= MEM_mallocN(sizeof(float)*RG.totelem, "formfactors");
389         else
390                 RG.formfactors= NULL;
391         
392 }
393
394 static void vecaddfac(float *vec, float *v1, float *v2, float fac)
395 {
396         vec[0]= v1[0] + fac*v2[0];
397         vec[1]= v1[1] + fac*v2[1];
398         vec[2]= v1[2] + fac*v2[2];
399
400 }
401
402 #if 0
403 /* unused now, doesnt work... */
404 static void filter_rad_values(void)
405 {
406         VlakRen *vlr=NULL;
407         VertRen *v1=NULL;
408         RadFace *rf;
409         float n1[3], n2[3], n3[4], n4[3], co[4];
410         int a;
411         
412         /* one filter pass */
413         for(a=0; a<R.totvert; a++) {
414                 if((a & 255)==0) v1= R.blove[a>>8]; else v1++;
415                 if(v1->accum>0.0) {
416                         v1->rad[0]= v1->rad[0]/v1->accum;
417                         v1->rad[1]= v1->rad[1]/v1->accum;
418                         v1->rad[2]= v1->rad[2]/v1->accum;
419                         v1->accum= 0.0;
420                 }
421         }
422         /* cosines in verts accumulate in faces */
423         for(a=0; a<R.totvlak; a++) {
424                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
425                 
426                 if(vlr->radface) {
427                         rf= vlr->radface;
428                         
429                         /* calculate cosines of angles, for weighted add (irregular faces) */
430                         VecSubf(n1, vlr->v2->co, vlr->v1->co);
431                         VecSubf(n2, vlr->v3->co, vlr->v2->co);
432                         Normalise(n1);
433                         Normalise(n2);
434         
435                         if(vlr->v4==NULL) {
436                                 VecSubf(n3, vlr->v1->co, vlr->v3->co);
437                                 Normalise(n3);
438                                 
439                                 co[0]= saacos(-n3[0]*n1[0]-n3[1]*n1[1]-n3[2]*n1[2])/M_PI;
440                                 co[1]= saacos(-n1[0]*n2[0]-n1[1]*n2[1]-n1[2]*n2[2])/M_PI;
441                                 co[2]= saacos(-n2[0]*n3[0]-n2[1]*n3[1]-n2[2]*n3[2])/M_PI;
442                                 co[0]= co[1]= co[2]= 1.0/3.0;
443                         }
444                         else {
445                                 VecSubf(n3, vlr->v4->co, vlr->v3->co);
446                                 VecSubf(n4, vlr->v1->co, vlr->v4->co);
447                                 Normalise(n3);
448                                 Normalise(n4);
449                                 
450                                 co[0]= saacos(-n4[0]*n1[0]-n4[1]*n1[1]-n4[2]*n1[2])/M_PI;
451                                 co[1]= saacos(-n1[0]*n2[0]-n1[1]*n2[1]-n1[2]*n2[2])/M_PI;
452                                 co[2]= saacos(-n2[0]*n3[0]-n2[1]*n3[1]-n2[2]*n3[2])/M_PI;
453                                 co[3]= saacos(-n3[0]*n4[0]-n3[1]*n4[1]-n3[2]*n4[2])/M_PI;
454                                 co[0]= co[1]= co[2]= co[3]= 1.0/4.0;
455                         }
456                         
457                         rf->totrad[0]= rf->totrad[1]= rf->totrad[2]= 0.0;
458
459                         vecaddfac(rf->totrad, rf->totrad, vlr->v1->rad, co[0]); 
460                         vecaddfac(rf->totrad, rf->totrad, vlr->v2->rad, co[1]); 
461                         vecaddfac(rf->totrad, rf->totrad, vlr->v3->rad, co[2]); 
462                         if(vlr->v4) {
463                                 vecaddfac(rf->totrad, rf->totrad, vlr->v4->rad, co[3]); 
464                         }
465                 }
466         }
467
468         /* accumulate vertexcolors again */
469         for(a=0; a<R.totvlak; a++) {
470                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
471                 
472                 if(vlr->radface) {
473                         rf= vlr->radface;
474
475                         vecaddfac(vlr->v1->rad, vlr->v1->rad, rf->totrad, rf->area); 
476                         vlr->v1->accum+= rf->area;
477                         vecaddfac(vlr->v2->rad, vlr->v2->rad, rf->totrad, rf->area); 
478                         vlr->v2->accum+= rf->area;
479                         vecaddfac(vlr->v3->rad, vlr->v3->rad, rf->totrad, rf->area); 
480                         vlr->v3->accum+= rf->area;
481                         if(vlr->v4) {
482                                 vecaddfac(vlr->v4->rad, vlr->v4->rad, rf->totrad, rf->area); 
483                                 vlr->v4->accum+= rf->area;
484                         }
485                 }
486         }
487
488 }
489 #endif
490
491 static void make_vertex_rad_values()
492 {
493         VertRen *v1=NULL;
494         VlakRen *vlr=NULL;
495         RadFace *rf;
496         int a;
497
498         RG.igamma= 1.0/RG.gamma;
499         RG.radfactor= RG.radfac*pow(64*64, RG.igamma)/128.0; /* compatible with radio-tool */
500
501         /* accumulate vertexcolors */
502         for(a=0; a<R.totvlak; a++) {
503                 if((a & 255)==0) vlr= R.blovl[a>>8]; else vlr++;
504                 
505                 if(vlr->radface) {
506                         rf= vlr->radface;
507                         
508                         /* apply correction */
509                         rf->totrad[0]= RG.radfactor*pow( rf->totrad[0], RG.igamma);
510                         rf->totrad[1]= RG.radfactor*pow( rf->totrad[1], RG.igamma);
511                         rf->totrad[2]= RG.radfactor*pow( rf->totrad[2], RG.igamma);
512                         
513                         /* correct rf->rad values for color */
514                         if(vlr->mat->r > 0.0) rf->totrad[0]/= vlr->mat->r;
515                         if(vlr->mat->g > 0.0) rf->totrad[1]/= vlr->mat->g;
516                         if(vlr->mat->b > 0.0) rf->totrad[2]/= vlr->mat->b;
517                         
518                         vecaddfac(vlr->v1->rad, vlr->v1->rad, rf->totrad, rf->area); 
519                         vlr->v1->accum+= rf->area;
520                         vecaddfac(vlr->v2->rad, vlr->v2->rad, rf->totrad, rf->area); 
521                         vlr->v2->accum+= rf->area;
522                         vecaddfac(vlr->v3->rad, vlr->v3->rad, rf->totrad, rf->area); 
523                         vlr->v3->accum+= rf->area;
524                         if(vlr->v4) {
525                                 vecaddfac(vlr->v4->rad, vlr->v4->rad, rf->totrad, rf->area); 
526                                 vlr->v4->accum+= rf->area;
527                         }
528                 }
529         }
530         
531         /* make vertex colors */
532         for(a=0; a<R.totvert; a++) {
533                 if((a & 255)==0) v1= R.blove[a>>8]; else v1++;
534                 if(v1->accum>0.0) {
535                         v1->rad[0]/= v1->accum;
536                         v1->rad[1]/= v1->accum;
537                         v1->rad[2]/= v1->accum;                 
538                 }
539         }
540
541 }
542
543 /* main call, extern */
544 void do_radio_render(void)
545 {
546         if(G.scene->radio==NULL) add_radio();
547         freeAllRad();   /* just in case radio-tool is still used */
548         
549         set_radglobal(); /* init the RG struct */
550
551         initradfaces();  /* add radface structs to render faces */
552         if(RG.totenergy>0.0) {
553
554                 initradiosity();        /* LUT's */
555                 inithemiwindows();      /* views, need RG.maxsize for clipping */
556         
557                 progressiverad_rr(); /* main radio loop */
558                 
559                 make_vertex_rad_values(); /* convert face energy to vertex ones */
560
561         }
562         
563         freeAllRad();   /* luts, hemis, sets vars at zero */
564 }
565
566 /* free call, after rendering, extern */
567 void end_radio_render(void)
568 {
569         if(radfaces) MEM_freeN(radfaces);
570         radfaces= NULL;
571 }
572