From 6bd0c8a916d18eecd106e754a3fb37eab080e369 Mon Sep 17 00:00:00 2001 From: William Johnson Date: Fri, 26 Jul 2024 16:58:19 -0700 Subject: [PATCH] Fix simple detection limit bugs. The right-clicked energy wasnt always being selected in Simple MDA tool, due to using the wrong/incorrect index to select gamma. Also, fixed a few bugs that could cause crash if no result was found. --- .../app_text/DetectionLimitSimple.xml | 1 + src/DetectionLimitCalc.cpp | 4 +- src/DetectionLimitSimple.cpp | 65 +++++++++++-------- src/DetectionLimitTool.cpp | 5 +- 4 files changed, 45 insertions(+), 30 deletions(-) diff --git a/InterSpec_resources/app_text/DetectionLimitSimple.xml b/InterSpec_resources/app_text/DetectionLimitSimple.xml index b291d56c..c1ead2f2 100644 --- a/InterSpec_resources/app_text/DetectionLimitSimple.xml +++ b/InterSpec_resources/app_text/DetectionLimitSimple.xml @@ -36,6 +36,7 @@ Current state of the Simple MDA tool. Detected: {1}, range [{2}, {3}], @{4} CL + Sorry, no result Excess counts: {1}, range [{2}, {3}], @{4} CL Activity < 0 {1} (observed significantly fewer counts than expected). Excess counts < 0 {2} (observed significantly fewer counts than expected). diff --git a/src/DetectionLimitCalc.cpp b/src/DetectionLimitCalc.cpp index e2f9c538..8b29a171 100644 --- a/src/DetectionLimitCalc.cpp +++ b/src/DetectionLimitCalc.cpp @@ -753,8 +753,8 @@ pair round_roi_to_channels( shared_ptrvalid() ) throw runtime_error( "mda_counts_calc: invalid energy calibration" ); - const float peak_region_lower_ch = cal->channel_for_energy( roi_lower_energy ); - const float peak_region_upper_ch = cal->channel_for_energy( roi_upper_energy ); + const double peak_region_lower_ch = std::max(0.0, cal->channel_for_energy( roi_lower_energy ) ); + const double peak_region_upper_ch = std::max(0.0, cal->channel_for_energy( roi_upper_energy ) ); //if( (peak_region_lower_ch - num_lower_side_channels) < 0.0 ) // throw runtime_error( "mda_counts_calc: lower energy goes off spectrum" ); diff --git a/src/DetectionLimitSimple.cpp b/src/DetectionLimitSimple.cpp index 97fb0e2b..10d190fd 100644 --- a/src/DetectionLimitSimple.cpp +++ b/src/DetectionLimitSimple.cpp @@ -824,6 +824,7 @@ void DetectionLimitSimple::setNuclide( const SandiaDecay::Nuclide *nuc, m_currentEnergy = energy; m_nucEnterController->setNuclideText( nuc ? nuc->symbol : string() ); + assert( m_currentNuclide == nuc ); m_currentNuclide = nuc; if( (age >= 0.0) && !m_nucEnterController->nuclideAgeStr().empty() ) @@ -998,6 +999,8 @@ void DetectionLimitSimple::handleNuclideChanged() ? det->peakResolutionSigma( static_cast(energyToSelect) ) : -1.0; + const float fwhm = std::max( 0.1f, PeakSearchGuiUtils::estimate_FWHM_of_foreground(energyToSelect) ); + size_t transition_index = 0; const SandiaDecay::Transition *transition = nullptr; PeakDef::SourceGammaType sourceGammaType = PeakDef::SourceGammaType::NormalGamma; @@ -1015,7 +1018,7 @@ void DetectionLimitSimple::handleNuclideChanged() } - double minDE = DBL_MAX; + double min_scale_delta_e = DBL_MAX; int currentIndex = -1; for( size_t i = 0; i < photons.size(); ++i ) @@ -1023,28 +1026,30 @@ void DetectionLimitSimple::handleNuclideChanged() double energy = photons[i].energy; const double intensity = photons[i].numPerSecond / dummy_activity; - const double dE = fabs( energyToSelect - energy ); - energy = floor(10000.0*energy + 0.5)/10000.0; - - if( dE < minDE ) - { - minDE = dE; - currentIndex = static_cast( i ); - } - - char text[128] = { '\0' }; - if( i < xrays.size() ) - { - snprintf( text, sizeof(text), "%.4f keV xray I=%.1e", energy, intensity ); - }else - { - snprintf( text, sizeof(text), "%.4f keV I=%.1e", energy, intensity ); - }//if( i < xrays.size() ) - if( intensity > std::numeric_limits::epsilon() ) { + const double delta_e = fabs( energyToSelect - energy ); + const double scale_delta_e = (0.1*fwhm + delta_e) / intensity; + + energy = floor(10000.0*energy + 0.5)/10000.0; + + char text[128] = { '\0' }; + if( i < xrays.size() ) + { + snprintf( text, sizeof(text), "%.4f keV xray I=%.1e", energy, intensity ); + }else + { + snprintf( text, sizeof(text), "%.4f keV I=%.1e", energy, intensity ); + }//if( i < xrays.size() ) + m_photoPeakEnergiesAndBr.push_back( make_pair(energy, intensity) ); m_photoPeakEnergy->addItem( text ); + + if( scale_delta_e < min_scale_delta_e ) + { + min_scale_delta_e = scale_delta_e; + currentIndex = static_cast( m_photoPeakEnergiesAndBr.size() - 1 ); + } }//if( intensity > 0.0 ) }//for each( const double energy, energies ) @@ -1068,7 +1073,7 @@ void DetectionLimitSimple::handleGammaChanged() const int gamma_index = m_photoPeakEnergy->currentIndex(); assert( gamma_index < static_cast(m_photoPeakEnergiesAndBr.size()) ); - if( (gamma_index >= 0) && (gamma_index < static_cast(m_photoPeakEnergiesAndBr.size())) ) + if( (gamma_index >= 0) && (gamma_index < static_cast(m_photoPeakEnergiesAndBr.size())) ) { m_currentEnergy = m_photoPeakEnergiesAndBr[gamma_index].first; }else @@ -1326,14 +1331,19 @@ void DetectionLimitSimple::updateSpectrumDecorationsAndResultText() const bool use_curie = use_curie_units(); const DetectorPeakResponse::EffGeometryType det_geom = drf ? drf->geometryType() : DetectorPeakResponse::EffGeometryType::FarField; - assert( (result->input.num_lower_side_channels != 0) + assert( !result + || (result->input.num_lower_side_channels != 0) || (result->input.num_upper_side_channels != 0) || (result->input.num_lower_side_channels == result->input.num_upper_side_channels) ); - const bool assertedIsBackground = ((result->input.num_lower_side_channels == 0) + const bool assertedIsBackground = (result + && (result->input.num_lower_side_channels == 0) && (result->input.num_upper_side_channels == 0)); - if( result->source_counts > result->decision_threshold ) + if( !result ) + { + result_txt = WString::tr("dls-det-no-result"); + }else if( result->source_counts > result->decision_threshold ) { assert( !assertedIsBackground ); @@ -1396,7 +1406,10 @@ void DetectionLimitSimple::updateSpectrumDecorationsAndResultText() WString full_result_txt; - if( assertedIsBackground ) + if( !result ) + { + full_result_txt = result_txt; + }else if( assertedIsBackground ) { full_result_txt = WString( "{1}" ); }else @@ -1405,14 +1418,14 @@ void DetectionLimitSimple::updateSpectrumDecorationsAndResultText() full_result_txt.arg(result_txt); }//if( assertedIsBackground ) / else - if( gammas_per_bq > 0.0 ) + if( result && (gammas_per_bq > 0.0) ) { const double detection_act = result->detection_limit / gammas_per_bq; const string act = PhysicalUnits::printToBestActivityUnits( detection_act, 2, use_curie ) + DetectorPeakResponse::det_eff_geom_type_postfix( det_geom ); full_result_txt.arg( WString::tr("dls-min-detectable-act").arg(act) ); - }else + }else if( result ) { const string counts = SpecUtils::printCompact(result->detection_limit, 4); full_result_txt.arg( WString::tr("dls-min-detectable-counts").arg( counts ) ); diff --git a/src/DetectionLimitTool.cpp b/src/DetectionLimitTool.cpp index c458d3b9..af4dc2a2 100644 --- a/src/DetectionLimitTool.cpp +++ b/src/DetectionLimitTool.cpp @@ -1488,7 +1488,8 @@ void DetectionLimitTool::update_spectrum_for_currie_result( D3SpectrumDisplayDiv // If number of of side channels is zero, then we expect the following // variables to all be zero - assert( (result->input.num_lower_side_channels != 0) + assert( !result + || (result->input.num_lower_side_channels != 0) || ((result->first_lower_continuum_channel == 0) && (result->last_lower_continuum_channel == 0) && (result->lower_continuum_counts_sum == 0) @@ -1497,7 +1498,7 @@ void DetectionLimitTool::update_spectrum_for_currie_result( D3SpectrumDisplayDiv && (result->upper_continuum_counts_sum == 0)) ); - const bool assertedNoSignal = (result->input.num_lower_side_channels == 0); + const bool assertedNoSignal = (input.num_lower_side_channels == 0); // We will prefer to set the highlight regions from the results, but if we dont have results // we'll do it from the input.