ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/tonemap.c
Revision: 3.50
Committed: Tue Apr 13 02:42:31 2021 UTC (3 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 3.49: +3 -5 lines
Log Message:
Added better termination criterion for linear tone curve

File Contents

# User Rev Content
1 greg 3.1 #ifndef lint
2 greg 3.50 static const char RCSid[] = "$Id: tonemap.c,v 3.49 2021/04/12 22:16:02 greg Exp $";
3 greg 3.1 #endif
4     /*
5     * Tone mapping functions.
6     * See tonemap.h for detailed function descriptions.
7 greg 3.9 * Added von Kries white-balance calculations 10/01 (GW).
8     *
9     * Externals declared in tonemap.h
10     */
11    
12 greg 3.10 #include "copyright.h"
13 greg 3.1
14     #include <stdio.h>
15 greg 3.34 #include <stdlib.h>
16 greg 3.1 #include <math.h>
17 greg 3.38 #include <string.h>
18 greg 3.1 #include "tmprivat.h"
19     #include "tmerrmsg.h"
20    
21     #define exp10(x) exp(M_LN10*(x))
22    
23 greg 3.4 /* our list of conversion packages */
24     struct tmPackage *tmPkg[TM_MAXPKG];
25     int tmNumPkgs = 0; /* number of registered packages */
26 greg 3.1
27 greg 3.23 /* luminance->brightness lookup */
28 greg 3.24 static TMbright *tmFloat2BrtLUT = NULL;
29 greg 3.23
30     #define tmCvLumLUfp(pf) tmFloat2BrtLUT[*(int32 *)(pf) >> 15]
31    
32 gregl 3.5
33 greg 3.16 TMstruct *
34     tmInit( /* initialize new tone mapping */
35     int flags,
36     RGBPRIMP monpri,
37     double gamval
38     )
39 greg 3.1 {
40     COLORMAT cmat;
41 greg 3.16 TMstruct *tmnew;
42     int i;
43 greg 3.1 /* allocate structure */
44 greg 3.16 tmnew = (TMstruct *)malloc(sizeof(TMstruct));
45 greg 3.1 if (tmnew == NULL)
46     return(NULL);
47    
48     tmnew->flags = flags & ~TM_F_UNIMPL;
49 greg 3.15 if (tmnew->flags & TM_F_BW)
50     tmnew->flags &= ~TM_F_MESOPIC;
51 greg 3.1 /* set monitor transform */
52     if (monpri == NULL || monpri == stdprims || tmnew->flags & TM_F_BW) {
53     tmnew->monpri = stdprims;
54     tmnew->clf[RED] = rgb2xyzmat[1][0];
55     tmnew->clf[GRN] = rgb2xyzmat[1][1];
56     tmnew->clf[BLU] = rgb2xyzmat[1][2];
57     } else {
58 greg 3.26 comprgb2xyzmat(cmat, tmnew->monpri=monpri);
59 greg 3.1 tmnew->clf[RED] = cmat[1][0];
60     tmnew->clf[GRN] = cmat[1][1];
61     tmnew->clf[BLU] = cmat[1][2];
62     }
63     /* set gamma value */
64     if (gamval < MINGAM)
65     tmnew->mongam = DEFGAM;
66     else
67     tmnew->mongam = gamval;
68 greg 3.4 /* set color divisors */
69     for (i = 0; i < 3; i++)
70 greg 3.40 tmnew->cdiv[i] = TM_BRES*pow(tmnew->clf[i], 1./tmnew->mongam);
71 greg 3.4
72 greg 3.1 /* set input transform */
73     tmnew->inppri = tmnew->monpri;
74     tmnew->cmat[0][0] = tmnew->cmat[1][1] = tmnew->cmat[2][2] =
75     tmnew->inpsf = WHTEFFICACY;
76     tmnew->cmat[0][1] = tmnew->cmat[0][2] = tmnew->cmat[1][0] =
77     tmnew->cmat[1][2] = tmnew->cmat[2][0] = tmnew->cmat[2][1] = 0.;
78 greg 3.18 tmnew->inpdat = NULL;
79 greg 3.9 tmnew->hbrmin = 10; tmnew->hbrmax = -10;
80 greg 3.1 tmnew->histo = NULL;
81 greg 3.9 tmnew->mbrmin = 10; tmnew->mbrmax = -10;
82 greg 3.1 tmnew->lumap = NULL;
83 greg 3.4 /* zero private data */
84     for (i = TM_MAXPKG; i--; )
85     tmnew->pd[i] = NULL;
86 greg 3.17 tmnew->lastError = TM_E_OK;
87     tmnew->lastFunc = "NoErr";
88 greg 3.16 /* return new TMstruct */
89     return(tmnew);
90 greg 3.1 }
91    
92    
93     int
94 greg 3.16 tmSetSpace( /* set input color space for conversions */
95     TMstruct *tms,
96     RGBPRIMP pri,
97 greg 3.18 double sf,
98     MEM_PTR dat
99 greg 3.16 )
100 greg 3.1 {
101 greg 3.17 static const char funcName[] = "tmSetSpace";
102 greg 3.16 int i, j;
103 greg 3.1 /* error check */
104 greg 3.16 if (tms == NULL)
105 greg 3.1 returnErr(TM_E_TMINVAL);
106     if (sf <= 1e-12)
107     returnErr(TM_E_ILLEGAL);
108     /* check if no change */
109 greg 3.18 if (pri == tms->inppri && FEQ(sf, tms->inpsf) && dat == tms->inpdat)
110 greg 3.1 returnOK;
111 greg 3.16 tms->inppri = pri; /* let's set it */
112     tms->inpsf = sf;
113 greg 3.18 tms->inpdat = dat;
114 greg 3.1
115 greg 3.16 if (tms->flags & TM_F_BW) { /* color doesn't matter */
116     tms->monpri = tms->inppri; /* eliminate xform */
117     if (tms->inppri == TM_XYZPRIM) {
118     tms->clf[CIEX] = tms->clf[CIEZ] = 0.;
119     tms->clf[CIEY] = 1.;
120 greg 3.1 } else {
121 greg 3.26 comprgb2xyzmat(tms->cmat, tms->monpri);
122 greg 3.16 tms->clf[RED] = tms->cmat[1][0];
123     tms->clf[GRN] = tms->cmat[1][1];
124     tms->clf[BLU] = tms->cmat[1][2];
125     }
126 greg 3.46 tms->cmat[0][0] = tms->cmat[1][1] = tms->cmat[2][2] = 1.;
127 greg 3.16 tms->cmat[0][1] = tms->cmat[0][2] = tms->cmat[1][0] =
128     tms->cmat[1][2] = tms->cmat[2][0] = tms->cmat[2][1] = 0.;
129 greg 3.1
130 greg 3.46 } else if (tms->inppri == TM_XYZPRIM) { /* input is XYZ */
131 greg 3.16 compxyz2rgbWBmat(tms->cmat, tms->monpri);
132 greg 3.1
133 greg 3.46 } else { /* input is RGB */
134 greg 3.16 if (tms->inppri != tms->monpri &&
135     PRIMEQ(tms->inppri, tms->monpri))
136     tms->inppri = tms->monpri; /* no xform */
137 greg 3.36 if (!comprgb2rgbWBmat(tms->cmat, tms->inppri, tms->monpri))
138     returnErr(TM_E_ILLEGAL);
139 greg 3.1 }
140     for (i = 0; i < 3; i++)
141     for (j = 0; j < 3; j++)
142 greg 3.16 tms->cmat[i][j] *= tms->inpsf;
143 greg 3.4 /* set color divisors */
144     for (i = 0; i < 3; i++)
145 greg 3.40 tms->cdiv[i] = TM_BRES*pow(tms->clf[i] < .001 ? .001 :
146 greg 3.36 tms->clf[i], 1./tms->mongam);
147 greg 3.4 /* notify packages */
148     for (i = tmNumPkgs; i--; )
149 greg 3.16 if (tms->pd[i] != NULL && tmPkg[i]->NewSpace != NULL)
150     (*tmPkg[i]->NewSpace)(tms);
151 greg 3.1 returnOK;
152     }
153    
154    
155     void
156 greg 3.23 tmClearHisto( /* clear current histogram */
157 greg 3.16 TMstruct *tms
158     )
159 greg 3.1 {
160 greg 3.16 if (tms == NULL || tms->histo == NULL)
161 greg 3.1 return;
162 greg 3.46 free(tms->histo);
163 greg 3.16 tms->histo = NULL;
164 greg 3.1 }
165    
166    
167 greg 3.23 TMbright
168     tmCvLuminance( /* convert a single luminance */
169     double lum
170     )
171     {
172     double d;
173    
174     #ifdef isfinite
175     if (!isfinite(lum) || lum <= TM_NOLUM)
176     #else
177     if (lum <= TM_NOLUM)
178     #endif
179     return(TM_NOBRT);
180     d = TM_BRTSCALE*log(lum);
181 greg 3.35 return((TMbright)(d + .5 - (d < 0.)));
182 greg 3.23 }
183    
184    
185     int
186     tmCvLums( /* convert luminances using lookup */
187     TMbright *ls,
188     float *scan,
189     int len
190     )
191     {
192     if (tmFloat2BrtLUT == NULL) { /* initialize lookup table */
193     int32 i;
194     tmFloat2BrtLUT = (TMbright *)malloc(sizeof(TMbright)*0x10000);
195     if (tmFloat2BrtLUT == NULL)
196     return(TM_E_NOMEM);
197     for (i = 0; i < 0x10000; i++) {
198     int32 l = (i<<1 | 1) << 14;
199     #ifndef isfinite
200     if ((l & 0x7f800000) == 0x7f800000)
201     tmFloat2BrtLUT[i] = TM_NOBRT;
202     else
203     #endif
204     tmFloat2BrtLUT[i] = tmCvLuminance(*(float *)&l);
205     }
206     }
207     if (len <= 0)
208     return(TM_E_OK);
209     if ((ls == NULL) | (scan == NULL))
210     return(TM_E_ILLEGAL);
211     while (len--) {
212     if (*scan <= TM_NOLUM) {
213     *ls++ = TM_NOBRT;
214     ++scan;
215     continue;
216     }
217     *ls++ = tmCvLumLUfp(scan++);
218     }
219     return(TM_E_OK);
220     }
221    
222    
223     int
224     tmCvGrays( /* convert float gray values */
225     TMstruct *tms,
226     TMbright *ls,
227     float *scan,
228     int len
229     )
230     {
231     static const char funcName[] = "tmCvGrays";
232     int i;
233    
234     if (tms == NULL)
235     returnErr(TM_E_TMINVAL);
236     if ((ls == NULL) | (scan == NULL) | (len < 0))
237     returnErr(TM_E_ILLEGAL);
238 greg 3.24 if (tmFloat2BrtLUT == NULL) /* initialize */
239     tmCvLums(NULL, NULL, 0);
240 greg 3.23 for (i = len; i--; ) {
241     float lum = tms->inpsf * scan[i];
242     if (lum <= TM_NOLUM)
243     ls[i] = TM_NOBRT;
244     else
245     ls[i] = tmCvLumLUfp(&lum);
246     }
247     returnOK;
248     }
249    
250    
251 greg 3.1 int
252 greg 3.16 tmCvColors( /* convert float colors */
253     TMstruct *tms,
254     TMbright *ls,
255 greg 3.37 uby8 *cs,
256 greg 3.16 COLOR *scan,
257     int len
258     )
259 greg 3.1 {
260 greg 3.17 static const char funcName[] = "tmCvColors";
261 greg 3.37 static uby8 gamtab[1024];
262 greg 3.20 static double curgam = .0;
263 greg 3.1 COLOR cmon;
264 greg 3.23 float lum, slum, d;
265 greg 3.16 int i;
266 greg 3.1
267 greg 3.16 if (tms == NULL)
268 greg 3.1 returnErr(TM_E_TMINVAL);
269 schorsch 3.13 if ((ls == NULL) | (scan == NULL) | (len < 0))
270 greg 3.1 returnErr(TM_E_ILLEGAL);
271 greg 3.24 if (tmFloat2BrtLUT == NULL) /* initialize */
272     tmCvLums(NULL, NULL, 0);
273 greg 3.21 if (cs != TM_NOCHROM && fabs(tms->mongam - curgam) > .02) {
274 greg 3.20 curgam = tms->mongam; /* (re)build table */
275     for (i = 1024; i--; )
276     gamtab[i] = (int)(256.*pow((i+.5)/1024., 1./curgam));
277     }
278 greg 3.1 for (i = len; i--; ) {
279 greg 3.16 if (tmNeedMatrix(tms)) { /* get monitor RGB */
280     colortrans(cmon, tms->cmat, scan[i]);
281 greg 3.15 } else {
282 greg 3.16 cmon[RED] = tms->inpsf*scan[i][RED];
283     cmon[GRN] = tms->inpsf*scan[i][GRN];
284     cmon[BLU] = tms->inpsf*scan[i][BLU];
285 greg 3.1 }
286 greg 3.23 #ifdef isfinite
287     if (!isfinite(cmon[RED]) || cmon[RED] < .0f) cmon[RED] = .0f;
288     if (!isfinite(cmon[GRN]) || cmon[GRN] < .0f) cmon[GRN] = .0f;
289     if (!isfinite(cmon[BLU]) || cmon[BLU] < .0f) cmon[BLU] = .0f;
290     #else
291     if (cmon[RED] < .0f) cmon[RED] = .0f;
292     if (cmon[GRN] < .0f) cmon[GRN] = .0f;
293     if (cmon[BLU] < .0f) cmon[BLU] = .0f;
294 greg 3.22 #endif
295 greg 3.1 /* world luminance */
296 greg 3.16 lum = tms->clf[RED]*cmon[RED] +
297     tms->clf[GRN]*cmon[GRN] +
298     tms->clf[BLU]*cmon[BLU] ;
299 greg 3.25 if (lum <= TM_NOLUM) { /* convert brightness */
300     lum = cmon[RED] = cmon[GRN] = cmon[BLU] = TM_NOLUM;
301 greg 3.23 ls[i] = TM_NOBRT;
302 greg 3.25 } else
303 greg 3.23 ls[i] = tmCvLumLUfp(&lum);
304 greg 3.1 if (cs == TM_NOCHROM) /* no color? */
305     continue;
306 greg 3.16 if (tms->flags & TM_F_MESOPIC && lum < LMESUPPER) {
307 greg 3.1 slum = scotlum(cmon); /* mesopic adj. */
308 greg 3.25 if (lum < LMESLOWER) {
309 greg 3.1 cmon[RED] = cmon[GRN] = cmon[BLU] = slum;
310 greg 3.25 } else {
311 greg 3.1 d = (lum - LMESLOWER)/(LMESUPPER - LMESLOWER);
312 greg 3.16 if (tms->flags & TM_F_BW)
313 greg 3.1 cmon[RED] = cmon[GRN] =
314     cmon[BLU] = d*lum;
315     else
316     scalecolor(cmon, d);
317 greg 3.23 d = (1.f-d)*slum;
318 greg 3.1 cmon[RED] += d;
319     cmon[GRN] += d;
320     cmon[BLU] += d;
321     }
322 greg 3.16 } else if (tms->flags & TM_F_BW) {
323 greg 3.1 cmon[RED] = cmon[GRN] = cmon[BLU] = lum;
324     }
325 greg 3.16 d = tms->clf[RED]*cmon[RED]/lum;
326 greg 3.23 cs[3*i ] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
327 greg 3.16 d = tms->clf[GRN]*cmon[GRN]/lum;
328 greg 3.23 cs[3*i+1] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
329 greg 3.16 d = tms->clf[BLU]*cmon[BLU]/lum;
330 greg 3.23 cs[3*i+2] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
331 greg 3.1 }
332     returnOK;
333     }
334    
335    
336 greg 3.9 int
337 greg 3.16 tmAddHisto( /* add values to histogram */
338     TMstruct *tms,
339     TMbright *ls,
340     int len,
341     int wt
342     )
343 greg 3.1 {
344 greg 3.17 static const char funcName[] = "tmAddHisto";
345 greg 3.9 int oldorig=0, oldlen, horig, hlen;
346 greg 3.16 int i, j;
347 greg 3.1
348 greg 3.16 if (tms == NULL)
349 greg 3.1 returnErr(TM_E_TMINVAL);
350 gwlarson 3.8 if (len < 0)
351     returnErr(TM_E_ILLEGAL);
352     if (len == 0)
353     returnOK;
354 greg 3.1 /* first, grow limits */
355 greg 3.16 if (tms->histo == NULL) {
356 greg 3.1 for (i = len; i-- && ls[i] < MINBRT; )
357     ;
358     if (i < 0)
359     returnOK;
360 greg 3.16 tms->hbrmin = tms->hbrmax = ls[i];
361 greg 3.1 oldlen = 0;
362     } else {
363 greg 3.31 oldorig = HISTI(tms->hbrmin);
364     oldlen = HISTI(tms->hbrmax) + 1 - oldorig;
365 greg 3.1 }
366     for (i = len; i--; ) {
367     if ((j = ls[i]) < MINBRT)
368     continue;
369 greg 3.16 if (j < tms->hbrmin)
370     tms->hbrmin = j;
371     else if (j > tms->hbrmax)
372     tms->hbrmax = j;
373 greg 3.1 }
374 greg 3.31 horig = HISTI(tms->hbrmin);
375     hlen = HISTI(tms->hbrmax) + 1 - horig;
376 greg 3.1 if (hlen > oldlen) { /* (re)allocate histogram */
377 greg 3.38 HIST_TYP *newhist = (HIST_TYP *)calloc(hlen, sizeof(HIST_TYP));
378 greg 3.1 if (newhist == NULL)
379     returnErr(TM_E_NOMEM);
380     if (oldlen) { /* copy and free old */
381 greg 3.46 memcpy(newhist+(oldorig-horig),
382     tms->histo, oldlen*sizeof(HIST_TYP));
383     free(tms->histo);
384 greg 3.1 }
385 greg 3.16 tms->histo = newhist;
386 greg 3.1 }
387     if (wt == 0)
388     returnOK;
389     for (i = len; i--; ) /* add in new counts */
390     if (ls[i] >= MINBRT)
391 greg 3.31 tms->histo[ HISTI(ls[i]) - horig ] += wt;
392 greg 3.1 returnOK;
393     }
394    
395    
396     static double
397 greg 3.16 htcontrs( /* human threshold contrast sensitivity, dL(La) */
398     double La
399     )
400 greg 3.1 {
401     double l10La, l10dL;
402     /* formula taken from Ferwerda et al. [SG96] */
403     if (La < 1.148e-4)
404     return(1.38e-3);
405     l10La = log10(La);
406     if (l10La < -1.44) /* rod response regime */
407     l10dL = pow(.405*l10La + 1.6, 2.18) - 2.86;
408     else if (l10La < -.0184)
409     l10dL = l10La - .395;
410     else if (l10La < 1.9) /* cone response regime */
411     l10dL = pow(.249*l10La + .65, 2.7) - .72;
412     else
413     l10dL = l10La - 1.255;
414    
415     return(exp10(l10dL));
416     }
417    
418    
419 greg 3.9 int
420 greg 3.16 tmFixedMapping( /* compute fixed, linear tone-mapping */
421     TMstruct *tms,
422     double expmult,
423 greg 3.46 double gamval,
424     double Lddyn
425 greg 3.16 )
426 greg 3.9 {
427 greg 3.17 static const char funcName[] = "tmFixedMapping";
428 greg 3.47 const int maxV = (1L<<(8*sizeof(TMAP_TYP))) - 1;
429     double minD;
430     int i;
431 greg 3.9
432 greg 3.16 if (!tmNewMap(tms))
433 greg 3.9 returnErr(TM_E_NOMEM);
434 greg 3.46 /* check arguments */
435     if (expmult <= .0) expmult = 1.;
436     if (gamval < MINGAM) gamval = tms->mongam;
437     if (Lddyn < MINLDDYN) Lddyn = DEFLDDYN;
438     minD = 1./Lddyn;
439 greg 3.28 for (i = tms->mbrmax-tms->mbrmin+1; i--; ) {
440 greg 3.46 double d;
441     d = expmult/tms->inpsf * tmLuminance(tms->mbrmin + i);
442 greg 3.48 if (d >= 2.*minD)
443     d -= minD;
444     else /* soft black crushing */
445     d *= d/(4.*minD);
446     d /= 1. - minD;
447 greg 3.46 d = TM_BRES*pow(d, 1./gamval);
448 greg 3.47 tms->lumap[i] = (d > maxV) ? maxV : (int)d;
449 greg 3.28 }
450 greg 3.9 returnOK;
451     }
452    
453    
454 greg 3.1 int
455 greg 3.16 tmComputeMapping( /* compute histogram tone-mapping */
456     TMstruct *tms,
457     double gamval,
458     double Lddyn,
459     double Ldmax
460     )
461 greg 3.1 {
462 greg 3.17 static const char funcName[] = "tmComputeMapping";
463 greg 3.38 HIST_TYP *histo;
464 greg 3.1 float *cumf;
465 greg 3.38 int brt0, histlen;
466 greg 3.39 HIST_TYP threshold, ceiling, trimmings, histot;
467 greg 3.44 double logLddyn, Ldmin, Lwavg, Tr, Lw, Ld;
468 greg 3.9 double sum;
469 greg 3.16 double d;
470     int i, j;
471 greg 3.1
472 greg 3.16 if (tms == NULL || tms->histo == NULL)
473 greg 3.1 returnErr(TM_E_TMINVAL);
474     /* check arguments */
475     if (Lddyn < MINLDDYN) Lddyn = DEFLDDYN;
476     if (Ldmax < MINLDMAX) Ldmax = DEFLDMAX;
477 greg 3.16 if (gamval < MINGAM) gamval = tms->mongam;
478 greg 3.1 /* compute handy values */
479     Ldmin = Ldmax/Lddyn;
480     logLddyn = log(Lddyn);
481 greg 3.31 i = HISTI(tms->hbrmin);
482 greg 3.32 brt0 = HISTV(i);
483 greg 3.31 histlen = HISTI(tms->hbrmax) + 1 - i;
484 greg 3.1 /* histogram total and mean */
485     histot = 0; sum = 0;
486     j = brt0 + histlen*HISTEP;
487     for (i = histlen; i--; ) {
488 greg 3.16 histot += tms->histo[i];
489 greg 3.30 sum += (double)(j -= HISTEP) * tms->histo[i];
490 greg 3.1 }
491 greg 3.29 if (!histot)
492 greg 3.1 returnErr(TM_E_TMFAIL);
493 greg 3.50 threshold = histot/40 + 1;
494 greg 3.1 Lwavg = tmLuminance( (double)sum / histot );
495 greg 3.29 /* use linear tone mapping? */
496 greg 3.30 if (tms->flags & TM_F_LINEAR || threshold < 4 ||
497     tms->hbrmax - tms->hbrmin < TM_BRTSCALE*logLddyn)
498 greg 3.29 goto linearmap;
499 greg 3.9 /* clamp histogram */
500 greg 3.38 histo = (HIST_TYP *)malloc(histlen*sizeof(HIST_TYP));
501 greg 3.9 cumf = (float *)malloc((histlen+2)*sizeof(float));
502 schorsch 3.13 if ((histo == NULL) | (cumf == NULL))
503 greg 3.9 returnErr(TM_E_NOMEM);
504     cumf[histlen+1] = 1.; /* guard for assignment code */
505 greg 3.42 /* make malleable copy */
506 greg 3.46 memcpy(histo, tms->histo, histlen*sizeof(HIST_TYP));
507 greg 3.9 do { /* iterate to solution */
508     sum = 0; /* cumulative probability */
509     for (i = 0; i < histlen; i++) {
510 greg 3.45 cumf[i] = sum/histot;
511     sum += (double)histo[i];
512 greg 3.9 }
513     cumf[histlen] = 1.;
514 greg 3.16 Tr = histot * (double)(tms->hbrmax - tms->hbrmin) /
515 greg 3.45 ((double)TM_BRTSCALE*histlen*logLddyn);
516 greg 3.9 ceiling = Tr + 1.;
517     trimmings = 0; /* clip to envelope */
518     for (i = histlen; i--; ) {
519 greg 3.16 if (tms->flags & TM_F_HCONTR) {
520 greg 3.9 Lw = tmLuminance(brt0 + i*HISTEP);
521     Ld = Ldmin * exp( logLddyn *
522     .5*(cumf[i]+cumf[i+1]) );
523     ceiling = Tr * (htcontrs(Ld) * Lw) /
524     (htcontrs(Lw) * Ld) + 1.;
525 greg 3.2 }
526 greg 3.9 if (histo[i] > ceiling) {
527     trimmings += histo[i] - ceiling;
528     histo[i] = ceiling;
529 greg 3.1 }
530     }
531 greg 3.9 /* check if we're out of data */
532     if ((histot -= trimmings) <= threshold) {
533 greg 3.46 free(histo);
534     free(cumf);
535 greg 3.9 goto linearmap;
536     }
537 greg 3.50 } while (40*trimmings > histot);
538 greg 3.30 /* allocate space for mapping */
539     if (!tmNewMap(tms))
540     returnErr(TM_E_NOMEM);
541 greg 3.9 /* assign tone-mapping */
542 greg 3.16 for (i = tms->mbrmax-tms->mbrmin+1; i--; ) {
543     j = d = (double)i/(tms->mbrmax-tms->mbrmin)*histlen;
544 greg 3.9 d -= (double)j;
545     Ld = Ldmin*exp(logLddyn*((1.-d)*cumf[j]+d*cumf[j+1]));
546     d = (Ld - Ldmin)/(Ldmax - Ldmin);
547 greg 3.40 tms->lumap[i] = TM_BRES*pow(d, 1./gamval);
548 greg 3.1 }
549 greg 3.46 free(histo); /* clean up and return */
550     free(cumf);
551 greg 3.1 returnOK;
552 greg 3.9 linearmap: /* linear tone-mapping */
553 greg 3.16 if (tms->flags & TM_F_HCONTR)
554 greg 3.44 d = htcontrs(sqrt(Ldmax*Ldmin)) / htcontrs(Lwavg);
555 greg 3.9 else
556 greg 3.44 d = Ldmax / tmLuminance(tms->hbrmax);
557 greg 3.46 return(tmFixedMapping(tms, tms->inpsf*d/Ldmax, gamval, Lddyn));
558 greg 3.1 }
559    
560    
561     int
562 greg 3.16 tmMapPixels( /* apply tone-mapping to pixel(s) */
563     TMstruct *tms,
564 greg 3.37 uby8 *ps,
565 greg 3.16 TMbright *ls,
566 greg 3.37 uby8 *cs,
567 greg 3.16 int len
568     )
569 greg 3.1 {
570 greg 3.17 static const char funcName[] = "tmMapPixels";
571 greg 3.41 TMbright lv;
572     TMAP_TYP li;
573     int pv;
574 greg 3.1
575 greg 3.16 if (tms == NULL || tms->lumap == NULL)
576 greg 3.1 returnErr(TM_E_TMINVAL);
577 schorsch 3.13 if ((ps == NULL) | (ls == NULL) | (len < 0))
578 greg 3.1 returnErr(TM_E_ILLEGAL);
579     while (len--) {
580 greg 3.41 if ((lv = *ls++) < tms->mbrmin) {
581 greg 3.9 li = 0;
582     } else {
583 greg 3.41 if (lv > tms->mbrmax)
584     lv = tms->mbrmax;
585     li = tms->lumap[lv - tms->mbrmin];
586 greg 3.9 }
587 greg 3.43 if (cs == TM_NOCHROM) {
588     #if !(TM_BRES & 0xff)
589     *ps++ = li>=TM_BRES ? 255 : li/(TM_BRES>>8);
590     #else
591     *ps++ = li>=TM_BRES ? 255 : (li<<8)/TM_BRES;
592     #endif
593     } else {
594 greg 3.16 pv = *cs++ * li / tms->cdiv[RED];
595 greg 3.1 *ps++ = pv>255 ? 255 : pv;
596 greg 3.16 pv = *cs++ * li / tms->cdiv[GRN];
597 greg 3.1 *ps++ = pv>255 ? 255 : pv;
598 greg 3.16 pv = *cs++ * li / tms->cdiv[BLU];
599 greg 3.1 *ps++ = pv>255 ? 255 : pv;
600     }
601     }
602     returnOK;
603     }
604    
605    
606 greg 3.16 TMstruct *
607 greg 3.42 tmDup( /* duplicate tone mapping */
608 greg 3.16 TMstruct *tms
609     )
610 greg 3.3 {
611     int len;
612 greg 3.16 int i;
613     TMstruct *tmnew;
614 greg 3.3
615 greg 3.16 if (tms == NULL) /* anything to duplicate? */
616 greg 3.3 return(NULL);
617 greg 3.16 tmnew = (TMstruct *)malloc(sizeof(TMstruct));
618 greg 3.3 if (tmnew == NULL)
619     return(NULL);
620 greg 3.16 *tmnew = *tms; /* copy everything */
621 greg 3.3 if (tmnew->histo != NULL) { /* duplicate histogram */
622 greg 3.31 len = HISTI(tmnew->hbrmax) + 1 - HISTI(tmnew->hbrmin);
623 greg 3.38 tmnew->histo = (HIST_TYP *)malloc(len*sizeof(HIST_TYP));
624 greg 3.3 if (tmnew->histo != NULL)
625 greg 3.46 memcpy(tmnew->histo, tms->histo, len*sizeof(HIST_TYP));
626 greg 3.3 }
627     if (tmnew->lumap != NULL) { /* duplicate luminance mapping */
628 gregl 3.6 len = tmnew->mbrmax-tmnew->mbrmin+1;
629 greg 3.38 tmnew->lumap = (TMAP_TYP *)malloc(len*sizeof(TMAP_TYP));
630 greg 3.3 if (tmnew->lumap != NULL)
631 greg 3.46 memcpy(tmnew->lumap, tms->lumap, len*sizeof(TMAP_TYP));
632 greg 3.3 }
633 greg 3.4 /* clear package data */
634     for (i = tmNumPkgs; i--; )
635     tmnew->pd[i] = NULL;
636 greg 3.16 /* return copy */
637     return(tmnew);
638 greg 3.1 }
639    
640    
641     void
642     tmDone(tms) /* done with tone mapping -- destroy it */
643 greg 3.16 TMstruct *tms;
644 greg 3.1 {
645 greg 3.16 int i;
646     /* NULL arg. is equiv. to tms */
647     if (tms == NULL)
648 greg 3.1 return;
649     /* free tables */
650     if (tms->histo != NULL)
651 greg 3.46 free(tms->histo);
652 greg 3.1 if (tms->lumap != NULL)
653 greg 3.46 free(tms->lumap);
654 greg 3.4 /* free private data */
655     for (i = tmNumPkgs; i--; )
656     if (tms->pd[i] != NULL)
657     (*tmPkg[i]->Free)(tms->pd[i]);
658 greg 3.46 free(tms); /* free basic structure */
659 greg 3.9 }
660    
661     /******************** Shared but Private library routines *********************/
662    
663 greg 3.37 uby8 tmMesofact[BMESUPPER-BMESLOWER];
664 greg 3.9
665     void
666     tmMkMesofact() /* build mesopic lookup factor table */
667     {
668 greg 3.16 int i;
669 greg 3.9
670     if (tmMesofact[BMESUPPER-BMESLOWER-1])
671     return;
672    
673     for (i = BMESLOWER; i < BMESUPPER; i++)
674     tmMesofact[i-BMESLOWER] = 256. *
675     (tmLuminance(i) - LMESLOWER) /
676     (LMESUPPER - LMESLOWER);
677     }
678    
679    
680     int
681 greg 3.33 tmNewMap( /* allocate new tone-mapping array */
682     TMstruct *tms
683     )
684     {
685     if (tms->lumap != NULL && (tms->mbrmax - tms->mbrmin) !=
686     (tms->hbrmax - tms->hbrmin)) {
687 greg 3.46 free(tms->lumap);
688 greg 3.33 tms->lumap = NULL;
689     }
690     tms->mbrmin = tms->hbrmin;
691     tms->mbrmax = tms->hbrmax;
692     if (tms->mbrmin > tms->mbrmax)
693 greg 3.47 return(0);
694 greg 3.33 if (tms->lumap == NULL)
695 greg 3.46 tms->lumap = (TMAP_TYP *)calloc(tms->mbrmax-tms->mbrmin+1,
696     sizeof(TMAP_TYP));
697     else
698     memset(tms->lumap, 0, (tms->mbrmax-tms->mbrmin+1)*sizeof(TMAP_TYP));
699    
700 greg 3.33 return(tms->lumap != NULL);
701     }
702    
703    
704     int
705 greg 3.16 tmErrorReturn( /* error return (with message) */
706 greg 3.17 const char *func,
707 greg 3.16 TMstruct *tms,
708     int err
709     )
710 greg 3.9 {
711 greg 3.17 if (tms != NULL) {
712     tms->lastFunc = func;
713     tms->lastError = err;
714     if (tms->flags & TM_F_NOSTDERR)
715     return(err);
716     }
717 greg 3.9 fputs(func, stderr);
718     fputs(": ", stderr);
719     fputs(tmErrorMessage[err], stderr);
720     fputs("!\n", stderr);
721     return(err);
722 greg 3.1 }