src/utils/Accumulator.hpp

/*
 * Copyright (c) 2012 The University of Notre Dame. All Rights Reserved.
 *
 * The University of Notre Dame grants you ("Licensee") a
 * non-exclusive, royalty free, license to use, modify and
 * redistribute this software in source and binary code form, provided
 * that the following conditions are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 *    notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in the
 *    documentation and/or other materials provided with the
 *    distribution.
 *
 * This software is provided "AS IS," without a warranty of any
 * kind. All express or implied conditions, representations and
 * warranties, including any implied warranty of merchantability,
 * fitness for a particular purpose or non-infringement, are hereby
 * excluded.  The University of Notre Dame and its licensors shall not
 * be liable for any damages suffered by licensee as a result of
 * using, modifying or distributing the software or its
 * derivatives. In no event will the University of Notre Dame or its
 * licensors be liable for any lost revenue, profit or data, or for
 * direct, indirect, special, consequential, incidental or punitive
 * damages, however caused and regardless of the theory of liability,
 * arising out of the use of or inability to use software, even if the
 * University of Notre Dame has been advised of the possibility of
 * such damages.
 *
 * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
 * research, please cite the appropriate papers when you publish your
 * work.  Good starting points are:
 *                                                                      
 * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).             
 * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
 * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
 * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
 * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
 */

#ifndef UTILS_ACCUMULATOR_HPP
#define UTILS_ACCUMULATOR_HPP

#include <cmath>
#include <cassert>
#include "math/Vector3.hpp"

namespace OpenMD {


  class BaseAccumulator {
  public:
    virtual void clear() = 0;
    /**
     * get the number of accumulated values
     */
    virtual size_t count()  {
      return Count_;
    }
    virtual ~BaseAccumulator() {};
  protected:
    size_t Count_;

  };   


  /** 
   * Basic Accumulator class for numbers. 
   */  
  class Accumulator : public BaseAccumulator {    

    typedef RealType ElementType;
    typedef RealType ResultType;

  public:
    
    Accumulator() : BaseAccumulator() {
      this->clear();
    }
    
    ~Accumulator() {};

    /**
     * Accumulate another value
     */
    virtual void add(ElementType const& val) {
      Count_++;
      Avg_  += (val       - Avg_ ) / Count_;
      Avg2_ += (val * val - Avg2_) / Count_;
      Val_   = val;
      if (Count_ <= 1) {
        Max_ = val;
        Min_ = val;
      } else {
        Max_ = val > Max_ ? val : Max_;
        Min_ = val < Min_ ? val : Min_;
      }
    }
    
    /**
     * reset the Accumulator to the empty state
     */
    void clear() {
      Count_ = 0;
      Avg_   = 0;
      Avg2_  = 0;
      Val_   = 0;
    }
    

    /**
     * return the most recently added value
     */
    void getLastValue(ElementType &ret)  {
      ret = Val_;
      return;
    }    

    /**
     * compute the Mean
     */
    void getAverage(ResultType &ret)  {
      assert(Count_ != 0);
      ret = Avg_;
      return;
    }

    /**
     * compute the Variance
     */
    void getVariance(ResultType &ret)  {
      assert(Count_ != 0);
      ret = (Avg2_ - Avg_  * Avg_);
      return;
    }
    
    /**
     * compute error of average value
     */
    void getStdDev(ResultType &ret)  {
      assert(Count_ != 0);
      RealType var;
      this->getVariance(var);
      ret = sqrt(var);
      return;
    }

    /**
     * return the largest value
     */
    void getMax(ElementType &ret)  {
      assert(Count_ != 0);
      ret = Max_;
      return;
    }

    /**
     * return the smallest value
     */
    void getMin(ElementType &ret)  {
      assert(Count_ != 0);
      ret = Max_;
      return;
    }

    /**
     * return the 95% confidence interval:
     *
     * That is returns c, such that we have 95% confidence that the
     * true mean is within 2c of the Average (x):
     *
     *   x - c <= true mean <= x + c
     *
     */
    void get95percentConfidenceInterval(ResultType &ret) {
      assert(Count_ != 0);
      RealType sd;
      this->getStdDev(sd);
      ret = 1.960 * sd / sqrt(RealType(Count_));
      return;
    }

  private:
    ElementType Val_;
    ResultType Avg_;
    ResultType Avg2_;
    ElementType Min_;
    ElementType Max_;
  };

  class VectorAccumulator : public BaseAccumulator {
    
    typedef Vector3d ElementType;
    typedef Vector3d ResultType;
    
  public:
    VectorAccumulator() : BaseAccumulator() {
      this->clear();
    }

    /**
     * Accumulate another value
     */
    void add(ElementType const& val) {
      Count_++;
      RealType len(0.0);
      for (unsigned int i =0; i < 3; i++) {
        Avg_[i]  += (val[i]       - Avg_[i] ) / Count_;
        Avg2_[i] += (val[i] * val[i] - Avg2_[i]) / Count_;
        Val_[i]   = val[i];
        len += val[i]*val[i];
      }
      len = sqrt(len);
      AvgLen_  += (len       - AvgLen_ ) / Count_;
      AvgLen2_ += (len * len - AvgLen2_) / Count_;

      if (Count_ <= 1) {
        Max_ = len;
        Min_ = len;
      } else {
        Max_ = len > Max_ ? len : Max_;
        Min_ = len < Min_ ? len : Min_;
      }
    }

    /**
     * reset the Accumulator to the empty state
     */
    void clear() {
      Count_ = 0;
      Avg_ = V3Zero;
      Avg2_ = V3Zero;
      Val_ = V3Zero;
      AvgLen_   = 0;
      AvgLen2_  = 0;
    }
    
    /**
     * return the most recently added value
     */
    void getLastValue(ElementType &ret) {
      ret = Val_;
      return;
    }
    
    /**
     * compute the Mean
     */
    void getAverage(ResultType &ret) {
      assert(Count_ != 0);
      ret = Avg_;
      return;
    }
    
    /**
     * compute the Variance
     */
    void getVariance(ResultType &ret) {
      assert(Count_ != 0);
      for (unsigned int i =0; i < 3; i++) {
        ret[i] = (Avg2_[i] - Avg_[i]  * Avg_[i]);
      }
      return;
    }
    
    /**
     * compute error of average value
     */
    void getStdDev(ResultType &ret) {
      assert(Count_ != 0);
      ResultType var;
      this->getVariance(var);
      ret[0] = sqrt(var[0]);
      ret[1] = sqrt(var[1]);
      ret[2] = sqrt(var[2]);
      return;
    }

    /**
     * return the 95% confidence interval:
     *
     * That is returns c, such that we have 95% confidence that the
     * true mean is within 2c of the Average (x):
     *
     *   x - c <= true mean <= x + c
     *
     */
    void get95percentConfidenceInterval(ResultType &ret) {
      assert(Count_ != 0);
      ResultType sd;
      this->getStdDev(sd);
      ret[0] = 1.960 * sd[0] / sqrt(RealType(Count_));
      ret[1] = 1.960 * sd[1] / sqrt(RealType(Count_));
      ret[2] = 1.960 * sd[2] / sqrt(RealType(Count_));
      return;
    }

    /**
     * return the largest length
     */
    void getMaxLength(RealType &ret) {
      assert(Count_ != 0);
      ret = Max_;
      return;
    }

    /**
     * return the smallest length
     */
    void getMinLength(RealType &ret) {
      assert(Count_ != 0);
      ret = Min_;
      return;
    }

    /**
     * return the largest length
     */
    void getAverageLength(RealType &ret) {
      assert(Count_ != 0);
      ret = AvgLen_;
      return;
    }

    /**
     * compute the Variance of the length
     */
    void getLengthVariance(RealType &ret) {
      assert(Count_ != 0);      
      ret= (AvgLen2_ - AvgLen_ * AvgLen_);
      return;
    }
    
    /**
     * compute error of average value
     */
    void getLengthStdDev(RealType &ret) {
      assert(Count_ != 0);
      RealType var;
      this->getLengthVariance(var);
      ret = sqrt(var);
      return;
    }

    /**
     * return the 95% confidence interval:
     *
     * That is returns c, such that we have 95% confidence that the
     * true mean is within 2c of the Average (x):
     *
     *   x - c <= true mean <= x + c
     *
     */
    void getLength95percentConfidenceInterval(ResultType &ret) {
      assert(Count_ != 0);
      RealType sd;
      this->getLengthStdDev(sd);
      ret = 1.960 * sd / sqrt(RealType(Count_));
      return;
    }


  private:
    ResultType Val_;
    ResultType Avg_;
    ResultType Avg2_;
    RealType AvgLen_;
    RealType AvgLen2_;
    RealType Min_;
    RealType Max_;

  };

  class MatrixAccumulator : public BaseAccumulator {
    
    typedef Mat3x3d ElementType;
    typedef Mat3x3d ResultType;
    
  public:
    MatrixAccumulator() : BaseAccumulator() {
      this->clear();
    }

    /**
     * Accumulate another value
     */
    void add(ElementType const& val) {
      Count_++;
      for (unsigned int i = 0; i < 3; i++) {
        for (unsigned int j = 0; j < 3; j++) {          
          Avg_(i,j)  += (val(i,j)       - Avg_(i,j) ) / Count_;
          Avg2_(i,j) += (val(i,j) * val(i,j) - Avg2_(i,j)) / Count_;
          Val_(i,j)   = val(i,j);
        }
      }
    }

    /**
     * reset the Accumulator to the empty state
     */
    void clear() {
      Count_ = 0;
      Avg_ *= 0.0;
      Avg2_ *= 0.0;
      Val_ *= 0.0;
    }
    
    /**
     * return the most recently added value
     */
    void getLastValue(ElementType &ret) {
      ret = Val_;
      return;
    }
    
    /**
     * compute the Mean
     */
    void getAverage(ResultType &ret) {
      assert(Count_ != 0);
      ret = Avg_;
      return;
    }

    /**
     * compute the Variance
     */
    void getVariance(ResultType &ret) {
      assert(Count_ != 0);
      for (unsigned int i = 0; i < 3; i++) {
        for (unsigned int j = 0; j < 3; j++) {          
          ret(i,j) = (Avg2_(i,j) - Avg_(i,j)  * Avg_(i,j));
        }
      }
      return;
    }
    
    /**
     * compute error of average value
     */
    void getStdDev(ResultType &ret) {
      assert(Count_ != 0);
      Mat3x3d var;
      this->getVariance(var);
      for (unsigned int i = 0; i < 3; i++) {
        for (unsigned int j = 0; j < 3; j++) {
          ret(i,j) = sqrt(var(i,j));  
        }
      }
      return;
    }
      
    /**
     * return the 95% confidence interval:
     *
     * That is returns c, such that we have 95% confidence that the
     * true mean is within 2c of the Average (x):
     *
     *   x - c <= true mean <= x + c
     *
     */
    void get95percentConfidenceInterval(ResultType &ret) {
      assert(Count_ != 0);
      Mat3x3d sd;
      this->getStdDev(sd);
      for (unsigned int i = 0; i < 3; i++) {
        for (unsigned int j = 0; j < 3; j++) {
          ret(i,j) = 1.960 * sd(i,j) / sqrt(RealType(Count_));
        }
      }
      return;
    }

  private:
    ElementType Val_;
    ResultType Avg_;
    ResultType Avg2_;
  };


} 

#endif
Revision:	2071
Committed:	Sat Mar 7 21:41:51 2015 UTC (10 years, 1 month ago) by gezelter
File size:	11339 byte(s)
Log Message:	Reducing the number of warnings when using g++ to compile.
#	Content
1	/*
2	* Copyright (c) 2012 The University of Notre Dame. All Rights Reserved.
3	*
4	* The University of Notre Dame grants you ("Licensee") a
5	* non-exclusive, royalty free, license to use, modify and
6	* redistribute this software in source and binary code form, provided
7	* that the following conditions are met:
8	*
9	* 1. Redistributions of source code must retain the above copyright
10	* notice, this list of conditions and the following disclaimer.
11	*
12	* 2. Redistributions in binary form must reproduce the above copyright
13	* notice, this list of conditions and the following disclaimer in the
14	* documentation and/or other materials provided with the
15	* distribution.
16	*
17	* This software is provided "AS IS," without a warranty of any
18	* kind. All express or implied conditions, representations and
19	* warranties, including any implied warranty of merchantability,
20	* fitness for a particular purpose or non-infringement, are hereby
21	* excluded. The University of Notre Dame and its licensors shall not
22	* be liable for any damages suffered by licensee as a result of
23	* using, modifying or distributing the software or its
24	* derivatives. In no event will the University of Notre Dame or its
25	* licensors be liable for any lost revenue, profit or data, or for
26	* direct, indirect, special, consequential, incidental or punitive
27	* damages, however caused and regardless of the theory of liability,
28	* arising out of the use of or inability to use software, even if the
29	* University of Notre Dame has been advised of the possibility of
30	* such damages.
31	*
32	* SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your
33	* research, please cite the appropriate papers when you publish your
34	* work. Good starting points are:
35	*
36	* [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
37	* [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
38	* [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).
39	* [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010).
40	* [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41	*/
42
43	#ifndef UTILS_ACCUMULATOR_HPP
44	#define UTILS_ACCUMULATOR_HPP
45
46	#include <cmath>
47	#include <cassert>
48	#include "math/Vector3.hpp"
49
50	namespace OpenMD {
51
52
53	class BaseAccumulator {
54	public:
55	virtual void clear() = 0;
56	/**
57	* get the number of accumulated values
58	*/
59	virtual size_t count() {
60	return Count_;
61	}
62	virtual ~BaseAccumulator() {};
63	protected:
64	size_t Count_;
65
66	};
67
68
69
70	/**
71	* Basic Accumulator class for numbers.
72	*/
73	class Accumulator : public BaseAccumulator {
74
75	typedef RealType ElementType;
76	typedef RealType ResultType;
77
78	public:
79
80	Accumulator() : BaseAccumulator() {
81	this->clear();
82	}
83
84	~Accumulator() {};
85
86	/**
87	* Accumulate another value
88	*/
89	virtual void add(ElementType const& val) {
90	Count_++;
91	Avg_ += (val - Avg_ ) / Count_;
92	Avg2_ += (val * val - Avg2_) / Count_;
93	Val_ = val;
94	if (Count_ <= 1) {
95	Max_ = val;
96	Min_ = val;
97	} else {
98	Max_ = val > Max_ ? val : Max_;
99	Min_ = val < Min_ ? val : Min_;
100	}
101	}
102
103	/**
104	* reset the Accumulator to the empty state
105	*/
106	void clear() {
107	Count_ = 0;
108	Avg_ = 0;
109	Avg2_ = 0;
110	Val_ = 0;
111	}
112
113
114	/**
115	* return the most recently added value
116	*/
117	void getLastValue(ElementType &ret) {
118	ret = Val_;
119	return;
120	}
121
122	/**
123	* compute the Mean
124	*/
125	void getAverage(ResultType &ret) {
126	assert(Count_ != 0);
127	ret = Avg_;
128	return;
129	}
130
131	/**
132	* compute the Variance
133	*/
134	void getVariance(ResultType &ret) {
135	assert(Count_ != 0);
136	ret = (Avg2_ - Avg_ * Avg_);
137	return;
138	}
139
140	/**
141	* compute error of average value
142	*/
143	void getStdDev(ResultType &ret) {
144	assert(Count_ != 0);
145	RealType var;
146	this->getVariance(var);
147	ret = sqrt(var);
148	return;
149	}
150
151	/**
152	* return the largest value
153	*/
154	void getMax(ElementType &ret) {
155	assert(Count_ != 0);
156	ret = Max_;
157	return;
158	}
159
160	/**
161	* return the smallest value
162	*/
163	void getMin(ElementType &ret) {
164	assert(Count_ != 0);
165	ret = Max_;
166	return;
167	}
168
169	/**
170	* return the 95% confidence interval:
171	*
172	* That is returns c, such that we have 95% confidence that the
173	* true mean is within 2c of the Average (x):
174	*
175	* x - c <= true mean <= x + c
176	*
177	*/
178	void get95percentConfidenceInterval(ResultType &ret) {
179	assert(Count_ != 0);
180	RealType sd;
181	this->getStdDev(sd);
182	ret = 1.960 * sd / sqrt(RealType(Count_));
183	return;
184	}
185
186	private:
187	ElementType Val_;
188	ResultType Avg_;
189	ResultType Avg2_;
190	ElementType Min_;
191	ElementType Max_;
192	};
193
194	class VectorAccumulator : public BaseAccumulator {
195
196	typedef Vector3d ElementType;
197	typedef Vector3d ResultType;
198
199	public:
200	VectorAccumulator() : BaseAccumulator() {
201	this->clear();
202	}
203
204	/**
205	* Accumulate another value
206	*/
207	void add(ElementType const& val) {
208	Count_++;
209	RealType len(0.0);
210	for (unsigned int i =0; i < 3; i++) {
211	Avg_[i] += (val[i] - Avg_[i] ) / Count_;
212	Avg2_[i] += (val[i] * val[i] - Avg2_[i]) / Count_;
213	Val_[i] = val[i];
214	len += val[i]*val[i];
215	}
216	len = sqrt(len);
217	AvgLen_ += (len - AvgLen_ ) / Count_;
218	AvgLen2_ += (len * len - AvgLen2_) / Count_;
219
220	if (Count_ <= 1) {
221	Max_ = len;
222	Min_ = len;
223	} else {
224	Max_ = len > Max_ ? len : Max_;
225	Min_ = len < Min_ ? len : Min_;
226	}
227	}
228
229	/**
230	* reset the Accumulator to the empty state
231	*/
232	void clear() {
233	Count_ = 0;
234	Avg_ = V3Zero;
235	Avg2_ = V3Zero;
236	Val_ = V3Zero;
237	AvgLen_ = 0;
238	AvgLen2_ = 0;
239	}
240
241	/**
242	* return the most recently added value
243	*/
244	void getLastValue(ElementType &ret) {
245	ret = Val_;
246	return;
247	}
248
249	/**
250	* compute the Mean
251	*/
252	void getAverage(ResultType &ret) {
253	assert(Count_ != 0);
254	ret = Avg_;
255	return;
256	}
257
258	/**
259	* compute the Variance
260	*/
261	void getVariance(ResultType &ret) {
262	assert(Count_ != 0);
263	for (unsigned int i =0; i < 3; i++) {
264	ret[i] = (Avg2_[i] - Avg_[i] * Avg_[i]);
265	}
266	return;
267	}
268
269	/**
270	* compute error of average value
271	*/
272	void getStdDev(ResultType &ret) {
273	assert(Count_ != 0);
274	ResultType var;
275	this->getVariance(var);
276	ret[0] = sqrt(var[0]);
277	ret[1] = sqrt(var[1]);
278	ret[2] = sqrt(var[2]);
279	return;
280	}
281
282	/**
283	* return the 95% confidence interval:
284	*
285	* That is returns c, such that we have 95% confidence that the
286	* true mean is within 2c of the Average (x):
287	*
288	* x - c <= true mean <= x + c
289	*
290	*/
291	void get95percentConfidenceInterval(ResultType &ret) {
292	assert(Count_ != 0);
293	ResultType sd;
294	this->getStdDev(sd);
295	ret[0] = 1.960 * sd[0] / sqrt(RealType(Count_));
296	ret[1] = 1.960 * sd[1] / sqrt(RealType(Count_));
297	ret[2] = 1.960 * sd[2] / sqrt(RealType(Count_));
298	return;
299	}
300
301	/**
302	* return the largest length
303	*/
304	void getMaxLength(RealType &ret) {
305	assert(Count_ != 0);
306	ret = Max_;
307	return;
308	}
309
310	/**
311	* return the smallest length
312	*/
313	void getMinLength(RealType &ret) {
314	assert(Count_ != 0);
315	ret = Min_;
316	return;
317	}
318
319	/**
320	* return the largest length
321	*/
322	void getAverageLength(RealType &ret) {
323	assert(Count_ != 0);
324	ret = AvgLen_;
325	return;
326	}
327
328	/**
329	* compute the Variance of the length
330	*/
331	void getLengthVariance(RealType &ret) {
332	assert(Count_ != 0);
333	ret= (AvgLen2_ - AvgLen_ * AvgLen_);
334	return;
335	}
336
337	/**
338	* compute error of average value
339	*/
340	void getLengthStdDev(RealType &ret) {
341	assert(Count_ != 0);
342	RealType var;
343	this->getLengthVariance(var);
344	ret = sqrt(var);
345	return;
346	}
347
348	/**
349	* return the 95% confidence interval:
350	*
351	* That is returns c, such that we have 95% confidence that the
352	* true mean is within 2c of the Average (x):
353	*
354	* x - c <= true mean <= x + c
355	*
356	*/
357	void getLength95percentConfidenceInterval(ResultType &ret) {
358	assert(Count_ != 0);
359	RealType sd;
360	this->getLengthStdDev(sd);
361	ret = 1.960 * sd / sqrt(RealType(Count_));
362	return;
363	}
364
365
366	private:
367	ResultType Val_;
368	ResultType Avg_;
369	ResultType Avg2_;
370	RealType AvgLen_;
371	RealType AvgLen2_;
372	RealType Min_;
373	RealType Max_;
374
375	};
376
377	class MatrixAccumulator : public BaseAccumulator {
378
379	typedef Mat3x3d ElementType;
380	typedef Mat3x3d ResultType;
381
382	public:
383	MatrixAccumulator() : BaseAccumulator() {
384	this->clear();
385	}
386
387	/**
388	* Accumulate another value
389	*/
390	void add(ElementType const& val) {
391	Count_++;
392	for (unsigned int i = 0; i < 3; i++) {
393	for (unsigned int j = 0; j < 3; j++) {
394	Avg_(i,j) += (val(i,j) - Avg_(i,j) ) / Count_;
395	Avg2_(i,j) += (val(i,j) * val(i,j) - Avg2_(i,j)) / Count_;
396	Val_(i,j) = val(i,j);
397	}
398	}
399	}
400
401	/**
402	* reset the Accumulator to the empty state
403	*/
404	void clear() {
405	Count_ = 0;
406	Avg_ *= 0.0;
407	Avg2_ *= 0.0;
408	Val_ *= 0.0;
409	}
410
411	/**
412	* return the most recently added value
413	*/
414	void getLastValue(ElementType &ret) {
415	ret = Val_;
416	return;
417	}
418
419	/**
420	* compute the Mean
421	*/
422	void getAverage(ResultType &ret) {
423	assert(Count_ != 0);
424	ret = Avg_;
425	return;
426	}
427
428	/**
429	* compute the Variance
430	*/
431	void getVariance(ResultType &ret) {
432	assert(Count_ != 0);
433	for (unsigned int i = 0; i < 3; i++) {
434	for (unsigned int j = 0; j < 3; j++) {
435	ret(i,j) = (Avg2_(i,j) - Avg_(i,j) * Avg_(i,j));
436	}
437	}
438	return;
439	}
440
441	/**
442	* compute error of average value
443	*/
444	void getStdDev(ResultType &ret) {
445	assert(Count_ != 0);
446	Mat3x3d var;
447	this->getVariance(var);
448	for (unsigned int i = 0; i < 3; i++) {
449	for (unsigned int j = 0; j < 3; j++) {
450	ret(i,j) = sqrt(var(i,j));
451	}
452	}
453	return;
454	}
455
456	/**
457	* return the 95% confidence interval:
458	*
459	* That is returns c, such that we have 95% confidence that the
460	* true mean is within 2c of the Average (x):
461	*
462	* x - c <= true mean <= x + c
463	*
464	*/
465	void get95percentConfidenceInterval(ResultType &ret) {
466	assert(Count_ != 0);
467	Mat3x3d sd;
468	this->getStdDev(sd);
469	for (unsigned int i = 0; i < 3; i++) {
470	for (unsigned int j = 0; j < 3; j++) {
471	ret(i,j) = 1.960 * sd(i,j) / sqrt(RealType(Count_));
472	}
473	}
474	return;
475	}
476
477	private:
478	ElementType Val_;
479	ResultType Avg_;
480	ResultType Avg2_;
481	};
482
483
484	}
485
486	#endif