Commit 4a3a09d7 authored by Philipp Schäfer's avatar Philipp Schäfer
Browse files

Atmospheric Ray Tracing

-----------------------
Rays
- removed Get... from get-functions
- added reflection order

Simulation/Engine
- Added function to add reflection and to copy period rays
- CWorker is now fully hidden in cpp file
parent baf383e8
......@@ -41,6 +41,7 @@ namespace ITAPropagationPathSim
std::vector<VistaVector3D> v3SamplingPoints;
std::vector<VistaVector3D> v3WavefrontNormals;
std::vector<double> dTimeStamps;
std::vector<unsigned int> iReflectionIndices;
public:
CRay(VistaVector3D v3SourcePos, VistaVector3D v3Direction)
......@@ -52,19 +53,22 @@ namespace ITAPropagationPathSim
virtual ~CRay() {}
public:
std::vector<VistaVector3D> GetSamplingPoints() const { return v3SamplingPoints; }
std::vector<VistaVector3D> GetWavefrontNormals() const { return v3WavefrontNormals; }
std::vector<double> GetTimeStamps() const { return dTimeStamps; }
std::vector<VistaVector3D> SamplingPoints() const { return v3SamplingPoints; }
std::vector<VistaVector3D> WavefrontNormals() const { return v3WavefrontNormals; }
std::vector<double> TimeStamps() const { return dTimeStamps; }
std::vector<unsigned int> RelectionIndices() const { return iReflectionIndices; }
unsigned int NumPoints() const { return v3SamplingPoints.size(); }
VistaVector3D GetSourcePoint() const { return v3SamplingPoints[0]; }
VistaVector3D GetInitialDirection() const { return v3WavefrontNormals[0]; }
VistaVector3D SourcePoint() const { return v3SamplingPoints[0]; }
VistaVector3D InitialDirection() const { return v3WavefrontNormals[0]; }
VistaVector3D GetLastPoint() const { return v3SamplingPoints[v3SamplingPoints.size()]; }
VistaVector3D GetLastWavefrontNormal() const { return v3WavefrontNormals[v3WavefrontNormals.size()]; }
double GetMaxTime() const { return dTimeStamps[dTimeStamps.size()]; }
VistaVector3D LastPoint() const { return v3SamplingPoints[v3SamplingPoints.size()]; }
VistaVector3D LastWavefrontNormal() const { return v3WavefrontNormals[v3WavefrontNormals.size()]; }
double MaxTime() const { return dTimeStamps[dTimeStamps.size()]; }
unsigned int ReflectionOrder() const { return iReflectionIndices.size(); }
bool SameDirection(const CRay& other) { return GetInitialDirection() == other.GetInitialDirection(); }
bool SameDirection(const CRay& other) { return InitialDirection() == other.InitialDirection(); }
void Append(const VistaVector3D& point, const VistaVector3D& wavefrontNormal, const double& timeStamp)
{
......@@ -72,6 +76,11 @@ namespace ITAPropagationPathSim
v3WavefrontNormals.push_back(wavefrontNormal);
dTimeStamps.push_back(timeStamp);
}
void AppendReflection(const VistaVector3D& point, const VistaVector3D& wavefrontNormal, const double& timeStamp)
{
Append(point, wavefrontNormal, timeStamp);
iReflectionIndices.push_back(NumPoints()-1);
}
};
}
}
......
......@@ -37,9 +37,6 @@ namespace ITAPropagationPathSim
class ITA_PROPAGATION_PATH_SIM_API CEngine
{
private:
class CWorker;
public:
SimulationSettings settings;
......@@ -57,7 +54,6 @@ namespace ITAPropagationPathSim
void Run(const ITAGeo::CStratifiedAtmosphere& atmosphere);
private:
void InitRays();
CRay TraceRay(const VistaVector3D& v3Direction);
};
}
}
......
......@@ -41,7 +41,7 @@ namespace ITAPropagationPathSim
double dTMax;
public:
AbortAtMaxTime(double tMax = 30) : dTMax(tMax) {}
bool AbortRequested(const CRay& ray) const { return ray.GetMaxTime() >= dTMax; }
bool AbortRequested(const CRay& ray) const { return ray.MaxTime() >= dTMax; }
};
enum ITA_PROPAGATION_PATH_SIM_API SolverMethod
......
......@@ -7,9 +7,12 @@
// STD
#include <cmath>
using namespace ITAPropagationPathSim::AtmosphericRayTracing;
using namespace ITAPropagationPathSim::AtmosphericRayTracing::Simulation;
class CEngine::CWorker
namespace ITAPropagationPathSim { namespace AtmosphericRayTracing { namespace Simulation {
class CWorker
{
private:
CRay& rRay;
......@@ -33,13 +36,47 @@ class CEngine::CWorker
rGround[Vista::Z] = 0;
dtGround = dt * timePortion;
}
void ExtendRayByOnePeriod()
{
std::vector<unsigned int> iReflectionIndices = rRay.RelectionIndices();
int reflectionOrder = rRay.ReflectionOrder();
if (reflectionOrder < 2)
ITA_EXCEPT_INVALID_PARAMETER("Cannot extend a ray with a reflection order below 2.");
int idxStartReflection = iReflectionIndices[reflectionOrder-2];
int idxEndReflection = iReflectionIndices[reflectionOrder - 1];
double tOffset = rRay.TimeStamps()[idxEndReflection] - rRay.TimeStamps()[idxStartReflection];
for (int idx = idxStartReflection+1; idx < idxEndReflection; idx++)
{
const double t = rRay.TimeStamps()[idx] + tOffset;
const VistaVector3D r = rRay.SamplingPoints()[idx];
const VistaVector3D n = rRay.WavefrontNormals()[idx];
rRay.Append(r, n, t);
if (rAbortCriterion.AbortRequested(rRay))
return;
}
const double t = rRay.TimeStamps()[idxEndReflection] + tOffset;
const VistaVector3D r = rRay.SamplingPoints()[idxEndReflection];
const VistaVector3D n = rRay.WavefrontNormals()[idxEndReflection];
rRay.AppendReflection(r, n, t);
}
void ExtendRayPeriodically()
{
while (!rAbortCriterion.AbortRequested(rRay))
ExtendRayByOnePeriod();
}
public:
void TraceRay(const ITAGeo::CStratifiedAtmosphere& atmosphere)
{
VistaVector3D r = rRay.GetLastPoint();
VistaVector3D r = rRay.LastPoint();
double rz = r[Vista::Z];
VistaVector3D n = rRay.GetLastWavefrontNormal();
double time = rRay.GetMaxTime();
VistaVector3D n = rRay.LastWavefrontNormal();
double time = rRay.MaxTime();
VistaVector3D s = ODESolver::NormalToSlowness(n, rz, atmosphere);
......@@ -50,14 +87,19 @@ class CEngine::CWorker
n = ODESolver::SlownessToNormal(s);
const double dt = mSolver.GetCurrentStepSize();
//if (ReflectionOccured(r[Vista::Z], rNew[Vista::Z]))
//{
// VistaVector3D rGround;
// double dtGround;
// GetReflectionPoint(r, r, dt, rGround, dtGround);
// rRay.AddReflection(rGround, n, time + dtGround);
//}
if (ReflectionOccured(r[Vista::Z], rNew[Vista::Z]))
{
VistaVector3D rGround;
double dtGround;
GetReflectionPoint(r, r, dt, rGround, dtGround);
rRay.AppendReflection(rGround, n, time + dtGround);
if (rRay.ReflectionOrder() >= 2) //Ray is periodic
{
ExtendRayPeriodically();
return;
}
}
r = rNew;
s = sNew;
......@@ -66,3 +108,21 @@ class CEngine::CWorker
}
}
};
}}}
void CEngine::Run(const ITAGeo::CStratifiedAtmosphere& atmosphere)
{
InitRays();
for (int idx = 0; idx < rays.size(); idx++)
{
CRay& ray = rays[idx];
CWorker(ray, settings, iAbortCriterion);
}
}
void CEngine::InitRays()
{
rays.clear();
for each (const VistaVector3D & v3Direction in v3RayDirections)
rays.push_back(CRay(v3SourcePosition, v3Direction));
}
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment