// Standard header files go here #include <iostream> #include <cmath> #include <cstdlib> #include <sstream> #include <vector> #include <list> #include <algorithm> // for std::sort #include <utility> // For std::pair #include <cassert> // Includes check for correct Boost version(s) #include <common/GGlobalDefines.hpp> //#include "GGlobalDefines.hpp" // Boost header files go here #include <boost/shared_ptr.hpp> #include <boost/bind.hpp> #include <boost/filesystem/operations.hpp> #include <boost/algorithm/string/trim.hpp> #include <boost/cast.hpp> #ifndef GOmegaPiIndividual_HPP_ #define GOmegaPiIndividual_HPP_ // For Microsoft-compatible compilers #if defined(_MSC_VER) && (_MSC_VER >= 1020) #pragma once #endif // Geneva header files go here #include <hap/GRandomT.hpp> #include <common/GCommonEnums.hpp> #include <common/GExceptions.hpp> #include <geneva/GConstrainedDoubleObject.hpp> #include <geneva/GConstrainedDoubleObjectCollection.hpp> #include <geneva/GDoubleGaussAdaptor.hpp> #include <geneva/GObjectExpectationChecksT.hpp> #include <geneva/GParameterObjectCollection.hpp> #include <geneva/GParameterSet.hpp> #ifdef GENEVATESTING #include <common/GUnitTestFrameworkT.hpp> #endif /* GENEVATESTING */ #include "PwaUtils/AbsStates.hh" #include "PwaUtils/DataUtils.hh" #include "Examples/MATpbarpToOmegaPi/OmegaPiData.hh" #include "Examples/MATpbarpToOmegaPi/OmegaPiEventList.hh" #include "Examples/MATpbarpToOmegaPi/pbarpToOmegaPi0States.hh" #include "Examples/MATpbarpToOmegaPi/OmegaPiLh.hh" namespace Gem { namespace Geneva { /************************************************************************************************/ /** * This individual searches for the minimum of a parabola of a given dimension, * It is part of a complete example that lets users adapt their optimization * problems more easily to the Geneva conventions. */ class GOmegaPiIndividual :public GParameterSet { /////////////////////////////////////////////////////////////////////// friend class boost::serialization::access; template<typename Archive> void serialize(Archive & ar, const unsigned int) { using boost::serialization::make_nvp; ar & BOOST_SERIALIZATION_BASE_OBJECT_NVP(GParameterSet); ar & BOOST_SERIALIZATION_NVP(_omegaPiLhPtr); //ar & BOOST_SERIALIZATION_NVP(_barpToOmegaPi0States); } /////////////////////////////////////////////////////////////////////// public: /********************************************************************************************/ /** * A simple constructor that initializes this object with a collection of bounded * double variables. * * @param dim The amount of variables * @param min The lower boundary of the variables * @param max The upper boundary of the variables */ GOmegaPiIndividual(boost::shared_ptr<pbarpToOmegaPi0States> theStates) : GParameterSet() , _omegaPiLhPtr(new OmegaPiLh(theStates)) , _barpToOmegaPi0States(theStates) { //_omegaPiLhPtr( new OmegaPiLh(theStates) ); // Set up a GConstrainedDoubleObjectCollection boost::shared_ptr<GConstrainedDoubleObjectCollection> gbdc_ptr(new GConstrainedDoubleObjectCollection()); // Create a suitable adaptor (sigma=0.1, sigma-adaption=0.5, min sigma=0, max sigma=0,5) boost::shared_ptr<GDoubleGaussAdaptor> gdga_ptr(new GDoubleGaussAdaptor(0.1, 0.5, 0., 0.5)); gdga_ptr->setAdaptionThreshold(1); // Adaption parameters are modified after each adaption gdga_ptr->setAdaptionProbability(0.05); // The likelihood for a parameter to be adapted // Add a GConstrainedDoubleObject object to the collection // gbdc_ptr->push_back(gbd_ptr); // gpoc_ptr->push_back(gbd_ptr); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaSinglet=_barpToOmegaPi0States->jpclsSinglet(); setFitParamVal(JPCLSOmegaSinglet, gbdc_ptr); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet0=_barpToOmegaPi0States->jpclsTriplet0(); setFitParamVal(JPCLSOmegaTriplet0, gbdc_ptr); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet1=_barpToOmegaPi0States->jpclsTriplet1(); setFitParamVal(JPCLSOmegaTriplet1, gbdc_ptr); for(std::size_t i=0; i<gbdc_ptr->size(); i++) { // Register the adaptor with GConstrainedDoubleObject objects gbdc_ptr->at(i)->addAdaptor(gdga_ptr); } // Add the collection to this object this->push_back(gbdc_ptr); } /********************************************************************************************/ /** * A standard copy constructor * * @param cp A copy of another GOmegaPiIndividual */ GOmegaPiIndividual(const GOmegaPiIndividual& cp) : GParameterSet(cp) ,_omegaPiLhPtr( new OmegaPiLh(cp.omegaPiLhPtr()) ) { _barpToOmegaPi0States=_omegaPiLhPtr->omegaPi0States(); } /********************************************************************************************/ /** * The standard destructor */ virtual ~GOmegaPiIndividual() { } /********************************************************************************************/ boost::shared_ptr<OmegaPiLh> omegaPiLhPtr() const {return _omegaPiLhPtr;} /** * A standard assignment operator * * @param cp A copy of another GOmegaPiIndividual object * @return A constant reference to this object */ const GOmegaPiIndividual& operator=(const GOmegaPiIndividual& cp){ GOmegaPiIndividual::load_(&cp); return *this; } /*******************************************************************************************/ /** * Checks for equality with another GOmegaPiIndividual object. * * NOTE: THIS FUNCTION IS OPTIONAL AND IS MAINLY USED IN CONJUNCTION WITH UNIT TESTS. * You do not need it if you do not intend to perform unit tests. * * @param cp A constant reference to another GOmegaPiIndividual object * @return A boolean indicating whether both objects are equal */ bool operator==(const GOmegaPiIndividual& cp) const { using namespace Gem::Common; // Means: The expectation of equality was fulfilled, if no error text was emitted (which converts to "true") return !checkRelationshipWith(cp, CE_EQUALITY, 0.,"GOmegaPiIndividual::operator==","cp", CE_SILENT); } /*******************************************************************************************/ /** * Checks for inequality with another GOmegaPiIndividual object. * * NOTE: THIS FUNCTION IS OPTIONAL AND IS MAINLY USED IN CONJUNCTION WITH UNIT TESTS. * You do not need it if you do not intend to perform unit tests. * * @param cp A constant reference to another GOmegaPiIndividual object * @return A boolean indicating whether both objects are inequal */ bool operator!=(const GOmegaPiIndividual& cp) const { using namespace Gem::Common; // Means: The expectation of inequality was fulfilled, if no error text was emitted (which converts to "true") return !checkRelationshipWith(cp, CE_INEQUALITY, 0.,"GOmegaPiIndividual::operator!=","cp", CE_SILENT); } bool getFitParams(OmegaPiData::fitParamVal& fitParmVal) { std::vector<double> theParms; // Extract the GDoubleCollection object. In a realistic scenario, you might want // to add error checks here upon first invocation. boost::shared_ptr<GConstrainedDoubleObjectCollection> vC = at<GConstrainedDoubleObjectCollection>(0); std::vector< boost::shared_ptr<const JPCLS> >::const_iterator itJPCLS; std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaSinglet=_barpToOmegaPi0States->jpclsSinglet(); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet0=_barpToOmegaPi0States->jpclsTriplet0(); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet1=_barpToOmegaPi0States->jpclsTriplet1(); std::vector<double> par; for(std::size_t i=0; i<vC->size(); i++){ double value = vC->at(i)->value(); par.push_back(value); } if (par.size()!= JPCLSOmegaSinglet.size()*2+JPCLSOmegaTriplet0.size()*2+JPCLSOmegaTriplet1.size()*2-3){ std::cout << "size of parameters wrong!!! par.size()=" << par.size() << "\tJPCLSOmegaSinglet.size()+JPCLSOmegaTriplet0.size()+JPCLSOmegaTriplet1.size()-3=" << JPCLSOmegaSinglet.size()*2+JPCLSOmegaTriplet0.size()*2+JPCLSOmegaTriplet1.size()*2-3 << std::endl; assert(0); } int counter=0; for ( itJPCLS=JPCLSOmegaSinglet.begin(); itJPCLS!=JPCLSOmegaSinglet.end(); ++itJPCLS){ //now fill the fitParameterMap double mag=par[counter]; counter++; double phi=0.; if (counter>1){ phi=par[counter]; counter++; } std::pair <double,double> tmpParameter=make_pair(mag,phi); fitParmVal.omegaProdSinglet[(*itJPCLS)]=tmpParameter; } for ( itJPCLS=JPCLSOmegaTriplet0.begin(); itJPCLS!=JPCLSOmegaTriplet0.end(); ++itJPCLS){ //now fill the fitParameterMap double mag=par[counter]; counter++; double phi=0.; if (counter>JPCLSOmegaSinglet.size()*2){ phi=par[counter]; counter++; } std::pair <double,double> tmpParameter=make_pair(mag,phi); fitParmVal.omegaProdTriplet0[(*itJPCLS)]=tmpParameter; } for ( itJPCLS=JPCLSOmegaTriplet1.begin(); itJPCLS!=JPCLSOmegaTriplet1.end(); ++itJPCLS){ //now fill the fitParameterMap double mag=par[counter]; counter++; double phi=0.; if (counter>JPCLSOmegaSinglet.size()*2+JPCLSOmegaTriplet0.size()*2-1){ phi=par[counter]; counter++; } std::pair <double,double> tmpParameter=make_pair(mag,phi); fitParmVal.omegaProdTriplet1[(*itJPCLS)]=tmpParameter; } return true; } void printFitParams(OmegaPiData::fitParamVal& fitParmVal) { std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaSinglet=_barpToOmegaPi0States->jpclsSinglet(); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet0=_barpToOmegaPi0States->jpclsTriplet0(); std::vector< boost::shared_ptr<const JPCLS> > JPCLSOmegaTriplet1=_barpToOmegaPi0States->jpclsTriplet1(); std::vector< boost::shared_ptr<const JPCLS> >::const_iterator itJPCLS; std::cout << "***fit parameter singlet states*** " <<std::endl; for ( itJPCLS=JPCLSOmegaSinglet.begin(); itJPCLS!=JPCLSOmegaSinglet.end(); ++itJPCLS){ std::cout << (*itJPCLS)->name()<< "\t"; std::pair<double, double> tmpParam=fitParmVal.omegaProdSinglet[(*itJPCLS)]; std::cout <<"\t mag:" << tmpParam.first <<"\t phi:" << tmpParam.second << std::endl; } std::cout << "***fit parameter triplet m=0 states*** " <<std::endl; for ( itJPCLS=JPCLSOmegaTriplet0.begin(); itJPCLS!=JPCLSOmegaTriplet0.end(); ++itJPCLS){ std::cout<< (*itJPCLS)->name()<< "\t"; std::pair<double, double> tmpParam=fitParmVal.omegaProdTriplet0[(*itJPCLS)]; std::cout <<"\t mag:" << tmpParam.first <<"\t phi:" << tmpParam.second << std::endl; } std::cout << "***fit parameter triplet m=1 states*** " <<std::endl; for ( itJPCLS=JPCLSOmegaTriplet1.begin(); itJPCLS!=JPCLSOmegaTriplet1.end(); ++itJPCLS){ std::cout<< (*itJPCLS)->name()<< "\t"; std::pair<double, double> tmpParam=fitParmVal.omegaProdTriplet1[(*itJPCLS)]; std::cout <<"\t mag:" << tmpParam.first <<"\t phi:" << tmpParam.second << std::endl; } std::cout << std::endl; return; } boost::shared_ptr<OmegaPiLh> getOmegaPiLhPtr() {return _omegaPiLhPtr;} /*******************************************************************************************/ /** * Checks whether a given expectation for the relationship between this object and another object * is fulfilled. * * NOTE: THIS FUNCTION IS OPTIONAL AND IS MAINLY USED IN CONJUNCTION WITH UNIT TESTS. * You do not need it if you do not intend to perform unit tests. * * @param cp A constant reference to another object, camouflaged as a GObject * @param e The expected outcome of the comparison * @param limit The maximum deviation for floating point values (important for similarity checks) * @param caller An identifier for the calling entity * @param y_name An identifier for the object that should be compared to this one * @param withMessages Whether or not information should be emitted in case of deviations from the expected outcome * @return A boost::optional<std::string> object that holds a descriptive string if expectations were not met */ boost::optional<std::string> checkRelationshipWith(const GObject& cp, const Gem::Common::expectation& e, const double& limit, const std::string& caller, const std::string& y_name, const bool& withMessages) const { using namespace Gem::Common; // Check that we are not accidently assigning this object to itself selfAssignmentCheck<GOmegaPiIndividual>(&cp); // Use this call instead when local data needs to be loaded: // const GOmegaPiIndividual *p_load = GObject::conversion_cast<GOmegaPiIndividual>(&cp); // Will hold possible deviations from the expectation, including explanations std::vector<boost::optional<std::string> > deviations; // Check our parent class'es data ... deviations.push_back(GParameterSet::checkRelationshipWith(cp, e, limit, "GOmegaPiIndividual", y_name, withMessages)); // Check local data like this: // deviations.push_back(checkExpectation(withMessages, "GOmegaPiIndividual", val_, p_load->val_, "val_", "p_load->val_", e , limit)); return evaluateDiscrepancies("GOmegaPiIndividual", caller, deviations, e); } #ifdef GENEVATESTING /*******************************************************************************************/ /** * Applies modifications to this object. This is needed for testing purposes * * @return A boolean which indicates whether modifications were made */ virtual bool modify_GUnitTests() { bool result; // Call the parent class'es function if(GParameterSet::modify_GUnitTests()) result = true; // Check that this individual actually contains data to be modified if(this->size() != 0) { adapt(); // Perform modifications result = true; } return result; } /*******************************************************************************************/ /** * Performs self tests that are expected to succeed. This is needed for testing purposes */ virtual void specificTestsNoFailureExpected_GUnitTests() { const boost::uint32_t NITERATIONS=100; // Call the parent class'es function GParameterSet::specificTestsNoFailureExpected_GUnitTests(); // Create an individual boost::shared_ptr<Gem::Geneva::GStartIndividual> p = boost::shared_ptr<Gem::Geneva::GStartIndividual>(new GStartIndividual(1000, -10, 10)); // Adapt a number of times and check that there were changes double oldfitness = p->fitness(); for(boost::uint32_t i=0; i<NITERATIONS; i++) { p->adapt(); double newfitness = p->fitness(); BOOST_CHECK_MESSAGE(newfitness != oldfitness, "Rare failures are normal for this test / " << i << "/" << NITERATIONS); oldfitness = newfitness; } } /*******************************************************************************************/ /** * Performs self tests that are expected to fail. This is needed for testing purposes */ virtual void specificTestsFailuresExpected_GUnitTests() { // Call the parent class'es function GParameterSet::specificTestsFailuresExpected_GUnitTests(); } #endif /* GENEVATESTING */ protected: /********************************************************************************************/ /** * Loads the data of another GOmegaPiIndividual, camouflaged as a GObject. * * @param cp A copy of another GOmegaPiIndividual, camouflaged as a GObject */ virtual void load_(const GObject* cp) { // Check that we are not accidently assigning this object to itself selfAssignmentCheck<GOmegaPiIndividual>(cp); // Use this call instead when local data needs to be loaded: // const GOmegaPiIndividual *p_load = GObject::conversion_cast<GOmegaPiIndividual>(cp); // Load our parent's data GParameterSet::load_(cp); // Load local data here like this: // myVar = p_load->myVar; } /********************************************************************************************/ /** * Creates a deep clone of this object * * @return A deep clone of this object, camouflaged as a GObject */ virtual GObject* clone_() const { return new GOmegaPiIndividual(*this); } /********************************************************************************************/ /** * The actual fitness calculation takes place here. * * @return The value of this object */ virtual double fitnessCalculation(){ double result = 0.; OmegaPiData::fitParamVal theFitParmValTmp; assert(getFitParams(theFitParmValTmp)); result=_omegaPiLhPtr->calcLogLh(theFitParmValTmp); return result; } private: /********************************************************************************************/ /** * The default constructor. Intentionally private and empty, as it is only needed for * serialization purposes. */ boost::shared_ptr<OmegaPiLh> _omegaPiLhPtr; boost::shared_ptr<pbarpToOmegaPi0States> _barpToOmegaPi0States; GOmegaPiIndividual() :GParameterSet() { _barpToOmegaPi0States = pbarpToOmegaPi0States::getStates(); } void setFitParamVal(std::vector< boost::shared_ptr<const JPCLS> > theJPCLSs, boost::shared_ptr<GConstrainedDoubleObjectCollection> theGbdc_ptr){ std::vector< boost::shared_ptr<const JPCLS> >::const_iterator itJPCLS; int counter=0; for ( itJPCLS=theJPCLSs.begin(); itJPCLS!=theJPCLSs.end(); ++itJPCLS){ //now fill the fitParameterMap boost::shared_ptr<GConstrainedDoubleObject> gbd_ptr(new GConstrainedDoubleObject(0., 1.) ); //JPCLS magnitude theGbdc_ptr->push_back(gbd_ptr); if (counter>0){ boost::shared_ptr<GConstrainedDoubleObject> gbd_ptr(new GConstrainedDoubleObject(-M_PI, M_PI) ); //JPCLS phi theGbdc_ptr->push_back(gbd_ptr); } counter++; } } /********************************************************************************************/ // You can add other variables here. Do not forget to serialize them if necessary // int myVar; }; /*************************************************************************************************/ } /* namespace Geneva */ } /* namespace Gem */ // Needed for serialization purposes #include <boost/serialization/export.hpp> BOOST_CLASS_EXPORT(Gem::Geneva::GOmegaPiIndividual) // Needed for testing purposes /*************************************************************************************************/ /////////////////////////////////////////////////////////////////////////////////////////////////// /*************************************************************************************************/ /** * As the Gem::Geneva::Gem::Geneva::GOmegaPiIndividual has a private default constructor, we need to provide a * specialization of the factory function that creates GStartProjectIndividual objects */ //template <> //boost::shared_ptr<Gem::Geneva::GOmegaPiIndividual> TFactory_GUnitTests<Gem::Geneva::GOmegaPiIndividual>() { // return boost::shared_ptr<Gem::Geneva::GOmegaPiIndividual>(new Gem::Geneva::GOmegaPiIndividual(1000,-10.,10.)); //} /*************************************************************************************************/ /////////////////////////////////////////////////////////////////////////////////////////////////// /*************************************************************************************************/ #endif /* GSTARTINDIVIDUAL_HPP_ */