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
00024 template<typename eT>
00025 inline
00026 void
00027 op_cor::direct_cor(Mat<eT>& out, const Mat<eT>& A, const u32 norm_type)
00028 {
00029 arma_extra_debug_sigprint();
00030
00031 if(A.is_vec())
00032 {
00033 out.set_size(1,1);
00034 out[0] = eT(1);
00035 }
00036 else
00037 {
00038 const u32 N = A.n_rows;
00039 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00040
00041 const Row<eT> acc = sum(A);
00042 const Row<eT> sd = stddev(A);
00043
00044 out = (trans(A) * A);
00045 out -= (trans(acc) * acc)/eT(N);
00046 out /= norm_val;
00047 out /= trans(sd) * sd;
00048 }
00049 }
00050
00051
00052
00053 template<typename T>
00054 inline
00055 void
00056 op_cor::direct_cor(Mat< std::complex<T> >& out, const Mat< std::complex<T> >& A, const u32 norm_type)
00057 {
00058 arma_extra_debug_sigprint();
00059
00060 typedef typename std::complex<T> eT;
00061
00062 if(A.is_vec())
00063 {
00064 out.set_size(1,1);
00065 out[0] = eT(1);
00066 }
00067 else
00068 {
00069 const u32 N = A.n_rows;
00070 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N);
00071
00072 const Row<eT> acc = sum(A);
00073 const Row<T> sd = stddev(A);
00074
00075 out = trans(conj(A)) * A;
00076 out -= (trans(conj(acc)) * acc)/eT(N);
00077 out /= norm_val;
00078
00079
00080 out /= conv_to< Mat<eT> >::from(trans(sd) * sd);
00081 }
00082 }
00083
00084
00085
00086 template<typename T1>
00087 inline
00088 void
00089 op_cor::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_cor>& in)
00090 {
00091 arma_extra_debug_sigprint();
00092
00093 typedef typename T1::elem_type eT;
00094
00095 const unwrap_check<T1> tmp(in.m, out);
00096 const Mat<eT>& A = tmp.M;
00097
00098 const u32 norm_type = in.aux_u32_a;
00099
00100 op_cor::direct_cor(out, A, norm_type);
00101 }
00102
00103
00104
00105