Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 template<typename eT>
00023 inline
00024 eT
00025 op_var::direct_var(const eT* const X, const u32 n_elem, const u32 norm_type)
00026 {
00027 arma_extra_debug_sigprint();
00028
00029 eT acc1 = eT(0);
00030
00031 for(u32 i=0; i<n_elem; ++i)
00032 {
00033 acc1 += X[i];
00034 }
00035
00036 const eT div_val = (n_elem > 0) ? eT(n_elem) : eT(1);
00037 acc1 /= div_val;
00038
00039 eT acc2 = eT(0);
00040 eT acc3 = eT(0);
00041
00042 for(u32 i=0; i<n_elem; ++i)
00043 {
00044 const eT tmp = acc1 - X[i];
00045
00046 acc2 += tmp*tmp;
00047 acc3 += tmp;
00048 }
00049
00050
00051 const eT norm_val = (norm_type == 0) ? ( (n_elem > 1) ? eT(n_elem-1) : eT(1) ) : eT(n_elem);
00052 const eT var_val = (acc2 - acc3*acc3/div_val) / norm_val;
00053
00054 return var_val;
00055 }
00056
00057
00058
00059
00060 template<typename T>
00061 inline
00062 T
00063 op_var::direct_var(const std::complex<T>* const X, const u32 n_elem, const u32 norm_type)
00064 {
00065 arma_extra_debug_sigprint();
00066
00067 typedef typename std::complex<T> eT;
00068
00069 eT acc1 = eT(0);
00070
00071 for(u32 i=0; i<n_elem; ++i)
00072 {
00073 acc1 += X[i];
00074 }
00075
00076 const T div_val = (n_elem > 0) ? T(n_elem) : T(1);
00077 acc1 /= div_val;
00078
00079 T acc2 = T(0);
00080 eT acc3 = eT(0);
00081
00082 for(u32 i=0; i<n_elem; ++i)
00083 {
00084 const eT tmp = acc1 - X[i];
00085
00086 acc2 += std::norm(tmp);
00087 acc3 += tmp;
00088 }
00089
00090 const T norm_val = (norm_type == 0) ? ( (n_elem > 1) ? T(n_elem-1) : T(1) ) : T(n_elem);
00091 const T var_val = (acc2 - std::norm(acc3)/div_val) / norm_val;
00092
00093 return var_val;
00094 }
00095
00096
00097
00098
00099 template<typename eT>
00100 inline
00101 typename get_pod_type<eT>::result
00102 op_var::direct_var(const subview_row<eT>& X, const u32 norm_type)
00103 {
00104 arma_extra_debug_sigprint();
00105
00106 const u32 n_elem = X.n_elem;
00107
00108 podarray<eT> tmp(n_elem);
00109
00110 eT* tmp_mem = tmp.memptr();
00111
00112 for(u32 i=0; i<n_elem; ++i)
00113 {
00114 tmp_mem[i] = X[i];
00115 }
00116
00117 return op_var::direct_var(tmp_mem, n_elem, norm_type);
00118 }
00119
00120
00121
00122
00123 template<typename eT>
00124 inline
00125 typename get_pod_type<eT>::result
00126 op_var::direct_var(const subview_col<eT>& X, const u32 norm_type)
00127 {
00128 arma_extra_debug_sigprint();
00129
00130 return op_var::direct_var(X.colptr(0), X.n_elem, norm_type);
00131 }
00132
00133
00134
00135
00136 template<typename eT>
00137 inline
00138 typename get_pod_type<eT>::result
00139 op_var::direct_var(const diagview<eT>& X, const u32 norm_type)
00140 {
00141 arma_extra_debug_sigprint();
00142
00143 const u32 n_elem = X.n_elem;
00144
00145 podarray<eT> tmp(n_elem);
00146
00147 eT* tmp_mem = tmp.memptr();
00148
00149 for(u32 i=0; i<n_elem; ++i)
00150 {
00151 tmp_mem[i] = X[i];
00152 }
00153
00154 return op_var::direct_var(tmp_mem, n_elem, norm_type);
00155 }
00156
00157
00158
00159
00160
00161
00162
00163 template<typename eT>
00164 inline
00165 void
00166 op_var::apply(Mat< typename get_pod_type<eT>::result >& out, const Mat<eT>& X, const u32 norm_type, const u32 dim)
00167 {
00168 arma_extra_debug_sigprint();
00169
00170 arma_debug_check( (X.n_elem == 0), "var(): given matrix has no elements" );
00171
00172 arma_debug_check( (norm_type > 1), "var(): incorrect usage. norm_type must be 0 or 1");
00173 arma_debug_check( (dim > 1), "var(): incorrect usage. dim must be 0 or 1" );
00174
00175
00176 if(dim == 0)
00177 {
00178 arma_extra_debug_print("op_var::apply(), dim = 0");
00179
00180 out.set_size(1, X.n_cols);
00181
00182 for(u32 col=0; col<X.n_cols; ++col)
00183 {
00184 out[col] = op_var::direct_var( X.colptr(col), X.n_rows, norm_type );
00185 }
00186 }
00187 else
00188 if(dim == 1)
00189 {
00190 arma_extra_debug_print("op_var::apply(), dim = 1");
00191
00192 const u32 n_rows = X.n_rows;
00193 const u32 n_cols = X.n_cols;
00194
00195 out.set_size(n_rows, 1);
00196
00197 podarray<eT> tmp(n_cols);
00198
00199 eT* tmp_mem = tmp.memptr();
00200
00201 for(u32 row=0; row<n_rows; ++row)
00202 {
00203 for(u32 col=0; col<n_cols; ++col)
00204 {
00205 tmp_mem[col] = X.at(row,col);
00206 }
00207
00208 out[row] = op_var::direct_var(tmp_mem, n_cols, norm_type);
00209 }
00210
00211 }
00212
00213 }
00214
00215
00216
00217