...
 
Commits (3)
  • Philipp Schäfer's avatar
    ART · c0cd7ae7
    Philipp Schäfer authored
    ++ documentation
    c0cd7ae7
  • Philipp Schäfer's avatar
    ART - AdaptiveRayGrid · 3a55471d
    Philipp Schäfer authored
    - advanced ray zooming now gets idxMin from ray directly
    - changed some variable names
    3a55471d
  • Philipp Schäfer's avatar
    ART - EigenraySearch · be667c2b
    Philipp Schäfer authored
    - added setting to (de)activate advanced
    - Engine now actually uses advanced ray zooming
    - spreading loss is now calculated for sound pressure and not intensity
    be667c2b
...@@ -58,7 +58,7 @@ namespace ITAPropagationPathSim ...@@ -58,7 +58,7 @@ namespace ITAPropagationPathSim
//! Simple adaptation method: Setting the limits to neighboring rays and shoots additional rays to double the angular resolution of the grid. //! Simple adaptation method: Setting the limits to neighboring rays and shoots additional rays to double the angular resolution of the grid.
void ZoomIntoRay(const std::shared_ptr<CRay> pRay); void ZoomIntoRay(const std::shared_ptr<CRay> pRay);
//! Advanced adaptation method: Further reduces the angular limits using additional information compared to simple method before doubling the angular resolution. //! Advanced adaptation method: Further reduces the angular limits using additional information compared to simple method before doubling the angular resolution.
void ZoomIntoRay(const std::shared_ptr<CRay> pRay, const int idxMinDist, const VistaVector3D& receiverPosition, const double& threshold); void ZoomIntoRay(const std::shared_ptr<CRay> pRay, const VistaVector3D& v3ReceiverPosition, const double& dThreshold);
//! Zooms into the given ray using a specific angular resolution by creating up to 8 new rays surrounding the given one. //! Zooms into the given ray using a specific angular resolution by creating up to 8 new rays surrounding the given one.
void ZoomIntoRay(const std::shared_ptr<CRay> pRay, const double& deltaTheta, const double& deltaPhi); void ZoomIntoRay(const std::shared_ptr<CRay> pRay, const double& deltaTheta, const double& deltaPhi);
......
...@@ -29,23 +29,24 @@ namespace ITAPropagationPathSim ...@@ -29,23 +29,24 @@ namespace ITAPropagationPathSim
namespace EigenraySearch { namespace EigenraySearch {
struct ITA_PROPAGATION_PATH_SIM_API RayAdaptationSettings { struct ITA_PROPAGATION_PATH_SIM_API RayAdaptationSettings {
struct { struct {
int maxNAdaptations = 15; //< Abort after X adaptations of the ray resolution int maxNAdaptations = 15; //< Abort after N adaptations of the ray resolution
double minAngleResolutionDeg = 0.001; //< [] double minAngleResolutionDeg = 0.001; //< Abort if angle between neighboring rays is below this value []
} abort; } abort;
struct { struct {
double maxReceiverRadius = 1; //< [m] double maxReceiverRadius = 1; //< Maximum radius of receiver sphere [m]
double maxSourceReceiverAngle = 1; //< [] double maxSourceReceiverAngle = 1; //< Maximum angle between vectors from source to receiver point and receiver sphere []
double maxAngleForGeomSpreading = 0.01; //< [] double maxAngleForGeomSpreading = 0.01; //< Angular resolution of rays which is required to calculate spreading loss []
} accuracy; } accuracy;
struct { struct {
double threshold = 0.1; //< [0 1] bool bActive = true; //< Switch to enable/disable advanced ray zooming
double threshold = 0.1; //< Threshold between 0 and 2 above which advanced ray zooming is performed (0 = always, 2 = never)
} advancedRayZooming; } advancedRayZooming;
}; };
struct ITA_PROPAGATION_PATH_SIM_API RayTracingAbortSettings { struct ITA_PROPAGATION_PATH_SIM_API RayTracingAbortSettings {
int maxReflectionOrder = 3; //< Maximum considered order of reflections int maxReflectionOrder = 1; //< Maximum considered order of reflections
double maxTime = 30; //< [s] double maxTime = 30; //< Maximum propagation time of rays [s]
}; };
struct ITA_PROPAGATION_PATH_SIM_API Settings { struct ITA_PROPAGATION_PATH_SIM_API Settings {
RayTracingAbortSettings rayTracing; RayTracingAbortSettings rayTracing;
......
...@@ -61,10 +61,10 @@ namespace ITAPropagationPathSim ...@@ -61,10 +61,10 @@ namespace ITAPropagationPathSim
struct ITA_PROPAGATION_PATH_SIM_API AdaptiveIntegrationSettings { struct ITA_PROPAGATION_PATH_SIM_API AdaptiveIntegrationSettings {
bool bActive = true; //!< If this is set to false, the adaptation is bypassed and the integration step size is therefore constant 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 dMaxError = 0.015; //< For errors above this threshold the time step is decreased.
double dUncriticalError = 0.005; double dUncriticalError = 0.005; //< An error below this limit allows to increase the time step again.
unsigned int iMaxAdaptationLevel = 31; //! Maximum times, the time step is halfed. Maximum valid value = 31. unsigned int iMaxAdaptationLevel = 31; //!< Maximum times, the time step is halfed. Maximum valid value = 31.
}; };
struct ITA_PROPAGATION_PATH_SIM_API Settings { struct ITA_PROPAGATION_PATH_SIM_API Settings {
......
...@@ -56,13 +56,17 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pRay) ...@@ -56,13 +56,17 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pRay)
DoubleRayResolution(); DoubleRayResolution();
} }
void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pRay, const int idxMinDist, const VistaVector3D& receiverPosition, const double& threshold) void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pRay, const VistaVector3D& v3ReceiverPosition, const double& dThreshold)
{ {
if (!Contains(pRay)) if (!Contains(pRay))
ITA_EXCEPT_INVALID_PARAMETER("Given ray is not part of the adapted grid."); ITA_EXCEPT_INVALID_PARAMETER("Given ray is not part of the adapted grid.");
const CRayReceiverData* receiverData = pRay->ReceiverDistanceData(v3ReceiverPosition);
if(!receiverData)
ITA_EXCEPT_INVALID_PARAMETER("Given ray has no data about given receiver.");
SetToSurroundingGrid(pRay); SetToSurroundingGrid(pRay);
SetAdvancedRayGridLimits(pRay, idxMinDist, receiverPosition, threshold); SetAdvancedRayGridLimits(pRay, receiverData->idxMinDist, v3ReceiverPosition, dThreshold);
DoubleRayResolution(); DoubleRayResolution();
} }
...@@ -118,37 +122,37 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pZoomRay, const d ...@@ -118,37 +122,37 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pZoomRay, const d
#pragma endregion #pragma endregion
#pragma region NEW RAY LIMITS - private #pragma region NEW RAY LIMITS - private
void CAdaptiveRayGrid::SetAdvancedRayGridLimits(const std::shared_ptr<CRay> pRay, const int idxMinDist, const VistaVector3D& receiverPosition, const double& threshold) void CAdaptiveRayGrid::SetAdvancedRayGridLimits(const std::shared_ptr<CRay> pRay, const int idxMinDist, const VistaVector3D& v3ReceiverPosition, const double& dThreshold)
{ {
if (NPhi() > 3) // There can be more than three rays if pRay is located at a pole if (NPhi() > 3) // There can be more than three rays if pRay is located at a pole
return; return;
std::vector<int> phiIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSameTheta(pRay), idxMinDist, receiverPosition, threshold ); std::vector<int> phiIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSameTheta(pRay), idxMinDist, v3ReceiverPosition, dThreshold );
std::vector<int> thetaIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSamePhi(pRay), idxMinDist, receiverPosition, threshold ); std::vector<int> thetaIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSamePhi(pRay), idxMinDist, v3ReceiverPosition, dThreshold );
FilterDirections(thetaIdxVec, phiIdxVec); FilterDirections(thetaIdxVec, phiIdxVec);
} }
std::vector<int> CAdaptiveRayGrid::FindAdvancedRayGridLimits1D(const std::vector<std::shared_ptr<CRay>>& rayVector, const int idxMinDist, const VistaVector3D& receiverPosition, const double& threshold) const std::vector<int> CAdaptiveRayGrid::FindAdvancedRayGridLimits1D(const std::vector<std::shared_ptr<CRay>>& vpRayVector, const int idxMinDist, const VistaVector3D& v3ReceiverPosition, const double& dThreshold) const
{ {
if (rayVector.size() == 1) if (vpRayVector.size() == 1)
return { 0 }; return { 0 };
if (rayVector.size() == 2) if (vpRayVector.size() == 2)
return {0, 1}; return {0, 1};
if (rayVector.size() != 3) if (vpRayVector.size() != 3)
ITA_EXCEPT_INVALID_PARAMETER("Expected a vector with one to three rays."); ITA_EXCEPT_INVALID_PARAMETER("Expected a vector with one to three rays.");
const std::shared_ptr<CRay>& pRayMin = rayVector[1]; const std::shared_ptr<CRay>& pRayMin = vpRayVector[1];
const std::shared_ptr<CRay>& pRay1 = rayVector[0]; const std::shared_ptr<CRay>& pRay1 = vpRayVector[0];
const std::shared_ptr<CRay>& pRay2 = rayVector[2]; const std::shared_ptr<CRay>& pRay2 = vpRayVector[2];
const double tMin = pRayMin->at(idxMinDist).timeStamp; const double tMin = pRayMin->at(idxMinDist).timeStamp;
VistaVector3D v1 = (receiverPosition - pRay1->AtTime(tMin).position).GetNormalized(); VistaVector3D v1 = (v3ReceiverPosition - pRay1->AtTime(tMin).position).GetNormalized();
VistaVector3D v2 = (receiverPosition - pRay2->AtTime(tMin).position).GetNormalized(); VistaVector3D v2 = (v3ReceiverPosition - pRay2->AtTime(tMin).position).GetNormalized();
VistaVector3D vMin = (receiverPosition - pRayMin->at(idxMinDist).position).GetNormalized(); VistaVector3D vMin = (v3ReceiverPosition - pRayMin->at(idxMinDist).position).GetNormalized();
const float vProd1 = v1.Dot(vMin); const float vProd1 = v1.Dot(vMin);
const float vProd2 = v2.Dot(vMin); const float vProd2 = v2.Dot(vMin);
if (std::abs(vProd1 - vProd2) < threshold) if (std::abs(vProd1 - vProd2) < dThreshold)
return {0, 1, 2}; return {0, 1, 2};
if (vProd1 < vProd2) if (vProd1 < vProd2)
return {0, 1}; return {0, 1};
......
...@@ -2,6 +2,7 @@ ...@@ -2,6 +2,7 @@
// STD // STD
#include <algorithm> #include <algorithm>
#include <cmath>
using namespace ITAPropagationPathSim::AtmosphericRayTracing; using namespace ITAPropagationPathSim::AtmosphericRayTracing;
...@@ -154,7 +155,11 @@ EigenraySearch::RayPtr EigenraySearch::CAdaptiveWorker::Run(const ITAGeo::CStrat ...@@ -154,7 +155,11 @@ EigenraySearch::RayPtr EigenraySearch::CAdaptiveWorker::Run(const ITAGeo::CStrat
m_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder); m_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder);
while (!EigenrayAccuracyReached()) while (!EigenrayAccuracyReached())
{ {
m_adaptiveRayGrid.ZoomIntoRay(m_pMinDistanceRay); if (m_rayAdaptationSettings.advancedRayZooming.bActive)
m_adaptiveRayGrid.ZoomIntoRay(m_pMinDistanceRay, VirtualReceiverPosition(), m_rayAdaptationSettings.advancedRayZooming.threshold);
else
m_adaptiveRayGrid.ZoomIntoRay(m_pMinDistanceRay);
m_simulationEngine.Run(atmosphere, m_adaptiveRayGrid.NewRaysOfLastAdaptation()); m_simulationEngine.Run(atmosphere, m_adaptiveRayGrid.NewRaysOfLastAdaptation());
m_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder); m_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder);
m_nAdaptations++; m_nAdaptations++;
...@@ -220,6 +225,6 @@ void EigenraySearch::CAdaptiveWorker::CalculateEigenraySpreadingLoss(const ITAGe ...@@ -220,6 +225,6 @@ void EigenraySearch::CAdaptiveWorker::CalculateEigenraySpreadingLoss(const ITAGe
const double cReceiver = atmosphere.SpeedOfSound(zReceiver); const double cReceiver = atmosphere.SpeedOfSound(zReceiver);
const double cRef = atmosphere.SpeedOfSound(zRef); const double cRef = atmosphere.SpeedOfSound(zRef);
m_pMinDistanceRay->SetSpreadingLoss( surfaceRef/surfaceReceiver * cReceiver/cRef ); m_pMinDistanceRay->SetSpreadingLoss( std::sqrt(surfaceRef/surfaceReceiver * cReceiver/cRef) );
} }
#pragma endregion #pragma endregion
\ No newline at end of file
...@@ -226,21 +226,27 @@ void Test1DThetaGrid() ...@@ -226,21 +226,27 @@ void Test1DThetaGrid()
void TestAdvancedRayZooming() void TestAdvancedRayZooming()
{ {
const auto sourcePos = VistaVector3D(0, 0, 1000); const auto sourcePos = VistaVector3D(0, 0, 1000);
const auto receiverPos = VistaVector3D(500, 0, 900);
CEquiangularRayDistribution rayGrid = CEquiangularRayDistribution(sourcePos, 7, 12); CEquiangularRayDistribution rayGrid = CEquiangularRayDistribution(sourcePos, 7, 12);
const float dEnd = 1000.0; const float dEnd = 1000.0;
const double tEnd = 10.0; const double tEnd = 10.0;
const int nPoints = 10; const int nPoints = 10;
for each (auto ray in rayGrid.UniqueRays()) for each (auto ray in rayGrid.UniqueRays())
for (int idxPoint = 1; idxPoint < nPoints; idxPoint++) {
ray->Append(rayGrid.SourcePosition() + ray->InitialDirection() * dEnd*idxPoint/nPoints, ray->InitialDirection(), tEnd*idxPoint/nPoints); for (int idxPoint = 1; idxPoint < nPoints; idxPoint++)
{
ray->Append(rayGrid.SourcePosition() + ray->InitialDirection() * dEnd * idxPoint / nPoints, ray->InitialDirection(), tEnd * idxPoint / nPoints);
ray->UpdateMinimumReceiverDistance(receiverPos);
}
ray->FinalizeMinimumReceiverDistances();
}
CAdaptiveRayGrid adaptiveRayGrid; CAdaptiveRayGrid adaptiveRayGrid;
std::shared_ptr<CRay> pRay; std::shared_ptr<CRay> pRay;
adaptiveRayGrid = CAdaptiveRayGrid(rayGrid); adaptiveRayGrid = CAdaptiveRayGrid(rayGrid);
const auto receiverPos = VistaVector3D(500, 0, 900);
pRay = rayGrid.At(3, 0); pRay = rayGrid.At(3, 0);
adaptiveRayGrid.ZoomIntoRay(pRay, 5, receiverPos, 0.1); adaptiveRayGrid.ZoomIntoRay(pRay, receiverPos, 0.1);
CheckAdaptedRayGrid(adaptiveRayGrid, 15, 15, 9); CheckAdaptedRayGrid(adaptiveRayGrid, 15, 15, 9);
} }
......
...@@ -137,16 +137,16 @@ void TestSourceReceiverAzimuthAbove324Deg() ...@@ -137,16 +137,16 @@ void TestSourceReceiverAzimuthAbove324Deg()
int main(int iNumInArgs, char* pcInArgs[]) int main(int iNumInArgs, char* pcInArgs[])
{ {
//Disable multi-threading for debugging purposes //Disable multi-threading for debugging purposes
omp_set_num_threads(1); omp_set_num_threads(1);
TestSourceReceiverAzimuthAbove324Deg();
const CStratifiedAtmosphere homAtmosphere = GetHomogeneousAtmosphere(); const CStratifiedAtmosphere homAtmosphere = GetHomogeneousAtmosphere();
const CStratifiedAtmosphere inhomAtmosphere = GetInhomogeneousAtmosphere(); const CStratifiedAtmosphere inhomAtmosphere = GetInhomogeneousAtmosphere();
TestReceiverNearGroundSourceElevation(homAtmosphere, "Homogeneous"); TestReceiverNearGroundSourceElevation(homAtmosphere, "Homogeneous");
TestReceiverNearGroundSourceElevation(inhomAtmosphere, "Inhomogeneous"); TestReceiverNearGroundSourceElevation(inhomAtmosphere, "Inhomogeneous");
TestReceiverNearGroundSourceAzimuth(inhomAtmosphere, "Inhomogeneous"); TestReceiverNearGroundSourceAzimuth(inhomAtmosphere, "Inhomogeneous");
TestSourceReceiverAzimuthAbove324Deg();
} }