1 //                                               -*- C++ -*-
2 /**
3  *  @brief Factory for MeixnerDistribution distribution
4  *
5  *  Copyright 2005-2021 Airbus-EDF-IMACS-ONERA-Phimeca
6  *
7  *  This library is free software: you can redistribute it and/or modify
8  *  it under the terms of the GNU Lesser General Public License as published by
9  *  the Free Software Foundation, either version 3 of the License, or
10  *  (at your option) any later version.
11  *
12  *  This library is distributed in the hope that it will be useful,
13  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  *  GNU Lesser General Public License for more details.
16  *
17  *  You should have received a copy of the GNU Lesser General Public License
18  *  along with this library.  If not, see <http://www.gnu.org/licenses/>.
19  *
20  */
21 #include "openturns/MeixnerDistributionFactory.hxx"
22 #include "openturns/SpecFunc.hxx"
23 #include "openturns/PersistentObjectFactory.hxx"
24 
25 BEGIN_NAMESPACE_OPENTURNS
26 
27 CLASSNAMEINIT(MeixnerDistributionFactory)
28 
29 static const Factory<MeixnerDistributionFactory> Factory_MeixnerDistributionFactory;
30 
31 /* Default constructor */
MeixnerDistributionFactory()32 MeixnerDistributionFactory::MeixnerDistributionFactory()
33   : DistributionFactoryImplementation()
34 {
35   // Nothing to do
36 }
37 
38 /* Virtual constructor */
clone() const39 MeixnerDistributionFactory * MeixnerDistributionFactory::clone() const
40 {
41   return new MeixnerDistributionFactory(*this);
42 }
43 
44 
45 /* Here is the interface that all derived class must implement */
46 
build(const Sample & sample) const47 Distribution MeixnerDistributionFactory::build(const Sample & sample) const
48 {
49   return buildAsMeixnerDistribution(sample).clone();
50 }
51 
build(const Point & parameters) const52 Distribution MeixnerDistributionFactory::build(const Point & parameters) const
53 {
54   return buildAsMeixnerDistribution(parameters).clone();
55 }
56 
build() const57 Distribution MeixnerDistributionFactory::build() const
58 {
59   return buildAsMeixnerDistribution().clone();
60 }
61 
buildAsMeixnerDistribution(const Sample & sample) const62 MeixnerDistribution MeixnerDistributionFactory::buildAsMeixnerDistribution(const Sample & sample) const
63 {
64   UnsignedInteger size = sample.getSize();
65   if (size < 4) throw InvalidArgumentException(HERE) << "Error: cannot build a MeixnerDistribution distribution from a sample of size < 4";
66   if (sample.getDimension() != 1) throw InvalidArgumentException(HERE) << "Error: can build a MeixnerDistribution distribution only from a sample of dimension 1, here dimension=" << sample.getDimension();
67   const Scalar gamma1 = sample.computeSkewness()[0];
68   const Scalar gamma2 = sample.computeKurtosis()[0];
69   const Scalar upperBound = 3.0 + 2.0 * gamma1 * gamma1;
70   if (gamma2 <= upperBound) throw InvalidArgumentException(HERE) << "Error: cannot estimate a MeixnerDistribution distribution if the sample kurtosis=" << gamma2 << " is not greater than 2*skewness^2+3=" << upperBound;
71   const Scalar m = sample.computeMean()[0];
72   const Scalar s2 = sample.computeVariance()[0];
73   const Scalar delta = 1.0 / (gamma2 - gamma1 * gamma1 - 3.0);
74   const Scalar beta = ((0.0 < gamma1) - (gamma1 < 0.0)) * std::acos(2.0 - delta * (gamma2 - 3.0));
75   const Scalar alpha = cbrt(s2 * (std::cos(beta) + 1.0));
76   const Scalar mu = m - alpha * delta * std::tan(0.5 * beta);
77   MeixnerDistribution result(alpha, beta, delta, mu);
78   result.setDescription(sample.getDescription());
79   return result;
80 }
81 
buildAsMeixnerDistribution(const Point & parameters) const82 MeixnerDistribution MeixnerDistributionFactory::buildAsMeixnerDistribution(const Point & parameters) const
83 {
84   try
85   {
86     MeixnerDistribution distribution;
87     distribution.setParameter(parameters);
88     return distribution;
89   }
90   catch (InvalidArgumentException &)
91   {
92     throw InvalidArgumentException(HERE) << "Error: cannot build a MeixnerDistribution distribution from the given parameters";
93   }
94 }
95 
buildAsMeixnerDistribution() const96 MeixnerDistribution MeixnerDistributionFactory::buildAsMeixnerDistribution() const
97 {
98   return MeixnerDistribution();
99 }
100 
101 END_NAMESPACE_OPENTURNS
102