Skip to content

Commit 7974509

Browse files
committed
Add ability to copy-paste DRF eff. fact coefficients into DRF select tool.
1 parent e9a819f commit 7974509

File tree

4 files changed

+153
-76
lines changed

4 files changed

+153
-76
lines changed

InterSpec/DetectorPeakResponse.h

Lines changed: 19 additions & 75 deletions
Original file line numberDiff line numberDiff line change
@@ -63,80 +63,15 @@ namespace rapidxml
6363
template<class Ch> class xml_document;
6464
}//namespace rapidxml
6565

66-
//Three expression evaluators have been tried; all have benefits and drawbacks.
67-
// CLHEP/Evaluator: only added 21 kb to final binary size, but slowest option
68-
// exprtk.hpp: added 5.3 Mb to final binary size, but 30x faster than CLHEP
69-
// muparserx: adds 430 kb to binary size, and about 10x faster than CLHEP
70-
//
71-
// Since muparserx is what TerminalWidget uses, it was decided on 20180408 to
72-
// only use muparserx.
73-
74-
//100000 evaluations of:
75-
// exp(-343.6330974237 + 269.1023287277*log(x) + -83.8077567526*log(x)^2
76-
// + 12.9980559362*log(x)^3 + -1.0068649823*log(x)^4 + 0.0311640084*log(x)^5)
77-
// muparser x took: cpu=0.143294s, wall=0.14341s (1.4 us/eval)
78-
// evaluator x took: cpu=1.10172s, wall=1.10209s (11 us/eval)
66+
/** TODO:
67+
- Add a "fixed" geometry option - then distances in like shielding/source display are not shown/used
68+
- Add a "setback" option
69+
*/
70+
7971

8072
//Forward declarations
81-
namespace mup
82-
{
83-
class Value;
84-
class ParserX;
85-
}
73+
struct FormulaWrapper;
8674

87-
struct FormulaWrapper
88-
{
89-
/** Constructor that takes an equation as a string, and creates a callable
90-
object to evaluate that equation.
91-
92-
\param fcnstr The function to use. Ex: "exp(-1.2 + 3*lox(x) + ...)"
93-
\param isMeV Only used to determine which energy value to use to test if the
94-
function is valid or not; If isMeV is true, uses 0.1, else 100.
95-
*/
96-
FormulaWrapper( const std::string &fcnstr, const bool isMev );
97-
~FormulaWrapper();
98-
99-
float efficiency( const float x );
100-
double operator()( const float x );
101-
102-
/** Finds the variable the user most likely intended to be the energy variable
103-
for detector response functions, by looking for arguments inside
104-
paranthesis.
105-
Returned answer is always lower case. Spaces, tabs, and newlines are all
106-
removed from input string before searching.
107-
Returns "x" if it doesnt find any other candidates.
108-
109-
e.x., assert(find_variable_name("exp(2.0*log(x) - 1.1*log(x)^2)") == "x");
110-
assert(find_variable_name("exp(2.0*log(E) - 1.1*log(E)^2)") == "e");
111-
assert(find_variable_name("3.2*exp(energy)") == "energy");
112-
assert(find_variable_name("") == "x");
113-
assert(find_variable_name("1-energy") == "x"); //no (...), so defaults to "x"
114-
*/
115-
static std::string find_variable_name( std::string eqn );
116-
117-
std::mutex m_mutex;
118-
std::string m_fcnstr;
119-
std::string m_var_name;
120-
121-
std::unique_ptr<mup::ParserX> m_parser;
122-
std::unique_ptr<mup::Value> m_value;
123-
};//struct FormulaWrapper
124-
125-
126-
/*
127-
//ToDo: Add table to database to track which DRF to use for a given detector
128-
// model, or serial number,
129-
struct DrfToUsePreference
130-
{
131-
Wt::Dbo::ptr<InterSpecUser> user;
132-
long int / <DetectorPeakResponse> //Dbo::weak_ptr<DetectorPeakResponse>
133-
enum MatchType{ MatchByDetectorModel, MatchBySerialNumber };
134-
MatchType match_type;
135-
string detector_dentifier;
136-
137-
//could add info such as detector name, or whatever
138-
}
139-
*/
14075

14176
class DetectorPeakResponse
14277
{
@@ -163,8 +98,6 @@ class DetectorPeakResponse
16398
std::cout << det.intrinsicEfficiency(700.0f) << std::endl; //0.219004
16499
std::cout << det.intrinsicEfficiency(800.0f) << std::endl; //0.197117
165100
*/
166-
167-
168101

169102
public:
170103
enum ResolutionFnctForm
@@ -317,6 +250,16 @@ class DetectorPeakResponse
317250
const float lowerEnergy,
318251
const float upperEnergy );
319252

253+
/** Makes a callable funtion from the provided mathematical formula
254+
255+
@param fcn Mathematical formula to evaluate. Example: "exp( 1.2 + 3.2*ln(x) + -2.1*ln(x)^2 )"
256+
@param isMeV Whether units are in MeV, or keV
257+
258+
Throws exception if formula is invalid.
259+
*/
260+
static std::function<float(float)> makeEfficiencyFunctionFromFormula( const std::string &formula,
261+
const bool isMeV );
262+
320263
//fromGadrasDefinition(...): accepts the Efficiency.csv and Detector.dat
321264
// files from GADRAS to define the detector.
322265
//Note that just the Detector.dat file contains enough info to define the
@@ -819,10 +762,11 @@ class DetectorPeakResponse
819762
if( (a.setsValue() || a.isSchema()) && m_efficiencyFormula.size() )
820763
{
821764
const bool isMeV = (m_efficiencyEnergyUnits > 10.0f);
765+
822766
try
823767
{
824-
std::shared_ptr<FormulaWrapper> expression = std::make_shared<FormulaWrapper>(m_efficiencyFormula,isMeV);
825-
m_efficiencyFcn = [expression](float a) -> float { return expression->efficiency(a); };
768+
m_efficiencyFcn
769+
= DetectorPeakResponse::makeEfficiencyFunctionFromFormula( m_efficiencyFormula, isMeV);
826770
}catch( std::exception & )
827771
{
828772
//In principle this shouldnt happen - in practice it might

InterSpec_resources/static_text/detector_edit_help.xml

Lines changed: 21 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -271,7 +271,27 @@ everywhere else.
271271
An example input formula could be:
272272
<div style="color:grey;margin-left:10px;">exp(-343.63 + 269.10*log(x) -83.80*ln(x)^2 + 13.00*pow(log(x),3) -1.01*log(x)^4)</div>
273273
</div>
274-
<div style="margin-top:10px;">
274+
275+
<div style="margin-top:20px;">
276+
As a shortcut, for efficiencies of the form
277+
<div style="color:grey;margin-left:10px;">exp(C<sub>0</sub> + C<sub>1</sub>*log(x) + C<sub>2</sub>*log(x)^2 + ... )</div>
278+
you can enter just the coefficients, and when you click the <b>Set</b> button, the equation will be tranformed for you. This is useful when copy-pasting the coefficients from a spreadsheet. The numbers can be seperated by commas, spaces, tabs, newlines, or semi-colons. For example, if you enter
279+
<code>
280+
<pre>
281+
-9.811E+00
282+
-7.99E-01
283+
4.618E-02
284+
2.81E-03
285+
-1.9432E-02
286+
0.000E+00
287+
0.000E+00
288+
</pre>
289+
</code>
290+
Then once you click the <b>Set</b> button, the text will be transformed into:
291+
<div style="color:grey;margin-left:10px;">exp( -9.811E+00 + -7.99E-01*log(x) + 4.618E-02*log(x)^2 + 2.81E-03*log(x)^3 + -1.9432E-02*log(x)^4 )</div>
292+
</div>
293+
294+
<div style="margin-top:20px;">
275295
After a formula is entered you may <u>click</u> the <b>Set</b> button to display the
276296
function on the chart to ensure the function is reasonable; the y-axis
277297
of this chart will always be set to display form 0 to 1.

src/DetectorPeakResponse.cpp

Lines changed: 62 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -244,6 +244,60 @@ namespace
244244
}//namespace
245245

246246

247+
/** Struct that a string formated equation for evaluation of efficiency.
248+
249+
Three expression evaluators have been tried; all have benefits and drawbacks.
250+
CLHEP/Evaluator: only added 21 kb to final binary size, but slowest option
251+
exprtk.hpp: added 5.3 Mb to final binary size, but 30x faster than CLHEP
252+
muparserx: adds 430 kb to binary size, and about 10x faster than CLHEP
253+
254+
Since muparserx is what TerminalWidget uses, it was decided on 20180408 to
255+
only use muparserx.
256+
257+
100000 evaluations of:
258+
exp(-343.6330974237 + 269.1023287277*log(x) + -83.8077567526*log(x)^2
259+
+ 12.9980559362*log(x)^3 + -1.0068649823*log(x)^4 + 0.0311640084*log(x)^5)
260+
muparser x took: cpu=0.143294s, wall=0.14341s (1.4 us/eval)
261+
evaluator x took: cpu=1.10172s, wall=1.10209s (11 us/eval)
262+
*/
263+
struct FormulaWrapper
264+
{
265+
/** Constructor that takes an equation as a string, and creates a callable
266+
object to evaluate that equation.
267+
268+
\param fcnstr The function to use. Ex: "exp(-1.2 + 3*lox(x) + ...)"
269+
\param isMeV Only used to determine which energy value to use to test if the
270+
function is valid or not; If isMeV is true, uses 0.1, else 100.
271+
*/
272+
FormulaWrapper( const std::string &fcnstr, const bool isMev );
273+
~FormulaWrapper();
274+
275+
float efficiency( const float x );
276+
double operator()( const float x );
277+
278+
/** Finds the variable the user most likely intended to be the energy variable
279+
for detector response functions, by looking for arguments inside
280+
paranthesis.
281+
Returned answer is always lower case. Spaces, tabs, and newlines are all
282+
removed from input string before searching.
283+
Returns "x" if it doesnt find any other candidates.
284+
285+
e.x., assert(find_variable_name("exp(2.0*log(x) - 1.1*log(x)^2)") == "x");
286+
assert(find_variable_name("exp(2.0*log(E) - 1.1*log(E)^2)") == "e");
287+
assert(find_variable_name("3.2*exp(energy)") == "energy");
288+
assert(find_variable_name("") == "x");
289+
assert(find_variable_name("1-energy") == "x"); //no (...), so defaults to "x"
290+
*/
291+
static std::string find_variable_name( std::string eqn );
292+
293+
std::mutex m_mutex;
294+
std::string m_fcnstr;
295+
std::string m_var_name;
296+
297+
std::unique_ptr<mup::ParserX> m_parser;
298+
std::unique_ptr<mup::Value> m_value;
299+
};//struct FormulaWrapper
300+
247301

248302
FormulaWrapper::FormulaWrapper( const std::string &fcnstr, const bool isMev )
249303
: m_fcnstr( fcnstr ), m_var_name( "x" )
@@ -821,6 +875,14 @@ void DetectorPeakResponse::setIntrinsicEfficiencyFormula( const string &fcnstr,
821875
}//void setIntrinsicEfficiencyFormula( const std::string &fcn )
822876

823877

878+
function<float(float)> DetectorPeakResponse::makeEfficiencyFunctionFromFormula(
879+
const std::string &formula,
880+
const bool isMeV )
881+
{
882+
shared_ptr<FormulaWrapper> expression = make_shared<FormulaWrapper>( formula, isMeV );
883+
return [expression](float a) -> float { return expression->efficiency(a); };
884+
}
885+
824886

825887
void DetectorPeakResponse::fromGadrasDefinition( std::istream &csvFile,
826888
std::istream &detDatFile )

src/DrfSelect.cpp

Lines changed: 51 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3532,6 +3532,57 @@ void DrfSelect::setDefineDetector()
35323532
// }//if(...)
35333533

35343534
string fcn = m_detectorManualFunctionText->text().narrow();
3535+
3536+
// Check if someone is copy-pasting in coefficiecnts from a spreadsheet
3537+
// We will allow space tab, comma, semicolon, newline seperators, but each entry
3538+
// must be
3539+
vector<string> potential_coefs;
3540+
SpecUtils::split( potential_coefs, fcn, ",\t; \n\r" );
3541+
vector<double> potential_coef_values;
3542+
size_t last_non_zero_coef = 0;
3543+
try
3544+
{
3545+
for( size_t i = 0; i < potential_coefs.size(); ++i )
3546+
{
3547+
const string &txt = potential_coefs[i];
3548+
size_t idx = 0;
3549+
const double val = stod( txt, &idx ); //throws exception if invalid float
3550+
3551+
// Dont accept inf, NaN, or allow trailing characters
3552+
if( IsInf(val) || IsNan(val) || (idx != txt.size()) )
3553+
throw exception();
3554+
potential_coef_values.push_back( val );
3555+
if( val != 0.0 )
3556+
last_non_zero_coef = i;
3557+
}//
3558+
}catch( std::exception & )
3559+
{
3560+
// If we're here, its an invalid double, or there where characters after the number
3561+
potential_coef_values.clear();
3562+
}//try catch to parse entered text as exactly a list of coefficecnts
3563+
3564+
if( (potential_coef_values.size() > 1) && (last_non_zero_coef > 1) )
3565+
{
3566+
assert( potential_coefs.size() == potential_coef_values.size() );
3567+
assert( !potential_coefs.empty() && (last_non_zero_coef < potential_coefs.size()) );
3568+
3569+
fcn = "";
3570+
3571+
for( size_t i = 0; i < potential_coefs.size(); ++i )
3572+
{
3573+
if( potential_coef_values[i] == 0.0 )
3574+
continue;
3575+
3576+
const string mfunc = (i > 0) ? "*log(x)" : "";
3577+
const string power = (i > 1) ? ("^" + to_string(i)) : string();
3578+
fcn += (fcn.empty() ? "exp( " : " + ") + potential_coefs[i] + mfunc + power;
3579+
}
3580+
fcn += fcn.empty() ? "" : " )";
3581+
3582+
m_detectorManualFunctionText->setText( WString::fromUTF8(fcn) );
3583+
}//if( is_all_coeffs )
3584+
3585+
35353586
const string name = m_detectorManualFunctionName->text().toUTF8();
35363587
string descr = m_detectorManualDescription->text().toUTF8();
35373588
if( descr.empty() )

0 commit comments

Comments
 (0)