Quite a large one this time... but now we have:
[blender.git] / source / blender / render / intern / source / rendercore.c
1 /**
2  * $Id$
3  *
4  * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
5  *
6  * This program is free software; you can redistribute it and/or
7  * modify it under the terms of the GNU General Public License
8  * as published by the Free Software Foundation; either version 2
9  * of the License, or (at your option) any later version. The Blender
10  * Foundation also sells licenses for use in proprietary software under
11  * the Blender License.  See http://www.blender.org/BL/ for information
12  * about this.
13  *
14  * This program is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17  * GNU General Public License for more details.
18  *
19  * You should have received a copy of the GNU General Public License
20  * along with this program; if not, write to the Free Software Foundation,
21  * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
22  *
23  * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
24  * All rights reserved.
25  *
26  * The Original Code is: all of this file.
27  *
28  * Contributor(s): Hos, Robert Wenzlaff.
29  *
30  * ***** END GPL/BL DUAL LICENSE BLOCK *****
31  */
32
33
34 /* External modules: */
35 #include "MEM_guardedalloc.h"
36 #include "BLI_arithb.h"
37 #include "IMB_imbuf_types.h"
38 #include "IMB_imbuf.h"
39 #include "MTC_matrixops.h"
40 #include "BKE_utildefines.h"
41
42 #include "DNA_mesh_types.h"
43 #include "DNA_meshdata_types.h"
44 #include "DNA_image_types.h"
45 #include "DNA_object_types.h"
46 #include "DNA_camera_types.h"
47 #include "DNA_lamp_types.h"
48 #include "DNA_texture_types.h"
49
50 #include "BKE_global.h"
51 #include "BKE_texture.h"
52
53 /* local include */
54 #include "RE_callbacks.h"
55 #include "old_zbuffer_types.h"
56 #include "render.h"
57 #include "render_intern.h"
58 #include "zbuf.h"               /* stuff like bgnaccumbuf, fillrect, ...*/
59 #include "pixelblending.h"
60 #include "shadbuf.h"
61 #include "renderHelp.h"
62
63 #include "jitter.h"
64 #include "texture.h"
65
66 /* system includes */
67 #include <math.h>
68 #include <string.h>
69 #include <stdlib.h>
70
71
72 /* own include */
73 #include "rendercore.h"
74 #include "rendercore_int.h"
75
76 #ifdef HAVE_CONFIG_H
77 #include <config.h>
78 #endif
79
80 /* globals for this file */
81 /* moved to renderData.c? Not yet... */
82 RE_Render R;
83 Osa O;
84
85 PixStrMain psmfirst;
86 int psmteller;
87
88 float fmask[256], centLut[16];
89 unsigned short usegamtab=0, *mask1[9], *mask2[9], *igamtab1, *igamtab2, *gamtab;
90 char cmask[256], *centmask;
91
92 /* functions */
93 /* comes from texture.c (only used here !) */
94 /*  extern void do_halo_tex(HaloRen *har, float xn, float yn, float *colf); */
95
96 void gamtabdit(unsigned short *in, char *out);
97 /*  int count_mask(unsigned short ); */
98 void scanlinehalo(unsigned int *rectz, unsigned int *rectt, short ys);
99 /*  void add_halo_flare(void); */
100 void edge_enhance(void);
101
102 /* Dither with gamma table? */
103 void gamtabdit(unsigned short *in, char *out)
104 /*  unsigned short *in; */
105 /*  char *out; */
106 {
107         static short rerr=0, gerr=0, berr=0;
108         unsigned int col;
109         char *cp;
110
111         cp= (char *)&col;
112         out[0]= in[0]>>8;
113
114         col= gamtab[in[2]]+berr;
115         if(col>65535) col= 65535;
116         out[1]= cp[2];
117         berr= cp[3];
118
119         col= gamtab[in[4]]+gerr;
120         if(col>65535) col= 65535;
121         out[2]= cp[2];
122         gerr= cp[3];
123
124         col= gamtab[in[6]]+rerr;
125         if(col>65535) col= 65535;
126         out[3]= cp[2];
127         rerr= cp[3];
128
129 }
130
131 float mistfactor(float *co)     /* dist en height, return alpha */
132 {
133         float fac, hi;
134         
135         fac= R.zcor - R.wrld.miststa;   /* R.zcor is calculated per pixel */
136
137         /* fac= -co[2]-R.wrld.miststa; */
138
139         if(fac>0.0) {
140                 if(fac< R.wrld.mistdist) {
141                         
142                         fac= (fac/(R.wrld.mistdist));
143                         
144                         if(R.wrld.mistype==0) fac*= fac;
145                         else if(R.wrld.mistype==1);
146                         else fac= sqrt(fac);
147                 }
148                 else fac= 1.0;
149         }
150         else fac= 0.0;
151         
152         /* height switched off mist */
153         if(R.wrld.misthi!=0.0 && fac!=0.0) {
154                 /* at height misthi the mist is completely gone */
155
156                 hi= R.viewinv[0][2]*co[0]+R.viewinv[1][2]*co[1]+R.viewinv[2][2]*co[2]+R.viewinv[3][2];
157                 
158                 if(hi>R.wrld.misthi) fac= 0.0;
159                 else if(hi>0.0) {
160                         hi= (R.wrld.misthi-hi)/R.wrld.misthi;
161                         fac*= hi*hi;
162                 }
163         }
164
165         return (1.0-fac)* (1-R.wrld.misi);      
166 }
167
168 void RE_sky(float *view, float *col)
169 {
170         float lo[3];
171
172         R.wrld.skytype |= WO_ZENUP;
173         
174         if(R.wrld.skytype & WO_SKYREAL) {
175         
176                 R.inprz= view[0]*R.grvec[0]+ view[1]*R.grvec[1]+ view[2]*R.grvec[2];
177
178                 if(R.inprz<0.0) R.wrld.skytype-= WO_ZENUP;
179                 R.inprz= fabs(R.inprz);
180         }
181         else if(R.wrld.skytype & WO_SKYPAPER) {
182                 R.inprz= 0.5+ 0.5*view[1];
183         }
184         else {
185                 R.inprz= fabs(0.5+ view[1]);
186         }
187
188         if(R.wrld.skytype & WO_SKYTEX) {
189                 VECCOPY(lo, view);
190                 if(R.wrld.skytype & WO_SKYREAL) {
191                         
192                         MTC_Mat3MulVecfl(R.imat, lo);
193         
194                         SWAP(float, lo[1],  lo[2]);
195                         
196                 }
197
198                 do_sky_tex(lo);
199                 
200         }
201
202         if(R.inprz>1.0) R.inprz= 1.0;
203         R.inprh= 1.0-R.inprz;
204
205         if(R.wrld.skytype & WO_SKYBLEND) {
206                 col[0]= (R.inprh*R.wrld.horr + R.inprz*R.wrld.zenr);
207                 col[1]= (R.inprh*R.wrld.horg + R.inprz*R.wrld.zeng);
208                 col[2]= (R.inprh*R.wrld.horb + R.inprz*R.wrld.zenb);
209         }
210         else {
211                 col[0]= R.wrld.horr;
212                 col[1]= R.wrld.horg;
213                 col[2]= R.wrld.horb;
214         }
215 }
216
217 void RE_sky_char(float *view, char *col)
218 {
219         float f, colf[3];
220         
221         RE_sky(view, colf);
222         f= 255.0*colf[0];
223         if(f<=0.0) col[0]= 0; else if(f>255.0) col[0]= 255;
224         else col[0]= (char)f;
225         f= 255.0*colf[1];
226         if(f<=0.0) col[1]= 0; else if(f>255.0) col[1]= 255;
227         else col[1]= (char)f;
228         f= 255.0*colf[2];
229         if(f<=0.0) col[2]= 0; else if(f>255.0) col[2]= 255;
230         else col[2]= (char)f;
231         col[3]= 1;      /* to prevent wrong optimalisation alphaover of flares */
232 }
233
234 /* ------------------------------------------------------------------------- */
235
236 void scanlinesky(char *rect, int y)
237 {
238         /* have to type this! set to :  addalphaUnder: char*, char*
239          * addalphaUnderGamma: ditto called with char *, uint* !!!
240          * unmangle this shit... */
241         void (*alphafunc)();
242         float fac, u, v, view[3];
243         int dx, x, ofs;
244         unsigned int col=0, *rt;
245         char *cp, *cp1;
246         
247         if(R.r.alphamode & R_ALPHAPREMUL) return;
248         
249         if(R.r.alphamode & R_ALPHAKEY) {
250                 
251                 cp= (char *)&col;
252                 cp[3]= 0;
253                 cp[0]= 255.0*R.wrld.horr;
254                 cp[1]= 255.0*R.wrld.horg;
255                 cp[2]= 255.0*R.wrld.horb;
256                 
257                 for(x=0; x<R.rectx; x++, rect+= 4) {
258                         if(rect[3]==0) {
259                                 *( ( unsigned int *)rect)= col;
260                         }
261                         else {
262                                 /* prevent  'col' to be in the image */
263                                 cp1= (char *)rect;
264                                 if( cp[0]==cp1[0] && cp[1]==cp1[1] && cp[2]==cp1[2] ) {
265
266                                         if(cp1[3]==255) cp1[3]= 254; 
267                                         else cp1[3]++; 
268                                 }
269
270                                 if(rect[3]!=255) { 
271                                         keyalpha(rect); 
272                                 } 
273                         }
274                 }
275                 return;
276         }
277
278         if(R.wrld.mode & WO_MIST) alphafunc= addalphaUnder;
279         else alphafunc= addalphaUnderGamma;
280
281
282         if(R.r.bufflag & 1) {
283                 if(R.backbuf->ok) {
284                         if(R.backbuf->ibuf==0) {
285                                 R.backbuf->ibuf= IMB_loadiffname(R.backbuf->name, IB_rect);
286                                 if(R.backbuf->ibuf==0) {
287                                         R.backbuf->ok= 0;
288                                         return;
289                                 }
290                         }
291                         /* which scanline/ */
292                         y= ((y+R.afmy+R.ystart)*R.backbuf->ibuf->y)/(R.recty);
293                         
294                         if(R.flag & R_SEC_FIELD) {
295                                 if((R.r.mode & R_ODDFIELD)==0) {
296                                         if( y<R.backbuf->ibuf->y) y++;
297                                 }
298                                 else {
299                                         if( y>0) y--;
300                                 }
301                         }
302                                 
303                         rt= (R.backbuf->ibuf->rect + y*R.backbuf->ibuf->x);
304                         
305                         /* at which location */
306                         fac= ((float)R.backbuf->ibuf->x)/(float)(2*R.afmx);
307                         ofs= (R.afmx+R.xstart)*fac;
308                         rt+= ofs;
309
310                         dx= (int) (65536.0*fac);
311                         
312                         ofs= 0;
313                         x= R.rectx;
314                         while( x-- ) {
315                                 if( rect[3] != 255) { 
316                                         if(rect[3]==0) *((unsigned int *)rect)= *rt;
317                                         else {
318                                                 alphafunc(rect, rt);
319                                         }
320                                 }
321                                 rect+= 4;
322                                 
323                                 ofs+= dx;
324                                 while( ofs > 65535 ) {
325                                         ofs-= 65536;
326                                         rt++;
327                                 }
328                         }
329
330                 }
331                 return;
332         }
333
334         if((R.wrld.skytype & (WO_SKYBLEND+WO_SKYTEX))==0) {
335                 for(x=0; x<R.rectx; x++, rect+= 4) {
336                         if(rect[3] != 255) { 
337                                 if(rect[3]==0) *((unsigned int *)rect)= R.wrld.fastcol;
338                                 else {
339                                         alphafunc(rect, &R.wrld.fastcol);
340                                 }
341                         }
342                 }
343         }
344         else {
345                 
346                 for(x=0; x<R.rectx; x++, rect+= 4) {
347                         if(rect[3] < 254) { 
348                                 if(R.wrld.skytype & WO_SKYPAPER) {
349                                         view[0]= (x+(R.xstart))/(float)R.afmx;
350                                         view[1]= (y+(R.ystart))/(float)R.afmy;
351                                         view[2]= 0.0;
352                                 }
353                                 else {
354                                         view[0]= (x+(R.xstart)+0.5);
355                 
356                                         if(R.flag & R_SEC_FIELD) {
357                                                 if(R.r.mode & R_ODDFIELD) view[1]= (y+R.ystart)*R.ycor;
358                                                 else view[1]= (y+R.ystart+1.0)*R.ycor;
359                                         }
360                                         else view[1]= (y+R.ystart+0.5)*R.ycor;
361                                         
362                                         view[2]= -R.viewfac;
363         
364                                         fac= Normalise(view);
365                                         if(R.wrld.skytype & WO_SKYTEX) {
366                                                 O.dxview= 1.0/fac;
367                                                 O.dyview= R.ycor/fac;
368                                         }
369                                 }
370                                 
371                                 if(R.r.mode & R_PANORAMA) {
372                                         float panoco, panosi;
373                                         panoco = getPanovCo();
374                                         panosi = getPanovSi();
375                                         u= view[0]; v= view[2];
376                                         
377                                         view[0]= panoco*u + panosi*v;
378                                         view[2]= -panosi*u + panoco*v;
379                                 }
380
381                                 RE_sky_char(view, (char *)&col);
382         
383                                 if(rect[3]==0) *((unsigned int *)rect)= col;
384                                 else alphafunc(rect, &col);
385                         }
386                 }
387         }       
388 }
389
390 /* ************************************** */
391
392
393 void spothalo(struct LampRen *lar, ShadeInput *shi, float *intens)
394 {
395         double a, b, c, disc, nray[3], npos[3];
396         float t0, t1 = 0.0, t2= 0.0, t3, haint;
397         float p1[3], p2[3], ladist, maxz = 0.0, maxy = 0.0;
398         int snijp, doclip=1, use_yco=0;
399         int ok1=0, ok2=0;
400         
401         *intens= 0.0;
402         haint= lar->haint;
403         
404         VECCOPY(npos, lar->sh_invcampos);       /* in initlamp calculated */
405         
406         /* rotate view */
407         VECCOPY(nray, shi->view);
408         MTC_Mat3MulVecd(lar->imat, nray);
409         
410         if(R.wrld.mode & WO_MIST) {
411                 /* patchy... */
412                 R.zcor= -lar->co[2];
413                 haint *= mistfactor(lar->co);
414                 if(haint==0.0) {
415                         return;
416                 }
417         }
418
419
420         /* rotate maxz */
421         if(shi->co[2]==0.0) doclip= 0;  /* for when halo at sky */
422         else {
423                 p1[0]= shi->co[0]-lar->co[0];
424                 p1[1]= shi->co[1]-lar->co[1];
425                 p1[2]= shi->co[2]-lar->co[2];
426         
427                 maxz= lar->imat[0][2]*p1[0]+lar->imat[1][2]*p1[1]+lar->imat[2][2]*p1[2];
428                 maxz*= lar->sh_zfac;
429                 maxy= lar->imat[0][1]*p1[0]+lar->imat[1][1]*p1[1]+lar->imat[2][1]*p1[2];
430
431                 if( fabs(nray[2]) <0.000001 ) use_yco= 1;
432         }
433         
434         /* scale z to make sure volume is normalized */ 
435         nray[2]*= lar->sh_zfac;
436         /* nray does not need normalization */
437         
438         ladist= lar->sh_zfac*lar->dist;
439         
440         /* solve */
441         a = nray[0] * nray[0] + nray[1] * nray[1] - nray[2]*nray[2];
442         b = nray[0] * npos[0] + nray[1] * npos[1] - nray[2]*npos[2];
443         c = npos[0] * npos[0] + npos[1] * npos[1] - npos[2]*npos[2];
444
445         snijp= 0;
446         if (fabs(a) < 0.00000001) {
447                 /*
448                  * Only one intersection point...
449                  */
450                 return;
451         }
452         else {
453                 disc = b*b - a*c;
454                 
455                 if(disc==0.0) {
456                         t1=t2= (-b)/ a;
457                         snijp= 2;
458                 }
459                 else if (disc > 0.0) {
460                         disc = sqrt(disc);
461                         t1 = (-b + disc) / a;
462                         t2 = (-b - disc) / a;
463                         snijp= 2;
464                 }
465         }
466         if(snijp==2) {
467                 /* sort */
468                 if(t1>t2) {
469                         a= t1; t1= t2; t2= a;
470                 }
471
472                 /* z of intersection points with diabolo */
473                 p1[2]= npos[2] + t1*nray[2];
474                 p2[2]= npos[2] + t2*nray[2];
475
476                 /* evaluate both points */
477                 if(p1[2]<=0.0) ok1= 1;
478                 if(p2[2]<=0.0 && t1!=t2) ok2= 1;
479                 
480                 /* at least 1 point with negative z */
481                 if(ok1==0 && ok2==0) return;
482                 
483                 /* intersction point with -ladist, the bottom of the cone */
484                 if(use_yco==0) {
485                         t3= (-ladist-npos[2])/nray[2];
486                                 
487                         /* de we have to replace one of the intersection points? */
488                         if(ok1) {
489                                 if(p1[2]<-ladist) t1= t3;
490                         }
491                         else {
492                                 ok1= 1;
493                                 t1= t3;
494                         }
495                         if(ok2) {
496                                 if(p2[2]<-ladist) t2= t3;
497                         }
498                         else {
499                                 ok2= 1;
500                                 t2= t3;
501                         }
502                 }
503                 else if(ok1==0 || ok2==0) return;
504                 
505                 /* at least 1 visible interesction point */
506                 if(t1<0.0 && t2<0.0) return;
507                 
508                 if(t1<0.0) t1= 0.0;
509                 if(t2<0.0) t2= 0.0;
510                 
511                 if(t1==t2) return;
512                 
513                 /* sort again to be sure */
514                 if(t1>t2) {
515                         a= t1; t1= t2; t2= a;
516                 }
517                 
518                 /* calculate t0: is the maximum visible z (when halo is intersected by face) */ 
519                 if(doclip) {
520                         if(use_yco==0) t0= (maxz-npos[2])/nray[2];
521                         else t0= (maxy-npos[1])/nray[1];
522
523                         if(t0<t1) return;
524                         if(t0<t2) t2= t0;
525                 }
526
527                 /* calc points */
528                 p1[0]= npos[0] + t1*nray[0];
529                 p1[1]= npos[1] + t1*nray[1];
530                 p1[2]= npos[2] + t1*nray[2];
531                 p2[0]= npos[0] + t2*nray[0];
532                 p2[1]= npos[1] + t2*nray[1];
533                 p2[2]= npos[2] + t2*nray[2];
534                 
535                         
536                 /* now we have 2 points, make three lengths with it */
537                 
538                 a= sqrt(p1[0]*p1[0]+p1[1]*p1[1]+p1[2]*p1[2]);
539                 b= sqrt(p2[0]*p2[0]+p2[1]*p2[1]+p2[2]*p2[2]);
540                 c= VecLenf(p1, p2);
541                 
542                 a/= ladist;
543                 a= sqrt(a);
544                 b/= ladist; 
545                 b= sqrt(b);
546                 c/= ladist;
547                 
548                 *intens= c*( (1.0-a)+(1.0-b) );
549
550                 /* WATCH IT: do not clip a,b en c at 1.0, this gives nasty little overflows
551                         at the edges (especially with narrow halos) */
552                 if(*intens<=0.0) return;
553
554                 /* soft area */
555                 /* not needed because t0 has been used for p1/p2 as well */
556                 /* if(doclip && t0<t2) { */
557                 /*      *intens *= (t0-t1)/(t2-t1); */
558                 /* } */
559                 
560                 *intens *= haint;
561                 
562                 if(lar->shb && lar->shb->shadhalostep) {
563                         *intens *= shadow_halo(lar, p1, p2);
564                 }
565                 
566         }
567 }
568
569 static void renderspothalo(ShadeInput *shi, float *col)
570 {
571         LampRen *lar;
572         float i;
573         int a;
574         
575         for(a=0; a<R.totlamp; a++) {
576                 lar= R.la[a];
577                 if(lar->type==LA_SPOT && (lar->mode & LA_HALO) && lar->haint>0) {
578         
579                         if(lar->org) {
580                                 lar->r= lar->org->r;
581                                 lar->g= lar->org->g;
582                                 lar->b= lar->org->b;
583                         }
584
585                         spothalo(lar, shi, &i);
586                         if(i>0.0) {
587                                 col[3]+= i;
588                                 col[0]+= i*lar->r;
589                                 col[1]+= i*lar->g;
590                                 col[2]+= i*lar->b;                      
591                         }
592                 }
593         }
594         /* clip alpha, is needed for unified 'alpha threshold' (vanillaRenderPipe.c) */
595         if(col[3]>1.0) col[3]= 1.0;
596 }
597
598 void render_lighting_halo(HaloRen *har, float *colf)
599 {
600         LampRen *lar;
601         float i, inp, inpr, rco[3], lv[3], lampdist, ld, t, *vn;
602         float ir, ig, ib, shadfac, soft;
603         int a;
604         
605         ir= ig= ib= 0.0;
606         VECCOPY(rco, har->co);
607         vn= har->no;
608         
609         O.dxco[0]= har->hasize;
610         O.dxco[1]= 0.0;
611         O.dxco[2]= 0.0;
612
613         O.dyco[0]= 0.0;
614         O.dyco[1]= har->hasize;
615         O.dyco[2]= 0.0;
616
617         for(a=0; a<R.totlamp; a++) {
618                 lar= R.la[a];
619
620                 /* test for lamplayer */
621                 if(lar->mode & LA_LAYER) if((lar->lay & har->lay)==0) continue;
622
623                 /* lampdist cacluation */
624                 if(lar->type==LA_SUN || lar->type==LA_HEMI) {
625                         VECCOPY(lv, lar->vec);
626                         lampdist= 1.0;
627                 }
628                 else {
629                         lv[0]= rco[0]-lar->co[0];
630                         lv[1]= rco[1]-lar->co[1];
631                         lv[2]= rco[2]-lar->co[2];
632                         ld= sqrt(lv[0]*lv[0]+lv[1]*lv[1]+lv[2]*lv[2]);
633                         lv[0]/= ld;
634                         lv[1]/= ld;
635                         lv[2]/= ld;
636                         
637                         /* ld is re-used further on (texco's) */
638                         
639                         if(lar->mode & LA_QUAD) {
640                                 t= 1.0;
641                                 if(lar->ld1>0.0)
642                                         t= lar->dist/(lar->dist+lar->ld1*ld);
643                                 if(lar->ld2>0.0)
644                                         t*= lar->distkw/(lar->distkw+lar->ld2*ld*ld);
645
646                                 lampdist= t;
647                         }
648                         else {
649                                 lampdist= (lar->dist/(lar->dist+ld));
650                         }
651
652                         if(lar->mode & LA_SPHERE) {
653                                 t= lar->dist - ld;
654                                 if(t<0.0) continue;
655                                 
656                                 t/= lar->dist;
657                                 lampdist*= (t);
658                         }
659                         
660                 }
661                 
662                 if(lar->mode & LA_TEXTURE) {
663                         ShadeInput shi;
664                         VECCOPY(shi.co, rco);
665                         shi.osatex= 0;
666                         do_lamp_tex(lar, lv, &shi);
667                 }
668                 
669                 if(lar->type==LA_SPOT) {
670
671                         if(lar->mode & LA_SQUARE) {
672                                 if(lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2]>0.0) {
673                                         float x, lvrot[3];
674                                         
675                                         /* rotate view to lampspace */
676                                         VECCOPY(lvrot, lv);
677                                         MTC_Mat3MulVecfl(lar->imat, lvrot);
678                                         
679                                         x= MAX2(fabs(lvrot[0]/lvrot[2]) , fabs(lvrot[1]/lvrot[2]));
680                                         /* 1.0/(sqrt(1+x*x)) is equivalent to cos(atan(x)) */
681
682                                         inpr= 1.0/(sqrt(1.0+x*x));
683                                 }
684                                 else inpr= 0.0;
685                         }
686                         else {
687                                 inpr= lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2];
688                         }
689
690                         t= lar->spotsi;
691                         if(inpr<t) continue;
692                         else {
693                                 t= inpr-t;
694                                 i= 1.0;
695                                 soft= 1.0;
696                                 if(t<lar->spotbl && lar->spotbl!=0.0) {
697                                         /* soft area */
698                                         i= t/lar->spotbl;
699                                         t= i*i;
700                                         soft= (3.0*t-2.0*t*i);
701                                         inpr*= soft;
702                                 }
703                                 if(lar->mode & LA_ONLYSHADOW) {
704                                         /* if(ma->mode & MA_SHADOW) { */
705                                                 /* dot product positive: front side face! */
706                                                 inp= vn[0]*lv[0] + vn[1]*lv[1] + vn[2]*lv[2];
707                                                 if(inp>0.0) {
708                                                         /* testshadowbuf==0.0 : 100% shadow */
709                                                         shadfac = testshadowbuf(lar->shb, rco, inp);
710                                                         if( shadfac>0.0 ) {
711                                                                 shadfac*= inp*soft*lar->energy;
712                                                                 ir -= shadfac;
713                                                                 ig -= shadfac;
714                                                                 ib -= shadfac;
715                                                                 
716                                                                 continue;
717                                                         }
718                                                 }
719                                         /* } */
720                                 }
721                                 lampdist*=inpr;
722                         }
723                         if(lar->mode & LA_ONLYSHADOW) continue;
724
725                 }
726
727                 /* dot product and  reflectivity*/
728                 
729                 inp= 1.0-fabs(vn[0]*lv[0] + vn[1]*lv[1] + vn[2]*lv[2]);
730                 
731                 /* inp= cos(0.5*M_PI-acos(inp)); */
732                 
733                 i= inp;
734                 
735                 if(lar->type==LA_HEMI) {
736                         i= 0.5*i+0.5;
737                 }
738                 if(i>0.0) {
739                         i*= lampdist;
740                         /* i*= lampdist*ma->ref; */
741                 }
742
743                 /* shadow  */
744                 if(i> -0.41) {                  /* heuristic valua! */
745                         shadfac= 1.0;
746                         if(lar->shb) {
747                                 /* if(ma->mode & MA_SHADOW) { */
748                                 shadfac = testshadowbuf(lar->shb, rco, inp);
749                                 if(shadfac==0.0) continue;
750                                 i*= shadfac;
751                                 /* } */
752                         }
753                 }
754                 
755                 if(i>0.0) {
756                         ir+= i*lar->r;
757                         ig+= i*lar->g;
758                         ib+= i*lar->b;
759                 }
760         }
761         
762         if(ir<0.0) ir= 0.0;
763         if(ig<0.0) ig= 0.0;
764         if(ib<0.0) ib= 0.0;
765
766         colf[0]*= ir;
767         colf[1]*= ig;
768         colf[2]*= ib;
769         
770 }
771
772
773 void RE_shadehalo(HaloRen *har, char *col, unsigned int zz, float dist, float xn, float yn, short flarec)
774 {
775         /* fill in in col */
776         extern float hashvectf[];
777         float t, zn, radist, ringf=0.0, linef=0.0, alpha, si, co, colf[4];
778         int colt, a;
779    
780         if(R.wrld.mode & WO_MIST) {
781        if(har->type & HA_ONLYSKY) {
782            /* stars have no mist */
783            alpha= har->alfa;
784        }
785        else {
786            /* patchy... */
787            R.zcor= -har->co[2];
788            alpha= mistfactor(har->co)*har->alfa;
789        }
790         }
791         else alpha= har->alfa;
792         
793         if(alpha==0.0) {
794                 col[0] = 0;
795                 col[1] = 0;
796                 col[2] = 0;
797                 col[3] = 0;
798
799 /*              *( (int *)col )=0; */
800
801                 return;
802         }
803
804         radist= sqrt(dist);
805
806         /* watch it: abused value: flarec was set to zero in pixstruct */
807         if(flarec) har->pixels+= (int)(har->rad-radist);
808
809         if(har->ringc) {
810                 float *rc, fac;
811                 int ofs;
812                 
813                 /* per ring an antialised circle */
814                 ofs= har->seed;
815                 
816                 for(a= har->ringc; a>0; a--, ofs+=2) {
817                         
818                         rc= hashvectf + (ofs % 768);
819                         
820                         fac= fabs( rc[1]*(har->rad*fabs(rc[0]) - radist) );
821                         
822                         if(fac< 1.0) {
823                                 ringf+= (1.0-fac);
824                         }
825                 }
826         }
827
828         if(har->type & HA_VECT) {
829                 dist= fabs( har->cos*(yn) - har->sin*(xn) )/har->rad;
830                 if(dist>1.0) dist= 1.0;
831                 if(har->tex) {
832                         zn= har->sin*xn - har->cos*yn;
833                         yn= har->cos*xn + har->sin*yn;
834                         xn= zn;
835                 }
836         }
837         else dist= dist/har->radsq;
838
839         if(har->type & HA_FLARECIRC) {
840                 
841                 dist= 0.5+fabs(dist-0.5);
842                 
843         }
844
845         if(har->hard>=30) {
846                 dist= sqrt(dist);
847                 if(har->hard>=40) {
848                         dist= sin(dist*M_PI_2);
849                         if(har->hard>=50) {
850                                 dist= sqrt(dist);
851                         }
852                 }
853         }
854         else if(har->hard<20) dist*=dist;
855
856         dist=(1.0-dist);
857         
858         if(har->linec) {
859                 float *rc, fac;
860                 int ofs;
861                 
862                 /* per starpoint an antialaised line */
863                 ofs= har->seed;
864                 
865                 for(a= har->linec; a>0; a--, ofs+=3) {
866                         
867                         rc= hashvectf + (ofs % 768);
868                         
869                         fac= fabs( (xn)*rc[0]+(yn)*rc[1]);
870                         
871                         if(fac< 1.0 ) {
872                                 linef+= (1.0-fac);
873                         }
874                 }
875                 
876                 linef*= dist;
877                 
878         }
879
880         if(har->starpoints) {
881                 float ster, hoek;
882                 /* rotation */
883                 hoek= atan2(yn, xn);
884                 hoek*= (1.0+0.25*har->starpoints);
885                 
886                 co= cos(hoek);
887                 si= sin(hoek);
888                 
889                 hoek= (co*xn+si*yn)*(co*yn-si*xn);
890                 
891                 ster= fabs(hoek);
892                 if(ster>1.0) {
893                         ster= (har->rad)/(ster);
894                         
895                         if(ster<1.0) dist*= sqrt(ster);
896                 }
897         }
898         
899         /* halo intersected? */
900         if(har->zs> zz-har->zd) {
901                 t= ((float)(zz-har->zs))/(float)har->zd;
902                 alpha*= sqrt(sqrt(t));
903         }
904
905         dist*= alpha;
906         ringf*= dist;
907         linef*= alpha;
908         
909         if(dist<0.003) {
910                 *( (int *)col )=0;
911                 return;
912         }
913
914         /* texture? */
915         if(har->tex) {
916                 colf[3]= dist;
917                 do_halo_tex(har, xn, yn, colf);
918                 
919                 /* dist== colf[3]; */
920                 
921                 colf[0]*= colf[3];
922                 colf[1]*= colf[3];
923                 colf[2]*= colf[3];
924                 
925         }
926         else {
927                 colf[0]= dist*har->r;
928                 colf[1]= dist*har->g;
929                 colf[2]= dist*har->b;
930                 
931                 if(har->type & HA_XALPHA) colf[3]= dist*dist;
932                 else colf[3]= dist;
933         }
934
935         if(har->mat && har->mat->mode & MA_HALO_SHADE) {
936                 /* we test for lights because of preview... */
937                 if(R.totlamp) render_lighting_halo(har, colf);
938         }
939
940         if(linef!=0.0) {
941                 Material *ma= har->mat;
942                 
943                 colf[0]+= 255.0*linef*ma->specr;
944                 colf[1]+= 255.0*linef*ma->specg;
945                 colf[2]+= 255.0*linef*ma->specb;
946                 
947                 if(har->type & HA_XALPHA) colf[3]+= linef*linef;
948                 else colf[3]+= linef;
949         }
950         if(ringf!=0.0) {
951                 Material *ma= har->mat;
952                 
953                 colf[0]+= 255.0*ringf*ma->mirr;
954                 colf[1]+= 255.0*ringf*ma->mirg;
955                 colf[2]+= 255.0*ringf*ma->mirb;
956                 
957                 if(har->type & HA_XALPHA) colf[3]+= ringf*ringf;
958                 else colf[3]+= ringf;
959         }
960
961         colt= 255.0*colf[3];
962         if(colt>254) col[3]= 255; else col[3]= colt;
963
964         colt= colf[2];
965         if(colt>254) col[2]= 255; else col[2]= colt;
966
967         colt= colf[1];
968         if(colt>254) col[1]= 255; else col[1]= colt;
969
970         colt= colf[0];
971         if(colt>254) col[0]= 255; else col[0]= colt;
972
973 }
974
975
976 unsigned int calchalo_z(HaloRen *har, unsigned int zz)
977 {
978
979         if(har->type & HA_ONLYSKY) {
980                 if(zz!=0x7FFFFFFF) zz= 0;
981         }
982         else {
983                 zz= (zz>>8);
984                 if(zz<0x800000) zz= (zz+0x7FFFFF);
985                 else zz= (zz-0x800000);
986         }
987         return zz;
988 }
989
990 void scanlinehaloPS(unsigned int *rectz, long *rectdelta, unsigned int *rectt, short ys)
991 {
992         HaloRen *har = NULL;
993         PixStr *ps;
994         float dist,xsq,ysq,xn,yn;
995         unsigned int a, *rz, *rt, zz;
996         long *rd;
997         int accol[4];
998         short minx,maxx,x,aantal, aantalm, flarec;
999         char col[4];
1000
1001         for(a=0;a<R.tothalo;a++) {
1002                 if((a & 255)==0) har= R.bloha[a>>8];
1003                 else har++;
1004
1005                 if( RE_local_test_break() ) break;  
1006
1007                 if(ys>har->maxy);
1008                 else if(ys<har->miny);
1009                 else {
1010                         minx= floor(har->xs-har->rad);
1011                         maxx= ceil(har->xs+har->rad);
1012                         
1013                         if(maxx<0);
1014                         else if(R.rectx<minx);
1015                         else {
1016                                 if(minx<0) minx= 0;
1017                                 if(maxx>=R.rectx) maxx= R.rectx-1;
1018
1019                                 rt= rectt+minx;
1020                                 rd= rectdelta+minx;
1021                                 rz= rectz+minx;
1022
1023                                 yn= (ys-har->ys)*R.ycor;
1024                                 ysq= yn*yn;
1025                                 for(x=minx; x<=maxx; x++) {
1026                                         
1027                                         flarec= har->flarec;    /* har->pixels is inly allowd to count once */
1028
1029                                         if( IS_A_POINTER_CODE(*rd)) {
1030                                                 xn= x-har->xs;
1031                                                 xsq= xn*xn;
1032                                                 dist= xsq+ysq;
1033                                                 if(dist<har->radsq) {
1034                                                         ps= (PixStr *) POINTER_FROM_CODE(*rd);
1035                                                         aantal= 0;
1036                                                         accol[0]=accol[1]=accol[2]=accol[3]= 0;
1037                                                         while(ps) {
1038                                                                 aantalm= count_mask(ps->mask);
1039                                                                 aantal+= aantalm;
1040
1041                                                                 zz= calchalo_z(har, ps->z);
1042                                                                 if(zz> har->zs) {
1043                                                                         *( (int *)col )= 0;
1044                                                                         RE_shadehalo(har, col, zz, dist, xn, yn, flarec);
1045                                                                         accol[0]+= aantalm*col[0];
1046                                                                         accol[1]+= aantalm*col[1];
1047                                                                         accol[2]+= aantalm*col[2];
1048                                                                         accol[3]+= aantalm*col[3];
1049                                                                         flarec= 0;
1050                                                                 }
1051
1052                                                                 ps= ps->next;
1053                                                         }
1054                                                         ps= (PixStr *) POINTER_FROM_CODE(*rd);
1055                                                         aantal= R.osa-aantal;
1056                                                         
1057                                                         zz= calchalo_z(har, *rz);
1058                                                         if(zz> har->zs) {
1059                                                                 *( (int *)col )= 0;
1060                                                                 RE_shadehalo(har, col, zz, dist, xn, yn, flarec);
1061                                                                 accol[0]+= aantal*col[0];
1062                                                                 accol[1]+= aantal*col[1];
1063                                                                 accol[2]+= aantal*col[2];
1064                                                                 accol[3]+= aantal*col[3];
1065                                                         }
1066
1067
1068                                                         col[0]= accol[0]/R.osa;
1069                                                         col[1]= accol[1]/R.osa;
1070                                                         col[2]= accol[2]/R.osa;
1071                                                         col[3]= accol[3]/R.osa;
1072
1073                                                         /* if(behind > (R.osa>>1)) addalphaUnder(rt,col); */
1074                                                         RE_addalphaAddfac((char *)rt, (char *)col, har->add);
1075                                                 }
1076                                         }
1077                                         else {
1078                                                 zz= calchalo_z(har, *rz);
1079                                                 if(zz> har->zs) {
1080                                                         xn= x- har->xs;
1081                                                         xsq= xn*xn;
1082                                                         dist= xsq+ysq;
1083                                                         if(dist<har->radsq) {
1084                                                                 RE_shadehalo(har, col, zz, dist, xn, yn, flarec);
1085                                                                 RE_addalphaAddfac((char *)rt, (char *)col, har->add);
1086                                                         }
1087                                                 }
1088                                         }
1089                                         rt++;
1090                                         rz++;
1091                                         rd++;
1092                                 }
1093                                 
1094                         }
1095                 }
1096         }
1097
1098 }
1099
1100 void scanlinehalo(unsigned int *rectz, unsigned int *rectt, short ys)
1101 {
1102         HaloRen *har = NULL;
1103         float dist,xsq,ysq,xn,yn;
1104         unsigned int a, *rz, *rt, zz;
1105         short minx,maxx, x;
1106         char col[4];
1107
1108         for(a=0; a<R.tothalo; a++) {
1109                 if((a & 255)==0) har= R.bloha[a>>8];
1110                 else har++;
1111
1112                 if(RE_local_test_break() ) break; 
1113
1114                 if(ys>har->maxy);
1115                 else if(ys<har->miny);
1116                 else {
1117                         minx= floor(har->xs-har->rad);
1118                         maxx= ceil(har->xs+har->rad);
1119                         
1120                         if(maxx<0);
1121                         else if(R.rectx<minx);
1122                         else {
1123                                 if(minx<0) minx= 0;
1124                                 if(maxx>=R.rectx) maxx= R.rectx-1;
1125
1126                                 rt= rectt+minx;
1127                                 rz= rectz+minx;
1128
1129                                 yn= (ys-har->ys)*R.ycor;
1130                                 ysq= yn*yn;
1131                                 for(x=minx; x<=maxx; x++) {
1132                                 
1133                                         zz= calchalo_z(har, *rz);
1134                                         if(zz> har->zs) {
1135                                                 xn= x- har->xs;
1136                                                 xsq= xn*xn;
1137                                                 dist= xsq+ysq;
1138                                                 if(dist<har->radsq) {
1139                                                         RE_shadehalo(har, col, zz, dist, xn, yn, har->flarec);
1140                                                         RE_addalphaAddfac((char *)rt, (char *)col, har->add);
1141                                                 }
1142                                         }
1143
1144                                         rt++;
1145                                         rz++;
1146                                 }
1147                         }
1148                 }
1149         }
1150
1151 }
1152
1153 void halovert()
1154 {
1155         HaloRen *har = NULL;
1156         float dist, xsq, ysq, xn, yn;
1157         unsigned int a, *rectz, *rz, *rectt, *rt, zz;
1158         int minx, maxx, miny, maxy, x, y;
1159         char col[4];
1160
1161
1162         for(a=0;a<R.tothalo;a++) {
1163                 if((a & 255)==0) har= R.bloha[a>>8];
1164                 else har++;
1165
1166                 if(RE_local_test_break() ) break; 
1167
1168                 if(har->maxy<0);
1169                 else if(R.recty<har->miny);
1170                 else {
1171                         minx= floor(har->xs-har->rad);
1172                         maxx= ceil(har->xs+har->rad);
1173                         
1174                         if(maxx<0);
1175                         else if(R.rectx<minx);
1176                         else {
1177                         
1178                                 miny= floor(har->ys-har->rad);
1179                                 maxy= ceil(har->ys+har->rad);
1180
1181                                 if(minx<0) minx= 0;
1182                                 if(maxx>=R.rectx) maxx= R.rectx-1;
1183                                 if(miny<0) miny= 0;
1184                                 if(maxy>R.recty) maxy= R.recty;
1185
1186                                 rectt= R.rectot+ R.rectx*miny;
1187                                 rectz= R.rectz+ R.rectx*miny;
1188
1189                                 for(y=miny;y<maxy;y++) {
1190
1191                                         rz= rectz+minx;
1192
1193                                         rt= (rectt+minx);
1194
1195                                         yn= (y - har->ys)*R.ycor;
1196                                         ysq= yn*yn;
1197                                         for(x=minx;x<=maxx;x++) {
1198                                                 
1199                                                 zz= calchalo_z(har, *rz);
1200                                                 
1201                                                 if(zz> har->zs) {
1202
1203                                                         xn= x - har->xs;
1204                                                         xsq= xn*xn;
1205                                                         dist= xsq+ysq;
1206                                                         if(dist<har->radsq) {
1207                                                                 RE_shadehalo(har, col, zz, dist, xn, yn, har->flarec);
1208                                                                 
1209                                                                 RE_addalphaAddfac((char *)rt, (char *)col, har->add);
1210                                                                 
1211                                                         }
1212                                                 }
1213                                                 rt++;
1214                                                 rz++;
1215                                         }
1216
1217                                         rectt+= R.rectx;
1218                                         rectz+= R.rectx;
1219                                         
1220                                         if(RE_local_test_break() ) break; 
1221                                 }
1222
1223                         }
1224                 }
1225         }
1226 }
1227
1228 /* ---------------- shaders ----------------------- */
1229
1230 static double Normalise_d(double *n)
1231 {
1232         double d;
1233         
1234         d= n[0]*n[0]+n[1]*n[1]+n[2]*n[2];
1235
1236         if(d>0.00000000000000001) {
1237                 d= sqrt(d);
1238
1239                 n[0]/=d; 
1240                 n[1]/=d; 
1241                 n[2]/=d;
1242         } else {
1243                 n[0]=n[1]=n[2]= 0.0;
1244                 d= 0.0;
1245         }
1246         return d;
1247 }
1248
1249
1250 /* Stoke's form factor. Need doubles here for extreme small area sizes */
1251 static float area_lamp_energy(float *co, float *vn, LampRen *lar)
1252 {
1253         double fac;
1254         double vec[4][3];       /* vectors of rendered co to vertices lamp */
1255         double cross[4][3];     /* cross products of this */
1256         double rad[4];          /* angles between vecs */
1257
1258         VECSUB(vec[0], co, lar->area[0]);
1259         VECSUB(vec[1], co, lar->area[1]);
1260         VECSUB(vec[2], co, lar->area[2]);
1261         VECSUB(vec[3], co, lar->area[3]);
1262         
1263         Normalise_d(vec[0]);
1264         Normalise_d(vec[1]);
1265         Normalise_d(vec[2]);
1266         Normalise_d(vec[3]);
1267
1268         /* cross product */
1269         CROSS(cross[0], vec[0], vec[1]);
1270         CROSS(cross[1], vec[1], vec[2]);
1271         CROSS(cross[2], vec[2], vec[3]);
1272         CROSS(cross[3], vec[3], vec[0]);
1273
1274         Normalise_d(cross[0]);
1275         Normalise_d(cross[1]);
1276         Normalise_d(cross[2]);
1277         Normalise_d(cross[3]);
1278
1279         /* angles */
1280         rad[0]= vec[0][0]*vec[1][0]+ vec[0][1]*vec[1][1]+ vec[0][2]*vec[1][2];
1281         rad[1]= vec[1][0]*vec[2][0]+ vec[1][1]*vec[2][1]+ vec[1][2]*vec[2][2];
1282         rad[2]= vec[2][0]*vec[3][0]+ vec[2][1]*vec[3][1]+ vec[2][2]*vec[3][2];
1283         rad[3]= vec[3][0]*vec[0][0]+ vec[3][1]*vec[0][1]+ vec[3][2]*vec[0][2];
1284
1285         rad[0]= acos(rad[0]);
1286         rad[1]= acos(rad[1]);
1287         rad[2]= acos(rad[2]);
1288         rad[3]= acos(rad[3]);
1289
1290         /* Stoke formula */
1291         VECMUL(cross[0], rad[0]);
1292         VECMUL(cross[1], rad[1]);
1293         VECMUL(cross[2], rad[2]);
1294         VECMUL(cross[3], rad[3]);
1295
1296         fac=  vn[0]*cross[0][0]+ vn[1]*cross[0][1]+ vn[2]*cross[0][2];
1297         fac+= vn[0]*cross[1][0]+ vn[1]*cross[1][1]+ vn[2]*cross[1][2];
1298         fac+= vn[0]*cross[2][0]+ vn[1]*cross[2][1]+ vn[2]*cross[2][2];
1299         fac+= vn[0]*cross[3][0]+ vn[1]*cross[3][1]+ vn[2]*cross[3][2];
1300
1301         if(fac<=0.0) return 0.0;
1302         return pow(fac*lar->areasize, lar->k);  // corrected for buttons size and lar->dist^2
1303 }
1304
1305
1306 float spec(float inp, int hard) 
1307 {
1308         float b1;
1309         
1310         if(inp>=1.0) return 1.0;
1311         else if (inp<=0.0) return 0.0;
1312         
1313         b1= inp*inp;
1314         /* avoid FPE */
1315         if(b1<0.01) b1= 0.01;   
1316         
1317         if((hard & 1)==0)  inp= 1.0;
1318         if(hard & 2)  inp*= b1;
1319         b1*= b1;
1320         if(hard & 4)  inp*= b1;
1321         b1*= b1;
1322         if(hard & 8)  inp*= b1;
1323         b1*= b1;
1324         if(hard & 16) inp*= b1;
1325         b1*= b1;
1326
1327         /* avoid FPE */
1328         if(b1<0.001) b1= 0.0;   
1329
1330         if(hard & 32) inp*= b1;
1331         b1*= b1;
1332         if(hard & 64) inp*=b1;
1333         b1*= b1;
1334         if(hard & 128) inp*=b1;
1335
1336         if(b1<0.001) b1= 0.0;   
1337
1338         if(hard & 256) {
1339                 b1*= b1;
1340                 inp*=b1;
1341         }
1342
1343         return inp;
1344 }
1345
1346 float Phong_Spec( float *n, float *l, float *v, int hard )
1347 {
1348         float h[3];
1349         float rslt;
1350
1351         h[0] = l[0] + v[0];
1352         h[1] = l[1] + v[1];
1353         h[2] = l[2] + v[2];
1354         Normalise(h);
1355
1356         rslt = h[0]*n[0] + h[1]*n[1] + h[2]*n[2];
1357
1358         if( rslt > 0.0 ) rslt= spec(rslt, hard);
1359         else rslt = 0.0;
1360
1361         return rslt;
1362 }
1363
1364
1365 /* reduced cook torrance spec (for off-specular peak) */
1366 float CookTorr_Spec(float *n, float *l, float *v, int hard)
1367 {
1368         float i, nh, nv, h[3];
1369
1370         h[0]= v[0]+l[0];
1371         h[1]= v[1]+l[1];
1372         h[2]= v[2]+l[2];
1373         Normalise(h);
1374
1375         nh= n[0]*h[0]+n[1]*h[1]+n[2]*h[2];
1376         if(nh<0.0) return 0.0;
1377         nv= n[0]*v[0]+n[1]*v[1]+n[2]*v[2];
1378         if(nv<0.0) nv= 0.0;
1379
1380         i= spec(nh, hard);
1381
1382         i= i/(0.1+nv);
1383         return i;
1384 }
1385
1386 /* Blinn spec */
1387 float Blinn_Spec(float *n, float *l, float *v, float refrac, float spec_power )
1388 {
1389         float i, nh, nv, nl, vh, h[3];
1390         float a, b, c, g=0.0, p, f, ang;
1391
1392         if(refrac < 1.0) return 0.0;
1393         if(spec_power == 0.0) return 0.0;
1394         
1395         /* conversion from 'hardness' (1-255) to 'spec_power' (50 maps at 0.1) */
1396         if(spec_power<100.0)
1397                 spec_power= sqrt(1.0/spec_power);
1398         else spec_power= 10.0/spec_power;
1399         
1400         h[0]= v[0]+l[0];
1401         h[1]= v[1]+l[1];
1402         h[2]= v[2]+l[2];
1403         Normalise(h);
1404
1405         nh= n[0]*h[0]+n[1]*h[1]+n[2]*h[2]; /* Dot product between surface normal and half-way vector */
1406
1407         if(nh<0.0) return 0.0;
1408
1409         nv= n[0]*v[0]+n[1]*v[1]+n[2]*v[2]; /* Dot product between surface normal and view vector */
1410
1411         if(nv<=0.0) nv= 0.01;
1412
1413         nl= n[0]*l[0]+n[1]*l[1]+n[2]*l[2]; /* Dot product between surface normal and light vector */
1414
1415         if(nl<=0.0) {
1416                 nl= 0.0;
1417                 return 0.0;
1418         }
1419
1420         vh= v[0]*h[0]+v[1]*h[1]+v[2]*h[2]; /* Dot product between view vector and half-way vector */
1421         if(vh<=0.0) vh= 0.01;
1422
1423         a = 1.0;
1424         b = (2.0*nh*nv)/vh;
1425         c = (2.0*nh*nl)/vh;
1426
1427         if( a < b && a < c ) g = a;
1428         else if( b < a && b < c ) g = b;
1429         else if( c < a && c < b ) g = c;
1430
1431         p = sqrt( (double)((refrac * refrac)+(vh*vh)-1.0) );
1432         f = (((p-vh)*(p-vh))/((p+vh)*(p+vh)))*(1+((((vh*(p+vh))-1.0)*((vh*(p+vh))-1.0))/(((vh*(p-vh))+1.0)*((vh*(p-vh))+1.0))));
1433         ang = (float)acos((double)(nh));
1434
1435         i= f * g * exp((double)(-(ang*ang) / (2.0*spec_power*spec_power)));
1436
1437         return i;
1438 }
1439
1440 /* cartoon render spec */
1441 float Toon_Spec( float *n, float *l, float *v, float size, float smooth )
1442 {
1443         float h[3];
1444         float ang;
1445         float rslt;
1446         
1447         h[0] = l[0] + v[0];
1448         h[1] = l[1] + v[1];
1449         h[2] = l[2] + v[2];
1450         Normalise(h);
1451         
1452         rslt = h[0]*n[0] + h[1]*n[1] + h[2]*n[2];
1453         
1454         ang = acos( rslt ); 
1455         
1456         if( ang < size ) rslt = 1.0;
1457         else if( ang >= (size + smooth) || smooth == 0.0 ) rslt = 0.0;
1458         else rslt = 1.0 - ((ang - size) / smooth);
1459         
1460         return rslt;
1461 }
1462
1463 /* cartoon render diffuse */
1464 float Toon_Diff( float *n, float *l, float *v, float size, float smooth )
1465 {
1466         float rslt, ang;
1467
1468         rslt = n[0]*l[0] + n[1]*l[1] + n[2]*l[2];
1469
1470         ang = acos( (double)(rslt) );
1471
1472         if( ang < size ) rslt = 1.0;
1473         else if( ang >= (size + smooth) || smooth == 0.0 ) rslt = 0.0;
1474         else rslt = 1.0 - ((ang - size) / smooth);
1475
1476         return rslt;
1477 }
1478
1479 /* Oren Nayar diffuse */
1480
1481 /* 'nl' is either dot product, or return value of area light */
1482 /* in latter case, only last multiplication uses 'nl' */
1483 float OrenNayar_Diff_i(float nl, float *n, float *l, float *v, float rough )
1484 {
1485         float i, nh, nv, vh, realnl, h[3];
1486         float a, b, t, A, B;
1487         float Lit_A, View_A, Lit_B[3], View_B[3];
1488         
1489         h[0]= v[0]+l[0];
1490         h[1]= v[1]+l[1];
1491         h[2]= v[2]+l[2];
1492         Normalise(h);
1493         
1494         nh= n[0]*h[0]+n[1]*h[1]+n[2]*h[2]; /* Dot product between surface normal and half-way vector */
1495         if(nh<0.0) nh = 0.0;
1496         
1497         nv= n[0]*v[0]+n[1]*v[1]+n[2]*v[2]; /* Dot product between surface normal and view vector */
1498         if(nv<=0.0) nv= 0.0;
1499         
1500         realnl= n[0]*l[0]+n[1]*l[1]+n[2]*l[2]; /* Dot product between surface normal and light vector */
1501         if(realnl<=0.0) return 0.0;
1502         if(nl<0.0) return 0.0;          /* value from area light */
1503         
1504         vh= v[0]*h[0]+v[1]*h[1]+v[2]*h[2]; /* Dot product between view vector and halfway vector */
1505         if(vh<=0.0) vh= 0.0;
1506         
1507         Lit_A = acos(realnl);
1508         View_A = acos( nv );
1509         
1510         Lit_B[0] = l[0] - (realnl * n[0]);
1511         Lit_B[1] = l[1] - (realnl * n[1]);
1512         Lit_B[2] = l[2] - (realnl * n[2]);
1513         Normalise( Lit_B );
1514         
1515         View_B[0] = v[0] - (nv * n[0]);
1516         View_B[1] = v[1] - (nv * n[1]);
1517         View_B[2] = v[2] - (nv * n[2]);
1518         Normalise( View_B );
1519         
1520         t = Lit_B[0]*View_B[0] + Lit_B[1]*View_B[1] + Lit_B[2]*View_B[2];
1521         if( t < 0 ) t = 0;
1522         
1523         if( Lit_A > View_A ) {
1524                 a = Lit_A;
1525                 b = View_A;
1526         }
1527         else {
1528                 a = View_A;
1529                 b = Lit_A;
1530         }
1531         
1532         A = 1 - (0.5 * ((rough * rough) / ((rough * rough) + 0.33)));
1533         B = 0.45 * ((rough * rough) / ((rough * rough) + 0.09));
1534         
1535         b*= 0.95;       /* prevent tangens from shooting to inf, 'nl' can be not a dot product here. */
1536                                 /* overflow only happens with extreme size area light, and higher roughness */
1537         i = nl * ( A + ( B * t * sin(a) * tan(b) ) );
1538         
1539         return i;
1540 }
1541
1542 /* Oren Nayar diffuse */
1543 float OrenNayar_Diff(float *n, float *l, float *v, float rough )
1544 {
1545         float nl= n[0]*l[0] + n[1]*l[1] + n[2]*l[2];
1546         return OrenNayar_Diff_i(nl, n, l, v, rough);
1547 }
1548
1549
1550 /* --------------------------------------------- */
1551
1552 void calc_R_ref(ShadeInput *shi)
1553 {
1554         float i;
1555
1556         /* shi->vn dot shi->view */
1557         i= -2*(shi->vn[0]*shi->view[0]+shi->vn[1]*shi->view[1]+shi->vn[2]*shi->view[2]);
1558
1559         shi->ref[0]= (shi->view[0]+i*shi->vn[0]);
1560         shi->ref[1]= (shi->view[1]+i*shi->vn[1]);
1561         shi->ref[2]= (shi->view[2]+i*shi->vn[2]);
1562         if(shi->osatex) {
1563                 if(shi->vlr->flag & R_SMOOTH) {
1564                         i= -2*( (shi->vn[0]+O.dxno[0])*(shi->view[0]+O.dxview) +
1565                                 (shi->vn[1]+O.dxno[1])*shi->view[1]+ (shi->vn[2]+O.dxno[2])*shi->view[2] );
1566
1567                         O.dxref[0]= shi->ref[0]- ( shi->view[0]+O.dxview+i*(shi->vn[0]+O.dxno[0]));
1568                         O.dxref[1]= shi->ref[1]- (shi->view[1]+ i*(shi->vn[1]+O.dxno[1]));
1569                         O.dxref[2]= shi->ref[2]- (shi->view[2]+ i*(shi->vn[2]+O.dxno[2]));
1570
1571                         i= -2*( (shi->vn[0]+O.dyno[0])*shi->view[0]+
1572                                 (shi->vn[1]+O.dyno[1])*(shi->view[1]+O.dyview)+ (shi->vn[2]+O.dyno[2])*shi->view[2] );
1573
1574                         O.dyref[0]= shi->ref[0]- (shi->view[0]+ i*(shi->vn[0]+O.dyno[0]));
1575                         O.dyref[1]= shi->ref[1]- (shi->view[1]+O.dyview+i*(shi->vn[1]+O.dyno[1]));
1576                         O.dyref[2]= shi->ref[2]- (shi->view[2]+ i*(shi->vn[2]+O.dyno[2]));
1577
1578                 }
1579                 else {
1580
1581                         i= -2*( shi->vn[0]*(shi->view[0]+O.dxview) +
1582                                 shi->vn[1]*shi->view[1]+ shi->vn[2]*shi->view[2] );
1583
1584                         O.dxref[0]= shi->ref[0]- (shi->view[0]+O.dxview+i*shi->vn[0]);
1585                         O.dxref[1]= shi->ref[1]- (shi->view[1]+ i*shi->vn[1]);
1586                         O.dxref[2]= shi->ref[2]- (shi->view[2]+ i*shi->vn[2]);
1587
1588                         i= -2*( shi->vn[0]*shi->view[0]+
1589                                 shi->vn[1]*(shi->view[1]+O.dyview)+ shi->vn[2]*shi->view[2] );
1590
1591                         O.dyref[0]= shi->ref[0]- (shi->view[0]+ i*shi->vn[0]);
1592                         O.dyref[1]= shi->ref[1]- (shi->view[1]+O.dyview+i*shi->vn[1]);
1593                         O.dyref[2]= shi->ref[2]- (shi->view[2]+ i*shi->vn[2]);
1594                 }
1595         }
1596
1597 }
1598
1599 /* mix of 'real' fresnel and allowing control. grad defines blending gradient */
1600 float fresnel_fac(float *view, float *vn, float grad, float fac)
1601 {
1602         float t1, t2;
1603         
1604         if(fac==0.0) return 1.0;
1605         
1606         t1= (view[0]*vn[0] + view[1]*vn[1] + view[2]*vn[2]);
1607         if(t1>0.0)  t2= 1.0+t1;
1608         else t2= 1.0-t1;
1609         
1610         t2= grad + (1.0-grad)*pow(t2, fac);
1611
1612         if(t2<0.0) return 0.0;
1613         else if(t2>1.0) return 1.0;
1614         return t2;
1615 }
1616
1617 void shade_color(ShadeInput *shi, ShadeResult *shr)
1618 {
1619         Material *ma= shi->matren;
1620
1621         if(ma->mode & (MA_VERTEXCOLP|MA_FACETEXTURE)) {
1622                 ma->r= shi->vcol[0];
1623                 ma->g= shi->vcol[1];
1624                 ma->b= shi->vcol[2];
1625         }
1626         
1627         ma->alpha= shi->mat->alpha;     // copy to render material, for fresnel and spectra
1628
1629         if(ma->texco) {
1630                 if(ma->mode & (MA_VERTEXCOLP|MA_FACETEXTURE)) {
1631                         shi->mat->r= shi->vcol[0];
1632                         shi->mat->g= shi->vcol[1];
1633                         shi->mat->b= shi->vcol[2];
1634                 }
1635                 do_material_tex(shi);
1636         }
1637
1638         if(ma->mode & (MA_ZTRA|MA_RAYTRANSP)) {
1639                 if(ma->fresnel_tra!=1.0) 
1640                         ma->alpha*= fresnel_fac(shi->view, shi->vn, ma->fresnel_tra_i, ma->fresnel_tra);
1641         }
1642
1643         shr->diff[0]= ma->r;
1644         shr->diff[1]= ma->g;
1645         shr->diff[2]= ma->b;
1646         shr->alpha= ma->alpha;
1647 }
1648
1649 /* ramp for at end of shade */
1650 void ramp_diffuse_result(float *diff, ShadeInput *shi)
1651 {
1652         Material *ma= shi->matren;
1653         float col[4], fac=0;
1654
1655         if(ma->ramp_col) {
1656                 if(ma->rampin_col==MA_RAMP_IN_RESULT) {
1657                         
1658                         fac= 0.3*diff[0] + 0.58*diff[1] + 0.12*diff[2];
1659                         do_colorband(ma->ramp_col, fac, col);
1660                         
1661                         /* blending method */
1662                         fac= col[3]*ma->rampfac_col;
1663                         switch(ma->rampblend_col) {
1664                         case MA_RAMP_BLEND:
1665                                 diff[0]= (1.0-fac)*diff[0] + fac*col[0];
1666                                 diff[1]= (1.0-fac)*diff[1] + fac*col[1];
1667                                 diff[2]= (1.0-fac)*diff[2] + fac*col[2];
1668                                 break;
1669                         case MA_RAMP_ADD:
1670                                 diff[0]+= fac*col[0];
1671                                 diff[1]+= fac*col[1];
1672                                 diff[2]+= fac*col[2];
1673                                 break;
1674                         case MA_RAMP_MULT:
1675                                 diff[0]*= (1.0-fac + fac*col[0]);
1676                                 diff[1]*= (1.0-fac + fac*col[1]);
1677                                 diff[2]*= (1.0-fac + fac*col[2]);
1678                                 break;
1679                         case MA_RAMP_SUB:
1680                                 diff[0]-= fac*col[0];
1681                                 diff[1]-= fac*col[1];
1682                                 diff[2]-= fac*col[2];
1683                                 break;
1684                         }
1685                 }
1686         }
1687 }
1688
1689 /* r,g,b denote energy, ramp is used with different values to make new material color */
1690 void add_to_diffuse(float *diff, ShadeInput *shi, float is, float r, float g, float b)
1691 {
1692         Material *ma= shi->matren;
1693         float col[4], fac=0;
1694         
1695         if(ma->ramp_col && (ma->mode & MA_RAMP_COL)) {
1696                 
1697                 /* MA_RAMP_IN_RESULT is exceptional */
1698                 if(ma->rampin_col==MA_RAMP_IN_RESULT) {
1699                         // normal add
1700                         diff[0] += r * ma->r;
1701                         diff[1] += g * ma->g;
1702                         diff[2] += b * ma->b;
1703                 }
1704                 else {
1705                         /* input */
1706                         switch(ma->rampin_col) {
1707                         case MA_RAMP_IN_ENERGY:
1708                                 fac= 0.3*r + 0.58*g + 0.12*b;
1709                                 break;
1710                         case MA_RAMP_IN_SHADER:
1711                                 fac= is;
1712                                 break;
1713                         case MA_RAMP_IN_NOR:
1714                                 fac= shi->view[0]*shi->vn[0] + shi->view[1]*shi->vn[1] + shi->view[2]*shi->vn[2];
1715                                 break;
1716                         }
1717         
1718                         do_colorband(ma->ramp_col, fac, col);
1719                         
1720                         /* blending method */
1721                         fac= col[3]*ma->rampfac_col;
1722                         switch(ma->rampblend_col) {
1723                         case MA_RAMP_BLEND:
1724                                 col[0]= (1.0-fac)*ma->r + fac*col[0];
1725                                 col[1]= (1.0-fac)*ma->g + fac*col[1];
1726                                 col[2]= (1.0-fac)*ma->b + fac*col[2];
1727                                 break;
1728                         case MA_RAMP_ADD:
1729                                 col[0]= ma->r + fac*col[0];
1730                                 col[1]= ma->g + fac*col[1];
1731                                 col[2]= ma->b + fac*col[2];
1732                                 break;
1733                         case MA_RAMP_MULT:
1734                                 col[0]= ma->r*(1.0-fac + fac*col[0]);
1735                                 col[1]= ma->g*(1.0-fac + fac*col[1]);
1736                                 col[2]= ma->b*(1.0-fac + fac*col[2]);
1737                                 break;
1738                         case MA_RAMP_SUB:
1739                                 col[0]= ma->r - fac*col[0];
1740                                 col[1]= ma->g - fac*col[1];
1741                                 col[2]= ma->b - fac*col[2];
1742                                 break;
1743                         }
1744                         
1745                         /* output to */
1746                         diff[0] += r * col[0];
1747                         diff[1] += g * col[1];
1748                         diff[2] += b * col[2];
1749                 }
1750         }
1751         else {
1752                 diff[0] += r * ma->r;
1753                 diff[1] += g * ma->g;
1754                 diff[2] += b * ma->b;
1755         }
1756 }
1757
1758 void ramp_spec_result(float *specr, float *specg, float *specb, ShadeInput *shi)
1759 {
1760         Material *ma= shi->matren;
1761         float col[4];
1762         float fac;
1763         
1764         if(ma->ramp_spec && (ma->rampin_spec==MA_RAMP_IN_RESULT)) {
1765                 fac= 0.3*(*specr) + 0.58*(*specg) + 0.12*(*specb);
1766                 do_colorband(ma->ramp_spec, fac, col);
1767                 
1768                 /* blending method */
1769                 fac= col[3]*ma->rampfac_spec;
1770                 switch(ma->rampblend_spec) {
1771                 case MA_RAMP_BLEND:
1772                         *specr= (1.0-fac)*(*specr) + fac*col[0];
1773                         *specg= (1.0-fac)*(*specg) + fac*col[1];
1774                         *specb= (1.0-fac)*(*specb) + fac*col[2];
1775                         break;
1776                 case MA_RAMP_ADD:
1777                         *specr += fac*col[0];
1778                         *specg += fac*col[1];
1779                         *specb += fac*col[2];
1780                         break;
1781                 case MA_RAMP_MULT:
1782                         *specr *= (1.0-fac + fac*col[0]);
1783                         *specg *= (1.0-fac + fac*col[1]);
1784                         *specb *= (1.0-fac + fac*col[2]);
1785                         break;
1786                 case MA_RAMP_SUB:
1787                         *specr -= fac*col[0];
1788                         *specg -= fac*col[1];
1789                         *specb -= fac*col[2];
1790                         break;
1791                 }
1792                 
1793         }
1794 }
1795
1796 /* is = dot product shade, t = spec energy */
1797 void do_specular_ramp(ShadeInput *shi, float is, float t, float *spec)
1798 {
1799         Material *ma= shi->matren;
1800         float col[4];
1801         float fac=0.0;
1802         
1803         /* MA_RAMP_IN_RESULT is exception */
1804         if(ma->ramp_spec && (ma->rampin_spec!=MA_RAMP_IN_RESULT)) {
1805                 
1806                 /* input */
1807                 switch(ma->rampin_spec) {
1808                 case MA_RAMP_IN_ENERGY:
1809                         fac= t;
1810                         break;
1811                 case MA_RAMP_IN_SHADER:
1812                         fac= is;
1813                         break;
1814                 case MA_RAMP_IN_NOR:
1815                         fac= shi->view[0]*shi->vn[0] + shi->view[1]*shi->vn[1] + shi->view[2]*shi->vn[2];
1816                         break;
1817                 }
1818                 
1819                 do_colorband(ma->ramp_spec, fac, col);
1820                 
1821                 /* blending method */
1822                 fac= col[3]*ma->rampfac_spec;
1823                 switch(ma->rampblend_spec) {
1824                 case MA_RAMP_BLEND:
1825                         spec[0]= (1.0-fac)*ma->specr + fac*col[0];
1826                         spec[1]= (1.0-fac)*ma->specg + fac*col[1];
1827                         spec[2]= (1.0-fac)*ma->specb + fac*col[2];
1828                         break;
1829                 case MA_RAMP_ADD:
1830                         spec[0]= ma->specr + fac*col[0];
1831                         spec[1]= ma->specg + fac*col[1];
1832                         spec[2]= ma->specb + fac*col[2];
1833                         break;
1834                 case MA_RAMP_MULT:
1835                         spec[0]= ma->specr*fac*col[0];
1836                         spec[1]= ma->specg*fac*col[1];
1837                         spec[2]= ma->specb*fac*col[2];
1838                         break;
1839                 case MA_RAMP_SUB:
1840                         spec[0]= ma->specr - fac*col[0];
1841                         spec[1]= ma->specg - fac*col[1];
1842                         spec[2]= ma->specb - fac*col[2];
1843                         break;
1844                 }
1845         }
1846         else {
1847                 spec[0]= ma->specr;
1848                 spec[1]= ma->specg;
1849                 spec[2]= ma->specb;
1850         }
1851 }
1852
1853
1854
1855 static void ambient_occlusion(World *wrld, ShadeInput *shi, ShadeResult *shr)
1856 {
1857         float f, shadfac[4];
1858         
1859         if((wrld->mode & WO_AMB_OCC) && (R.r.mode & R_RAYTRACE) && shi->matren->amb!=0.0) {
1860                 ray_ao(shi, wrld, shadfac);
1861
1862                 if(wrld->aocolor==WO_AOPLAIN) {
1863                         if (wrld->aomix==WO_AOADDSUB) shadfac[3] = 2.0*shadfac[3]-1.0;
1864                         else if (wrld->aomix==WO_AOSUB) shadfac[3] = shadfac[3]-1.0;
1865
1866                         f= wrld->aoenergy*shadfac[3]*shi->matren->amb;
1867                         add_to_diffuse(shr->diff, shi, f, f, f, f);
1868                         //shr->diff[0] += f;
1869                         //shr->diff[1] += f;
1870                         //shr->diff[2] += f;
1871                 }
1872                 else {
1873                         if (wrld->aomix==WO_AOADDSUB) {
1874                                 shadfac[0] = 2.0*shadfac[0]-1.0;
1875                                 shadfac[1] = 2.0*shadfac[1]-1.0;
1876                                 shadfac[2] = 2.0*shadfac[2]-1.0;
1877                         }
1878                         else if (wrld->aomix==WO_AOSUB) {
1879                                 shadfac[0] = shadfac[0]-1.0;
1880                                 shadfac[1] = shadfac[1]-1.0;
1881                                 shadfac[2] = shadfac[2]-1.0;
1882                         }
1883                         f= wrld->aoenergy*shi->matren->amb;
1884                         add_to_diffuse(shr->diff, shi, f, f*shadfac[0], f*shadfac[1], f*shadfac[2]);
1885                         //shr->diff[0] += f*shadfac[0];
1886                         //shr->diff[1] += f*shadfac[1];
1887                         //shr->diff[2] += f*shadfac[2];
1888                 }
1889         }
1890 }
1891
1892 void shade_lamp_loop(ShadeInput *shi, ShadeResult *shr)
1893 {
1894         LampRen *lar;
1895         Material *ma;
1896         float i, inp, inpr, is, t, lv[3], lampdist, ld = 0;
1897         float lvrot[3], *vn, *view, shadfac[4], soft;   // shadfac = rgba
1898         int a;
1899
1900         vn= shi->vn;
1901         view= shi->view;
1902         ma= shi->matren;
1903         
1904         memset(shr, 0, sizeof(ShadeResult));
1905         
1906         /* separate loop */
1907         if(ma->mode & MA_ONLYSHADOW) {
1908                 float ir;
1909                 
1910                 shadfac[3]= ir= 0.0;
1911                 for(a=0; a<R.totlamp; a++) {
1912                         lar= R.la[a];
1913                         
1914                         if(lar->mode & LA_LAYER) if((lar->lay & shi->vlr->lay)==0) continue;
1915                         
1916                         lv[0]= shi->co[0]-lar->co[0];
1917                         lv[1]= shi->co[1]-lar->co[1];
1918                         lv[2]= shi->co[2]-lar->co[2];
1919
1920                         if(lar->shb) {
1921                                 /* only test within spotbundel */
1922                                 Normalise(lv);
1923                                 inpr= lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2];
1924                                 if(inpr>lar->spotsi) {
1925                                         
1926                                         inp= vn[0]*lv[0] + vn[1]*lv[1] + vn[2]*lv[2];
1927                                         
1928                                         i = testshadowbuf(lar->shb, shi->co, inp);
1929                                         
1930                                         t= inpr - lar->spotsi;
1931                                         if(t<lar->spotbl && lar->spotbl!=0.0) {
1932                                                 t/= lar->spotbl;
1933                                                 t*= t;
1934                                                 i= t*i+(1.0-t);
1935                                         }
1936                                         
1937                                         shadfac[3]+= i;
1938                                         ir+= 1.0;
1939                                 }
1940                                 else {
1941                                         shadfac[3]+= 1.0;
1942                                         ir+= 1.0;
1943                                 }
1944                         }
1945                         else if(lar->mode & LA_SHAD_RAY) {
1946                                 float shad[4];
1947                                 
1948                                 /* single sided? */
1949                                 if( shi->vlr->n[0]*lv[0] + shi->vlr->n[1]*lv[1] + shi->vlr->n[2]*lv[2] > -0.01) {
1950                                         ray_shadow(shi, lar, shad);
1951                                         shadfac[3]+= shad[3];
1952                                         ir+= 1.0;
1953                                 }
1954                         }
1955
1956                 }
1957                 if(ir>0.0) shadfac[3]/= ir;
1958                 shr->alpha= (shi->mat->alpha)*(1.0-shadfac[3]);
1959                 
1960                 return;
1961         }
1962                 
1963         if(ma->mode & (MA_VERTEXCOLP|MA_FACETEXTURE)) {
1964                 ma->r= shi->vcol[0];
1965                 ma->g= shi->vcol[1];
1966                 ma->b= shi->vcol[2];
1967         }
1968         
1969         ma->alpha= shi->mat->alpha;     // copy to render material, for fresnel and spectra
1970         
1971         /* envmap hack, always reset */
1972         shi->refcol[0]= shi->refcol[1]= shi->refcol[2]= shi->refcol[3]= 0.0;
1973
1974         if(ma->texco) {
1975                 if(ma->mode & (MA_VERTEXCOLP|MA_FACETEXTURE)) {
1976                         shi->mat->r= shi->vcol[0];
1977                         shi->mat->g= shi->vcol[1];
1978                         shi->mat->b= shi->vcol[2];
1979                 }
1980                 do_material_tex(shi);
1981         }
1982         
1983         if(ma->mode & MA_SHLESS) {
1984                 shr->diff[0]= ma->r;
1985                 shr->diff[1]= ma->g;
1986                 shr->diff[2]= ma->b;
1987                 shr->alpha= ma->alpha;
1988                 return;
1989         }
1990
1991         if( (ma->mode & (MA_VERTEXCOL+MA_VERTEXCOLP))== MA_VERTEXCOL ) {
1992                 // add_to_diffuse(shr->diff, shi, 1.0, ma->emit+shi->vcol[0], ma->emit+shi->vcol[1], ma->emit+shi->vcol[2]);
1993                 shr->diff[0]= ma->r*(ma->emit+shi->vcol[0]);
1994                 shr->diff[1]= ma->g*(ma->emit+shi->vcol[1]);
1995                 shr->diff[2]= ma->b*(ma->emit+shi->vcol[2]);
1996         }
1997         else {
1998                 // add_to_diffuse(shr->diff, shi, 1.0, ma->emit, ma->emit, ma->emit);
1999                 shr->diff[0]= ma->r*ma->emit;
2000                 shr->diff[1]= ma->g*ma->emit;
2001                 shr->diff[2]= ma->b*ma->emit;
2002         }
2003         
2004         ambient_occlusion(&R.wrld, shi, shr);
2005
2006         for(a=0; a<R.totlamp; a++) {
2007                 lar= R.la[a];
2008
2009                 /* test for lamp layer */
2010                 if(lar->mode & LA_LAYER) if((lar->lay & shi->vlr->lay)==0) continue;
2011                 
2012                 /* lampdist calculation */
2013                 if(lar->type==LA_SUN || lar->type==LA_HEMI) {
2014                         VECCOPY(lv, lar->vec);
2015                         lampdist= 1.0;
2016                 }
2017                 else {
2018                         lv[0]= shi->co[0]-lar->co[0];
2019                         lv[1]= shi->co[1]-lar->co[1];
2020                         lv[2]= shi->co[2]-lar->co[2];
2021                         ld= sqrt(lv[0]*lv[0]+lv[1]*lv[1]+lv[2]*lv[2]);
2022                         lv[0]/= ld;
2023                         lv[1]/= ld;
2024                         lv[2]/= ld;
2025                         
2026                         /* ld is re-used further on (texco's) */
2027                         if(lar->type==LA_AREA) {
2028                                 lampdist= 1.0;
2029                         }
2030                         else {
2031                                 if(lar->mode & LA_QUAD) {
2032                                         t= 1.0;
2033                                         if(lar->ld1>0.0)
2034                                                 t= lar->dist/(lar->dist+lar->ld1*ld);
2035                                         if(lar->ld2>0.0)
2036                                                 t*= lar->distkw/(lar->distkw+lar->ld2*ld*ld);
2037         
2038                                         lampdist= t;
2039                                 }
2040                                 else {
2041                                         lampdist= (lar->dist/(lar->dist+ld));
2042                                 }
2043         
2044                                 if(lar->mode & LA_SPHERE) {
2045                                         t= lar->dist - ld;
2046                                         if(t<0.0) continue;
2047                                         
2048                                         t/= lar->dist;
2049                                         lampdist*= (t);
2050                                 }
2051                         }
2052                 }
2053
2054                 if(lar->mode & LA_TEXTURE)  do_lamp_tex(lar, lv, shi);
2055
2056                 /* init transp shadow */
2057                 shadfac[3]= 1.0;
2058                 if(ma->mode & MA_SHADOW_TRA) shadfac[0]= shadfac[1]= shadfac[2]= 1.0;
2059
2060                 if(lar->type==LA_SPOT) {
2061                         
2062                         if(lar->mode & LA_SQUARE) {
2063                                 if(lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2]>0.0) {
2064                                         float x;
2065                                         
2066                                         /* rotate view to lampspace */
2067                                         VECCOPY(lvrot, lv);
2068                                         MTC_Mat3MulVecfl(lar->imat, lvrot);
2069                                         
2070                                         x= MAX2(fabs(lvrot[0]/lvrot[2]) , fabs(lvrot[1]/lvrot[2]));
2071                                         /* 1.0/(sqrt(1+x*x)) is equivalent to cos(atan(x)) */
2072
2073                                         inpr= 1.0/(sqrt(1+x*x));
2074                                 }
2075                                 else inpr= 0.0;
2076                         }
2077                         else {
2078                                 inpr= lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2];
2079                         }
2080
2081                         t= lar->spotsi;
2082                         if(inpr<t) continue;
2083                         else {
2084                                 t= inpr-t;
2085                                 i= 1.0;
2086                                 soft= 1.0;
2087                                 if(t<lar->spotbl && lar->spotbl!=0.0) {
2088                                         /* soft area */
2089                                         i= t/lar->spotbl;
2090                                         t= i*i;
2091                                         soft= (3.0*t-2.0*t*i);
2092                                         inpr*= soft;
2093                                 }
2094                                 lampdist*=inpr;
2095                         }
2096
2097                         if(lar->mode & LA_OSATEX) {
2098                                 shi->osatex= 1; /* signal for multitex() */
2099                                 
2100                                 O.dxlv[0]= lv[0] - (shi->co[0]-lar->co[0]+O.dxco[0])/ld;
2101                                 O.dxlv[1]= lv[1] - (shi->co[1]-lar->co[1]+O.dxco[1])/ld;
2102                                 O.dxlv[2]= lv[2] - (shi->co[2]-lar->co[2]+O.dxco[2])/ld;
2103
2104                                 O.dylv[0]= lv[0] - (shi->co[0]-lar->co[0]+O.dyco[0])/ld;
2105                                 O.dylv[1]= lv[1] - (shi->co[1]-lar->co[1]+O.dyco[1])/ld;
2106                                 O.dylv[2]= lv[2] - (shi->co[2]-lar->co[2]+O.dyco[2])/ld;
2107                         }
2108                         
2109                 }
2110
2111                 /* dot product and reflectivity */
2112                 /* inp = dotproduct, is = shader result, i = lamp energy (with shadow) */
2113                 
2114                 inp= vn[0]*lv[0] + vn[1]*lv[1] + vn[2]*lv[2];
2115                                 
2116                 if(lar->mode & LA_NO_DIFF) {
2117                         is= 0.0;        // skip shaders
2118                 }
2119                 else if(lar->type==LA_HEMI) {
2120                         is= 0.5*inp + 0.5;
2121                 }
2122                 else {
2123                 
2124                         if(lar->type==LA_AREA) {
2125                                 /* single sided */
2126                                 if(lv[0]*lar->vec[0]+lv[1]*lar->vec[1]+lv[2]*lar->vec[2]>0.0)
2127                                         inp= area_lamp_energy(shi->co, shi->vn, lar);
2128                                 else inp= 0.0;
2129                         }
2130                         
2131                         /* diffuse shaders (oren nayer gets inp from area light) */
2132                         if(ma->diff_shader==MA_DIFF_ORENNAYAR) is= OrenNayar_Diff_i(inp, vn, lv, view, ma->roughness);
2133                         else if(ma->diff_shader==MA_DIFF_TOON) is= Toon_Diff(vn, lv, view, ma->param[0], ma->param[1]);
2134                         else is= inp;   // Lambert
2135                 }
2136                 
2137                 i= is;
2138                 if(i>0.0) {
2139                         i*= lampdist*ma->ref;
2140                 }
2141
2142                 /* shadow and spec, (lampdist==0 outside spot) */
2143                 if(lampdist> 0.0) {
2144                         
2145                         if(i>0.0 && (R.r.mode & R_SHADOW)) {
2146                                 if(ma->mode & MA_SHADOW) {
2147                                         if(lar->type==LA_HEMI); // no shadow
2148                                         else {
2149                                                 if(lar->shb) {
2150                                                         shadfac[3] = testshadowbuf(lar->shb, shi->co, inp);
2151                                                 }
2152                                                 else if(lar->mode & LA_SHAD_RAY) {
2153                                                         // this extra 0.001 prevents boundary cases (shadow on smooth sphere)
2154                                                         if((shi->vlr->n[0]*lv[0] + shi->vlr->n[1]*lv[1] + shi->vlr->n[2]*lv[2]) > -0.001) 
2155                                                                 ray_shadow(shi, lar, shadfac);
2156                                                         else shadfac[3]= 0.0;
2157                                                 }
2158         
2159                                                 /* warning, here it skips the loop */
2160                                                 if(lar->mode & LA_ONLYSHADOW) {
2161                                                         
2162                                                         shadfac[3]= i*lar->energy*(1.0-shadfac[3]);
2163                                                         shr->diff[0] -= shadfac[3]*ma->r;
2164                                                         shr->diff[1] -= shadfac[3]*ma->g;
2165                                                         shr->diff[2] -= shadfac[3]*ma->b;
2166                                                         
2167                                                         continue;
2168                                                 }
2169                                                 
2170                                                 if(shadfac[3]==0.0) continue;
2171         
2172                                                 i*= shadfac[3];
2173                                         }
2174                                 }
2175                         }
2176                 
2177                         /* specularity */
2178                         if(shadfac[3]>0.0 && ma->spec!=0.0 && !(lar->mode & LA_NO_SPEC)) {
2179                                 
2180                                 if(lar->type==LA_HEMI) {
2181                                         /* hemi uses no spec shaders (yet) */
2182                                         
2183                                         lv[0]+= view[0];
2184                                         lv[1]+= view[1];
2185                                         lv[2]+= view[2];
2186                                         
2187                                         Normalise(lv);
2188                                         
2189                                         t= vn[0]*lv[0]+vn[1]*lv[1]+vn[2]*lv[2];
2190                                         
2191                                         if(lar->type==LA_HEMI) {
2192                                                 t= 0.5*t+0.5;
2193                                         }
2194                                         
2195                                         t= shadfac[3]*ma->spec*spec(t, ma->har);
2196                                         shr->spec[0]+= t*(lar->r * ma->specr);
2197                                         shr->spec[1]+= t*(lar->g * ma->specg);
2198                                         shr->spec[2]+= t*(lar->b * ma->specb);
2199                                 }
2200                                 else {
2201                                         /* specular shaders */
2202                                         float specfac;
2203
2204                                         if(ma->spec_shader==MA_SPEC_PHONG) 
2205                                                 specfac= Phong_Spec(vn, lv, view, ma->har);
2206                                         else if(ma->spec_shader==MA_SPEC_COOKTORR) 
2207                                                 specfac= CookTorr_Spec(vn, lv, view, ma->har);
2208                                         else if(ma->spec_shader==MA_SPEC_BLINN) 
2209                                                 specfac= Blinn_Spec(vn, lv, view, ma->refrac, (float)ma->har);
2210                                         else 
2211                                                 specfac= Toon_Spec(vn, lv, view, ma->param[2], ma->param[3]);
2212                                 
2213                                         /* area lamp correction */
2214                                         if(lar->type==LA_AREA) specfac*= inp;
2215                                         
2216                                         t= shadfac[3]*ma->spec*lampdist*specfac;
2217                                         
2218                                         if(ma->mode & MA_RAMP_SPEC) {
2219                                                 float spec[3];
2220                                                 do_specular_ramp(shi, specfac, t, spec);
2221                                                 shr->spec[0]+= t*(lar->r * spec[0]);
2222                                                 shr->spec[1]+= t*(lar->g * spec[1]);
2223                                                 shr->spec[2]+= t*(lar->b * spec[2]);
2224                                         }
2225                                         else {
2226                                                 shr->spec[0]+= t*(lar->r * ma->specr);
2227                                                 shr->spec[1]+= t*(lar->g * ma->specg);
2228                                                 shr->spec[2]+= t*(lar->b * ma->specb);
2229                                         }
2230                                 }
2231                         }
2232                 }
2233                 
2234                 /* in case 'no diffuse' we still do most calculus, spec can be in shadow */
2235                 if(i>0.0 && !(lar->mode & LA_NO_DIFF)) {
2236                         if(ma->mode & MA_SHADOW_TRA) {
2237                                 add_to_diffuse(shr->diff, shi, is, i*shadfac[0]*lar->r, i*shadfac[1]*lar->g, i*shadfac[2]*lar->b);
2238                                 //shr->diff[0]+= i*shadfac[0]*lar->r;
2239                                 //shr->diff[1]+= i*shadfac[1]*lar->g;
2240                                 //shr->diff[2]+= i*shadfac[2]*lar->b;
2241                         }
2242                         else {
2243                                 add_to_diffuse(shr->diff, shi, is, i*lar->r, i*lar->g, i*lar->b);
2244                                 //shr->diff[0]+= i*lar->r;
2245                                 //shr->diff[1]+= i*lar->g;
2246                                 //shr->diff[2]+= i*lar->b;
2247                         }
2248                 }
2249         }
2250
2251         if(ma->mode & (MA_ZTRA|MA_RAYTRANSP)) {
2252                 if(ma->fresnel_tra!=1.0) 
2253                         ma->alpha*= fresnel_fac(shi->view, shi->vn, ma->fresnel_tra_i, ma->fresnel_tra);
2254
2255                 if(ma->spectra!=0.0) {
2256
2257                         t = MAX3(shr->spec[0], shr->spec[1], shr->spec[2]);
2258                         t *= ma->spectra;
2259                         if(t>1.0) t= 1.0;
2260                         ma->alpha= (1.0-t)*ma->alpha+t;
2261                 }
2262         }
2263
2264         shr->alpha= ma->alpha;
2265
2266         if(shr->spec[0]<0.0) shr->spec[0]= 0.0;
2267         if(shr->spec[1]<0.0) shr->spec[1]= 0.0;
2268         if(shr->spec[2]<0.0) shr->spec[2]= 0.0;
2269
2270         shr->diff[0]+= ma->r*ma->amb*shi->rad[0];
2271         //shr->diff[0]*= ma->r;
2272         shr->diff[0]+= ma->ambr;
2273         if(shr->diff[0]<0.0) shr->diff[0]= 0.0;
2274         
2275         shr->diff[1]+= ma->g*ma->amb*shi->rad[1];
2276         //shr->diff[1]*= ma->g;
2277         shr->diff[1]+= ma->ambg;
2278         if(shr->diff[1]<0.0) shr->diff[1]= 0.0;
2279         
2280         shr->diff[2]+= ma->b*ma->amb*shi->rad[2];
2281         //shr->diff[2]*= ma->b;
2282         shr->diff[2]+= ma->ambb;
2283         if(shr->diff[2]<0.0) shr->diff[2]= 0.0;
2284         
2285         if(ma->mode & MA_RAMP_COL) ramp_diffuse_result(shr->diff, shi);
2286         if(ma->mode & MA_RAMP_SPEC) ramp_spec_result(shr->spec, shr->spec+1, shr->spec+2, shi);
2287         
2288         /* refcol is for envmap only */
2289         if(shi->refcol[0]!=0.0) {
2290                 shr->diff[0]= ma->mirr*shi->refcol[1] + (1.0 - ma->mirr*shi->refcol[0])*shr->diff[0];
2291                 shr->diff[1]= ma->mirg*shi->refcol[2] + (1.0 - ma->mirg*shi->refcol[0])*shr->diff[1];
2292                 shr->diff[2]= ma->mirb*shi->refcol[3] + (1.0 - ma->mirb*shi->refcol[0])*shr->diff[2];
2293         }
2294
2295 }
2296
2297 void shade_input_set_coords(ShadeInput *shi, float u, float v, int i1, int i2, int i3)
2298 {
2299         VertRen *v1, *v2, *v3;
2300         VlakRen *vlr= shi->vlr;
2301         float l, dl;
2302         short texco= shi->matren->texco;
2303         int mode= shi->matren->mode;
2304         char p1, p2, p3;
2305         
2306         /* for rendering of quads, the following values are used to denote vertices:
2307            0 1 2        scanline tria & first half quad, and ray tria
2308            0 2 3    scanline 2nd half quad
2309            0 1 3    raytracer first half quad
2310            2 1 3    raytracer 2nd half quad
2311         */
2312
2313         if(i1==0) {
2314                 v1= vlr->v1;
2315                 p1= ME_FLIPV1;
2316         } else {
2317                 v1= vlr->v3;
2318                 p1= ME_FLIPV3;
2319         }
2320
2321         if(i2==1) {
2322                 v2= vlr->v2;
2323                 p2= ME_FLIPV2;
2324         } else {
2325                 v2= vlr->v3;
2326                 p2= ME_FLIPV3;
2327         }
2328         
2329         if(i3==2) {
2330                 v3= vlr->v3;
2331                 p3= ME_FLIPV3;
2332         } else {
2333                 v3= vlr->v4;
2334                 p3= ME_FLIPV4;
2335         }
2336
2337         /* calculate U and V, for scanline (normal u and v are -1 to 0) */
2338         if(u==1.0) {
2339                 /* exception case for wire render of edge */
2340                 if(vlr->v2==vlr->v3);
2341                 else if( (vlr->flag & R_SMOOTH) || (texco & NEED_UV)) {
2342                         float detsh, t00, t10, t01, t11;
2343                         
2344                         if(vlr->snproj==0) {
2345                                 t00= v3->co[0]-v1->co[0]; t01= v3->co[1]-v1->co[1];
2346                                 t10= v3->co[0]-v2->co[0]; t11= v3->co[1]-v2->co[1];
2347                         }
2348                         else if(vlr->snproj==1) {
2349                                 t00= v3->co[0]-v1->co[0]; t01= v3->co[2]-v1->co[2];
2350                                 t10= v3->co[0]-v2->co[0]; t11= v3->co[2]-v2->co[2];
2351                         }
2352                         else {
2353                                 t00= v3->co[1]-v1->co[1]; t01= v3->co[2]-v1->co[2];
2354                                 t10= v3->co[1]-v2->co[1]; t11= v3->co[2]-v2->co[2];
2355                         }
2356                         
2357                         detsh= t00*t11-t10*t01;
2358                         t00/= detsh; t01/=detsh; 
2359                         t10/=detsh; t11/=detsh;
2360                 
2361                         if(vlr->snproj==0) {
2362                                 u= (shi->co[0]-v3->co[0])*t11-(shi->co[1]-v3->co[1])*t10;
2363                                 v= (shi->co[1]-v3->co[1])*t00-(shi->co[0]-v3->co[0])*t01;
2364                                 if(shi->osatex) {
2365                                         O.dxuv[0]=  O.dxco[0]*t11- O.dxco[1]*t10;
2366                                         O.dxuv[1]=  O.dxco[1]*t00- O.dxco[0]*t01;
2367                                         O.dyuv[0]=  O.dyco[0]*t11- O.dyco[1]*t10;
2368                                         O.dyuv[1]=  O.dyco[1]*t00- O.dyco[0]*t01;
2369                                 }
2370                         }
2371                         else if(vlr->snproj==1) {
2372                                 u= (shi->co[0]-v3->co[0])*t11-(shi->co[2]-v3->co[2])*t10;
2373                                 v= (shi->co[2]-v3->co[2])*t00-(shi->co[0]-v3->co[0])*t01;
2374                                 if(shi->osatex) {
2375                                         O.dxuv[0]=  O.dxco[0]*t11- O.dxco[2]*t10;
2376                                         O.dxuv[1]=  O.dxco[2]*t00- O.dxco[0]*t01;
2377                                         O.dyuv[0]=  O.dyco[0]*t11- O.dyco[2]*t10;
2378                                         O.dyuv[1]=  O.dyco[2]*t00- O.dyco[0]*t01;
2379                                 }
2380                         }
2381                         else {
2382                                 u= (shi->co[1]-v3->co[1])*t11-(shi->co[2]-v3->co[2])*t10;
2383                                 v= (shi->co[2]-v3->co[2])*t00-(shi->co[1]-v3->co[1])*t01;
2384                                 if(shi->osatex) {
2385                                         O.dxuv[0]=  O.dxco[1]*t11- O.dxco[2]*t10;
2386                                         O.dxuv[1]=  O.dxco[2]*t00- O.dxco[1]*t01;
2387                                         O.dyuv[0]=  O.dyco[1]*t11- O.dyco[2]*t10;
2388                                         O.dyuv[1]=  O.dyco[2]*t00- O.dyco[1]*t01;
2389                                 }
2390                         }
2391                 }
2392         
2393         }
2394         l= 1.0+u+v;
2395         
2396         /* calculate punos (vertexnormals) */
2397         if(vlr->flag & R_SMOOTH) { 
2398                 float n1[3], n2[3], n3[3];
2399                 
2400                 if(vlr->puno & p1) {
2401                         n1[0]= -v1->n[0]; n1[1]= -v1->n[1]; n1[2]= -v1->n[2];
2402                 } else {
2403                         n1[0]= v1->n[0]; n1[1]= v1->n[1]; n1[2]= v1->n[2];
2404                 }
2405                 if(vlr->puno & p2) {
2406                         n2[0]= -v2->n[0]; n2[1]= -v2->n[1]; n2[2]= -v2->n[2];
2407                 } else {
2408                         n2[0]= v2->n[0]; n2[1]= v2->n[1]; n2[2]= v2->n[2];
2409                 }
2410                 
2411                 if(vlr->puno & p3) {
2412                         n3[0]= -v3->n[0]; n3[1]= -v3->n[1]; n3[2]= -v3->n[2];
2413                 } else {
2414                         n3[0]= v3->n[0]; n3[1]= v3->n[1]; n3[2]= v3->n[2];
2415                 }
2416
2417                 shi->vn[0]= l*n3[0]-u*n1[0]-v*n2[0];
2418                 shi->vn[1]= l*n3[1]-u*n1[1]-v*n2[1];
2419                 shi->vn[2]= l*n3[2]-u*n1[2]-v*n2[2];
2420
2421                 Normalise(shi->vn);
2422
2423                 if(shi->osatex && (texco & (TEXCO_NORM|TEXCO_REFL)) ) {
2424                         dl= O.dxuv[0]+O.dxuv[1];
2425                         O.dxno[0]= dl*n3[0]-O.dxuv[0]*n1[0]-O.dxuv[1]*n2[0];
2426                         O.dxno[1]= dl*n3[1]-O.dxuv[0]*n1[1]-O.dxuv[1]*n2[1];
2427                         O.dxno[2]= dl*n3[2]-O.dxuv[0]*n1[2]-O.dxuv[1]*n2[2];
2428                         dl= O.dyuv[0]+O.dyuv[1];
2429                         O.dyno[0]= dl*n3[0]-O.dyuv[0]*n1[0]-O.dyuv[1]*n2[0];
2430                         O.dyno[1]= dl*n3[1]-O.dyuv[0]*n1[1]-O.dyuv[1]*n2[1];
2431                         O.dyno[2]= dl*n3[2]-O.dyuv[0]*n1[2]-O.dyuv[1]*n2[2];
2432
2433                 }
2434         }
2435         else {
2436                 VECCOPY(shi->vn, vlr->n);
2437         }
2438
2439         /* texture coordinates. O.dxuv O.dyuv have been set */
2440         if(texco & NEED_UV) {
2441                 if(texco & TEXCO_ORCO) {
2442                         if(v1->orco) {
2443                                 float *o1, *o2, *o3;
2444                                 
2445                                 o1= v1->orco;
2446                                 o2= v2->orco;
2447                                 o3= v3->orco;
2448                                 
2449                                 shi->lo[0]= l*o3[0]-u*o1[0]-v*o2[0];
2450                                 shi->lo[1]= l*o3[1]-u*o1[1]-v*o2[1];
2451                                 shi->lo[2]= l*o3[2]-u*o1[2]-v*o2[2];
2452         
2453                                 if(shi->osatex) {
2454                                         dl= O.dxuv[0]+O.dxuv[1];
2455                                         O.dxlo[0]= dl*o3[0]-O.dxuv[0]*o1[0]-O.dxuv[1]*o2[0];
2456                                         O.dxlo[1]= dl*o3[1]-O.dxuv[0]*o1[1]-O.dxuv[1]*o2[1];
2457                                         O.dxlo[2]= dl*o3[2]-O.dxuv[0]*o1[2]-O.dxuv[1]*o2[2];
2458                                         dl= O.dyuv[0]+O.dyuv[1];
2459                                         O.dylo[0]= dl*o3[0]-O.dyuv[0]*o1[0]-O.dyuv[1]*o2[0];
2460                                         O.dylo[1]= dl*o3[1]-O.dyuv[0]*o1[1]-O.dyuv[1]*o2[1];
2461                                         O.dylo[2]= dl*o3[2]-O.dyuv[0]*o1[2]-O.dyuv[1]*o2[2];
2462                                 }
2463                         }
2464                 }
2465                 
2466                 if(texco & TEXCO_GLOB) {
2467                         VECCOPY(shi->gl, shi->co);
2468                         MTC_Mat4MulVecfl(R.viewinv, shi->gl);
2469                         if(shi->osatex) {
2470                                 VECCOPY(O.dxgl, O.dxco);
2471                                 MTC_Mat3MulVecfl(R.imat, O.dxco);
2472                                 VECCOPY(O.dygl, O.dyco);
2473                                 MTC_Mat3MulVecfl(R.imat, O.dyco);
2474                         }
2475                 }
2476                 if((texco & TEXCO_UV) || (mode & (MA_VERTEXCOL|MA_FACETEXTURE)))  {
2477                         int j1=i1, j2=i2, j3=i3;
2478                         
2479                         /* to prevent storing new tfaces or vcols, we check a split runtime */
2480                         /*              4---3           4---3 */
2481                         /*              |\ 1|   or  |1 /| */
2482                         /*              |0\ |           |/ 0| */
2483                         /*              1---2           1---2   0 = orig face, 1 = new face */
2484                         
2485                         /* Update vert nums to point to correct verts of original face */
2486                         if(vlr->flag & R_DIVIDE_24) {  
2487                                 if(vlr->flag & R_FACE_SPLIT) {
2488                                         j1++; j2++; j3++;
2489                                 }
2490                                 else {
2491                                         j3++;
2492                                 }
2493                         }
2494                         else if(vlr->flag & R_FACE_SPLIT) {
2495                                 j2++; j3++; 
2496                         }
2497                         
2498                         if(mode & MA_VERTEXCOL) {
2499                                 
2500                                 if(vlr->vcol) {
2501                                         char *cp1, *cp2, *cp3;
2502                                         
2503                                         cp1= (char *)(vlr->vcol+j1);
2504                                         cp2= (char *)(vlr->vcol+j2);
2505                                         cp3= (char *)(vlr->vcol+j3);
2506
2507                                         shi->vcol[0]= (l*cp3[3]-u*cp1[3]-v*cp2[3])/255.0;
2508                                         shi->vcol[1]= (l*cp3[2]-u*cp1[2]-v*cp2[2])/255.0;
2509                                         shi->vcol[2]= (l*cp3[1]-u*cp1[1]-v*cp2[1])/255.0;
2510                                         
2511                                 }
2512                                 else {
2513                                         shi->vcol[0]= 0.0;
2514                                         shi->vcol[1]= 0.0;
2515                                         shi->vcol[2]= 0.0;
2516                                 }
2517                         }
2518                         if(vlr->tface) {
2519                                 float *uv1, *uv2, *uv3;
2520                                 
2521                                 uv1= vlr->tface->uv[j1];
2522                                 uv2= vlr->tface->uv[j2];
2523                                 uv3= vlr->tface->uv[j3];
2524                                 
2525                                 shi->uv[0]= -1.0 + 2.0*(l*uv3[0]-u*uv1[0]-v*uv2[0]);
2526                                 shi->uv[1]= -1.0 + 2.0*(l*uv3[1]-u*uv1[1]-v*uv2[1]);
2527                                 shi->uv[2]= 0.0;        // texture.c assumes there are 3 coords
2528                                 
2529                                 if(shi->osatex) {
2530                                         float duv[2];
2531                                         
2532                                         dl= O.dxuv[0]+O.dxuv[1];
2533                                         duv[0]= O.dxuv[0]; 
2534                                         duv[1]= O.dxuv[1];
2535                                         
2536                                         O.dxuv[0]= 2.0*(dl*uv3[0]-duv[0]*uv1[0]-duv[1]*uv2[0]);
2537                                         O.dxuv[1]= 2.0*(dl*uv3[1]-duv[0]*uv1[1]-duv[1]*uv2[1]);
2538         
2539                                         dl= O.dyuv[0]+O.dyuv[1];
2540                                         duv[0]= O.dyuv[0]; 
2541                                         duv[1]= O.dyuv[1];
2542         
2543                                         O.dyuv[0]= 2.0*(dl*uv3[0]-duv[0]*uv1[0]-duv[1]*uv2[0]);
2544                                         O.dyuv[1]= 2.0*(dl*uv3[1]-duv[0]*uv1[1]-duv[1]*uv2[1]);
2545                                 }
2546                                 
2547                                 if(mode & MA_FACETEXTURE) {
2548                                         if((mode & MA_VERTEXCOL)==0) {
2549                                                 shi->vcol[0]= 1.0;
2550                                                 shi->vcol[1]= 1.0;
2551                                                 shi->vcol[2]= 1.0;
2552                                         }
2553                                         if(vlr->tface) render_realtime_texture(shi);
2554                                 }
2555                         }
2556                         else {
2557                                 shi->uv[0]= 2.0*(u+.5);
2558                                 shi->uv[1]= 2.0*(v+.5);
2559                                 shi->uv[2]= 0.0;        // texture.c assumes there are 3 coords
2560                         }
2561                 }
2562                 if(texco & TEXCO_NORM) {
2563                         shi->orn[0]= shi->vn[0];
2564                         shi->orn[1]= shi->vn[1];
2565                         shi->orn[2]= shi->vn[2];
2566                 }
2567                                 
2568                 if(mode & MA_RADIO) {
2569                         shi->rad[0]= (l*v3->rad[0] - u*v1->rad[0] - v*v2->rad[0]);
2570                         shi->rad[1]= (l*v3->rad[1] - u*v1->rad[1] - v*v2->rad[1]);
2571                         shi->rad[2]= (l*v3->rad[2] - u*v1->rad[2] - v*v2->rad[2]);
2572                 }
2573                 else {
2574                         shi->rad[0]= shi->rad[1]= shi->rad[2]= 0.0;
2575                 }
2576                         
2577                 if(texco & TEXCO_REFL) {
2578                         /* mirror reflection colour textures (and envmap) */
2579                         calc_R_ref(shi);
2580                 }
2581                 
2582         }
2583         else {
2584                 shi->rad[0]= shi->rad[1]= shi->rad[2]= 0.0;
2585         }
2586 }
2587
2588   /* x,y: window coordinate from 0 to rectx,y */
2589   /* return pointer to rendered face */
2590   
2591 float bluroffsx, bluroffsy;     // set in initrender.c (ton)
2592
2593 void *shadepixel(float x, float y, int vlaknr, int mask, float *col)
2594 {
2595         ShadeResult shr;
2596         ShadeInput shi;
2597         VlakRen *vlr=NULL;
2598         
2599         if(vlaknr< 0) { /* error */
2600                 return NULL;
2601         }
2602         /* currently in use for dithering soft shadow */
2603         shi.xs= x;
2604         shi.ys= y;
2605         
2606         /* mask is used to indicate amount of samples (ray shad/mir and AO) */
2607         shi.mask= mask;
2608         
2609         if(vlaknr==0) { /* sky */
2610                 col[0]= 0.0; col[1]= 0.0; col[2]= 0.0; col[3]= 0.0;
2611         }
2612         else if( (vlaknr & 0x7FFFFF) <= R.totvlak) {
2613                 VertRen *v1, *v2, *v3;
2614                 float alpha, fac, dvlak, deler;
2615                 
2616                 vlr= RE_findOrAddVlak( (vlaknr-1) & 0x7FFFFF);
2617                 
2618                 shi.mat= vlr->mat;
2619                 shi.matren= shi.mat->ren;
2620                 shi.vlr= vlr;
2621                 shi.osatex= (shi.matren->texco & TEXCO_OSA);
2622
2623                 v1= vlr->v1;
2624                 dvlak= v1->co[0]*vlr->n[0]+v1->co[1]*vlr->n[1]+v1->co[2]*vlr->n[2];
2625
2626                 /* COXYZ AND VIEW VECTOR  */
2627                 shi.view[0]= (x+(R.xstart)+bluroffsx +0.5);
2628
2629                 if(R.flag & R_SEC_FIELD) {
2630                         if(R.r.mode & R_ODDFIELD) shi.view[1]= (y+R.ystart)*R.ycor;
2631                         else shi.view[1]= (y+R.ystart+1.0)*R.ycor;
2632                 }
2633                 else shi.view[1]= (y+R.ystart+bluroffsy+0.5)*R.ycor;
2634                 
2635                 shi.view[2]= -R.viewfac;
2636
2637                 if(R.r.mode & R_PANORAMA) {
2638                         float panoco, panosi, u, v;
2639                         panoco = getPanovCo();
2640                         panosi = getPanovSi();
2641
2642                         u= shi.view[0]; v= shi.view[2];
2643                         shi.view[0]= panoco*u + panosi*v;
2644                         shi.view[2]= -panosi*u + panoco*v;
2645                 }
2646
2647                 deler= vlr->n[0]*shi.view[0] + vlr->n[1]*shi.view[1] + vlr->n[2]*shi.view[2];
2648                 if (deler!=0.0) fac= R.zcor= dvlak/deler;
2649                 else fac= R.zcor= 0.0;
2650                 
2651                 shi.co[0]= fac*shi.view[0];
2652                 shi.co[1]= fac*shi.view[1];
2653                 shi.co[2]= fac*shi.view[2];
2654                 
2655                 /* pixel dx/dy for render coord */
2656                 if(shi.osatex || (R.r.mode & R_SHADOW) ) {
2657                         float u= dvlak/(deler-vlr->n[0]);
2658                         float v= dvlak/(deler- R.ycor*vlr->n[1]);
2659
2660                         O.dxco[0]= shi.co[0]- (shi.view[0]-1.0)*u;
2661                         O.dxco[1]= shi.co[1]- (shi.view[1])*u;
2662                         O.dxco[2]= shi.co[2]- (shi.view[2])*u;
2663
2664                         O.dyco[0]= shi.co[0]- (shi.view[0])*v;
2665                         O.dyco[1]= shi.co[1]- (shi.view[1]-1.0*R.ycor)*v;
2666                         O.dyco[2]= shi.co[2]- (shi.view[2])*v;
2667
2668                 }
2669
2670                 fac= Normalise(shi.view);
2671                 R.zcor*= fac;   /* for mist */
2672                 
2673                 if(shi.osatex) {
2674                         if( (shi.matren->texco & TEXCO_REFL) ) {
2675                                 O.dxview= -1.0/fac;
2676                                 O.dyview= -R.ycor/fac;
2677                         }
2678                 }
2679                 
2680                 /* calcuate normals, texture coords, vertex colors, etc */
2681                 if(vlaknr & 0x800000)
2682                         shade_input_set_coords(&shi, 1.0, 1.0, 0, 2, 3);
2683                 else 
2684                         shade_input_set_coords(&shi, 1.0, 1.0, 0, 1, 2);
2685
2686                 /* this only avalailable for scanline */
2687                 if(shi.matren->texco & TEXCO_WINDOW) {
2688                         shi.winco[0]= (x+(R.xstart))/(float)R.afmx;
2689                         shi.winco[1]= (y+(R.ystart))/(float)R.afmy;
2690                 }
2691                 /* after this the u and v AND O.dxuv and O.dyuv are incorrect */
2692                 if(shi.matren->texco & TEXCO_STICKY) {
2693                         if(v1->sticky) {
2694                                 extern float Zmulx, Zmuly;
2695                                 float *o1, *o2, *o3, hox, hoy, l, dl, u, v;
2696                                 float s00, s01, s10, s11, detsh;
2697                                 
2698                                 if(vlaknr & 0x800000) {
2699                                         v2= vlr->v3; v3= vlr->v4;
2700                                 } else {
2701                                         v2= vlr->v2; v3= vlr->v3;
2702                                 }
2703                                 
2704                                 s00= v3->ho[0]/v3->ho[3] - v1->ho[0]/v1->ho[3];
2705                                 s01= v3->ho[1]/v3->ho[3] - v1->ho[1]/v1->ho[3];
2706                                 s10= v3->ho[0]/v3->ho[3] - v2->ho[0]/v2->ho[3];
2707                                 s11= v3->ho[1]/v3->ho[3] - v2->ho[1]/v2->ho[3];
2708                                 
2709                                 detsh= s00*s11-s10*s01;
2710                                 s00/= detsh; s01/=detsh; 
2711                                 s10/=detsh; s11/=detsh;
2712         
2713                                 /* recalc u and v again */
2714                                 hox= x/Zmulx -1.0;
2715                                 hoy= y/Zmuly -1.0;
2716                                 u= (hox - v3->ho[0]/v3->ho[3])*s11 - (hoy - v3->ho[1]/v3->ho[3])*s10;
2717                                 v= (hoy - v3->ho[1]/v3->ho[3])*s00 - (hox - v3->ho[0]/v3->ho[3])*s01;
2718                                 l= 1.0+u+v;
2719                                 
2720                                 o1= v1->sticky;
2721                                 o2= v2->sticky;
2722                                 o3= v3->sticky;
2723                                 
2724                                 shi.sticky[0]= l*o3[0]-u*o1[0]-v*o2[0];
2725                                 shi.sticky[1]= l*o3[1]-u*o1[1]-v*o2[1];
2726         
2727                                 if(shi.osatex) {
2728                                         O.dxuv[0]=  s11/Zmulx;
2729                                         O.dxuv[1]=  - s01/Zmulx;
2730                                         O.dyuv[0]=  - s10/Zmuly;
2731                                         O.dyuv[1]=  s00/Zmuly;
2732                                         
2733                                         dl= O.dxuv[0]+O.dxuv[1];
2734                                         O.dxsticky[0]= dl*o3[0]-O.dxuv[0]*o1[0]-O.dxuv[1]*o2[0];
2735                                         O.dxsticky[1]= dl*o3[1]-O.dxuv[0]*o1[1]-O.dxuv[1]*o2[1];
2736                                         dl= O.dyuv[0]+O.dyuv[1];
2737                                         O.dysticky[0]= dl*o3[0]-O.dyuv[0]*o1[0]-O.dyuv[1]*o2[0];
2738                                         O.dysticky[1]= dl*o3[1]-O.dyuv[0]*o1[1]-O.dyuv[1]*o2[1];
2739                                 }
2740                         }
2741                 }
2742                 
2743                 /* ------  main shading loop */
2744                 shade_lamp_loop(&shi, &shr);
2745
2746                 if(shi.matren->translucency!=0.0) {
2747                         ShadeResult shr_t;
2748                         
2749                         VecMulf(shi.vn, -1.0);
2750                         VecMulf(shi.vlr->n, -1.0);
2751                         shade_lamp_loop(&shi, &shr_t);
2752                         shr.diff[0]+= shi.matren->translucency*shr_t.diff[0];
2753                         shr.diff[1]+= shi.matren->translucency*shr_t.diff[1];
2754                         shr.diff[2]+= shi.matren->translucency*shr_t.diff[2];
2755                         VecMulf(shi.vn, -1.0);
2756                         VecMulf(shi.vlr->n, -1.0);
2757                 }
2758                 
2759                 if(R.r.mode & R_RAYTRACE) {
2760                         if(shi.matren->ray_mirror!=0.0 || (shi.mat->mode & MA_RAYTRANSP && shr.alpha!=1.0)) {
2761                                 ray_trace(&shi, &shr);
2762                         }
2763                 }
2764                 else {
2765                         // doesnt look 'correct', but is better for preview, plus envmaps dont raytrace this
2766                         if(shi.mat->mode & MA_RAYTRANSP) shr.alpha= 1.0;
2767                 }
2768                 
2769                 VecAddf(col, shr.diff, shr.spec);
2770                 
2771                 /* exposure correction */
2772                 if(R.wrld.exp!=0.0 || R.wrld.range!=1.0) {
2773                         if((shi.matren->mode & MA_SHLESS)==0) {
2774                                 col[0]= R.wrld.linfac*(1.0-exp( col[0]*R.wrld.logfac) );
2775                                 col[1]= R.wrld.linfac*(1.0-exp( col[1]*R.wrld.logfac) );
2776                                 col[2]= R.wrld.linfac*(1.0-exp( col[2]*R.wrld.logfac) );
2777                         }
2778                 }
2779                 
2780                 /* MIST */
2781                 if( (R.wrld.mode & WO_MIST) && (shi.matren->mode & MA_NOMIST)==0 ){
2782                         alpha= mistfactor(shi.co);
2783                 }
2784                 else alpha= 1.0;
2785
2786                 if(shr.alpha!=1.0 || alpha!=1.0) {
2787                         fac= alpha*(shr.alpha);
2788                         
2789                         col[3]= fac;
2790                         col[0]*= fac;
2791                         col[1]*= fac;
2792                         col[2]*= fac;
2793                 }
2794                 else col[3]= 1.0;
2795         }
2796         
2797         if(R.flag & R_LAMPHALO) {
2798                 if(vlaknr<=0) { /* calc view vector and put shi.co at far */
2799                 
2800                         shi.view[0]= (x+(R.xstart)+0.5);
2801
2802                         if(R.flag & R_SEC_FIELD) {
2803                                 if(R.r.mode & R_ODDFIELD) shi.view[1]= (y+R.ystart)*R.ycor;
2804                                 else shi.view[1]= (y+R.ystart+1.0)*R.ycor;
2805                         }
2806                         else shi.view[1]= (y+R.ystart+0.5)*R.ycor;
2807                         
2808                         shi.view[2]= -R.viewfac;
2809                         
2810                         if(R.r.mode & R_PANORAMA) {
2811                                 float u,v, panoco, panosi;
2812                                 panoco = getPanovCo();
2813                                 panosi = getPanovSi();
2814                                 
2815                                 u= shi.view[0]; v= shi.view[2];
2816                                 shi.view[0]= panoco*u + panosi*v;
2817                                 shi.view[2]= -panosi*u + panoco*v;
2818                         }
2819
2820                         shi.co[2]= 0.0;
2821                         
2822                 }
2823                 renderspothalo(&shi, col);
2824         }
2825         
2826         return vlr;
2827 }
2828
2829 void shadepixel_short(float x, float y, int vlaknr, int mask, unsigned short *shortcol)
2830 {
2831         float colf[4];
2832         
2833         shadepixel(x, y, vlaknr, mask, colf);
2834         
2835         if(colf[0]<=0.0) shortcol[0]= 0; else if(colf[0]>=1.0) shortcol[0]= 65535;
2836         else shortcol[0]= 65535.0*colf[0];
2837         if(colf[1]<=0.0) shortcol[1]= 0; else if(colf[1]>=1.0) shortcol[1]= 65535;
2838         else shortcol[1]= 65535.0*colf[1];
2839         if(colf[2]<=0.0) shortcol[2]= 0; else if(colf[2]>=1.0) shortcol[2]= 65535;
2840         else shortcol[2]= 65535.0*colf[2];
2841         if(colf[3]<=0.0) shortcol[3]= 0; else if(colf[3]>=1.0) shortcol[3]= 65535;
2842         else shortcol[3]= 65535.0*colf[3];
2843
2844         if(usegamtab) {
2845                 shortcol[0]= igamtab2[ shortcol[0] ];
2846                 shortcol[1]= igamtab2[ shortcol[1] ];
2847                 shortcol[2]= igamtab2[ shortcol[2] ];
2848         }
2849 }
2850
2851 PixStr *addpsmain()
2852 {
2853         PixStrMain *psm;
2854
2855         psm= &psmfirst;
2856
2857         while(psm->next) {
2858                 psm= psm->next;
2859         }
2860
2861         psm->next= (PixStrMain *)MEM_mallocN(sizeof(PixStrMain),"pixstrMain");
2862
2863         psm= psm->next;
2864         psm->next=0;
2865         psm->ps= (PixStr *)MEM_mallocN(4096*sizeof(PixStr),"pixstr");
2866         psmteller= 0;
2867