ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/pcond3.c
(Generate patch)

Comparing ray/src/px/pcond3.c (file contents):
Revision 3.4 by greg, Wed Jan 8 17:57:01 1997 UTC vs.
Revision 3.7 by greg, Sat Jan 11 10:39:25 1997 UTC

# Line 1 | Line 1
1 < /* Copyright (c) 1996 Regents of the University of California */
1 > /* Copyright (c) 1997 Regents of the University of California */
2  
3   #ifndef lint
4   static char SCCSid[] = "$SunId$ LBL";
# Line 11 | Line 11 | static char SCCSid[] = "$SunId$ LBL";
11   #include "pcond.h"
12  
13  
14 < #define CVRATIO         0.025           /* fraction of pixels allowed > env. */
14 > #define CVRATIO         0.025           /* fraction of samples allowed > env. */
15  
16   #define exp10(x)        exp(2.302585093*(x))
17  
18 < int     modhist[HISTRES];               /* modified histogram */
18 > float   modhist[HISTRES];               /* modified histogram */
19 > double  mhistot;                        /* modified histogram total */
20   float   cumf[HISTRES+1];                /* cumulative distribution function */
21  
22  
23   mkcumf()                        /* make cumulative distribution function */
24   {
25          register int    i;
26 <        register long   sum;
26 >        register double sum;
27  
28 <        cumf[0] = 0.;
29 <        sum = modhist[0];
30 <        for (i = 1; i < HISTRES; i++) {
31 <                cumf[i] = (double)sum/histot;
28 >        mhistot = 0.;           /* compute modified total */
29 >        for (i = 0; i < HISTRES; i++)
30 >                mhistot += modhist[i];
31 >
32 >        sum = 0.;               /* compute cumulative function */
33 >        for (i = 0; i < HISTRES; i++) {
34 >                cumf[i] = sum/mhistot;
35                  sum += modhist[i];
36          }
37          cumf[HISTRES] = 1.;
# Line 96 | Line 100 | double Lw;
100  
101  
102   int
99 shiftdir(bw)            /* compute shift direction for histogram */
100 double  bw;
101 {
102        if (what2do&DO_HSENS && cf(bw) - (bw - bwmin)/(Bldmax - bwmin))
103                return(1);
104        return(-1);
105 }
106
107
108 int
103   mkbrmap()                       /* make dynamic range map */
104   {
111        int     hdiffs[HISTRES], above, below;
105          double  T, b, s;
106 <        int     maxd, maxi, sd;
106 >        double  ceiling, trimmings;
107          register int    i;
108                                          /* copy initial histogram */
109 <        for (i = 0; i < HISTRES; i++)
117 <                modhist[i] = bwhist[i];
118 <        T = histot * (bwmax - bwmin) / HISTRES;
109 >        bcopy((char *)bwhist, (char *)modhist, sizeof(modhist));
110          s = (bwmax - bwmin)/HISTRES;
111                                          /* loop until satisfactory */
112 <        for ( ; ; ) {
113 <                mkcumf();               /* sync brightness mapping */
114 <                above = below = 0;      /* compute visibility overflow */
112 >        do {
113 >                mkcumf();                       /* sync brightness mapping */
114 >                if (mhistot <= histot*CVRATIO)
115 >                        return(-1);             /* no compression needed! */
116 >                T = mhistot * (bwmax - bwmin) / HISTRES;
117 >                trimmings = 0.;                 /* clip to envelope */
118                  for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s) {
119 <                        hdiffs[i] = modhist[i] - (int)(T*clampf(Lb(b)) + .5);
120 <                        if (hdiffs[i] > 0) above += hdiffs[i];
121 <                        else below -= hdiffs[i];
119 >                        ceiling = T*clampf(Lb(b));
120 >                        if (modhist[i] > ceiling) {
121 >                                trimmings += modhist[i] - ceiling;
122 >                                modhist[i] = ceiling;
123 >                        }
124                  }
125 <                if (above <= histot*CVRATIO)
126 <                        break;          /* close enough */
127 <                if (above-below >= 0)
132 <                        return(-1);     /* Houston, we have a problem.... */
133 <                /* original looped here as well (BEGIN_L2) */
134 <                maxd = 0;               /* find largest overvis */
135 <                for (i = 0; i < HISTRES; i++)
136 <                        if (hdiffs[i] > maxd)
137 <                                maxd = hdiffs[maxi=i];
138 <                /* broke loop here when (maxd == 0) (BREAK_L2) */
139 <                for (sd = shiftdir((maxi+.5)/HISTRES*(bwmax-bwmin)+bwmin);
140 <                                hdiffs[maxi] == maxd; sd = -sd)
141 <                        for (i = maxi+sd; i >= 0 & i < HISTRES; i += sd)
142 <                                if (hdiffs[i] < 0) {
143 <                                        if (hdiffs[i] <= -maxd) {
144 <                                                modhist[i] += maxd;
145 <                                                modhist[maxi] -= maxd;
146 <                                                hdiffs[i] += maxd;
147 <                                                hdiffs[maxi] = 0;
148 <                                        } else {
149 <                                                modhist[maxi] += hdiffs[i];
150 <                                                modhist[i] -= hdiffs[i];
151 <                                                hdiffs[maxi] += hdiffs[i];
152 <                                                hdiffs[i] = 0;
153 <                                        }
154 <                                        break;
155 <                                }
156 <                /* (END_L2) */
157 <        }
158 <        return(0);
125 >        } while (trimmings > histot*CVRATIO);
126 >
127 >        return(0);                      /* we got it */
128   }
129  
130  
# Line 216 | Line 185 | FILE   *fp;
185   {
186          double  b, s;
187          register int    i;
188 <        double  wlum, sf;
188 >        double  wlum, sf, dlum;
189  
190          sf = scalef*inpexp;
191          if (lumf == cielum) sf *= WHTEFFICACY;
192          s = (bwmax - bwmin)/HISTRES;
193          for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s) {
194                  wlum = Lb(b);
195 <                if (what2do&DO_LINEAR)
196 <                        fprintf(fp, "%e %e\n", wlum, sf*wlum);
197 <                else
195 >                if (what2do&DO_LINEAR) {
196 >                        dlum = sf*wlum;
197 >                        if (dlum > ldmax) dlum = ldmax;
198 >                        else if (dlum < ldmin) dlum = ldmin;
199 >                        fprintf(fp, "%e %e\n", wlum, dlum);
200 >                } else
201                          fprintf(fp, "%e %e\n", wlum, Lb(BLw(wlum)));
202          }
203   }
232
233
234 #ifdef DEBUG
235 doplots()                       /* generate debugging plots */
236 {
237        double  T, b, s;
238        FILE    *fp;
239        char    fname[128];
240        register int    i;
241
242        T = histot * (bwmax - bwmin) / HISTRES;
243        s = (bwmax - bwmin)/HISTRES;
244
245        sprintf(fname, "%s_hist.plt", infn);
246        if ((fp = fopen(fname, "w")) == NULL)
247                syserror(fname);
248        fputs("include=curve.plt\n", fp);
249        fputs("title=\"Brightness Frequency Distribution\"\n", fp);
250        fprintf(fp, "subtitle=%s\n", infn);
251        fputs("ymin=0\n", fp);
252        fputs("xlabel=\"Perceptual Brightness B(Lw)\"\n", fp);
253        fputs("ylabel=\"Frequency Count\"\n", fp);
254        fputs("Alabel=\"Histogram\"\n", fp);
255        fputs("Alintype=0\n", fp);
256        fputs("Blabel=\"Envelope\"\n", fp);
257        fputs("Bsymsize=0\n", fp);
258        fputs("Adata=\n", fp);
259        for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s)
260                fprintf(fp, "\t%f %d\n", b, modhist[i]);
261        fputs(";\nBdata=\n", fp);
262        for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s)
263                fprintf(fp, "\t%f %f\n", b, T*clampf(Lb(b)));
264        fputs(";\n", fp);
265        fclose(fp);
266 }
267 #endif

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines