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