eu07/maszyna
1
0
Fork 0
mirror of https://github.com/MaSzyna-EU07/maszyna.git synced 2026-03-22 15:05:03 +01:00
Code Issues Packages Projects Releases Wiki Activity

Reorganize source files into logical subdirectories

Browse Source 
Co-authored-by: Hirek193 <23196899+Hirek193@users.noreply.github.com>
This commit is contained in:
copilot-swe-agent[bot]
2026-03-14 19:01:57 +00:00
parent f981f81d55
commit 0531086bb9
221 changed files with 131 additions and 108 deletions
Download Patch File Download Diff File Expand all files Collapse all files

315
environment/moon.cpp Normal file
View File

@@ -0,0 +1,315 @@
#include "stdafx.h"
#include "moon.h"
#include "Globals.h"
#include "mtable.h"
#include "utilities.h"
#include "simulationtime.h"
//////////////////////////////////////////////////////////////////////////////////////////
// cSun -- class responsible for dynamic calculation of position and intensity of the Sun,
cMoon::cMoon() {
setLocation( 19.00f, 52.00f ); // default location roughly in centre of Poland
m_observer.press = 1013.0; // surface pressure, millibars
m_observer.temp = 15.0; // ambient dry-bulb temperature, degrees C
}
cMoon::~cMoon() {
}
void
cMoon::init() {
m_observer.timezone = -1.0 * simulation::Time.zone_bias();
// NOTE: we're calculating phase just once, because it's unlikely simulation will last a few days,
// plus a sudden texture change would be pretty jarring
phase();
}
void
cMoon::update( bool const Includephase ) {
m_observer.temp = Global.AirTemperature;
move();
if( Includephase ) { phase(); }
glm::vec3 position( 0.f, 0.f, -1.f );
position = glm::rotateX( position, glm::radians( static_cast<float>( m_body.elevref ) ) );
position = glm::rotateY( position, glm::radians( static_cast<float>( -m_body.hrang ) ) );
m_position = glm::normalize( position );
}
void
cMoon::render() {
//m7t
}
glm::vec3
cMoon::getDirection() const {
return m_position;
}
float
cMoon::getAngle() const {
return (float)m_body.elevref;
}
float cMoon::getIntensity() {
irradiance();
// NOTE: we don't have irradiance model for the moon so we cheat here
// calculating intensity of the sun instead, and returning 15% of the value,
// which roughly matches how much sunlight is reflected by the moon
// We alter the intensity further based on current phase of the moon
auto const phasefactor = 1.0f - std::abs( m_phase - 29.53f * 0.5f ) / ( 29.53f * 0.5f );
return static_cast<float>( ( m_body.etr/ 1399.0 ) * phasefactor * 0.15 ); // arbitrary scaling factor taken from etrn value
}
void cMoon::setLocation( float const Longitude, float const Latitude ) {
// convert fraction from geographical base of 6o minutes
m_observer.longitude = (int)Longitude + (Longitude - (int)(Longitude)) * 100.0 / 60.0;
m_observer.latitude = (int)Latitude + (Latitude - (int)(Latitude)) * 100.0 / 60.0 ;
}
// sets current time, overriding one acquired from the system clock
void cMoon::setTime( int const Hour, int const Minute, int const Second ) {
m_observer.hour = clamp( Hour, -1, 23 );
m_observer.minute = clamp( Minute, -1, 59 );
m_observer.second = clamp( Second, -1, 59 );
}
void cMoon::setTemperature( float const Temperature ) {
m_observer.temp = Temperature;
}
void cMoon::setPressure( float const Pressure ) {
m_observer.press = Pressure;
}
void cMoon::move() {
static double radtodeg = 57.295779513; // converts from radians to degrees
static double degtorad = 0.0174532925; // converts from degrees to radians
SYSTEMTIME localtime = simulation::Time.data(); // time for the calculation
if( m_observer.hour >= 0 ) { localtime.wHour = m_observer.hour; }
if( m_observer.minute >= 0 ) { localtime.wMinute = m_observer.minute; }
if( m_observer.second >= 0 ) { localtime.wSecond = m_observer.second; }
double localut =
localtime.wHour
+ localtime.wMinute / 60.0 // too low resolution, noticeable skips
+ localtime.wSecond / 3600.0; // good enough in normal circumstances
/*
+ localtime.wMilliseconds / 3600000.0; // for really smooth movement
*/
double daynumber
= 367 * localtime.wYear
- 7 * ( localtime.wYear + ( localtime.wMonth + 9 ) / 12 ) / 4
+ 275 * localtime.wMonth / 9
+ localtime.wDay
- 730530
+ ( localut / 24.0 );
// Universal Coordinated (Greenwich standard) time
m_observer.utime = localut - m_observer.timezone;
// obliquity of the ecliptic
m_body.oblecl = clamp_circular( 23.4393 - 3.563e-7 * daynumber );
// moon parameters
double longascnode = clamp_circular( 125.1228 - 0.0529538083 * daynumber ); // N, degrees
double const inclination = 5.1454; // i, degrees
double const mndistance = 60.2666; // a, in earth radii
// argument of perigee
double const perigeearg = clamp_circular( 318.0634 + 0.1643573223 * daynumber ); // w, degrees
// mean anomaly
m_body.mnanom = clamp_circular( 115.3654 + 13.0649929509 * daynumber ); // M, degrees
// eccentricity
double const e = 0.054900;
// eccentric anomaly
double E0 = m_body.mnanom + radtodeg * e * std::sin( degtorad * m_body.mnanom ) * ( 1.0 + e * std::cos( degtorad * m_body.mnanom ) );
double E1 = E0 - ( E0 - radtodeg * e * std::sin( degtorad * E0 ) - m_body.mnanom ) / ( 1.0 - e * std::cos( degtorad * E0 ) );
while( std::abs( E0 - E1 ) > 0.005 ) { // arbitrary precision tolerance threshold
E0 = E1;
E1 = E0 - ( E0 - radtodeg * e * std::sin( degtorad * E0 ) - m_body.mnanom ) / ( 1.0 - e * std::cos( degtorad * E0 ) );
}
double const E = E1;
// lunar orbit plane rectangular coordinates
double const xv = mndistance * ( std::cos( degtorad * E ) - e );
double const yv = mndistance * std::sin( degtorad * E ) * std::sqrt( 1.0 - e*e );
// distance
m_body.distance = std::sqrt( xv*xv + yv*yv ); // r
// true anomaly
m_body.tranom = clamp_circular( radtodeg * std::atan2( yv, xv ) ); // v
// ecliptic rectangular coordinates
double const vpluswinrad = degtorad * ( m_body.tranom + perigeearg );
double const xeclip = m_body.distance * ( std::cos( degtorad * longascnode ) * std::cos( vpluswinrad ) - std::sin( degtorad * longascnode ) * std::sin( vpluswinrad ) * std::cos( degtorad * inclination ) );
double const yeclip = m_body.distance * ( std::sin( degtorad * longascnode ) * std::cos( vpluswinrad ) + std::cos( degtorad * longascnode ) * std::sin( vpluswinrad ) * std::cos( degtorad * inclination ) );
double const zeclip = m_body.distance * std::sin( vpluswinrad ) * std::sin( degtorad * inclination );
// ecliptic coordinates
double ecliplat = radtodeg * std::atan2( zeclip, std::sqrt( xeclip*xeclip + yeclip*yeclip ) );
m_body.eclong = clamp_circular( radtodeg * std::atan2( yeclip, xeclip ) );
// distance
m_body.distance = std::sqrt( xeclip*xeclip + yeclip*yeclip + zeclip*zeclip );
// perturbations
// NOTE: perturbation calculation can be safely disabled if we don't mind error of 1-2 degrees
// Sun's mean anomaly: Ms (already computed)
double const sunmnanom = clamp_circular( 356.0470 + 0.9856002585 * daynumber ); // M
// Sun's mean longitude: Ls (already computed)
double const sunphlong = clamp_circular( 282.9404 + 4.70935e-5 * daynumber );
double const sunmnlong = clamp_circular( sunphlong + sunmnanom ); // L = w + M
// Moon's mean anomaly: Mm (already computed)
// Moon's mean longitude: Lm = N + w + M (for the Moon)
m_body.mnlong = clamp_circular( longascnode + perigeearg + m_body.mnanom );
// Moon's mean elongation: D = Lm - Ls
double const mnelong = clamp_circular( m_body.mnlong - sunmnlong );
// Moon's argument of latitude: F = Lm - N
double const arglat = clamp_circular( m_body.mnlong - longascnode );
// longitude perturbations
double const pertevection = -1.274 * std::sin( degtorad * ( m_body.mnanom - 2.0 * mnelong ) ); // Evection
double const pertvariation = +0.658 * std::sin( degtorad * ( 2.0 * mnelong ) ); // Variation
double const pertyearlyeqt = -0.186 * std::sin( degtorad * sunmnanom ); // Yearly equation
// latitude perturbations
double const pertlat = -0.173 * std::sin( degtorad * ( arglat - 2.0 * mnelong ) );
m_body.eclong += pertevection + pertvariation + pertyearlyeqt;
ecliplat += pertlat;
// declination
m_body.declin = radtodeg * std::asin( std::sin (m_body.oblecl * degtorad) * std::sin(m_body.eclong * degtorad) );
// right ascension
double top = std::cos( degtorad * m_body.oblecl ) * std::sin( degtorad * m_body.eclong );
double bottom = std::cos( degtorad * m_body.eclong );
m_body.rascen = clamp_circular( radtodeg * std::atan2( top, bottom ) );
// Greenwich mean sidereal time
m_observer.gmst = 6.697375 + 0.0657098242 * daynumber + m_observer.utime;
m_observer.gmst -= 24.0 * (int)( m_observer.gmst / 24.0 );
if( m_observer.gmst < 0.0 ) m_observer.gmst += 24.0;
// local mean sidereal time
m_observer.lmst = m_observer.gmst * 15.0 + m_observer.longitude;
m_observer.lmst -= 360.0 * (int)( m_observer.lmst / 360.0 );
if( m_observer.lmst < 0.0 ) m_observer.lmst += 360.0;
// hour angle
m_body.hrang = m_observer.lmst - m_body.rascen;
if( m_body.hrang < -180.0 ) m_body.hrang += 360.0; // (force it between -180 and 180 degrees)
else if( m_body.hrang > 180.0 ) m_body.hrang -= 360.0;
double cz; // cosine of the solar zenith angle
double tdatcd = std::cos( degtorad * m_body.declin );
double tdatch = std::cos( degtorad * m_body.hrang );
double tdatcl = std::cos( degtorad * m_observer.latitude );
double tdatsd = std::sin( degtorad * m_body.declin );
double tdatsl = std::sin( degtorad * m_observer.latitude );
cz = tdatsd * tdatsl + tdatcd * tdatcl * tdatch;
// (watch out for the roundoff errors)
if( fabs (cz) > 1.0 ) { cz >= 0.0 ? cz = 1.0 : cz = -1.0; }
m_body.zenetr = std::acos( cz ) * radtodeg;
m_body.elevetr = 90.0 - m_body.zenetr;
refract();
}
void cMoon::refract() {
static double degtorad = 0.0174532925; // converts from degrees to radians
double prestemp; // temporary pressure/temperature correction
double refcor; // temporary refraction correction
double tanelev; // tangent of the solar elevation angle
// if the sun is near zenith, the algorithm bombs; refraction near 0.
if( m_body.elevetr > 85.0 )
refcor = 0.0;
else {
tanelev = tan( degtorad * m_body.elevetr );
if( m_body.elevetr >= 5.0 )
refcor = 58.1 / tanelev
- 0.07 / pow( tanelev, 3 )
+ 0.000086 / pow( tanelev, 5 );
else if( m_body.elevetr >= -0.575 )
refcor = 1735.0
+ m_body.elevetr * ( -518.2 + m_body.elevetr *
( 103.4 + m_body.elevetr * ( -12.79 + m_body.elevetr * 0.711 ) ) );
else
refcor = -20.774 / tanelev;
prestemp = ( m_observer.press * 283.0 ) / ( 1013.0 * ( 273.0 + m_observer.temp ) );
refcor *= prestemp / 3600.0;
}
// refracted solar elevation angle
m_body.elevref = m_body.elevetr + refcor;
// refracted solar zenith angle
m_body.zenref = 90.0 - m_body.elevref;
}
void cMoon::irradiance() {
static double radtodeg = 57.295779513; // converts from radians to degrees
static double degtorad = 0.0174532925; // converts from degrees to radians
m_body.dayang = ( simulation::Time.year_day() - 1 ) * 360.0 / 365.0;
double sd = sin( degtorad * m_body.dayang ); // sine of the day angle
double cd = cos( degtorad * m_body.dayang ); // cosine of the day angle or delination
m_body.erv = 1.000110 + 0.034221*cd + 0.001280*sd;
double d2 = 2.0 * m_body.dayang;
double c2 = cos( degtorad * d2 );
double s2 = sin( degtorad * d2 );
m_body.erv += 0.000719*c2 + 0.000077*s2;
double solcon = 1367.0; // Solar constant, 1367 W/sq m
m_body.coszen = cos( degtorad * m_body.zenref );
if( m_body.coszen > 0.0 ) {
m_body.etrn = solcon * m_body.erv;
m_body.etr = m_body.etrn * m_body.coszen;
}
else {
m_body.etrn = 0.0;
m_body.etr = 0.0;
}
}
void
cMoon::phase() {
SYSTEMTIME lt = simulation::Time.data();
if ((lt.wMonth==5)&&(lt.wDay==4)) //May the forth be with you!
m_phase = 50;
else {
// calculate moon's age in days from new moon
float ip = normalize( ( simulation::Time.julian_day() - 2451550.1f ) / 29.530588853f );
m_phase = ip * 29.53f;
}
}
// normalize values to range 0...1
float
cMoon::normalize( const float Value ) const {
float value = Value - floor( Value );
if( value < 0.f ) { ++value; }
return value;
}

99
environment/moon.h Normal file
View File

@@ -0,0 +1,99 @@
#pragma once
#include "winheaders.h"
// TODO: sun and moon share code as celestial bodies, we could make a base class out of it
class cMoon {
public:
// types:
// methods:
void init();
void update( bool const Includephase = false );
void render();
// returns vector pointing at the sun
glm::vec3 getDirection() const;
// returns current elevation above horizon
float getAngle() const;
// returns current intensity of the sun
float getIntensity();
// returns current phase of the moon
float getPhase() const { return m_phase; }
// sets current time, overriding one acquired from the system clock
void setTime( int const Hour, int const Minute, int const Second );
// sets current geographic location
void setLocation( float const Longitude, float const Latitude );
// sets ambient temperature in degrees C.
void setTemperature( float const Temperature );
// sets surface pressure in milibars
void setPressure( float const Pressure );
// constructors:
cMoon();
// deconstructor:
~cMoon();
// members:
protected:
// types:
// methods:
// calculates sun position on the sky given specified time and location
void move();
// calculates position adjustment due to refraction
void refract();
// calculates light intensity at current moment
void irradiance();
void phase();
// helper, normalize values to range 0...1
float normalize( const float Value ) const;
// members:
struct celestialbody { // main planet parameters
double dayang; // day angle (daynum*360/year-length) degrees
double phlong; // longitude of perihelion
double mnlong; // mean longitude, degrees
double mnanom; // mean anomaly, degrees
double tranom; // true anomaly, degrees
double eclong; // ecliptic longitude, degrees.
double oblecl; // obliquity of ecliptic.
double declin; // declination--zenith angle of solar noon at equator, degrees NORTH.
double rascen; // right ascension, degrees
double hrang; // hour angle--hour of sun from solar noon, degrees WEST
double zenetr; // solar zenith angle, no atmospheric correction (= ETR)
double zenref; // solar zenith angle, deg. from zenith, refracted
double coszen; // cosine of refraction corrected solar zenith angle
double elevetr; // solar elevation, no atmospheric correction (= ETR)
double elevref; // solar elevation angle, deg. from horizon, refracted.
double distance; // distance from earth in AUs
double erv; // earth radius vector (multiplied to solar constant)
double etr; // extraterrestrial (top-of-atmosphere) W/sq m global horizontal solar irradiance
double etrn; // extraterrestrial (top-of-atmosphere) W/sq m direct normal solar irradiance
};
struct observer { // weather, time and position data in observer's location
double latitude; // latitude, degrees north (south negative)
double longitude; // longitude, degrees east (west negative)
int hour{ -1 }; // current time, used for calculation of utime. if set to -1, time for
int minute{ -1 };// calculation will be obtained from the local clock
int second{ -1 };
double utime; // universal (Greenwich) standard time
double timezone; // time zone, east (west negative). USA: Mountain = -7, Central = -6, etc.
double gmst; // Greenwich mean sidereal time, hours
double lmst; // local mean sidereal time, degrees
double temp; // ambient dry-bulb temperature, degrees C, used for refraction correction
double press; // surface pressure, millibars, used for refraction correction and ampress
};
celestialbody m_body;
observer m_observer;
glm::vec3 m_position;
float m_phase;
};

27
environment/sky.cpp Normal file
View File

@@ -0,0 +1,27 @@
/*
This Source Code Form is subject to the
terms of the Mozilla Public License, v.
2.0. If a copy of the MPL was not
distributed with this file, You can
obtain one at
http://mozilla.org/MPL/2.0/.
*/
#include "stdafx.h"
#include "sky.h"
#include "Globals.h"
#include "MdlMngr.h"
//---------------------------------------------------------------------------
//GLfloat lightPos[4] = {0.0f, 0.0f, 0.0f, 1.0f};
void TSky::Init() {
if( ( Global.asSky != "1" )
&& ( Global.asSky != "0" ) ) {
mdCloud = TModelsManager::GetModel( Global.asSky );
}
};
//---------------------------------------------------------------------------

29
environment/sky.h Normal file
View File

@@ -0,0 +1,29 @@
/*
This Source Code Form is subject to the
terms of the Mozilla Public License, v.
2.0. If a copy of the MPL was not
distributed with this file, You can
obtain one at
http://mozilla.org/MPL/2.0/.
*/
//---------------------------------------------------------------------------
#pragma once
#include "Classes.h"
class TSky {
friend opengl_renderer;
friend opengl33_renderer;
public:
TSky() = default;
void Init();
private:
TModel3d *mdCloud { nullptr };
};
//---------------------------------------------------------------------------

318
environment/skydome.cpp Normal file
View File

@@ -0,0 +1,318 @@
#include "stdafx.h"
#include "skydome.h"
#include "color.h"
#include "utilities.h"
#include "simulationenvironment.h"
#include "Globals.h"
// sky gradient based on "A practical analytic model for daylight"
// by A. J. Preetham Peter Shirley Brian Smits (University of Utah)
float CSkyDome::m_distributionluminance[ 5 ][ 2 ] = { // Perez distributions
{ 0.17872f , -1.66303f }, // a = darkening or brightening of the horizon
{ -0.35540f , 0.42750f }, // b = luminance gradient near the horizon,
{ -0.02266f , 5.32505f }, // c = relative intensity of the circumsolar region
{ 0.12064f , -2.57705f }, // d = width of the circumsolar region
{ -0.06696f , 0.37027f } // e = relative backscattered light
};
float CSkyDome::m_distributionxcomp[ 5 ][ 2 ] = {
{ -0.01925f , -0.25922f },
{ -0.06651f , 0.00081f },
{ -0.00041f , 0.21247f },
{ -0.06409f , -0.89887f },
{ -0.00325f , 0.04517f }
};
float CSkyDome::m_distributionycomp[ 5 ][ 2 ] = {
{ -0.01669f , -0.26078f },
{ -0.09495f , 0.00921f },
{ -0.00792f , 0.21023f },
{ -0.04405f , -1.65369f },
{ -0.01092f , 0.05291f }
};
float CSkyDome::m_zenithxmatrix[ 3 ][ 4 ] = {
{ 0.00165f, -0.00375f, 0.00209f, 0.00000f },
{ -0.02903f, 0.06377f, -0.03202f, 0.00394f },
{ 0.11693f, -0.21196f, 0.06052f, 0.25886f }
};
float CSkyDome::m_zenithymatrix[ 3 ][ 4 ] = {
{ 0.00275f, -0.00610f, 0.00317f, 0.00000f },
{ -0.04214f, 0.08970f, -0.04153f, 0.00516f },
{ 0.15346f, -0.26756f, 0.06670f, 0.26688f }
};
//******************************************************************************//
CSkyDome::CSkyDome (int const Tesselation) :
m_tesselation( Tesselation ) {
// SetSunPosition( Math3D::vector3(75.0f, 0.0f, 0.0f) );
SetTurbidity( Global.fTurbidity );
SetExposure( true, 10.0f );
SetOvercastFactor( 0.05f );
Generate();
}
CSkyDome::~CSkyDome() {
}
//******************************************************************************//
void CSkyDome::Generate() {
// radius of dome
float const radius = 1.0f;
float const offset = 0.1f * radius; // horizontal offset, a cheap way to prevent a gap between ground and horizon
// create geometry chunk
int const latitudes = m_tesselation / 2 / 2; // half-sphere only
int const longitudes = m_tesselation;
std::uint16_t index = 0;
for( int i = 0; i <= latitudes; ++i ) {
float const latitude = M_PI * ( -0.5f + (float)( i ) / latitudes / 2 ); // half-sphere only
float const z = std::sin( latitude );
float const zr = std::cos( latitude );
for( int j = 0; j <= longitudes; ++j ) {
float const longitude = 2.0 * M_PI * (float)( j ) / longitudes;
float const x = std::cos( longitude );
float const y = std::sin( longitude );
/*
m_vertices.emplace_back( float3( x * zr, y * zr - offset, z ) * radius );
// we aren't using normals, but the code is left here in case it's ever needed
// m_normals.emplace_back( float3( x * zr, -y * zr, -z ) );
*/
// cartesian to opengl swap: -x, -z, -y
m_vertices.emplace_back( glm::vec3( -x * zr, -z - offset, -y * zr ) * radius );
m_colours.emplace_back( glm::vec3( 0.75f, 0.75f, 0.75f ) ); // placeholder
if( (i == 0) || (j == 0) ) {
// initial edge of the dome, don't start indices yet
++index;
}
else {
// indices for two triangles, formed between current and previous latitude
m_indices.emplace_back( index - 1 - (longitudes + 1) );
m_indices.emplace_back( index - 1 );
m_indices.emplace_back( index );
m_indices.emplace_back( index );
m_indices.emplace_back( index - ( longitudes + 1 ) );
m_indices.emplace_back( index - 1 - ( longitudes + 1 ) );
++index;
}
}
}
}
void CSkyDome::Update( glm::vec3 const &Sun ) {
SetTurbidity(Global.fTurbidity);
if( true == SetSunPosition( Sun ) ) {
// build colors if there's a change in sun position
RebuildColors();
}
}
// render skydome to screen
bool CSkyDome::SetSunPosition( glm::vec3 const &Direction ) {
if( Direction == m_sundirection ) {
return false;
}
m_sundirection = Direction;
m_thetasun = std::acos( m_sundirection.y );
m_phisun = std::atan2( m_sundirection.z, m_sundirection.x );
return true;
}
void CSkyDome::SetTurbidity( float const Turbidity ) {
m_turbidity = clamp( Turbidity, 1.f, 128.f );
}
void CSkyDome::SetExposure( bool const Linearexposure, float const Expfactor ) {
m_linearexpcontrol = Linearexposure;
m_expfactor = 1.0f / clamp( Expfactor, 1.0f, std::numeric_limits<float>::infinity() );
}
void CSkyDome::SetOvercastFactor( float const Overcast ) {
m_overcast = clamp( Overcast, 0.0f, 1.0f ) * 0.75f; // going above 0.65 makes the model go pretty bad, appearance-wise
}
void CSkyDome::GetPerez( float *Perez, float Distribution[ 5 ][ 2 ], const float Turbidity ) {
Perez[ 0 ] = Distribution[ 0 ][ 0 ] * Turbidity + Distribution[ 0 ][ 1 ];
Perez[ 1 ] = Distribution[ 1 ][ 0 ] * Turbidity + Distribution[ 1 ][ 1 ];
Perez[ 2 ] = Distribution[ 2 ][ 0 ] * Turbidity + Distribution[ 2 ][ 1 ];
Perez[ 3 ] = Distribution[ 3 ][ 0 ] * Turbidity + Distribution[ 3 ][ 1 ];
Perez[ 4 ] = Distribution[ 4 ][ 0 ] * Turbidity + Distribution[ 4 ][ 1 ];
}
float CSkyDome::GetZenith( float Zenithmatrix[ 3 ][ 4 ], const float Theta, const float Turbidity ) {
const float theta2 = Theta*Theta;
const float theta3 = Theta*theta2;
return ( Zenithmatrix[0][0] * theta3 + Zenithmatrix[0][1] * theta2 + Zenithmatrix[0][2] * Theta + Zenithmatrix[0][3]) * Turbidity * Turbidity +
( Zenithmatrix[1][0] * theta3 + Zenithmatrix[1][1] * theta2 + Zenithmatrix[1][2] * Theta + Zenithmatrix[1][3]) * Turbidity +
( Zenithmatrix[2][0] * theta3 + Zenithmatrix[2][1] * theta2 + Zenithmatrix[2][2] * Theta + Zenithmatrix[2][3]);
}
float CSkyDome::PerezFunctionO1( float Perezcoeffs[ 5 ], const float Thetasun, const float Zenithval ) {
const float val = ( 1.0f + Perezcoeffs[ 0 ] * std::exp( Perezcoeffs[ 1 ] ) ) *
( 1.0f + Perezcoeffs[ 2 ] * std::exp( Perezcoeffs[ 3 ] * Thetasun ) + Perezcoeffs[ 4 ] * std::pow( std::cos( Thetasun ), 2 ) );
return Zenithval / val;
}
float CSkyDome::PerezFunctionO2( float Perezcoeffs[ 5 ], const float Icostheta, const float Gamma, const float Cosgamma2, const float Zenithval ) {
// iCosTheta = 1.0f / cosf(theta)
// cosGamma2 = SQR( cosf( gamma ) )
return Zenithval * ( 1.0f + Perezcoeffs[ 0 ] * std::exp( Perezcoeffs[ 1 ] * Icostheta ) ) *
( 1.0f + Perezcoeffs[ 2 ] * std::exp( Perezcoeffs[ 3 ] * Gamma ) + Perezcoeffs[ 4 ] * Cosgamma2 );
}
void CSkyDome::RebuildColors() {
float twilightfactor = clamp( -simulation::Environment.sun().getAngle(), 0.0f, 18.0f ) / 18.0f;
auto gammacorrection = interpolate( glm::vec3( 0.45f ), glm::vec3( 1.0f ), twilightfactor );
// get zenith luminance
float const chi = ( (4.0f / 9.0f) - (m_turbidity / 120.0f) ) * ( M_PI - (2.0f * m_thetasun) );
float zenithluminance = ( (4.0453f * m_turbidity) - 4.9710f ) * std::tan( chi ) - (0.2155f * m_turbidity) + 2.4192f;
// get x / y zenith
float zenithx = GetZenith( m_zenithxmatrix, m_thetasun, m_turbidity );
float zenithy = GetZenith( m_zenithymatrix, m_thetasun, m_turbidity );
// get perez function parametrs
float perezluminance[5], perezx[5], perezy[5];
GetPerez( perezluminance, m_distributionluminance, m_turbidity );
GetPerez( perezx, m_distributionxcomp, m_turbidity );
GetPerez( perezy, m_distributionycomp, m_turbidity );
// make some precalculation
zenithx = PerezFunctionO1( perezx, m_thetasun, zenithx );
zenithy = PerezFunctionO1( perezy, m_thetasun, zenithy );
zenithluminance = PerezFunctionO1( perezluminance, m_thetasun, zenithluminance );
// start with fresh average for the new pass
glm::vec3 averagecolor, averagehorizoncolor;
// trough all vertices
glm::vec3 vertex;
glm::vec3 color, colorconverter, shiftedcolor;
for ( unsigned int i = 0; i < m_vertices.size(); ++i ) {
// grab it
vertex = glm::normalize( m_vertices[ i ] );
// angle between sun and vertex
const float gamma = std::acos( glm::dot( vertex, m_sundirection ) );
// warning : major hack!!! .. i had to do something with values under horizon
//vertex.y = Clamp<float>( vertex.y, 0.05f, 1.0f );
if ( vertex.y < 0.05f ) vertex.y = 0.05f;
// from paper:
// const float theta = arccos( vertex.y );
// const float iCosTheta = 1.0f / cosf( theta );
// optimized:
// iCosTheta =
// = 1.0f / cosf( arccos( vertex.y ) );
// = 1.0f / vertex.y;
float const icostheta = 1.0f / vertex.y;
float const cosgamma2 = std::pow( std::cos( gamma ), 2 );
// Compute x,y values
float const x = PerezFunctionO2( perezx, icostheta, gamma, cosgamma2, zenithx );
float const y = PerezFunctionO2( perezy, icostheta, gamma, cosgamma2, zenithy );
// luminance(Y) for clear & overcast sky
float const yclear = std::max( 0.01f, PerezFunctionO2( perezluminance, icostheta, gamma, cosgamma2, zenithluminance ) );
float const yover = std::max( 0.01f, zenithluminance * ( 1.0f + 2.0f * vertex.y ) / 3.0f );
float const Y = interpolate( yclear, yover, m_overcast );
float const X = (x / y) * Y;
float const Z = ((1.0f - x - y) / y) * Y;
colorconverter = glm::vec3( X, Y, Z );
color = colors::XYZtoRGB( colorconverter );
colorconverter = colors::RGBtoHSV(color);
if ( m_linearexpcontrol ) {
// linear scale
colorconverter.z *= m_expfactor;
} else {
// exp scale
colorconverter.z = 1.0f - std::exp( -m_expfactor * colorconverter.z );
}
if( colorconverter.z > 0.85f ) {
colorconverter.z = 0.85f + ( colorconverter.z - 0.85f ) * 0.35f;
}
colorconverter.y = clamp( colorconverter.y * Global.m_skysaturationcorrection, 0.0f, 1.0f );
// desaturate sky colour, based on overcast level
if( colorconverter.y > 0.0f ) {
colorconverter.y *= ( 1.0f - 0.5f * m_overcast );
}
// override the hue, based on sun height above the horizon. crude way to deal with model shortcomings
// correction begins when the sun is higher than 10 degrees above the horizon, and fully in effect at 10+15 degrees
float const degreesabovehorizon = 90.0f - m_thetasun * ( 180.0f / M_PI );
auto const sunbasedphase = clamp( (1.0f / 15.0f) * ( degreesabovehorizon - 10.0f ), 0.0f, 1.0f );
// correction is applied in linear manner from the bottom, becomes fully in effect for vertices with y = 0.50
auto const heightbasedphase = clamp( vertex.y * 2.0f, 0.0f, 1.0f );
// this height-based factor is reduced the farther the sun is up in the sky
float const shiftfactor = clamp( interpolate(heightbasedphase, sunbasedphase, sunbasedphase), 0.0f, 1.0f );
// h = 210 makes for 'typical' sky tone
shiftedcolor = glm::vec3( 210.0f, colorconverter.y, colorconverter.z );
shiftedcolor = colors::HSVtoRGB( shiftedcolor );
color = colors::HSVtoRGB(colorconverter);
color = interpolate( color, shiftedcolor, shiftfactor * Global.m_skyhuecorrection );
// crude correction for the times where the model breaks (late night)
// TODO: use proper night sky calculation for these times instead
if( ( color.x <= 0.05f )
&& ( color.y <= 0.05f ) ) {
// darken the sky as the sun goes deeper below the horizon
// 15:50:75 is picture-based night sky colour. it may not be accurate but looks 'right enough'
color.z = 0.75f * std::max( color.z + m_sundirection.y, 0.075f );
color.x = 0.20f * color.z;
color.y = 0.65f * color.z;
}
// simple gradient, darkening towards the top
color *= clamp( ( 1.0f - vertex.y ), 0.f, 1.f );
//color *= ( 0.25f - vertex.y );
// gamma correction
color = glm::pow( color, gammacorrection );
m_colours[ i ] = color;
averagecolor += color;
if( ( m_vertices.size() - i ) <= ( m_tesselation * 10 + 10 ) ) {
// calculate horizon colour from the bottom band of tris
averagehorizoncolor += color;
}
}
m_averagecolour = glm::max( glm::vec3(), averagecolor / static_cast<float>( m_vertices.size() ) );
m_averagehorizoncolour = glm::max( glm::vec3(), averagehorizoncolor / static_cast<float>( m_tesselation * 10 + 10 ) );
m_dirty = true;
}
//******************************************************************************//

67
environment/skydome.h Normal file
View File

@@ -0,0 +1,67 @@
#pragma once
// sky gradient based on "A practical analytic model for daylight"
// by A. J. Preetham Peter Shirley Brian Smits (University of Utah)
class CSkyDome {
public:
CSkyDome( int const Tesselation = 54 );
~CSkyDome();
void Generate();
void RebuildColors();
bool SetSunPosition( glm::vec3 const &Direction );
void SetTurbidity( const float Turbidity = 5.0f );
void SetExposure( const bool Linearexposure, const float Expfactor );
void SetOvercastFactor( const float Overcast = 0.0f );
// update skydome
void Update( glm::vec3 const &Sun );
// retrieves average colour of the sky dome
glm::vec3 GetAverageColor() { return m_averagecolour * 8.f / 6.f; }
glm::vec3 GetAverageHorizonColor() { return m_averagehorizoncolour; }
std::vector<glm::vec3> const & vertices() const {
return m_vertices; }
std::vector<glm::vec3> & colors() {
return m_colours; }
std::vector<std::uint16_t> const & indices() const {
return m_indices; }
auto const & is_dirty() const { return m_dirty; }
auto & is_dirty() { return m_dirty; }
private:
// shading parametrs
glm::vec3 m_sundirection;
float m_thetasun, m_phisun;
float m_turbidity;
bool m_linearexpcontrol;
float m_expfactor;
float m_overcast;
glm::vec3 m_averagecolour;
glm::vec3 m_averagehorizoncolour;
// data
int const m_tesselation;
std::vector<glm::vec3> m_vertices;
std::vector<std::uint16_t> m_indices;
// std::vector<float3> m_normals;
std::vector<glm::vec3> m_colours;
bool m_dirty { true }; // indicates sync state between simulation and gpu sides
static float m_distributionluminance[ 5 ][ 2 ];
static float m_distributionxcomp[ 5 ][ 2 ];
static float m_distributionycomp[ 5 ][ 2 ];
static float m_zenithxmatrix[ 3 ][ 4 ];
static float m_zenithymatrix[ 3 ][ 4 ];
// coloring
void GetPerez( float *Perez, float Distribution[ 5 ][ 2 ], const float Turbidity );
float GetZenith( float Zenithmatrix[ 3 ][ 4 ], const float Theta, const float Turbidity );
float PerezFunctionO1( float Perezcoeffs[ 5 ], const float Thetasun, const float Zenithval );
float PerezFunctionO2( float Perezcoeffs[ 5 ], const float Icostheta, const float Gamma, const float Cosgamma2, const float Zenithval );
};

13
environment/stars.cpp Normal file
View File

@@ -0,0 +1,13 @@
#include "stdafx.h"
#include "stars.h"
#include "Globals.h"
#include "MdlMngr.h"
//////////////////////////////////////////////////////////////////////////////////////////
// cStars -- simple starfield model, simulating appearance of starry sky
void
cStars::init() {
m_stars = TModelsManager::GetModel( "skydome_stars.t3d", false );
}

34
environment/stars.h Normal file
View File

@@ -0,0 +1,34 @@
#pragma once
#include "Classes.h"
//////////////////////////////////////////////////////////////////////////////////////////
// cStars -- simple starfield model, simulating appearance of starry sky
class cStars {
friend opengl_renderer;
friend opengl33_renderer;
public:
// types:
// methods:
void init();
// constructors:
cStars() = default;
// deconstructor:
// members:
private:
// types:
// methods:
// members:
float m_longitude{ 19.0f }; // geograpic coordinates hardcoded roughly to Poland location, for the time being
float m_latitude{ 52.0f };
TModel3d *m_stars { nullptr };
};

264
environment/sun.cpp Normal file
View File

@@ -0,0 +1,264 @@
#include "stdafx.h"
#include "sun.h"
#include "Globals.h"
#include "mtable.h"
#include "utilities.h"
#include "simulationtime.h"
//////////////////////////////////////////////////////////////////////////////////////////
// cSun -- class responsible for dynamic calculation of position and intensity of the Sun,
cSun::cSun() {
setLocation( 19.00f, 52.00f ); // default location roughly in centre of Poland
m_observer.press = 1013.0; // surface pressure, millibars
m_observer.temp = 15.0; // ambient dry-bulb temperature, degrees C
}
cSun::~cSun()
{
}
void
cSun::init() {
m_observer.timezone = -1.0 * simulation::Time.zone_bias();
}
void
cSun::update() {
m_observer.temp = Global.AirTemperature;
move();
glm::vec3 position( 0.f, 0.f, -1.f );
position = glm::rotateX( position, glm::radians( static_cast<float>( m_body.elevref ) ) );
position = glm::rotateY( position, glm::radians( static_cast<float>( -m_body.hrang ) ) );
m_position = glm::normalize( position );
}
void
cSun::render() {
//m7t
}
/*
glm::vec3
cSun::getPosition() {
return m_position * 1000.f * Global.fDistanceFactor;
}
*/
glm::vec3
cSun::getDirection() const {
return m_position;
}
float
cSun::getAngle() const {
return (float)m_body.elevref;
}
// return current hour angle
double
cSun::getHourAngle() const {
return m_body.hrang;
}
float cSun::getIntensity() {
irradiance();
return (float)( m_body.etr/ 1399.0 ); // arbitrary scaling factor taken from etrn value
}
void cSun::setLocation( float const Longitude, float const Latitude ) {
// convert fraction from geographical base of 6o minutes
m_observer.longitude = (int)Longitude + (Longitude - (int)(Longitude)) * 100.0 / 60.0;
m_observer.latitude = (int)Latitude + (Latitude - (int)(Latitude)) * 100.0 / 60.0 ;
}
// sets current time, overriding one acquired from the system clock
void cSun::setTime( int const Hour, int const Minute, int const Second ) {
m_observer.hour = clamp( Hour, -1, 23 );
m_observer.minute = clamp( Minute, -1, 59 );
m_observer.second = clamp( Second, -1, 59 );
}
void cSun::setTemperature( float const Temperature ) {
m_observer.temp = Temperature;
}
void cSun::setPressure( float const Pressure ) {
m_observer.press = Pressure;
}
void cSun::move() {
static double radtodeg = 57.295779513; // converts from radians to degrees
static double degtorad = 0.0174532925; // converts from degrees to radians
SYSTEMTIME localtime = simulation::Time.data(); // time for the calculation
if( m_observer.hour >= 0 ) { localtime.wHour = m_observer.hour; }
if( m_observer.minute >= 0 ) { localtime.wMinute = m_observer.minute; }
if( m_observer.second >= 0 ) { localtime.wSecond = m_observer.second; }
double localut =
localtime.wHour
+ localtime.wMinute / 60.0 // too low resolution, noticeable skips
+ localtime.wSecond / 3600.0; // good enough in normal circumstances
/*
+ localtime.wMilliseconds / 3600000.0; // for really smooth movement
*/
double daynumber =
367 * localtime.wYear
- 7 * ( localtime.wYear + ( localtime.wMonth + 9 ) / 12 ) / 4
+ 275 * localtime.wMonth / 9
+ localtime.wDay
- 730530
+ ( localut / 24.0 );
// Universal Coordinated (Greenwich standard) time
m_observer.utime = localut - m_observer.timezone;
// perihelion longitude
m_body.phlong = 282.9404 + 4.70935e-5 * daynumber; // w
// orbit eccentricity
double const e = 0.016709 - 1.151e-9 * daynumber;
// mean anomaly
m_body.mnanom = clamp_circular( 356.0470 + 0.9856002585 * daynumber ); // M
// obliquity of the ecliptic
m_body.oblecl = 23.4393 - 3.563e-7 * daynumber;
// mean longitude
m_body.mnlong = clamp_circular( m_body.phlong + m_body.mnanom ); // L = w + M
// eccentric anomaly
double const E = m_body.mnanom + radtodeg * e * std::sin( degtorad * m_body.mnanom ) * ( 1.0 + e * std::cos( degtorad * m_body.mnanom ) );
// ecliptic plane rectangular coordinates
double const xv = std::cos( degtorad * E ) - e;
double const yv = std::sin( degtorad * E ) * std::sqrt( 1.0 - e*e );
// distance
m_body.distance = std::sqrt( xv*xv + yv*yv ); // r
// true anomaly
m_body.tranom = radtodeg * std::atan2( yv, xv ); // v
// ecliptic longitude
m_body.eclong = clamp_circular( m_body.tranom + m_body.phlong ); // lon = v + w
/*
// ecliptic rectangular coordinates
double const x = m_body.distance * std::cos( degtorad * m_body.eclong );
double const y = m_body.distance * std::sin( degtorad * m_body.eclong );
double const z = 0.0;
// equatorial rectangular coordinates
double const xequat = x;
double const yequat = y * std::cos( degtorad * m_body.oblecl ) - 0.0 * std::sin( degtorad * m_body.oblecl );
double const zequat = y * std::sin( degtorad * m_body.oblecl ) + 0.0 * std::cos( degtorad * m_body.oblecl );
// declination
m_body.declin = radtodeg * std::atan2( zequat, std::sqrt( xequat*xequat + yequat*yequat ) );
// right ascension
m_body.rascen = radtodeg * std::atan2( yequat, xequat );
*/
// declination
m_body.declin = radtodeg * std::asin( std::sin( m_body.oblecl * degtorad ) * std::sin( m_body.eclong * degtorad ) );
// right ascension
double top = std::cos( degtorad * m_body.oblecl ) * std::sin( degtorad * m_body.eclong );
double bottom = std::cos( degtorad * m_body.eclong );
m_body.rascen = clamp_circular( radtodeg * std::atan2( top, bottom ) );
// Greenwich mean sidereal time (hours)
m_observer.gmst = clamp_circular( 6.697375 + 0.0657098242 * daynumber + m_observer.utime, 24.0 );
// local mean sidereal time (deg)
m_observer.lmst = clamp_circular( m_observer.gmst * 15.0 + m_observer.longitude );
// hour angle
m_body.hrang = m_observer.lmst - m_body.rascen;
if( m_body.hrang < -180.0 ) m_body.hrang += 360.0; // (force it between -180 and 180 degrees)
else if( m_body.hrang > 180.0 ) m_body.hrang -= 360.0;
double cz; // cosine of the solar zenith angle
double tdatcd = std::cos( degtorad * m_body.declin );
double tdatch = std::cos( degtorad * m_body.hrang );
double tdatcl = std::cos( degtorad * m_observer.latitude );
double tdatsd = std::sin( degtorad * m_body.declin );
double tdatsl = std::sin( degtorad * m_observer.latitude );
cz = tdatsd * tdatsl + tdatcd * tdatcl * tdatch;
// (watch out for the roundoff errors)
if( fabs( cz ) > 1.0 ) { cz >= 0.0 ? cz = 1.0 : cz = -1.0; }
m_body.zenetr = std::acos( cz ) * radtodeg;
m_body.elevetr = 90.0 - m_body.zenetr;
refract();
}
void cSun::refract() {
static double raddeg = 0.0174532925; // converts from degrees to radians
double prestemp; // temporary pressure/temperature correction
double refcor; // temporary refraction correction
double tanelev; // tangent of the solar elevation angle
// if the sun is near zenith, the algorithm bombs; refraction near 0.
if( m_body.elevetr > 85.0 )
refcor = 0.0;
else {
tanelev = tan( raddeg * m_body.elevetr );
if( m_body.elevetr >= 5.0 )
refcor = 58.1 / tanelev
- 0.07 / pow( tanelev, 3 )
+ 0.000086 / pow( tanelev, 5 );
else if( m_body.elevetr >= -0.575 )
refcor = 1735.0
+ m_body.elevetr * ( -518.2 + m_body.elevetr *
( 103.4 + m_body.elevetr * ( -12.79 + m_body.elevetr * 0.711 ) ) );
else
refcor = -20.774 / tanelev;
prestemp = ( m_observer.press * 283.0 ) / ( 1013.0 * ( 273.0 + m_observer.temp ) );
refcor *= prestemp / 3600.0;
}
// refracted solar elevation angle
m_body.elevref = m_body.elevetr + refcor;
// refracted solar zenith angle
m_body.zenref = 90.0 - m_body.elevref;
}
void cSun::irradiance() {
m_body.dayang = ( simulation::Time.year_day() - 1 ) * 360.0 / 365.0;
double sd = std::sin( glm::radians( m_body.dayang ) ); // sine of the day angle
double cd = std::cos( glm::radians( m_body.dayang ) ); // cosine of the day angle or delination
m_body.erv = 1.000110 + 0.034221*cd + 0.001280*sd;
double d2 = 2.0 * m_body.dayang;
double c2 = std::cos( glm::radians( d2 ) );
double s2 = std::sin( glm::radians( d2 ) );
m_body.erv += 0.000719*c2 + 0.000077*s2;
double solcon = 1367.0; // Solar constant, 1367 W/sq m
m_body.coszen = std::cos( glm::radians( m_body.zenref ) );
if( m_body.coszen > 0.0 ) {
m_body.etrn = solcon * m_body.erv;
m_body.etr = m_body.etrn * m_body.coszen;
}
else {
m_body.etrn = 0.0;
m_body.etr = 0.0;
}
}

101
environment/sun.h Normal file
View File

@@ -0,0 +1,101 @@
#pragma once
#include "winheaders.h"
//////////////////////////////////////////////////////////////////////////////////////////
// cSun -- class responsible for dynamic calculation of position and intensity of the Sun,
// given current weather, time and geographic location.
class cSun {
public:
// types:
// methods:
void init();
void update();
void render();
/*
// returns location of the sun in the 3d scene
glm::vec3 getPosition();
*/
// returns vector pointing at the sun
glm::vec3 getDirection() const;
// returns current elevation above horizon
float getAngle() const;
// return current hour angle
double getHourAngle() const;
// returns current intensity of the sun
float getIntensity();
// sets current time, overriding one acquired from the system clock
void setTime( int const Hour, int const Minute, int const Second );
// sets current geographic location
void setLocation( float const Longitude, float const Latitude );
// sets ambient temperature in degrees C.
void setTemperature( float const Temperature );
// sets surface pressure in milibars
void setPressure( float const Pressure );
// constructors:
cSun();
// deconstructor:
~cSun();
// members:
protected:
// types:
// methods:
// calculates sun position on the sky given specified time and location
void move();
// calculates position adjustment due to refraction
void refract();
// calculates light intensity at current moment
void irradiance();
// members:
struct celestialbody { // main planet parameters
double dayang; // day angle (daynum*360/year-length) degrees
double phlong; // longitude of perihelion
double mnlong; // mean longitude, degrees
double mnanom; // mean anomaly, degrees
double tranom; // true anomaly, degrees
double eclong; // ecliptic longitude, degrees.
double oblecl; // obliquity of ecliptic.
double declin; // declination--zenith angle of solar noon at equator, degrees NORTH.
double rascen; // right ascension, degrees
double hrang; // hour angle--hour of sun from solar noon, degrees WEST
double zenetr; // solar zenith angle, no atmospheric correction (= ETR)
double zenref; // solar zenith angle, deg. from zenith, refracted
double coszen; // cosine of refraction corrected solar zenith angle
double elevetr; // solar elevation, no atmospheric correction (= ETR)
double elevref; // solar elevation angle, deg. from horizon, refracted.
double distance; // distance from earth in AUs
double erv; // earth radius vector (multiplied to solar constant)
double etr; // extraterrestrial (top-of-atmosphere) W/sq m global horizontal solar irradiance
double etrn; // extraterrestrial (top-of-atmosphere) W/sq m direct normal solar irradiance
};
struct observer { // weather, time and position data in observer's location
double latitude; // latitude, degrees north (south negative)
double longitude; // longitude, degrees east (west negative)
int hour{ -1 }; // current time, used for calculation of utime. if set to -1, time for
int minute{ -1 };// calculation will be obtained from the local clock
int second{ -1 };
double utime; // universal (Greenwich) standard time
double timezone; // time zone, east (west negative). USA: Mountain = -7, Central = -6, etc.
double gmst; // Greenwich mean sidereal time, hours
double lmst; // local mean sidereal time, degrees
double temp; // ambient dry-bulb temperature, degrees C, used for refraction correction
double press; // surface pressure, millibars, used for refraction correction and ampress
};
celestialbody m_body;
observer m_observer;
glm::vec3 m_position;
};