diff --git a/Code/Mantid/Framework/MDEvents/inc/MantidMDEvents/MDHistoWorkspace.h b/Code/Mantid/Framework/MDEvents/inc/MantidMDEvents/MDHistoWorkspace.h index 07be4fed3cd9d1f89971133bd5fbd04412cafa69..451e70ea3aa1ab3e434cea0b93d38e0e15c3f983 100644 --- a/Code/Mantid/Framework/MDEvents/inc/MantidMDEvents/MDHistoWorkspace.h +++ b/Code/Mantid/Framework/MDEvents/inc/MantidMDEvents/MDHistoWorkspace.h @@ -66,22 +66,30 @@ namespace MDEvents // -------------------------------------------------------------------------------------------- MDHistoWorkspace & operator+=(const MDHistoWorkspace & b); - MDHistoWorkspace & operator+=(const Mantid::DataObjects::WorkspaceSingleValue & b); + void add(const MDHistoWorkspace & b); + void add(const signal_t signal, const signal_t error); MDHistoWorkspace & operator-=(const MDHistoWorkspace & b); - MDHistoWorkspace & operator-=(const Mantid::DataObjects::WorkspaceSingleValue & b); + void subtract(const MDHistoWorkspace & b); + void subtract(const signal_t signal, const signal_t error); MDHistoWorkspace & operator*=(const MDHistoWorkspace & b); - MDHistoWorkspace & operator*=(const Mantid::DataObjects::WorkspaceSingleValue & b); + void multiply(const MDHistoWorkspace & b); + void multiply(const signal_t signal, const signal_t error); MDHistoWorkspace & operator/=(const MDHistoWorkspace & b); - MDHistoWorkspace & operator/=(const Mantid::DataObjects::WorkspaceSingleValue & b); + void divide(const MDHistoWorkspace & b); + void divide(const signal_t signal, const signal_t error); MDHistoWorkspace & operator&=(const MDHistoWorkspace & b); MDHistoWorkspace & operator|=(const MDHistoWorkspace & b); MDHistoWorkspace & operator^=(const MDHistoWorkspace & b); void operatorNot(); + void log(); + void log10(); + void exp(); + void power(double exponent); diff --git a/Code/Mantid/Framework/MDEvents/src/MDHistoWorkspace.cpp b/Code/Mantid/Framework/MDEvents/src/MDHistoWorkspace.cpp index 4b7858a2cfd369e05ad4b2aad72008748b6dab71..1316aae60c6867e71b6dd5be44cdcffdd4014e37 100644 --- a/Code/Mantid/Framework/MDEvents/src/MDHistoWorkspace.cpp +++ b/Code/Mantid/Framework/MDEvents/src/MDHistoWorkspace.cpp @@ -372,33 +372,42 @@ namespace MDEvents * @return *this after operation */ MDHistoWorkspace & MDHistoWorkspace::operator+=(const MDHistoWorkspace & b) { - checkWorkspaceSize(b, "+="); + add(b); + return *this; + } + + //---------------------------------------------------------------------------------------------- + /** Perform the += operation, element-by-element, for two MDHistoWorkspace's + * + * @param b :: workspace on the RHS of the operation + * */ + void MDHistoWorkspace::add(const MDHistoWorkspace & b) + { + checkWorkspaceSize(b, "add"); for (size_t i=0; i<m_length; ++i) { m_signals[i] += b.m_signals[i]; m_errorsSquared[i] += b.m_errorsSquared[i]; } - return *this; } //---------------------------------------------------------------------------------------------- /** Perform the += operation with a scalar as the RHS argument * - * @param b :: WorkspaceSingleValue (signal and error) as the RHS argument - * @return *this after operation */ - MDHistoWorkspace & MDHistoWorkspace::operator+=(const Mantid::DataObjects::WorkspaceSingleValue & b) + * @param signal :: signal to apply + * @param error :: error (not squared) to apply + * */ + void MDHistoWorkspace::add(const signal_t signal, const signal_t error) { - signal_t signal = b.dataY(0)[0]; - signal_t errorSquared = b.dataE(0)[0]; - errorSquared *= errorSquared; + signal_t errorSquared = error * error; for (size_t i=0; i<m_length; ++i) { m_signals[i] += signal; m_errorsSquared[i] += errorSquared; } - return *this; } + //---------------------------------------------------------------------------------------------- /** Perform the -= operation, element-by-element, for two MDHistoWorkspace's * @@ -406,31 +415,39 @@ namespace MDEvents * @return *this after operation */ MDHistoWorkspace & MDHistoWorkspace::operator-=(const MDHistoWorkspace & b) { - checkWorkspaceSize(b, "-="); + subtract(b); + return *this; + } + + //---------------------------------------------------------------------------------------------- + /** Perform the -= operation, element-by-element, for two MDHistoWorkspace's + * + * @param b :: workspace on the RHS of the operation + * */ + void MDHistoWorkspace::subtract(const MDHistoWorkspace & b) + { + checkWorkspaceSize(b, "subtract"); for (size_t i=0; i<m_length; ++i) { m_signals[i] -= b.m_signals[i]; m_errorsSquared[i] += b.m_errorsSquared[i]; } - return *this; } //---------------------------------------------------------------------------------------------- /** Perform the -= operation with a scalar as the RHS argument * - * @param b :: WorkspaceSingleValue (signal and error) as the RHS argument - * @return *this after operation */ - MDHistoWorkspace & MDHistoWorkspace::operator-=(const Mantid::DataObjects::WorkspaceSingleValue & b) + * @param signal :: signal to apply + * @param error :: error (not squared) to apply + * */ + void MDHistoWorkspace::subtract(const signal_t signal, const signal_t error) { - signal_t signal = b.dataY(0)[0]; - signal_t errorSquared = b.dataE(0)[0]; - errorSquared *= errorSquared; + signal_t errorSquared = error * error; for (size_t i=0; i<m_length; ++i) { m_signals[i] -= signal; m_errorsSquared[i] += errorSquared; } - return *this; } //---------------------------------------------------------------------------------------------- @@ -443,7 +460,21 @@ namespace MDEvents * @return *this after operation */ MDHistoWorkspace & MDHistoWorkspace::operator*=(const MDHistoWorkspace & b_ws) { - checkWorkspaceSize(b_ws, "*="); + multiply(b_ws); + return *this; + } + + //---------------------------------------------------------------------------------------------- + /** Perform the *= operation, element-by-element, for two MDHistoWorkspace's + * + * Error propagation of \f$ f = a * b \f$ is given by: + * \f$ df^2 = f^2 * (da^2 / a^2 + db^2 / b^2) \f$ + * + * @param b :: workspace on the RHS of the operation + * */ + void MDHistoWorkspace::multiply(const MDHistoWorkspace & b_ws) + { + checkWorkspaceSize(b_ws, "multiply"); for (size_t i=0; i<m_length; ++i) { signal_t a = m_signals[i]; @@ -458,7 +489,6 @@ namespace MDEvents m_signals[i] = f; m_errorsSquared[i] = df2; } - return *this; } //---------------------------------------------------------------------------------------------- @@ -467,12 +497,13 @@ namespace MDEvents * Error propagation of \f$ f = a * b \f$ is given by: * \f$ df^2 = f^2 * (da^2 / a^2 + db^2 / b^2) \f$ * - * @param b_ws :: WorkspaceSingleValue (signal and error) as the RHS argument + * @param signal :: signal to apply + * @param error :: error (not squared) to apply * @return *this after operation */ - MDHistoWorkspace & MDHistoWorkspace::operator*=(const Mantid::DataObjects::WorkspaceSingleValue & b_ws) + void MDHistoWorkspace:: multiply(const signal_t signal, const signal_t error) { - signal_t b = b_ws.dataY(0)[0]; - signal_t db2 = b_ws.dataE(0)[0] * b_ws.dataE(0)[0]; + signal_t b = signal; + signal_t db2 = error * error; signal_t db2_relative = db2 / (b*b); for (size_t i=0; i<m_length; ++i) { @@ -485,7 +516,6 @@ namespace MDEvents m_signals[i] = f; m_errorsSquared[i] = df2; } - return *this; } //---------------------------------------------------------------------------------------------- @@ -498,7 +528,22 @@ namespace MDEvents * @return *this after operation */ MDHistoWorkspace & MDHistoWorkspace::operator/=(const MDHistoWorkspace & b_ws) { - checkWorkspaceSize(b_ws, "/="); + divide(b_ws); + return *this; + } + + + //---------------------------------------------------------------------------------------------- + /** Perform the /= operation, element-by-element, for two MDHistoWorkspace's + * + * Error propagation of \f$ f = a / b \f$ is given by: + * \f$ df^2 = f^2 * (da^2 / a^2 + db^2 / b^2) \f$ + * + * @param b_ws :: workspace on the RHS of the operation + **/ + void MDHistoWorkspace::divide(const MDHistoWorkspace & b_ws) + { + checkWorkspaceSize(b_ws, "divide"); for (size_t i=0; i<m_length; ++i) { signal_t a = m_signals[i]; @@ -513,9 +558,10 @@ namespace MDEvents m_signals[i] = f; m_errorsSquared[i] = df2; } - return *this; } + + //---------------------------------------------------------------------------------------------- /** Perform the /= operation with a scalar as the RHS argument * @@ -523,11 +569,11 @@ namespace MDEvents * \f$ df^2 = f^2 * (da^2 / a^2 + db^2 / b^2) \f$ * * @param b_ws :: WorkspaceSingleValue (signal and error) as the RHS argument - * @return *this after operation */ - MDHistoWorkspace & MDHistoWorkspace::operator/=(const Mantid::DataObjects::WorkspaceSingleValue & b_ws) + **/ + void MDHistoWorkspace::divide(const signal_t signal, const signal_t error) { - signal_t b = b_ws.dataY(0)[0]; - signal_t db2 = b_ws.dataE(0)[0] * b_ws.dataE(0)[0]; + signal_t b = signal; + signal_t db2 = error * error; signal_t db2_relative = db2 / (b*b); for (size_t i=0; i<m_length; ++i) { @@ -540,7 +586,6 @@ namespace MDEvents m_signals[i] = f; m_errorsSquared[i] = df2; } - return *this; } @@ -613,6 +658,76 @@ namespace MDEvents } } + //---------------------------------------------------------------------------------------------- + /** Perform the natural logarithm on each signal in the workspace. + * + * Error propagation of \f$ f = ln(a) \f$ is given by: + * \f$ df^2 = a^2 / da^2 \f$ + */ + void MDHistoWorkspace::log() + { + for (size_t i=0; i<m_length; ++i) + { + signal_t a = m_signals[i]; + signal_t da2 = m_errorsSquared[i]; + m_signals[i] = std::log(a); + m_errorsSquared[i] = da2 / (a*a); + } + } + + //---------------------------------------------------------------------------------------------- + /** Perform the base-10 logarithm on each signal in the workspace. + * + * Error propagation of \f$ f = ln(a) \f$ is given by: + * \f$ df^2 = a^2 / da^2 \f$ + */ + void MDHistoWorkspace::log10() + { + for (size_t i=0; i<m_length; ++i) + { + signal_t a = m_signals[i]; + signal_t da2 = m_errorsSquared[i]; + m_signals[i] = std::log10(a); + m_errorsSquared[i] = da2 / (a*a); + } + } + + //---------------------------------------------------------------------------------------------- + /** Perform the exp() function on each signal in the workspace. + * + * Error propagation of \f$ f = exp(a) \f$ is given by: + * \f$ df^2 = f^2 * da^2 \f$ + */ + void MDHistoWorkspace::exp() + { + for (size_t i=0; i<m_length; ++i) + { + signal_t f = std::exp(m_signals[i]); + signal_t da2 = m_errorsSquared[i]; + m_signals[i] = f; + m_errorsSquared[i] = f*f * da2; + } + } + + //---------------------------------------------------------------------------------------------- + /** Perform the power function (signal^exponent) on each signal S in the workspace. + * + * Error propagation of \f$ f = a^b \f$ is given by: + * \f$ df^2 = f^2 * b^2 * (da^2 / a^2) \f$ + */ + void MDHistoWorkspace::power(double exponent) + { + double exponent_squared = exponent * exponent; + for (size_t i=0; i<m_length; ++i) + { + signal_t a = m_signals[i]; + signal_t f = std::pow(a, exponent); + signal_t da2 = m_errorsSquared[i]; + m_signals[i] = f; + m_errorsSquared[i] = f*f * exponent_squared * da2 / (a*a); + } + } + } // namespace Mantid } // namespace MDEvents diff --git a/Code/Mantid/Framework/MDEvents/test/MDHistoWorkspaceTest.h b/Code/Mantid/Framework/MDEvents/test/MDHistoWorkspaceTest.h index 54eb7e0e1508f9bc9be68a7d743f7abb72ba791f..94eefc77fd6b50fb3207abefd2da9b7f5b5def50 100644 --- a/Code/Mantid/Framework/MDEvents/test/MDHistoWorkspaceTest.h +++ b/Code/Mantid/Framework/MDEvents/test/MDHistoWorkspaceTest.h @@ -411,8 +411,7 @@ public: void test_plus_scalar() { MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.0, 2, 5, 10.0, 2.5 /*errorSquared*/); - WorkspaceSingleValue b(3.0, sqrt(3.5) ); - *a += b; + a->add(3.0, sqrt(3.5)); checkWorkspace(a, 5.0, 6.0); } @@ -428,8 +427,7 @@ public: void test_minus_scalar() { MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(3.0, 2, 5, 10.0, 2.5 /*errorSquared*/); - WorkspaceSingleValue b(2.0, sqrt(3.5) ); - *a -= b; + a->subtract(2.0, sqrt(3.5)); checkWorkspace(a, 1.0, 6.0); } @@ -446,13 +444,12 @@ public: void test_times_scalar() { MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.0, 2, 5, 10.0, 2.0 /*errorSquared*/); - WorkspaceSingleValue b(3.0, sqrt(3.0) ); - *a *= b; + a->multiply(3.0, sqrt(3.0)); checkWorkspace(a, 6.0, 36. * (.5 + 1./3.)); // Scalar without error MDHistoWorkspace_sptr d = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.0, 2, 5, 10.0, 2.0 /*errorSquared*/); WorkspaceSingleValue e(3.0, 0); - *d *= e; + d->multiply(3.0, 0); checkWorkspace(d, 6.0, 9 * 2.0); } @@ -469,11 +466,42 @@ public: void test_divide_scalar() { MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(3.0, 2, 5, 10.0, 3.0 /*errorSquared*/); - WorkspaceSingleValue b(2.0, sqrt(2.0) ); - *a /= b; + a->divide(2.0, sqrt(2.0)); checkWorkspace(a, 1.5, 1.5 * 1.5 * (.5 + 1./3.)); } + //-------------------------------------------------------------------------------------- + void test_exp() + { + MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.0, 2, 5, 10.0, 3.0); + a->exp(); + checkWorkspace(a, std::exp(2.0), std::exp(2.0)*std::exp(2.0) * 3.0 ); + } + + //-------------------------------------------------------------------------------------- + void test_log() + { + MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.71828, 2, 5, 10.0, 3.0); + a->log(); + checkWorkspace(a, 1.0, 3.0/(2.71828*2.71828)); + } + + //-------------------------------------------------------------------------------------- + void test_log10() + { + MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(10.0, 2, 5, 10.0, 3.0); + a->log10(); + checkWorkspace(a, 1.0, 3./100.); + } + + //-------------------------------------------------------------------------------------- + void test_power() + { + MDHistoWorkspace_sptr a = MDEventsTestHelper::makeFakeMDHistoWorkspace(2.0, 2, 5, 10.0, 3.0); + a->power(2.); + checkWorkspace(a, 4.0, 16*4*3./4.); + } + //-------------------------------------------------------------------------------------- void test_boolean_and() {