--- ray/src/common/tonemap.c 2006/03/23 22:01:43 3.22 +++ ray/src/common/tonemap.c 2008/07/03 03:30:21 3.31 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: tonemap.c,v 3.22 2006/03/23 22:01:43 greg Exp $"; +static const char RCSid[] = "$Id: tonemap.c,v 3.31 2008/07/03 03:30:21 greg Exp $"; #endif /* * Tone mapping functions. @@ -22,7 +22,12 @@ static const char RCSid[] = "$Id: tonemap.c,v 3.22 200 struct tmPackage *tmPkg[TM_MAXPKG]; int tmNumPkgs = 0; /* number of registered packages */ + /* luminance->brightness lookup */ +static TMbright *tmFloat2BrtLUT = NULL; +#define tmCvLumLUfp(pf) tmFloat2BrtLUT[*(int32 *)(pf) >> 15] + + TMstruct * tmInit( /* initialize new tone mapping */ int flags, @@ -48,7 +53,7 @@ double gamval tmnew->clf[GRN] = rgb2xyzmat[1][1]; tmnew->clf[BLU] = rgb2xyzmat[1][2]; } else { - comprgb2xyzWBmat(cmat, tmnew->monpri=monpri); + comprgb2xyzmat(cmat, tmnew->monpri=monpri); tmnew->clf[RED] = cmat[1][0]; tmnew->clf[GRN] = cmat[1][1]; tmnew->clf[BLU] = cmat[1][2]; @@ -111,7 +116,7 @@ MEM_PTR dat tms->clf[CIEX] = tms->clf[CIEZ] = 0.; tms->clf[CIEY] = 1.; } else { - comprgb2xyzWBmat(tms->cmat, tms->monpri); + comprgb2xyzmat(tms->cmat, tms->monpri); tms->clf[RED] = tms->cmat[1][0]; tms->clf[GRN] = tms->cmat[1][1]; tms->clf[BLU] = tms->cmat[1][2]; @@ -149,7 +154,7 @@ MEM_PTR dat void -tmClearHisto( /* clear current histogram */ +tmClearHisto( /* clear current histogram */ TMstruct *tms ) { @@ -160,7 +165,93 @@ TMstruct *tms } +TMbright +tmCvLuminance( /* convert a single luminance */ +double lum +) +{ + double d; + +#ifdef isfinite + if (!isfinite(lum) || lum <= TM_NOLUM) +#else + if (lum <= TM_NOLUM) +#endif + return(TM_NOBRT); + d = TM_BRTSCALE*log(lum); + if (d > 0.) + return((TMbright)(d+.5)); + return((TMbright)(d-.5)); +} + + int +tmCvLums( /* convert luminances using lookup */ +TMbright *ls, +float *scan, +int len +) +{ + if (tmFloat2BrtLUT == NULL) { /* initialize lookup table */ + int32 i; + tmFloat2BrtLUT = (TMbright *)malloc(sizeof(TMbright)*0x10000); + if (tmFloat2BrtLUT == NULL) + return(TM_E_NOMEM); + for (i = 0; i < 0x10000; i++) { + int32 l = (i<<1 | 1) << 14; +#ifndef isfinite + if ((l & 0x7f800000) == 0x7f800000) + tmFloat2BrtLUT[i] = TM_NOBRT; + else +#endif + tmFloat2BrtLUT[i] = tmCvLuminance(*(float *)&l); + } + } + if (len <= 0) + return(TM_E_OK); + if ((ls == NULL) | (scan == NULL)) + return(TM_E_ILLEGAL); + while (len--) { + if (*scan <= TM_NOLUM) { + *ls++ = TM_NOBRT; + ++scan; + continue; + } + *ls++ = tmCvLumLUfp(scan++); + } + return(TM_E_OK); +} + + +int +tmCvGrays( /* convert float gray values */ +TMstruct *tms, +TMbright *ls, +float *scan, +int len +) +{ + static const char funcName[] = "tmCvGrays"; + int i; + + if (tms == NULL) + returnErr(TM_E_TMINVAL); + if ((ls == NULL) | (scan == NULL) | (len < 0)) + returnErr(TM_E_ILLEGAL); + if (tmFloat2BrtLUT == NULL) /* initialize */ + tmCvLums(NULL, NULL, 0); + for (i = len; i--; ) { + float lum = tms->inpsf * scan[i]; + if (lum <= TM_NOLUM) + ls[i] = TM_NOBRT; + else + ls[i] = tmCvLumLUfp(&lum); + } + returnOK; +} + + +int tmCvColors( /* convert float colors */ TMstruct *tms, TMbright *ls, @@ -170,18 +261,18 @@ int len ) { static const char funcName[] = "tmCvColors"; - static COLOR csmall = {.5*MINLUM, .5*MINLUM, .5*MINLUM}; static BYTE gamtab[1024]; static double curgam = .0; COLOR cmon; - double lum, slum; - double d; + float lum, slum, d; int i; if (tms == NULL) returnErr(TM_E_TMINVAL); if ((ls == NULL) | (scan == NULL) | (len < 0)) returnErr(TM_E_ILLEGAL); + if (tmFloat2BrtLUT == NULL) /* initialize */ + tmCvLums(NULL, NULL, 0); if (cs != TM_NOCHROM && fabs(tms->mongam - curgam) > .02) { curgam = tms->mongam; /* (re)build table */ for (i = 1024; i--; ) @@ -195,41 +286,38 @@ int len cmon[GRN] = tms->inpsf*scan[i][GRN]; cmon[BLU] = tms->inpsf*scan[i][BLU]; } -#ifdef isnan - if (isnan(cmon[RED]) || isinf(cmon[RED])) cmon[RED] = .0f; - if (isnan(cmon[GRN]) || isinf(cmon[GRN])) cmon[GRN] = .0f; - if (isnan(cmon[BLU]) || isinf(cmon[BLU])) cmon[BLU] = .0f; +#ifdef isfinite + if (!isfinite(cmon[RED]) || cmon[RED] < .0f) cmon[RED] = .0f; + if (!isfinite(cmon[GRN]) || cmon[GRN] < .0f) cmon[GRN] = .0f; + if (!isfinite(cmon[BLU]) || cmon[BLU] < .0f) cmon[BLU] = .0f; +#else + if (cmon[RED] < .0f) cmon[RED] = .0f; + if (cmon[GRN] < .0f) cmon[GRN] = .0f; + if (cmon[BLU] < .0f) cmon[BLU] = .0f; #endif /* world luminance */ lum = tms->clf[RED]*cmon[RED] + tms->clf[GRN]*cmon[GRN] + tms->clf[BLU]*cmon[BLU] ; - /* check range */ - if (clipgamut(cmon, lum, CGAMUT_LOWER, csmall, cwhite)) - lum = tms->clf[RED]*cmon[RED] + - tms->clf[GRN]*cmon[GRN] + - tms->clf[BLU]*cmon[BLU] ; - if (lum < MINLUM) { - ls[i] = MINBRT-1; /* bogus value */ - lum = MINLUM; - } else { - d = TM_BRTSCALE*log(lum); /* encode it */ - ls[i] = d>0. ? (int)(d+.5) : (int)(d-.5); - } + if (lum <= TM_NOLUM) { /* convert brightness */ + lum = cmon[RED] = cmon[GRN] = cmon[BLU] = TM_NOLUM; + ls[i] = TM_NOBRT; + } else + ls[i] = tmCvLumLUfp(&lum); if (cs == TM_NOCHROM) /* no color? */ continue; if (tms->flags & TM_F_MESOPIC && lum < LMESUPPER) { slum = scotlum(cmon); /* mesopic adj. */ - if (lum < LMESLOWER) + if (lum < LMESLOWER) { cmon[RED] = cmon[GRN] = cmon[BLU] = slum; - else { + } else { d = (lum - LMESLOWER)/(LMESUPPER - LMESLOWER); if (tms->flags & TM_F_BW) cmon[RED] = cmon[GRN] = cmon[BLU] = d*lum; else scalecolor(cmon, d); - d = (1.-d)*slum; + d = (1.f-d)*slum; cmon[RED] += d; cmon[GRN] += d; cmon[BLU] += d; @@ -238,60 +326,17 @@ int len cmon[RED] = cmon[GRN] = cmon[BLU] = lum; } d = tms->clf[RED]*cmon[RED]/lum; - cs[3*i ] = d>=.999 ? 255 : gamtab[(int)(1024.*d)]; + cs[3*i ] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)]; d = tms->clf[GRN]*cmon[GRN]/lum; - cs[3*i+1] = d>=.999 ? 255 : gamtab[(int)(1024.*d)]; + cs[3*i+1] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)]; d = tms->clf[BLU]*cmon[BLU]/lum; - cs[3*i+2] = d>=.999 ? 255 : gamtab[(int)(1024.*d)]; + cs[3*i+2] = d>=.999f ? 255 : gamtab[(int)(1024.f*d)]; } returnOK; } -TMbright -tmCvLuminance( /* convert a single luminance */ -double lum -) -{ - double d; - - if (lum <= TM_NOLUM) - return(TM_NOBRT); - d = TM_BRTSCALE*log(lum); - if (d > 0.) - return((TMbright)(d+.5)); - return((TMbright)(d-.5)); -} - - int -tmCvGrays( /* convert float gray values */ -TMstruct *tms, -TMbright *ls, -float *scan, -int len -) -{ - static const char funcName[] = "tmCvGrays"; - double d; - int i; - - if (tms == NULL) - returnErr(TM_E_TMINVAL); - if ((ls == NULL) | (scan == NULL) | (len < 0)) - returnErr(TM_E_ILLEGAL); - for (i = len; i--; ) - if (scan[i] <= TM_NOLUM) { - ls[i] = TM_NOBRT; /* bogus value */ - } else { - d = TM_BRTSCALE*log(scan[i]); /* encode it */ - ls[i] = d>0. ? (int)(d+.5) : (int)(d-.5); - } - returnOK; -} - - -int tmAddHisto( /* add values to histogram */ TMstruct *tms, TMbright *ls, @@ -318,8 +363,8 @@ int wt tms->hbrmin = tms->hbrmax = ls[i]; oldlen = 0; } else { - oldorig = (tms->hbrmin-MINBRT)/HISTEP; - oldlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - oldorig; + oldorig = HISTI(tms->hbrmin); + oldlen = HISTI(tms->hbrmax) + 1 - oldorig; } for (i = len; i--; ) { if ((j = ls[i]) < MINBRT) @@ -329,8 +374,8 @@ int wt else if (j > tms->hbrmax) tms->hbrmax = j; } - horig = (tms->hbrmin-MINBRT)/HISTEP; - hlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - horig; + horig = HISTI(tms->hbrmin); + hlen = HISTI(tms->hbrmax) + 1 - horig; if (hlen > oldlen) { /* (re)allocate histogram */ int *newhist = (int *)calloc(hlen, sizeof(int)); if (newhist == NULL) @@ -346,7 +391,7 @@ int wt returnOK; for (i = len; i--; ) /* add in new counts */ if (ls[i] >= MINBRT) - tms->histo[ (ls[i]-MINBRT)/HISTEP - horig ] += wt; + tms->histo[ HISTI(ls[i]) - horig ] += wt; returnOK; } @@ -413,10 +458,12 @@ double gamval if (gamval < MINGAM) gamval = tms->mongam; d = log(expmult/tms->inpsf); - for (i = tms->mbrmax-tms->mbrmin+1; i--; ) - tms->lumap[i] = 256. * exp( + for (i = tms->mbrmax-tms->mbrmin+1; i--; ) { + double val = 256. * exp( ( d + (tms->mbrmin+i)*(1./TM_BRTSCALE) ) - / gamval ); + / gamval); + tms->lumap[i] = val >= (double)0xffff ? 0xffff : (int)val; + } returnOK; } @@ -449,25 +496,23 @@ double Ldmax Ldmin = Ldmax/Lddyn; logLddyn = log(Lddyn); Ldavg = sqrt(Ldmax*Ldmin); - i = (tms->hbrmin-MINBRT)/HISTEP; + i = HISTI(tms->hbrmin); brt0 = MINBRT + HISTEP/2 + i*HISTEP; - histlen = (tms->hbrmax-MINBRT)/HISTEP + 1 - i; + histlen = HISTI(tms->hbrmax) + 1 - i; /* histogram total and mean */ histot = 0; sum = 0; j = brt0 + histlen*HISTEP; for (i = histlen; i--; ) { histot += tms->histo[i]; - sum += (j -= HISTEP) * tms->histo[i]; + sum += (double)(j -= HISTEP) * tms->histo[i]; } threshold = histot*0.005 + .5; - if (threshold < 4) + if (!histot) returnErr(TM_E_TMFAIL); Lwavg = tmLuminance( (double)sum / histot ); - /* allocate space for mapping */ - if (!tmNewMap(tms)) - returnErr(TM_E_NOMEM); /* use linear tone mapping? */ - if (tms->flags & TM_F_LINEAR) + if (tms->flags & TM_F_LINEAR || threshold < 4 || + tms->hbrmax - tms->hbrmin < TM_BRTSCALE*logLddyn) goto linearmap; /* clamp histogram */ histo = (int *)malloc(histlen*sizeof(int)); @@ -508,6 +553,9 @@ double Ldmax goto linearmap; } } while (trimmings > threshold); + /* allocate space for mapping */ + if (!tmNewMap(tms)) + returnErr(TM_E_NOMEM); /* assign tone-mapping */ for (i = tms->mbrmax-tms->mbrmin+1; i--; ) { j = d = (double)i/(tms->mbrmax-tms->mbrmin)*histlen; @@ -585,8 +633,7 @@ TMstruct *tms return(NULL); *tmnew = *tms; /* copy everything */ if (tmnew->histo != NULL) { /* duplicate histogram */ - len = (tmnew->hbrmax-MINBRT)/HISTEP + 1 - - (tmnew->hbrmin-MINBRT)/HISTEP; + len = HISTI(tmnew->hbrmax) + 1 - HISTI(tmnew->hbrmin); tmnew->histo = (int *)malloc(len*sizeof(int)); if (tmnew->histo != NULL) for (i = len; i--; )