0.5.0
Loading...
Searching...
No Matches
NAV::Ambiguity Namespace Reference

Namespaces

namespace  internal

Functions

template<typename DerivedA, typename DerivedQ, typename DerivedL, typename DerivedD>
double calcChi2 (const Eigen::MatrixBase< DerivedA > &a, const Eigen::MatrixBase< DerivedQ > &Q, const Eigen::MatrixBase< DerivedL > &L_LTDL_Q, const Eigen::MatrixBase< DerivedD > &D_LTDL_Q, const Eigen::Index &ncand=2, double factor=1.5)
 Calculates $ \chi^2 $, the size of the ellipsoidal region.
double criticalValueFailureRateLookup (double Pf, double Pf_ILS, size_t nAmb)
 Look-up table for the critical value µ, depending on failure rate Pf_ILS.
template<typename DerivedA, typename DerivedQ>
std::optional< std::tuple< typename DerivedQ::PlainObject, typename DerivedQ::PlainObject, typename DerivedQ::PlainObject, typename DerivedA::PlainObject, typename DerivedA::PlainObject > > decorrelate_ztrafo (const Eigen::MatrixBase< DerivedA > &a, const Eigen::MatrixBase< DerivedQ > &Q)
 Decorrelates the ambiguities.
template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool differenceTest (const Eigen::MatrixBase< DerivedAFix1 > &aFix1, const Eigen::MatrixBase< DerivedAFix2 > &aFix2, const Eigen::MatrixBase< DerivedAFloat > &aFloat, const Eigen::MatrixBase< DerivedQ > &Qa, double criticalValue)
 Difference test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::DifferenceTest.
bool differenceTest (double sqNormFix1, double sqNormFix2, double criticalValue)
 Difference test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::DifferenceTest.
template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ, typename DerivedD>
bool fixedFailureRateRatioTest (double Pf, const Eigen::MatrixBase< DerivedAFix1 > &aFix1, const Eigen::MatrixBase< DerivedAFix2 > &aFix2, const Eigen::MatrixBase< DerivedAFloat > &aFloat, const Eigen::MatrixBase< DerivedQ > &Qa, const Eigen::MatrixBase< DerivedD > &D_LTDL_Q)
 Ratio test with a fixed-failure rate for ambiguity resolution, see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestFailureRate.
template<typename DerivedD>
bool fixedFailureRateRatioTest (double Pf, double sqNormFix1, double sqNormFix2, size_t nAmb, const Eigen::MatrixBase< DerivedD > &D_LTDL_Q, bool validateBootstrappedSuccessRate=true)
 Ratio test with a fixed-failure rate for ambiguity resolution, see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestFailureRate.
template<typename DerivedA, typename DerivedQ, typename DerivedL, typename DerivedD>
std::pair< Eigen::MatrixXd, Eigen::VectorXd > integerLeastSquaresSearch (const Eigen::MatrixBase< DerivedA > &a, const Eigen::MatrixBase< DerivedQ > &Q, const Eigen::MatrixBase< DerivedL > &L_LTDL_Q, const Eigen::MatrixBase< DerivedD > &D_LTDL_Q, const Eigen::Index &numCandidates=2)
 Performs an integer least-squares to search for integer candidates for the ambiguities.
template<typename DerivedA, typename DerivedL, typename DerivedD>
std::pair< Eigen::MatrixXd, Eigen::VectorXd > integerLeastSquaresSearchAndShrink (const Eigen::MatrixBase< DerivedA > &zh, const Eigen::MatrixBase< DerivedL > &L_LTDL_Q, const Eigen::MatrixBase< DerivedD > &D_LTDL_Q, const Eigen::Index &numCandidates=2)
 Performs an integer least-squares to search for integer candidates for the ambiguities using the search-and-shrink technique (MLAMBDA)
template<typename DerivedA, typename DerivedQ>
DerivedA::PlainObject integerSearchBootstrapping (const Eigen::MatrixBase< DerivedA > &a, const Eigen::MatrixBase< DerivedQ > &Qa)
 Searches for the integer ambiguities by integer bootstrapping.
template<typename Derived>
Derived::PlainObject integerSearchRounding (const Eigen::MatrixBase< Derived > &a)
 Searches for the integer ambiguities by integer rounding.
template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool projectorTest (const Eigen::MatrixBase< DerivedAFix1 > &aFix1, const Eigen::MatrixBase< DerivedAFix2 > &aFix2, const Eigen::MatrixBase< DerivedAFloat > &aFloat, const Eigen::MatrixBase< DerivedQ > &Qa, double mu)
 Projector test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::ProjectorTest.
template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool ratioTest (const Eigen::MatrixBase< DerivedAFix1 > &aFix1, const Eigen::MatrixBase< DerivedAFix2 > &aFix2, const Eigen::MatrixBase< DerivedAFloat > &aFloat, const Eigen::MatrixBase< DerivedQ > &Qa, double mu)
 Ratio test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestCriticalValue.
bool ratioTest (double sqNormFix1, double sqNormFix2, double mu)
 Ratio test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestCriticalValue.
template<typename Derived>
double successRateBootstrapping (const Eigen::MatrixBase< Derived > &D_LTDL_Q)
 Calculates the bootstrapped success rate.

Variables

std::unordered_map< double, std::vector< std::vector< double > > > lookupTableFailureProbability

Function Documentation

◆ calcChi2()

template<typename DerivedA, typename DerivedQ, typename DerivedL, typename DerivedD>
double NAV::Ambiguity::calcChi2 ( const Eigen::MatrixBase< DerivedA > & a,
const Eigen::MatrixBase< DerivedQ > & Q,
const Eigen::MatrixBase< DerivedL > & L_LTDL_Q,
const Eigen::MatrixBase< DerivedD > & D_LTDL_Q,
const Eigen::Index & ncand = 2,
double factor = 1.5 )

Calculates $ \chi^2 $, the size of the ellipsoidal region.

Parameters
[in]aFloat ambiguity vector [cycles]
[in]QVariance/covariance matrix of the ambiguities
[in]L_LTDL_QLower-triangular matrix from the L^T * D * L decomposition of Q
[in]D_LTDL_QDiagonal matrix from the L^T * D * L decomposition of Q
[in]ncandRequested number of candidates (default = 2)
[in]factorMultiplication factor for the volume of the resulting search ellipsoid (default = 1.5)
Returns
Size of the ellipsoidal region $ \chi^2 $

Definition at line 129 of file EllipsoidalRegion.hpp.

◆ criticalValueFailureRateLookup()

double NAV::Ambiguity::criticalValueFailureRateLookup ( double Pf,
double Pf_ILS,
size_t nAmb )

Look-up table for the critical value µ, depending on failure rate Pf_ILS.

Parameters
PfFailure probability to not be higher
Pf_ILSUpper bound of the failure probability
nAmbAmount of ambiguities
Returns
Critical value µ

Definition at line 108 of file Validate.cpp.

◆ decorrelate_ztrafo()

template<typename DerivedA, typename DerivedQ>
std::optional< std::tuple< typename DerivedQ::PlainObject, typename DerivedQ::PlainObject, typename DerivedQ::PlainObject, typename DerivedA::PlainObject, typename DerivedA::PlainObject > > NAV::Ambiguity::decorrelate_ztrafo ( const Eigen::MatrixBase< DerivedA > & a,
const Eigen::MatrixBase< DerivedQ > & Q )
nodiscard

Decorrelates the ambiguities.

Parameters
[in]aFloat ambiguity vector [cycles]
[in]QVariance/covariance matrix of the ambiguities
Returns
[Qz, Z, L, D, z] L, D are a (L^T * D * L) decomposition of Q_z
Note
See [11] Chang 2005, Reduction algorithm

Definition at line 94 of file Decorrelation.hpp.

◆ differenceTest() [1/2]

template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool NAV::Ambiguity::differenceTest ( const Eigen::MatrixBase< DerivedAFix1 > & aFix1,
const Eigen::MatrixBase< DerivedAFix2 > & aFix2,
const Eigen::MatrixBase< DerivedAFloat > & aFloat,
const Eigen::MatrixBase< DerivedQ > & Qa,
double criticalValue )

Difference test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::DifferenceTest.

Parameters
[in]aFix1Best integer ambiguities solution
[in]aFix2Second best integer ambiguities solution
[in]aFloatFloat ambiguities
[in]QaCovariance matrix of the ambiguities
[in]criticalValueCritical value to check against (empirically determined value)
Returns
True if the test passes and the best integer solution should be accepted
Note
See [52] Verhagen 2006, eq. 31

Definition at line 51 of file Validate.hpp.

◆ differenceTest() [2/2]

bool NAV::Ambiguity::differenceTest ( double sqNormFix1,
double sqNormFix2,
double criticalValue )

Difference test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::DifferenceTest.

Parameters
[in]sqNormFix1Squared norm of the best integer ambiguities solution
[in]sqNormFix2Squared norm of the second best integer ambiguities solution
[in]criticalValueCritical value to check against (empirically determined value)
Returns
True if the test passes and the best integer solution should be accepted
Note
See [52] Verhagen 2006, eq. 31

Definition at line 98 of file Validate.cpp.

◆ fixedFailureRateRatioTest() [1/2]

template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ, typename DerivedD>
bool NAV::Ambiguity::fixedFailureRateRatioTest ( double Pf,
const Eigen::MatrixBase< DerivedAFix1 > & aFix1,
const Eigen::MatrixBase< DerivedAFix2 > & aFix2,
const Eigen::MatrixBase< DerivedAFloat > & aFloat,
const Eigen::MatrixBase< DerivedQ > & Qa,
const Eigen::MatrixBase< DerivedD > & D_LTDL_Q )

Ratio test with a fixed-failure rate for ambiguity resolution, see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestFailureRate.

Parameters
[in]PfFixed-failure rate
[in]aFix1Best integer ambiguities solution
[in]aFix2Second best integer ambiguities solution
[in]aFloatFloat ambiguities
[in]QaCovariance matrix of the ambiguities
[in]D_LTDL_QDiagonal entries from the L^T * D * L decomposition of Q
Returns
True if the test passes and the best integer solution should be accepted
Note
See [53] Verhagen 2013

Definition at line 123 of file Validate.hpp.

◆ fixedFailureRateRatioTest() [2/2]

template<typename DerivedD>
bool NAV::Ambiguity::fixedFailureRateRatioTest ( double Pf,
double sqNormFix1,
double sqNormFix2,
size_t nAmb,
const Eigen::MatrixBase< DerivedD > & D_LTDL_Q,
bool validateBootstrappedSuccessRate = true )

Ratio test with a fixed-failure rate for ambiguity resolution, see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestFailureRate.

Parameters
[in]PfFixed-failure rate
[in]sqNormFix1Squared norm of the best integer ambiguities solution
[in]sqNormFix2Squared norm of the second best integer ambiguities solution
[in]nAmbNumber of ambiguities
[in]D_LTDL_QDiagonal entries from the L^T * D * L decomposition of Q
[in]validateBootstrappedSuccessRateWhether to check the bootstrapped success rate before doing the ratio test
Returns
True if the test passes and the best integer solution should be accepted
Note
See [53] Verhagen 2013

Definition at line 149 of file Validate.hpp.

◆ integerLeastSquaresSearch()

template<typename DerivedA, typename DerivedQ, typename DerivedL, typename DerivedD>
std::pair< Eigen::MatrixXd, Eigen::VectorXd > NAV::Ambiguity::integerLeastSquaresSearch ( const Eigen::MatrixBase< DerivedA > & a,
const Eigen::MatrixBase< DerivedQ > & Q,
const Eigen::MatrixBase< DerivedL > & L_LTDL_Q,
const Eigen::MatrixBase< DerivedD > & D_LTDL_Q,
const Eigen::Index & numCandidates = 2 )

Performs an integer least-squares to search for integer candidates for the ambiguities.

Parameters
[in]aDecorrelated float ambiguity vector [cycles]
[in]QVariance/covariance matrix of the ambiguities
[in]L_LTDL_QLower-triangular matrix from the L^T * D * L decomposition of Q
[in]D_LTDL_QDiagonal entries from the L^T * D * L decomposition of Q
[in]numCandidatesRequested number of candidates (default = 2)
Returns
Pair of: First a list of integer candidates (last column is the most likely), Vector with squared norm of the integer candidates
Note
See [13] de Jonge 1996, Algorithm FI71

Definition at line 79 of file Search.hpp.

◆ integerLeastSquaresSearchAndShrink()

template<typename DerivedA, typename DerivedL, typename DerivedD>
std::pair< Eigen::MatrixXd, Eigen::VectorXd > NAV::Ambiguity::integerLeastSquaresSearchAndShrink ( const Eigen::MatrixBase< DerivedA > & zh,
const Eigen::MatrixBase< DerivedL > & L_LTDL_Q,
const Eigen::MatrixBase< DerivedD > & D_LTDL_Q,
const Eigen::Index & numCandidates = 2 )

Performs an integer least-squares to search for integer candidates for the ambiguities using the search-and-shrink technique (MLAMBDA)

Parameters
[in]zhDecorrelated float ambiguity vector [cycles]
[in]L_LTDL_QLower-triangular matrix from the L^T * D * L decomposition of Q
[in]D_LTDL_QDiagonal entries from the L^T * D * L decomposition of Q
[in]numCandidatesRequested number of candidates (default = 2)
Returns
Pair of: First a list of integer candidates (last column is the most likely), Vector with squared norm of the integer candidates
Note
See [11] Chang 2005, Argorithm 3.3
See [13] de Jonge 1996, ch. 4.8
See ssearch.m file from TUDelft Matlab code

Definition at line 243 of file Search.hpp.

◆ integerSearchBootstrapping()

template<typename DerivedA, typename DerivedQ>
DerivedA::PlainObject NAV::Ambiguity::integerSearchBootstrapping ( const Eigen::MatrixBase< DerivedA > & a,
const Eigen::MatrixBase< DerivedQ > & Qa )

Searches for the integer ambiguities by integer bootstrapping.

Parameters
[in]aDecorrelated float ambiguity vector [cycles]
[in]QaVariance/covariance matrix of the ambiguities
Returns
Integer ambiguity vector
Note
See [48] Springer Handbook GNSS, ch. 23.2.3

Definition at line 48 of file Search.hpp.

◆ integerSearchRounding()

template<typename Derived>
Derived::PlainObject NAV::Ambiguity::integerSearchRounding ( const Eigen::MatrixBase< Derived > & a)

Searches for the integer ambiguities by integer rounding.

Parameters
[in]aFloat ambiguity vector [cycles]
Returns
Integer ambiguity vector
Note
See [48] Springer Handbook GNSS, ch. 23.2.2

Definition at line 34 of file Search.hpp.

◆ projectorTest()

template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool NAV::Ambiguity::projectorTest ( const Eigen::MatrixBase< DerivedAFix1 > & aFix1,
const Eigen::MatrixBase< DerivedAFix2 > & aFix2,
const Eigen::MatrixBase< DerivedAFloat > & aFloat,
const Eigen::MatrixBase< DerivedQ > & Qa,
double mu )

Projector test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::ProjectorTest.

Parameters
[in]aFix1Best integer ambiguities solution
[in]aFix2Second best integer ambiguities solution
[in]aFloatFloat ambiguities
[in]QaCovariance matrix of the ambiguities
[in]muCritical value to check against, with range (0, 1]
Returns
True if the test passes and the best integer solution should be accepted
Note
See [52] Verhagen 2006, eq. 35

Definition at line 172 of file Validate.hpp.

◆ ratioTest() [1/2]

template<typename DerivedAFix1, typename DerivedAFix2, typename DerivedAFloat, typename DerivedQ>
bool NAV::Ambiguity::ratioTest ( const Eigen::MatrixBase< DerivedAFix1 > & aFix1,
const Eigen::MatrixBase< DerivedAFix2 > & aFix2,
const Eigen::MatrixBase< DerivedAFloat > & aFloat,
const Eigen::MatrixBase< DerivedQ > & Qa,
double mu )

Ratio test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestCriticalValue.

Parameters
[in]aFix1Best integer ambiguities solution
[in]aFix2Second best integer ambiguities solution
[in]aFloatFloat ambiguities
[in]QaCovariance matrix of the ambiguities
[in]muCritical value to check against, with range (0, 1]
Returns
True if the test passes and the best integer solution should be accepted
Note
See [48] Springer Handbook GNSS, ch. 23.6.4, eq. 23.78 or [52] Verhagen 2006, eq. 28, 29

Definition at line 83 of file Validate.hpp.

◆ ratioTest() [2/2]

bool NAV::Ambiguity::ratioTest ( double sqNormFix1,
double sqNormFix2,
double mu )

Ratio test for the ambiguity resolution see NAV::AmbiguityResolutionParameters::ValidationAlgorithm::RatioTestCriticalValue.

Parameters
[in]sqNormFix1Squared norm of the best integer ambiguities solution
[in]sqNormFix2Squared norm of the second best integer ambiguities solution
[in]muCritical value to check against, with range (0, 1]
Returns
True if the test passes and the best integer solution should be accepted
Note
See [48] Springer Handbook GNSS, ch. 23.6.4, eq. 23.78 or [52] Verhagen 2006, eq. 28, 29

Definition at line 103 of file Validate.cpp.

◆ successRateBootstrapping()

template<typename Derived>
double NAV::Ambiguity::successRateBootstrapping ( const Eigen::MatrixBase< Derived > & D_LTDL_Q)

Calculates the bootstrapped success rate.

Parameters
[in]D_LTDL_QDiagonal entries from the L^T * D * L decomposition of Q
Returns
The bootstrapped success rate
Note
See [48] Springer Handbook GNSS, ch. 23.2.4, eq. 23.28

Definition at line 29 of file Validate.hpp.

Variable Documentation

◆ lookupTableFailureProbability

std::unordered_map<double, std::vector<std::vector<double> > > NAV::Ambiguity::lookupTableFailureProbability

Definition at line 27 of file Validate.cpp.