...
 
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
//! 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);
//! 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.
void ZoomIntoRay(const std::shared_ptr<CRay> pRay, const double& deltaTheta, const double& deltaPhi);
......
......@@ -29,23 +29,24 @@ namespace ITAPropagationPathSim
namespace EigenraySearch {
struct ITA_PROPAGATION_PATH_SIM_API RayAdaptationSettings {
struct {
int maxNAdaptations = 15; //< Abort after X adaptations of the ray resolution
double minAngleResolutionDeg = 0.001; //< []
int maxNAdaptations = 15; //< Abort after N adaptations of the ray resolution
double minAngleResolutionDeg = 0.001; //< Abort if angle between neighboring rays is below this value []
} abort;
struct {
double maxReceiverRadius = 1; //< [m]
double maxSourceReceiverAngle = 1; //< []
double maxAngleForGeomSpreading = 0.01; //< []
double maxReceiverRadius = 1; //< Maximum radius of receiver sphere [m]
double maxSourceReceiverAngle = 1; //< Maximum angle between vectors from source to receiver point and receiver sphere []
double maxAngleForGeomSpreading = 0.01; //< Angular resolution of rays which is required to calculate spreading loss []
} accuracy;
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;
};
struct ITA_PROPAGATION_PATH_SIM_API RayTracingAbortSettings {
int maxReflectionOrder = 3; //< Maximum considered order of reflections
double maxTime = 30; //< [s]
int maxReflectionOrder = 1; //< Maximum considered order of reflections
double maxTime = 30; //< Maximum propagation time of rays [s]
};
struct ITA_PROPAGATION_PATH_SIM_API Settings {
RayTracingAbortSettings rayTracing;
......
......@@ -62,9 +62,9 @@ namespace ITAPropagationPathSim
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
double dMaxError = 0.015;
double dUncriticalError = 0.005;
unsigned int iMaxAdaptationLevel = 31; //! Maximum times, the time step is halfed. Maximum valid value = 31.
double dMaxError = 0.015; //< For errors above this threshold the time step is decreased.
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.
};
struct ITA_PROPAGATION_PATH_SIM_API Settings {
......
......@@ -56,13 +56,17 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pRay)
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))
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);
SetAdvancedRayGridLimits(pRay, idxMinDist, receiverPosition, threshold);
SetAdvancedRayGridLimits(pRay, receiverData->idxMinDist, v3ReceiverPosition, dThreshold);
DoubleRayResolution();
}
......@@ -118,37 +122,37 @@ void CAdaptiveRayGrid::ZoomIntoRay(const std::shared_ptr<CRay> pZoomRay, const d
#pragma endregion
#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
return;
std::vector<int> phiIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSameTheta(pRay), idxMinDist, receiverPosition, threshold );
std::vector<int> thetaIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSamePhi(pRay), idxMinDist, receiverPosition, threshold );
std::vector<int> phiIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSameTheta(pRay), idxMinDist, v3ReceiverPosition, dThreshold );
std::vector<int> thetaIdxVec = FindAdvancedRayGridLimits1D( GetRaysWithSamePhi(pRay), idxMinDist, v3ReceiverPosition, dThreshold );
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 };
if (rayVector.size() == 2)
if (vpRayVector.size() == 2)
return {0, 1};
if (rayVector.size() != 3)
if (vpRayVector.size() != 3)
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>& pRay1 = rayVector[0];
const std::shared_ptr<CRay>& pRay2 = rayVector[2];
const std::shared_ptr<CRay>& pRayMin = vpRayVector[1];
const std::shared_ptr<CRay>& pRay1 = vpRayVector[0];
const std::shared_ptr<CRay>& pRay2 = vpRayVector[2];
const double tMin = pRayMin->at(idxMinDist).timeStamp;
VistaVector3D v1 = (receiverPosition - pRay1->AtTime(tMin).position).GetNormalized();
VistaVector3D v2 = (receiverPosition - pRay2->AtTime(tMin).position).GetNormalized();
VistaVector3D vMin = (receiverPosition - pRayMin->at(idxMinDist).position).GetNormalized();
VistaVector3D v1 = (v3ReceiverPosition - pRay1->AtTime(tMin).position).GetNormalized();
VistaVector3D v2 = (v3ReceiverPosition - pRay2->AtTime(tMin).position).GetNormalized();
VistaVector3D vMin = (v3ReceiverPosition - pRayMin->at(idxMinDist).position).GetNormalized();
const float vProd1 = v1.Dot(vMin);
const float vProd2 = v2.Dot(vMin);
if (std::abs(vProd1 - vProd2) < threshold)
if (std::abs(vProd1 - vProd2) < dThreshold)
return {0, 1, 2};
if (vProd1 < vProd2)
return {0, 1};
......
......@@ -2,6 +2,7 @@
// STD
#include <algorithm>
#include <cmath>
using namespace ITAPropagationPathSim::AtmosphericRayTracing;
......@@ -154,7 +155,11 @@ EigenraySearch::RayPtr EigenraySearch::CAdaptiveWorker::Run(const ITAGeo::CStrat
m_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder);
while (!EigenrayAccuracyReached())
{
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_pMinDistanceRay = FindMinimumDistanceRay(m_adaptiveRayGrid.UniqueRays(), m_iActiveReflexionOrder);
m_nAdaptations++;
......@@ -220,6 +225,6 @@ void EigenraySearch::CAdaptiveWorker::CalculateEigenraySpreadingLoss(const ITAGe
const double cReceiver = atmosphere.SpeedOfSound(zReceiver);
const double cRef = atmosphere.SpeedOfSound(zRef);
m_pMinDistanceRay->SetSpreadingLoss( surfaceRef/surfaceReceiver * cReceiver/cRef );
m_pMinDistanceRay->SetSpreadingLoss( std::sqrt(surfaceRef/surfaceReceiver * cReceiver/cRef) );
}
#pragma endregion
\ No newline at end of file
......@@ -226,21 +226,27 @@ void Test1DThetaGrid()
void TestAdvancedRayZooming()
{
const auto sourcePos = VistaVector3D(0, 0, 1000);
const auto receiverPos = VistaVector3D(500, 0, 900);
CEquiangularRayDistribution rayGrid = CEquiangularRayDistribution(sourcePos, 7, 12);
const float dEnd = 1000.0;
const double tEnd = 10.0;
const int nPoints = 10;
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);
{
ray->Append(rayGrid.SourcePosition() + ray->InitialDirection() * dEnd * idxPoint / nPoints, ray->InitialDirection(), tEnd * idxPoint / nPoints);
ray->UpdateMinimumReceiverDistance(receiverPos);
}
ray->FinalizeMinimumReceiverDistances();
}
CAdaptiveRayGrid adaptiveRayGrid;
std::shared_ptr<CRay> pRay;
adaptiveRayGrid = CAdaptiveRayGrid(rayGrid);
const auto receiverPos = VistaVector3D(500, 0, 900);
pRay = rayGrid.At(3, 0);
adaptiveRayGrid.ZoomIntoRay(pRay, 5, receiverPos, 0.1);
adaptiveRayGrid.ZoomIntoRay(pRay, receiverPos, 0.1);
CheckAdaptedRayGrid(adaptiveRayGrid, 15, 15, 9);
}
......
......@@ -140,13 +140,13 @@ int main(int iNumInArgs, char* pcInArgs[])
//Disable multi-threading for debugging purposes
omp_set_num_threads(1);
TestSourceReceiverAzimuthAbove324Deg();
const CStratifiedAtmosphere homAtmosphere = GetHomogeneousAtmosphere();
const CStratifiedAtmosphere inhomAtmosphere = GetInhomogeneousAtmosphere();
TestReceiverNearGroundSourceElevation(homAtmosphere, "Homogeneous");
TestReceiverNearGroundSourceElevation(inhomAtmosphere, "Inhomogeneous");
TestReceiverNearGroundSourceAzimuth(inhomAtmosphere, "Inhomogeneous");
TestSourceReceiverAzimuthAbove324Deg();
}