1 // -*- C++ -*-
2 /**
3 * @brief This class build a stationary covariance model using a time grid and a collection of covariance matrices
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/UserDefinedStationaryCovarianceModel.hxx"
22 #include "openturns/PersistentObjectFactory.hxx"
23 #include "openturns/Exception.hxx"
24
25 BEGIN_NAMESPACE_OPENTURNS
26
27 CLASSNAMEINIT(UserDefinedStationaryCovarianceModel)
28
29 static const Factory<UserDefinedStationaryCovarianceModel> Factory_UserDefinedStationaryCovarianceModel;
30
31 /* Constructor with parameters */
UserDefinedStationaryCovarianceModel()32 UserDefinedStationaryCovarianceModel::UserDefinedStationaryCovarianceModel()
33 : CovarianceModelImplementation()
34 , covarianceCollection_(0)
35 , mesh_()
36 , nearestNeighbour_()
37 {
38 inputDimension_ = 1;
39 outputDimension_ = 0;
40 isStationary_ = true;
41 }
42
43 // Classical constructor
44 // For a stationary model, we need N covariance matrices with N the number of time stamps in the time grid
UserDefinedStationaryCovarianceModel(const RegularGrid & mesh,const SquareMatrixCollection & covarianceFunction)45 UserDefinedStationaryCovarianceModel::UserDefinedStationaryCovarianceModel(const RegularGrid & mesh,
46 const SquareMatrixCollection & covarianceFunction)
47 : CovarianceModelImplementation()
48 , covarianceCollection_(0)
49 , mesh_(mesh)
50 , nearestNeighbour_(mesh)
51 {
52 isStationary_ = true;
53 const UnsignedInteger size = mesh.getVerticesNumber();
54 if (size != covarianceFunction.getSize())
55 throw InvalidArgumentException(HERE) << "Error: for a non stationary covariance model, sizes are incoherent"
56 << " mesh size = " << size << "covariance function size = " << covarianceFunction.getSize();
57 inputDimension_ = mesh.getDimension();
58 covarianceCollection_ = SquareMatrixCollection(size);
59 // put the first element
60 covarianceCollection_[0] = covarianceFunction[0];
61 outputDimension_ = covarianceCollection_[0].getDimension();
62 // put the next elements if dimension is ok
63 for (UnsignedInteger k = 1; k < size; ++k)
64 {
65 if (covarianceFunction[k].getDimension() != outputDimension_)
66 throw InvalidArgumentException(HERE) << " Error with dimension; the covariance matrices should be of same dimension";
67 covarianceCollection_[k] = covarianceFunction[k];
68 }
69 }
70
71 /* Virtual constructor */
clone() const72 UserDefinedStationaryCovarianceModel * UserDefinedStationaryCovarianceModel::clone() const
73 {
74 return new UserDefinedStationaryCovarianceModel(*this);
75 }
76
computeAsScalar(const Point & tau) const77 Scalar UserDefinedStationaryCovarianceModel::computeAsScalar(const Point &tau) const
78 {
79 if (outputDimension_ > 1)
80 throw InvalidArgumentException(HERE) << "Error : UserDefinedStationaryCovarianceModel::computeAsScalar(tau) should be only used if output dimension is 1. Here, output dimension = " << outputDimension_;
81 if (tau.getDimension() != inputDimension_)
82 throw InvalidArgumentException(HERE) << "In UserDefinedStationaryCovarianceModel::computeStandardRepresentative: expected a shift of dimension=" << getInputDimension() << ", got dimension=" << tau.getDimension();
83 // We filter the collection & return corresponding values
84 if (mesh_.getN() == 1)
85 return covarianceCollection_[0](0, 0);
86 UnsignedInteger index;
87 if (tau[0] < 0.0)
88 index = nearestNeighbour_.query(tau * (-1.0));
89 else
90 index = nearestNeighbour_.query(tau);
91 // index
92 return covarianceCollection_[index](0, 0);
93 }
94
computeAsScalar(const Collection<Scalar>::const_iterator & s_begin,const Collection<Scalar>::const_iterator & t_begin) const95 Scalar UserDefinedStationaryCovarianceModel::computeAsScalar(const Collection<Scalar>::const_iterator &s_begin,
96 const Collection<Scalar>::const_iterator &t_begin) const
97 {
98 if (outputDimension_ != 1)
99 throw InvalidArgumentException(HERE) << "Error : UserDefinedStationaryCovarianceModel::computeAsScalar(it, it) should be only used if output dimension is 1. Here, output dimension = " << outputDimension_;
100
101 // Work on iterators
102 // Unfortunately there is no other solution then generating the tau point
103 Collection<Scalar>::const_iterator s_it = s_begin;
104 Collection<Scalar>::const_iterator t_it = t_begin;
105 Point tau(inputDimension_, 0.0);
106 for (UnsignedInteger i = 0; i < inputDimension_; ++i, ++s_it, ++t_it)
107 {
108 tau[i] = *s_it - *t_it;
109 }
110 return computeAsScalar(tau);
111 }
112
113 /* Computation of the covariance function */
operator ()(const Point & tau) const114 SquareMatrix UserDefinedStationaryCovarianceModel::operator()(const Point & tau) const
115 {
116 if (tau.getDimension() != inputDimension_) throw InvalidArgumentException(HERE) << "Error: expected a shift of dimension=" << inputDimension_ << ", got dimension=" << tau.getDimension();
117 // If the grid size is one, return the covariance function
118 // else find in the grid the nearest instant values with nonnegative first component
119 if (mesh_.getN() == 1) return covarianceCollection_[0];
120
121 if (tau[0] < 0.0) return covarianceCollection_[nearestNeighbour_.query(tau * (-1.0))];
122 return covarianceCollection_[nearestNeighbour_.query(tau)];
123 }
124
discretize(const Mesh & mesh) const125 CovarianceMatrix UserDefinedStationaryCovarianceModel::discretize(const Mesh & mesh) const
126 {
127 return discretize(RegularGrid(mesh));
128 }
129
discretize(const Sample & vertices) const130 CovarianceMatrix UserDefinedStationaryCovarianceModel::discretize(const Sample & vertices) const
131 {
132 return discretize(Mesh(vertices));
133 }
134
135 /* Mesh accessor */
getTimeGrid() const136 RegularGrid UserDefinedStationaryCovarianceModel::getTimeGrid() const
137 {
138 return mesh_;
139 }
140
141 /* String converter */
__repr__() const142 String UserDefinedStationaryCovarianceModel::__repr__() const
143 {
144 OSS oss(true);
145 oss << "class=" << UserDefinedStationaryCovarianceModel::GetClassName()
146 << " mesh=" << mesh_.__repr__()
147 << " covarianceCollection=" << covarianceCollection_;
148 return oss;
149 }
150
151 /* String converter */
__str__(const String &) const152 String UserDefinedStationaryCovarianceModel::__str__(const String & ) const
153 {
154 return __repr__();
155 }
156
157 /* Method save() stores the object through the StorageManager */
save(Advocate & adv) const158 void UserDefinedStationaryCovarianceModel::save(Advocate & adv) const
159 {
160 CovarianceModelImplementation::save(adv);
161 adv.saveAttribute( "covarianceCollection_", covarianceCollection_);
162 adv.saveAttribute( "mesh_", mesh_);
163 }
164
165 /* Method load() reloads the object from the StorageManager */
load(Advocate & adv)166 void UserDefinedStationaryCovarianceModel::load(Advocate & adv)
167 {
168 CovarianceModelImplementation::load(adv);
169 adv.loadAttribute( "covarianceCollection_", covarianceCollection_);
170 adv.loadAttribute( "mesh_", mesh_);
171 }
172
173 END_NAMESPACE_OPENTURNS
174