glue_cov_meat.hpp
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
00023 template<typename eT>
00024 inline
00025 void
00026 glue_cov::direct_cov(Mat<eT>& out, const Mat<eT>& A, const Mat<eT>& B, const u32 norm_type)
00027 {
00028 arma_extra_debug_sigprint();
00029
00030 if(A.is_vec() && B.is_vec())
00031 {
00032 arma_debug_check( (A.n_elem != B.n_elem), "cov(): the number of elements in A and B must match" );
00033
00034 const eT* A_ptr = A.memptr();
00035 const eT* B_ptr = B.memptr();
00036
00037 eT A_acc = eT(0);
00038 eT B_acc = eT(0);
00039 eT out_acc = eT(0);
00040
00041 const u32 N = A.n_elem;
00042
00043 for(u32 i=0; i<N; ++i)
00044 {
00045 const eT A_tmp = A_ptr[i];
00046 const eT B_tmp = B_ptr[i];
00047
00048 A_acc += A_tmp;
00049 B_acc += B_tmp;
00050
00051 out_acc += A_tmp * B_tmp;
00052 }
00053
00054 out_acc -= (A_acc * B_acc)/eT(N);
00055
00056 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00057
00058 out.set_size(1,1);
00059 out[0] = out_acc/norm_val;
00060 }
00061 else
00062 {
00063 arma_debug_assert_same_size(A, B, "cov()");
00064
00065 const u32 N = A.n_rows;
00066 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00067
00068 out = trans(A) * B;
00069 out -= (trans(sum(A)) * sum(B))/eT(N);
00070 out /= norm_val;
00071 }
00072 }
00073
00074
00075
00076 template<typename T>
00077 inline
00078 void
00079 glue_cov::direct_cov(Mat< std::complex<T> >& out, const Mat< std::complex<T> >& A, const Mat< std::complex<T> >& B, const u32 norm_type)
00080 {
00081 arma_extra_debug_sigprint();
00082
00083 typedef typename std::complex<T> eT;
00084
00085 if(A.is_vec() && B.is_vec())
00086 {
00087 arma_debug_check( (A.n_elem != B.n_elem), "cov(): the number of elements in A and B must match" );
00088
00089 const eT* A_ptr = A.memptr();
00090 const eT* B_ptr = B.memptr();
00091
00092 eT A_acc = eT(0);
00093 eT B_acc = eT(0);
00094 eT out_acc = eT(0);
00095
00096 const u32 N = A.n_elem;
00097
00098 for(u32 i=0; i<N; ++i)
00099 {
00100 const eT A_tmp = A_ptr[i];
00101 const eT B_tmp = B_ptr[i];
00102
00103 A_acc += A_tmp;
00104 B_acc += B_tmp;
00105
00106 out_acc += std::conj(A_tmp) * B_tmp;
00107 }
00108
00109 out_acc -= (std::conj(A_acc) * B_acc)/eT(N);
00110
00111 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00112
00113 out.set_size(1,1);
00114 out[0] = out_acc/norm_val;
00115 }
00116 else
00117 {
00118 arma_debug_assert_same_size(A, B, "cov()");
00119
00120 const u32 N = A.n_rows;
00121 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00122
00123 out = trans(conj(A)) * B;
00124 out -= (trans(conj(sum(A))) * sum(B))/eT(N);
00125 out /= norm_val;
00126 }
00127 }
00128
00129
00130
00131 template<typename T1, typename T2>
00132 inline
00133 void
00134 glue_cov::apply(Mat<typename T1::elem_type>& out, const Glue<T1,T2,glue_cov>& X)
00135 {
00136 arma_extra_debug_sigprint();
00137
00138 typedef typename T1::elem_type eT;
00139
00140 const unwrap_check<T1> A_tmp(X.A, out);
00141 const unwrap_check<T2> B_tmp(X.B, out);
00142
00143 const Mat<eT>& A = A_tmp.M;
00144 const Mat<eT>& B = B_tmp.M;
00145
00146 const u32 norm_type = X.aux_u32;
00147
00148 if(&A != &B)
00149 {
00150 glue_cov::direct_cov(out, A, B, norm_type);
00151 }
00152 else
00153 {
00154 op_cov::direct_cov(out, A, norm_type);
00155 }
00156
00157 }
00158
00159
00160
00161