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

Comparing ray/src/common/tonemap.c (file contents):
Revision 3.7 by gregl, Mon Dec 1 09:52:00 1997 UTC vs.
Revision 3.16 by greg, Fri Jan 7 20:33:02 2005 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1997 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   * Tone mapping functions.
6   * See tonemap.h for detailed function descriptions.
7 + * Added von Kries white-balance calculations 10/01 (GW).
8 + *
9 + * Externals declared in tonemap.h
10   */
11  
12 + #include "copyright.h"
13 +
14   #include        <stdio.h>
15   #include        <math.h>
16   #include        "tmprivat.h"
# Line 16 | Line 18 | static char SCCSid[] = "$SunId$ LBL";
18  
19   #define exp10(x)        exp(M_LN10*(x))
20  
19 struct tmStruct *tmTop = NULL;          /* current tone mapping stack */
20
21                                          /* our list of conversion packages */
22   struct tmPackage        *tmPkg[TM_MAXPKG];
23   int     tmNumPkgs = 0;                  /* number of registered packages */
# Line 26 | Line 26 | int    tmLastError;                    /* last error incurred by library *
26   char    *tmLastFunction;                /* error-generating function name */
27  
28  
29 < int
30 < tmErrorReturn(func, err)                /* error return (with message) */
31 < char    *func;
32 < int     err;
29 > TMstruct *
30 > tmInit(                                 /* initialize new tone mapping */
31 > int     flags,
32 > RGBPRIMP        monpri,
33 > double  gamval
34 > )
35   {
34        tmLastFunction = func;
35        tmLastError = err;
36        if (tmTop != NULL && tmTop->flags & TM_F_NOSTDERR)
37                return(err);
38        fputs(func, stderr);
39        fputs(": ", stderr);
40        fputs(tmErrorMessage[err], stderr);
41        fputs("!\n", stderr);
42        return(err);
43 }
44
45
46 struct tmStruct *
47 tmInit(flags, monpri, gamval)           /* initialize new tone mapping */
48 int     flags;
49 RGBPRIMP        monpri;
50 double  gamval;
51 {
52        static char     funcName[] = "tmInit";
36          COLORMAT        cmat;
37 <        register struct tmStruct        *tmnew;
38 <        register int    i;
37 >        TMstruct        *tmnew;
38 >        int     i;
39                                                  /* allocate structure */
40 <        tmnew = (struct tmStruct *)malloc(sizeof(struct tmStruct));
40 >        tmnew = (TMstruct *)malloc(sizeof(TMstruct));
41          if (tmnew == NULL)
42                  return(NULL);
43  
44          tmnew->flags = flags & ~TM_F_UNIMPL;
45 +        if (tmnew->flags & TM_F_BW)
46 +                tmnew->flags &= ~TM_F_MESOPIC;
47                                                  /* set monitor transform */
48          if (monpri == NULL || monpri == stdprims || tmnew->flags & TM_F_BW) {
49                  tmnew->monpri = stdprims;
# Line 66 | Line 51 | double gamval;
51                  tmnew->clf[GRN] = rgb2xyzmat[1][1];
52                  tmnew->clf[BLU] = rgb2xyzmat[1][2];
53          } else {
54 <                comprgb2xyzmat(cmat, tmnew->monpri=monpri);
54 >                comprgb2xyzWBmat(cmat, tmnew->monpri=monpri);
55                  tmnew->clf[RED] = cmat[1][0];
56                  tmnew->clf[GRN] = cmat[1][1];
57                  tmnew->clf[BLU] = cmat[1][2];
# Line 86 | Line 71 | double gamval;
71                          tmnew->inpsf = WHTEFFICACY;
72          tmnew->cmat[0][1] = tmnew->cmat[0][2] = tmnew->cmat[1][0] =
73          tmnew->cmat[1][2] = tmnew->cmat[2][0] = tmnew->cmat[2][1] = 0.;
74 <        tmnew->hbrmin = tmnew->hbrmax = 0;
74 >        tmnew->hbrmin = 10; tmnew->hbrmax = -10;
75          tmnew->histo = NULL;
76 <        tmnew->mbrmin = tmnew->mbrmax = 0;
76 >        tmnew->mbrmin = 10; tmnew->mbrmax = -10;
77          tmnew->lumap = NULL;
78                                                  /* zero private data */
79          for (i = TM_MAXPKG; i--; )
80                  tmnew->pd[i] = NULL;
81 <                                                /* make tmnew current */
82 <        tmnew->tmprev = tmTop;
98 <        return(tmTop = tmnew);
81 >                                                /* return new TMstruct */
82 >        return(tmnew);
83   }
84  
85  
86   int
87 < tmSetSpace(pri, sf)             /* set input color space for conversions */
88 < RGBPRIMP        pri;
89 < double  sf;
87 > tmSetSpace(                     /* set input color space for conversions */
88 > TMstruct        *tms,
89 > RGBPRIMP        pri,
90 > double  sf
91 > )
92   {
93          static char     funcName[] = "tmSetSpace";
94 <        register int    i, j;
94 >        int     i, j;
95                                                  /* error check */
96 <        if (tmTop == NULL)
96 >        if (tms == NULL)
97                  returnErr(TM_E_TMINVAL);
98          if (sf <= 1e-12)
99                  returnErr(TM_E_ILLEGAL);
100                                                  /* check if no change */
101 <        if (pri == tmTop->inppri && FEQ(sf, tmTop->inpsf))
101 >        if (pri == tms->inppri && FEQ(sf, tms->inpsf))
102                  returnOK;
103 <        tmTop->inppri = pri;                    /* let's set it */
104 <        tmTop->inpsf = sf;
103 >        tms->inppri = pri;                      /* let's set it */
104 >        tms->inpsf = sf;
105  
106 <        if (tmTop->flags & TM_F_BW) {           /* color doesn't matter */
107 <                tmTop->monpri = tmTop->inppri;          /* eliminate xform */
108 <                if (tmTop->inppri == TM_XYZPRIM) {
109 <                        tmTop->clf[CIEX] = tmTop->clf[CIEZ] = 0.;
110 <                        tmTop->clf[CIEY] = 1.;
106 >        if (tms->flags & TM_F_BW) {             /* color doesn't matter */
107 >                tms->monpri = tms->inppri;              /* eliminate xform */
108 >                if (tms->inppri == TM_XYZPRIM) {
109 >                        tms->clf[CIEX] = tms->clf[CIEZ] = 0.;
110 >                        tms->clf[CIEY] = 1.;
111                  } else {
112 <                        comprgb2xyzmat(tmTop->cmat, tmTop->monpri);
113 <                        tmTop->clf[RED] = tmTop->cmat[1][0];
114 <                        tmTop->clf[GRN] = tmTop->cmat[1][1];
115 <                        tmTop->clf[BLU] = tmTop->cmat[1][2];
112 >                        comprgb2xyzWBmat(tms->cmat, tms->monpri);
113 >                        tms->clf[RED] = tms->cmat[1][0];
114 >                        tms->clf[GRN] = tms->cmat[1][1];
115 >                        tms->clf[BLU] = tms->cmat[1][2];
116                  }
117 <                tmTop->cmat[0][0] = tmTop->cmat[1][1] = tmTop->cmat[2][2] =
118 <                                tmTop->inpsf;
119 <                tmTop->cmat[0][1] = tmTop->cmat[0][2] = tmTop->cmat[1][0] =
120 <                tmTop->cmat[1][2] = tmTop->cmat[2][0] = tmTop->cmat[2][1] = 0.;
117 >                tms->cmat[0][0] = tms->cmat[1][1] = tms->cmat[2][2] =
118 >                                tms->inpsf;
119 >                tms->cmat[0][1] = tms->cmat[0][2] = tms->cmat[1][0] =
120 >                tms->cmat[1][2] = tms->cmat[2][0] = tms->cmat[2][1] = 0.;
121  
122 <        } else if (tmTop->inppri == TM_XYZPRIM) /* input is XYZ */
123 <                compxyz2rgbmat(tmTop->cmat, tmTop->monpri);
122 >        } else if (tms->inppri == TM_XYZPRIM)   /* input is XYZ */
123 >                compxyz2rgbWBmat(tms->cmat, tms->monpri);
124  
125          else {                                  /* input is RGB */
126 <                if (tmTop->inppri != tmTop->monpri &&
127 <                                PRIMEQ(tmTop->inppri, tmTop->monpri))
128 <                        tmTop->inppri = tmTop->monpri;  /* no xform */
129 <                comprgb2rgbmat(tmTop->cmat, tmTop->inppri, tmTop->monpri);
126 >                if (tms->inppri != tms->monpri &&
127 >                                PRIMEQ(tms->inppri, tms->monpri))
128 >                        tms->inppri = tms->monpri;      /* no xform */
129 >                comprgb2rgbWBmat(tms->cmat, tms->inppri, tms->monpri);
130          }
131          for (i = 0; i < 3; i++)
132                  for (j = 0; j < 3; j++)
133 <                        tmTop->cmat[i][j] *= tmTop->inpsf;
133 >                        tms->cmat[i][j] *= tms->inpsf;
134                                                  /* set color divisors */
135          for (i = 0; i < 3; i++)
136 <                if (tmTop->clf[i] > .001)
137 <                        tmTop->cdiv[i] =
138 <                                256.*pow(tmTop->clf[i], 1./tmTop->mongam);
136 >                if (tms->clf[i] > .001)
137 >                        tms->cdiv[i] =
138 >                                256.*pow(tms->clf[i], 1./tms->mongam);
139                  else
140 <                        tmTop->cdiv[i] = 1;
140 >                        tms->cdiv[i] = 1;
141                                                  /* notify packages */
142          for (i = tmNumPkgs; i--; )
143 <                if (tmTop->pd[i] != NULL && tmPkg[i]->NewSpace != NULL)
144 <                        (*tmPkg[i]->NewSpace)(tmTop);
143 >                if (tms->pd[i] != NULL && tmPkg[i]->NewSpace != NULL)
144 >                        (*tmPkg[i]->NewSpace)(tms);
145          returnOK;
146   }
147  
148  
149   void
150 < tmClearHisto()                  /* clear current histogram */
150 > tmClearHisto(                   /* clear current histogram */
151 > TMstruct        *tms
152 > )
153   {
154 <        if (tmTop == NULL || tmTop->histo == NULL)
154 >        if (tms == NULL || tms->histo == NULL)
155                  return;
156 <        free((MEM_PTR)tmTop->histo);
157 <        tmTop->histo = NULL;
156 >        free((MEM_PTR)tms->histo);
157 >        tms->histo = NULL;
158   }
159  
160  
161   int
162 < tmCvColors(ls, cs, scan, len)           /* convert float colors */
163 < TMbright        *ls;
164 < BYTE    *cs;
165 < COLOR   *scan;
166 < int     len;
162 > tmCvColors(                             /* convert float colors */
163 > TMstruct        *tms,
164 > TMbright        *ls,
165 > BYTE    *cs,
166 > COLOR   *scan,
167 > int     len
168 > )
169   {
170          static char     funcName[] = "tmCvColors";
171 <        static COLOR    csmall = {1e-6, 1e-6, 1e-6};
171 >        static COLOR    csmall = {.5*MINLUM, .5*MINLUM, .5*MINLUM};
172          COLOR   cmon;
173          double  lum, slum;
174 <        register double d;
175 <        register int    i;
174 >        double  d;
175 >        int     i;
176  
177 <        if (tmTop == NULL)
177 >        if (tms == NULL)
178                  returnErr(TM_E_TMINVAL);
179 <        if (ls == NULL | scan == NULL | len <= 0)
179 >        if ((ls == NULL) | (scan == NULL) | (len < 0))
180                  returnErr(TM_E_ILLEGAL);
181          for (i = len; i--; ) {
182 <                if (tmNeedMatrix(tmTop))                /* get monitor RGB */
183 <                        colortrans(cmon, tmTop->cmat, scan[i]);
184 <                else {
185 <                        cmon[RED] = tmTop->inpsf*scan[i][RED];
186 <                        cmon[GRN] = tmTop->inpsf*scan[i][GRN];
187 <                        cmon[BLU] = tmTop->inpsf*scan[i][BLU];
182 >                if (tmNeedMatrix(tms)) {                /* get monitor RGB */
183 >                        colortrans(cmon, tms->cmat, scan[i]);
184 >                } else {
185 >                        cmon[RED] = tms->inpsf*scan[i][RED];
186 >                        cmon[GRN] = tms->inpsf*scan[i][GRN];
187 >                        cmon[BLU] = tms->inpsf*scan[i][BLU];
188                  }
189                                                          /* world luminance */
190 <                lum =   tmTop->clf[RED]*cmon[RED] +
191 <                        tmTop->clf[GRN]*cmon[GRN] +
192 <                        tmTop->clf[BLU]*cmon[BLU] ;
190 >                lum =   tms->clf[RED]*cmon[RED] +
191 >                        tms->clf[GRN]*cmon[GRN] +
192 >                        tms->clf[BLU]*cmon[BLU] ;
193                                                          /* check range */
194                  if (clipgamut(cmon, lum, CGAMUT_LOWER, csmall, cwhite))
195 <                        lum =   tmTop->clf[RED]*cmon[RED] +
196 <                                tmTop->clf[GRN]*cmon[GRN] +
197 <                                tmTop->clf[BLU]*cmon[BLU] ;
195 >                        lum =   tms->clf[RED]*cmon[RED] +
196 >                                tms->clf[GRN]*cmon[GRN] +
197 >                                tms->clf[BLU]*cmon[BLU] ;
198                  if (lum < MINLUM) {
199                          ls[i] = MINBRT-1;               /* bogus value */
200                          lum = MINLUM;
# Line 214 | Line 204 | int    len;
204                  }
205                  if (cs == TM_NOCHROM)                   /* no color? */
206                          continue;
207 <                if (tmTop->flags & TM_F_MESOPIC && lum < LMESUPPER) {
207 >                if (tms->flags & TM_F_MESOPIC && lum < LMESUPPER) {
208                          slum = scotlum(cmon);           /* mesopic adj. */
209                          if (lum < LMESLOWER)
210                                  cmon[RED] = cmon[GRN] = cmon[BLU] = slum;
211                          else {
212                                  d = (lum - LMESLOWER)/(LMESUPPER - LMESLOWER);
213 <                                if (tmTop->flags & TM_F_BW)
213 >                                if (tms->flags & TM_F_BW)
214                                          cmon[RED] = cmon[GRN] =
215                                                          cmon[BLU] = d*lum;
216                                  else
# Line 230 | Line 220 | int    len;
220                                  cmon[GRN] += d;
221                                  cmon[BLU] += d;
222                          }
223 <                } else if (tmTop->flags & TM_F_BW) {
223 >                } else if (tms->flags & TM_F_BW) {
224                          cmon[RED] = cmon[GRN] = cmon[BLU] = lum;
225                  }
226 <                d = tmTop->clf[RED]*cmon[RED]/lum;
226 >                d = tms->clf[RED]*cmon[RED]/lum;
227                  cs[3*i  ] = d>=.999 ? 255 :
228 <                                (int)(256.*pow(d, 1./tmTop->mongam));
229 <                d = tmTop->clf[GRN]*cmon[GRN]/lum;
228 >                                (int)(256.*pow(d, 1./tms->mongam));
229 >                d = tms->clf[GRN]*cmon[GRN]/lum;
230                  cs[3*i+1] = d>=.999 ? 255 :
231 <                                (int)(256.*pow(d, 1./tmTop->mongam));
232 <                d = tmTop->clf[BLU]*cmon[BLU]/lum;
231 >                                (int)(256.*pow(d, 1./tms->mongam));
232 >                d = tms->clf[BLU]*cmon[BLU]/lum;
233                  cs[3*i+2] = d>=.999 ? 255 :
234 <                                (int)(256.*pow(d, 1./tmTop->mongam));
234 >                                (int)(256.*pow(d, 1./tms->mongam));
235          }
236          returnOK;
237   }
238  
239  
240   int
241 < tmAddHisto(ls, len, wt)                 /* add values to histogram */
242 < register TMbright       *ls;
243 < int     len;
244 < int     wt;
241 > tmCvGrays(                              /* convert float gray values */
242 > TMstruct        *tms,
243 > TMbright        *ls,
244 > float   *scan,
245 > int     len
246 > )
247   {
248 <        static char     funcName[] = "tmAddHisto";
249 <        int     oldorig, oldlen, horig, hlen;
250 <        register int    i, j;
248 >        static char     funcName[] = "tmCvGrays";
249 >        double  d;
250 >        int     i;
251  
252 <        if (len <= 0)
252 >        if (tms == NULL)
253 >                returnErr(TM_E_TMINVAL);
254 >        if ((ls == NULL) | (scan == NULL) | (len < 0))
255                  returnErr(TM_E_ILLEGAL);
256 <        if (tmTop == NULL)
256 >        for (i = len; i--; )
257 >                if (scan[i] <= TM_NOLUM) {
258 >                        ls[i] = TM_NOBRT;               /* bogus value */
259 >                } else {
260 >                        d = TM_BRTSCALE*log(scan[i]);   /* encode it */
261 >                        ls[i] = d>0. ? (int)(d+.5) : (int)(d-.5);
262 >                }
263 >        returnOK;
264 > }
265 >
266 >
267 > int
268 > tmAddHisto(                             /* add values to histogram */
269 > TMstruct        *tms,
270 > TMbright        *ls,
271 > int     len,
272 > int     wt
273 > )
274 > {
275 >        static char     funcName[] = "tmAddHisto";
276 >        int     oldorig=0, oldlen, horig, hlen;
277 >        int     i, j;
278 >
279 >        if (tms == NULL)
280                  returnErr(TM_E_TMINVAL);
281 +        if (len < 0)
282 +                returnErr(TM_E_ILLEGAL);
283 +        if (len == 0)
284 +                returnOK;
285                                                  /* first, grow limits */
286 <        if (tmTop->histo == NULL) {
286 >        if (tms->histo == NULL) {
287                  for (i = len; i-- && ls[i] < MINBRT; )
288                          ;
289                  if (i < 0)
290                          returnOK;
291 <                tmTop->hbrmin = tmTop->hbrmax = ls[i];
291 >                tms->hbrmin = tms->hbrmax = ls[i];
292                  oldlen = 0;
293          } else {
294 <                oldorig = (tmTop->hbrmin-MINBRT)/HISTEP;
295 <                oldlen = (tmTop->hbrmax-MINBRT)/HISTEP + 1 - oldorig;
294 >                oldorig = (tms->hbrmin-MINBRT)/HISTEP;
295 >                oldlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - oldorig;
296          }
297          for (i = len; i--; ) {
298                  if ((j = ls[i]) < MINBRT)
299                          continue;
300 <                if (j < tmTop->hbrmin)
301 <                        tmTop->hbrmin = j;
302 <                else if (j > tmTop->hbrmax)
303 <                        tmTop->hbrmax = j;
300 >                if (j < tms->hbrmin)
301 >                        tms->hbrmin = j;
302 >                else if (j > tms->hbrmax)
303 >                        tms->hbrmax = j;
304          }
305 <        horig = (tmTop->hbrmin-MINBRT)/HISTEP;
306 <        hlen = (tmTop->hbrmax-MINBRT)/HISTEP + 1 - horig;
305 >        horig = (tms->hbrmin-MINBRT)/HISTEP;
306 >        hlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - horig;
307          if (hlen > oldlen) {                    /* (re)allocate histogram */
308 <                register int    *newhist = (int *)calloc(hlen, sizeof(int));
308 >                int     *newhist = (int *)calloc(hlen, sizeof(int));
309                  if (newhist == NULL)
310                          returnErr(TM_E_NOMEM);
311                  if (oldlen) {                   /* copy and free old */
312                          for (i = oldlen, j = i+oldorig-horig; i; )
313 <                                newhist[--j] = tmTop->histo[--i];
314 <                        free((MEM_PTR)tmTop->histo);
313 >                                newhist[--j] = tms->histo[--i];
314 >                        free((MEM_PTR)tms->histo);
315                  }
316 <                tmTop->histo = newhist;
316 >                tms->histo = newhist;
317          }
318          if (wt == 0)
319                  returnOK;
320          for (i = len; i--; )                    /* add in new counts */
321                  if (ls[i] >= MINBRT)
322 <                        tmTop->histo[ (ls[i]-MINBRT)/HISTEP - horig ] += wt;
322 >                        tms->histo[ (ls[i]-MINBRT)/HISTEP - horig ] += wt;
323          returnOK;
324   }
325  
326  
327   static double
328 < htcontrs(La)            /* human threshold contrast sensitivity, dL(La) */
329 < double  La;
328 > htcontrs(               /* human threshold contrast sensitivity, dL(La) */
329 > double  La
330 > )
331   {
332          double  l10La, l10dL;
333                                  /* formula taken from Ferwerda et al. [SG96] */
# Line 325 | Line 347 | double La;
347   }
348  
349  
350 + static int
351 + tmNewMap(                       /* allocate new tone-mapping array */
352 + TMstruct        *tms
353 + )
354 + {
355 +        if (tms->lumap != NULL && (tms->mbrmax - tms->mbrmin) !=
356 +                                        (tms->hbrmax - tms->hbrmin)) {
357 +                free((MEM_PTR)tms->lumap);
358 +                tms->lumap = NULL;
359 +        }
360 +        tms->mbrmin = tms->hbrmin;
361 +        tms->mbrmax = tms->hbrmax;
362 +        if (tms->mbrmin > tms->mbrmax)
363 +                return 0;
364 +        if (tms->lumap == NULL)
365 +                tms->lumap = (unsigned short *)malloc(sizeof(unsigned short)*
366 +                                        (tms->mbrmax-tms->mbrmin+1));
367 +        return(tms->lumap != NULL);
368 + }
369 +
370 +
371   int
372 < tmComputeMapping(gamval, Lddyn, Ldmax)
373 < double  gamval;
374 < double  Lddyn;
375 < double  Ldmax;
372 > tmFixedMapping(                 /* compute fixed, linear tone-mapping */
373 > TMstruct        *tms,
374 > double  expmult,
375 > double  gamval
376 > )
377   {
378 +        static char     funcName[] = "tmFixedMapping";
379 +        double          d;
380 +        int     i;
381 +        
382 +        if (!tmNewMap(tms))
383 +                returnErr(TM_E_NOMEM);
384 +        if (expmult <= .0)
385 +                expmult = 1.;
386 +        if (gamval < MINGAM)
387 +                gamval = tms->mongam;
388 +        d = log(expmult/tms->inpsf);
389 +        for (i = tms->mbrmax-tms->mbrmin+1; i--; )
390 +                tms->lumap[i] = 256. * exp(
391 +                        ( d + (tms->mbrmin+i)*(1./TM_BRTSCALE) )
392 +                        / gamval );
393 +        returnOK;
394 + }
395 +
396 +
397 + int
398 + tmComputeMapping(                       /* compute histogram tone-mapping */
399 + TMstruct        *tms,
400 + double  gamval,
401 + double  Lddyn,
402 + double  Ldmax
403 + )
404 + {
405          static char     funcName[] = "tmComputeMapping";
406          int     *histo;
407          float   *cumf;
408 <        int     brt0, histlen, histot, threshold, ceiling, trimmings;
408 >        int     brt0, histlen, threshold, ceiling, trimmings;
409          double  logLddyn, Ldmin, Ldavg, Lwavg, Tr, Lw, Ld;
410 <        int4    sum;
411 <        register double d;
412 <        register int    i, j;
410 >        int32   histot;
411 >        double  sum;
412 >        double  d;
413 >        int     i, j;
414  
415 <        if (tmTop == NULL || tmTop->histo == NULL)
415 >        if (tms == NULL || tms->histo == NULL)
416                  returnErr(TM_E_TMINVAL);
417                                          /* check arguments */
418          if (Lddyn < MINLDDYN) Lddyn = DEFLDDYN;
419          if (Ldmax < MINLDMAX) Ldmax = DEFLDMAX;
420 <        if (gamval < MINGAM) gamval = tmTop->mongam;
420 >        if (gamval < MINGAM) gamval = tms->mongam;
421                                          /* compute handy values */
422          Ldmin = Ldmax/Lddyn;
423          logLddyn = log(Lddyn);
424          Ldavg = sqrt(Ldmax*Ldmin);
425 <        i = (tmTop->hbrmin-MINBRT)/HISTEP;
425 >        i = (tms->hbrmin-MINBRT)/HISTEP;
426          brt0 = MINBRT + HISTEP/2 + i*HISTEP;
427 <        histlen = (tmTop->hbrmax-MINBRT)/HISTEP + 1 - i;
427 >        histlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - i;
428                                          /* histogram total and mean */
429          histot = 0; sum = 0;
430          j = brt0 + histlen*HISTEP;
431          for (i = histlen; i--; ) {
432 <                histot += tmTop->histo[i];
433 <                sum += (j -= HISTEP) * tmTop->histo[i];
432 >                histot += tms->histo[i];
433 >                sum += (j -= HISTEP) * tms->histo[i];
434          }
435 <        threshold = histot*.025 + .5;
435 >        threshold = histot*0.005 + .5;
436          if (threshold < 4)
437                  returnErr(TM_E_TMFAIL);
438          Lwavg = tmLuminance( (double)sum / histot );
439 <        if (!(tmTop->flags & TM_F_LINEAR)) {    /* clamp histogram */
440 <                histo = (int *)malloc(histlen*sizeof(int));
441 <                cumf = (float *)malloc((histlen+1)*sizeof(float));
442 <                if (histo == NULL | cumf == NULL)
443 <                        returnErr(TM_E_NOMEM);
444 <                for (i = histlen; i--; )        /* make malleable copy */
445 <                        histo[i] = tmTop->histo[i];
446 <                do {                            /* iterate to solution */
447 <                        sum = 0;                /* cumulative probability */
448 <                        for (i = 0; i < histlen; i++) {
449 <                                cumf[i] = (double)sum/histot;
450 <                                sum += histo[i];
439 >                                        /* allocate space for mapping */
440 >        if (!tmNewMap(tms))
441 >                returnErr(TM_E_NOMEM);
442 >                                        /* use linear tone mapping? */
443 >        if (tms->flags & TM_F_LINEAR)
444 >                goto linearmap;
445 >                                        /* clamp histogram */
446 >        histo = (int *)malloc(histlen*sizeof(int));
447 >        cumf = (float *)malloc((histlen+2)*sizeof(float));
448 >        if ((histo == NULL) | (cumf == NULL))
449 >                returnErr(TM_E_NOMEM);
450 >        cumf[histlen+1] = 1.;           /* guard for assignment code */
451 >        for (i = histlen; i--; )        /* make malleable copy */
452 >                histo[i] = tms->histo[i];
453 >        do {                            /* iterate to solution */
454 >                sum = 0;                /* cumulative probability */
455 >                for (i = 0; i < histlen; i++) {
456 >                        cumf[i] = (double)sum/histot;
457 >                        sum += histo[i];
458 >                }
459 >                cumf[histlen] = 1.;
460 >                Tr = histot * (double)(tms->hbrmax - tms->hbrmin) /
461 >                        ((double)histlen*TM_BRTSCALE) / logLddyn;
462 >                ceiling = Tr + 1.;
463 >                trimmings = 0;          /* clip to envelope */
464 >                for (i = histlen; i--; ) {
465 >                        if (tms->flags & TM_F_HCONTR) {
466 >                                Lw = tmLuminance(brt0 + i*HISTEP);
467 >                                Ld = Ldmin * exp( logLddyn *
468 >                                        .5*(cumf[i]+cumf[i+1]) );
469 >                                ceiling = Tr * (htcontrs(Ld) * Lw) /
470 >                                        (htcontrs(Lw) * Ld) + 1.;
471                          }
472 <                        cumf[i] = 1.;
473 <                        Tr = histot * (double)(tmTop->hbrmax - tmTop->hbrmin) /
474 <                                ((double)histlen*TM_BRTSCALE) / logLddyn;
383 <                        ceiling = Tr + 1.;
384 <                        trimmings = 0;                  /* clip to envelope */
385 <                        for (i = histlen; i--; ) {
386 <                                if (tmTop->flags & TM_F_HCONTR) {
387 <                                        Lw = tmLuminance(brt0 + i*HISTEP);
388 <                                        Ld = Ldmin * exp( logLddyn *
389 <                                                .5*(cumf[i]+cumf[i+1]) );
390 <                                        ceiling = Tr * (htcontrs(Ld) * Lw) /
391 <                                                (htcontrs(Lw) * Ld) + 1.;
392 <                                }
393 <                                if (histo[i] > ceiling) {
394 <                                        trimmings += histo[i] - ceiling;
395 <                                        histo[i] = ceiling;
396 <                                }
472 >                        if (histo[i] > ceiling) {
473 >                                trimmings += histo[i] - ceiling;
474 >                                histo[i] = ceiling;
475                          }
398                } while ((histot -= trimmings) > threshold &&
399                                                trimmings > threshold);
400        }
401                                                /* allocate luminance map */
402        if (tmTop->lumap != NULL)
403                free((MEM_PTR)tmTop->lumap);
404        tmTop->mbrmin = tmTop->hbrmin;
405        tmTop->mbrmax = tmTop->hbrmax;
406        tmTop->lumap = (unsigned short *)malloc(
407                (tmTop->mbrmax-tmTop->mbrmin+1)*sizeof(unsigned short) );
408        if (tmTop->lumap == NULL)
409                returnErr(TM_E_NOMEM);
410        if (tmTop->flags & TM_F_LINEAR || histot <= threshold) {
411                                                /* linear tone mapping */
412                if (tmTop->flags & TM_F_HCONTR)
413                        d = htcontrs(Ldavg) / htcontrs(Lwavg);
414                else
415                        d = Ldavg / Lwavg;
416                d = log(d/Ldmax);
417                for (i = tmTop->mbrmax-tmTop->mbrmin+1; i--; )
418                        tmTop->lumap[i] = 256. * exp(
419                                ( d + (tmTop->mbrmin+i)/(double)TM_BRTSCALE )
420                                / gamval );
421        } else {
422                                                /* histogram adjustment */
423                for (i = tmTop->mbrmax-tmTop->mbrmin+1; i--; ) {
424                        j = d = (double)i/(tmTop->mbrmax-tmTop->mbrmin)*histlen;
425                        d -= (double)j;
426                        Ld = Ldmin*exp(logLddyn*((1.-d)*cumf[j]+d*cumf[j+1]));
427                        d = (Ld - Ldmin)/(Ldmax - Ldmin);
428                        tmTop->lumap[i] = 256.*pow(d, 1./gamval);
476                  }
477 +                                        /* check if we're out of data */
478 +                if ((histot -= trimmings) <= threshold) {
479 +                        free((MEM_PTR)histo);
480 +                        free((MEM_PTR)cumf);
481 +                        goto linearmap;
482 +                }
483 +        } while (trimmings > threshold);
484 +                                        /* assign tone-mapping */
485 +        for (i = tms->mbrmax-tms->mbrmin+1; i--; ) {
486 +                j = d = (double)i/(tms->mbrmax-tms->mbrmin)*histlen;
487 +                d -= (double)j;
488 +                Ld = Ldmin*exp(logLddyn*((1.-d)*cumf[j]+d*cumf[j+1]));
489 +                d = (Ld - Ldmin)/(Ldmax - Ldmin);
490 +                tms->lumap[i] = 256.*pow(d, 1./gamval);
491          }
492 <        if (!(tmTop->flags & TM_F_LINEAR)) {
493 <                free((MEM_PTR)histo);
433 <                free((MEM_PTR)cumf);
434 <        }
492 >        free((MEM_PTR)histo);           /* clean up and return */
493 >        free((MEM_PTR)cumf);
494          returnOK;
495 + linearmap:                              /* linear tone-mapping */
496 +        if (tms->flags & TM_F_HCONTR)
497 +                d = htcontrs(Ldavg) / htcontrs(Lwavg);
498 +        else
499 +                d = Ldavg / Lwavg;
500 +        return(tmFixedMapping(tms, tms->inpsf*d/Ldmax, gamval));
501   }
502  
503  
504   int
505 < tmMapPixels(ps, ls, cs, len)
506 < register BYTE   *ps;
507 < TMbright        *ls;
508 < register BYTE   *cs;
509 < int     len;
505 > tmMapPixels(                    /* apply tone-mapping to pixel(s) */
506 > TMstruct        *tms,
507 > BYTE    *ps,
508 > TMbright        *ls,
509 > BYTE    *cs,
510 > int     len
511 > )
512   {
513          static char     funcName[] = "tmMapPixels";
514 <        register int4   li, pv;
514 >        int32   li, pv;
515  
516 <        if (tmTop == NULL || tmTop->lumap == NULL)
516 >        if (tms == NULL || tms->lumap == NULL)
517                  returnErr(TM_E_TMINVAL);
518 <        if (ps == NULL | ls == NULL | len <= 0)
518 >        if ((ps == NULL) | (ls == NULL) | (len < 0))
519                  returnErr(TM_E_ILLEGAL);
520          while (len--) {
521 <                if ((li = *ls++) < tmTop->mbrmin)
522 <                        li = tmTop->mbrmin;
523 <                else if (li > tmTop->mbrmax)
524 <                        li = tmTop->mbrmax;
525 <                li = tmTop->lumap[li - tmTop->mbrmin];
521 >                if ((li = *ls++) < tms->mbrmin) {
522 >                        li = 0;
523 >                } else {
524 >                        if (li > tms->mbrmax)
525 >                                li = tms->mbrmax;
526 >                        li = tms->lumap[li - tms->mbrmin];
527 >                }
528                  if (cs == TM_NOCHROM)
529                          *ps++ = li>255 ? 255 : li;
530                  else {
531 <                        pv = *cs++ * li / tmTop->cdiv[RED];
531 >                        pv = *cs++ * li / tms->cdiv[RED];
532                          *ps++ = pv>255 ? 255 : pv;
533 <                        pv = *cs++ * li / tmTop->cdiv[GRN];
533 >                        pv = *cs++ * li / tms->cdiv[GRN];
534                          *ps++ = pv>255 ? 255 : pv;
535 <                        pv = *cs++ * li / tmTop->cdiv[BLU];
535 >                        pv = *cs++ * li / tms->cdiv[BLU];
536                          *ps++ = pv>255 ? 255 : pv;
537                  }
538          }
# Line 471 | Line 540 | int    len;
540   }
541  
542  
474 struct tmStruct *
475 tmPop()                         /* pop top tone mapping off stack */
476 {
477        register struct tmStruct        *tms;
543  
479        if ((tms = tmTop) != NULL)
480                tmTop = tms->tmprev;
481        return(tms);
482 }
544  
545 <
546 < int
547 < tmPull(tms)                     /* pull a tone mapping from stack */
548 < register struct tmStruct        *tms;
545 > TMstruct *
546 > tmDup(                          /* duplicate top tone mapping */
547 > TMstruct        *tms
548 > )
549   {
489        register struct tmStruct        *tms2;
490                                        /* special cases first */
491        if (tms == NULL | tmTop == NULL)
492                return(0);
493        if (tms == tmTop) {
494                tmTop = tms->tmprev;
495                tms->tmprev = NULL;
496                return(1);
497        }
498        for (tms2 = tmTop; tms2->tmprev != NULL; tms2 = tms2->tmprev)
499                if (tms == tms2->tmprev) {      /* remove it */
500                        tms2->tmprev = tms->tmprev;
501                        tms->tmprev = NULL;
502                        return(1);
503                }
504        return(0);                      /* not found on stack */
505 }
506
507
508 struct tmStruct *
509 tmDup()                         /* duplicate top tone mapping */
510 {
550          int     len;
551 <        register int    i;
552 <        register struct tmStruct        *tmnew;
551 >        int     i;
552 >        TMstruct        *tmnew;
553  
554 <        if (tmTop == NULL)              /* anything to duplicate? */
554 >        if (tms == NULL)                /* anything to duplicate? */
555                  return(NULL);
556 <        tmnew = (struct tmStruct *)malloc(sizeof(struct tmStruct));
556 >        tmnew = (TMstruct *)malloc(sizeof(TMstruct));
557          if (tmnew == NULL)
558                  return(NULL);
559 <        *tmnew = *tmTop;                /* copy everything */
559 >        *tmnew = *tms;          /* copy everything */
560          if (tmnew->histo != NULL) {     /* duplicate histogram */
561                  len = (tmnew->hbrmax-MINBRT)/HISTEP + 1 -
562                                  (tmnew->hbrmin-MINBRT)/HISTEP;
563                  tmnew->histo = (int *)malloc(len*sizeof(int));
564                  if (tmnew->histo != NULL)
565                          for (i = len; i--; )
566 <                                tmnew->histo[i] = tmTop->histo[i];
566 >                                tmnew->histo[i] = tms->histo[i];
567          }
568          if (tmnew->lumap != NULL) {     /* duplicate luminance mapping */
569                  len = tmnew->mbrmax-tmnew->mbrmin+1;
# Line 532 | Line 571 | tmDup()                                /* duplicate top tone mapping */
571                                                  len*sizeof(unsigned short) );
572                  if (tmnew->lumap != NULL)
573                          for (i = len; i--; )
574 <                                tmnew->lumap[i] = tmTop->lumap[i];
574 >                                tmnew->lumap[i] = tms->lumap[i];
575          }
576                                          /* clear package data */
577          for (i = tmNumPkgs; i--; )
578                  tmnew->pd[i] = NULL;
579 <        tmnew->tmprev = tmTop;          /* make copy current */
580 <        return(tmTop = tmnew);
579 >                                        /* return copy */
580 >        return(tmnew);
581   }
582  
583  
545 int
546 tmPush(tms)                     /* push tone mapping on top of stack */
547 register struct tmStruct        *tms;
548 {
549        static char     funcName[] = "tmPush";
550                                        /* check validity */
551        if (tms == NULL)
552                returnErr(TM_E_ILLEGAL);
553        if (tms == tmTop)               /* check necessity */
554                returnOK;
555                                        /* pull if already in stack */
556        (void)tmPull(tms);
557                                        /* push it on top */
558        tms->tmprev = tmTop;
559        tmTop = tms;
560        returnOK;
561 }
562
563
584   void
585   tmDone(tms)                     /* done with tone mapping -- destroy it */
586 < register struct tmStruct        *tms;
586 > TMstruct        *tms;
587   {
588 <        register int    i;
589 <                                        /* NULL arg. is equiv. to tmTop */
590 <        if (tms == NULL && (tms = tmTop) == NULL)
588 >        int     i;
589 >                                        /* NULL arg. is equiv. to tms */
590 >        if (tms == NULL)
591                  return;
572                                        /* take out of stack if present */
573        (void)tmPull(tms);
592                                          /* free tables */
593          if (tms->histo != NULL)
594                  free((MEM_PTR)tms->histo);
# Line 581 | Line 599 | register struct tmStruct       *tms;
599                  if (tms->pd[i] != NULL)
600                          (*tmPkg[i]->Free)(tms->pd[i]);
601          free((MEM_PTR)tms);             /* free basic structure */
602 + }
603 +
604 + /******************** Shared but Private library routines *********************/
605 +
606 + BYTE    tmMesofact[BMESUPPER-BMESLOWER];
607 +
608 + void
609 + tmMkMesofact()                          /* build mesopic lookup factor table */
610 + {
611 +        int     i;
612 +
613 +        if (tmMesofact[BMESUPPER-BMESLOWER-1])
614 +                return;
615 +
616 +        for (i = BMESLOWER; i < BMESUPPER; i++)
617 +                tmMesofact[i-BMESLOWER] = 256. *
618 +                                (tmLuminance(i) - LMESLOWER) /
619 +                                (LMESUPPER - LMESLOWER);
620 + }
621 +
622 +
623 + int
624 + tmErrorReturn(                          /* error return (with message) */
625 + char    *func,
626 + TMstruct        *tms,
627 + int     err
628 + )
629 + {
630 +        tmLastFunction = func;
631 +        tmLastError = err;
632 +        if (tms != NULL && tms->flags & TM_F_NOSTDERR)
633 +                return(err);
634 +        fputs(func, stderr);
635 +        fputs(": ", stderr);
636 +        fputs(tmErrorMessage[err], stderr);
637 +        fputs("!\n", stderr);
638 +        return(err);
639   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines