annotate armadillo-3.900.4/include/armadillo_bits/spop_mean_meat.hpp @ 84:55a047986812 tip

Update library URI so as not to be document-local
author Chris Cannam
date Wed, 22 Apr 2020 14:21:57 +0100
parents 1ec0e2823891
children
rev   line source
Chris@49 1 // Copyright (C) 2012 Ryan Curtin
Chris@49 2 // Copyright (C) 2012 Conrad Sanderson
Chris@49 3 //
Chris@49 4 // This Source Code Form is subject to the terms of the Mozilla Public
Chris@49 5 // License, v. 2.0. If a copy of the MPL was not distributed with this
Chris@49 6 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
Chris@49 7
Chris@49 8
Chris@49 9 //! \addtogroup spop_mean
Chris@49 10 //! @{
Chris@49 11
Chris@49 12
Chris@49 13
Chris@49 14 template<typename T1>
Chris@49 15 inline
Chris@49 16 void
Chris@49 17 spop_mean::apply(SpMat<typename T1::elem_type>& out, const SpOp<T1, spop_mean>& in)
Chris@49 18 {
Chris@49 19 arma_extra_debug_sigprint();
Chris@49 20
Chris@49 21 typedef typename T1::elem_type eT;
Chris@49 22
Chris@49 23 const uword dim = in.aux_uword_a;
Chris@49 24 arma_debug_check((dim > 1), "mean(): incorrect usage. dim must be 0 or 1");
Chris@49 25
Chris@49 26 SpProxy<T1> p(in.m);
Chris@49 27
Chris@49 28 if(p.is_alias(out) == false)
Chris@49 29 {
Chris@49 30 spop_mean::apply_noalias(out, p, dim);
Chris@49 31 }
Chris@49 32 else
Chris@49 33 {
Chris@49 34 SpMat<eT> tmp;
Chris@49 35
Chris@49 36 spop_mean::apply_noalias(tmp, p, dim);
Chris@49 37
Chris@49 38 out.steal_mem(tmp);
Chris@49 39 }
Chris@49 40 }
Chris@49 41
Chris@49 42
Chris@49 43
Chris@49 44 template<typename T1>
Chris@49 45 inline
Chris@49 46 void
Chris@49 47 spop_mean::apply_noalias
Chris@49 48 (
Chris@49 49 SpMat<typename T1::elem_type>& out_ref,
Chris@49 50 const SpProxy<T1>& p,
Chris@49 51 const uword dim
Chris@49 52 )
Chris@49 53 {
Chris@49 54 arma_extra_debug_sigprint();
Chris@49 55
Chris@49 56 typedef typename T1::elem_type eT;
Chris@49 57
Chris@49 58 const uword p_n_rows = p.get_n_rows();
Chris@49 59 const uword p_n_cols = p.get_n_cols();
Chris@49 60
Chris@49 61 if (dim == 0)
Chris@49 62 {
Chris@49 63 arma_extra_debug_print("spop_mean::apply_noalias(), dim = 0");
Chris@49 64
Chris@49 65 out_ref.set_size((p_n_rows > 0) ? 1 : 0, p_n_cols);
Chris@49 66
Chris@49 67 if(p_n_rows > 0)
Chris@49 68 {
Chris@49 69 for(uword col = 0; col < p_n_cols; ++col)
Chris@49 70 {
Chris@49 71 // Do we have to use an iterator or can we use memory directly?
Chris@49 72 if(SpProxy<T1>::must_use_iterator == true)
Chris@49 73 {
Chris@49 74 typename SpProxy<T1>::const_iterator_type it = p.begin_col(col);
Chris@49 75 typename SpProxy<T1>::const_iterator_type end = p.begin_col(col + 1);
Chris@49 76
Chris@49 77 const uword n_zero = p.get_n_rows() - (end.pos() - it.pos());
Chris@49 78
Chris@49 79 out_ref.at(col) = spop_mean::iterator_mean(it, end, n_zero, eT(0));
Chris@49 80 }
Chris@49 81 else
Chris@49 82 {
Chris@49 83 out_ref.at(col) = spop_mean::direct_mean
Chris@49 84 (
Chris@49 85 &p.get_values()[p.get_col_ptrs()[col]],
Chris@49 86 p.get_col_ptrs()[col + 1] - p.get_col_ptrs()[col],
Chris@49 87 p.get_n_rows()
Chris@49 88 );
Chris@49 89 }
Chris@49 90 }
Chris@49 91 }
Chris@49 92 }
Chris@49 93 else if (dim == 1)
Chris@49 94 {
Chris@49 95 arma_extra_debug_print("spop_mean::apply_noalias(), dim = 1");
Chris@49 96
Chris@49 97 out_ref.set_size(p_n_rows, (p_n_cols > 0) ? 1 : 0);
Chris@49 98
Chris@49 99 if(p_n_cols > 0)
Chris@49 100 {
Chris@49 101 for(uword row = 0; row < p_n_rows; ++row)
Chris@49 102 {
Chris@49 103 // We must use an iterator regardless of how it is stored.
Chris@49 104 typename SpProxy<T1>::const_row_iterator_type it = p.begin_row(row);
Chris@49 105 typename SpProxy<T1>::const_row_iterator_type end = p.end_row(row);
Chris@49 106
Chris@49 107 const uword n_zero = p.get_n_cols() - (end.pos() - it.pos());
Chris@49 108
Chris@49 109 out_ref.at(row) = spop_mean::iterator_mean(it, end, n_zero, eT(0));
Chris@49 110 }
Chris@49 111 }
Chris@49 112 }
Chris@49 113 }
Chris@49 114
Chris@49 115
Chris@49 116
Chris@49 117 template<typename eT>
Chris@49 118 inline
Chris@49 119 eT
Chris@49 120 spop_mean::direct_mean
Chris@49 121 (
Chris@49 122 const eT* const X,
Chris@49 123 const uword length,
Chris@49 124 const uword N
Chris@49 125 )
Chris@49 126 {
Chris@49 127 arma_extra_debug_sigprint();
Chris@49 128
Chris@49 129 typedef typename get_pod_type<eT>::result T;
Chris@49 130
Chris@49 131 const eT result = arrayops::accumulate(X, length) / T(N);
Chris@49 132
Chris@49 133 return arma_isfinite(result) ? result : spop_mean::direct_mean_robust(X, length, N);
Chris@49 134 }
Chris@49 135
Chris@49 136
Chris@49 137
Chris@49 138 template<typename eT>
Chris@49 139 inline
Chris@49 140 eT
Chris@49 141 spop_mean::direct_mean_robust
Chris@49 142 (
Chris@49 143 const eT* const X,
Chris@49 144 const uword length,
Chris@49 145 const uword N
Chris@49 146 )
Chris@49 147 {
Chris@49 148 arma_extra_debug_sigprint();
Chris@49 149
Chris@49 150 typedef typename get_pod_type<eT>::result T;
Chris@49 151
Chris@49 152 uword i, j;
Chris@49 153
Chris@49 154 eT r_mean = eT(0);
Chris@49 155
Chris@49 156 const uword diff = (N - length); // number of zeros
Chris@49 157
Chris@49 158 for(i = 0, j = 1; j < length; i += 2, j += 2)
Chris@49 159 {
Chris@49 160 const eT Xi = X[i];
Chris@49 161 const eT Xj = X[j];
Chris@49 162
Chris@49 163 r_mean += (Xi - r_mean) / T(diff + j);
Chris@49 164 r_mean += (Xj - r_mean) / T(diff + j + 1);
Chris@49 165 }
Chris@49 166
Chris@49 167 if(i < length)
Chris@49 168 {
Chris@49 169 const eT Xi = X[i];
Chris@49 170
Chris@49 171 r_mean += (Xi - r_mean) / T(diff + i + 1);
Chris@49 172 }
Chris@49 173
Chris@49 174 return r_mean;
Chris@49 175 }
Chris@49 176
Chris@49 177
Chris@49 178
Chris@49 179 template<typename T1>
Chris@49 180 inline
Chris@49 181 typename T1::elem_type
Chris@49 182 spop_mean::mean_all(const SpBase<typename T1::elem_type, T1>& X)
Chris@49 183 {
Chris@49 184 arma_extra_debug_sigprint();
Chris@49 185
Chris@49 186 SpProxy<T1> p(X.get_ref());
Chris@49 187
Chris@49 188 if (SpProxy<T1>::must_use_iterator == true)
Chris@49 189 {
Chris@49 190 typename SpProxy<T1>::const_iterator_type it = p.begin();
Chris@49 191 typename SpProxy<T1>::const_iterator_type end = p.end();
Chris@49 192
Chris@49 193 return spop_mean::iterator_mean(it, end, p.get_n_elem() - p.get_n_nonzero(), typename T1::elem_type(0));
Chris@49 194 }
Chris@49 195 else // must_use_iterator == false; that is, we can directly access the values array
Chris@49 196 {
Chris@49 197 return spop_mean::direct_mean(p.get_values(), p.get_n_nonzero(), p.get_n_elem());
Chris@49 198 }
Chris@49 199 }
Chris@49 200
Chris@49 201
Chris@49 202
Chris@49 203 template<typename T1, typename eT>
Chris@49 204 inline
Chris@49 205 eT
Chris@49 206 spop_mean::iterator_mean(T1& it, const T1& end, const uword n_zero, const eT junk)
Chris@49 207 {
Chris@49 208 arma_extra_debug_sigprint();
Chris@49 209 arma_ignore(junk);
Chris@49 210
Chris@49 211 typedef typename get_pod_type<eT>::result T;
Chris@49 212
Chris@49 213 eT sum = eT(0);
Chris@49 214
Chris@49 215 T1 backup_it(it); // in case we have to use robust iterator_mean
Chris@49 216
Chris@49 217 const uword it_begin_pos = it.pos();
Chris@49 218
Chris@49 219 while (it != end)
Chris@49 220 {
Chris@49 221 sum += (*it);
Chris@49 222 ++it;
Chris@49 223 }
Chris@49 224
Chris@49 225 const eT result = sum / T(n_zero + (it.pos() - it_begin_pos));
Chris@49 226
Chris@49 227 return arma_isfinite(result) ? result : spop_mean::iterator_mean_robust(backup_it, end, n_zero, eT(0));
Chris@49 228 }
Chris@49 229
Chris@49 230
Chris@49 231
Chris@49 232 template<typename T1, typename eT>
Chris@49 233 inline
Chris@49 234 eT
Chris@49 235 spop_mean::iterator_mean_robust(T1& it, const T1& end, const uword n_zero, const eT junk)
Chris@49 236 {
Chris@49 237 arma_extra_debug_sigprint();
Chris@49 238 arma_ignore(junk);
Chris@49 239
Chris@49 240 typedef typename get_pod_type<eT>::result T;
Chris@49 241
Chris@49 242 eT r_mean = eT(0);
Chris@49 243
Chris@49 244 const uword it_begin_pos = it.pos();
Chris@49 245
Chris@49 246 while (it != end)
Chris@49 247 {
Chris@49 248 r_mean += ((*it - r_mean) / T(n_zero + (it.pos() - it_begin_pos) + 1));
Chris@49 249 ++it;
Chris@49 250 }
Chris@49 251
Chris@49 252 return r_mean;
Chris@49 253 }
Chris@49 254
Chris@49 255
Chris@49 256
Chris@49 257 //! @}