workbench/libs/lcms2/src/cmssm.c

   1 //---------------------------------------------------------------------------------
   2 //
   3 //  Little Color Management System
   4 //  Copyright (c) 1998-2011 Marti Maria Saguer
   5 //
   6 // Permission is hereby granted, free of charge, to any person obtaining
   7 // a copy of this software and associated documentation files (the "Software"),
   8 // to deal in the Software without restriction, including without limitation
   9 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
  10 // and/or sell copies of the Software, and to permit persons to whom the Software
  11 // is furnished to do so, subject to the following conditions:
  12 //
  13 // The above copyright notice and this permission notice shall be included in
  14 // all copies or substantial portions of the Software.
  15 //
  16 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
  17 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO
  18 // THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
  19 // NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
  20 // LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
  21 // OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
  22 // WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
  23 //
  24 //---------------------------------------------------------------------------------
  25 //
  26
  27 #include "lcms2_internal.h"
  28
  29
  30 // ------------------------------------------------------------------------
  31
  32 // Gamut boundary description by using Jan Morovic's Segment maxima method
  33 // Many thanks to Jan for allowing me to use his algorithm.
  34
  35 // r = C*
  36 // alpha = Hab
  37 // theta = L*
  38
  39 #define SECTORS 16      // number of divisions in alpha and theta
  40
  41 // Spherical coordinates
  42 typedef struct {
  43
  44     cmsFloat64Number r;
  45     cmsFloat64Number alpha;
  46     cmsFloat64Number theta;
  47
  48 } cmsSpherical;
  49
  50 typedef  enum {
  51         GP_EMPTY,
  52         GP_SPECIFIED,
  53         GP_MODELED
  54
  55     } GDBPointType;
  56
  57
  58 typedef struct {
  59
  60     GDBPointType Type;
  61     cmsSpherical p;         // Keep also alpha & theta of maximum
  62
  63 } cmsGDBPoint;
  64
  65
  66 typedef struct {
  67
  68     cmsContext ContextID;
  69     cmsGDBPoint Gamut[SECTORS][SECTORS];
  70
  71 } cmsGDB;
  72
  73
  74 // A line using the parametric form
  75 // P = a + t*u
  76 typedef struct {
  77
  78     cmsVEC3 a;
  79     cmsVEC3 u;
  80
  81 } cmsLine;
  82
  83
  84 // A plane using the parametric form
  85 // Q = b + r*v + s*w
  86 typedef struct {
  87
  88     cmsVEC3 b;
  89     cmsVEC3 v;
  90     cmsVEC3 w;
  91
  92 } cmsPlane;
  93
  94
  95
  96 // --------------------------------------------------------------------------------------------
  97
  98 // ATAN2() which always returns degree positive numbers
  99
 100 static
 101 cmsFloat64Number _cmsAtan2(cmsFloat64Number y, cmsFloat64Number x)
 102 {
 103     cmsFloat64Number a;
 104
 105     // Deal with undefined case
 106     if (x == 0.0 && y == 0.0) return 0;
 107
 108     a = (atan2(y, x) * 180.0) / M_PI;
 109
 110     while (a < 0) {
 111         a += 360;
 112     }
 113
 114     return a;
 115 }
 116
 117 // Convert to spherical coordinates
 118 static
 119 void ToSpherical(cmsSpherical* sp, const cmsVEC3* v)
 120 {
 121
 122     cmsFloat64Number L, a, b;
 123
 124     L = v ->n[VX];
 125     a = v ->n[VY];
 126     b = v ->n[VZ];
 127
 128     sp ->r = sqrt( L*L + a*a + b*b );
 129
 130    if (sp ->r == 0) {
 131         sp ->alpha = sp ->theta = 0;
 132         return;
 133     }
 134
 135     sp ->alpha = _cmsAtan2(a, b);
 136     sp ->theta = _cmsAtan2(sqrt(a*a + b*b), L);
 137 }
 138
 139
 140 // Convert to cartesian from spherical
 141 static
 142 void ToCartesian(cmsVEC3* v, const cmsSpherical* sp)
 143 {
 144     cmsFloat64Number sin_alpha;
 145     cmsFloat64Number cos_alpha;
 146     cmsFloat64Number sin_theta;
 147     cmsFloat64Number cos_theta;
 148     cmsFloat64Number L, a, b;
 149
 150     sin_alpha = sin((M_PI * sp ->alpha) / 180.0);
 151     cos_alpha = cos((M_PI * sp ->alpha) / 180.0);
 152     sin_theta = sin((M_PI * sp ->theta) / 180.0);
 153     cos_theta = cos((M_PI * sp ->theta) / 180.0);
 154
 155     a = sp ->r * sin_theta * sin_alpha;
 156     b = sp ->r * sin_theta * cos_alpha;
 157     L = sp ->r * cos_theta;
 158
 159     v ->n[VX] = L;
 160     v ->n[VY] = a;
 161     v ->n[VZ] = b;
 162 }
 163
 164
 165 // Quantize sector of a spherical coordinate. Saturate 360, 180 to last sector
 166 // The limits are the centers of each sector, so
 167 static
 168 void QuantizeToSector(const cmsSpherical* sp, int* alpha, int* theta)
 169 {
 170     *alpha = (int) floor(((sp->alpha * (SECTORS)) / 360.0) );
 171     *theta = (int) floor(((sp->theta * (SECTORS)) / 180.0) );
 172
 173     if (*alpha >= SECTORS)
 174         *alpha = SECTORS-1;
 175     if (*theta >= SECTORS)
 176         *theta = SECTORS-1;
 177 }
 178
 179
 180 // Line determined by 2 points
 181 static
 182 void LineOf2Points(cmsLine* line, cmsVEC3* a, cmsVEC3* b)
 183 {
 184
 185     _cmsVEC3init(&line ->a, a ->n[VX], a ->n[VY], a ->n[VZ]);
 186     _cmsVEC3init(&line ->u, b ->n[VX] - a ->n[VX],
 187                             b ->n[VY] - a ->n[VY],
 188                             b ->n[VZ] - a ->n[VZ]);
 189 }
 190
 191
 192 // Evaluate parametric line
 193 static
 194 void GetPointOfLine(cmsVEC3* p, const cmsLine* line, cmsFloat64Number t)
 195 {
 196     p ->n[VX] = line ->a.n[VX] + t * line->u.n[VX];
 197     p ->n[VY] = line ->a.n[VY] + t * line->u.n[VY];
 198     p ->n[VZ] = line ->a.n[VZ] + t * line->u.n[VZ];
 199 }
 200
 201
 202
 203 /*
 204     Closest point in sector line1 to sector line2 (both are defined as 0 <=t <= 1)
 205     http://softsurfer.com/Archive/algorithm_0106/algorithm_0106.htm
 206
 207     Copyright 2001, softSurfer (www.softsurfer.com)
 208     This code may be freely used and modified for any purpose
 209     providing that this copyright notice is included with it.
 210     SoftSurfer makes no warranty for this code, and cannot be held
 211     liable for any real or imagined damage resulting from its use.
 212     Users of this code must verify correctness for their application.
 213
 214 */
 215
 216 static
 217 cmsBool ClosestLineToLine(cmsVEC3* r, const cmsLine* line1, const cmsLine* line2)
 218 {
 219     cmsFloat64Number a, b, c, d, e, D;
 220     cmsFloat64Number sc, sN, sD;
 221     cmsFloat64Number tc, tN, tD;
 222     cmsVEC3 w0;
 223
 224     _cmsVEC3minus(&w0, &line1 ->a, &line2 ->a);
 225
 226     a  = _cmsVEC3dot(&line1 ->u, &line1 ->u);
 227     b  = _cmsVEC3dot(&line1 ->u, &line2 ->u);
 228     c  = _cmsVEC3dot(&line2 ->u, &line2 ->u);
 229     d  = _cmsVEC3dot(&line1 ->u, &w0);
 230     e  = _cmsVEC3dot(&line2 ->u, &w0);
 231
 232     D  = a*c - b * b;      // Denominator
 233     sD = tD = D;           // default sD = D >= 0
 234
 235     if (D <  MATRIX_DET_TOLERANCE) {   // the lines are almost parallel
 236
 237         sN = 0.0;        // force using point P0 on segment S1
 238         sD = 1.0;        // to prevent possible division by 0.0 later
 239         tN = e;
 240         tD = c;
 241     }
 242     else {                // get the closest points on the infinite lines
 243
 244         sN = (b*e - c*d);
 245         tN = (a*e - b*d);
 246
 247         if (sN < 0.0) {       // sc < 0 => the s=0 edge is visible
 248
 249             sN = 0.0;
 250             tN = e;
 251             tD = c;
 252         }
 253         else if (sN > sD) {   // sc > 1 => the s=1 edge is visible
 254             sN = sD;
 255             tN = e + b;
 256             tD = c;
 257         }
 258     }
 259
 260     if (tN < 0.0) {           // tc < 0 => the t=0 edge is visible
 261
 262         tN = 0.0;
 263         // recompute sc for this edge
 264         if (-d < 0.0)
 265             sN = 0.0;
 266         else if (-d > a)
 267             sN = sD;
 268         else {
 269             sN = -d;
 270             sD = a;
 271         }
 272     }
 273     else if (tN > tD) {      // tc > 1 => the t=1 edge is visible
 274
 275         tN = tD;
 276
 277         // recompute sc for this edge
 278         if ((-d + b) < 0.0)
 279             sN = 0;
 280         else if ((-d + b) > a)
 281             sN = sD;
 282         else {
 283             sN = (-d + b);
 284             sD = a;
 285         }
 286     }
 287     // finally do the division to get sc and tc
 288     sc = (fabs(sN) < MATRIX_DET_TOLERANCE ? 0.0 : sN / sD);
 289     tc = (fabs(tN) < MATRIX_DET_TOLERANCE ? 0.0 : tN / tD);
 290
 291     GetPointOfLine(r, line1, sc);
 292     return TRUE;
 293 }
 294
 295
 296
 297 // ------------------------------------------------------------------ Wrapper
 298
 299
 300 // Allocate & free structure
 301 cmsHANDLE  CMSEXPORT cmsGBDAlloc(cmsContext ContextID)
 302 {
 303     cmsGDB* gbd = (cmsGDB*) _cmsMallocZero(ContextID, sizeof(cmsGDB));
 304     if (gbd == NULL) return NULL;
 305
 306     gbd -> ContextID = ContextID;
 307
 308     return (cmsHANDLE) gbd;
 309 }
 310
 311
 312 void CMSEXPORT cmsGBDFree(cmsHANDLE hGBD)
 313 {
 314     cmsGDB* gbd = (cmsGDB*) hGBD;
 315     if (hGBD != NULL)
 316         _cmsFree(gbd->ContextID, (void*) gbd);
 317 }
 318
 319
 320 // Auxiliar to retrieve a pointer to the segmentr containing the Lab value
 321 static
 322 cmsGDBPoint* GetPoint(cmsGDB* gbd, const cmsCIELab* Lab, cmsSpherical* sp)
 323 {
 324     cmsVEC3 v;
 325     int alpha, theta;
 326
 327     // Housekeeping
 328     _cmsAssert(gbd != NULL);
 329     _cmsAssert(Lab != NULL);
 330     _cmsAssert(sp != NULL);
 331
 332     // Center L* by substracting half of its domain, that's 50
 333     _cmsVEC3init(&v, Lab ->L - 50.0, Lab ->a, Lab ->b);
 334
 335     // Convert to spherical coordinates
 336     ToSpherical(sp, &v);
 337
 338     if (sp ->r < 0 || sp ->alpha < 0 || sp->theta < 0) {
 339          cmsSignalError(gbd ->ContextID, cmsERROR_RANGE, "spherical value out of range");
 340          return NULL;
 341     }
 342
 343     // On which sector it falls?
 344     QuantizeToSector(sp, &alpha, &theta);
 345
 346     if (alpha < 0 || theta < 0 || alpha >= SECTORS || theta >= SECTORS) {
 347          cmsSignalError(gbd ->ContextID, cmsERROR_RANGE, " quadrant out of range");
 348          return NULL;
 349     }
 350
 351     // Get pointer to the sector
 352     return &gbd ->Gamut[theta][alpha];
 353 }
 354
 355 // Add a point to gamut descriptor. Point to add is in Lab color space.
 356 // GBD is centered on a=b=0 and L*=50
 357 cmsBool CMSEXPORT cmsGDBAddPoint(cmsHANDLE hGBD, const cmsCIELab* Lab)
 358 {
 359     cmsGDB* gbd = (cmsGDB*) hGBD;
 360     cmsGDBPoint* ptr;
 361     cmsSpherical sp;
 362
 363
 364     // Get pointer to the sector
 365     ptr = GetPoint(gbd, Lab, &sp);
 366     if (ptr == NULL) return FALSE;
 367
 368     // If no samples at this sector, add it
 369     if (ptr ->Type == GP_EMPTY) {
 370
 371         ptr -> Type = GP_SPECIFIED;
 372         ptr -> p    = sp;
 373     }
 374     else {
 375
 376
 377         // Substitute only if radius is greater
 378         if (sp.r > ptr -> p.r) {
 379
 380                 ptr -> Type = GP_SPECIFIED;
 381                 ptr -> p    = sp;
 382         }
 383     }
 384
 385     return TRUE;
 386 }
 387
 388 // Check if a given point falls inside gamut
 389 cmsBool CMSEXPORT cmsGDBCheckPoint(cmsHANDLE hGBD, const cmsCIELab* Lab)
 390 {
 391     cmsGDB* gbd = (cmsGDB*) hGBD;
 392     cmsGDBPoint* ptr;
 393     cmsSpherical sp;
 394
 395     // Get pointer to the sector
 396     ptr = GetPoint(gbd, Lab, &sp);
 397     if (ptr == NULL) return FALSE;
 398
 399     // If no samples at this sector, return no data
 400     if (ptr ->Type == GP_EMPTY) return FALSE;
 401
 402     // In gamut only if radius is greater
 403
 404     return (sp.r <= ptr -> p.r);
 405 }
 406
 407 // -----------------------------------------------------------------------------------------------------------------------
 408
 409 // Find near sectors. The list of sectors found is returned on Close[].
 410 // The function returns the number of sectors as well.
 411
 412 // 24   9  10  11  12
 413 // 23   8   1   2  13
 414 // 22   7   *   3  14
 415 // 21   6   5   4  15
 416 // 20  19  18  17  16
 417 //
 418 // Those are the relative movements
 419 // {-2,-2}, {-1, -2}, {0, -2}, {+1, -2}, {+2,  -2},
 420 // {-2,-1}, {-1, -1}, {0, -1}, {+1, -1}, {+2,  -1},
 421 // {-2, 0}, {-1,  0}, {0,  0}, {+1,  0}, {+2,   0},
 422 // {-2,+1}, {-1, +1}, {0, +1}, {+1,  +1}, {+2,  +1},
 423 // {-2,+2}, {-1, +2}, {0, +2}, {+1,  +2}, {+2,  +2}};
 424
 425
 426 static
 427 const struct _spiral {
 428
 429     int AdvX, AdvY;
 430
 431     } Spiral[] = { {0,  -1}, {+1, -1}, {+1,  0}, {+1, +1}, {0,  +1}, {-1, +1},
 432                    {-1,  0}, {-1, -1}, {-1, -2}, {0,  -2}, {+1, -2}, {+2, -2},
 433                    {+2, -1}, {+2,  0}, {+2, +1}, {+2, +2}, {+1, +2}, {0,  +2},
 434                    {-1, +2}, {-2, +2}, {-2, +1}, {-2, 0},  {-2, -1}, {-2, -2} };
 435
 436 #define NSTEPS (sizeof(Spiral) / sizeof(struct _spiral))
 437
 438 static
 439 int FindNearSectors(cmsGDB* gbd, int alpha, int theta, cmsGDBPoint* Close[])
 440 {
 441     int nSectors = 0;
 442     int a, t;
 443     cmsUInt32Number i;
 444     cmsGDBPoint* pt;
 445
 446     for (i=0; i < NSTEPS; i++) {
 447
 448         a = alpha + Spiral[i].AdvX;
 449         t = theta + Spiral[i].AdvY;
 450
 451         // Cycle at the end
 452         a %= SECTORS;
 453         t %= SECTORS;
 454
 455         // Cycle at the begin
 456         if (a < 0) a = SECTORS + a;
 457         if (t < 0) t = SECTORS + t;
 458
 459         pt = &gbd ->Gamut[t][a];
 460
 461         if (pt -> Type != GP_EMPTY) {
 462
 463             Close[nSectors++] = pt;
 464         }
 465     }
 466
 467     return nSectors;
 468 }
 469
 470
 471 // Interpolate a missing sector. Method identifies whatever this is top, bottom or mid
 472 static
 473 cmsBool InterpolateMissingSector(cmsGDB* gbd, int alpha, int theta)
 474 {
 475     cmsSpherical sp;
 476     cmsVEC3 Lab;
 477     cmsVEC3 Centre;
 478     cmsLine ray;
 479     int nCloseSectors;
 480     cmsGDBPoint* Close[NSTEPS + 1];
 481     cmsSpherical closel, templ;
 482     cmsLine edge;
 483     int k, m;
 484
 485     // Is that point already specified?
 486     if (gbd ->Gamut[theta][alpha].Type != GP_EMPTY) return TRUE;
 487
 488     // Fill close points
 489     nCloseSectors = FindNearSectors(gbd, alpha, theta, Close);
 490
 491
 492     // Find a central point on the sector
 493     sp.alpha = (cmsFloat64Number) ((alpha + 0.5) * 360.0) / (SECTORS);
 494     sp.theta = (cmsFloat64Number) ((theta + 0.5) * 180.0) / (SECTORS);
 495     sp.r     = 50.0;
 496
 497     // Convert to Cartesian
 498     ToCartesian(&Lab, &sp);
 499
 500     // Create a ray line from centre to this point
 501     _cmsVEC3init(&Centre, 50.0, 0, 0);
 502     LineOf2Points(&ray, &Lab, &Centre);
 503
 504     // For all close sectors
 505     closel.r = 0.0;
 506     closel.alpha = 0;
 507     closel.theta = 0;
 508
 509     for (k=0; k < nCloseSectors; k++) {
 510
 511         for(m = k+1; m < nCloseSectors; m++) {
 512
 513             cmsVEC3 temp, a1, a2;
 514
 515             // A line from sector to sector
 516             ToCartesian(&a1, &Close[k]->p);
 517             ToCartesian(&a2, &Close[m]->p);
 518
 519             LineOf2Points(&edge, &a1, &a2);
 520
 521             // Find a line
 522             ClosestLineToLine(&temp, &ray, &edge);
 523
 524             // Convert to spherical
 525             ToSpherical(&templ, &temp);
 526
 527
 528             if ( templ.r > closel.r &&
 529                  templ.theta >= (theta*180.0/SECTORS) &&
 530                  templ.theta <= ((theta+1)*180.0/SECTORS) &&
 531                  templ.alpha >= (alpha*360.0/SECTORS) &&
 532                  templ.alpha <= ((alpha+1)*360.0/SECTORS)) {
 533
 534                 closel = templ;
 535             }
 536         }
 537     }
 538
 539     gbd ->Gamut[theta][alpha].p = closel;
 540     gbd ->Gamut[theta][alpha].Type = GP_MODELED;
 541
 542     return TRUE;
 543
 544 }
 545
 546
 547 // Interpolate missing parts. The algorithm fist computes slices at
 548 // theta=0 and theta=Max.
 549 cmsBool CMSEXPORT cmsGDBCompute(cmsHANDLE hGBD, cmsUInt32Number dwFlags)
 550 {
 551     int alpha, theta;
 552     cmsGDB* gbd = (cmsGDB*) hGBD;
 553
 554     _cmsAssert(hGBD != NULL);
 555
 556     // Interpolate black
 557     for (alpha = 0; alpha < SECTORS; alpha++) {
 558
 559         if (!InterpolateMissingSector(gbd, alpha, 0)) return FALSE;
 560     }
 561
 562     // Interpolate white
 563     for (alpha = 0; alpha < SECTORS; alpha++) {
 564
 565         if (!InterpolateMissingSector(gbd, alpha, SECTORS-1)) return FALSE;
 566     }
 567
 568
 569     // Interpolate Mid
 570     for (theta = 1; theta < SECTORS; theta++) {
 571         for (alpha = 0; alpha < SECTORS; alpha++) {
 572
 573             if (!InterpolateMissingSector(gbd, alpha, theta)) return FALSE;
 574         }
 575     }
 576
 577     // Done
 578     return TRUE;
 579
 580     cmsUNUSED_PARAMETER(dwFlags);
 581 }
 582
 583
 584
 585
 586 // --------------------------------------------------------------------------------------------------------
 587
 588 // Great for debug, but not suitable for real use
 589
 590 #if 0
 591 cmsBool cmsGBDdumpVRML(cmsHANDLE hGBD, const char* fname)
 592 {
 593     FILE* fp;
 594     int   i, j;
 595     cmsGDB* gbd = (cmsGDB*) hGBD;
 596     cmsGDBPoint* pt;
 597
 598     fp = fopen (fname, "wt");
 599     if (fp == NULL)
 600         return FALSE;
 601
 602     fprintf (fp, "#VRML V2.0 utf8\n");
 603
 604     // set the viewing orientation and distance
 605     fprintf (fp, "DEF CamTest Group {\n");
 606     fprintf (fp, "\tchildren [\n");
 607     fprintf (fp, "\t\tDEF Cameras Group {\n");
 608     fprintf (fp, "\t\t\tchildren [\n");
 609     fprintf (fp, "\t\t\t\tDEF DefaultView Viewpoint {\n");
 610     fprintf (fp, "\t\t\t\t\tposition 0 0 340\n");
 611     fprintf (fp, "\t\t\t\t\torientation 0 0 1 0\n");
 612     fprintf (fp, "\t\t\t\t\tdescription \"default view\"\n");
 613     fprintf (fp, "\t\t\t\t}\n");
 614     fprintf (fp, "\t\t\t]\n");
 615     fprintf (fp, "\t\t},\n");
 616     fprintf (fp, "\t]\n");
 617     fprintf (fp, "}\n");
 618
 619     // Output the background stuff
 620     fprintf (fp, "Background {\n");
 621     fprintf (fp, "\tskyColor [\n");
 622     fprintf (fp, "\t\t.5 .5 .5\n");
 623     fprintf (fp, "\t]\n");
 624     fprintf (fp, "}\n");
 625
 626     // Output the shape stuff
 627     fprintf (fp, "Transform {\n");
 628     fprintf (fp, "\tscale .3 .3 .3\n");
 629     fprintf (fp, "\tchildren [\n");
 630
 631     // Draw the axes as a shape:
 632     fprintf (fp, "\t\tShape {\n");
 633     fprintf (fp, "\t\t\tappearance Appearance {\n");
 634     fprintf (fp, "\t\t\t\tmaterial Material {\n");
 635     fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
 636     fprintf (fp, "\t\t\t\t\temissiveColor 1.0 1.0 1.0\n");
 637     fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
 638     fprintf (fp, "\t\t\t\t}\n");
 639     fprintf (fp, "\t\t\t}\n");
 640     fprintf (fp, "\t\t\tgeometry IndexedLineSet {\n");
 641     fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
 642     fprintf (fp, "\t\t\t\t\tpoint [\n");
 643     fprintf (fp, "\t\t\t\t\t0.0 0.0 0.0,\n");
 644     fprintf (fp, "\t\t\t\t\t%f 0.0 0.0,\n",  255.0);
 645     fprintf (fp, "\t\t\t\t\t0.0 %f 0.0,\n",  255.0);
 646     fprintf (fp, "\t\t\t\t\t0.0 0.0 %f]\n",  255.0);
 647     fprintf (fp, "\t\t\t\t}\n");
 648     fprintf (fp, "\t\t\t\tcoordIndex [\n");
 649     fprintf (fp, "\t\t\t\t\t0, 1, -1\n");
 650     fprintf (fp, "\t\t\t\t\t0, 2, -1\n");
 651     fprintf (fp, "\t\t\t\t\t0, 3, -1]\n");
 652     fprintf (fp, "\t\t\t}\n");
 653     fprintf (fp, "\t\t}\n");
 654
 655
 656     fprintf (fp, "\t\tShape {\n");
 657     fprintf (fp, "\t\t\tappearance Appearance {\n");
 658     fprintf (fp, "\t\t\t\tmaterial Material {\n");
 659     fprintf (fp, "\t\t\t\t\tdiffuseColor 0 0.8 0\n");
 660     fprintf (fp, "\t\t\t\t\temissiveColor 1 1 1\n");
 661     fprintf (fp, "\t\t\t\t\tshininess 0.8\n");
 662     fprintf (fp, "\t\t\t\t}\n");
 663     fprintf (fp, "\t\t\t}\n");
 664     fprintf (fp, "\t\t\tgeometry PointSet {\n");
 665
 666     // fill in the points here
 667     fprintf (fp, "\t\t\t\tcoord Coordinate {\n");
 668     fprintf (fp, "\t\t\t\t\tpoint [\n");
 669
 670     // We need to transverse all gamut hull.
 671     for (i=0; i < SECTORS; i++)
 672         for (j=0; j < SECTORS; j++) {
 673
 674             cmsVEC3 v;
 675
 676             pt = &gbd ->Gamut[i][j];
 677             ToCartesian(&v, &pt ->p);
 678
 679             fprintf (fp, "\t\t\t\t\t%g %g %g", v.n[0]+50, v.n[1], v.n[2]);
 680
 681             if ((j == SECTORS - 1) && (i == SECTORS - 1))
 682                 fprintf (fp, "]\n");
 683             else
 684                 fprintf (fp, ",\n");
 685
 686         }
 687
 688         fprintf (fp, "\t\t\t\t}\n");
 689
 690
 691
 692     // fill in the face colors
 693     fprintf (fp, "\t\t\t\tcolor Color {\n");
 694     fprintf (fp, "\t\t\t\t\tcolor [\n");
 695
 696     for (i=0; i < SECTORS; i++)
 697         for (j=0; j < SECTORS; j++) {
 698
 699            cmsVEC3 v;
 700
 701             pt = &gbd ->Gamut[i][j];
 702
 703
 704             ToCartesian(&v, &pt ->p);
 705
 706
 707         if (pt ->Type == GP_EMPTY)
 708             fprintf (fp, "\t\t\t\t\t%g %g %g", 0.0, 0.0, 0.0);
 709         else
 710             if (pt ->Type == GP_MODELED)
 711                 fprintf (fp, "\t\t\t\t\t%g %g %g", 1.0, .5, .5);
 712             else {
 713                 fprintf (fp, "\t\t\t\t\t%g %g %g", 1.0, 1.0, 1.0);
 714
 715             }
 716
 717         if ((j == SECTORS - 1) && (i == SECTORS - 1))
 718                 fprintf (fp, "]\n");
 719             else
 720                 fprintf (fp, ",\n");
 721     }
 722     fprintf (fp, "\t\t\t}\n");
 723
 724
 725     fprintf (fp, "\t\t\t}\n");
 726     fprintf (fp, "\t\t}\n");
 727     fprintf (fp, "\t]\n");
 728     fprintf (fp, "}\n");
 729
 730     fclose (fp);
 731
 732     return TRUE;
 733 }
 734 #endif
 735