VectorRule.h
1//-*-C++-*-
2/***************************************************************************
3 *
4 * Copyright (C) 2005 - 2015 by Willem van Straten
5 * Licensed under the Academic Free License version 2.1
6 *
7 ***************************************************************************/
8
9// psrchive/More/MEAL/MEAL/VectorRule.h
10
11#ifndef __MEAL_VectorRule_H
12#define __MEAL_VectorRule_H
13
14#include "MEAL/ProjectGradient.h"
15#include "MEAL/Composite.h"
16#include "stringtok.h"
17
18namespace MEAL {
19
21 template<class T>
22 class VectorRule : public T
23 {
24
25 public:
26
27 typedef typename T::Result Result;
28
30 VectorRule (Composite* policy = 0);
31
33 VectorRule (const VectorRule& v);
34
37
40
42 void push_back (T* model);
43
45 void assign (T* model);
46
48 void erase (unsigned index);
49
52
55
57 unsigned size () const { return model.size(); }
58
60 void set_index (unsigned index);
61
63 unsigned get_index () const { return model_index; }
64
66 std::string get_name () const
67 { return "VectorRule<" + std::string(T::Name)+ ">"; }
68
69 protected:
70
72 void calculate (Result& result, std::vector<Result>* grad);
73
74 private:
75
77 std::vector< Project<T> > model;
78
80 unsigned model_index;
81
84
85 };
86
87}
88
89template<class T>
91{
92 if (policy)
93 composite = policy;
94 else
95 composite = new Composite (this);
96
97 model_index = 0;
98}
99
101template<class T>
103{
104 composite = new Composite (this);
105 operator = (v);
106}
107
108template<class T>
111{
112 if (this == &copy)
113 return *this;
114
115 unsigned nmodel = copy.model.size();
116 for (unsigned imodel=0; imodel < nmodel; imodel++)
117 push_back (copy.model[imodel]);
118
119 model_index = copy.model_index;
120
121 return *this;
122}
123
124template<class T>
126{
127 if (T::very_verbose)
128 std::cerr << get_name() + "push_back" << std::endl;
129
130 model_index = model.size();
131 model.push_back (Project<T>(x));
132 composite->map (model.back());
133}
134
135template<class T>
137{
138 if (T::very_verbose)
139 std::cerr << get_name() + "assign" << std::endl;
140
141 if (!model.size())
142 {
143 push_back (x);
144 return;
145 }
146
147 if (model[model_index])
148 {
149 if (T::very_verbose)
150 std::cerr << get_name() + "assign unmap old" << std::endl;
151
152 composite->unmap (model[model_index]);
153 }
154
155 model[model_index] = x;
156
157 if (!x)
158 return;
159
160 if (T::very_verbose)
161 std::cerr << get_name() + "assign map new" << std::endl;
162
163 composite->map (model[model_index]);
164}
165
166template<class T>
167void MEAL::VectorRule<T>::erase (unsigned index)
168{
169 if (T::very_verbose)
170 std::cerr << get_name() + "erase" << std::endl;
171
172 if (!model[index])
173 return;
174
175 if (T::very_verbose)
176 std::cerr << get_name() + "erase unmap" << std::endl;
177
178 composite->unmap (model[index]);
179
180 if (T::very_verbose)
181 std::cerr << get_name() + "erase vector::erase" << std::endl;
182
183 model.erase( model.begin() + index );
184}
185
186template<class T>
188{
189 if (model_index >= model.size())
190 throw Error (InvalidRange, "MEAL::"+get_name()+"::get_current",
191 "index=%d >= nmodel=%d", model_index, model.size());
192
193 if (!model[model_index])
194 throw Error (InvalidState, "MEAL::"+get_name()+"::get_current",
195 "index=%d not set; nmodel=%d", model_index);
196
197 return model[model_index];
198}
199
200template<class T>
202{
203 if (model_index >= model.size())
204 throw Error (InvalidRange, "MEAL::"+get_name()+"::get_projection",
205 "index=%d >= nmodel=%d", model_index, model.size());
206
207 if (!model[model_index])
208 throw Error (InvalidState, "MEAL::"+get_name()+"::get_projection",
209 "index=%d not set; nmodel=%d", model_index);
210
211 return model[model_index];
212}
213
214template<class T>
216{
217 if (index == model_index)
218 return;
219
220 if (index >= model.size())
221 throw Error (InvalidRange, "MEAL::"+get_name()+"::set_index",
222 "index=%d >= nmodel=%d", index, model.size());
223
224 model_index = index;
225 composite->attribute_changed( Function::Evaluation );
226}
227
228template<class T>
230 std::vector<Result>* grad)
231{
232 unsigned nmodel = model.size();
233 if (T::very_verbose)
234 std::cerr << get_name() + "calculate nmodel=" << nmodel << std::endl;
235
236 if (nmodel == 0)
237 throw Error (InvalidState, "MEAL::"+get_name()+"::calculate", "nmodel = 0");
238
239 // the gradient of each component
240 std::vector<Result> comp_gradient;
241
242 // the pointer to the above array, if grad != 0
243 std::vector<Result>* comp_gradient_ptr = 0;
244
245 if (grad)
246 comp_gradient_ptr = &comp_gradient;
247
248 if (this->get_verbose())
249 std::cerr << get_name() + "calculate evaluate "
250 << model[model_index]->get_name() << std::endl;
251
252 try {
253
254 // evaluate the model and its gradient
255 result = model[model_index]->evaluate (comp_gradient_ptr);
256
257 if (this->get_verbose())
258 std::cerr << get_name() + "calculate "
259 << model[model_index]->get_name()
260 << " result=" << result << std::endl;
261
262 }
263 catch (Error& error) {
264 error += get_name() + "calculate";
265 throw error << " model=" << model[model_index]->get_name();
266 }
267
268 if (grad)
269 {
270 if (model[model_index]->get_nparam() != comp_gradient.size())
271 throw Error (InvalidState, (get_name() + "calculate").c_str(),
272 "model[%d]=%s.get_nparam=%d != gradient.size=%d",
273 model_index, model[model_index]->get_name().c_str(),
274 model[model_index]->get_nparam(), comp_gradient.size());
275
276 unsigned nparam = this->get_nparam();
277
278 grad->resize (nparam);
279 for (unsigned iparam=0; iparam<nparam; iparam++)
280 (*grad)[iparam] = 0.0;
281
282 // re-map the elements of the component gradient into the Composite space
283 ProjectGradient (model[model_index], comp_gradient, *grad);
284 }
285
286 if (T::very_verbose)
287 {
288 std::cerr << get_name() + "calculate model index=" << model_index
289 << " result\n " << result << std::endl;
290 if (grad)
291 {
292 std::cerr << get_name() + "calculate gradient" << std::endl;
293 for (unsigned i=0; i<grad->size(); i++)
294 std::cerr << " " << i << ":" << this->get_infit(i)
295 << "=" << (*grad)[i] << std::endl;
296 }
297 }
298}
299
300#endif
301
Parameter policy for composite functions.
Definition Composite.h:20
virtual void copy(const Function *model)
Does the work for operator =.
Definition Function.C:58
@ Evaluation
Function evaluation, as returned by the evaluate method.
Definition Function.h:171
Template combines a reference to a Component and its Projection.
Definition Projection.h:60
VectorRule implements a vector of Functions with a current index.
Definition VectorRule.h:23
unsigned size() const
Get the size of the array.
Definition VectorRule.h:57
void push_back(T *model)
Add an element to the array.
Definition VectorRule.h:125
std::string get_name() const
Return the name of the class.
Definition VectorRule.h:66
void assign(T *model)
Assign the current element of the array.
Definition VectorRule.h:136
void erase(unsigned index)
Erase an element from the array.
Definition VectorRule.h:167
~VectorRule()
Destructor.
Definition VectorRule.h:39
T * get_current()
Get the current element of the array.
Definition VectorRule.h:187
Project< T > & get_projection()
Get the projection of the current element of the array.
Definition VectorRule.h:201
VectorRule(Composite *policy=0)
Default constructor.
Definition VectorRule.h:90
unsigned get_index() const
Get the current element of the array.
Definition VectorRule.h:63
void calculate(Result &result, std::vector< Result > *grad)
Return the result and its gradient.
Definition VectorRule.h:229
void set_index(unsigned index)
Set the current element of the array.
Definition VectorRule.h:215
VectorRule & operator=(const VectorRule &copy)
Assignment operator.
Definition VectorRule.h:110
Namespace in which all modeling and calibration related code is declared.
Definition ExampleComplex2.h:16

Generated using doxygen 1.14.0