simgear/scene/sky/oursun.cxx

   1 // oursun.hxx -- model earth's sun
   2 //
   3 // Written by Durk Talsma. Originally started October 1997, for distribution
   4 // with the FlightGear project. Version 2 was written in August and
   5 // September 1998. This code is based upon algorithms and data kindly
   6 // provided by Mr. Paul Schlyter. (pausch@saaf.se).
   7 //
   8 // Separated out rendering pieces and converted to ssg by Curt Olson,
   9 // March 2000
  10 // This library is free software; you can redistribute it and/or
  11 // modify it under the terms of the GNU Library General Public
  12 // License as published by the Free Software Foundation; either
  13 // version 2 of the License, or (at your option) any later version.
  14 //
  15 // This library is distributed in the hope that it will be useful,
  16 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  18 // Library General Public License for more details.
  19 //
  20 // You should have received a copy of the GNU General Public License
  21 // along with this program; if not, write to the Free Software
  22 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
  23 //
  24 // $Id$
  25
  26
  27 #ifdef HAVE_CONFIG_H
  28 #  include <simgear_config.h>
  29 #endif
  30
  31 #include <simgear/compiler.h>
  32
  33 #include <osg/AlphaFunc>
  34 #include <osg/BlendFunc>
  35 #include <osg/Fog>
  36 #include <osg/Geode>
  37 #include <osg/Geometry>
  38 #include <osg/Material>
  39 #include <osg/ShadeModel>
  40 #include <osg/TexEnv>
  41 #include <osg/Texture2D>
  42 #include <osgDB/ReadFile>
  43
  44 #include <simgear/screen/colors.hxx>
  45 #include <simgear/scene/model/model.hxx>
  46 #include "oursun.hxx"
  47
  48 // Constructor
  49 SGSun::SGSun( void ) {
  50     prev_sun_angle = -9999.0;
  51     visibility = -9999.0;
  52 }
  53
  54
  55 // Destructor
  56 SGSun::~SGSun( void ) {
  57 }
  58
  59
  60 // initialize the sun object and connect it into our scene graph root
  61 osg::Node*
  62 SGSun::build( SGPath path, double sun_size, SGPropertyNode *property_tree_Node ) {
  63
  64     env_node = property_tree_Node;
  65
  66     SGPath ihalopath = path, ohalopath = path;
  67
  68     // build the ssg scene graph sub tree for the sky and connected
  69     // into the provide scene graph branch
  70     sun_transform = new osg::MatrixTransform;
  71     osg::StateSet* stateSet = sun_transform->getOrCreateStateSet();
  72
  73     osg::TexEnv* texEnv = new osg::TexEnv;
  74     texEnv->setMode(osg::TexEnv::MODULATE);
  75     stateSet->setTextureAttribute(0, texEnv, osg::StateAttribute::ON);
  76
  77     osg::Material* material = new osg::Material;
  78     material->setColorMode(osg::Material::AMBIENT_AND_DIFFUSE);
  79     material->setEmission(osg::Material::FRONT_AND_BACK,
  80                           osg::Vec4(0, 0, 0, 1));
  81     material->setSpecular(osg::Material::FRONT_AND_BACK,
  82                           osg::Vec4(0, 0, 0, 1));
  83     stateSet->setAttribute(material);
  84
  85     osg::ShadeModel* shadeModel = new osg::ShadeModel;
  86     shadeModel->setMode(osg::ShadeModel::FLAT);
  87     stateSet->setAttributeAndModes(shadeModel);
  88
  89     osg::AlphaFunc* alphaFunc = new osg::AlphaFunc;
  90     alphaFunc->setFunction(osg::AlphaFunc::GREATER);
  91     alphaFunc->setReferenceValue(0.01);
  92     stateSet->setAttributeAndModes(alphaFunc);
  93
  94     osg::BlendFunc* blendFunc = new osg::BlendFunc;
  95     blendFunc->setSource(osg::BlendFunc::SRC_ALPHA);
  96     blendFunc->setDestination(osg::BlendFunc::ONE_MINUS_SRC_ALPHA);
  97     stateSet->setAttributeAndModes(blendFunc);
  98
  99     stateSet->setMode(GL_FOG, osg::StateAttribute::OFF);
 100     stateSet->setMode(GL_LIGHTING, osg::StateAttribute::OFF);
 101     stateSet->setMode(GL_CULL_FACE, osg::StateAttribute::OFF);
 102     stateSet->setMode(GL_DEPTH_TEST, osg::StateAttribute::OFF);
 103
 104
 105     osg::Geode* geode = new osg::Geode;
 106     stateSet = geode->getOrCreateStateSet();
 107
 108     stateSet->setRenderBinDetails(-8, "RenderBin");
 109
 110     // set up the sun-state
 111     path.append( "sun.rgba" );
 112     osg::Texture2D* texture = SGLoadTexture2D(path);
 113     stateSet->setTextureAttributeAndModes(0, texture);
 114
 115     // Build scenegraph
 116     sun_cl = new osg::Vec4Array;
 117     sun_cl->push_back(osg::Vec4(1, 1, 1, 1));
 118
 119     osg::Vec3Array* sun_vl = new osg::Vec3Array;
 120     sun_vl->push_back(osg::Vec3(-sun_size, 0, -sun_size));
 121     sun_vl->push_back(osg::Vec3(sun_size, 0, -sun_size));
 122     sun_vl->push_back(osg::Vec3(-sun_size, 0, sun_size));
 123     sun_vl->push_back(osg::Vec3(sun_size, 0, sun_size));
 124
 125     osg::Vec2Array* sun_tl = new osg::Vec2Array;
 126     sun_tl->push_back(osg::Vec2(0, 0));
 127     sun_tl->push_back(osg::Vec2(1, 0));
 128     sun_tl->push_back(osg::Vec2(0, 1));
 129     sun_tl->push_back(osg::Vec2(1, 1));
 130
 131     osg::Geometry* geometry = new osg::Geometry;
 132     geometry->setUseDisplayList(false);
 133     geometry->setVertexArray(sun_vl);
 134     geometry->setColorArray(sun_cl.get());
 135     geometry->setColorBinding(osg::Geometry::BIND_OVERALL);
 136     geometry->setNormalBinding(osg::Geometry::BIND_OFF);
 137     geometry->setTexCoordArray(0, sun_tl);
 138     geometry->addPrimitiveSet(new osg::DrawArrays(GL_TRIANGLE_STRIP, 0, 4));
 139     geode->addDrawable(geometry);
 140
 141     sun_transform->addChild( geode );
 142
 143
 144     // set up the inner-halo state
 145     geode = new osg::Geode;
 146     stateSet = geode->getOrCreateStateSet();
 147     stateSet->setRenderBinDetails(-7, "RenderBin");
 148
 149     ihalopath.append( "inner_halo.rgba" );
 150     texture = SGLoadTexture2D(path);
 151     stateSet->setTextureAttributeAndModes(0, texture);
 152
 153     // Build ssg structure
 154     ihalo_cl = new osg::Vec4Array;
 155     ihalo_cl->push_back(osg::Vec4(1, 1, 1, 1));
 156
 157     float ihalo_size = sun_size * 2.0;
 158     osg::Vec3Array* ihalo_vl = new osg::Vec3Array;
 159     ihalo_vl->push_back(osg::Vec3(-ihalo_size, 0, -ihalo_size));
 160     ihalo_vl->push_back(osg::Vec3(ihalo_size, 0, -ihalo_size));
 161     ihalo_vl->push_back(osg::Vec3(-ihalo_size, 0, ihalo_size));
 162     ihalo_vl->push_back(osg::Vec3(ihalo_size, 0, ihalo_size));
 163
 164     osg::Vec2Array* ihalo_tl = new osg::Vec2Array;
 165     ihalo_tl->push_back(osg::Vec2(0, 0));
 166     ihalo_tl->push_back(osg::Vec2(1, 0));
 167     ihalo_tl->push_back(osg::Vec2(0, 1));
 168     ihalo_tl->push_back(osg::Vec2(1, 1));
 169
 170     geometry = new osg::Geometry;
 171     geometry->setUseDisplayList(false);
 172     geometry->setVertexArray(ihalo_vl);
 173     geometry->setColorArray(ihalo_cl.get());
 174     geometry->setColorBinding(osg::Geometry::BIND_OVERALL);
 175     geometry->setNormalBinding(osg::Geometry::BIND_OFF);
 176     geometry->setTexCoordArray(0, ihalo_tl);
 177     geometry->addPrimitiveSet(new osg::DrawArrays(GL_TRIANGLE_STRIP, 0, 4));
 178     geode->addDrawable(geometry);
 179
 180     sun_transform->addChild( geode );
 181
 182
 183     // set up the outer halo state
 184
 185     geode = new osg::Geode;
 186     stateSet = geode->getOrCreateStateSet();
 187     stateSet->setRenderBinDetails(-6, "RenderBin");
 188
 189     ohalopath.append( "outer_halo.rgba" );
 190     texture = SGLoadTexture2D(path);
 191     stateSet->setTextureAttributeAndModes(0, texture);
 192
 193     // Build ssg structure
 194     ohalo_cl = new osg::Vec4Array;
 195     ohalo_cl->push_back(osg::Vec4(1, 1, 1, 1));
 196
 197     double ohalo_size = sun_size * 7.0;
 198     osg::Vec3Array* ohalo_vl = new osg::Vec3Array;
 199     ohalo_vl->push_back(osg::Vec3(-ohalo_size, 0, -ohalo_size));
 200     ohalo_vl->push_back(osg::Vec3(ohalo_size, 0, -ohalo_size));
 201     ohalo_vl->push_back(osg::Vec3(-ohalo_size, 0, ohalo_size));
 202     ohalo_vl->push_back(osg::Vec3(ohalo_size, 0, ohalo_size));
 203
 204     osg::Vec2Array* ohalo_tl = new osg::Vec2Array;
 205     ohalo_tl->push_back(osg::Vec2(0, 0));
 206     ohalo_tl->push_back(osg::Vec2(1, 0));
 207     ohalo_tl->push_back(osg::Vec2(0, 1));
 208     ohalo_tl->push_back(osg::Vec2(1, 1));
 209
 210     geometry = new osg::Geometry;
 211     geometry->setUseDisplayList(false);
 212     geometry->setVertexArray(ihalo_vl);
 213     geometry->setColorArray(ihalo_cl.get());
 214     geometry->setColorBinding(osg::Geometry::BIND_OVERALL);
 215     geometry->setNormalBinding(osg::Geometry::BIND_OFF);
 216     geometry->setTexCoordArray(0, ihalo_tl);
 217     geometry->addPrimitiveSet(new osg::DrawArrays(GL_TRIANGLE_STRIP, 0, 4));
 218     geode->addDrawable(geometry);
 219
 220     sun_transform->addChild( geode );
 221
 222
 223     // force a repaint of the sun colors with arbitrary defaults
 224     repaint( 0.0, 1.0 );
 225
 226     return sun_transform.get();
 227 }
 228
 229
 230 // repaint the sun colors based on current value of sun_angle in
 231 // degrees relative to verticle
 232 // 0 degrees = high noon
 233 // 90 degrees = sun rise/set
 234 // 180 degrees = darkest midnight
 235 bool SGSun::repaint( double sun_angle, double new_visibility ) {
 236
 237         if ( visibility != new_visibility ) {
 238                 visibility = new_visibility;
 239
 240                 static const double sqrt_m_log01 = sqrt( -log( 0.01 ) );
 241                 sun_exp2_punch_through = sqrt_m_log01 / ( visibility * 15 );
 242         }
 243
 244         if ( prev_sun_angle != sun_angle ) {
 245                 prev_sun_angle = sun_angle;
 246
 247                 // determine how much aerosols are in the air (rough guess)
 248                 double aerosol_factor;
 249                 if ( visibility < 100 ){
 250                         aerosol_factor = 8000;
 251                 }
 252                 else {
 253                         aerosol_factor = 80.5 / log( visibility / 100 );
 254                 }
 255
 256                 // get environmental data from property tree or use defaults
 257                 double rel_humidity, density_avg;
 258
 259                 if ( !env_node )
 260                 {
 261                         rel_humidity = 0.5;
 262                         density_avg = 0.7;
 263                 }
 264                 else
 265                 {
 266                         rel_humidity = env_node->getFloatValue( "relative-humidity" );
 267                         density_avg =  env_node->getFloatValue( "atmosphere/density-tropo-avg" );
 268                 }
 269
 270                 // ok, now let's go and generate the sun color
 271                 osg::Vec4 i_halo_color, o_halo_color, sun_color;
 272
 273                 // Some comments:
 274                 // When the sunangle changes, light has to travel a longer distance through the atmosphere.
 275                 // So it's scattered more due to raleigh scattering, which affects blue more than green light.
 276                 // Red is almost not scattered and effectively only get's touched when the sun is near the horizon.
 277                 // Visability also affects suncolor inasmuch as more particles are in the air that cause more scattering.
 278                 // We base our calculation on the halo's color, which is most scattered.
 279
 280                 // Red - is almost not scattered
 281                 // Lambda is 700 nm
 282
 283                 double red_scat_f = ( aerosol_factor * path_distance * density_avg ) / 5E+07;
 284                 sun_color[0] = 1 - red_scat_f;
 285                 i_halo_color[0] = 1 - ( 1.1 * red_scat_f );
 286                 o_halo_color[0] = 1 - ( 1.4 * red_scat_f );
 287
 288                 // Green - 546.1 nm
 289                 double green_scat_f = ( aerosol_factor * path_distance * density_avg ) / 8.8938E+06;
 290                 sun_color[1] = 1 - green_scat_f;
 291                 i_halo_color[1] = 1 - ( 1.1 * green_scat_f );
 292                 o_halo_color[1] = 1 - ( 1.4 * green_scat_f );
 293
 294                 // Blue - 435.8 nm
 295                 double blue_scat_f = ( aerosol_factor * path_distance * density_avg ) / 3.607E+06;
 296                 sun_color[2] = 1 - blue_scat_f;
 297                 i_halo_color[2] = 1 - ( 1.1 * blue_scat_f );
 298                 o_halo_color[2] = 1 - ( 1.4 * blue_scat_f );
 299
 300                 // Alpha
 301                 sun_color[3] = 1;
 302                 i_halo_color[3] = 1;
 303
 304                 o_halo_color[3] = blue_scat_f;
 305                 if ( ( new_visibility < 10000 ) &&  ( blue_scat_f > 1 )){
 306                         o_halo_color[3] = 2 - blue_scat_f;
 307                 }
 308
 309
 310                 // Now that we have the color calculated
 311                 // let's consider the saturation which is produced by mie scattering
 312                 double saturation = 1 - ( rel_humidity / 200 );
 313                 sun_color[1] += (( 1 - saturation ) * ( 1 - sun_color[1] ));
 314                 sun_color[2] += (( 1 - saturation ) * ( 1 - sun_color[2] ));
 315
 316                 i_halo_color[1] += (( 1 - saturation ) * ( 1 - i_halo_color[1] ));
 317                 i_halo_color[2] += (( 1 - saturation ) * ( 1 - i_halo_color[2] ));
 318
 319                 o_halo_color[1] += (( 1 - saturation ) * ( 1 - o_halo_color[1] ));
 320                 o_halo_color[2] += (( 1 - saturation ) * ( 1 - o_halo_color[2] ));
 321
 322                 // just to make sure we're in the limits
 323                 if ( sun_color[0] < 0 ) sun_color[0] = 0;
 324                 else if ( sun_color[0] > 1) sun_color[0] = 1;
 325                 if ( i_halo_color[0] < 0 ) i_halo_color[0] = 0;
 326                 else if ( i_halo_color[0] > 1) i_halo_color[0] = 1;
 327                 if ( o_halo_color[0] < 0 ) o_halo_color[0] = 0;
 328                 else if ( o_halo_color[0] > 1) o_halo_color[0] = 1;
 329
 330                 if ( sun_color[1] < 0 ) sun_color[1] = 0;
 331                 else if ( sun_color[1] > 1) sun_color[1] = 1;
 332                 if ( i_halo_color[1] < 0 ) i_halo_color[1] = 0;
 333                 else if ( i_halo_color[1] > 1) i_halo_color[1] = 1;
 334                 if ( o_halo_color[1] < 0 ) o_halo_color[1] = 0;
 335                 else if ( o_halo_color[1] > 1) o_halo_color[1] = 1;
 336
 337                 if ( sun_color[2] < 0 ) sun_color[2] = 0;
 338                 else if ( sun_color[2] > 1) sun_color[2] = 1;
 339                 if ( i_halo_color[2] < 0 ) i_halo_color[2] = 0;
 340                 else if ( i_halo_color[2] > 1) i_halo_color[2] = 1;
 341                 if ( o_halo_color[2] < 0 ) o_halo_color[2] = 0;
 342                 else if ( o_halo_color[2] > 1) o_halo_color[2] = 1;
 343                 if ( o_halo_color[3] < 0 ) o_halo_color[2] = 0;
 344                 else if ( o_halo_color[3] > 1) o_halo_color[3] = 1;
 345
 346
 347                 gamma_correct_rgb( i_halo_color._v );
 348                 gamma_correct_rgb( o_halo_color._v );
 349                 gamma_correct_rgb( sun_color._v );
 350
 351                 (*sun_cl)[0] = sun_color;
 352                 sun_cl->dirty();
 353                 (*ihalo_cl)[0] = i_halo_color;
 354                 ihalo_cl->dirty();
 355                 (*ohalo_cl)[0] = o_halo_color;
 356                 ohalo_cl->dirty();
 357     }
 358
 359     return true;
 360 }
 361
 362
 363 // reposition the sun at the specified right ascension and
 364 // declination, offset by our current position (p) so that it appears
 365 // fixed at a great distance from the viewer.  Also add in an optional
 366 // rotation (i.e. for the current time of day.)
 367 // Then calculate stuff needed for the sun-coloring
 368 bool SGSun::reposition( const SGVec3f& p, double angle,
 369                         double rightAscension, double declination,
 370                         double sun_dist, double lat, double alt_asl, double sun_angle)
 371 {
 372     // GST - GMT sidereal time
 373     osg::Matrix T1, T2, GST, RA, DEC;
 374
 375     T1.makeTranslate(p.osg());
 376     GST.makeRotate(SGD_DEGREES_TO_RADIANS*angle, osg::Vec3(0, 0, -1));
 377
 378     // xglRotatef( ((SGD_RADIANS_TO_DEGREES * rightAscension)- 90.0),
 379     //             0.0, 0.0, 1.0);
 380     RA.makeRotate(rightAscension - 90*SGD_DEGREES_TO_RADIANS, osg::Vec3(0, 0, 1));
 381
 382     // xglRotatef((SGD_RADIANS_TO_DEGREES * declination), 1.0, 0.0, 0.0);
 383     DEC.makeRotate(declination, osg::Vec3(1, 0, 0));
 384
 385     // xglTranslatef(0,sun_dist);
 386     T2.makeTranslate(osg::Vec3(0, sun_dist, 0));
 387
 388     sun_transform->setMatrix(T2*DEC*RA*GST*T1);
 389
 390     // Suncolor related things:
 391     if ( prev_sun_angle != sun_angle ) {
 392       if ( sun_angle == 0 ) sun_angle = 0.1;
 393          const double r_earth_pole = 6356752.314;
 394          const double r_tropo_pole = 6356752.314 + 8000;
 395          const double epsilon_earth2 = 6.694380066E-3;
 396          const double epsilon_tropo2 = 9.170014946E-3;
 397
 398          double r_tropo = r_tropo_pole / sqrt ( 1 - ( epsilon_tropo2 * pow ( cos( lat ), 2 )));
 399          double r_earth = r_earth_pole / sqrt ( 1 - ( epsilon_earth2 * pow ( cos( lat ), 2 )));
 400
 401          double position_radius = r_earth + alt_asl;
 402
 403          double gamma =  SG_PI - sun_angle;
 404          double sin_beta =  ( position_radius * sin ( gamma )  ) / r_tropo;
 405          double alpha =  SG_PI - gamma - asin( sin_beta );
 406
 407          // OK, now let's calculate the distance the light travels
 408          path_distance = sqrt( pow( position_radius, 2 ) + pow( r_tropo, 2 )
 409                         - ( 2 * position_radius * r_tropo * cos( alpha ) ));
 410
 411          double alt_half = sqrt( pow ( r_tropo, 2 ) + pow( path_distance / 2, 2 ) - r_tropo * path_distance * cos( asin( sin_beta )) ) - r_earth;
 412
 413          if ( alt_half < 0.0 ) alt_half = 0.0;
 414
 415          // Push the data to the property tree, so it can be used in the enviromental code
 416          if ( env_node ){
 417             env_node->setDoubleValue( "atmosphere/altitude-troposphere-top", r_tropo - r_earth );
 418             env_node->setDoubleValue( "atmosphere/altitude-half-to-sun", alt_half );
 419       }
 420     }
 421
 422     return true;
 423 }
 424
 425 SGVec4f
 426 SGSun::get_color()
 427 {
 428     return SGVec4f((*sun_cl)[0][0], (*sun_cl)[0][1], (*sun_cl)[0][2], (*sun_cl)[0][3]);
 429 }