ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/tonemap.c
Revision: 3.54
Committed: Sat Jan 15 16:57:46 2022 UTC (2 years, 3 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R4, HEAD
Changes since 3.53: +2 -2 lines
Log Message:
refactor: eliminated use of MEM_PTR and (char *) for untyped memory

File Contents

# User Rev Content
1 greg 3.1 #ifndef lint
2 greg 3.54 static const char RCSid[] = "$Id: tonemap.c,v 3.53 2022/01/07 23:01:01 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 greg 3.54 void *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.51 tms->hbrmin = 10; tms->hbrmax = -10;
164 greg 3.16 tms->histo = NULL;
165 greg 3.1 }
166    
167    
168 greg 3.23 TMbright
169     tmCvLuminance( /* convert a single luminance */
170     double lum
171     )
172     {
173     double d;
174    
175     #ifdef isfinite
176     if (!isfinite(lum) || lum <= TM_NOLUM)
177     #else
178     if (lum <= TM_NOLUM)
179     #endif
180     return(TM_NOBRT);
181     d = TM_BRTSCALE*log(lum);
182 greg 3.35 return((TMbright)(d + .5 - (d < 0.)));
183 greg 3.23 }
184    
185    
186     int
187     tmCvLums( /* convert luminances using lookup */
188     TMbright *ls,
189     float *scan,
190     int len
191     )
192     {
193     if (tmFloat2BrtLUT == NULL) { /* initialize lookup table */
194     int32 i;
195     tmFloat2BrtLUT = (TMbright *)malloc(sizeof(TMbright)*0x10000);
196     if (tmFloat2BrtLUT == NULL)
197     return(TM_E_NOMEM);
198     for (i = 0; i < 0x10000; i++) {
199     int32 l = (i<<1 | 1) << 14;
200     #ifndef isfinite
201     if ((l & 0x7f800000) == 0x7f800000)
202     tmFloat2BrtLUT[i] = TM_NOBRT;
203     else
204     #endif
205     tmFloat2BrtLUT[i] = tmCvLuminance(*(float *)&l);
206     }
207     }
208     if (len <= 0)
209     return(TM_E_OK);
210     if ((ls == NULL) | (scan == NULL))
211     return(TM_E_ILLEGAL);
212     while (len--) {
213     if (*scan <= TM_NOLUM) {
214     *ls++ = TM_NOBRT;
215     ++scan;
216     continue;
217     }
218     *ls++ = tmCvLumLUfp(scan++);
219     }
220     return(TM_E_OK);
221     }
222    
223    
224     int
225     tmCvGrays( /* convert float gray values */
226     TMstruct *tms,
227     TMbright *ls,
228     float *scan,
229     int len
230     )
231     {
232     static const char funcName[] = "tmCvGrays";
233     int i;
234    
235     if (tms == NULL)
236     returnErr(TM_E_TMINVAL);
237     if ((ls == NULL) | (scan == NULL) | (len < 0))
238     returnErr(TM_E_ILLEGAL);
239 greg 3.24 if (tmFloat2BrtLUT == NULL) /* initialize */
240     tmCvLums(NULL, NULL, 0);
241 greg 3.23 for (i = len; i--; ) {
242     float lum = tms->inpsf * scan[i];
243     if (lum <= TM_NOLUM)
244     ls[i] = TM_NOBRT;
245     else
246     ls[i] = tmCvLumLUfp(&lum);
247     }
248     returnOK;
249     }
250    
251    
252 greg 3.1 int
253 greg 3.16 tmCvColors( /* convert float colors */
254     TMstruct *tms,
255     TMbright *ls,
256 greg 3.37 uby8 *cs,
257 greg 3.16 COLOR *scan,
258     int len
259     )
260 greg 3.1 {
261 greg 3.17 static const char funcName[] = "tmCvColors";
262 greg 3.37 static uby8 gamtab[1024];
263 greg 3.20 static double curgam = .0;
264 greg 3.1 COLOR cmon;
265 greg 3.23 float lum, slum, d;
266 greg 3.16 int i;
267 greg 3.1
268 greg 3.16 if (tms == NULL)
269 greg 3.1 returnErr(TM_E_TMINVAL);
270 schorsch 3.13 if ((ls == NULL) | (scan == NULL) | (len < 0))
271 greg 3.1 returnErr(TM_E_ILLEGAL);
272 greg 3.24 if (tmFloat2BrtLUT == NULL) /* initialize */
273     tmCvLums(NULL, NULL, 0);
274 greg 3.21 if (cs != TM_NOCHROM && fabs(tms->mongam - curgam) > .02) {
275 greg 3.20 curgam = tms->mongam; /* (re)build table */
276     for (i = 1024; i--; )
277     gamtab[i] = (int)(256.*pow((i+.5)/1024., 1./curgam));
278     }
279 greg 3.1 for (i = len; i--; ) {
280 greg 3.16 if (tmNeedMatrix(tms)) { /* get monitor RGB */
281     colortrans(cmon, tms->cmat, scan[i]);
282 greg 3.15 } else {
283 greg 3.16 cmon[RED] = tms->inpsf*scan[i][RED];
284     cmon[GRN] = tms->inpsf*scan[i][GRN];
285     cmon[BLU] = tms->inpsf*scan[i][BLU];
286 greg 3.1 }
287 greg 3.23 #ifdef isfinite
288     if (!isfinite(cmon[RED]) || cmon[RED] < .0f) cmon[RED] = .0f;
289     if (!isfinite(cmon[GRN]) || cmon[GRN] < .0f) cmon[GRN] = .0f;
290     if (!isfinite(cmon[BLU]) || cmon[BLU] < .0f) cmon[BLU] = .0f;
291     #else
292     if (cmon[RED] < .0f) cmon[RED] = .0f;
293     if (cmon[GRN] < .0f) cmon[GRN] = .0f;
294     if (cmon[BLU] < .0f) cmon[BLU] = .0f;
295 greg 3.22 #endif
296 greg 3.1 /* world luminance */
297 greg 3.16 lum = tms->clf[RED]*cmon[RED] +
298     tms->clf[GRN]*cmon[GRN] +
299     tms->clf[BLU]*cmon[BLU] ;
300 greg 3.25 if (lum <= TM_NOLUM) { /* convert brightness */
301     lum = cmon[RED] = cmon[GRN] = cmon[BLU] = TM_NOLUM;
302 greg 3.23 ls[i] = TM_NOBRT;
303 greg 3.25 } else
304 greg 3.23 ls[i] = tmCvLumLUfp(&lum);
305 greg 3.1 if (cs == TM_NOCHROM) /* no color? */
306     continue;
307 greg 3.16 if (tms->flags & TM_F_MESOPIC && lum < LMESUPPER) {
308 greg 3.1 slum = scotlum(cmon); /* mesopic adj. */
309 greg 3.25 if (lum < LMESLOWER) {
310 greg 3.1 cmon[RED] = cmon[GRN] = cmon[BLU] = slum;
311 greg 3.25 } else {
312 greg 3.1 d = (lum - LMESLOWER)/(LMESUPPER - LMESLOWER);
313 greg 3.16 if (tms->flags & TM_F_BW)
314 greg 3.1 cmon[RED] = cmon[GRN] =
315     cmon[BLU] = d*lum;
316     else
317     scalecolor(cmon, d);
318 greg 3.23 d = (1.f-d)*slum;
319 greg 3.1 cmon[RED] += d;
320     cmon[GRN] += d;
321     cmon[BLU] += d;
322     }
323 greg 3.16 } else if (tms->flags & TM_F_BW) {
324 greg 3.53 int j = 3;
325     while (j--) cs[3*i+j] = tms->cdiv[j]/(TM_BRES>>8);
326     continue;
327 greg 3.1 }
328 greg 3.16 d = tms->clf[RED]*cmon[RED]/lum;
329 greg 3.23 cs[3*i ] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
330 greg 3.16 d = tms->clf[GRN]*cmon[GRN]/lum;
331 greg 3.23 cs[3*i+1] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
332 greg 3.16 d = tms->clf[BLU]*cmon[BLU]/lum;
333 greg 3.23 cs[3*i+2] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)];
334 greg 3.1 }
335     returnOK;
336     }
337    
338    
339 greg 3.9 int
340 greg 3.16 tmAddHisto( /* add values to histogram */
341     TMstruct *tms,
342     TMbright *ls,
343     int len,
344     int wt
345     )
346 greg 3.1 {
347 greg 3.17 static const char funcName[] = "tmAddHisto";
348 greg 3.9 int oldorig=0, oldlen, horig, hlen;
349 greg 3.51 int i;
350 greg 3.1
351 greg 3.16 if (tms == NULL)
352 greg 3.1 returnErr(TM_E_TMINVAL);
353 gwlarson 3.8 if (len < 0)
354     returnErr(TM_E_ILLEGAL);
355     if (len == 0)
356     returnOK;
357 greg 3.1 /* first, grow limits */
358 greg 3.16 if (tms->histo == NULL) {
359 greg 3.1 for (i = len; i-- && ls[i] < MINBRT; )
360     ;
361     if (i < 0)
362     returnOK;
363 greg 3.16 tms->hbrmin = tms->hbrmax = ls[i];
364 greg 3.1 oldlen = 0;
365     } else {
366 greg 3.31 oldorig = HISTI(tms->hbrmin);
367     oldlen = HISTI(tms->hbrmax) + 1 - oldorig;
368 greg 3.1 }
369     for (i = len; i--; ) {
370 greg 3.51 if (ls[i] < MINBRT)
371 greg 3.1 continue;
372 greg 3.51 if (ls[i] < tms->hbrmin)
373     tms->hbrmin = ls[i];
374     else if (ls[i] > tms->hbrmax)
375     tms->hbrmax = ls[i];
376 greg 3.1 }
377 greg 3.31 horig = HISTI(tms->hbrmin);
378     hlen = HISTI(tms->hbrmax) + 1 - horig;
379 greg 3.1 if (hlen > oldlen) { /* (re)allocate histogram */
380 greg 3.38 HIST_TYP *newhist = (HIST_TYP *)calloc(hlen, sizeof(HIST_TYP));
381 greg 3.1 if (newhist == NULL)
382     returnErr(TM_E_NOMEM);
383     if (oldlen) { /* copy and free old */
384 greg 3.46 memcpy(newhist+(oldorig-horig),
385     tms->histo, oldlen*sizeof(HIST_TYP));
386     free(tms->histo);
387 greg 3.1 }
388 greg 3.16 tms->histo = newhist;
389 greg 3.1 }
390     if (wt == 0)
391     returnOK;
392     for (i = len; i--; ) /* add in new counts */
393     if (ls[i] >= MINBRT)
394 greg 3.31 tms->histo[ HISTI(ls[i]) - horig ] += wt;
395 greg 3.1 returnOK;
396     }
397    
398    
399     static double
400 greg 3.16 htcontrs( /* human threshold contrast sensitivity, dL(La) */
401     double La
402     )
403 greg 3.1 {
404     double l10La, l10dL;
405     /* formula taken from Ferwerda et al. [SG96] */
406     if (La < 1.148e-4)
407     return(1.38e-3);
408     l10La = log10(La);
409     if (l10La < -1.44) /* rod response regime */
410     l10dL = pow(.405*l10La + 1.6, 2.18) - 2.86;
411     else if (l10La < -.0184)
412     l10dL = l10La - .395;
413     else if (l10La < 1.9) /* cone response regime */
414     l10dL = pow(.249*l10La + .65, 2.7) - .72;
415     else
416     l10dL = l10La - 1.255;
417    
418     return(exp10(l10dL));
419     }
420    
421    
422 greg 3.9 int
423 greg 3.16 tmFixedMapping( /* compute fixed, linear tone-mapping */
424     TMstruct *tms,
425     double expmult,
426 greg 3.46 double gamval,
427     double Lddyn
428 greg 3.16 )
429 greg 3.9 {
430 greg 3.17 static const char funcName[] = "tmFixedMapping";
431 greg 3.47 const int maxV = (1L<<(8*sizeof(TMAP_TYP))) - 1;
432     double minD;
433     int i;
434 greg 3.9
435 greg 3.16 if (!tmNewMap(tms))
436 greg 3.9 returnErr(TM_E_NOMEM);
437 greg 3.46 /* check arguments */
438     if (expmult <= .0) expmult = 1.;
439     if (gamval < MINGAM) gamval = tms->mongam;
440     if (Lddyn < MINLDDYN) Lddyn = DEFLDDYN;
441     minD = 1./Lddyn;
442 greg 3.28 for (i = tms->mbrmax-tms->mbrmin+1; i--; ) {
443 greg 3.46 double d;
444     d = expmult/tms->inpsf * tmLuminance(tms->mbrmin + i);
445 greg 3.48 if (d >= 2.*minD)
446     d -= minD;
447     else /* soft black crushing */
448     d *= d/(4.*minD);
449     d /= 1. - minD;
450 greg 3.46 d = TM_BRES*pow(d, 1./gamval);
451 greg 3.47 tms->lumap[i] = (d > maxV) ? maxV : (int)d;
452 greg 3.28 }
453 greg 3.9 returnOK;
454     }
455    
456    
457 greg 3.1 int
458 greg 3.16 tmComputeMapping( /* compute histogram tone-mapping */
459     TMstruct *tms,
460     double gamval,
461     double Lddyn,
462     double Ldmax
463     )
464 greg 3.1 {
465 greg 3.17 static const char funcName[] = "tmComputeMapping";
466 greg 3.38 HIST_TYP *histo;
467 greg 3.1 float *cumf;
468 greg 3.38 int brt0, histlen;
469 greg 3.39 HIST_TYP threshold, ceiling, trimmings, histot;
470 greg 3.44 double logLddyn, Ldmin, Lwavg, Tr, Lw, Ld;
471 greg 3.9 double sum;
472 greg 3.16 double d;
473     int i, j;
474 greg 3.1
475 greg 3.16 if (tms == NULL || tms->histo == NULL)
476 greg 3.1 returnErr(TM_E_TMINVAL);
477     /* check arguments */
478     if (Lddyn < MINLDDYN) Lddyn = DEFLDDYN;
479     if (Ldmax < MINLDMAX) Ldmax = DEFLDMAX;
480 greg 3.16 if (gamval < MINGAM) gamval = tms->mongam;
481 greg 3.1 /* compute handy values */
482     Ldmin = Ldmax/Lddyn;
483     logLddyn = log(Lddyn);
484 greg 3.31 i = HISTI(tms->hbrmin);
485 greg 3.32 brt0 = HISTV(i);
486 greg 3.31 histlen = HISTI(tms->hbrmax) + 1 - i;
487 greg 3.1 /* histogram total and mean */
488     histot = 0; sum = 0;
489     j = brt0 + histlen*HISTEP;
490     for (i = histlen; i--; ) {
491 greg 3.16 histot += tms->histo[i];
492 greg 3.30 sum += (double)(j -= HISTEP) * tms->histo[i];
493 greg 3.1 }
494 greg 3.29 if (!histot)
495 greg 3.1 returnErr(TM_E_TMFAIL);
496 greg 3.52 threshold = histot/500 + 1;
497 greg 3.1 Lwavg = tmLuminance( (double)sum / histot );
498 greg 3.29 /* use linear tone mapping? */
499 greg 3.52 if (tms->flags & TM_F_LINEAR ||
500 greg 3.30 tms->hbrmax - tms->hbrmin < TM_BRTSCALE*logLddyn)
501 greg 3.29 goto linearmap;
502 greg 3.9 /* clamp histogram */
503 greg 3.38 histo = (HIST_TYP *)malloc(histlen*sizeof(HIST_TYP));
504 greg 3.9 cumf = (float *)malloc((histlen+2)*sizeof(float));
505 schorsch 3.13 if ((histo == NULL) | (cumf == NULL))
506 greg 3.9 returnErr(TM_E_NOMEM);
507     cumf[histlen+1] = 1.; /* guard for assignment code */
508 greg 3.42 /* make malleable copy */
509 greg 3.46 memcpy(histo, tms->histo, histlen*sizeof(HIST_TYP));
510 greg 3.9 do { /* iterate to solution */
511     sum = 0; /* cumulative probability */
512     for (i = 0; i < histlen; i++) {
513 greg 3.45 cumf[i] = sum/histot;
514     sum += (double)histo[i];
515 greg 3.9 }
516     cumf[histlen] = 1.;
517 greg 3.16 Tr = histot * (double)(tms->hbrmax - tms->hbrmin) /
518 greg 3.45 ((double)TM_BRTSCALE*histlen*logLddyn);
519 greg 3.9 ceiling = Tr + 1.;
520     trimmings = 0; /* clip to envelope */
521     for (i = histlen; i--; ) {
522 greg 3.16 if (tms->flags & TM_F_HCONTR) {
523 greg 3.9 Lw = tmLuminance(brt0 + i*HISTEP);
524     Ld = Ldmin * exp( logLddyn *
525     .5*(cumf[i]+cumf[i+1]) );
526     ceiling = Tr * (htcontrs(Ld) * Lw) /
527     (htcontrs(Lw) * Ld) + 1.;
528 greg 3.2 }
529 greg 3.9 if (histo[i] > ceiling) {
530     trimmings += histo[i] - ceiling;
531     histo[i] = ceiling;
532 greg 3.1 }
533     }
534 greg 3.9 /* check if we're out of data */
535     if ((histot -= trimmings) <= threshold) {
536 greg 3.46 free(histo);
537     free(cumf);
538 greg 3.9 goto linearmap;
539     }
540 greg 3.50 } while (40*trimmings > histot);
541 greg 3.30 /* allocate space for mapping */
542     if (!tmNewMap(tms))
543     returnErr(TM_E_NOMEM);
544 greg 3.9 /* assign tone-mapping */
545 greg 3.16 for (i = tms->mbrmax-tms->mbrmin+1; i--; ) {
546     j = d = (double)i/(tms->mbrmax-tms->mbrmin)*histlen;
547 greg 3.9 d -= (double)j;
548     Ld = Ldmin*exp(logLddyn*((1.-d)*cumf[j]+d*cumf[j+1]));
549     d = (Ld - Ldmin)/(Ldmax - Ldmin);
550 greg 3.40 tms->lumap[i] = TM_BRES*pow(d, 1./gamval);
551 greg 3.1 }
552 greg 3.46 free(histo); /* clean up and return */
553     free(cumf);
554 greg 3.1 returnOK;
555 greg 3.9 linearmap: /* linear tone-mapping */
556 greg 3.16 if (tms->flags & TM_F_HCONTR)
557 greg 3.44 d = htcontrs(sqrt(Ldmax*Ldmin)) / htcontrs(Lwavg);
558 greg 3.9 else
559 greg 3.44 d = Ldmax / tmLuminance(tms->hbrmax);
560 greg 3.46 return(tmFixedMapping(tms, tms->inpsf*d/Ldmax, gamval, Lddyn));
561 greg 3.1 }
562    
563    
564     int
565 greg 3.16 tmMapPixels( /* apply tone-mapping to pixel(s) */
566     TMstruct *tms,
567 greg 3.37 uby8 *ps,
568 greg 3.16 TMbright *ls,
569 greg 3.37 uby8 *cs,
570 greg 3.16 int len
571     )
572 greg 3.1 {
573 greg 3.17 static const char funcName[] = "tmMapPixels";
574 greg 3.41 TMbright lv;
575     TMAP_TYP li;
576     int pv;
577 greg 3.1
578 greg 3.16 if (tms == NULL || tms->lumap == NULL)
579 greg 3.1 returnErr(TM_E_TMINVAL);
580 schorsch 3.13 if ((ps == NULL) | (ls == NULL) | (len < 0))
581 greg 3.1 returnErr(TM_E_ILLEGAL);
582     while (len--) {
583 greg 3.41 if ((lv = *ls++) < tms->mbrmin) {
584 greg 3.9 li = 0;
585     } else {
586 greg 3.41 if (lv > tms->mbrmax)
587     lv = tms->mbrmax;
588     li = tms->lumap[lv - tms->mbrmin];
589 greg 3.9 }
590 greg 3.43 if (cs == TM_NOCHROM) {
591     #if !(TM_BRES & 0xff)
592     *ps++ = li>=TM_BRES ? 255 : li/(TM_BRES>>8);
593     #else
594     *ps++ = li>=TM_BRES ? 255 : (li<<8)/TM_BRES;
595     #endif
596     } else {
597 greg 3.16 pv = *cs++ * li / tms->cdiv[RED];
598 greg 3.1 *ps++ = pv>255 ? 255 : pv;
599 greg 3.16 pv = *cs++ * li / tms->cdiv[GRN];
600 greg 3.1 *ps++ = pv>255 ? 255 : pv;
601 greg 3.16 pv = *cs++ * li / tms->cdiv[BLU];
602 greg 3.1 *ps++ = pv>255 ? 255 : pv;
603     }
604     }
605     returnOK;
606     }
607    
608    
609 greg 3.16 TMstruct *
610 greg 3.42 tmDup( /* duplicate tone mapping */
611 greg 3.16 TMstruct *tms
612     )
613 greg 3.3 {
614     int len;
615 greg 3.16 int i;
616     TMstruct *tmnew;
617 greg 3.3
618 greg 3.16 if (tms == NULL) /* anything to duplicate? */
619 greg 3.3 return(NULL);
620 greg 3.16 tmnew = (TMstruct *)malloc(sizeof(TMstruct));
621 greg 3.3 if (tmnew == NULL)
622     return(NULL);
623 greg 3.16 *tmnew = *tms; /* copy everything */
624 greg 3.3 if (tmnew->histo != NULL) { /* duplicate histogram */
625 greg 3.31 len = HISTI(tmnew->hbrmax) + 1 - HISTI(tmnew->hbrmin);
626 greg 3.38 tmnew->histo = (HIST_TYP *)malloc(len*sizeof(HIST_TYP));
627 greg 3.3 if (tmnew->histo != NULL)
628 greg 3.46 memcpy(tmnew->histo, tms->histo, len*sizeof(HIST_TYP));
629 greg 3.3 }
630     if (tmnew->lumap != NULL) { /* duplicate luminance mapping */
631 gregl 3.6 len = tmnew->mbrmax-tmnew->mbrmin+1;
632 greg 3.38 tmnew->lumap = (TMAP_TYP *)malloc(len*sizeof(TMAP_TYP));
633 greg 3.3 if (tmnew->lumap != NULL)
634 greg 3.46 memcpy(tmnew->lumap, tms->lumap, len*sizeof(TMAP_TYP));
635 greg 3.3 }
636 greg 3.4 /* clear package data */
637     for (i = tmNumPkgs; i--; )
638     tmnew->pd[i] = NULL;
639 greg 3.16 /* return copy */
640     return(tmnew);
641 greg 3.1 }
642    
643    
644     void
645     tmDone(tms) /* done with tone mapping -- destroy it */
646 greg 3.16 TMstruct *tms;
647 greg 3.1 {
648 greg 3.16 int i;
649     /* NULL arg. is equiv. to tms */
650     if (tms == NULL)
651 greg 3.1 return;
652     /* free tables */
653     if (tms->histo != NULL)
654 greg 3.46 free(tms->histo);
655 greg 3.1 if (tms->lumap != NULL)
656 greg 3.46 free(tms->lumap);
657 greg 3.4 /* free private data */
658     for (i = tmNumPkgs; i--; )
659     if (tms->pd[i] != NULL)
660     (*tmPkg[i]->Free)(tms->pd[i]);
661 greg 3.46 free(tms); /* free basic structure */
662 greg 3.9 }
663    
664     /******************** Shared but Private library routines *********************/
665    
666 greg 3.37 uby8 tmMesofact[BMESUPPER-BMESLOWER];
667 greg 3.9
668     void
669     tmMkMesofact() /* build mesopic lookup factor table */
670     {
671 greg 3.16 int i;
672 greg 3.9
673     if (tmMesofact[BMESUPPER-BMESLOWER-1])
674     return;
675    
676     for (i = BMESLOWER; i < BMESUPPER; i++)
677     tmMesofact[i-BMESLOWER] = 256. *
678     (tmLuminance(i) - LMESLOWER) /
679     (LMESUPPER - LMESLOWER);
680     }
681    
682    
683     int
684 greg 3.33 tmNewMap( /* allocate new tone-mapping array */
685     TMstruct *tms
686     )
687     {
688     if (tms->lumap != NULL && (tms->mbrmax - tms->mbrmin) !=
689     (tms->hbrmax - tms->hbrmin)) {
690 greg 3.46 free(tms->lumap);
691 greg 3.33 tms->lumap = NULL;
692     }
693     tms->mbrmin = tms->hbrmin;
694     tms->mbrmax = tms->hbrmax;
695     if (tms->mbrmin > tms->mbrmax)
696 greg 3.47 return(0);
697 greg 3.33 if (tms->lumap == NULL)
698 greg 3.46 tms->lumap = (TMAP_TYP *)calloc(tms->mbrmax-tms->mbrmin+1,
699     sizeof(TMAP_TYP));
700     else
701     memset(tms->lumap, 0, (tms->mbrmax-tms->mbrmin+1)*sizeof(TMAP_TYP));
702    
703 greg 3.33 return(tms->lumap != NULL);
704     }
705    
706    
707     int
708 greg 3.16 tmErrorReturn( /* error return (with message) */
709 greg 3.17 const char *func,
710 greg 3.16 TMstruct *tms,
711     int err
712     )
713 greg 3.9 {
714 greg 3.17 if (tms != NULL) {
715     tms->lastFunc = func;
716     tms->lastError = err;
717     if (tms->flags & TM_F_NOSTDERR)
718     return(err);
719     }
720 greg 3.9 fputs(func, stderr);
721     fputs(": ", stderr);
722     fputs(tmErrorMessage[err], stderr);
723     fputs("!\n", stderr);
724     return(err);
725 greg 3.1 }