diff --git a/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.cxx b/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.cxx index 19cc9effc912d939037fb0b97453f570ae58cdd8..c51529a7eaaf7b5d89cbc71ed3944a82be0d3455 100755 --- a/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.cxx +++ b/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration */ /**************************************************************************************** @@ -33,7 +33,7 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M return modes; } - std::vector< Mixture > mixture[5]; + std::array<std::vector< Mixture >,5> mixture; fillMixture( multiComponentState, mixture ); @@ -49,7 +49,6 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M double para_start = para_startMap.begin()->second; modes[i] = findMode( para_start, i, mixture, log ) ; - // Calculate the FWHM and return this back so that it can be used to correct the covariance matrix if( para_start != modes[i] ){ // mode calculation was successful now calulate FWHM @@ -63,7 +62,7 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M double upperbound =modes[i] + 1.5 * currentWidth; while(true){ if( pdf( upperbound, i, mixture ) > pdfVal*0.5 ){ - upperbound += currentWidth; + upperbound += currentWidth; } else { break; } @@ -77,7 +76,7 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M double lowerbound =modes[i] - 1.5 * currentWidth; while(true){ if( pdf( lowerbound, i, mixture ) > pdfVal*0.5 ){ - lowerbound -= currentWidth; + lowerbound -= currentWidth; } else { break; } @@ -92,7 +91,6 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M double FWHM = highX - lowX; if( ATH_UNLIKELY( log.level() <= MSG::DEBUG ) ) log << "PDFval, high val, low val [ " << pdfVal << ", " << pdf( highX, i, mixture ) << ", " << pdf( lowX, i, mixture) << "]" << endmsg; - if( FWHM <= 0 ) { if( ATH_UNLIKELY( log.level() <= MSG::DEBUG ) ) log << i << " Width is neagtive? " << highX << " " << lowX << " " << modes[i] << endmsg; @@ -117,15 +115,13 @@ Amg::VectorX Trk::MultiComponentStateModeCalculator::calculateMode( const Trk::M } } } - - para_startMap.clear(); } return modes; } -void Trk::MultiComponentStateModeCalculator::fillMixture( const Trk::MultiComponentState& multiComponentState, std::vector< Mixture >* mixture ) +void Trk::MultiComponentStateModeCalculator::fillMixture( const Trk::MultiComponentState& multiComponentState, std::array<std::vector< Mixture >,5>& mixture ) { for (int i=0; i<5; i++){ @@ -160,16 +156,14 @@ void Trk::MultiComponentStateModeCalculator::fillMixture( const Trk::MultiCompon mean -= 2 * M_PI; } } - Mixture mix(weight, mean, sigma ); - mixture[i].push_back( mix ); } } return; } -double Trk::MultiComponentStateModeCalculator::findMode( double xStart, int i, std::vector< Mixture >* mixture, MsgStream &log ) +double Trk::MultiComponentStateModeCalculator::findMode( double xStart, int i, std::array<std::vector< Mixture >,5>& mixture, MsgStream &log ) { int iteration(0); @@ -181,7 +175,6 @@ double Trk::MultiComponentStateModeCalculator::findMode( double xStart, int i, s while( iteration < 20 && tolerance > 1.e-8 ){ double previousMode(mode); - double d2pdfVal = d2pdf( previousMode, i, mixture); if (d2pdfVal != 0.0) { @@ -218,7 +211,7 @@ double Trk::MultiComponentStateModeCalculator::findMode( double xStart, int i, s } -double Trk::MultiComponentStateModeCalculator::findModeGlobal(double mean,int i, std::vector< Mixture >* mixture) +double Trk::MultiComponentStateModeCalculator::findModeGlobal(double mean,int i, std::array<std::vector< Mixture >,5>& mixture) { double start(-1); @@ -246,7 +239,7 @@ double Trk::MultiComponentStateModeCalculator::findModeGlobal(double mean,int i, return mode; } -double Trk::MultiComponentStateModeCalculator::pdf( double x, int i, std::vector< Mixture >* mixture ) +double Trk::MultiComponentStateModeCalculator::pdf( double x, int i, std::array<std::vector< Mixture >,5>& mixture ) { double pdf(0.); @@ -260,7 +253,7 @@ double Trk::MultiComponentStateModeCalculator::pdf( double x, int i, std::vector } -double Trk::MultiComponentStateModeCalculator::d1pdf( double x, int i,std::vector< Mixture >* mixture ) +double Trk::MultiComponentStateModeCalculator::d1pdf( double x, int i,std::array<std::vector< Mixture >,5>& mixture ) { double result(0.); @@ -279,7 +272,7 @@ double Trk::MultiComponentStateModeCalculator::d1pdf( double x, int i,std::vecto } -double Trk::MultiComponentStateModeCalculator::d2pdf( double x, int i, std::vector< Mixture >* mixture ) +double Trk::MultiComponentStateModeCalculator::d2pdf( double x, int i, std::array<std::vector< Mixture >,5>& mixture ) { double result(0.); @@ -313,7 +306,7 @@ double Trk::MultiComponentStateModeCalculator::gaus( double x, double mean, doub } -double Trk::MultiComponentStateModeCalculator::width( int i, std::vector< Mixture >* mixture ) +double Trk::MultiComponentStateModeCalculator::width( int i, std::array<std::vector< Mixture >,5>& mixture ) { double pdf(0.); @@ -329,7 +322,7 @@ double Trk::MultiComponentStateModeCalculator::width( int i, std::vector< Mixtur -double Trk::MultiComponentStateModeCalculator::findRoot(double &result, double xlo, double xhi, double value, double i, std::vector< Mixture >* mixture, MsgStream &log) +double Trk::MultiComponentStateModeCalculator::findRoot(double &result, double xlo, double xhi, double value, double i, std::array<std::vector< Mixture >,5>& mixture, MsgStream &log) { // Do the root finding using the Brent-Decker method. Returns a boolean status and // loads 'result' with our best guess at the root if true. @@ -338,8 +331,8 @@ double Trk::MultiComponentStateModeCalculator::findRoot(double &result, double x double a(xlo),b(xhi); - double fa= pdf(a, i, mixture) - value; - double fb= pdf(b, i, mixture) - value; + double fa= pdf(a,i,mixture) - value; + double fb= pdf(b,i,mixture) - value; if(fb*fa > 0) { diff --git a/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.h b/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.h index dfd78eda76c595af3f47b582fc37baeb321211cf..2d73b82d51bb6df13c364f97bea171a522411368 100755 --- a/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.h +++ b/Tracking/TrkFitter/TrkGaussianSumFilter/src/MultiComponentStateModeCalculator.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration */ /*********************************************************************************** @@ -8,15 +8,15 @@ begin : Thursday 6th July 2006 author : atkinson email : Tom.Atkinson@cern.ch -description : Class to calculate the mode (q/p) of a gaussian mixture +description : Class to calculate the mode (q/p) of a gaussian mixtureArray ***********************************************************************************/ #ifndef Trk_MultiComponentStateModeCalculator_H #define Trk_MultiComponentStateModeCalculator_H - #include "GaudiKernel/MsgStream.h" #include "EventPrimitives/EventPrimitives.h" +#include <array> namespace Trk{ @@ -24,56 +24,52 @@ class MultiComponentState; namespace MultiComponentStateModeCalculator{ - struct Mixture{ // Default constructor - Mixture() - : + Mixture(): weight(0.), mean(0.), sigma(0.) {} // Constructor with arguments - Mixture( double aWeight, double aMean, double aSigma ) - : + Mixture( double aWeight, double aMean, double aSigma ): weight( aWeight ), mean( aMean ), sigma( aSigma ) {} - double weight; double mean; double sigma; }; - //!< IMultiComponentStateModeCalculator interface method to calculate mode + //!< IMultiComponentStateModeCalculator interface method to calculate mode Amg::VectorX calculateMode( const MultiComponentState&, MsgStream &log ) ; //!< Private method to extract the weight, mean and sigma values from the multi-component state - void fillMixture( const MultiComponentState&, std::vector< Mixture >* ) ; + void fillMixture( const MultiComponentState&, std::array<std::vector< Mixture >,5>& ) ; //!< Private method to find the mode using the Newton-Raphson method based on a starting guess - double findMode( double, int, std::vector< Mixture >* mixture, MsgStream &log ) ; + double findMode( double, int, std::array<std::vector< Mixture >,5>& mixture, MsgStream &log ) ; //!< Private method to determine the pdf of the cashed mixture at a given value - double pdf( double, int, std::vector< Mixture >* mixture ) ; + double pdf( double, int, std::array<std::vector< Mixture >,5>& mixture ) ; //!< Private method to determine the first order derivative of the pdf at a given value - double d1pdf( double, int, std::vector< Mixture >* mixture ) ; + double d1pdf( double, int, std::array<std::vector< Mixture >,5>& mixture ) ; //!< Private method to determine the second order derivative of the pdf at a given value - double d2pdf( double, int, std::vector< Mixture >* mixture) ; + double d2pdf( double, int, std::array<std::vector< Mixture >,5>& mixture) ; //!< Private method to determine the value of the a gaussian distribution at a given value double gaus( double x, double mean, double sigma ) ; - double findModeGlobal(double, int,std::vector< Mixture >* mixture) ; + double findModeGlobal(double, int,std::array<std::vector< Mixture >,5>& mixture) ; - double width( int i, std::vector< Mixture >* mixture) ; + double width( int i, std::array<std::vector< Mixture >,5>& mixture) ; - double findRoot(double &result, double xlo, double xhi, double value, double i,std::vector< Mixture >* mixture, MsgStream &log) ; + double findRoot(double &result, double xlo, double xhi, double value, double i,std::array<std::vector< Mixture >,5>& mixture, MsgStream &log) ; } // end MultiComponentStateModeCalculator namespace