SkyDomeImage.cpp

Go to the documentation of this file.

/*
 * SkyDomeImage.cpp
 *
 * Part of Fly! Legacy project
 *
 * Copyright 2003 Chris Wallace
 *
 * Fly! Legacy is free software; you can redistribute it and/or modify
 *   it under the terms of the GNU General Public License as published by
 *   the Free Software Foundation; either version 2 of the License, or
 *   (at your option) any later version.
 *
 * Fly! Legacy is distributed in the hope that it will be useful,
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *   GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 *   along with Fly! Legacy; if not, write to the Free Software
 *    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 *
 */

#include "../Include/Globals.h"
#include "../Include/Sky.h"
#include "../Include/Utility.h"
#include "../Include/Ui.h"
#include <plib/ssg.h>

using namespace std;

//
// SSG pre-draw callback disables z-buffering and fog
//
static int skydomeimage_predraw (ssgEntity *e)
{
  ssgLeaf *f = (ssgLeaf*) e;
  if (f->hasState()) f->getState()->apply();

  if (globals->settings->skydomeWireframe) {
    glPushAttrib (GL_DEPTH_BUFFER_BIT | GL_FOG_BIT | GL_POLYGON_BIT);
    glPolygonMode( GL_FRONT_AND_BACK, GL_LINE );
  } else {
    glPushAttrib (GL_DEPTH_BUFFER_BIT | GL_FOG_BIT);
  }
  glDisable (GL_DEPTH_TEST);
  glDisable (GL_FOG);

    return true;
}


//
// SSG post-draw callback restores original state of
//   z-buffering and fog
//
static int skydomeimage_postdraw (ssgEntity *e)
{
  glPopAttrib ();

  return true;
}


/*
 * Perez sky dome
 */

//
// Perez function
//
float PerezFunction (float theta, float gamma, const SPerezCoefficients &p)
{
  // Check for indeterminate case of cos(theta)=0, where exp(p.B/cos(theta))=0,
  //   assuming p.B is negative
  float f1;
  double epsilon = 1.0e-10;
  if (cos(theta) < epsilon) {
    f1 = 1;
  } else {
    f1 = 1 + (p.A * exp (p.B / cos(theta)));
  }
  float cosgamma = cos(gamma);
  float cos2gamma = cosgamma * cosgamma;
  float f2 = 1 + (p.C * exp (p.D * gamma)) + (p.E * cos2gamma);

  return f1 * f2;
}

           
//
// Calculate Perez distribution coefficients
//
// Inputs:
//  T   Turbidity
//  parm  Coefficient parameters
//
// Outputs:
//  coeff Coefficients
//
void CalcPerezCoefficients (float T,
                const SCoefficientParameters &parm,
                SPerezCoefficients &coeff)
{
  coeff.A = parm.A[1] * T + parm.A[0];
  coeff.B = parm.B[1] * T + parm.B[0];
  coeff.C = parm.C[1] * T + parm.C[0];
  coeff.D = parm.D[1] * T + parm.D[0];
  coeff.E = parm.E[1] * T + parm.E[0];
}

          
//
// Calculate the Y (luminance component) for the zenith
//
// Input parameters:
//  theta Solar zenith angle in radians
//  T   Turbidity
//
// Return code:
//  float Y (luminance)
//typedef struct {
//  float   A, B, C, D;
//} SLuminanceZenith;
//
static float ZenithLuminance (float theta, float T, const SLuminanceZenith &p)
{
  float chi = ((4.0 / 9.0) - (T / 120.0)) * (SG_PI - 2.0 * theta);
  float Yz = ((p.A * T) + p.B) * tan(chi) + (p.C * T) + p.D;
  Yz *= 1000.0f;

  return Yz;
}


//
// Calculate chrominance (x or y) component for the zenith
//
// Input parameters:
//  theta Solar zenith angle in radians
//  T   Turbidity
//  p   Zenith chrominance parameter matrix
//
// Return value:
//  float Chrominance component
//
static float ZenithChrominance (float theta, float T, const SChromaZenith &p)
{
  float rc = 0.0;
  float T2 = T * T;
  float theta2 = theta * theta;
  float theta3 = theta * theta * theta;
  rc = 
    (p.t2[3] * theta3 + p.t2[2] * theta2 + p.t2[1] * theta + p.t2[0]) * T2 +
    (p.t1[3] * theta3 + p.t1[2] * theta2 + p.t1[1] * theta + p.t1[0]) * T +
    (p.t0[3] * theta3 + p.t0[2] * theta2 + p.t0[1] * theta + p.t0[0]);

  return rc;
}


//
// Convert CIE Yxy colour value to RGB.  This is done in two steps, first to
//   CIE XYZ format, then from XYZ to RGB.  White point D65 is used.
//
// Input parameters:
//  Y CIE Yxy luminance Y value   [0..100]
//  x CIE Yxy chrominance x value   [0..1]
//  y CIE Yxy chrominance y value   [0..1]
//
static void YxyToRGB (float Y, float x, float y, float* r, float* g, float* b)
{
  // First convert to CIE XYZ format
  float X = (x / y) * Y;
  float Z = ((1 - x - y) / y) * Y;

  // Convert XYZ to RGB using D65 white point
  float factor = 1.0f - exp(-(1.0f / 15000.0f));
  float R = factor * ( 3.240479 * X + -1.537150 * Y + -0.498535 * Z);
  float G = factor * (-0.969256 * X +  1.875992 * Y +  0.041556 * Z);
  float B = factor * ( 0.055648 * X + -0.204043 * Y +  1.057311 * Z);

  if (R < 0.0f) R = 0.0f;
  if (R > 1.0f) R = 1.0f;

  if (G < 0.0f) G = 0.0f;
  if (G > 1.0f) G = 1.0f;
  
  if (B < 0.0f) B = 0.0f;
  if (B > 1.0f) B = 1.0f;

  // Gamma-adjust
  float gamma = 1.8;
  float invgamma = 1.0f / gamma;
  *r = pow (R, invgamma);
  *g = pow (G, invgamma);
  *b = pow (B, invgamma);
}


static void InitPerezParameters (SPerezParameters* p,
                 CIniFile *inisky,
                 const char* name)
{
  // Initialize to defaults
  memset (p, 0, sizeof(SPerezParameters));

  // Initialize name
  strcpy (p->name, name);

  // Initialize zenith parameters
  inisky->Get (name, "zenith_Y_A", &(p->Y_zenith.A));
  inisky->Get (name, "zenith_Y_B", &(p->Y_zenith.B));
  inisky->Get (name, "zenith_Y_C", &(p->Y_zenith.C));
  inisky->Get (name, "zenith_Y_D", &(p->Y_zenith.D));

  inisky->Get (name, "zenith_x_theta3_T2", &(p->x_zenith.t2[3]));
  inisky->Get (name, "zenith_x_theta3_T1", &(p->x_zenith.t1[3]));
  inisky->Get (name, "zenith_x_theta3_T0", &(p->x_zenith.t0[3]));
  inisky->Get (name, "zenith_x_theta2_T2", &(p->x_zenith.t2[2]));
  inisky->Get (name, "zenith_x_theta2_T1", &(p->x_zenith.t1[2]));
  inisky->Get (name, "zenith_x_theta2_T0", &(p->x_zenith.t0[2]));
  inisky->Get (name, "zenith_x_theta1_T2", &(p->x_zenith.t2[1]));
  inisky->Get (name, "zenith_x_theta1_T1", &(p->x_zenith.t1[1]));
  inisky->Get (name, "zenith_x_theta1_T0", &(p->x_zenith.t0[1]));
  inisky->Get (name, "zenith_x_theta0_T2", &(p->x_zenith.t2[0]));
  inisky->Get (name, "zenith_x_theta0_T1", &(p->x_zenith.t1[0]));
  inisky->Get (name, "zenith_x_theta0_T0", &(p->x_zenith.t0[0]));

  inisky->Get (name, "zenith_y_theta3_T2", &(p->y_zenith.t2[3]));
  inisky->Get (name, "zenith_y_theta3_T1", &(p->y_zenith.t1[3]));
  inisky->Get (name, "zenith_y_theta3_T0", &(p->y_zenith.t0[3]));
  inisky->Get (name, "zenith_y_theta2_T2", &(p->y_zenith.t2[2]));
  inisky->Get (name, "zenith_y_theta2_T1", &(p->y_zenith.t1[2]));
  inisky->Get (name, "zenith_y_theta2_T0", &(p->y_zenith.t0[2]));
  inisky->Get (name, "zenith_y_theta1_T2", &(p->y_zenith.t2[1]));
  inisky->Get (name, "zenith_y_theta1_T1", &(p->y_zenith.t1[1]));
  inisky->Get (name, "zenith_y_theta1_T0", &(p->y_zenith.t0[1]));
  inisky->Get (name, "zenith_y_theta0_T2", &(p->y_zenith.t2[0]));
  inisky->Get (name, "zenith_y_theta0_T1", &(p->y_zenith.t1[0]));
  inisky->Get (name, "zenith_y_theta0_T0", &(p->y_zenith.t0[0]));

  // Initialize curve parameters
  inisky->Get (name, "perez_Y_A_T1", &(p->Y_curve.A[1]));
  inisky->Get (name, "perez_Y_A_T0", &(p->Y_curve.A[0]));
  inisky->Get (name, "perez_Y_B_T1", &(p->Y_curve.B[1]));
  inisky->Get (name, "perez_Y_B_T0", &(p->Y_curve.B[0]));
  inisky->Get (name, "perez_Y_C_T1", &(p->Y_curve.C[1]));
  inisky->Get (name, "perez_Y_C_T0", &(p->Y_curve.C[0]));
  inisky->Get (name, "perez_Y_D_T1", &(p->Y_curve.D[1]));
  inisky->Get (name, "perez_Y_D_T0", &(p->Y_curve.D[0]));
  inisky->Get (name, "perez_Y_E_T1", &(p->Y_curve.E[1]));
  inisky->Get (name, "perez_Y_E_T0", &(p->Y_curve.E[0]));

  inisky->Get (name, "perez_x_A_T1", &(p->x_curve.A[1]));
  inisky->Get (name, "perez_x_A_T0", &(p->x_curve.A[0]));
  inisky->Get (name, "perez_x_B_T1", &(p->x_curve.B[1]));
  inisky->Get (name, "perez_x_B_T0", &(p->x_curve.B[0]));
  inisky->Get (name, "perez_x_C_T1", &(p->x_curve.C[1]));
  inisky->Get (name, "perez_x_C_T0", &(p->x_curve.C[0]));
  inisky->Get (name, "perez_x_D_T1", &(p->x_curve.D[1]));
  inisky->Get (name, "perez_x_D_T0", &(p->x_curve.D[0]));
  inisky->Get (name, "perez_x_E_T1", &(p->x_curve.E[1]));
  inisky->Get (name, "perez_x_E_T0", &(p->x_curve.E[0]));

  inisky->Get (name, "perez_y_A_T1", &(p->y_curve.A[1]));
  inisky->Get (name, "perez_y_A_T0", &(p->y_curve.A[0]));
  inisky->Get (name, "perez_y_B_T1", &(p->y_curve.B[1]));
  inisky->Get (name, "perez_y_B_T0", &(p->y_curve.B[0]));
  inisky->Get (name, "perez_y_C_T1", &(p->y_curve.C[1]));
  inisky->Get (name, "perez_y_C_T0", &(p->y_curve.C[0]));
  inisky->Get (name, "perez_y_D_T1", &(p->y_curve.D[1]));
  inisky->Get (name, "perez_y_D_T0", &(p->y_curve.D[0]));
  inisky->Get (name, "perez_y_E_T1", &(p->y_curve.E[1]));
  inisky->Get (name, "perez_y_E_T0", &(p->y_curve.E[0]));
}


SPerezParameters DefaultPerezParameters =
{
  "Default",

  // Y distribution
  {
    {-1.4630, +0.1787}, // A
    {+0.4275, -0.3554}, // B
    {+5.3251, -0.0227}, // C
    {-2.5771, +0.1206}, // D
    {+0.3703, -0.0670}, // E
  },

  // x distribution
  {
    {-0.2592, -0.0193}, // A
    {+0.0008, -0.0665}, // B
    {+0.2125, -0.0004}, // C
    {-0.8989, -0.0641}, // D
    {+0.0452, -0.0033}, // E
  },

  // y distribution
  {
    {-0.2608, -0.0167}, // A
    {+0.0092, -0.0950}, // B
    {+0.2102, -0.0079}, // C
    {-1.6537, -0.0441}, // D
    {+0.0529, -0.0109}, // E
  },

  // Y zenith
  {
    +4.0453,    // A
    -4.9710,    // B
    -0.2155,    // C
    +2.4192,    // D
  },

  // x zenith
  {
    //   1       theta    theta^2   theta^3
    {+0.00000, +0.00209, -0.00375, +0.00166}, // T^2
    {+0.00394, -0.03202, +0.06377, -0.02903}, // T
    {+0.25886, +0.06052, -0.21196, +0.11693}, // 1
  },

  // y zenith
  {
    //   1       theta    theta^2   theta^3
    {+0.00000, +0.00317, -0.00610, +0.00275}, // T^2
    {+0.00516, -0.04153, +0.08970, -0.04214}, // T
    {+0.26688, +0.06670, -0.26756, +0.15346}, // 1
  }
};

CSkyDomeImage::CSkyDomeImage (double distance)
{
  int i, j;
  
  // Add default Perez parameters
  SPerezParameters *pDefault = new SPerezParameters;
  memcpy (pDefault, &DefaultPerezParameters, sizeof(SPerezParameters));
  perez[DefaultPerezParameters.name] = pDefault;

  //
  // Load Perez parameter sets from SKIES.INI
  //
  CIniFile *inisky = new CIniFile;
  if (inisky->Load ("System\\Skies.ini") != 0) {
    // SKIES.INI could be read; parse parameters from INI settings
    int nSkies = inisky->GetNumSections();
    for (i=0; i<nSkies; i++) {
      SPerezParameters *p = new SPerezParameters;
      memcpy (p, &DefaultPerezParameters, sizeof(SPerezParameters));
      char* skyname = inisky->GetSectionName (i);
      InitPerezParameters (p, inisky, skyname);
      perez[skyname] = p;
    }
  }
  delete inisky;

  // Get default sky name from INI settings
  perezCurrent = &DefaultPerezParameters;
  char skyname[64];
  strcpy (skyname, "");
  GetIniString ("Graphics", "skyDefaultName", skyname, 64);
  if (strlen (skyname) > 0) {
    // Parameter was specified, search available skies for the default
    map<string,SPerezParameters*>::iterator i = perez.find(skyname);
    if (i != perez.end()) {
      perezCurrent = i->second;
    }
  }

  //
  // Generate raw vertex data
  //
  // The sky dome vertex raw data is made up of a single zenith vertex, plus
  //   SKYDOME_STACKS "rings" of SKYDOME_SLICES vertices in each ring.  The top of
  //   the sky dome is constructed of a triangle fan between the zenith
  //   vertex and each of the vertices in the top ring.  The remainder of the
  //   dome is constructed of triangle strips spanning each successive vertex ring.
  //

  // Zenith vertex is easy
  sgSetVec3 (vtx_zenith, 0, 0, distance);

  // Iterate over all vertices of all rings
  double dPhi = (SGD_PI / 2) / SKYDOME_STACKS;
  double dTheta = (2.0 * SGD_PI) / SKYDOME_SLICES;
  for (i=0; i<SKYDOME_STACKS; i++) {
    // Calculate elevation angle phi for all vertices in this ring
    double phi = (SGD_PI / 2) - ((i+1) * dPhi);
    double ring_radius = distance * cos(phi);
    double ring_elevation = distance * sin(phi);

    for (j=0; j<SKYDOME_SLICES; j++) {
      // Calculate azimuth angle theta
      double theta = j * dTheta;
      sgSetVec3 (vtx_ring[i][j],
                 cos(theta) * ring_radius,
                 sin(theta) * ring_radius,
                 ring_elevation);
    }
  }

  //
  // Create top-level SSG object to be the parent of all leaf objects
  //
  top = new ssgTransform;
  top->setName ("SkyDomeImage");

  //  Define a dummy colour for construction of the skydome.
  //  The Repaint method will set the proper colours for each vertex
  sgVec4 c;
  sgSetVec4 (c, 0.39, 0.50, 0.75, 1.0);

  // Define a "special" colour that can be used to highlight particular
  //   areas of the skydome for testing
  sgVec4 special;
  sgSetVec4 (special, 1.0, 0.0, 0.0, 1.0);

  //
  // Create SSG leaf for zenith triangle fan
  //
  ssgVertexArray *va = new ssgVertexArray(SKYDOME_SLICES + 2);
  ssgColourArray *rgba = new ssgColourArray(SKYDOME_SLICES + 2);

  // Generate zenith disk
  va->add (vtx_zenith);
  rgba->add (c);
  for (i=0; i<SKYDOME_SLICES; i++) {
    va->add (vtx_ring[0][i]);
    rgba->add (c);
  }
  va->add (vtx_ring[0][0]);
  rgba->add (c);

  // Create SSG state data to be applied to zenith fan
  ssgSimpleState *state = new ssgSimpleState ();
  state->setShadeModel (GL_SMOOTH);
  state->disable (GL_LIGHTING);
  state->disable (GL_CULL_FACE);
  state->disable (GL_TEXTURE_2D);
  state->enable (GL_BLEND);
  state->disable (GL_ALPHA_TEST);
  state->enable (GL_COLOR_MATERIAL);
  state->setColourMaterial (GL_AMBIENT_AND_DIFFUSE);
  state->setMaterial (GL_EMISSION, 0, 0, 0, 1);
  state->setMaterial (GL_SPECULAR, 0, 0, 0, 1);

  // Create SSG leaf, assign name, state and callbacks
  ssgVtxTable *vtZenith = new ssgVtxTable (GL_TRIANGLE_FAN, va, NULL, NULL, rgba);
  vtZenith->setName ("Zenith");
  vtZenith->setState (state);
  vtZenith->setCallback (SSG_CALLBACK_PREDRAW,  skydomeimage_predraw);
  vtZenith->setCallback (SSG_CALLBACK_POSTDRAW, skydomeimage_postdraw);

  // Attach leaf to top-level transform
  top->addKid (vtZenith);

  //
  // Create SSG leaf objects for each stack
  //
  for (i=0; i<SKYDOME_STACKS-1; i++) {
    // Allocate vertex and colour arrays
    va = new ssgVertexArray((SKYDOME_SLICES * 2) + 2);
    rgba = new ssgColourArray((SKYDOME_SLICES * 2) + 2);

    // Generate arrays
    for (j=0; j<SKYDOME_SLICES; j++) {
      va->add (vtx_ring[i][j]);
      rgba->add (c);
      va->add (vtx_ring[i+1][j]);
      rgba->add (c);
    }
    va->add (vtx_ring[i][0]);
    rgba->add (c);
    va->add (vtx_ring[i+1][0]);
    rgba->add (c);

    // Create SSG state data to be applied to this slice
    ssgSimpleState *state = new ssgSimpleState ();
    state->setShadeModel (GL_SMOOTH);
    state->disable (GL_LIGHTING);
    state->disable (GL_CULL_FACE);
    state->disable (GL_TEXTURE_2D);
    state->enable (GL_BLEND);
    state->disable (GL_ALPHA_TEST);
    state->enable (GL_COLOR_MATERIAL);
    state->setColourMaterial (GL_AMBIENT_AND_DIFFUSE);
    state->setMaterial (GL_EMISSION, 0, 0, 0, 1);
    state->setMaterial (GL_SPECULAR, 0, 0, 0, 1);

    // Create SSG leaf, assign name, state and callbacks
    ssgVtxTable *vtab = new ssgVtxTable (GL_TRIANGLE_STRIP, va, NULL, NULL, rgba);
    char name[64];
    sprintf (name, "Stack %d", i);
    vtab->setName (name);
    vtab->setState (state);
    vtab->setCallback (SSG_CALLBACK_PREDRAW,  skydomeimage_predraw);
    vtab->setCallback (SSG_CALLBACK_POSTDRAW, skydomeimage_postdraw);

    // Attach leaf to top-level transform
    top->addKid (vtab);
  }

  // Initialize repaint solar zenith angle
  prevTheta = 0;
}


CSkyDomeImage::~CSkyDomeImage (void)
{
  map<string,SPerezParameters*>::iterator i;
  for (i=perez.begin(); i!=perez.end(); i++) {
    SPerezParameters *p = i->second;
    delete p;
  }
}


//
// Return the top-level SSG entity for the sky dome
//
ssgEntity*      CSkyDomeImage::GetSSGEntity (void)
{
  return top;
}


//
// Repaint the sky dome
//
// Input parameters:
//  solTheta    Zenith angle of the sun, radians
//  solPhi      Azimuth angle of the sun, in radians W of S
//
void CSkyDomeImage::Repaint (double solTheta, double solPhi)
{
  int i, j;
  double thetaThreshold = DegToRad (0.01);

  if (fabs(solTheta - prevTheta) < thetaThreshold) {
    // Sun angle has not changed significantly; nothing to do
    return;
  }
  prevTheta = solTheta;

  //
  // Perez sky shading is intended for day/dusk only.  Set dome colour to
  //   all black, these values will be overwritten with Perez-derived
  //   colours if the sun is above astronomical night
  float nightTheta = DegToRad (108.0f);   // -18 deg solar elevation

  //
  // Step 1 : Calculate new sky dome colours
  //

  // TEMP Assume a reasonable value for turbidity, should be based on visibility

  if (solTheta > nightTheta) {
    // Sun is below level of astronomical night...set all vertices to black
    sgSetVec3 (rgb_zenith, 0, 0, 0);
    for (i=0; i<SKYDOME_STACKS; i++) {
      for (j=0; j<SKYDOME_SLICES; j++) {
        sgSetVec3 (rgb_dome[i][j], 0, 0, 0);
      }
    }
  } else {
    // Sun is above level of astronomical night...use Perez equations to compute
    //   zenith and dome vertex colours
    float turbidity = 4.0;
    float r = 0, g = 0, b = 0;

    // Calculate distribution coefficients
    SPerezCoefficients cY, cx, cy;
    CalcPerezCoefficients (turbidity, perezCurrent->Y_curve, cY);
    CalcPerezCoefficients (turbidity, perezCurrent->x_curve, cx);
    CalcPerezCoefficients (turbidity, perezCurrent->y_curve, cy);

    // Calculate zenith luminance and chrominance
    float Yz = ZenithLuminance (solTheta, turbidity, perezCurrent->Y_zenith);
    float xz = ZenithChrominance (solTheta, turbidity, perezCurrent->x_zenith);
    float yz = ZenithChrominance (solTheta, turbidity, perezCurrent->y_zenith);

    // Calculate Perez function zenith values for later use
    float Ypz = PerezFunction (0, solTheta, cY);
    float xpz = PerezFunction (0, solTheta, cx);
    float ypz = PerezFunction (0, solTheta, cy);

    // Convert to RGB
    YxyToRGB (Yz, xz, yz, &r, &g, &b);

    // Set zenith colour
    sgSetVec3 (rgb_zenith, r, g, b);

// DEBUG
//  char debug[80];
//  sprintf (debug, "Yz = %7.3f xz = %7.3f yz = %7.3f   rz=%7.3f gz=%7.3f bz=%7.3f",
//    Yz, xz, yz, rz, gz, bz);
//  DrawNoticeToUser (debug, 1);

    // Pre-calculate theta and phi steps per dome stack/slice
    float dTheta = (SG_PI / 2.0f) / SKYDOME_STACKS;
    float dPhi = (SG_PI * 2.0f) / SKYDOME_SLICES;

    // Since sky dome 0th vertex is aligned with sun, the sun vector by
    //   definition is at phi=0, regardless of the actual phi angle
    //   cos(0) = 1, sin(0) = 0
    sgVec3 vSun;
    sgSetVec3 (vSun, sin(solTheta), 0, cos(solTheta));

    // Iterate over all dome vertices
    for (int i=0; i<SKYDOME_STACKS; i++) {
      float domeTheta = (i+1) * dTheta;
      for (int j=0; j<SKYDOME_SLICES; j++) {
        float domePhi = j * dPhi;

        // Calcualate gamma angle for this point on the dome
        sgVec3 vDome;
        sgSetVec3 (vDome,
                   cos(domePhi) * sin(domeTheta),
                   sin(domePhi) * sin(domeTheta),
                   cos(domeTheta));
        float gamma = DegToRad (sgAngleBetweenVec3 (vSun, vDome));

        // Calculate dome luminance and chrominance
        float Yd = Yz * PerezFunction (domeTheta, gamma, cY) / Ypz;
        float xd = xz * PerezFunction (domeTheta, gamma, cx) / xpz;
        float yd = yz * PerezFunction (domeTheta, gamma, cy) / ypz;

        // Convert to RGB
        YxyToRGB (Yd, xd, yd, &r, &g, &b);

        sgSetVec3 (rgb_dome[i][j], r, g, b);
      }
    }
  }

  //
  // Step 2 : Apply calculated colours to sky dome
  //

  // Apply colour to zenith vertex
  ssgVtxTable* vt = (ssgVtxTable*) top->getKid (0);
  float *c = vt->getColour (0);
  sgSetVec4 (c, rgb_zenith[0], rgb_zenith[1], rgb_zenith[2], 1.0);

  // Apply colours to bottom row of zenith
  for (j=0; j<SKYDOME_SLICES; j++) {
    c = vt->getColour (j+1);
    sgSetVec4 (c, rgb_dome[0][j][0], rgb_dome[0][j][1], rgb_dome[0][j][2], 1.0);
  }
  c = vt->getColour (SKYDOME_SLICES+1);
  sgSetVec4 (c, rgb_dome[0][0][0], rgb_dome[0][0][1], rgb_dome[0][0][2], 1.0);

  // Iterate over stacks
  for (i=0; i<SKYDOME_STACKS-1; i++) {
    ssgVtxTable *vt = (ssgVtxTable *) top->getKid (i+1);

    // Apply colours top vertices of this stack
    for (j=0; j<SKYDOME_SLICES; j++) {
      c = vt->getColour (j * 2);
      sgSetVec4 (c, rgb_dome[i][j][0], rgb_dome[i][j][1], rgb_dome[i][j][2], 1.0);
    }
    c = vt->getColour (SKYDOME_SLICES * 2);
    sgSetVec4 (c, rgb_dome[i][0][0], rgb_dome[i][0][1], rgb_dome[i][0][2], 1.0);

    // Apply colours to bottom vertices of this stack
    for (j=0; j<=SKYDOME_SLICES; j++) {
      c = vt->getColour ((j * 2) + 1);
      sgSetVec4 (c, rgb_dome[i+1][j][0], rgb_dome[i+1][j][1], rgb_dome[i+1][j][2], 1.0);
    }
    c = vt->getColour ((SKYDOME_SLICES * 2) + 1);
    sgSetVec4 (c, rgb_dome[i+1][0][0], rgb_dome[i+1][0][1], rgb_dome[i+1][0][2], 1.0);
  }

/*
  // TEST code to colour sky dome vertices
  sgVec4 test;
  sgSetVec4 (test, 1.0, 0.0, 0.0, 1.0);

  int ring = 0;
  ssgVtxTable *zenith = NULL;
  ssgVtxTable *topStack = NULL;
  ssgVtxTable *botStack = NULL;

  if (ring == 0) {
    zenith = vtZenith;
  }

  if (ring < SKYDOME_STACKS-1) {
    topStack = vtStack[ring];
  }

  if (ring > 0) {
    botStack = vtStack[ring-1];
  }

  for (int j=0; j<SKYDOME_SLICES; j++) {
    // Apply to zenith fan vertices
    if (vtZenith != NULL) {
      c = vtZenith->getColour (j+1);
      sgCopyVec4 (c, test);
    }

    // Apply to top vertices of stack
    if (topStack != NULL) {
      c = topStack->getColour (j * 2);
      sgCopyVec4 (c, test);

      if (j==0) {
        c = topStack->getColour (SKYDOME_SLICES * 2);
        sgCopyVec4 (c, test);
      }
    }
  }
*/
}


//
// Reposition the sky dome.  The only repositioning required is rotation of the
//   dome so that the 0th vertex of each ring is aligned with the solar azimuth
//
// Input parameters:
//  solAzimuth    Azimuth angle of the sun in radians W of S
//
void CSkyDomeImage::Reposition (double solAzimuth)
{
  // Make rotation matrix for sun azimuth angle
  sgVec3 axis;
  sgSetVec3 (axis, 0, 0, -1);
  sgMat4 SPIN;
  sgMakeRotMat4 (SPIN, RadToDeg (solAzimuth) - 90, axis);

  // Make top-level transform
  sgMat4 T;
  sgCopyMat4 (T, SPIN);
  top->setTransform (T);

/*
  *** Obsolete sky dome repositioning code to be deleted

  sgMat4 LON, LAT, SPIN;
  sgVec3 axis;

  // Rotate for eye position latitude and longitude
  sgSetVec3 (axis, 0, 0, 1);
  sgMakeRotMat4 (LON, lon, axis);

  sgSetVec3 (axis, 0, 1, 0);
  sgMakeRotMat4 (LAT, 90 - lat, axis);

  // Rotate for sun azimuth.  This ensures that the 0th vertex of each ring of the
  //   sky dome mesh is aligned with the current sun position, making it easier
  //   to apply complex sky dome shading algorithms in the future
  sgSetVec3 (axis, 0, 0, 1);
  sgMakeRotMat4 (SPIN, spin, axis);

  sgMat4 T;
  sgMakeIdentMat4 (T);
  sgPreMultMat4 (T, LON);
  sgPreMultMat4 (T, LAT);
  sgPreMultMat4 (T, SPIN);

  sgCoord skypos;
  sgSetCoord (&skypos, T);

  top->setTransform (&skypos);
*/
}


//
// Dump sky dome attributes to a file for debugging
//
// This is called from the parent CSkyManager class
//
void CSkyDomeImage::Print (FILE *f)
{

  // Display Perez parameter sets
  map<string,SPerezParameters*>::iterator i;
  for (i=perez.begin(); i!=perez.end(); i++) {
    SPerezParameters* p = i->second;

    fprintf (f, "Perez parameter set : %s\n", p->name);
    fprintf (f, "  Y (Luminance) curve parameters:\n");
    fprintf (f, "    A = %10.8f T + %10.8f\n", p->Y_curve.A[1], p->Y_curve.A[0]);
    fprintf (f, "    B = %10.8f T + %10.8f\n", p->Y_curve.B[1], p->Y_curve.B[0]);
    fprintf (f, "    C = %10.8f T + %10.8f\n", p->Y_curve.C[1], p->Y_curve.C[0]);
    fprintf (f, "    D = %10.8f T + %10.8f\n", p->Y_curve.D[1], p->Y_curve.D[0]);
    fprintf (f, "    E = %10.8f T + %10.8f\n", p->Y_curve.E[1], p->Y_curve.E[0]);

    fprintf (f, "  x (Chrominance) curve parameters:\n");
    fprintf (f, "    A = %10.8f T + %10.8f\n", p->x_curve.A[1], p->x_curve.A[0]);
    fprintf (f, "    B = %10.8f T + %10.8f\n", p->x_curve.B[1], p->x_curve.B[0]);
    fprintf (f, "    C = %10.8f T + %10.8f\n", p->x_curve.C[1], p->x_curve.C[0]);
    fprintf (f, "    D = %10.8f T + %10.8f\n", p->x_curve.D[1], p->x_curve.D[0]);
    fprintf (f, "    E = %10.8f T + %10.8f\n", p->x_curve.E[1], p->x_curve.E[0]);

    fprintf (f, "  y (Chrominance) curve parameters:\n");
    fprintf (f, "    A = %10.8f T + %10.8f\n", p->y_curve.A[1], p->y_curve.A[0]);
    fprintf (f, "    B = %10.8f T + %10.8f\n", p->y_curve.B[1], p->y_curve.B[0]);
    fprintf (f, "    C = %10.8f T + %10.8f\n", p->y_curve.C[1], p->y_curve.C[0]);
    fprintf (f, "    D = %10.8f T + %10.8f\n", p->y_curve.D[1], p->y_curve.D[0]);
    fprintf (f, "    E = %10.8f T + %10.8f\n", p->y_curve.E[1], p->y_curve.E[0]);

    fprintf (f, "  Y (Luminance) zenith parameters:\n");
    fprintf (f, "    A = %10.8f\n", p->Y_zenith.A);
    fprintf (f, "    B = %10.8f\n", p->Y_zenith.B);
    fprintf (f, "    C = %10.8f\n", p->Y_zenith.C);
    fprintf (f, "    D = %10.8f\n", p->Y_zenith.D);

    fprintf (f, "  x (Chrominance) zenith parameters:\n");
    fprintf (f, "      theta^3 (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->x_zenith.t2[3], p->x_zenith.t1[3], p->x_zenith.t0[3]);
    fprintf (f, "    + theta^2 (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->x_zenith.t2[2], p->x_zenith.t1[2], p->x_zenith.t0[2]);
    fprintf (f, "    + theta   (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->x_zenith.t2[1], p->x_zenith.t1[1], p->x_zenith.t0[1]);
    fprintf (f, "    +         (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->x_zenith.t2[0], p->x_zenith.t1[0], p->x_zenith.t0[0]);

    fprintf (f, "  y (Chrominance) zenith parameters:\n");
    fprintf (f, "      theta^3 (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->y_zenith.t2[3], p->y_zenith.t1[3], p->y_zenith.t0[3]);
    fprintf (f, "    + theta^2 (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->y_zenith.t2[2], p->y_zenith.t1[2], p->y_zenith.t0[2]);
    fprintf (f, "    + theta   (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->y_zenith.t2[1], p->y_zenith.t1[1], p->y_zenith.t0[1]);
    fprintf (f, "    +         (%10.8f T^2 + %10.8f T + %10.8f\n",
      p->y_zenith.t2[0], p->y_zenith.t1[0], p->y_zenith.t0[0]);

    fprintf (f, "\n");
  }
}

Documentation generated by