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.1 by greg, Thu Oct 3 16:52:50 1996 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 BotMesopic      5.62e-3         /* top of scotopic range */
17 #define TopMesopic      5.62            /* bottom of photopic range */
18
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 60 | Line 61 | double Lw;
61                  return(Bldmin);
62          if (b >= bwmax-FTINY)
63                  return(Bldmax);
64 <        return(Bldmin + cf(Bl(Lw))*(Bldmax-Bldmin));
64 >        return(Bldmin + cf(b)*(Bldmax-Bldmin));
65   }
66  
67  
# Line 99 | Line 100 | double Lw;
100  
101  
102   int
102 shiftdir(bw)            /* compute shift direction for histogram */
103 double  bw;
104 {
105        if (what2do&DO_HSENS && cf(bw) - (bw - bwmin)/(Bldmax - bwmin))
106                return(1);
107        return(-1);
108 }
109
110
111 int
103   mkbrmap()                       /* make dynamic range map */
104   {
114        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++)
120 <                modhist[i] = bwhist[i];
121 <        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)
135 <                        return(-1);     /* Houston, we have a problem.... */
136 <                /* original looped here as well (BEGIN_L2) */
137 <                maxd = 0;               /* find largest overvis */
138 <                for (i = 0; i < HISTRES; i++)
139 <                        if (hdiffs[i] > maxd)
140 <                                maxd = hdiffs[maxi=i];
141 <                /* broke loop here when (maxd == 0) (BREAK_L2) */
142 <                for (sd = shiftdir((maxi+.5)/HISTRES*(bwmax-bwmin)+bwmin);
143 <                                hdiffs[maxi] == maxd; sd = -sd)
144 <                        for (i = maxi+sd; i >= 0 & i < HISTRES; i += sd)
145 <                                if (hdiffs[i] < 0) {
146 <                                        if (hdiffs[i] <= -maxd) {
147 <                                                modhist[i] += maxd;
148 <                                                modhist[maxi] -= maxd;
149 <                                                hdiffs[i] += maxd;
150 <                                                hdiffs[maxi] = 0;
151 <                                        } else {
152 <                                                modhist[maxi] += hdiffs[i];
153 <                                                modhist[i] -= hdiffs[i];
154 <                                                hdiffs[maxi] += hdiffs[i];
155 <                                                hdiffs[i] = 0;
156 <                                        }
157 <                                        break;
158 <                                }
159 <                /* (END_L2) */
160 <        }
161 <        return(0);
125 >        } while (trimmings > histot*CVRATIO);
126 >
127 >        return(0);                      /* we got it */
128   }
129  
130  
# Line 180 | Line 146 | int    xres;
146                          incolor = (Lw - BotMesopic) /
147                                          (TopMesopic - BotMesopic);
148                  if (incolor < 1.-FTINY) {
149 +                        b = (1.-incolor)*slum(scan[i])*inpexp/SWNORM;
150 +                        if (lumf == rgblum) b /= WHTEFFICACY;
151 +                        setcolor(ctmp, b, b, b);
152                          if (incolor <= FTINY)
153                                  setcolor(scan[i], 0., 0., 0.);
154                          else
155                                  scalecolor(scan[i], incolor);
187                        b = (1.-incolor)*slum(scan[i])*inpexp/SWNORM;
188                        if (lumf == rgblum) b /= WHTEFFICACY;
189                        setcolor(ctmp, b, b, b);
156                          addcolor(scan[i], ctmp);
157                  }
158          }
# Line 214 | Line 180 | int    xres;
180   }
181  
182  
183 < #ifdef DEBUG
184 < doplots()                       /* generate debugging plots */
183 > putmapping(fp)                  /* put out mapping function */
184 > FILE    *fp;
185   {
186 <        double  T, b, s;
221 <        FILE    *fp;
222 <        char    fname[128];
186 >        double  b, s;
187          register int    i;
188 +        double  wlum, sf, dlum;
189  
190 <        T = histot * (bwmax - bwmin) / HISTRES;
190 >        sf = scalef*inpexp;
191 >        if (lumf == cielum) sf *= WHTEFFICACY;
192          s = (bwmax - bwmin)/HISTRES;
193 <
194 <        sprintf(fname, "%s_hist.plt", infn);
195 <        if ((fp = fopen(fname, "w")) == NULL)
196 <                syserror(fname);
197 <        fputs("include=curve.plt\n", fp);
198 <        fputs("title=\"Brightness Frequency Distribution\"\n", fp);
199 <        fprintf(fp, "subtitle=%s\n", infn);
200 <        fputs("ymin=0\n", fp);
201 <        fputs("xlabel=\"Perceptual Brightness B(Lw)\"\n", fp);
202 <        fputs("ylabel=\"Frequency Count\"\n", fp);
237 <        fputs("Alabel=\"Histogram\"\n", fp);
238 <        fputs("Alintype=0\n", fp);
239 <        fputs("Blabel=\"Envelope\"\n", fp);
240 <        fputs("Bsymsize=0\n", fp);
241 <        fputs("Adata=\n", fp);
242 <        for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s)
243 <                fprintf(fp, "\t%f %d\n", b, modhist[i]);
244 <        fputs(";\nBdata=\n", fp);
245 <        for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s)
246 <                fprintf(fp, "\t%f %f\n", b, T*clampf(Lb(b)));
247 <        fputs(";\n", fp);
248 <        fclose(fp);
193 >        for (i = 0, b = bwmin + .5*s; i < HISTRES; i++, b += s) {
194 >                wlum = Lb(b);
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   }
250 #endif

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines