Atmospheric Ray Tracing - Simulation Engine

- moved alot of the API to the cpp file

include/ITAPropagationPathSim/AtmosphericRayTracing/SimulationEngine.h

@@ -16,8 +16,8 @@
 *
 */
 
-#ifndef IW_ITA_PROPAGATIONPATHSIM_ART_ENGINE
-#define IW_ITA_PROPAGATIONPATHSIM_ART_ENGINE
+#ifndef IW_ITA_PROPAGATIONPATHSIM_ART_SIMULATIONENGINE
+#define IW_ITA_PROPAGATIONPATHSIM_ART_SIMULATIONENGINE
 
 #include <ITAPropagationPathSim/Definitions.h>
 
@@ -29,70 +29,67 @@
 #include <ITAPropagationPathSim/AtmosphericRayTracing/Rays.h>
 #include <ITAGeo/Atmosphere/StratifiedAtmosphere.h>
 
-
 namespace ITAPropagationPathSim
 {
 	namespace AtmosphericRayTracing
 	{
-		struct SimulationSettings {
-			ODESolver::SolverMethod solverMethod = ODESolver::EULER;
-			double dTMax = 30;
-			double dIntegrationTimeStep = 0.1;
-			IAbortCriterion* pAbortCriterion = nullptr; //!< Pointer to external abort criterion. If this is a NULL-pointer, it is ignored.
-		};
-
-		class ITA_PROPAGATION_PATH_SIM_API IAbortCriterion
-		{
-			virtual bool Abort(const CRay& ray) const = 0;
-		};
-
-		class ITA_PROPAGATION_PATH_SIM_API CSimulationEngine
-		{
-		private:
-			class CIntegrationStepAdaptor
+		namespace Simulation {
+			class ITA_PROPAGATION_PATH_SIM_API IAbortCriterion
+			{
+			public:
+				//! Returns true if the abort criterion for tracing the given ray is reached
+				virtual bool Check(const CRay& ray) const = 0; //TODO: Rename to Reached() ?
+			};
+			class ITA_PROPAGATION_PATH_SIM_API AbortAtMaxTime : public IAbortCriterion
 			{
+			private:
+				double dTMax;
 			public:
-				bool bActive = true;	//!< If this is set to false, this class is bypassed and the integration step size is therefore not adapted
+				AbortAtMaxTime(double tMax = 30) : dTMax(tMax) {}
+				bool Check(const CRay& ray) const { return ray.GetMaxTime() >= dTMax; }
+			};
 
+			struct AdaptiveIntegrationSettings {
+				bool bActive = true;	//!< If this is set to false, the adaptation is bypassed and the integration step size is therefore constant
 				double dMaxError = 0.015;
 				double dUncriticalError = 0.005;
+			};
 
-				CIntegrationStepAdaptor() {}
-				//CIntegrationStepAdaptor(double maxError, double uncriticalError) : dMaxError(maxError), dUncriticalError(uncriticalError) {}
+			struct SimulationSettings {
+				ODESolver::SolverMethod solverMethod = ODESolver::EULER;
+				double dIntegrationTimeStep = 0.1;
+				AdaptiveIntegrationSettings AdaptiveIntegration;
 			};
-			CIntegrationStepAdaptor integrationStepAdaptor;
 
-			class CWorker
+
+
+			class ITA_PROPAGATION_PATH_SIM_API CEngine
 			{
 			private:
-				CRay& rRay;
-				const SimulationSettings& rSettings;
-				const CIntegrationStepAdaptor& rIntegrationStepAdaptor;
-			public:
-				CWorker(CRay& ray, const SimulationSettings& settings, const CIntegrationStepAdaptor& adaptor) : rRay(ray), rSettings(settings), rIntegrationStepAdaptor(adaptor){}
+				class CWorker;
+				class CIntegrationStepAdaptor;
 
-				void TraceRay();
-			};
+			public:
+				SimulationSettings settings;
 
-		public:
-			SimulationSettings settings;
+				VistaVector3D v3SourcePosition;
+				std::vector<VistaVector3D> v3RayDirections;
+			private:
+				std::vector<CRay> rays;
+				const IAbortCriterion& iAbortCriterion; //!< Reference to externally defined abort criterion.
 
-			VistaVector3D v3SourcePosition;
-			std::vector<VistaVector3D> v3RayDirections;
-		private:
-			std::vector<CRay> rays;
-			
-		public:
-			CSimulationEngine() {}
-			virtual ~CSimulationEngine() {}
+			public:
+				CEngine(const IAbortCriterion& abortCriterion = AbortAtMaxTime()): iAbortCriterion(abortCriterion) {}
+				virtual ~CEngine() {}
 
-		public:
-			void Run();
-		private:
-			void InitRays();
-			CRay TraceRay(const VistaVector3D& v3Direction);
-		};
+			public:
+				void Run(const ITAGeo::CStratifiedAtmosphere& atmosphere);
+			private:
+				void InitRays();
+				CRay TraceRay(const VistaVector3D& v3Direction);
+			};
+		}
 	}
 }
 
-#endif //IW_ITA_PROPAGATIONPATHSIM_ART_ENGINE
\ No newline at end of file
+#endif //IW_ITA_PROPAGATIONPATHSIM_ART_SIMULATIONENGINE
\ No newline at end of file

src/ITAPropagationPathSim/AtmosphericRayTracing/SimulationEngine.cpp

 #include <ITAPropagationPathSim/AtmosphericRayTracing/SimulationEngine.h>
 
 // ITA includes
+#include <ITAException.h>
 #include <ITAPropagationPathSim/AtmosphericRayTracing/ODESolver.h>
 
 // STD
-//#include <cmath>
+#include <cmath>
 
-using namespace ITAPropagationPathSim::AtmosphericRayTracing;
+using namespace ITAPropagationPathSim::AtmosphericRayTracing::Simulation;
 
+class CEngine::CIntegrationStepAdaptor
+{
+private:
+	const AdaptiveIntegrationSettings settings;
+public:
+	CIntegrationStepAdaptor(const AdaptiveIntegrationSettings& settings) : settings(settings) {}
+public:
+	//TODO: Functionality of adapting dt according to settings, current ray values (and atmosphere?)
+};
+
+
+class CEngine::CWorker
+ {
+ private:
+	 CRay& rRay;
+	 const SimulationSettings& rSettings;
+	 const IAbortCriterion& iAbortCriterion; //!< Reference to externally defined abort criterion.
+	 CIntegrationStepAdaptor integrationStepAdaptor;
+ public:
+	 CWorker(CRay& ray, const SimulationSettings& simSettings, const IAbortCriterion& abortCriterion) :
+		 rRay(ray), rSettings(simSettings), iAbortCriterion(abortCriterion), integrationStepAdaptor(CIntegrationStepAdaptor(simSettings.AdaptiveIntegration)) {}
+ private:
+ public:
+	 void TraceRay(const ITAGeo::CStratifiedAtmosphere& atmosphere)
+	 {
+		 VistaVector3D r = rRay.GetLastPoint();
+		 double rz = r[Vista::Z];
+		 VistaVector3D n = rRay.GetLastWavefrontNormal();
+		 double time = rRay.GetMaxTime();
+
+		 VistaVector3D s = ODESolver::NormalToSlowness(n, rz, atmosphere);
+
+		 double dt = rSettings.dIntegrationTimeStep;
+
+		 while ( !iAbortCriterion.Check(rRay) )
+		 {
+			 switch (rSettings.solverMethod)
+			 {
+			 case ODESolver::EULER:
+				 ODESolver::Euler(r, s, dt, atmosphere);
+				 break;
+			 case ODESolver::RUNGE_KUTTA:
+				 ODESolver::RungeKutta(r, s, dt, atmosphere);
+				 break;
+			 default:
+				 ITA_EXCEPT_INVALID_PARAMETER("Unknown solver method.");
+			 }
+			 time += dt;
+
+			 rz = abs(r[Vista::Z]);
+			 n = ODESolver::SlownessToNormal(s, rz, atmosphere);
+			 rRay.Append(r, n, time);
+		 }
+	 }
+ };