diff --git a/InterSpec/RelEffChart.h b/InterSpec/RelEffChart.h index 5eadf3dc..63b8a1e3 100644 --- a/InterSpec/RelEffChart.h +++ b/InterSpec/RelEffChart.h @@ -30,12 +30,14 @@ class RelEffChart : public Wt::WContainerWidget void setData( const double live_time, const std::vector &fit_peaks, const std::vector &rel_acts, - const std::string &jsRelEffEqn ); + const std::string &jsRelEffEqn, + const Wt::WString &chi2_title_str ); /** Set data from an "manual" relative efficiency fit. */ void setData( const std::vector &peaks, const std::map> &relActsColors, - std::string relEffEqn ); + std::string relEffEqn, + const Wt::WString &chi2_title_str ); void setLineColor( const Wt::WColor &color ); void setTextColor( const Wt::WColor &color ); diff --git a/InterSpec_resources/RelEffPlot.js b/InterSpec_resources/RelEffPlot.js index e818f6d7..951cda8a 100644 --- a/InterSpec_resources/RelEffPlot.js +++ b/InterSpec_resources/RelEffPlot.js @@ -169,12 +169,17 @@ RelEffPlot.prototype.setYAxisTitle = function( title, dontCallResize ){ }//RelEffPlot.prototype.setYAxisTitle -RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { +RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn, chi2_txt) { const self = this; if( !Array.isArray(data_vals) || (data_vals.length === 0) ) data_vals = null; + console.log( 'chi2_txt:', chi2_txt ); + this.chi2Txt = chi2_txt; + if( (typeof this.chi2Txt !== "string") || (this.chi2Txt.length === 0) ) + this.chi2Txt = null; + const parentWidth = this.chart.clientWidth; const parentHeight = this.chart.clientHeight; @@ -185,6 +190,26 @@ RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { const yaxistitleBB = this.yaxistitle ? this.yaxistitle.node().getBBox() : null; const ytitleh = yaxistitleBB ? yaxistitleBB.height + titlePad : 0; + if( this.chi2Txt ) + { + if( !this.chartInfoTitle ){ + this.chartInfoTitle = this.svg.append("text") + .attr("class", "ChartInfoTitle") + .attr("y",0) // right at the top + .style("text-anchor", "end") // Align text to the end (right) + .attr("dominant-baseline", "hanging") // Align text to the top of the text box + ; + } + + this.chartInfoTitle.attr("x", parentWidth - 10) + .text(this.chi2Txt); + }else if( this.chartInfoTitle ) + { + this.chartInfoTitle.remove(); + this.chartInfoTitle = null; + } + const chi2_txt_pad = this.chartInfoTitle ? 0.5*this.chartInfoTitle.node().getBBox().height : 0; + // We will perform an initial estimate of the chart are - using the axis we have. // Later on we will put current values into the axis, and then get the size, and // update things @@ -195,7 +220,7 @@ RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { let ytickw = yticks.empty() ? 37 : Math.max( 37, d3.max(yticks[0], function(t){ return d3.select(t).node().getBBox().width; }) ); //If we have no labels, we'll get like 7 px, so will require at least 37 px, which is the nominal value we should get let chartAreaWidth = parentWidth - this.options.margins.left - this.options.margins.right - ytitleh - ytickw; - let chartAreaHeight = parentHeight - this.options.margins.top - this.options.margins.bottom - xtitleh - xtickh; + let chartAreaHeight = parentHeight - this.options.margins.top - this.options.margins.bottom - xtitleh - xtickh - chi2_txt_pad; this.data_vals = data_vals; this.fit_eqn = fit_eqn; @@ -260,7 +285,7 @@ RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { ytickw = yticks.empty() ? 7 : d3.max(yticks[0], function(t){ return d3.select(t).node().getBBox().width; }); chartAreaWidth = parentWidth - this.options.margins.left - this.options.margins.right - ytitleh - ytickw; - chartAreaHeight = parentHeight - this.options.margins.top - this.options.margins.bottom - xtitleh - xtickh; + chartAreaHeight = parentHeight - this.options.margins.top - this.options.margins.bottom - xtitleh - xtickh - chi2_txt_pad; this.xScale.range([0, chartAreaWidth]); this.yScale.range([chartAreaHeight, 0]); @@ -279,13 +304,13 @@ RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { if( this.yaxistitle ) this.yaxistitle .attr("y", this.options.margins.left ) - .attr("x",this.options.margins.top - 0.5*chartAreaHeight ); + .attr("x",this.options.margins.top + chi2_txt_pad - 0.5*chartAreaHeight ); this.chartArea.selectAll('.xAxis') .attr("transform", "translate(0," + chartAreaHeight + ")"); this.chartArea .attr("transform", "translate(" + (this.options.margins.left + ytickw + ytitleh) - + "," + this.options.margins.top + ")"); + + "," + (this.options.margins.top + chi2_txt_pad) + ")"); if( fit_eqn ){ // Define the line @@ -405,5 +430,5 @@ RelEffPlot.prototype.setRelEffData = function (data_vals, fit_eqn) { RelEffPlot.prototype.handleResize = function () { - this.setRelEffData(this.data_vals, this.fit_eqn); + this.setRelEffData(this.data_vals, this.fit_eqn, this.chi2_txt ); };//RelEffPlot.prototype.handleResize diff --git a/InterSpec_resources/app_text/RelActManualGui.xml b/InterSpec_resources/app_text/RelActManualGui.xml index 838c4302..e654b550 100644 --- a/InterSpec_resources/app_text/RelActManualGui.xml +++ b/InterSpec_resources/app_text/RelActManualGui.xml @@ -95,4 +95,7 @@ You may want to consider adding "Add. Uncert" for Uranium problems. You are adding "Add. Uncert" - interpret computed uncertainties with care. You are currently adding "Add. Uncert", which will cause computed uncertainties to not be correct. + χ²/dof = {1}/{2}{3}{4} + (pval={1}) + (1-pval={1}) diff --git a/src/PeakSearchGuiUtils.cpp b/src/PeakSearchGuiUtils.cpp index 59206fc3..8485e9dd 100644 --- a/src/PeakSearchGuiUtils.cpp +++ b/src/PeakSearchGuiUtils.cpp @@ -927,7 +927,8 @@ class PeakSelectorWindow : public AuxWindow m_chartPanel->setHidden( (max_num_peaks < 2) ); string relEffEqn = ""; - m_rel_eff_chart->setData( peaks, relActsColors, relEffEqn ); + WString title_chi2_info; + m_rel_eff_chart->setData( peaks, relActsColors, relEffEqn, title_chi2_info ); }//void refreshRelEffChart() #endif // USE_REL_ACT_TOOL diff --git a/src/RelActAutoGui.cpp b/src/RelActAutoGui.cpp index 64db4fed..56e9e8c0 100644 --- a/src/RelActAutoGui.cpp +++ b/src/RelActAutoGui.cpp @@ -4664,8 +4664,66 @@ void RelActAutoGui::updateFromCalc( std::shared_ptrm_rel_eff_coefficients ); const double live_time = answer->m_foreground ? answer->m_foreground->live_time() : 1.0f; + + + WString chi2_title("χ²/dof = {1}/{2}{3}"); + chi2_title.arg( SpecUtils::printCompact(answer->m_chi2, 3) ) + .arg( static_cast(answer->m_dof) ); + + // If we have U or Pu, we'll give the enrichment, or if we have two nuclides we'll + // give their ratio + set isotopes; + for( const auto &relact : answer->m_rel_activities ) + { + if( relact.nuclide ) + isotopes.insert( relact.nuclide ); + } + + const SandiaDecay::SandiaDecayDataBase *db = DecayDataBaseServer::database(); + assert( db ); + const SandiaDecay::Nuclide * const u235 = db->nuclide( "U235" ); + const SandiaDecay::Nuclide * const u238 = db->nuclide( "U238" ); + const SandiaDecay::Nuclide * const pu239 = db->nuclide( "Pu239" ); + const SandiaDecay::Nuclide * const pu240 = db->nuclide( "Pu240" ); + assert( u235 && u238 && pu239 && pu240 ); + + if( (isotopes.count(u235) && isotopes.count(u238) && !isotopes.count(pu239)) + || (isotopes.count(pu239) && isotopes.count(pu240) && !isotopes.count(u235)) ) + { + const SandiaDecay::Nuclide * const iso = isotopes.count(u235) ? u235 : pu239; + const double nominal = answer->mass_enrichment_fraction( iso ); + // TODO: add errors + string enrich = ", " + SpecUtils::printCompact(100.0*nominal, 4) + // + "±" + SpecUtils::printCompact(100.0*error, 4) + + "% " + iso->symbol; + chi2_title.arg( enrich ); + }else if( isotopes.size() == 2 ) + { + const vector &rel_acts = answer->m_rel_activities; + const int num_index = (rel_acts[0].rel_activity > rel_acts[1].rel_activity) ? 1 : 0; + const int denom_index = (num_index ? 0 : 1); + const SandiaDecay::Nuclide * const num_nuc = rel_acts[num_index].nuclide; + const SandiaDecay::Nuclide * const den_nuc = rel_acts[denom_index].nuclide; + assert( num_nuc && den_nuc ); + if( num_nuc && den_nuc ) + { + const double ratio = answer->activity_ratio( num_nuc, den_nuc ); + // TODO: add errors + string ratio_txt = ", act(" + num_nuc->symbol + "/" + den_nuc->symbol + ")=" + + SpecUtils::printCompact(ratio, 4); + + chi2_title.arg( ratio_txt ); + }else + { + chi2_title.arg( "" ); + } + }else + { + chi2_title.arg( "" ); + } - m_rel_eff_chart->setData( live_time, answer->m_fit_peaks, answer->m_rel_activities, rel_eff_eqn_js ); + m_rel_eff_chart->setData( live_time, answer->m_fit_peaks, answer->m_rel_activities, + rel_eff_eqn_js, chi2_title ); bool any_nucs_updated = false; diff --git a/src/RelActManualGui.cpp b/src/RelActManualGui.cpp index 95e55a16..cdf0d01d 100644 --- a/src/RelActManualGui.cpp +++ b/src/RelActManualGui.cpp @@ -25,6 +25,8 @@ #include "rapidxml/rapidxml.hpp" +#include + #include #include #include @@ -1098,7 +1100,7 @@ void RelActManualGui::render( Wt::WFlags flags ) void RelActManualGui::calculateSolution() { m_currentSolution.reset(); - m_chart->setData( vector{}, {}, "" ); + m_chart->setData( vector{}, {}, "", {} ); m_results->clear(); // TODO: should do the actual computation not on the GUI thread! @@ -1584,7 +1586,7 @@ void RelActManualGui::updateGuiWithResults( shared_ptrsetData( vector{}, {}, "" ); + m_chart->setData( vector{}, {}, "", {} ); break; } @@ -1640,8 +1642,79 @@ void RelActManualGui::updateGuiWithResults( shared_ptr 0.99 ) + pval_str = WString::tr("ramg-1-pval").arg( SpecUtils::printCompact(1.0 - prob, 3) ); + else + pval_str = WString::tr("ramg-pval").arg( SpecUtils::printCompact(prob, 3) ); + }catch( std::exception &e ) + { + pval_str = ""; + } + + WString chi2_title = WString::tr("ramg-chart-info-title"); + chi2_title.arg( SpecUtils::printCompact(solution.m_chi2, 3) ) + .arg( static_cast(solution.m_dof) ) + .arg( pval_str ); + + // If we have U or Pu, we'll give the enrichment, or if we have two nuclides we'll + // give their ratio. If we have U and Pu, we wont give enrichment. + set isotopes; + for( const auto &relact : solution.m_rel_activities ) + isotopes.insert( relact.m_isotope ); + if( (isotopes.count("U235") && isotopes.count("U238") && !isotopes.count("Pu239")) + || (isotopes.count("Pu239") && isotopes.count("Pu240") && !isotopes.count("U235")) ) + { + string enrich; + const string iso = isotopes.count("U235") ? "U235" : "Pu239"; + + try + { + const double nominal = solution.mass_fraction(iso); + const double plus = solution.mass_fraction( iso, 1.0 ); + const double minus = solution.mass_fraction( iso, -1.0 ); + const double error = 0.5*( fabs(plus - nominal) + fabs(nominal - minus) ); + enrich = ", " + PhysicalUnits::printValueWithUncertainty(100.0*nominal, 100.0*error, 4) + + "% "+ iso; + }catch( std::exception & ) + { + // We dont have the covariance matrix required to get mass fraction errors + // We'll try leaving out the errors + try + { + const double nominal = solution.mass_fraction(iso); + enrich = ", " + SpecUtils::printCompact(100.0*nominal, 4) + "% "+ iso; + }catch( std::exception & ) + { + // I really dont think we should ever get here + assert( 0 ); + } + }//try / catch + + chi2_title.arg( enrich ); + }else if( isotopes.size() == 2 ) + { + const vector &rel_acts = solution.m_rel_activities; + const int num_index = (rel_acts[0].m_rel_activity > rel_acts[1].m_rel_activity) ? 1 : 0; + const int denom_index = (num_index ? 0 : 1); + const string num_nuc = rel_acts[num_index].m_isotope; + const string den_nuc = rel_acts[denom_index].m_isotope; + const double ratio = solution.activity_ratio(num_nuc, den_nuc); + const double uncert = solution.activity_ratio_uncert(num_nuc, den_nuc); + + string ratio_txt = ", act(" + num_nuc + "/" + den_nuc + ")=" + + PhysicalUnits::printValueWithUncertainty(ratio, uncert, 4); + chi2_title.arg( ratio_txt ); + }else + { + chi2_title.arg( "" ); + } - m_chart->setData( solution.m_input_peak, relActsColors, relEffEqn ); + m_chart->setData( solution.m_input_peak, relActsColors, relEffEqn, chi2_title ); // Now update the text diff --git a/src/RelEffChart.cpp b/src/RelEffChart.cpp index e2e7824f..cc89c9bf 100644 --- a/src/RelEffChart.cpp +++ b/src/RelEffChart.cpp @@ -136,7 +136,8 @@ void RelEffChart::render( Wt::WFlags flags ) void RelEffChart::setData( const double live_time, const vector &fit_peaks, const std::vector &rel_acts, - const std::string &relEffEqn ) + const std::string &relEffEqn, + const Wt::WString &chi2_title_str ) { std::vector peaks; map> relActsColors; @@ -189,13 +190,14 @@ void RelEffChart::setData( const double live_time, peaks.push_back( peak ); }//for( const PeakDef &p : fit_peaks ) - setData( peaks, relActsColors, relEffEqn ); + setData( peaks, relActsColors, relEffEqn, chi2_title_str ); }//void setData( std::vector m_fit_peaks, std::string relEffEqn ) void RelEffChart::setData( const std::vector &peaks, const map> &relActsColors, - string relEffEqn ) + string relEffEqn, + const Wt::WString &chi2_title_str ) { char buffer[512] = { '\0' }; @@ -288,7 +290,10 @@ void RelEffChart::setData( const std::vector if( njson_entries < 1 ) rel_eff_plot_values.str("null"); - const string js = m_jsgraph + ".setRelEffData(" + rel_eff_plot_values.str() + "," + relEffEqn + ");"; + const string js = m_jsgraph + ".setRelEffData(" + rel_eff_plot_values.str() + "," + + relEffEqn + "," + + (chi2_title_str.empty() ? string("null") : chi2_title_str.jsStringLiteral()) + + ");"; if( isRendered() ) doJavaScript( js ); else