ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/pictool.c
Revision: 2.6
Committed: Fri Jul 19 17:37:56 2019 UTC (4 years, 9 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.5: +2 -2 lines
Log Message:
Moved declarations and definitions for header.c from resolu.h to rtio.h

File Contents

# User Rev Content
1 greg 2.2 #ifndef lint
2 greg 2.6 static const char RCSid[] = "$Id: pictool.c,v 2.5 2018/01/24 04:39:52 greg Exp $";
3 greg 2.2 #endif
4 greg 2.1 #include "pictool.h"
5     #include "g3sphere.h"
6     #include <math.h>
7 greg 2.6 #include "rtio.h"
8 greg 2.1
9     #define index strchr
10    
11    
12    
13     int pict_init(pict* p)
14     {
15     p->pic = NULL;
16     p->use_lut = 0;
17     p->lut = NULL;
18     p->pjitter = 0;
19     p->comment[0] = '\0';
20     #ifdef PICT_GLARE
21     p->pinfo = NULL;
22     p->glareinfo = g3fl_create(PICT_GLSIZE);
23     #endif
24     p->valid_view = 1;
25     p->view = stdview;
26     if (!pict_update_view(p))
27     return 0;
28     p->resol.rt = PIXSTANDARD;
29     return pict_resize(p,512,512);
30     }
31    
32     void pict_new_gli(pict* p)
33     {
34     int i;
35    
36     g3Float* gli;
37     gli = g3fl_append_new(p->glareinfo);
38     for(i=0;i<PICT_GLSIZE;i++)
39     gli[i] = 0.0;
40     }
41    
42     void pict_set_comment(pict* p,const char* c)
43     {
44     strcpy(p->comment,c);
45     }
46    
47     pict* pict_get_copy(pict* p)
48     {
49     pict* res;
50    
51     if (p->use_lut) {
52     fprintf(stderr,"pict_get_copy: can't copy in lut mode");
53     return NULL;
54     }
55     res = pict_create();
56     pict_resize(res,pict_get_xsize(p),pict_get_ysize(p));
57     res->view = p->view;
58     res->resol = p->resol;
59     res->pjitter = p->pjitter;
60     return res;
61     }
62    
63     pict* pict_create()
64     {
65     pict* p = (pict*) malloc(sizeof(pict));
66    
67     if (!pict_init(p)) {
68     pict_free(p);
69     return NULL;
70     }
71     return p;
72     }
73    
74     void pict_free(pict* p)
75     {
76     if (p->pic)
77     free(p->pic);
78     #ifdef PICT_GLARE
79     if (p->pinfo)
80     free(p->pinfo);
81     g3fl_free(p->glareinfo);
82     p->pinfo = NULL;
83     #endif
84     p->pic = NULL;
85     p->lut = NULL;
86     free(p);
87     }
88     static int exp_headline(char* s,void* exparg) {
89     double* exposure;
90 greg 2.3 if(strstr(s, EXPOSSTR) != NULL ) {
91     fprintf(stderr,"EXP AND tab\n");
92     }
93 greg 2.1
94     exposure = (double*) exparg;
95     if (isexpos(s)) {
96     (*exposure) *= exposval(s);
97     }
98     return 0;
99     }
100    
101     static double read_exposure(FILE* fp)
102     {
103     double exposure = 1;
104    
105    
106     if ((getheader(fp, exp_headline, &exposure)) < 0) {
107     fprintf(stderr,"error reading header\n");
108     return 0;
109     }
110     return exposure;
111     }
112    
113     int pict_update_view(pict* p)
114     {
115     char* ret;
116     if ((ret = setview(&(p->view)))) {
117     fprintf(stderr,"pict_update_view: %s\n",ret);
118     return 0;
119     }
120     return 1;
121     }
122    
123     int pict_set_pj(pict* p,double pj)
124     {
125     if ((pj < 0) || (pj > 1))
126     fprintf(stderr,"pict_set_pj: warning jitter factor out of range\n");
127     p->pjitter = pj;
128     return 1;
129     }
130    
131     int pict_setvup(pict* p,FVECT vup)
132     {
133     VCOPY(p->view.vup,vup);
134     return pict_update_view(p);
135     }
136    
137     int pict_trans_cam(pict* p,g3ATrans t)
138     {
139     g3Vec o = g3v_create();
140     g3Vec d = g3v_create();
141     g3Vec u = g3v_create();
142    
143     fprintf(stderr,"vdpre %f %f %f\n",p->view.vdir[0],p->view.vdir[1],p->view.vdir[2]);
144     fprintf(stderr,"vupre %f %f %f\n",p->view.vup[0],p->view.vup[1],p->view.vup[2]);
145     g3v_fromrad(o,p->view.vp);
146     g3v_fromrad(d,p->view.vdir);
147     g3v_fromrad(u,p->view.vup);
148    
149     g3v_add(d,d,o);
150     g3v_add(u,u,o);
151     g3at_apply(t,o);
152     g3at_apply(t,u);
153     g3at_apply(t,d);
154     g3v_sub(d,d,o);
155     g3v_sub(u,u,o);
156     g3v_torad(p->view.vp,o);
157     g3v_torad(p->view.vdir,d);
158     g3v_torad(p->view.vup,u);
159     fprintf(stderr,"vd %f %f %f\n",p->view.vdir[0],p->view.vdir[1],p->view.vdir[2]);
160     fprintf(stderr,"vu %f %f %f\n",p->view.vup[0],p->view.vup[1],p->view.vup[2]);
161    
162     g3v_free(o);
163     g3v_free(d);
164     g3v_free(u);
165     return pict_update_view(p);
166     }
167    
168     int pict_resize(pict* p,int x,int y)
169     {
170     p->resol.xr = x;
171     p->resol.yr = y;
172    
173     if (!(p->pic = (COLOR*) realloc(p->pic,(p->resol.xr)*(p->resol.yr)*sizeof(COLOR)))) {
174     fprintf(stderr,"Can't alloc picture\n");
175     return 0;
176     }
177     if (pict_lut_used(p)) {
178     if (!(p->lut = (COLOR*) realloc(p->lut,(p->resol.xr)*(p->resol.yr)*sizeof(COLOR)))) {
179     fprintf(stderr,"Can't alloc LUT\n");
180     return 0;
181     }
182     }
183     #ifdef PICT_GLARE
184     if (!(p->pinfo = (pixinfo*)
185     realloc(p->pinfo,(p->resol.xr)*p->resol.yr*sizeof(pixinfo)))) {
186     fprintf(stderr,"Can't alloc picture info\n");
187     return 0;
188     }
189     /* fprintf(stderr,"resize %d %d %d\n",p->resol.xr,p->resol.yr,sizeof(pixinfo));*/
190     #endif
191     return 1;
192     }
193    
194    
195     int pict_zero(pict* p)
196     {
197     int x,y,i;
198    
199     for(x=0;x<pict_get_xsize(p);x++)
200     for(y=0;y<pict_get_ysize(p);y++)
201     for(i=0;i<3;i++)
202     pict_get_color(p,x,y)[i] = 0;
203     return 1;
204     }
205    
206     struct hinfo {
207     VIEW *hv;
208     int ok;
209     double exposure;
210     };
211    
212    
213     static int
214     gethinfo( /* get view from header */
215     char *s,
216     void *v
217     )
218     {
219 greg 2.3 if(strstr(s, EXPOSSTR) != NULL && strstr(s, "\t") != NULL) {
220 greg 2.4
221     fprintf(stderr,"error: header contains invalid exposure entry!!!!\n");
222     fprintf(stderr,"check exposure and correct header setting !\n");
223    
224    
225 greg 2.3 fprintf(stderr,"stopping !!!!\n");
226     exit(1);
227     }
228 greg 2.1 if (isview(s) && sscanview(((struct hinfo*)v)->hv, s) > 0) {
229     ((struct hinfo*)v)->ok++;
230     } else if (isexpos(s)) {
231     /*fprintf(stderr,"readexp %f %f\n",((struct hinfo*)v)->exposure,exposval(s));*/
232     ((struct hinfo*)v)->exposure *= exposval(s);
233     }
234     return(0);
235     }
236    
237     int pict_read_fp(pict* p,FILE* fp)
238     {
239     struct hinfo hi;
240     int x,y,yy;
241    
242    
243     hi.hv = &(p->view);
244     hi.ok = 0;
245     hi.exposure = 1;
246     getheader(fp, gethinfo, &hi);
247     /*fprintf(stderr,"expscale %f\n",hi.exposure);*/
248    
249     if (!(pict_update_view(p)) || !(hi.ok))
250     p->valid_view = 0;
251     /* printf("dir %f %f %f\n",p->view.vdir[0],p->view.vdir[1],p->view.vdir[2]); */
252     if (fgetsresolu(&(p->resol),fp) < 0) {
253     fprintf(stderr,"No resolution given\n");
254     return 0;
255     }
256     if (!(p->resol.rt & YMAJOR)) {
257     x = p->resol.xr;
258     p->resol.xr = p->resol.yr;
259     p->resol.yr = x;
260     p->resol.rt |= YMAJOR;
261     }
262     if (!pict_resize(p,p->resol.xr,p->resol.yr))
263     return 0;
264     for(yy=0;yy<p->resol.yr;yy++) {
265     y = (p->resol.rt & YDECR) ? p->resol.yr - 1 - yy : yy;
266     if (freadscan((p->pic + y*p->resol.xr),p->resol.xr,fp) < 0) {
267     fprintf(stderr,"error reading line %d\n",y);
268     return 0;
269     }
270     for(x=0;x<p->resol.xr;x++) {
271     scalecolor(pict_get_color(p,x,y),1.0/hi.exposure);
272     }
273     }
274     #ifdef PICT_GLARE
275     pict_update_evalglare_caches(p);
276     #endif
277     return 1;
278     }
279    
280     #ifdef PICT_GLARE
281     void pict_update_evalglare_caches(pict* p)
282     {
283     int x,y,yy;
284     float lum;
285     for(yy=0;yy<p->resol.yr;yy++) {
286     y = (p->resol.rt & YDECR) ? p->resol.yr - 1 - yy : yy;
287     for(x=0;x<p->resol.xr;x++) {
288     pict_get_omega(p,x,y) = pict_get_sangle(p,x,y);
289     pict_get_dir(p,x,y,pict_get_cached_dir(p,x,y));
290     pict_get_gsn(p,x,y) = 0;
291     pict_get_pgs(p,x,y) = 0;
292     /* make picture grey */
293     lum=luminance(pict_get_color(p,x,y))/179.0;
294     pict_get_color(p,x,y)[RED]=lum;
295     pict_get_color(p,x,y)[GRN]=lum;
296     pict_get_color(p,x,y)[BLU]=lum;
297     }
298     }
299     }
300    
301     #endif
302    
303     int pict_read(pict* p,char* fn)
304     {
305     FILE* fp;
306    
307     if (!(fp = fopen(fn,"rb"))) {
308     fprintf(stderr,"Can't open %s\n",fn);
309     return 0;
310     }
311     if (!(pict_read_fp(p,fp)))
312     return 0;
313     fclose(fp);
314     return 1;
315     }
316    
317     static void ch_pct(char* line)
318     {
319     char* pct;
320    
321     while ((pct = index(line,',')))
322     *pct = '.';
323     }
324    
325     char* nextline(FILE* fp)
326     {
327     static char* line = NULL;
328     static int lsize = 500;
329     int done;
330     long fpos;
331    
332     if (!line) {
333     if (!(line = (char*) malloc((lsize + 1)*sizeof(char)))) {
334     fprintf(stderr,"allocation problem in nextline\n");
335     return NULL;
336     }
337     }
338     done = 0;
339     fpos = ftell(fp);
340     while(!done) {
341     if (!fgets(line,lsize,fp))
342     return NULL;
343     if (!index(line,'\n')) {
344     lsize *= 2;
345     if (!(line = (char*) realloc(line,(lsize + 1)*sizeof(char)))) {
346     fprintf(stderr,"allocation problem in nextline\n");
347     return NULL;
348     }
349     fseek(fp,fpos,SEEK_SET);
350     } else {
351     done = 1;
352     }
353     }
354     return line;
355     }
356    
357     static int read_tt_txt(pict* p,int bufsel,int channel,char* fn)
358     {
359     FILE* fp;
360     char* line = NULL;
361     int x,y,vsize,hsize;
362     float* col;
363     float c;
364     char* curr;
365     char** endp;
366     COLOR* picbuf;
367    
368     if (!(fp = fopen(fn,"r"))) {
369     fprintf(stderr,"pict_read_tt_txt: Can't open file %s\n",fn);
370     return 0;
371     }
372     if (!(line = nextline(fp))) {
373     fprintf(stderr,"pict_read_tt_txt: reading from file failed\n");
374     return 0;
375     }
376     if (strncmp(line,"float",strlen("float"))) {
377     fprintf(stderr,"pict_read_tt_txt: Unexpected format\n");
378     return 0;
379     }
380     if (!(line = nextline(fp))) {
381     fprintf(stderr,"pict_read_tt_txt: reading from file failed\n");
382     return 0;
383     }
384     sscanf(line,"%d %d %d %d",&x,&vsize,&y,&hsize);
385     vsize = vsize - x;
386     hsize = hsize - y;
387     if (!pict_resize(p,vsize,hsize))
388     return 0;
389     endp = malloc(sizeof(char*));
390     picbuf = (bufsel == PICT_LUT) ? p->lut : p->pic;
391     for(x=0;x<vsize;x++) {
392     if (!(line = nextline(fp))) {
393     fprintf(stderr,"pict_read_tt_txt: too few lines %d of %d\n",x,vsize);
394     free(endp);
395     return 0;
396     }
397     ch_pct(line);
398     curr = line;
399     for(y=0;y<hsize;y++) {
400     col = pict_colbuf(p,picbuf,x,y);
401     c = strtod(curr,endp);
402     if (channel & PICT_CH1)
403     col[0] = c;
404     if (channel & PICT_CH2)
405     col[1] = c;
406     if (channel & PICT_CH3)
407     col[2] = c;
408     if (curr == (*endp)) {
409     fprintf(stderr,"pict_read_tt_txt: too few values in line %d\n",x);
410     free(endp);
411     return 0;
412     }
413     curr = (*endp);
414     }
415     }
416     free(endp);
417     return 1;
418     }
419    
420     int pict_read_tt_txt(pict* p,char* fn)
421     {
422     return read_tt_txt(p,PICT_VIEW,PICT_ALLCH,fn);
423     }
424    
425     int pict_setup_lut(pict* p,char* theta_fn,char* phi_fn,char* omega_fn)
426     {
427     int xs,ys;
428    
429     p->use_lut = 1;
430     if (!read_tt_txt(p,PICT_LUT,PICT_CH1,theta_fn)) {
431     fprintf(stderr,"pict_setup_lut: error reading theta file\n");
432     p->use_lut = 0;
433     return 0;
434     }
435     xs = pict_get_xsize(p);
436     ys = pict_get_ysize(p);
437     if (!read_tt_txt(p,PICT_LUT,PICT_CH2,phi_fn)) {
438     fprintf(stderr,"pict_setup_lut: error reading phi file\n");
439     p->use_lut = 0;
440     return 0;
441     }
442     if ((xs != pict_get_xsize(p)) || (ys != pict_get_ysize(p))) {
443     fprintf(stderr,"pict_setup_lut: different resolutions for lut files\n");
444     p->use_lut = 0;
445     return 0;
446     }
447     if (!read_tt_txt(p,PICT_LUT,PICT_CH3,omega_fn)) {
448     fprintf(stderr,"pict_setup_lut: error reading omega file\n");
449     p->use_lut = 0;
450     return 0;
451     }
452     if ((xs != pict_get_xsize(p)) || (ys != pict_get_ysize(p))) {
453     fprintf(stderr,"pict_setup_lut: different resolutions for lut files\n");
454     p->use_lut = 0;
455     return 0;
456     }
457     return 1;
458     }
459    
460     int pict_update_lut(pict* p)
461     {
462     g3ATrans transf;
463     g3Vec d;
464     g3Vec u;
465     int x,y;
466    
467     d = g3v_fromrad(g3v_create(),p->view.vdir);
468     u = g3v_fromrad(g3v_create(),p->view.vup);
469     g3at_ctrans(transf,d,u);
470     for(x=0;x<p->resol.xr;x++) {
471     for(y=0;y<p->resol.yr;y++) {
472     d[G3S_RAD] = 1.0;
473     d[G3S_THETA] = DEG2RAD(pict_colbuf(p,p->lut,x,y)[0]);
474     d[G3S_PHI] = DEG2RAD(pict_colbuf(p,p->lut,x,y)[1]);
475     g3s_sphtocc(d,d);
476     g3at_apply(transf,d);
477     g3s_cctosph(d,d);
478     pict_colbuf(p,p->lut,x,y)[0] = RAD2DEG(d[G3S_THETA]);
479     pict_colbuf(p,p->lut,x,y)[1] = RAD2DEG(d[G3S_PHI]);
480     }
481     }
482     g3v_free(d);
483     g3v_free(u);
484     return 1;
485     }
486    
487    
488     int pict_write_dir(pict* p,FILE* fp)
489     {
490     int x,y;
491     FVECT dir,orig;
492    
493     for(y=0;y<p->resol.yr;y++) {
494     for(x=0;x<p->resol.xr;x++) {
495     if (!(pict_get_ray(p,x,y,orig,dir))) {
496     fprintf(fp,"%f %f %f 0.0 0.0 0.0\n",orig[0],orig[1],orig[2]);
497     } else {
498     fprintf(fp,"%f %f %f %f %f %f\n",orig[0],orig[1],orig[2],dir[0],dir[1],dir[2]);
499     }
500     }
501     }
502     return 1;
503     }
504    
505     int pict_read_from_list(pict* p,FILE* fp)
506     {
507     int x,y,sz;
508     char ln[1024];
509     double val;
510    
511     sz = 0;
512     for(y=0;y<p->resol.yr;y++) {
513     for(x=0;x<p->resol.xr;x++) {
514     if (!(fgets(ln,1023,fp))) {
515     fprintf(stderr,"pict_read_from_list: Read only %d values\n",sz);
516     return 0;
517     }
518     if (sscanf(ln,"%lg",&val) != 1) {
519     fprintf(stderr,"pict_read_from_list: errline %d\n",sz);
520     } else {
521     sz++;
522     setcolor(pict_get_color(p,x,y),val,val,val);
523     }
524     }
525     }
526     return 1;
527     }
528    
529     int pict_write(pict* p,char* fn)
530     {
531     FILE* fp;
532    
533     if (!(fp = fopen(fn,"wb"))) {
534     fprintf(stderr,"Can't open %s for writing\n",fn);
535     return 0;
536     }
537     if (!(pict_write_fp(p,fp))) {
538     fclose(fp);
539     return 0;
540     }
541     fclose(fp);
542     return 1;
543     }
544    
545     int pict_write_fp(pict* p,FILE* fp)
546     {
547     char res[100];
548     int y,yy;
549     newheader("RADIANCE",fp);
550     fputnow(fp);
551     if (strlen(p->comment) > 0)
552     fputs(p->comment,fp);
553     fputs(VIEWSTR, fp);
554     fprintview(&(p->view), fp);
555     fprintf(fp,"\n");
556     fputformat(COLRFMT,fp);
557     fprintf(fp,"\n");
558    
559     resolu2str(res,&(p->resol));
560     fprintf(fp,"%s",res);
561     for(y=0;y<p->resol.yr;y++) {
562     yy = (p->resol.rt & YDECR) ? p->resol.yr - 1 - y : y;
563     if (fwritescan((p->pic + yy*p->resol.xr),p->resol.xr,fp) < 0) {
564     fprintf(stderr,"writing pic-file failed\n");
565     return 0;
566     }
567     }
568     return 1;
569     }
570    
571     int pict_get_dir(pict* p,int x,int y,FVECT dir)
572     {
573     FVECT orig;
574     if (pict_lut_used(p)) {
575     if (!(p->lut)) {
576     fprintf(stderr,"pict_get_dir: no lut defined\n");
577     return 0;
578     }
579     dir[G3S_RAD] = 1.0;
580     dir[G3S_THETA] = DEG2RAD(pict_colbuf(p,p->lut,x,y)[0]);
581     dir[G3S_PHI] = DEG2RAD(pict_colbuf(p,p->lut,x,y)[1]);
582     g3s_sphtocc(dir,dir);
583     return 1;
584     }
585     return pict_get_ray(p,x,y,orig,dir);
586     }
587    
588     int pict_get_dir_ic(pict* p,double x,double y,FVECT dir)
589     {
590     FVECT orig;
591     if (pict_lut_used(p)) {
592     int res[2];
593     loc2pix(res,&(p->resol),x,y);
594     return pict_get_dir(p,res[0],res[1],dir);
595     }
596     return pict_get_ray_ic(p,x,y,orig,dir);
597     }
598    
599     int pict_get_ray(pict* p,int x,int y,FVECT orig,FVECT dir)
600     {
601     RREAL pc[2];
602    
603     pix2loc(pc,&(p->resol),x,p->resol.yr - 1 - y);
604     if (pict_lut_used(p)) {
605     if (viewray(orig,dir,&(p->view),pc[0],pc[1]) <= 0)
606     return 0;
607     return pict_get_dir(p,x,y,dir);
608     }
609     return pict_get_ray_ic(p,pc[0],pc[1],orig,dir);
610     }
611    
612     int pict_get_ray_ic(pict* p,double x,double y,FVECT orig,FVECT dir)
613     {
614     if (pict_lut_used(p)) {
615     int res[2];
616     loc2pix(res,&(p->resol),x,y);
617     return pict_get_ray(p,res[0],res[1],orig,dir);
618     }
619     x = (x + (0.5*p->pjitter)*(0.5-gb_random())/p->resol.xr);
620     y = (y + (0.5*p->pjitter)*(0.5-gb_random())/p->resol.yr);
621     if (viewray(orig,dir,&(p->view),x,y) < 0)
622     return 0;
623     return 1;
624     }
625    
626    
627     static int splane_normal(pict* p,double e1x,double e1y,
628     double e2x,double e2y,FVECT n)
629     {
630     FVECT e1 = {0,0,0};
631     FVECT e2 = {0,0,0};
632    
633     pict_get_dir_ic(p,e1x,e1y,e1);
634     pict_get_dir_ic(p,e2x,e2y,n);
635     VSUB(e2,n,e1);
636     VCROSS(n,e1,e2);
637     if (normalize(n) == 0.0) {
638     return 0;
639     }
640     return 1;
641     }
642    
643     double pict_get_sangle(pict* p,int x,int y)
644     {
645     int i;
646     double ang,a,hpx,hpy;
647     RREAL pc[2];
648     FVECT n[4] = {{0,0,0},{0,0,0},{0,0,0},{0,0,0}};
649    
650     if (pict_lut_used(p)) {
651     return pict_colbuf(p,p->lut,x,y)[2];
652     }
653     pix2loc(pc,&(p->resol),x,y);
654     hpx = (0.5/p->resol.xr);
655     hpy = (0.5/p->resol.yr);
656     pc[0] -= hpx;
657     pc[1] -= hpy;
658    
659     i = splane_normal(p,pc[0],pc[1],pc[0],(pc[1]+2.0*hpy),n[0]);
660     i &= splane_normal(p,pc[0],(pc[1]+2.0*hpy),(pc[0]+2.0*hpx),(pc[1]+2.0*hpy),n[1]);
661     i &= splane_normal(p,(pc[0]+2.0*hpx),(pc[1]+2.0*hpy),(pc[0]+2.0*hpx),pc[1],n[2]);
662     i &= splane_normal(p,(pc[0]+2.0*hpx),pc[1],pc[0],pc[1],n[3]);
663    
664     if (!i) {
665     return 0;
666     }
667     ang = 0;
668     for(i=0;i<4;i++) {
669     a = DOT(n[i],n[(i+1)%4]);
670     a = acos(a);
671     a = fabs(a);
672     ang += M_PI - a;
673     }
674     ang = ang - 2.0*M_PI;
675     if ((ang > (2.0*M_PI)) || ang < 0) {
676     fprintf(stderr,"Normal error in pict_get_sangle %f %d %d\n",ang,x,y);
677     return -1.0;
678     }
679     return ang;
680     }
681    
682     int pict_locate(pict* p,FVECT pt,int* x,int* y)
683     {
684     FVECT pp;
685     int res[2];
686    
687     if (pict_lut_used(p)) {
688     fprintf(stderr,"pict_locate: Not implemented for lut mode\n");
689     return 0;
690     }
691    
692 greg 2.5 if (viewloc(pp,&(p->view),pt) <= 0)
693 greg 2.1 return 0;
694     loc2pix(res,&(p->resol),pp[0],pp[1]);
695     *x = res[0];
696     *y = res[1];
697     return 1;
698     }
699    
700     int pict_get_sigma(pict* p,double x,double y,FVECT vdir,FVECT vup,double* s)
701     {
702     FVECT pvdir;
703    
704     if (!pict_get_dir(p,x,y,pvdir))
705     return 0;
706     if (normalize(vdir) == 0.0) {
707     fprintf(stderr,"get_sigma: view dir is zero\n");
708     return 0;
709     }
710     if (normalize(vup) == 0.0) {
711     fprintf(stderr,"get_sigma: view up is zero\n");
712     return 0;
713     }
714     *s = acos(DOT(vdir,pvdir));
715     return 1;
716     }
717    
718     int pict_get_tau(pict* p,double x,double y,FVECT vdir,FVECT vup,double* t)
719     {
720     FVECT hv,pvdir;
721     double s;
722     int i;
723     if (!pict_get_sigma(p,x,y,vdir,vup,&s))
724     return 0;
725     if (!pict_get_dir(p,x,y,pvdir))
726     return 0;
727     VCOPY(hv,pvdir);
728     normalize(hv);
729     for(i=0;i<3;i++){
730     hv[i] /= cos(s);
731     }
732     VSUB(hv,hv,vdir);
733     normalize(hv);
734     *t = acos(fdot(vup,hv));
735     return 1;
736     }
737    
738     static int ortho_coord(FVECT vdir,FVECT vup,FVECT vv,FVECT hv)
739     {
740     if (normalize(vdir) == 0.0) {
741     fprintf(stderr,"ortho_coord: view direction is zero\n");
742     return 0;
743     }
744     if (normalize(vup) == 0.0) {
745     fprintf(stderr,"ortho_coord: view direction is zero\n");
746     return 0;
747     }
748     fcross(hv,vdir,vup);
749     fcross(vv,vdir,hv);
750     return 1;
751     }
752    
753     int pict_get_vangle(pict* p,double x,double y,FVECT vdir,FVECT vup,double* a)
754     {
755     FVECT hv,vv,pvdir;
756     if (!ortho_coord(vdir,vup,vv,hv)) {
757     return 0;
758     }
759     if (!pict_get_dir(p,x,y,pvdir))
760     return 0;
761     *a = acos(fdot(vv,pvdir)) - (M_PI/2.0) ;
762     return 1;
763     }
764    
765     int pict_get_hangle(pict* p,double x,double y,FVECT vdir,FVECT vup,double* a)
766     {
767     FVECT hv,vv,pvdir;
768     if (!ortho_coord(vdir,vup,vv,hv)) {
769     return 0;
770     }
771     if (!pict_get_dir(p,x,y,pvdir)) {
772     return 0;
773     }
774     *a = ((M_PI/2.0) - acos(fdot(hv,pvdir)));
775     return 1;
776     }
777    
778     int pict_setorigin(pict* p,FVECT orig)
779     {
780     VCOPY(p->view.vp,orig);
781     return pict_update_view(p);
782     }
783    
784     int pict_setdir(pict* p,FVECT dir)
785     {
786     VCOPY(p->view.vdir,dir);
787     return pict_update_view(p);
788     }
789    
790     int pict_sethview(pict* p,double ha)
791     {
792     p->view.horiz = ha;
793     return pict_update_view(p);
794     }
795    
796     int pict_setvview(pict* p,double va)
797     {
798     p->view.vert = va;
799     return pict_update_view(p);
800     }
801    
802     int pict_setvtype(pict* p,char vt)
803     {
804     p->view.type = vt;
805     return pict_update_view(p);
806     }
807    
808     #ifdef TEST_PICTTOOL
809    
810     char* progname = "pictool";
811    
812     int main(int argc,char** argv)
813     {
814     pict* p = pict_create();
815     FVECT dir = {0,-1,0};
816     FVECT vdir,edir;
817     FVECT up = {0,0,1};
818     int x,y,i;
819     double sang,h,v;
820    
821     if (argc < 2) {
822     fprintf(stderr,"usage: %s <pic-file>\n",argv[0]);
823     return EXIT_FAILURE;
824     }
825     fprintf(stderr,"readdd\n");
826     pict_read(p,argv[1]);
827     pict_set_comment(p,"pictool test output\n");
828     /*
829     pict_get_dir(p,100,100,edir); // YDECR
830     for(y=0;y<p->resol.yr;y++) {
831     for(x=0;x<p->resol.xr;x++) {
832     pict_get_hangle(p,x,y,dir,up,&h);
833     pict_get_vangle(p,x,y,dir,up,&v);
834     pict_get_dir(p,x,y,vdir);
835     if (acos(DOT(vdir,edir)) < DEG2RAD(3))
836     setcolor(pict_get_color(p,x,y),h+M_PI/2,v+M_PI/2,1);
837     else
838     setcolor(pict_get_color(p,x,y),h+M_PI/2,v+M_PI/2,0);
839     }
840     }
841     */
842     printf("resolution: %dx%d\n",pict_get_xsize(p),pict_get_ysize(p));
843     sang = 0.0;
844    
845     pict_new_gli(p);
846     for(x=0;x<pict_get_xsize(p);x++)
847     for(y=0;y<pict_get_ysize(p);y++) {
848     pict_get_dir(p,x,y,dir);
849     if (DOT(dir,p->view.vdir) >= 0.0) {
850     pict_get_av_omega(p,0) += pict_get_omega(p,x,y);
851     }
852     sang = pict_get_sangle(p,x,y);
853     for(i=0;i<3;i++)
854     if (sang <= 0 && pict_is_validpixel(p,x,y)) {
855     pict_get_color(p,x,y)[i] = 1;
856     } else {
857     pict_get_color(p,x,y)[i] =0;
858     }
859     }
860     printf("s %f\n",pict_get_av_omega(p,0));
861     pict_write(p,"test.pic");
862     while(1) {
863     printf("pos\n");
864     if (scanf("%d %d",&x,&y) != 2)
865     continue;
866     pict_get_dir(p,x,y,dir);
867     printf("dir: \t\t %f %f %f\n",dir[0],dir[1],dir[2]);
868     printf("solid angle: \t %f \n",pict_get_sangle(p,x,y));
869     }
870     return EXIT_SUCCESS;
871     }
872    
873     #endif