$search
00001 // Copyright (C) 2009-2011 NICTA (www.nicta.com.au) 00002 // Copyright (C) 2009-2011 Conrad Sanderson 00003 // Copyright (C) 2009-2010 Dimitrios Bouzas 00004 // 00005 // This file is part of the Armadillo C++ library. 00006 // It is provided without any warranty of fitness 00007 // for any purpose. You can redistribute this file 00008 // and/or modify it under the terms of the GNU 00009 // Lesser General Public License (LGPL) as published 00010 // by the Free Software Foundation, either version 3 00011 // of the License or (at your option) any later version. 00012 // (see http://www.opensource.org/licenses for more info) 00013 00014 00015 00018 00019 00020 00021 template<typename eT> 00022 inline 00023 void 00024 op_cov::direct_cov(Mat<eT>& out, const Mat<eT>& A, const uword norm_type) 00025 { 00026 arma_extra_debug_sigprint(); 00027 00028 if(A.is_vec()) 00029 { 00030 if(A.n_rows == 1) 00031 { 00032 out = var(trans(A), norm_type); 00033 } 00034 else 00035 { 00036 out = var(A, norm_type); 00037 } 00038 } 00039 else 00040 { 00041 const uword N = A.n_rows; 00042 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N); 00043 00044 const Row<eT> acc = sum(A); 00045 00046 out = trans(A) * A; 00047 out -= (trans(acc) * acc)/eT(N); 00048 out /= norm_val; 00049 } 00050 } 00051 00052 00053 00054 template<typename T> 00055 inline 00056 void 00057 op_cov::direct_cov(Mat< std::complex<T> >& out, const Mat< std::complex<T> >& A, const uword norm_type) 00058 { 00059 arma_extra_debug_sigprint(); 00060 00061 typedef typename std::complex<T> eT; 00062 00063 if(A.is_vec()) 00064 { 00065 if(A.n_rows == 1) 00066 { 00067 const Mat<T> tmp_mat = var(trans(A), norm_type); 00068 out.set_size(1,1); 00069 out[0] = tmp_mat[0]; 00070 } 00071 else 00072 { 00073 const Mat<T> tmp_mat = var(A, norm_type); 00074 out.set_size(1,1); 00075 out[0] = tmp_mat[0]; 00076 } 00077 } 00078 else 00079 { 00080 const uword N = A.n_rows; 00081 const eT norm_val = (norm_type == 0) ? ( (N > 1) ? eT(N-1) : eT(1) ) : eT(N); 00082 00083 const Row<eT> acc = sum(A); 00084 00085 out = trans(A) * A; // out = strans(conj(A)) * A; 00086 out -= (trans(acc) * acc)/eT(N); // out -= (strans(conj(acc)) * acc)/eT(N); 00087 out /= norm_val; 00088 } 00089 } 00090 00091 00092 00093 template<typename T1> 00094 inline 00095 void 00096 op_cov::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_cov>& in) 00097 { 00098 arma_extra_debug_sigprint(); 00099 00100 typedef typename T1::elem_type eT; 00101 00102 const unwrap_check<T1> tmp(in.m, out); 00103 const Mat<eT>& A = tmp.M; 00104 00105 const uword norm_type = in.aux_uword_a; 00106 00107 op_cov::direct_cov(out, A, norm_type); 00108 } 00109 00110 00111