Commit f42b94b7d0c63c8902de6160eb6a3fd01da4bf57
1 parent
f9535133
Re-implement model parameter boundaries storing the factor boundaries (#3072)
The re-implementation should fix a precision problem that occurs on OpenSolaris.
Showing
6 changed files
with
324 additions
and
222 deletions
NEWS
| 1 | 1 | New Features and Important Changes in GammaLib 1.7.0 |
| 2 | 2 | |
| 3 | -26 June 2020 | |
| 3 | +29 June 2020 | |
| 4 | 4 | |
| 5 | 5 | |
| 6 | 6 | 1. Introduction |
| ... | ... | @@ -46,6 +46,10 @@ The following methods have been added: |
| 46 | 46 | - GSkyMap::extract() - Extract section of sky map |
| 47 | 47 | - GSkyRegionCircle::operator==() - Equality operator |
| 48 | 48 | - GSkyRegionCircle::operator!=() - Non-equality operator |
| 49 | +- GModelPar::has_factor_min() - Test if parameter has factor minimum | |
| 50 | +- GModelPar::has_factor_max() - Test if parameter has factor maximum | |
| 51 | +- GModelPar::remove_factor_min - Remove factor minimum | |
| 52 | +- GModelPar::remove_factor_max - Remove factor maximum | |
| 49 | 53 | - GModelSpatial::flux() - Extract flux from sky region |
| 50 | 54 | - GModelSpatialDiffuseCube::mc_cone() - Return MC simulation cone |
| 51 | 55 | - GModelSpatialDiffuseMap::mc_cone() - Return MC simulation cone | ... | ... |
include/GOptimizerPar.hpp
| ... | ... | @@ -79,6 +79,10 @@ |
| 79 | 79 | * value = m_factor_value * m_scale |
| 80 | 80 | * error = m_factor_error * m_scale |
| 81 | 81 | * gradient = m_factor_gradient * m_scale |
| 82 | + * min = m_factor_min * m_scale for m_scale > 0 | |
| 83 | + * m_factor_max * m_scale otherwise | |
| 84 | + * max = m_factor_max * m_scale for m_scale > 0 | |
| 85 | + * m_factor_min * m_scale otherwise | |
| 82 | 86 | * |
| 83 | 87 | * The @p factor and @p scale terms can be set and retrieved using the |
| 84 | 88 | * factor_value(), factor_error(), factor_gradient(), factor_min(), |
| ... | ... | @@ -115,8 +119,8 @@ public: |
| 115 | 119 | const double& factor_value(void) const; |
| 116 | 120 | const double& factor_error(void) const; |
| 117 | 121 | const double& factor_gradient(void) const; |
| 118 | - double factor_min(void) const; | |
| 119 | - double factor_max(void) const; | |
| 122 | + const double& factor_min(void) const; | |
| 123 | + const double& factor_max(void) const; | |
| 120 | 124 | const double& scale(void) const; |
| 121 | 125 | void factor_value(const double& value); |
| 122 | 126 | void factor_error(const double& error); |
| ... | ... | @@ -129,9 +133,13 @@ public: |
| 129 | 133 | // Boundary methods |
| 130 | 134 | bool has_min(void) const; |
| 131 | 135 | bool has_max(void) const; |
| 136 | + bool has_factor_min(void) const; | |
| 137 | + bool has_factor_max(void) const; | |
| 132 | 138 | bool has_range(void) const; |
| 133 | 139 | void remove_min(void); |
| 134 | 140 | void remove_max(void); |
| 141 | + void remove_factor_min(void); | |
| 142 | + void remove_factor_max(void); | |
| 135 | 143 | void remove_range(void); |
| 136 | 144 | |
| 137 | 145 | // Property methods |
| ... | ... | @@ -162,15 +170,15 @@ protected: |
| 162 | 170 | // Proteced data members |
| 163 | 171 | std::string m_name; //!< Parameter name |
| 164 | 172 | std::string m_unit; //!< Parameter unit |
| 165 | - double m_factor_value; //!< Parameter value factor | |
| 166 | - double m_factor_error; //!< Uncertainty in parameter value factor | |
| 167 | - double m_min_value; //!< Parameter minimum value | |
| 168 | - double m_max_value; //!< Parameter maximum value | |
| 169 | - mutable double m_factor_gradient; //!< Function gradient factor | |
| 173 | + double m_factor_value; //!< Parameter factor value | |
| 174 | + double m_factor_error; //!< Uncertainty in parameter factor value | |
| 175 | + double m_factor_min; //!< Parameter minimum factor value | |
| 176 | + double m_factor_max; //!< Parameter maximum factor value | |
| 177 | + mutable double m_factor_gradient; //!< Function factor gradient | |
| 170 | 178 | double m_scale; //!< Parameter scaling (true = factor * scale) |
| 171 | 179 | bool m_free; //!< Parameter is free |
| 172 | - bool m_has_min; //!< Parameter has minimum boundary | |
| 173 | - bool m_has_max; //!< Parameter has maximum boundary | |
| 180 | + bool m_has_factor_min; //!< Parameter has minimum factor boundary | |
| 181 | + bool m_has_factor_max; //!< Parameter has maximum factor boundary | |
| 174 | 182 | bool m_has_grad; //!< Parameter has analytic gradient |
| 175 | 183 | }; |
| 176 | 184 | |
| ... | ... | @@ -247,7 +255,7 @@ double GOptimizerPar::gradient(void) const |
| 247 | 255 | inline |
| 248 | 256 | double GOptimizerPar::min(void) const |
| 249 | 257 | { |
| 250 | - return (m_min_value); | |
| 258 | + return ((m_scale < 0) ? m_factor_max * m_scale : m_factor_min * m_scale); | |
| 251 | 259 | } |
| 252 | 260 | |
| 253 | 261 | |
| ... | ... | @@ -261,16 +269,16 @@ double GOptimizerPar::min(void) const |
| 261 | 269 | inline |
| 262 | 270 | double GOptimizerPar::max(void) const |
| 263 | 271 | { |
| 264 | - return (m_max_value); | |
| 272 | + return ((m_scale < 0) ? m_factor_min * m_scale : m_factor_max * m_scale); | |
| 265 | 273 | } |
| 266 | 274 | |
| 267 | 275 | |
| 268 | 276 | /***********************************************************************//** |
| 269 | - * @brief Return parameter value factor | |
| 277 | + * @brief Return parameter factor value | |
| 270 | 278 | * |
| 271 | - * @return Parameter value factor. | |
| 279 | + * @return Parameter factor value. | |
| 272 | 280 | * |
| 273 | - * Returns the parameter value factor. | |
| 281 | + * Returns the parameter factor value. | |
| 274 | 282 | ***************************************************************************/ |
| 275 | 283 | inline |
| 276 | 284 | const double& GOptimizerPar::factor_value(void) const |
| ... | ... | @@ -280,11 +288,11 @@ const double& GOptimizerPar::factor_value(void) const |
| 280 | 288 | |
| 281 | 289 | |
| 282 | 290 | /***********************************************************************//** |
| 283 | - * @brief Return parameter error factor | |
| 291 | + * @brief Return parameter factor error | |
| 284 | 292 | * |
| 285 | - * @return Parameter error factor. | |
| 293 | + * @return Parameter factor error. | |
| 286 | 294 | * |
| 287 | - * Returns the parameter error factor. | |
| 295 | + * Returns the parameter factor error. | |
| 288 | 296 | ***************************************************************************/ |
| 289 | 297 | inline |
| 290 | 298 | const double& GOptimizerPar::factor_error(void) const |
| ... | ... | @@ -294,11 +302,11 @@ const double& GOptimizerPar::factor_error(void) const |
| 294 | 302 | |
| 295 | 303 | |
| 296 | 304 | /***********************************************************************//** |
| 297 | - * @brief Return parameter gradient factor | |
| 305 | + * @brief Return parameter factor gradient | |
| 298 | 306 | * |
| 299 | - * @return Parameter gradient factor. | |
| 307 | + * @return Parameter factor gradient. | |
| 300 | 308 | * |
| 301 | - * Returns the parameter gradient factor. | |
| 309 | + * Returns the parameter factor gradient. | |
| 302 | 310 | ***************************************************************************/ |
| 303 | 311 | inline |
| 304 | 312 | const double& GOptimizerPar::factor_gradient(void) const |
| ... | ... | @@ -308,6 +316,34 @@ const double& GOptimizerPar::factor_gradient(void) const |
| 308 | 316 | |
| 309 | 317 | |
| 310 | 318 | /***********************************************************************//** |
| 319 | + * @brief Return parameter minimum factor boundary | |
| 320 | + * | |
| 321 | + * @return Minimum parameter factor boundary. | |
| 322 | + * | |
| 323 | + * Returns the minimum parameter factor boundary. | |
| 324 | + ***************************************************************************/ | |
| 325 | +inline | |
| 326 | +const double& GOptimizerPar::factor_min(void) const | |
| 327 | +{ | |
| 328 | + return (m_factor_min); | |
| 329 | +} | |
| 330 | + | |
| 331 | + | |
| 332 | +/***********************************************************************//** | |
| 333 | + * @brief Return parameter maximum factor boundary | |
| 334 | + * | |
| 335 | + * @return Maximum parameter factor boundary. | |
| 336 | + * | |
| 337 | + * Returns the maximum parameter factor boundary. | |
| 338 | + ***************************************************************************/ | |
| 339 | +inline | |
| 340 | +const double& GOptimizerPar::factor_max(void) const | |
| 341 | +{ | |
| 342 | + return (m_factor_max); | |
| 343 | +} | |
| 344 | + | |
| 345 | + | |
| 346 | +/***********************************************************************//** | |
| 311 | 347 | * @brief Set minimum and maximum parameter boundaries |
| 312 | 348 | * |
| 313 | 349 | * @param[in] min Parameter minimum. |
| ... | ... | @@ -339,11 +375,11 @@ const double& GOptimizerPar::scale(void) const |
| 339 | 375 | |
| 340 | 376 | |
| 341 | 377 | /***********************************************************************//** |
| 342 | - * @brief Set parameter error factor | |
| 378 | + * @brief Set parameter factor error | |
| 343 | 379 | * |
| 344 | - * @param[in] error Parameter error factor. | |
| 380 | + * @param[in] error Parameter factor error. | |
| 345 | 381 | * |
| 346 | - * Sets the parameter error factor. | |
| 382 | + * Sets the parameter factor error. | |
| 347 | 383 | ***************************************************************************/ |
| 348 | 384 | inline |
| 349 | 385 | void GOptimizerPar::factor_error(const double& error) |
| ... | ... | @@ -354,11 +390,11 @@ void GOptimizerPar::factor_error(const double& error) |
| 354 | 390 | |
| 355 | 391 | |
| 356 | 392 | /***********************************************************************//** |
| 357 | - * @brief Set parameter gradient factor | |
| 393 | + * @brief Set parameter factor gradient | |
| 358 | 394 | * |
| 359 | - * @param[in] gradient Parameter gradient factor. | |
| 395 | + * @param[in] gradient Parameter factor gradient. | |
| 360 | 396 | * |
| 361 | - * Sets the parameter gradient factor. | |
| 397 | + * Sets the parameter factor gradient. | |
| 362 | 398 | ***************************************************************************/ |
| 363 | 399 | inline |
| 364 | 400 | void GOptimizerPar::factor_gradient(const double& gradient) const |
| ... | ... | @@ -378,7 +414,7 @@ void GOptimizerPar::factor_gradient(const double& gradient) const |
| 378 | 414 | inline |
| 379 | 415 | bool GOptimizerPar::has_min(void) const |
| 380 | 416 | { |
| 381 | - return m_has_min; | |
| 417 | + return ((m_scale < 0) ? m_has_factor_max : m_has_factor_min); | |
| 382 | 418 | } |
| 383 | 419 | |
| 384 | 420 | |
| ... | ... | @@ -392,7 +428,35 @@ bool GOptimizerPar::has_min(void) const |
| 392 | 428 | inline |
| 393 | 429 | bool GOptimizerPar::has_max(void) const |
| 394 | 430 | { |
| 395 | - return m_has_max; | |
| 431 | + return ((m_scale < 0) ? m_has_factor_min : m_has_factor_max); | |
| 432 | +} | |
| 433 | + | |
| 434 | + | |
| 435 | +/***********************************************************************//** | |
| 436 | + * @brief Signal if parameter has minimum factor boundary | |
| 437 | + * | |
| 438 | + * @return True if parameter has minimum factor boundary, false otherwise. | |
| 439 | + * | |
| 440 | + * Signals if the parameter has a minimum factor boundary. | |
| 441 | + ***************************************************************************/ | |
| 442 | +inline | |
| 443 | +bool GOptimizerPar::has_factor_min(void) const | |
| 444 | +{ | |
| 445 | + return m_has_factor_min; | |
| 446 | +} | |
| 447 | + | |
| 448 | + | |
| 449 | +/***********************************************************************//** | |
| 450 | + * @brief Signal if parameter has maximum factor boundary | |
| 451 | + * | |
| 452 | + * @return True if parameter has maximum factor boundary, false otherwise. | |
| 453 | + * | |
| 454 | + * Signals if the parameter has a maximum factor boundary. | |
| 455 | + ***************************************************************************/ | |
| 456 | +inline | |
| 457 | +bool GOptimizerPar::has_factor_max(void) const | |
| 458 | +{ | |
| 459 | + return m_has_factor_max; | |
| 396 | 460 | } |
| 397 | 461 | |
| 398 | 462 | |
| ... | ... | @@ -407,7 +471,7 @@ bool GOptimizerPar::has_max(void) const |
| 407 | 471 | inline |
| 408 | 472 | bool GOptimizerPar::has_range(void) const |
| 409 | 473 | { |
| 410 | - return (m_has_min && m_has_max); | |
| 474 | + return (m_has_factor_min && m_has_factor_max); | |
| 411 | 475 | } |
| 412 | 476 | |
| 413 | 477 | |
| ... | ... | @@ -419,7 +483,14 @@ bool GOptimizerPar::has_range(void) const |
| 419 | 483 | inline |
| 420 | 484 | void GOptimizerPar::remove_min(void) |
| 421 | 485 | { |
| 422 | - m_has_min = false; | |
| 486 | + if (m_scale < 0) { | |
| 487 | + m_has_factor_max = false; | |
| 488 | + m_factor_max = 0.0; | |
| 489 | + } | |
| 490 | + else { | |
| 491 | + m_has_factor_min = false; | |
| 492 | + m_factor_min = 0.0; | |
| 493 | + } | |
| 423 | 494 | return; |
| 424 | 495 | } |
| 425 | 496 | |
| ... | ... | @@ -432,7 +503,42 @@ void GOptimizerPar::remove_min(void) |
| 432 | 503 | inline |
| 433 | 504 | void GOptimizerPar::remove_max(void) |
| 434 | 505 | { |
| 435 | - m_has_max = false; | |
| 506 | + if (m_scale < 0) { | |
| 507 | + m_has_factor_min = false; | |
| 508 | + m_factor_min = 0.0; | |
| 509 | + } | |
| 510 | + else { | |
| 511 | + m_has_factor_max = false; | |
| 512 | + m_factor_max = 0.0; | |
| 513 | + } | |
| 514 | + return; | |
| 515 | +} | |
| 516 | + | |
| 517 | + | |
| 518 | +/***********************************************************************//** | |
| 519 | + * @brief Removes minimum factor boundary | |
| 520 | + * | |
| 521 | + * Removes minimum factor boundary from the parameter. | |
| 522 | + ***************************************************************************/ | |
| 523 | +inline | |
| 524 | +void GOptimizerPar::remove_factor_min(void) | |
| 525 | +{ | |
| 526 | + m_has_factor_min = false; | |
| 527 | + m_factor_min = 0.0; | |
| 528 | + return; | |
| 529 | +} | |
| 530 | + | |
| 531 | + | |
| 532 | +/***********************************************************************//** | |
| 533 | + * @brief Removes maximum factor boundary | |
| 534 | + * | |
| 535 | + * Removes maximum factor boundary from the parameter. | |
| 536 | + ***************************************************************************/ | |
| 537 | +inline | |
| 538 | +void GOptimizerPar::remove_factor_max(void) | |
| 539 | +{ | |
| 540 | + m_has_factor_max = false; | |
| 541 | + m_factor_max = 0.0; | |
| 436 | 542 | return; |
| 437 | 543 | } |
| 438 | 544 | |
| ... | ... | @@ -445,8 +551,8 @@ void GOptimizerPar::remove_max(void) |
| 445 | 551 | inline |
| 446 | 552 | void GOptimizerPar::remove_range(void) |
| 447 | 553 | { |
| 448 | - m_has_min = false; | |
| 449 | - m_has_max = false; | |
| 554 | + m_has_factor_min = false; | |
| 555 | + m_has_factor_max = false; | |
| 450 | 556 | return; |
| 451 | 557 | } |
| 452 | 558 | ... | ... |
pyext/GOptimizerPar.i
| ... | ... | @@ -61,8 +61,8 @@ public: |
| 61 | 61 | const double& factor_value(void) const; |
| 62 | 62 | const double& factor_error(void) const; |
| 63 | 63 | const double& factor_gradient(void) const; |
| 64 | - double factor_min(void) const; | |
| 65 | - double factor_max(void) const; | |
| 64 | + const double& factor_min(void) const; | |
| 65 | + const double& factor_max(void) const; | |
| 66 | 66 | const double& scale(void) const; |
| 67 | 67 | void factor_value(const double& value); |
| 68 | 68 | void factor_error(const double& error); |
| ... | ... | @@ -75,9 +75,13 @@ public: |
| 75 | 75 | // Boundary methods |
| 76 | 76 | bool has_min(void) const; |
| 77 | 77 | bool has_max(void) const; |
| 78 | + bool has_factor_min(void) const; | |
| 79 | + bool has_factor_max(void) const; | |
| 78 | 80 | bool has_range(void) const; |
| 79 | 81 | void remove_min(void); |
| 80 | 82 | void remove_max(void); |
| 83 | + void remove_factor_min(void); | |
| 84 | + void remove_factor_max(void); | |
| 81 | 85 | void remove_range(void); |
| 82 | 86 | |
| 83 | 87 | // Property methods | ... | ... |
src/model/GModelPar.cpp
| ... | ... | @@ -270,7 +270,7 @@ void GModelPar::read(const GXmlElement& xml) |
| 270 | 270 | factor_min(gammalib::todouble(arg)); |
| 271 | 271 | } |
| 272 | 272 | else { |
| 273 | - remove_min(); | |
| 273 | + remove_factor_min(); | |
| 274 | 274 | } |
| 275 | 275 | |
| 276 | 276 | // Get max |
| ... | ... | @@ -279,7 +279,7 @@ void GModelPar::read(const GXmlElement& xml) |
| 279 | 279 | factor_max(gammalib::todouble(arg)); |
| 280 | 280 | } |
| 281 | 281 | else { |
| 282 | - remove_max(); | |
| 282 | + remove_factor_max(); | |
| 283 | 283 | } |
| 284 | 284 | |
| 285 | 285 | // Get free |
| ... | ... | @@ -293,7 +293,7 @@ void GModelPar::read(const GXmlElement& xml) |
| 293 | 293 | |
| 294 | 294 | // If there is a minimum and maximum, make sure that the maximum is |
| 295 | 295 | // not smaller than the minimum |
| 296 | - if (m_has_min && m_has_max) { | |
| 296 | + if (has_min() && has_max()) { | |
| 297 | 297 | if (min() > max()) { |
| 298 | 298 | std::string msg = "The model parameter \""+m_name+ |
| 299 | 299 | "\" in the XML document has a minimum boundary "+ |
| ... | ... | @@ -305,7 +305,7 @@ void GModelPar::read(const GXmlElement& xml) |
| 305 | 305 | } |
| 306 | 306 | |
| 307 | 307 | // If there is a minimum, make sure that the value is not below it |
| 308 | - if (m_has_min && value() < min()) { | |
| 308 | + if (has_min() && (value() < min())) { | |
| 309 | 309 | std::string msg = "The model parameter \""+m_name+ |
| 310 | 310 | "\" in the XML document has a value "+ |
| 311 | 311 | gammalib::str(value())+ |
| ... | ... | @@ -315,7 +315,7 @@ void GModelPar::read(const GXmlElement& xml) |
| 315 | 315 | } |
| 316 | 316 | |
| 317 | 317 | // If there is a maximum, make sure that the value is not above it |
| 318 | - if (m_has_max && value() > max()) { | |
| 318 | + if (has_max() && (value() > max())) { | |
| 319 | 319 | std::string msg = "The model parameter \""+m_name+ |
| 320 | 320 | "\" in the XML document has a value "+ |
| 321 | 321 | gammalib::str(value())+ |
| ... | ... | @@ -365,12 +365,12 @@ void GModelPar::write(GXmlElement& xml) const |
| 365 | 365 | xml.attribute("scale", gammalib::str(m_scale)); |
| 366 | 366 | |
| 367 | 367 | // Set minimum |
| 368 | - if (has_min()) { | |
| 368 | + if (has_factor_min()) { | |
| 369 | 369 | xml.attribute("min", gammalib::str(factor_min())); |
| 370 | 370 | } |
| 371 | 371 | |
| 372 | 372 | // Set maximum |
| 373 | - if (has_max()) { | |
| 373 | + if (has_factor_max()) { | |
| 374 | 374 | xml.attribute("max", gammalib::str(factor_max())); |
| 375 | 375 | } |
| 376 | 376 | ... | ... |
src/opt/GOptimizerPar.cpp
| ... | ... | @@ -321,11 +321,17 @@ void GOptimizerPar::min(const double& min) |
| 321 | 321 | throw GException::invalid_argument(G_MIN, msg); |
| 322 | 322 | } |
| 323 | 323 | |
| 324 | - // Store minimum | |
| 325 | - m_min_value = min; | |
| 324 | + // If scale is negative then set maximum factor boundary | |
| 325 | + if (m_scale < 0) { | |
| 326 | + m_factor_max = min / m_scale; | |
| 327 | + m_has_factor_max = true; | |
| 328 | + } | |
| 326 | 329 | |
| 327 | - // Flag that minimum was set | |
| 328 | - m_has_min = true; | |
| 330 | + // ... otherwise set minimum factor boundary | |
| 331 | + else { | |
| 332 | + m_factor_min = min / m_scale; | |
| 333 | + m_has_factor_min = true; | |
| 334 | + } | |
| 329 | 335 | |
| 330 | 336 | // Return |
| 331 | 337 | return; |
| ... | ... | @@ -352,11 +358,17 @@ void GOptimizerPar::max(const double& max) |
| 352 | 358 | throw GException::invalid_argument(G_MAX, msg); |
| 353 | 359 | } |
| 354 | 360 | |
| 355 | - // Store maximum | |
| 356 | - m_max_value = max; | |
| 361 | + // If scale is negative then set minimum factor boundary | |
| 362 | + if (m_scale < 0) { | |
| 363 | + m_factor_min = max / m_scale; | |
| 364 | + m_has_factor_min = true; | |
| 365 | + } | |
| 357 | 366 | |
| 358 | - // Flag that maximum was set | |
| 359 | - m_has_max = true; | |
| 367 | + // ... otherwise set maximum factor boundary | |
| 368 | + else { | |
| 369 | + m_factor_max = max / m_scale; | |
| 370 | + m_has_factor_max = true; | |
| 371 | + } | |
| 360 | 372 | |
| 361 | 373 | // Return |
| 362 | 374 | return; |
| ... | ... | @@ -364,54 +376,6 @@ void GOptimizerPar::max(const double& max) |
| 364 | 376 | |
| 365 | 377 | |
| 366 | 378 | /***********************************************************************//** |
| 367 | - * @brief Return minimum parameter boundary factor | |
| 368 | - * | |
| 369 | - * @return Minimum parameter boundary factor. | |
| 370 | - * | |
| 371 | - * Returns the minimum parameter boundary factor. | |
| 372 | - * | |
| 373 | - * For positive scale factors the minimum boundary factor is the minimum | |
| 374 | - * boundary value divided by the scale factor, for negative scale factors | |
| 375 | - * it is the maximum boundary value divided by the scale factor. | |
| 376 | - ***************************************************************************/ | |
| 377 | -double GOptimizerPar::factor_min(void) const | |
| 378 | -{ | |
| 379 | - // Determine the minimum boundary factor. If the scale factor is negative | |
| 380 | - // the maximum boundary value divided by the scale factor is returned, | |
| 381 | - // otherwise the minimum boundary value divided by the scale factor is | |
| 382 | - // returned | |
| 383 | - double min = (m_scale < 0) ? m_max_value/m_scale : m_min_value/m_scale; | |
| 384 | - | |
| 385 | - // Return minimum boundary factor | |
| 386 | - return (min); | |
| 387 | -} | |
| 388 | - | |
| 389 | - | |
| 390 | -/***********************************************************************//** | |
| 391 | - * @brief Return maximum parameter boundary factor | |
| 392 | - * | |
| 393 | - * @return Maximum parameter boundary factor. | |
| 394 | - * | |
| 395 | - * Returns the maximum parameter boundary factor. | |
| 396 | - * | |
| 397 | - * For positive scale factors the maximum boundary factor is the maximum | |
| 398 | - * boundary value divided by the scale factor, for negative scale factors | |
| 399 | - * it is the minimum boundary value divided by the scale factor. | |
| 400 | - ***************************************************************************/ | |
| 401 | -double GOptimizerPar::factor_max(void) const | |
| 402 | -{ | |
| 403 | - // Determine the minimum boundary factor. If the scale factor is negative | |
| 404 | - // the minimum boundary value divided by the scale factor is returned, | |
| 405 | - // otherwise the maximum boundary value divided by the scale factor is | |
| 406 | - // returned | |
| 407 | - double max = (m_scale < 0) ? m_min_value/m_scale : m_max_value/m_scale; | |
| 408 | - | |
| 409 | - // Return maximum boundary factor | |
| 410 | - return (max); | |
| 411 | -} | |
| 412 | - | |
| 413 | - | |
| 414 | -/***********************************************************************//** | |
| 415 | 379 | * @brief Set parameter value factor |
| 416 | 380 | * |
| 417 | 381 | * @param[in] value Value factor. |
| ... | ... | @@ -426,7 +390,7 @@ void GOptimizerPar::factor_value(const double& value) |
| 426 | 390 | { |
| 427 | 391 | // If there is a minimum boundary and if value is below this boundary |
| 428 | 392 | // then throw an exception |
| 429 | - if (m_has_min && (value*m_scale < min())) { | |
| 393 | + if (has_factor_min() && (value < factor_min())) { | |
| 430 | 394 | std::string msg = "Specified value factor "+gammalib::str(value)+ |
| 431 | 395 | " is smaller than the minimum boundary "+ |
| 432 | 396 | gammalib::str(factor_min())+"."; |
| ... | ... | @@ -435,7 +399,7 @@ void GOptimizerPar::factor_value(const double& value) |
| 435 | 399 | |
| 436 | 400 | // If there is a maximum boundary and if value is above this boundary |
| 437 | 401 | // then throw an exception |
| 438 | - if (m_has_max && (value*m_scale > max())) { | |
| 402 | + if (has_factor_max() && (value > factor_max())) { | |
| 439 | 403 | std::string msg = "Specified value factor "+gammalib::str(value)+ |
| 440 | 404 | " is larger than the maximum boundary "+ |
| 441 | 405 | gammalib::str(factor_max())+"."; |
| ... | ... | @@ -472,16 +436,9 @@ void GOptimizerPar::factor_min(const double& min) |
| 472 | 436 | throw GException::invalid_argument(G_FACTOR_VALUE, msg); |
| 473 | 437 | } |
| 474 | 438 | |
| 475 | - // For negative scale factors, the minimum is stored as the maximum | |
| 476 | - // value after multiplication with the scale factor. For positive | |
| 477 | - // scale factors it is stored as minimum value after multiplication | |
| 478 | - // with the scale factor. | |
| 479 | - if (m_scale < 0) { | |
| 480 | - this->max(min * m_scale); | |
| 481 | - } | |
| 482 | - else { | |
| 483 | - this->min(min * m_scale); | |
| 484 | - } | |
| 439 | + // Set value and flag | |
| 440 | + m_factor_min = min; | |
| 441 | + m_has_factor_min = true; | |
| 485 | 442 | |
| 486 | 443 | // Return |
| 487 | 444 | return; |
| ... | ... | @@ -510,16 +467,9 @@ void GOptimizerPar::factor_max(const double& max) |
| 510 | 467 | throw GException::invalid_argument(G_FACTOR_VALUE, msg); |
| 511 | 468 | } |
| 512 | 469 | |
| 513 | - // For negative scale factors, the maximum is stored as the minimum | |
| 514 | - // value after multiplication with the scale factor. For positive | |
| 515 | - // scale factors it is stored as maximum value after multiplication | |
| 516 | - // with the scale factor. | |
| 517 | - if (m_scale < 0) { | |
| 518 | - this->min(max * m_scale); | |
| 519 | - } | |
| 520 | - else { | |
| 521 | - this->max(max * m_scale); | |
| 522 | - } | |
| 470 | + // Set value and flag | |
| 471 | + m_factor_max = max; | |
| 472 | + m_has_factor_max = true; | |
| 523 | 473 | |
| 524 | 474 | // Return |
| 525 | 475 | return; |
| ... | ... | @@ -578,7 +528,25 @@ void GOptimizerPar::scale(const double& scale) |
| 578 | 528 | // Set values, error, gradient, min and max |
| 579 | 529 | m_factor_value *= rescale; |
| 580 | 530 | m_factor_error *= rescale; |
| 581 | - m_factor_gradient *= rescale; | |
| 531 | + m_factor_gradient /= rescale; | |
| 532 | + m_factor_min *= rescale; | |
| 533 | + m_factor_max *= rescale; | |
| 534 | + | |
| 535 | + // If re-scaling changes sign then swap minimum and maximum boundary | |
| 536 | + if (rescale < 0.0) { | |
| 537 | + | |
| 538 | + // Store maximum in swap space | |
| 539 | + bool bswap = m_has_factor_max; | |
| 540 | + double dswap = m_factor_max; | |
| 541 | + | |
| 542 | + // Set maximum | |
| 543 | + m_has_factor_max = m_has_factor_min; | |
| 544 | + m_factor_max = m_factor_min; | |
| 545 | + | |
| 546 | + // Set minimum | |
| 547 | + m_has_factor_min = bswap; | |
| 548 | + m_factor_min = dswap; | |
| 549 | + } | |
| 582 | 550 | |
| 583 | 551 | // Return |
| 584 | 552 | return; |
| ... | ... | @@ -603,16 +571,8 @@ void GOptimizerPar::autoscale(void) |
| 603 | 571 | // Continue only if the value is non-zero |
| 604 | 572 | if (value != 0.0) { |
| 605 | 573 | |
| 606 | - // Get the renormalization factor | |
| 607 | - double renormalization = m_scale / value; | |
| 608 | - | |
| 609 | - // Set the new scale factor to the actual value | |
| 610 | - m_scale = value; | |
| 611 | - | |
| 612 | - // Renormalize values, error, gradient, min and max | |
| 613 | - m_factor_value *= renormalization; | |
| 614 | - m_factor_error *= renormalization; | |
| 615 | - m_factor_gradient /= renormalization; | |
| 574 | + // Re-scale to value | |
| 575 | + this->scale(value); | |
| 616 | 576 | |
| 617 | 577 | } // endif: value was non-zero |
| 618 | 578 | |
| ... | ... | @@ -647,13 +607,13 @@ std::string GOptimizerPar::print(const GChatter& chatter) const |
| 647 | 607 | } |
| 648 | 608 | |
| 649 | 609 | // Append parameter limites if they exist |
| 650 | - if (m_has_min && m_has_max) { | |
| 610 | + if (has_min() && has_max()) { | |
| 651 | 611 | result.append(" ["+gammalib::str(min()) + ","+gammalib::str(max())+"]"); |
| 652 | 612 | } |
| 653 | - else if (m_has_min) { | |
| 613 | + else if (has_min()) { | |
| 654 | 614 | result.append(" ["+gammalib::str(min()) + ",infty["); |
| 655 | 615 | } |
| 656 | - else if (m_has_max) { | |
| 616 | + else if (has_max()) { | |
| 657 | 617 | result.append(" ]-infty,"+gammalib::str(max())+"]"); |
| 658 | 618 | } |
| 659 | 619 | |
| ... | ... | @@ -703,12 +663,12 @@ void GOptimizerPar::init_members(void) |
| 703 | 663 | m_factor_value = 0.0; |
| 704 | 664 | m_factor_error = 0.0; |
| 705 | 665 | m_factor_gradient = 0.0; |
| 706 | - m_min_value = 0.0; | |
| 707 | - m_max_value = 0.0; | |
| 666 | + m_factor_min = 0.0; | |
| 667 | + m_factor_max = 0.0; | |
| 708 | 668 | m_scale = 1.0; |
| 709 | 669 | m_free = true; |
| 710 | - m_has_min = false; | |
| 711 | - m_has_max = false; | |
| 670 | + m_has_factor_min = false; | |
| 671 | + m_has_factor_max = false; | |
| 712 | 672 | m_has_grad = false; |
| 713 | 673 | |
| 714 | 674 | // Return |
| ... | ... | @@ -729,12 +689,12 @@ void GOptimizerPar::copy_members(const GOptimizerPar& par) |
| 729 | 689 | m_factor_value = par.m_factor_value; |
| 730 | 690 | m_factor_error = par.m_factor_error; |
| 731 | 691 | m_factor_gradient = par.m_factor_gradient; |
| 732 | - m_min_value = par.m_min_value; | |
| 733 | - m_max_value = par.m_max_value; | |
| 692 | + m_factor_min = par.m_factor_min; | |
| 693 | + m_factor_max = par.m_factor_max; | |
| 734 | 694 | m_scale = par.m_scale; |
| 735 | 695 | m_free = par.m_free; |
| 736 | - m_has_min = par.m_has_min; | |
| 737 | - m_has_max = par.m_has_max; | |
| 696 | + m_has_factor_min = par.m_has_factor_min; | |
| 697 | + m_has_factor_max = par.m_has_factor_max; | |
| 738 | 698 | m_has_grad = par.m_has_grad; |
| 739 | 699 | |
| 740 | 700 | // Return | ... | ... |
test/test_GModel.cpp
| ... | ... | @@ -252,37 +252,31 @@ TestGModel* TestGModel::clone(void) const |
| 252 | 252 | void TestGModel::test_model_par(void) |
| 253 | 253 | { |
| 254 | 254 | // Test void constructor |
| 255 | - test_try("Test void constructor"); | |
| 256 | - try { | |
| 257 | - GModelPar par; | |
| 258 | - test_try_success(); | |
| 259 | - } | |
| 260 | - catch (std::exception &e) { | |
| 261 | - test_try_failure(e); | |
| 262 | - } | |
| 255 | + GModelPar par1; | |
| 256 | + test_value(par1.value(), 0.0); | |
| 257 | + test_value(par1.error(), 0.0); | |
| 258 | + test_value(par1.gradient(), 0.0); | |
| 259 | + test_value(par1.min(), 0.0); | |
| 260 | + test_value(par1.max(), 0.0); | |
| 261 | + test_assert(!par1.has_min(), "Parameter shall have no minimum."); | |
| 262 | + test_assert(!par1.has_max(), "Parameter shall have no maximum."); | |
| 263 | + test_assert(!par1.has_range(), "Parameter shall have no range."); | |
| 263 | 264 | |
| 264 | - // Test parameter constructor (value version) | |
| 265 | - test_try("Test parameter constructor"); | |
| 266 | - try { | |
| 267 | - GModelPar par("Test parameter", 47.0); | |
| 268 | - test_value(par.value(), 47.0); | |
| 269 | - test_value(par.error(), 0.0); | |
| 270 | - test_value(par.gradient(), 0.0); | |
| 271 | - test_value(par.min(), 0.0); | |
| 272 | - test_value(par.max(), 0.0); | |
| 273 | - test_assert(!par.has_min(), "Parameter shall have no minimum."); | |
| 274 | - test_assert(!par.has_max(), "Parameter shall have no maximum."); | |
| 275 | - test_assert(!par.has_range(), "Parameter shall have no range."); | |
| 276 | - test_try_success(); | |
| 277 | - } | |
| 278 | - catch (std::exception &e) { | |
| 279 | - test_try_failure(e); | |
| 280 | - } | |
| 265 | + // Test parameter value constructor (value version) | |
| 266 | + GModelPar par2("Test parameter", 47.0); | |
| 267 | + test_value(par2.value(), 47.0); | |
| 268 | + test_value(par2.error(), 0.0); | |
| 269 | + test_value(par2.gradient(), 0.0); | |
| 270 | + test_value(par2.min(), 0.0); | |
| 271 | + test_value(par2.max(), 0.0); | |
| 272 | + test_assert(!par2.has_min(), "Parameter shall have no minimum."); | |
| 273 | + test_assert(!par2.has_max(), "Parameter shall have no maximum."); | |
| 274 | + test_assert(!par2.has_range(), "Parameter shall have no range."); | |
| 281 | 275 | |
| 282 | 276 | // Test invalid parameter constructor (factor & scale version) |
| 283 | 277 | test_try("Test invalid parameter constructor"); |
| 284 | 278 | try { |
| 285 | - GModelPar par("Test parameter", 1.0, 0.0); | |
| 279 | + GModelPar par3("Test parameter", 1.0, 0.0); | |
| 286 | 280 | test_try_failure("Parameter constructor with zero scale factor" |
| 287 | 281 | " shall throw an exception."); |
| 288 | 282 | } |
| ... | ... | @@ -294,57 +288,50 @@ void TestGModel::test_model_par(void) |
| 294 | 288 | } |
| 295 | 289 | |
| 296 | 290 | // Test valid parameter constructor (factor & scale version) |
| 297 | - test_try("Test invalid parameter constructor"); | |
| 298 | - try { | |
| 299 | - GModelPar par("Test parameter", 47.01, 2.0); | |
| 300 | - par.factor_error(2.003); | |
| 301 | - par.factor_gradient(51.0); | |
| 302 | - par.unit("MeV"); | |
| 303 | - par.free(); | |
| 304 | - test_value(par.value(), 94.02); | |
| 305 | - test_value(par.error(), 4.006); | |
| 306 | - test_value(par.gradient(), 25.5); | |
| 307 | - test_value(par.factor_value(), 47.01); | |
| 308 | - test_value(par.factor_error(), 2.003); | |
| 309 | - test_value(par.factor_gradient(), 51.0); | |
| 310 | - test_value(par.scale(), 2.0); | |
| 311 | - test_assert(par.name() == "Test parameter", "Parameter name"); | |
| 312 | - test_assert(par.unit() == "MeV", "Parameter unit"); | |
| 313 | - test_assert(par.is_free(), "Parameter freezing"); | |
| 314 | - test_assert(par.print() == | |
| 315 | - " Test parameter ...........: 94.02 +/- 4.006 MeV" | |
| 316 | - " (free,scale=2)", "Parameter printing"); | |
| 317 | - GModelPar par2("Another test parameter", 3.14, 3.0); | |
| 318 | - test_value(par2.value(), 9.42); | |
| 319 | - test_value(par2.factor_value(), 3.14); | |
| 320 | - test_value(par2.scale(), 3.0); | |
| 321 | - test_assert(par2.name() == "Another test parameter", "Parameter name"); | |
| 322 | - test_assert(par2.is_free(), "Parameter freezing"); | |
| 323 | - test_assert(par2.print() == | |
| 324 | - " Another test parameter ...: 9.42 +/- 0 (free,scale=3)", | |
| 325 | - "Parameter printing"); | |
| 326 | - test_try_success(); | |
| 327 | - } | |
| 328 | - catch (std::exception &e) { | |
| 329 | - test_try_failure(e); | |
| 330 | - } | |
| 291 | + GModelPar par4("Test parameter", 47.01, 2.0); | |
| 292 | + par4.factor_error(2.003); | |
| 293 | + par4.factor_gradient(51.0); | |
| 294 | + par4.unit("MeV"); | |
| 295 | + par4.free(); | |
| 296 | + test_value(par4.value(), 94.02); | |
| 297 | + test_value(par4.error(), 4.006); | |
| 298 | + test_value(par4.gradient(), 25.5); | |
| 299 | + test_value(par4.factor_value(), 47.01); | |
| 300 | + test_value(par4.factor_error(), 2.003); | |
| 301 | + test_value(par4.factor_gradient(), 51.0); | |
| 302 | + test_value(par4.scale(), 2.0); | |
| 303 | + test_assert(par4.name() == "Test parameter", "Parameter name"); | |
| 304 | + test_assert(par4.unit() == "MeV", "Parameter unit"); | |
| 305 | + test_assert(par4.is_free(), "Parameter freezing"); | |
| 306 | + test_assert(par4.print() == | |
| 307 | + " Test parameter ...........: 94.02 +/- 4.006 MeV" | |
| 308 | + " (free,scale=2)", "Parameter printing"); | |
| 309 | + GModelPar par5("Another test parameter", 3.14, 3.0); | |
| 310 | + test_value(par5.value(), 9.42); | |
| 311 | + test_value(par5.factor_value(), 3.14); | |
| 312 | + test_value(par5.scale(), 3.0); | |
| 313 | + test_assert(par5.name() == "Another test parameter", "Parameter name"); | |
| 314 | + test_assert(par5.is_free(), "Parameter freezing"); | |
| 315 | + test_assert(par5.print() == | |
| 316 | + " Another test parameter ...: 9.42 +/- 0 (free,scale=3)", | |
| 317 | + "Parameter printing"); | |
| 331 | 318 | |
| 332 | 319 | // Test boundary handling 1 |
| 333 | 320 | test_try("Test boundary handling (1/4)"); |
| 334 | 321 | try { |
| 335 | - GModelPar par("Test boundary", 1.0); | |
| 336 | - test_assert(!par.has_min(), "Parameter shall have no minimum."); | |
| 337 | - test_assert(!par.has_max(), "Parameter shall have no maximum."); | |
| 338 | - test_assert(!par.has_range(), "Parameter shall have no range."); | |
| 339 | - par.min(0.5); | |
| 340 | - test_assert(par.has_min(), "Parameter shall have minimum."); | |
| 341 | - test_assert(!par.has_max(), "Parameter shall have no maximum."); | |
| 342 | - test_assert(!par.has_range(), "Parameter shall have no range."); | |
| 343 | - par.max(2.0); | |
| 344 | - test_assert(par.has_min(), "Parameter shall have minimum."); | |
| 345 | - test_assert(par.has_max(), "Parameter shall have maximum."); | |
| 346 | - test_assert(par.has_range(), "Parameter shall have range."); | |
| 347 | - par.value(5.0); | |
| 322 | + GModelPar par6("Test boundary", 1.0); | |
| 323 | + test_assert(!par6.has_min(), "Parameter shall have no minimum."); | |
| 324 | + test_assert(!par6.has_max(), "Parameter shall have no maximum."); | |
| 325 | + test_assert(!par6.has_range(), "Parameter shall have no range."); | |
| 326 | + par6.min(0.5); | |
| 327 | + test_assert(par6.has_min(), "Parameter shall have minimum."); | |
| 328 | + test_assert(!par6.has_max(), "Parameter shall have no maximum."); | |
| 329 | + test_assert(!par6.has_range(), "Parameter shall have no range."); | |
| 330 | + par6.max(2.0); | |
| 331 | + test_assert(par6.has_min(), "Parameter shall have minimum."); | |
| 332 | + test_assert(par6.has_max(), "Parameter shall have maximum."); | |
| 333 | + test_assert(par6.has_range(), "Parameter shall have range."); | |
| 334 | + par6.value(5.0); | |
| 348 | 335 | test_try_failure("Setting a parameter outside boundaries shall" |
| 349 | 336 | " generate an exception."); |
| 350 | 337 | } |
| ... | ... | @@ -358,8 +345,8 @@ void TestGModel::test_model_par(void) |
| 358 | 345 | // Test boundary handling 2/4 |
| 359 | 346 | test_try("Test boundary handling (2/4)"); |
| 360 | 347 | try { |
| 361 | - GModelPar par("Test boundary", 1.0); | |
| 362 | - par.min(2.0); | |
| 348 | + GModelPar par7("Test boundary", 1.0); | |
| 349 | + par7.min(2.0); | |
| 363 | 350 | test_try_failure("Setting the minimum boundary that is larger" |
| 364 | 351 | " than the parameter value shall generate an" |
| 365 | 352 | " exception."); |
| ... | ... | @@ -374,8 +361,8 @@ void TestGModel::test_model_par(void) |
| 374 | 361 | // Test boundary handling 3/4 |
| 375 | 362 | test_try("Test boundary handling (3/4)"); |
| 376 | 363 | try { |
| 377 | - GModelPar par("Test boundary", 1.0); | |
| 378 | - par.max(0.5); | |
| 364 | + GModelPar par8("Test boundary", 1.0); | |
| 365 | + par8.max(0.5); | |
| 379 | 366 | test_try_failure("Setting the maximum boundary that is smaller" |
| 380 | 367 | " than the parameter value shall generate an" |
| 381 | 368 | " exception."); |
| ... | ... | @@ -390,8 +377,8 @@ void TestGModel::test_model_par(void) |
| 390 | 377 | // Test boundary handling 4/4 |
| 391 | 378 | test_try("Test boundary handling (4/4)"); |
| 392 | 379 | try { |
| 393 | - GModelPar par("Test boundary", 1.0); | |
| 394 | - par.range(2.0, 0.5); | |
| 380 | + GModelPar par9("Test boundary", 1.0); | |
| 381 | + par9.range(2.0, 0.5); | |
| 395 | 382 | test_try_failure("Setting the minimum boundary that is larger" |
| 396 | 383 | " than the maximum boundary shall generate an" |
| 397 | 384 | " exception."); |
| ... | ... | @@ -403,6 +390,47 @@ void TestGModel::test_model_par(void) |
| 403 | 390 | test_try_failure(e); |
| 404 | 391 | } |
| 405 | 392 | |
| 393 | + // Test boundary handling | |
| 394 | + GModelPar par10("Test parameter", 3.0, 3.0); | |
| 395 | + par10.min(1.0); | |
| 396 | + test_assert(par10.has_min(), "1Parameter shall have minimum."); | |
| 397 | + test_assert(!par10.has_max(), "1Parameter shall have no maximum."); | |
| 398 | + test_assert(par10.has_factor_min(), "1Parameter shall have factor minimum."); | |
| 399 | + test_assert(!par10.has_factor_max(), "1Parameter shall have no factor maximum."); | |
| 400 | + test_assert(!par10.has_range(), "1Parameter shall have no range."); | |
| 401 | + test_value(par10.min(), 1.0); | |
| 402 | + test_value(par10.factor_min(), 1.0/3.0); | |
| 403 | + test_value(par10.factor_max(), 0.0); | |
| 404 | + par10.scale(-1.0); | |
| 405 | + test_assert(par10.has_min(), "2Parameter shall have minimum."); | |
| 406 | + test_assert(!par10.has_max(), "2Parameter shall have no maximum."); | |
| 407 | + test_assert(!par10.has_factor_min(), "2Parameter shall have factor minimum."); | |
| 408 | + test_assert(par10.has_factor_max(), "2Parameter shall have no factor maximum."); | |
| 409 | + test_assert(!par10.has_range(), "2Parameter shall have no range."); | |
| 410 | + test_value(par10.min(), 1.0); | |
| 411 | + test_value(par10.factor_min(), 0.0); | |
| 412 | + test_value(par10.factor_max(), -1.0); | |
| 413 | + // | |
| 414 | + par10.remove_min(); | |
| 415 | + par10.max(10.0); | |
| 416 | + test_assert(!par10.has_min(), "3Parameter shall have no minimum."); | |
| 417 | + test_assert(par10.has_max(), "3Parameter shall have maximum."); | |
| 418 | + test_assert(par10.has_factor_min(), "3Parameter shall have no factor minimum."); | |
| 419 | + test_assert(!par10.has_factor_max(), "3Parameter shall have factor maximum."); | |
| 420 | + test_assert(!par10.has_range(), "3Parameter shall have no range."); | |
| 421 | + test_value(par10.max(), 10.0); | |
| 422 | + test_value(par10.factor_min(), -10.0); | |
| 423 | + test_value(par10.factor_max(), 0.0); | |
| 424 | + par10.scale(3.0); | |
| 425 | + test_assert(!par10.has_min(), "4Parameter shall have no minimum."); | |
| 426 | + test_assert(par10.has_max(), "4Parameter shall have maximum."); | |
| 427 | + test_assert(!par10.has_factor_min(), "4Parameter shall have no factor minimum."); | |
| 428 | + test_assert(par10.has_factor_max(), "4Parameter shall have factor maximum."); | |
| 429 | + test_assert(!par10.has_range(), "4Parameter shall have no range."); | |
| 430 | + test_value(par10.max(), 10.0); | |
| 431 | + test_value(par10.factor_min(), 0.0); | |
| 432 | + test_value(par10.factor_max(), 10.0/3.0); | |
| 433 | + | |
| 406 | 434 | // Test autoscaling |
| 407 | 435 | GModelPar par("Test parameter", 3.0); |
| 408 | 436 | par.scale(1.0); | ... | ... |