@@ -52,6 +52,43 @@ class Mean : public Lazy<typename Rcpp::traits::storage_type<RTYPE>::type, Mean<
5252 const VEC_TYPE& object ;
5353};
5454
55+ template <bool NA, typename T>
56+ class Mean <CPLXSXP,NA,T> : public Lazy<Rcomplex, Mean<CPLXSXP,NA,T> > {
57+ public:
58+ typedef typename Rcpp::VectorBase<CPLXSXP,NA,T> VEC_TYPE;
59+
60+ Mean (const VEC_TYPE& object_) : object(object_) {}
61+
62+ Rcomplex get () const {
63+ ComplexVector input = object;
64+ int n = input.size (); // double pass (as in summary.c)
65+ long double s = 0.0 , si = 0.0 ;
66+ for (int i=0 ; i<n; i++) {
67+ Rcomplex z = input[i];
68+ s += z.r ;
69+ si += z.i ;
70+ }
71+ s /= n;
72+ si /= n;
73+ if (R_FINITE ((double )s) && R_FINITE ((double )si)) {
74+ long double t = 0.0 , ti = 0.0 ;
75+ for (int i = 0 ; i < n; i++) {
76+ Rcomplex z = input[i];
77+ t += z.r - s;
78+ ti += z.i - si;
79+ }
80+ s += t/n;
81+ si += ti/n;
82+ }
83+ Rcomplex z;
84+ z.r = s;
85+ z.i = si;
86+ return z;
87+ }
88+ private:
89+ const VEC_TYPE& object ;
90+ };
91+
5592} // sugar
5693
5794template <bool NA, typename T>
@@ -64,6 +101,11 @@ inline sugar::Mean<INTSXP,NA,T> mean(const VectorBase<INTSXP,NA,T>& t) {
64101 return sugar::Mean<INTSXP,NA,T>(t);
65102}
66103
104+ template <bool NA, typename T>
105+ inline sugar::Mean<CPLXSXP,NA,T> mean (const VectorBase<CPLXSXP,NA,T>& t) {
106+ return sugar::Mean<CPLXSXP,NA,T>(t);
107+ }
108+
67109} // Rcpp
68110#endif
69111
0 commit comments