src/tiff/libtiff/tif_color.c

   1 /* $Id$ */
   2
   3 /*
   4  * Copyright (c) 1988-1997 Sam Leffler
   5  * Copyright (c) 1991-1997 Silicon Graphics, Inc.
   6  *
   7  * Permission to use, copy, modify, distribute, and sell this software and
   8  * its documentation for any purpose is hereby granted without fee, provided
   9  * that (i) the above copyright notices and this permission notice appear in
  10  * all copies of the software and related documentation, and (ii) the names of
  11  * Sam Leffler and Silicon Graphics may not be used in any advertising or
  12  * publicity relating to the software without the specific, prior written
  13  * permission of Sam Leffler and Silicon Graphics.
  14  *
  15  * THE SOFTWARE IS PROVIDED "AS-IS" AND WITHOUT WARRANTY OF ANY KIND,
  16  * EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
  17  * WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
  18  *
  19  * IN NO EVENT SHALL SAM LEFFLER OR SILICON GRAPHICS BE LIABLE FOR
  20  * ANY SPECIAL, INCIDENTAL, INDIRECT OR CONSEQUENTIAL DAMAGES OF ANY KIND,
  21  * OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS,
  22  * WHETHER OR NOT ADVISED OF THE POSSIBILITY OF DAMAGE, AND ON ANY THEORY OF
  23  * LIABILITY, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE
  24  * OF THIS SOFTWARE.
  25  */
  26
  27 /*
  28  * CIE L*a*b* to CIE XYZ and CIE XYZ to RGB conversion routines are taken
  29  * from the VIPS library (http://www.vips.ecs.soton.ac.uk) with
  30  * the permission of John Cupitt, the VIPS author.
  31  */
  32
  33 /*
  34  * TIFF Library.
  35  *
  36  * Color space conversion routines.
  37  */
  38
  39 #include "tiffiop.h"
  40 #include <math.h>
  41
  42 /*
  43  * Convert color value from the CIE L*a*b* 1976 space to CIE XYZ.
  44  */
  45 void
  46 TIFFCIELabToXYZ(TIFFCIELabToRGB *cielab, uint32 l, int32 a, int32 b,
  47                 float *X, float *Y, float *Z)
  48 {
  49         float L = (float)l * 100.0F / 255.0F;
  50         float cby, tmp;
  51
  52         if( L < 8.856F ) {
  53                 *Y = (L * cielab->Y0) / 903.292F;
  54                 cby = 7.787F * (*Y / cielab->Y0) + 16.0F / 116.0F;
  55         } else {
  56                 cby = (L + 16.0F) / 116.0F;
  57                 *Y = cielab->Y0 * cby * cby * cby;
  58         }
  59
  60         tmp = (float)a / 500.0F + cby;
  61         if( tmp < 0.2069F )
  62                 *X = cielab->X0 * (tmp - 0.13793F) / 7.787F;
  63         else
  64                 *X = cielab->X0 * tmp * tmp * tmp;
  65
  66         tmp = cby - (float)b / 200.0F;
  67         if( tmp < 0.2069F )
  68                 *Z = cielab->Z0 * (tmp - 0.13793F) / 7.787F;
  69         else
  70                 *Z = cielab->Z0 * tmp * tmp * tmp;
  71 }
  72
  73 #define RINT(R) ((uint32)((R)>0?((R)+0.5):((R)-0.5)))
  74 /*
  75  * Convert color value from the XYZ space to RGB.
  76  */
  77 void
  78 TIFFXYZToRGB(TIFFCIELabToRGB *cielab, float X, float Y, float Z,
  79              uint32 *r, uint32 *g, uint32 *b)
  80 {
  81         int i;
  82         float Yr, Yg, Yb;
  83         float *matrix = &cielab->display.d_mat[0][0];
  84
  85         /* Multiply through the matrix to get luminosity values. */
  86         Yr =  matrix[0] * X + matrix[1] * Y + matrix[2] * Z;
  87         Yg =  matrix[3] * X + matrix[4] * Y + matrix[5] * Z;
  88         Yb =  matrix[6] * X + matrix[7] * Y + matrix[8] * Z;
  89
  90         /* Clip input */
  91         Yr = TIFFmax(Yr, cielab->display.d_Y0R);
  92         Yg = TIFFmax(Yg, cielab->display.d_Y0G);
  93         Yb = TIFFmax(Yb, cielab->display.d_Y0B);
  94
  95         /* Avoid overflow in case of wrong input values */
  96         Yr = TIFFmin(Yr, cielab->display.d_YCR);
  97         Yg = TIFFmin(Yg, cielab->display.d_YCG);
  98         Yb = TIFFmin(Yb, cielab->display.d_YCB);
  99
 100         /* Turn luminosity to colour value. */
 101         i = (int)((Yr - cielab->display.d_Y0R) / cielab->rstep);
 102         i = TIFFmin(cielab->range, i);
 103         *r = RINT(cielab->Yr2r[i]);
 104
 105         i = (int)((Yg - cielab->display.d_Y0G) / cielab->gstep);
 106         i = TIFFmin(cielab->range, i);
 107         *g = RINT(cielab->Yg2g[i]);
 108
 109         i = (int)((Yb - cielab->display.d_Y0B) / cielab->bstep);
 110         i = TIFFmin(cielab->range, i);
 111         *b = RINT(cielab->Yb2b[i]);
 112
 113         /* Clip output. */
 114         *r = TIFFmin(*r, cielab->display.d_Vrwr);
 115         *g = TIFFmin(*g, cielab->display.d_Vrwg);
 116         *b = TIFFmin(*b, cielab->display.d_Vrwb);
 117 }
 118 #undef RINT
 119
 120 /*
 121  * Allocate conversion state structures and make look_up tables for
 122  * the Yr,Yb,Yg <=> r,g,b conversions.
 123  */
 124 int
 125 TIFFCIELabToRGBInit(TIFFCIELabToRGB* cielab,
 126                     const TIFFDisplay *display, float *refWhite)
 127 {
 128         int i;
 129         double gamma;
 130
 131         cielab->range = CIELABTORGB_TABLE_RANGE;
 132
 133         _TIFFmemcpy(&cielab->display, display, sizeof(TIFFDisplay));
 134
 135         /* Red */
 136         gamma = 1.0 / cielab->display.d_gammaR ;
 137         cielab->rstep =
 138                 (cielab->display.d_YCR - cielab->display.d_Y0R) / cielab->range;
 139         for(i = 0; i <= cielab->range; i++) {
 140                 cielab->Yr2r[i] = cielab->display.d_Vrwr
 141                     * ((float)pow((double)i / cielab->range, gamma));
 142         }
 143
 144         /* Green */
 145         gamma = 1.0 / cielab->display.d_gammaG ;
 146         cielab->gstep =
 147             (cielab->display.d_YCR - cielab->display.d_Y0R) / cielab->range;
 148         for(i = 0; i <= cielab->range; i++) {
 149                 cielab->Yg2g[i] = cielab->display.d_Vrwg
 150                     * ((float)pow((double)i / cielab->range, gamma));
 151         }
 152
 153         /* Blue */
 154         gamma = 1.0 / cielab->display.d_gammaB ;
 155         cielab->bstep =
 156             (cielab->display.d_YCR - cielab->display.d_Y0R) / cielab->range;
 157         for(i = 0; i <= cielab->range; i++) {
 158                 cielab->Yb2b[i] = cielab->display.d_Vrwb
 159                     * ((float)pow((double)i / cielab->range, gamma));
 160         }
 161
 162         /* Init reference white point */
 163         cielab->X0 = refWhite[0];
 164         cielab->Y0 = refWhite[1];
 165         cielab->Z0 = refWhite[2];
 166
 167         return 0;
 168 }
 169
 170 /*
 171  * Convert color value from the YCbCr space to CIE XYZ.
 172  * The colorspace conversion algorithm comes from the IJG v5a code;
 173  * see below for more information on how it works.
 174  */
 175 #define SHIFT                   16
 176 #define FIX(x)                  ((int32)((x) * (1L<<SHIFT) + 0.5))
 177 #define ONE_HALF                ((int32)(1<<(SHIFT-1)))
 178 #define Code2V(c, RB, RW, CR)   ((((c)-(int32)(RB))*(float)(CR))/(float)(((RW)-(RB)) ? ((RW)-(RB)) : 1))
 179 #define CLAMP(f,min,max)        ((f)<(min)?(min):(f)>(max)?(max):(f))
 180 #define HICLAMP(f,max)          ((f)>(max)?(max):(f))
 181
 182 void
 183 TIFFYCbCrtoRGB(TIFFYCbCrToRGB *ycbcr, uint32 Y, int32 Cb, int32 Cr,
 184                uint32 *r, uint32 *g, uint32 *b)
 185 {
 186         int32 i;
 187
 188         /* XXX: Only 8-bit YCbCr input supported for now */
 189         Y = HICLAMP(Y, 255), Cb = CLAMP(Cb, 0, 255), Cr = CLAMP(Cr, 0, 255);
 190
 191         i = ycbcr->Y_tab[Y] + ycbcr->Cr_r_tab[Cr];
 192         *r = CLAMP(i, 0, 255);
 193         i = ycbcr->Y_tab[Y]
 194             + (int)((ycbcr->Cb_g_tab[Cb] + ycbcr->Cr_g_tab[Cr]) >> SHIFT);
 195         *g = CLAMP(i, 0, 255);
 196         i = ycbcr->Y_tab[Y] + ycbcr->Cb_b_tab[Cb];
 197         *b = CLAMP(i, 0, 255);
 198 }
 199
 200 /*
 201  * Initialize the YCbCr->RGB conversion tables.  The conversion
 202  * is done according to the 6.0 spec:
 203  *
 204  *    R = Y + Cr*(2 - 2*LumaRed)
 205  *    B = Y + Cb*(2 - 2*LumaBlue)
 206  *    G =   Y
 207  *        - LumaBlue*Cb*(2-2*LumaBlue)/LumaGreen
 208  *        - LumaRed*Cr*(2-2*LumaRed)/LumaGreen
 209  *
 210  * To avoid floating point arithmetic the fractional constants that
 211  * come out of the equations are represented as fixed point values
 212  * in the range 0...2^16.  We also eliminate multiplications by
 213  * pre-calculating possible values indexed by Cb and Cr (this code
 214  * assumes conversion is being done for 8-bit samples).
 215  */
 216 int
 217 TIFFYCbCrToRGBInit(TIFFYCbCrToRGB* ycbcr, float *luma, float *refBlackWhite)
 218 {
 219     TIFFRGBValue* clamptab;
 220     int i;
 221
 222 #define LumaRed     luma[0]
 223 #define LumaGreen   luma[1]
 224 #define LumaBlue    luma[2]
 225
 226     clamptab = (TIFFRGBValue*)(
 227         (uint8*) ycbcr+TIFFroundup_32(sizeof (TIFFYCbCrToRGB), sizeof (long)));
 228     _TIFFmemset(clamptab, 0, 256);              /* v < 0 => 0 */
 229     ycbcr->clamptab = (clamptab += 256);
 230     for (i = 0; i < 256; i++)
 231         clamptab[i] = (TIFFRGBValue) i;
 232     _TIFFmemset(clamptab+256, 255, 2*256);      /* v > 255 => 255 */
 233     ycbcr->Cr_r_tab = (int*) (clamptab + 3*256);
 234     ycbcr->Cb_b_tab = ycbcr->Cr_r_tab + 256;
 235     ycbcr->Cr_g_tab = (int32*) (ycbcr->Cb_b_tab + 256);
 236     ycbcr->Cb_g_tab = ycbcr->Cr_g_tab + 256;
 237     ycbcr->Y_tab = ycbcr->Cb_g_tab + 256;
 238
 239     { float f1 = 2-2*LumaRed;           int32 D1 = FIX(f1);
 240       float f2 = LumaRed*f1/LumaGreen;  int32 D2 = -FIX(f2);
 241       float f3 = 2-2*LumaBlue;          int32 D3 = FIX(f3);
 242       float f4 = LumaBlue*f3/LumaGreen; int32 D4 = -FIX(f4);
 243       int x;
 244
 245 #undef LumaBlue
 246 #undef LumaGreen
 247 #undef LumaRed
 248
 249       /*
 250        * i is the actual input pixel value in the range 0..255
 251        * Cb and Cr values are in the range -128..127 (actually
 252        * they are in a range defined by the ReferenceBlackWhite
 253        * tag) so there is some range shifting to do here when
 254        * constructing tables indexed by the raw pixel data.
 255        */
 256       for (i = 0, x = -128; i < 256; i++, x++) {
 257             int32 Cr = (int32)Code2V(x, refBlackWhite[4] - 128.0F,
 258                             refBlackWhite[5] - 128.0F, 127);
 259             int32 Cb = (int32)Code2V(x, refBlackWhite[2] - 128.0F,
 260                             refBlackWhite[3] - 128.0F, 127);
 261
 262             ycbcr->Cr_r_tab[i] = (int32)((D1*Cr + ONE_HALF)>>SHIFT);
 263             ycbcr->Cb_b_tab[i] = (int32)((D3*Cb + ONE_HALF)>>SHIFT);
 264             ycbcr->Cr_g_tab[i] = D2*Cr;
 265             ycbcr->Cb_g_tab[i] = D4*Cb + ONE_HALF;
 266             ycbcr->Y_tab[i] =
 267                     (int32)Code2V(x + 128, refBlackWhite[0], refBlackWhite[1], 255);
 268       }
 269     }
 270
 271     return 0;
 272 }
 273 #undef  HICLAMP
 274 #undef  CLAMP
 275 #undef  Code2V
 276 #undef  SHIFT
 277 #undef  ONE_HALF
 278 #undef  FIX
 279
 280 /* vim: set ts=8 sts=8 sw=8 noet: */
 281 /*
 282  * Local Variables:
 283  * mode: c
 284  * c-basic-offset: 8
 285  * fill-column: 78
 286  * End:
 287  */