Commit 31f3d4d4 authored by Armin Erraji's avatar Armin Erraji

Implementation of Maekawa and UTD(iwp).

parent 64eb7942
......@@ -82,7 +82,7 @@ void benchmark_single_wedge()
oReceiver->qOrient = VistaQuaternion( 0, 1, 0, 0 );
oReceiver->vPos.SetValues( -.4f, .6f, .1f );
auto pWedge = std::make_shared< ITAGeo::CITADiffractionWedgeApertureBase >();
auto pWedge = std::make_shared< ITAGeo::CITADiffractionOuterWedgeAperture >();
pWedge->v3MainWedgeFaceNormal = VistaVector3D( 1.0f, 1.0f, .0f ).GetNormalized();
pWedge->v3OppositeWedgeFaceNormal = VistaVector3D( -1.0f, 1.0f, .0f ).GetNormalized();
pWedge->v3VertextStart.SetValues( .0f, 1.0f, -1.0f );
......
......@@ -70,27 +70,27 @@ namespace ITAPropagationModels
* @param[in] v3SourcePos Source position
* @param[in] v3TargetPos Target position
* @param[in] pApex Aperture point on wedge
* @param[out] dPropagationLengthDirect Propagation path length without occlusion [m]
* @param[out] dPropagationLengthDetour Propagation path detour length with occlusion [m]
* @param[out] dPropagationLengthFreeField Propagation path length without occlusion [m]
* @param[out] dPropagationLengthOverApex Propagation path detour length with occlusion [m]
*
* @sa IsMaekawaApplicable
*
* @note Throws ITAException on error
*
*/
void ITA_PROPAGATION_MODELS_API GetDirectLengthAndDetourLength( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, double& dPropagationLengthDirect, double& dPropagationLengthDetour );
void ITA_PROPAGATION_MODELS_API GetDirectLengthAndDetourLength( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, double& dPropagationLengthFreeField, double& dPropagationLengthOverApex );
//! Calculates the detour-and-direct-path difference of a propagation path for the Maekawa method
/**
* @param[in] dPropagationLengthDirect Direct path length [m] (sometimes referred to as 'd')
* @param[in] dPropagationLengthDetour Detour path length [m] (sometimes referred to as 'A+B')
* @param[in] dPropagationLengthFreeField Direct path length [m] (sometimes referred to as 'd')
* @param[in] dPropagationLengthOverApex Detour path length [m] (sometimes referred to as 'A+B')
*
* @return Kirchhoff delta (difference) [m] (\delta, sometimes referred to as 'A+B-d')
*
* @note Throws ITAException on error
*
*/
double ITA_PROPAGATION_MODELS_API KirchhoffDelta( const double dPropagationLengthDirect, const double dPropagationLengthDetour );
double ITA_PROPAGATION_MODELS_API KirchhoffDelta( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex );
//! Calculates the detour-and-direct-path inverse difference factor of a propagation path for the Maekawa method
/**
......@@ -99,10 +99,10 @@ namespace ITAPropagationModels
* @note Throws ITAException on error
*
*/
inline double KirchhoffDeltaInverseFactor( const double dPropagationLengthDirect, const double dPropagationLengthDetour )
inline double KirchhoffDeltaInverseFactor( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex )
{
const double dKirchhoffDelta = Maekawa::KirchhoffDelta( dPropagationLengthDirect, dPropagationLengthDetour );
const double dKirchhoffDeltaInverseFactor = ( dPropagationLengthDirect + dKirchhoffDelta ) / dPropagationLengthDirect;
const double dKirchhoffDelta = Maekawa::KirchhoffDelta( dPropagationLengthFreeField, dPropagationLengthOverApex );
const double dKirchhoffDeltaInverseFactor = ( dPropagationLengthFreeField + dKirchhoffDelta ) / dPropagationLengthFreeField;
return dKirchhoffDeltaInverseFactor;
};
......@@ -111,13 +111,13 @@ namespace ITAPropagationModels
* This function calculates the Maekawa diffraction filter coefficient based on detour length around a wedge obstacle
* and the direct path length with and without any obstacles.
*
* @param[in] dPropagationLengthDirect Direct path length [m]
* @param[in] dPropagationLengthDetour Detour path length [m]
* @param[in] dPropagationLengthFreeField Direct path length [m]
* @param[in] dPropagationLengthOverApex Detour path length [m]
* @param[in] fFrequency Frequency [Hz]
* @param[in] fSpeedOfSound Speed of sound [m/s]
*
*/
double ITA_PROPAGATION_MODELS_API CalculateDiffractionFilterCoefficient( const double dPropagationLengthDirect, const double dPropagationLengthDetour, const float fFrequency, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
double ITA_PROPAGATION_MODELS_API CalculateDiffractionFilterCoefficient( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex, const float fFrequency, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
//! Calculate diffraction factor based on the Maekawa detour method
/**
......@@ -125,13 +125,13 @@ namespace ITAPropagationModels
* and the direct path length with and without any obstacles.
* The diffraction factor contains the influence of the wedge without the influence of the propagation path length.
*
* @param[in] dPropagationLengthDirect Direct path length [m]
* @param[in] dPropagationLengthDetour Detour path length [m]
* @param[in] dPropagationLengthFreeField Direct path length [m]
* @param[in] dPropagationLengthOverApex Detour path length [m]
* @param[in] fFrequency Frequency [Hz]
* @param[in] fSpeedOfSound Speed of sound [m/s]
*
*/
double ITA_PROPAGATION_MODELS_API CalculateDiffractionFactor(const double dPropagationLengthDirect, const double dPropagationLengthDetour, const float fFrequency, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
double ITA_PROPAGATION_MODELS_API CalculateDiffractionFactor(const double dPropagationLengthFreeField, const double dPropagationLengthOverApex, const float fFrequency, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
//! Calculate diffraction filter (frequency domain) based on the Maekawa detour method
......@@ -139,13 +139,13 @@ namespace ITAPropagationModels
* This function calculates the Maekawa diffraction filter based on detour length around a wedge obstacle
* and the direct path length with and withont any obstacles.
*
* @param[in] dPropagationLengthDirect Direct path length [m]
* @param[in] dPropagationLengthDetour Detour path length [m]
* @param[in] dPropagationLengthFreeField Direct path length [m]
* @param[in] dPropagationLengthOverApex Detour path length [m]
* @param[out] oTransferFunction Destination filter of simulation (zero-phase frequency domain magnitude spectrum)
* @param[in] fSpeedOfSound Speed of sound [m/s]
*
*/
void ITA_PROPAGATION_MODELS_API CalculateDiffractionFilter( const double dPropagationLengthDirect, const double dPropagationLengthDetour, ITABase::CHDFTSpectrum& oTransferFunction, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
void ITA_PROPAGATION_MODELS_API CalculateDiffractionFilter( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex, ITABase::CHDFTSpectrum& oTransferFunction, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F);
}
}
......
......@@ -40,6 +40,9 @@ namespace ITAPropagationModels
{
namespace UTD
{
using namespace std;
//! UTD diffraction coefficient calculation methods (only approximation implemented yet!)
enum UTD_DIFFRACTION_ALGORITHM
{
......@@ -60,10 +63,24 @@ namespace ITAPropagationModels
* @param[in] v3TargetPos target position (outside wedge)
* @param[in] pApex Aperture point including wedge parameters (angles, etc)
* @param[in] fWaveNumber Wave number ( k )
* @param[out] cfCoeff Complex-valued diffraction coefficient ( D( n, k,\rho, r, \phi_i, \alpha_i, \alpha_d ) )
* @param[out] cfFactor Complex-valued diffraction coefficient ( D( n, k,\rho, r, \phi_i, \alpha_i, \alpha_d ) )
* @param[in] iMethod Approximation algorithm
*/
ITA_PROPAGATION_MODELS_API bool CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fWaveNumber, std::complex< float >& cfCoeff, const int iMethod = UTD_APPROX_KAWAI_KOUYOUMJIAN );
bool ITA_PROPAGATION_MODELS_API CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fFrequency, complex< float >& cfCoeff, const int iMethod = UTD_APPROX_KAWAI_KOUYOUMJIAN, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F );
//! Calculate diffraction factor based on the Maekawa detour method
/**
* This function calculates the Maekawa diffraction factor based on detour length around a wedge obstacle
* and the direct path length with and without any obstacles.
* The diffraction factor contains the influence of the wedge without the influence of the propagation path length.
*
* @param[in] dPropagationLengthFreeField Direct path length [m]
* @param[in] dPropagationLengthOverApex Detour path length [m]
* @param[in] fFrequency Frequency [Hz]
* @param[in] fSpeedOfSound Speed of sound [m/s]
*
*/
bool ITA_PROPAGATION_MODELS_API CalculateDiffractionFactor(const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fFrequency, complex< float >& cfFactor, const int iMethod = UTD_APPROX_KAWAI_KOUYOUMJIAN, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F );
//! Calculates the UTD diffraction filter in frequency domain for a given constellation
/**
......@@ -78,7 +95,7 @@ namespace ITAPropagationModels
* @param[in] iMethod Approximation algorithm
* @param[in] fSpeedOfSound Sound speed
*/
ITA_PROPAGATION_MODELS_API void CalculateDiffractionFilter( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, ITABase::CHDFTSpectrum& oTF, const int iMethod = UTD_APPROX_KAWAI_KOUYOUMJIAN, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F );
void ITA_PROPAGATION_MODELS_API CalculateDiffractionFilter( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, ITABase::CHDFTSpectrum& oTF, const int iMethod = UTD_APPROX_KAWAI_KOUYOUMJIAN, const float fSpeedOfSound = ITAConstants::SPEED_OF_SOUND_F );
}
}
......
......@@ -121,7 +121,7 @@ void CFilterEngine::ApplyDiffractionModel(ITAGeo::CPropagationPathList & oPathLi
//TODO:Change diffraction calculation, just for testing
double dDirectLength, dDetourLength;
Maekawa::GetDirectLengthAndDetourLength(pLastAnchor->v3InteractionPoint, pNextAnchor->v3InteractionPoint, pApex, dDirectLength, dDetourLength);
dDiffractionFactor= Maekawa::CalculateDiffractionFactor(dDirectLength, dDetourLength, fFrequency, ITAConstants::SPEED_OF_SOUND_F);
dDiffractionFactor = Maekawa::CalculateDiffractionFactor(dDirectLength, dDetourLength, fFrequency, ITAConstants::SPEED_OF_SOUND_F);
//Absorption coefficient alpha = 1 - |D|^2
vfAbsorptionCoeffs.push_back(1.0f - pow(dDiffractionFactor,2));
......
......@@ -12,27 +12,27 @@ bool Maekawa::IsApplicable( const VistaVector3D& v3SourcePos, const VistaVector3
return pApex->IsOccluding( v3SourcePos, v3TargetPos );
}
void Maekawa::GetDirectLengthAndDetourLength( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, double& dPropagationLengthDirect, double& dPropagationLengthDetour )
void Maekawa::GetDirectLengthAndDetourLength( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, double& dPropagationLengthFreeField, double& dPropagationLengthOverApex )
{
dPropagationLengthDetour = ( pApex->v3InteractionPoint - v3SourcePos ).GetLength() + ( v3TargetPos - pApex->v3InteractionPoint ).GetLength();
dPropagationLengthDirect = ( v3TargetPos - v3SourcePos ).GetLength();
dPropagationLengthOverApex = ( pApex->v3InteractionPoint - v3SourcePos ).GetLength() + ( v3TargetPos - pApex->v3InteractionPoint ).GetLength();
dPropagationLengthFreeField = ( v3TargetPos - v3SourcePos ).GetLength();
assert( dPropagationLengthDetour > 0.0f );
assert( dPropagationLengthDirect > 0.0f );
assert( dPropagationLengthOverApex > 0.0f );
assert( dPropagationLengthFreeField > 0.0f );
}
double Maekawa::KirchhoffDelta( const double dPropagationLengthDirect, const double dPropagationLengthDetour )
double Maekawa::KirchhoffDelta( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex )
{
if( dPropagationLengthDirect <= 0.0f )
if( dPropagationLengthFreeField <= 0.0f )
ITA_EXCEPT1( INVALID_PARAMETER, "Diffraction direct path length can not be zero or negative" );
if( dPropagationLengthDetour <= 0.0f )
if( dPropagationLengthOverApex <= 0.0f )
ITA_EXCEPT1( INVALID_PARAMETER, "Diffraction detour path length can not be zero or negative" )
if( dPropagationLengthDetour < dPropagationLengthDirect )
if( dPropagationLengthOverApex < dPropagationLengthFreeField )
ITA_EXCEPT1( INVALID_PARAMETER, "Diffraction detour path length can not be smaller than direct path" );
const double dKirchhoffDelta = dPropagationLengthDetour - dPropagationLengthDirect; // Difference in path length, no magic here
const double dKirchhoffDelta = dPropagationLengthOverApex - dPropagationLengthFreeField; // Difference in path length, no magic here
return dKirchhoffDelta;
}
......@@ -60,25 +60,25 @@ void Maekawa::CalculateDiffractionFilter( const double dDirect, const double dDe
}
}
double Maekawa::CalculateDiffractionFilterCoefficient( const double dPropagationLengthDirect, const double dPropagationLengthDetour, const float fFrequency, const float fSpeedOfSound /*= ITAConstants::SPEED_OF_SOUND_F */ )
double Maekawa::CalculateDiffractionFilterCoefficient( const double dPropagationLengthFreeField, const double dPropagationLengthOverApex, const float fFrequency, const float fSpeedOfSound /*= ITAConstants::SPEED_OF_SOUND_F */ )
{
//Difference between detour length and direct length
const double dDifference = dPropagationLengthDetour - dPropagationLengthDirect;
const double dDifference = dPropagationLengthOverApex - dPropagationLengthFreeField;
//Fresnel number N according to Maekawa //TODO: Equation number
const double dN = 2.0 * dDifference / fSpeedOfSound * fFrequency;
//Calculation of the diffraction filter coefficient according to Kurze, Anderson
const double dTemp = sqrt(ITAConstants::TWO_PI_D * dN);
const double dFilterCoefficient = pow(pow(10.0, 5.0 / 20.0)* dTemp / tanh(dTemp), -1) / dPropagationLengthDirect;
const double dFilterCoefficient = pow(pow(10.0, 5.0 / 20.0)* dTemp / tanh(dTemp), -1) / dPropagationLengthFreeField;
return dFilterCoefficient;
}
double Maekawa::CalculateDiffractionFactor(const double dPropagationLengthDirect, const double dPropagationLengthDetour, const float fFrequency, const float fSpeedOfSound /*= ITAConstants::SPEED_OF_SOUND_F */)
double Maekawa::CalculateDiffractionFactor(const double dPropagationLengthFreeField, const double dPropagationLengthOverApex, const float fFrequency, const float fSpeedOfSound /*= ITAConstants::SPEED_OF_SOUND_F */)
{
//Return the diffraction factor without the influence of the propagation path length
return dPropagationLengthDetour * CalculateDiffractionFilterCoefficient(dPropagationLengthDirect, dPropagationLengthDetour, fFrequency, fSpeedOfSound);
return dPropagationLengthOverApex * CalculateDiffractionFilterCoefficient(dPropagationLengthFreeField, dPropagationLengthOverApex, fFrequency, fSpeedOfSound);
}
\ No newline at end of file
......@@ -8,7 +8,7 @@
using namespace ITAConstants;
using namespace ITAPropagationModels;
using namespace std;
//! Solves the Fresnel integral for given value, in literature often referred to as F( X )
std::complex< float > UTDHelperFunctionFresnelIntegralSolution( const float fX );
......@@ -30,7 +30,7 @@ float UTDHelperFunction_a_plus( const float fBeta, const float n );
float UTDHelperFunction_a_minus( const float fBeta, const float n );
void UTD::CalculateDiffractionFilter( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, ITABase::CHDFTSpectrum& oTF, const int iMethod /* = UTD_APPROX_KAWAI_KOUYOUMJIAN */, const float fSpeedOfSound /* = ITAConstants::SPEED_OF_SOUND_F */ )
void UTD::CalculateDiffractionFilter( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, ITABase::CHDFTSpectrum& oTF, const int iMethod /* = UTD_APPROX_KAWAI_KOUYOUMJIAN */, const float fSpeedOfSound /* = ITAConstants::SPEED_OF_SOUND_F */ )
{
std::complex< float > cfCoeff;
oTF.SetCoeffRI( 0, 1.0f, 0.0f ); // DC
......@@ -38,13 +38,12 @@ void UTD::CalculateDiffractionFilter( const VistaVector3D& v3SourcePos, const Vi
{
const float fFrequency = float( i ) * oTF.GetFrequencyResolution();
assert( fFrequency > 0.0f );
const float fWaveNumberRad = 2.0f * PI_F * fFrequency / fSpeedOfSound;
UTD::CalculateDiffractionCoefficient( v3SourcePos, v3TargetPos, pApex, fWaveNumberRad, cfCoeff, iMethod );
UTD::CalculateDiffractionCoefficient( v3SourcePos, v3TargetPos, pApex, fFrequency, cfCoeff, iMethod );
oTF.SetCoeff( i, cfCoeff );
}
}
bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, std::shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fWaveNumber, std::complex< float >& cfCoeff, const int iMethod /* = UTD_APPROX_KAWAI_KOUYOUMJIAN */ )
bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fFrequency, complex< float >& cfCoeff, const int iMethod /* = UTD_APPROX_KAWAI_KOUYOUMJIAN */, const float fSpeedOfSound /* = ITAConstants::SPEED_OF_SOUND_F */)
{
if( ITAGeoUtils::IsPointInsideDiffrationWedge( pApex->v3VertextStart, pApex->v3VertextEnd, pApex->v3MainWedgeFaceNormal, pApex->v3OppositeWedgeFaceNormal, v3SourcePos ) )
return false;
......@@ -61,7 +60,7 @@ bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, con
const float theta_i = pApex->GetIncidenceWaveAngleRad( v3SourcePos );
assert( theta_i != 0.0f );
const float k = fWaveNumber;
const float k = TWO_PI_F_L * fFrequency / fSpeedOfSound;;
assert( k > 0.0f );
const std::complex< float > j( 0.0f, 1.0f ); // Complex value "j" (engineering) or "i" (math, physics)
......@@ -69,6 +68,7 @@ bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, con
const std::complex< float > cfNominator = -1.0f * exp( -j * PI_F / 4.0f );
const std::complex< float > cfDenominator = 2.0f * n * sqrt( 2.0f * k * PI_F ) * sin( theta_i );
//TODO: GetWedgeMainFaceElevationToPointRad() doesn't return the correct angle
const float alpha_i = ITAGeoUtils::GetWedgeMainFaceElevationToPointRad( pApex->v3VertextStart, pApex->v3VertextEnd, pApex->v3MainWedgeFaceNormal, v3SourcePos );
const float alpha_d = ITAGeoUtils::GetWedgeMainFaceElevationToPointRad( pApex->v3VertextStart, pApex->v3VertextEnd, pApex->v3MainWedgeFaceNormal, v3TargetPos );
......@@ -106,13 +106,13 @@ bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, con
{
// Kawai approx, equation 22
const float fX1 = k * L * UTDHelperFunction_a_plus( fAlphaDiff, n );
const float fX1 = k * L * UTDHelperFunction_a_plus(fAlphaSum, n );
cfTerm1 = 1.0f / tan( fCotArg1 ) * UTDHelperFunctionKawaiApproximationFresnelIntegral( fX1 );
const float fX2 = k * L * UTDHelperFunction_a_minus( fAlphaDiff, n );
const float fX2 = k * L * UTDHelperFunction_a_minus(fAlphaSum, n );
cfTerm2 = 1.0f / tan( fCotArg2 ) * UTDHelperFunctionKawaiApproximationFresnelIntegral( fX2 );
const float fX3 = k * L * UTDHelperFunction_a_plus( fAlphaSum, n );
const float fX3 = k * L * UTDHelperFunction_a_plus(fAlphaDiff, n );
cfTerm3 = 1.0f / tan( fCotArg3 ) * UTDHelperFunctionKawaiApproximationFresnelIntegral( fX3 );
const float fX4 = k * L * UTDHelperFunction_a_minus( fAlphaSum, n );
const float fX4 = k * L * UTDHelperFunction_a_minus(fAlphaDiff, n );
cfTerm4 = 1.0f / tan( fCotArg4 ) * UTDHelperFunctionKawaiApproximationFresnelIntegral( fX4 );
}
}
......@@ -135,12 +135,23 @@ bool UTD::CalculateDiffractionCoefficient( const VistaVector3D& v3SourcePos, con
}
// Diffraction coefficient formula
cfCoeff = ( cfNominator / cfDenominator ) * ( cfTerm1 + cfTerm2 + cfTerm3 + cfTerm4 ); // D( .. )
complex<float> cfD = ( cfNominator / cfDenominator ) * ( cfTerm1 + cfTerm2 + cfTerm3 + cfTerm4 ); // D( .. )
// Amplitude scaling factor
float fA = sqrt(rho / (r*(rho + r)))*sin(theta_i);
//Complete diffraction coefficient with A and influence of rho
cfCoeff = cfD * fA /rho;
return true;
}
std::complex< float > UTDHelperFunctionFresnelIntegralSolution( const float fX )
bool CalculateDiffractionFactor(const VistaVector3D& v3SourcePos, const VistaVector3D& v3TargetPos, shared_ptr< const ITAGeo::CITADiffractionWedgeApertureBase > pApex, const float fFrequency, complex< float >& cfFactor, const int iMethod /* = UTD_APPROX_KAWAI_KOUYOUMJIAN*/, const float fSpeedOfSound /* = ITAConstants::SPEED_OF_SOUND_F */)
{
ITA_EXCEPT_NOT_IMPLEMENTED;
}
complex< float > UTDHelperFunctionFresnelIntegralSolution( const float fX )
{
const std::complex< float > j( 0.0f, 1.0f ); // Complex helper value, "j" or "i"
const std::complex< float > cfIntegralResult( 0.0f, 0.0f ); // @todo jst: solve.
......@@ -151,27 +162,28 @@ std::complex< float > UTDHelperFunctionFresnelIntegralSolution( const float fX )
//return cfResult;
}
std::complex< float > UTDHelperFunctionKawaiApproximationFresnelIntegral( const float fX )
complex< float > UTDHelperFunctionKawaiApproximationFresnelIntegral( const float fX )
{
const std::complex< float > j( 0.0f, 1.0f ); // Complex helper value, "j" or "i"
const complex< float > j( 0.0f, 1.0f ); // Complex helper value, "j" or "i"
// Some sanity
assert( fX >= 0.0f );
if( fX < 0.0f )
return 0.0f;
//Phase term is always the same
const complex< float > cfTermPhase = exp(j * PI_F / 4.0f *(1 - sqrt(fX / (fX + 1.4f))));
if( fX < 0.8f )
{
const float fTerm1 = sqrt( PI_F * fX );
const float fTerm2 = 1 - sqrt( fX ) / ( 0.7f * sqrt( fX ) + 1.2f );
const std::complex< float > cfTerm3 = exp( j * PI_F / 4.0f * sqrt( fX / ( fX + 1.4f ) ) );
return ( fTerm1 * fTerm2 * cfTerm3 );
const float fTermAmplitude1 = sqrt( PI_F * fX );
const float fTermAmplitude2 = 1 - sqrt( fX ) / ( 0.7f * sqrt( fX ) + 1.2f );
return ( fTermAmplitude1 * fTermAmplitude2 * cfTermPhase );
}
else
{
const float fTerm1 = 1 - 0.8f / ( ( fX + 1.25f ) * ( fX + 1.25f ) );
const std::complex< float > cfTerm2 = exp( j * PI_F / 4.0f * sqrt( fX / ( fX + 1.4f ) ) );
return ( fTerm1 * cfTerm2 );
const float fTermAmplitude = 1 - 0.8f / ( ( fX + 1.25f ) * ( fX + 1.25f ) );
return ( fTermAmplitude * cfTermPhase);
}
}
......
......@@ -45,14 +45,14 @@ int main(int, char**)
auto pReflection = make_shared< CSpecularReflection >(VistaVector3D(0.0f, 2.0f, 0.0f));
auto pDiffraction = make_shared< CITADiffractionWedgeApertureBase >();
auto pDiffraction = make_shared< CITADiffractionOuterWedgeAperture >();
pDiffraction->v3AperturePoint = VistaVector3D(0.0f, -2.0f, 0.0f);
pDiffraction->v3VertextStart = VistaVector3D(0.0f, -2.0f, 1.0f);
pDiffraction->v3VertextEnd = VistaVector3D(0.0f, -2.0f, -1.0f);
pDiffraction->v3MainWedgeFaceNormal = VistaVector3D(-1.0f, 0.0f, 0.0f);
pDiffraction->v3OppositeWedgeFaceNormal = VistaVector3D(0.0f, 1.0f, 0.0f);
auto pW = std::make_shared< CITADiffractionWedgeApertureBase >();
auto pW = std::make_shared< CITADiffractionOuterWedgeAperture >();
pW->v3AperturePoint.SetValues(0.0f, 1.0f, 0.0f);
pW->v3MainWedgeFaceNormal.SetValues(-1.0f, 1.0f, .0f);
pW->v3MainWedgeFaceNormal.Normalize();
......@@ -76,11 +76,11 @@ int main(int, char**)
oPathDiffraction.push_back(pReceiver);
CPropagationPathList oPathList;
oPathList.push_back(oPathDirect);
/*oPathList.push_back(oPathDirect);
oPathList.push_back(oPathReflection);
oPathList.push_back(oPathDiffraction);
//oPathList.Load("PropagationPathListExample.json");
*/
oPathList.Load("PropagationPathListExample.json");
CFilterEngine oFilter;
......
......@@ -45,7 +45,7 @@ int main( int, char** )
auto pS = std::make_shared< CEmitter >( VistaVector3D( -2.0f, 0, 0 ) );
auto pR = std::make_shared< CSensor >( VistaVector3D( 2.0f, 0, 0 ) );
auto pW = std::make_shared< CITADiffractionWedgeApertureBase >();
auto pW = std::make_shared< CITADiffractionOuterWedgeAperture >();
pW->v3AperturePoint.SetValues( 0.0f, 1.0f, 0.0f );
pW->v3MainWedgeFaceNormal.SetValues( -1.0f, 1.0f, .0f );
pW->v3MainWedgeFaceNormal.Normalize();
......
......@@ -16,10 +16,11 @@ c = 344.0;
maekawa_ir = maekawa_ir_cpp;
maekawa_ir.freqData = ita_diffraction_maekawa( w, s, r, maekawa_ir_cpp.freqVector, c );
utd_ir = maekawa_ir_cpp;
utd_ir.freqData = ita_diffraction_utd( w, s, r, maekawa_ir_cpp.freqVector, c );
%% Combined result
maekawa_eval = ita_merge( maekawa_ir, maekawa_ir_cpp );
maekawa_eval = ita_merge( maekawa_ir, maekawa_ir_cpp,utd_ir );
%% Plot result
......
......@@ -48,7 +48,7 @@ int main( int, char** )
auto pS = std::make_shared< CEmitter >( VistaVector3D( -2.0f, 0.0f, 0.0f ) );
auto pR = std::make_shared< CSensor >( VistaVector3D( 2.0f, 0.0f, 0.0f ) );
auto pW = std::make_shared< CITADiffractionWedgeApertureBase >();
auto pW = std::make_shared< CITADiffractionOuterWedgeAperture >();
pW->v3AperturePoint.SetValues( 0.0f, 1.0f, 0.0f );
pW->v3MainWedgeFaceNormal.SetValues( -1.0f, 1.0f, .0f );
pW->v3MainWedgeFaceNormal.Normalize();
......
......@@ -47,14 +47,14 @@ int main( int, char** )
auto pS = std::make_shared< CEmitter >( VistaVector3D( -2.0f, 0.0f, 0.0f ) );
auto pR = std::make_shared< CSensor >( VistaVector3D( 2.0f, 0.0f, 0.0f ) );
auto pW = std::make_shared< CITADiffractionWedgeApertureBase >();
auto pW = std::make_shared< CITADiffractionOuterWedgeAperture >();
pW->v3AperturePoint.SetValues( 0.0f, 1.0f, 0.0f );
pW->v3MainWedgeFaceNormal.SetValues( -1.0f, 1.0f, .0f );
pW->v3MainWedgeFaceNormal.SetValues(-1.0f, 1.0f, .0f);
pW->v3MainWedgeFaceNormal.Normalize();
pW->v3OppositeWedgeFaceNormal.SetValues( 1.0f, 1.0f, .0f );
pW->v3OppositeWedgeFaceNormal.SetValues(1.0f, 1.0f, .0f);
pW->v3OppositeWedgeFaceNormal.Normalize();
pW->v3VertextStart.SetValues( .0f, 1.0f, -1.0f );
pW->v3VertextEnd.SetValues( .0f, 1.0f, 1.0f );
pW->v3VertextStart.SetValues(.0f, 1.0f, 1.0f);
pW->v3VertextEnd.SetValues(.0f, 1.0f, -1.0f);
CPropagationPath oPropPath;
oPropPath.push_back( pS );
......
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