From b181925b6913e744818947cbc1fb4bf23fb72735 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Thu, 18 Oct 2018 12:02:39 -0600 Subject: [PATCH 1/7] Changed triggering to in Node to be more sensitive, and corrected math in terms of distance allowable away from node for triggering --- glasscore/glasslib/include/Node.h | 8 +++ glasscore/glasslib/src/Node.cpp | 108 +++++++++++++++++------------- 2 files changed, 68 insertions(+), 48 deletions(-) diff --git a/glasscore/glasslib/include/Node.h b/glasscore/glasslib/include/Node.h index 3dc17a9b..99b07f34 100644 --- a/glasscore/glasslib/include/Node.h +++ b/glasscore/glasslib/include/Node.h @@ -398,6 +398,14 @@ class CNode { * NOTE: = sqrt(2)/2) */ static constexpr double k_dGridPointVsResolutionRatio = 0.7071; + + /** + * \brief A factor for allowing distance of triggering relative to the distance + * between detection nodes. A factor of one would just allow a residual change + * corresponding to a center of a node cubiod. 3 allows it to span across a + * kitty corner node in a cubiod plus halfway to the next. + */ + static constexpr double k_residualDistanceAllowanceFactor = 2.; }; } // namespace glasscore #endif // NODE_H diff --git a/glasscore/glasslib/src/Node.cpp b/glasscore/glasslib/src/Node.cpp index c9f3aeba..7f7c0ace 100644 --- a/glasscore/glasslib/src/Node.cpp +++ b/glasscore/glasslib/src/Node.cpp @@ -28,14 +28,14 @@ constexpr double CNode::k_dGridPointVsResolutionRatio; // site Link sorting function // Compares site links using travel times bool sortSiteLink(const SiteLink &lhs, const SiteLink &rhs) { - double travelTime1 = std::get< LINK_TT1>(lhs); + double travelTime1 = std::get < LINK_TT1 > (lhs); if (travelTime1 < 0) { - travelTime1 = std::get< LINK_TT2>(lhs); + travelTime1 = std::get < LINK_TT2 > (lhs); } - double travelTime2 = std::get< LINK_TT1>(rhs); + double travelTime2 = std::get < LINK_TT1 > (rhs); if (travelTime2 < 0) { - travelTime2 = std::get< LINK_TT2>(rhs); + travelTime2 = std::get < LINK_TT2 > (rhs); } // compare @@ -67,7 +67,7 @@ CNode::~CNode() { // ---------------------------------------------------------clear void CNode::clear() { - std::lock_guard nodeGuard(m_NodeMutex); + std::lock_guard < std::recursive_mutex > nodeGuard(m_NodeMutex); clearSiteLinks(); @@ -88,11 +88,11 @@ void CNode::clearSiteLinks() { } // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // remove any links that sites have TO this node for (auto &link : m_vSiteLinkList) { - std::shared_ptr aSite = std::get< LINK_PTR>(link); + std::shared_ptr aSite = std::get < LINK_PTR > (link); aSite->removeNode(getID()); } @@ -103,7 +103,7 @@ void CNode::clearSiteLinks() { // ---------------------------------------------------------initialize bool CNode::initialize(std::string name, double lat, double lon, double z, double resolution, double maxDepth) { - std::lock_guard nodeGuard(m_NodeMutex); + std::lock_guard < std::recursive_mutex > nodeGuard(m_NodeMutex); clear(); @@ -137,7 +137,7 @@ bool CNode::linkSite(std::shared_ptr site, std::shared_ptr node, } // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // Link node to site using traveltime // NOTE: No validation on travel times or distance @@ -172,7 +172,7 @@ bool CNode::unlinkSite(std::shared_ptr site) { std::shared_ptr foundSite; for (const auto &link : m_vSiteLinkList) { // get the site - std::shared_ptr currentSite = std::get< LINK_PTR>(link); + std::shared_ptr currentSite = std::get < LINK_PTR > (link); // if the new station would be before // the current station @@ -226,7 +226,7 @@ bool CNode::unlinkLastSite() { lastSite->removeNode(getID()); // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // unlink last site from node m_vSiteLinkList.pop_back(); @@ -239,7 +239,7 @@ bool CNode::unlinkLastSite() { // ---------------------------------------------------------nucleate std::shared_ptr CNode::nucleate(double tOrigin) { - std::lock_guard nodeGuard(m_NodeMutex); + std::lock_guard < std::recursive_mutex > nodeGuard(m_NodeMutex); // nullchecks // check web @@ -267,10 +267,10 @@ std::shared_ptr CNode::nucleate(double tOrigin) { double dSum = 0.0; int nCount = 0; - std::vector> vPick; + std::vector < std::shared_ptr < CPick >> vPick; // lock mutex for this scope (iterating through the site links) - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // search through each site linked to this node for (const auto &link : m_vSiteLinkList) { @@ -281,7 +281,7 @@ std::shared_ptr CNode::nucleate(double tOrigin) { std::shared_ptr pickBest; // get shared pointer to site - std::shared_ptr site = std::get< LINK_PTR>(link); + std::shared_ptr site = std::get < LINK_PTR > (link); // Ignore if station out of service if (!site->getUse()) { @@ -292,9 +292,9 @@ std::shared_ptr CNode::nucleate(double tOrigin) { } // get traveltime(s) to site - double travelTime1 = std::get< LINK_TT1>(link); - double travelTime2 = std::get< LINK_TT2>(link); - double distDeg = std::get< LINK_DIST>(link); + double travelTime1 = std::get < LINK_TT1 > (link); + double travelTime2 = std::get < LINK_TT2 > (link); + double distDeg = std::get < LINK_DIST > (link); // the minimum and maximum time windows for picks double min = 0.0; @@ -522,29 +522,41 @@ double CNode::getBestSignificance(double tObservedTT, double travelTime1, if (tRes1 > 0) { // calculate the PDF based on the number of SDs we are from mean, by // allowing for WEb Resolution slop converted to seconds - dSig1 = glass3::util::GlassMath::sig( - std::max( - 0.0, - (tRes1 - - travelTime1 / distDeg - * (m_dResolution - + k_dDepthShellResolutionKm) - * glass3::util::Geo::k_KmToDegrees - * k_dGridPointVsResolutionRatio)), - CGlass::k_dNucleationSecondsPerSigma); + dSig1 = + glass3::util::GlassMath::sig( + std::max( + 0.0, + (tRes1 + - (travelTime1 / distDeg) + * (std::sqrt( + 3. + * (m_dResolution + * m_dResolution) + + (k_dDepthShellResolutionKm + * k_dDepthShellResolutionKm)) + * .5) + * glass3::util::Geo::k_KmToDegrees + * k_residualDistanceAllowanceFactor)), + CGlass::k_dNucleationSecondsPerSigma); } double dSig2 = 0; if (tRes2 > 0) { - dSig2 = glass3::util::GlassMath::sig( - std::max( - 0.0, - (tRes2 - - travelTime2 / distDeg - * (m_dResolution - + k_dDepthShellResolutionKm) - * glass3::util::Geo::k_KmToDegrees - * k_dGridPointVsResolutionRatio)), - CGlass::k_dNucleationSecondsPerSigma); + dSig2 = + glass3::util::GlassMath::sig( + std::max( + 0.0, + (tRes2 + - (travelTime2 / distDeg) + * (std::sqrt( + 3. + * (m_dResolution + * m_dResolution) + + (k_dDepthShellResolutionKm + * k_dDepthShellResolutionKm)) + * .5) + * glass3::util::Geo::k_KmToDegrees + * k_residualDistanceAllowanceFactor)), + CGlass::k_dNucleationSecondsPerSigma); } // return the higher of the two significances @@ -562,14 +574,14 @@ std::shared_ptr CNode::getSite(std::string siteID) { } // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // NOTE: could be made more efficient (faster) // if we had a std::map // for all sites for (const auto &link : m_vSiteLinkList) { // get the site - auto aSite = std::get< LINK_PTR>(link); + auto aSite = std::get < LINK_PTR > (link); if (aSite->getSCNL() == siteID) { // found @@ -588,10 +600,10 @@ std::shared_ptr CNode::getLastSite() { } // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); SiteLink lastLink = m_vSiteLinkList[m_vSiteLinkList.size() - 1]; - std::shared_ptr lastSite = std::get< LINK_PTR>(lastLink); + std::shared_ptr lastSite = std::get < LINK_PTR > (lastLink); // found return (lastSite); @@ -600,7 +612,7 @@ std::shared_ptr CNode::getLastSite() { // ---------------------------------------------------------sortSiteLinks void CNode::sortSiteLinks() { // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); // sort sites sort(m_vSiteLinkList.begin(), m_vSiteLinkList.end(), sortSiteLink); @@ -609,13 +621,13 @@ void CNode::sortSiteLinks() { // ---------------------------------------------------------getSitesString std::string CNode::getSitesString() { // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); std::string siteString = ""; // write to station file for (const auto &link : m_vSiteLinkList) { // get the site - std::shared_ptr currentSite = std::get< LINK_PTR>(link); + std::shared_ptr currentSite = std::get < LINK_PTR > (link); double lat, lon, r; currentSite->getGeo().getGeographic(&lat, &lon, &r); @@ -631,7 +643,7 @@ std::string CNode::getSitesString() { // ---------------------------------------------------------getSiteLinksCount int CNode::getSiteLinksCount() const { // lock mutex for this scope - std::lock_guard guard(m_SiteLinkListMutex); + std::lock_guard < std::mutex > guard(m_SiteLinkListMutex); return (m_vSiteLinkList.size()); } @@ -680,13 +692,13 @@ glass3::util::Geo CNode::getGeo() const { // ---------------------------------------------------------getWeb CWeb * CNode::getWeb() const { - std::lock_guard nodeGuard(m_NodeMutex); + std::lock_guard < std::recursive_mutex > nodeGuard(m_NodeMutex); return (m_pWeb); } // ---------------------------------------------------------setWeb void CNode::setWeb(CWeb* web) { - std::lock_guard nodeGuard(m_NodeMutex); + std::lock_guard < std::recursive_mutex > nodeGuard(m_NodeMutex); m_pWeb = web; } From b3a9ed60eba1f1b9e11a7ec18595faaf65775de5 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Thu, 18 Oct 2018 13:22:40 -0600 Subject: [PATCH 2/7] Further tweaking constants for triggering purposes --- glasscore/glasslib/include/Glass.h | 2 +- glasscore/glasslib/include/Hypo.h | 2 +- glasscore/glasslib/include/Node.h | 4 ++-- 3 files changed, 4 insertions(+), 4 deletions(-) diff --git a/glasscore/glasslib/include/Glass.h b/glasscore/glasslib/include/Glass.h index 47908388..eee631cc 100644 --- a/glasscore/glasslib/include/Glass.h +++ b/glasscore/glasslib/include/Glass.h @@ -474,7 +474,7 @@ class CGlass { * \brief The seconds per sigma value used for nucleation, intentionally * a looser sigma value than association */ - static constexpr double k_dNucleationSecondsPerSigma = 3.0; + static constexpr double k_dNucleationSecondsPerSigma = 10.0; /** * \brief The maximum allowed depth diff --git a/glasscore/glasslib/include/Hypo.h b/glasscore/glasslib/include/Hypo.h index e731b87b..10e2deef 100644 --- a/glasscore/glasslib/include/Hypo.h +++ b/glasscore/glasslib/include/Hypo.h @@ -516,7 +516,7 @@ class CHypo { * step size in seconds, default 0.5 seconds * \return Returns a double value containing the final baysian fit. */ - double anneal(int nIter = 250, double dStart = 100.0, double dStop = 1.0, + double anneal(int nIter = 5000, double dStart = 100.0, double dStop = 1.0, double tStart = 5., double tStop = .5); /** diff --git a/glasscore/glasslib/include/Node.h b/glasscore/glasslib/include/Node.h index 99b07f34..9ebdbaea 100644 --- a/glasscore/glasslib/include/Node.h +++ b/glasscore/glasslib/include/Node.h @@ -391,7 +391,7 @@ class CNode { * \brief The depth shell resolution in kilometers used in nucleation * NOTE: this could be calculated from the grid in the future. */ - static constexpr double k_dDepthShellResolutionKm = 10.0; + static constexpr double k_dDepthShellResolutionKm = 100.0; /** * \brief The the ratio between the grid points and the grid resolution @@ -405,7 +405,7 @@ class CNode { * corresponding to a center of a node cubiod. 3 allows it to span across a * kitty corner node in a cubiod plus halfway to the next. */ - static constexpr double k_residualDistanceAllowanceFactor = 2.; + static constexpr double k_residualDistanceAllowanceFactor = 4.; }; } // namespace glasscore #endif // NODE_H From 18490834dc5b5dc742c21b6fbba76bd4f3497337 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Thu, 18 Oct 2018 15:49:58 -0600 Subject: [PATCH 3/7] Changing parameters for triggering --- glasscore/glasslib/include/Glass.h | 2 +- glasscore/glasslib/include/Node.h | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/glasscore/glasslib/include/Glass.h b/glasscore/glasslib/include/Glass.h index eee631cc..3743fb65 100644 --- a/glasscore/glasslib/include/Glass.h +++ b/glasscore/glasslib/include/Glass.h @@ -474,7 +474,7 @@ class CGlass { * \brief The seconds per sigma value used for nucleation, intentionally * a looser sigma value than association */ - static constexpr double k_dNucleationSecondsPerSigma = 10.0; + static constexpr double k_dNucleationSecondsPerSigma = 5.0; /** * \brief The maximum allowed depth diff --git a/glasscore/glasslib/include/Node.h b/glasscore/glasslib/include/Node.h index 9ebdbaea..1b27ecf6 100644 --- a/glasscore/glasslib/include/Node.h +++ b/glasscore/glasslib/include/Node.h @@ -405,7 +405,7 @@ class CNode { * corresponding to a center of a node cubiod. 3 allows it to span across a * kitty corner node in a cubiod plus halfway to the next. */ - static constexpr double k_residualDistanceAllowanceFactor = 4.; + static constexpr double k_residualDistanceAllowanceFactor = 2.; }; } // namespace glasscore #endif // NODE_H From 887c71b78f4e58e7f0b5808b8b04877a45796977 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Fri, 19 Oct 2018 10:52:31 -0600 Subject: [PATCH 4/7] Fixed comment for getBestSig in Node.cpp --- glasscore/glasslib/src/Node.cpp | 23 ++++++----------------- 1 file changed, 6 insertions(+), 17 deletions(-) diff --git a/glasscore/glasslib/src/Node.cpp b/glasscore/glasslib/src/Node.cpp index 7f7c0ace..026d4448 100644 --- a/glasscore/glasslib/src/Node.cpp +++ b/glasscore/glasslib/src/Node.cpp @@ -503,25 +503,14 @@ double CNode::getBestSignificance(double tObservedTT, double travelTime1, // compute significances using residuals and web resolution // should trigger be a looser cutoff than location cutoff - // CNode::nucleate() will ignore anything with a residual > 2.15 sigma. - // since we are using the web resolution as our starting point, that says we - // should be able to pull in anything that is 2.15 times the web resolution - // away, which should let us go plenty far (anything over sqrt(2)/2 * - // resolution should theoretically be closer to another node than it is to - // us, but we're already reducing the slop by up to 8x over what it was) - // even though we have not accounted for pick error. I think this really - // should be some combination of: location error(converted to seconds) - - // sqrt(2)/2 * m_dSurfaceResolution* glass3::util::Geo::k_KmToDegrees * - // tt1.rayparam(dtdx) + 0.5 * m_dDepthResolution*tt1.dtdz tt error - - // tt1.residual_PDF_SD for sigma - observed difference between theoretical - // and actual for tt1 pick error - estimated picking error from pick data. - // but that's more complicated than what we are prepared to deal with at - // this time, so let's go with what's below, and refine it empirically: - // where we subtract location error from tRes1, and then + // The significance is defined in a way that allows for picks to still be + // significant even if an event is not directly on a node. This is done in + // the form of a residual allowance which calculates the maximum off grid + // distance assuming the nodes form a cuboid and multiply but compute + // slowness at that region, then multiplies by a factor (2) for slop. double dSig1 = 0; if (tRes1 > 0) { - // calculate the PDF based on the number of SDs we are from mean, by - // allowing for WEb Resolution slop converted to seconds + dSig1 = glass3::util::GlassMath::sig( std::max( From ac0c001f049d89cf27dc507874a9a9d5547ecc0d Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Fri, 19 Oct 2018 15:34:45 -0600 Subject: [PATCH 5/7] Changed version and disclaimer --- DISCLAIMER.md | 16 ++++++++-------- cmake/version.cmake | 2 +- 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/DISCLAIMER.md b/DISCLAIMER.md index 29afbc12..e951aaff 100644 --- a/DISCLAIMER.md +++ b/DISCLAIMER.md @@ -1,11 +1,11 @@ Disclaimer ========== -This software has been approved for release by the U.S. Geological Survey -(USGS). Although the software has been subjected to rigorous review, the USGS -reserves the right to update the software as needed pursuant to further analysis -and review. No warranty, expressed or implied, is made by the USGS or the U.S. -Government as to the functionality of the software and related material nor shall -the fact of release constitute any such warranty. Furthermore, the software is -released on condition that neither the USGS nor the U.S. Government shall be -held liable for any damages resulting from its authorized or unauthorized use. +This software is preliminary or provisional and is subject to revision. It is +being provided to meet the need for timely best science. The software has not +received final approval by the U.S. Geological Survey (USGS). No warranty, +expressed or implied, is made by the USGS or the U.S. Government as to the +functionality of the software and related material nor shall the fact of release +constitute any such warranty. The software is provided on the condition that +neither the USGS nor the U.S. Government shall be held liable for any damages +resulting from the authorized or unauthorized use of the software. diff --git a/cmake/version.cmake b/cmake/version.cmake index 2cdd670b..7a6b018e 100644 --- a/cmake/version.cmake +++ b/cmake/version.cmake @@ -1,4 +1,4 @@ # version.cmake - a CMake script that defines the overall project version set (PROJECT_VERSION_MAJOR 1) set (PROJECT_VERSION_MINOR 0) -set (PROJECT_VERSION_PATCH 0) +set (PROJECT_VERSION_PATCH 1) From 3c7604d0c965d3de450f408bf519daaab11256f0 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Mon, 22 Oct 2018 08:49:50 -0600 Subject: [PATCH 6/7] Fixed extra whitespace between comment blocks --- glasscore/glasslib/src/Node.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/glasscore/glasslib/src/Node.cpp b/glasscore/glasslib/src/Node.cpp index 026d4448..bea85d30 100644 --- a/glasscore/glasslib/src/Node.cpp +++ b/glasscore/glasslib/src/Node.cpp @@ -502,7 +502,6 @@ double CNode::getBestSignificance(double tObservedTT, double travelTime1, // compute significances using residuals and web resolution // should trigger be a looser cutoff than location cutoff - // The significance is defined in a way that allows for picks to still be // significant even if an event is not directly on a node. This is done in // the form of a residual allowance which calculates the maximum off grid From e103f387f16b5c62d6bd8ca2a0d2ddf73e04a042 Mon Sep 17 00:00:00 2001 From: William L Yeck Date: Mon, 22 Oct 2018 09:05:27 -0600 Subject: [PATCH 7/7] another comment tweak --- glasscore/glasslib/src/Node.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/glasscore/glasslib/src/Node.cpp b/glasscore/glasslib/src/Node.cpp index bea85d30..d8327658 100644 --- a/glasscore/glasslib/src/Node.cpp +++ b/glasscore/glasslib/src/Node.cpp @@ -502,6 +502,7 @@ double CNode::getBestSignificance(double tObservedTT, double travelTime1, // compute significances using residuals and web resolution // should trigger be a looser cutoff than location cutoff + // // The significance is defined in a way that allows for picks to still be // significant even if an event is not directly on a node. This is done in // the form of a residual allowance which calculates the maximum off grid