MNE-CPP  0.1.9
A Framework for Electrophysiology
rtcov.cpp
Go to the documentation of this file.
1 //=============================================================================================================
37 //=============================================================================================================
38 // INCLUDES
39 //=============================================================================================================
40 
41 #include "rtcov.h"
42 
43 
44 //=============================================================================================================
45 // QT INCLUDES
46 //=============================================================================================================
47 
48 #include <QDebug>
49 #include <QtConcurrent>
50 
51 //=============================================================================================================
52 // USED NAMESPACES
53 //=============================================================================================================
54 
55 using namespace RTPROCESSINGLIB;
56 using namespace FIFFLIB;
57 using namespace Eigen;
58 
59 //=============================================================================================================
60 // DEFINE MEMBER METHODS
61 //=============================================================================================================
62 
63 RtCov::RtCov(QSharedPointer<FIFFLIB::FiffInfo> pFiffInfo)
64 : m_fiffInfo(*pFiffInfo)
65 , m_iSamples(0)
66 {
67 }
68 
69 //=============================================================================================================
70 
71 FiffCov RtCov::estimateCovariance(const Eigen::MatrixXd& matData,
72  int iNewMaxSamples)
73 {
74  if(m_fiffInfo.chs.isEmpty()) {
75  qWarning() << "[RtCov::estimateCovariance] FiffInfo was not set. Returning empty covariance estimation.";
76  return FiffCov();
77  }
78 
79  m_lData.append(matData);
80  m_iSamples += matData.cols();
81 
82  if(m_iSamples < iNewMaxSamples) {
83  return FiffCov();
84  }
85 
86  QFuture<RtCovComputeResult> result = QtConcurrent::mappedReduced(m_lData,
87  compute,
88  reduce);
89 
90  result.waitForFinished();
91 
92  RtCovComputeResult finalResult = result.result();
93 
94  //Final computation
95  FiffCov computedCov;
96  computedCov.data = finalResult.matData;
97 
98  QStringList exclude;
99  for(int i = 0; i<m_fiffInfo.chs.size(); i++) {
100  if(m_fiffInfo.chs.at(i).kind != FIFFV_MEG_CH &&
101  m_fiffInfo.chs.at(i).kind != FIFFV_EEG_CH) {
102  exclude << m_fiffInfo.chs.at(i).ch_name;
103  }
104  }
105  bool doProj = true;
106 
107  if(m_iSamples > 0) {
108  finalResult.mu /= (float)m_iSamples;
109  computedCov.data.array() -= m_iSamples * (finalResult.mu * finalResult.mu.transpose()).array();
110  computedCov.data.array() /= (m_iSamples - 1);
111 
112  computedCov.kind = FIFFV_MNE_NOISE_COV;
113  computedCov.diag = false;
114  computedCov.dim = computedCov.data.rows();
115 
116  //ToDo do picks
117  computedCov.names = m_fiffInfo.ch_names;
118  computedCov.projs = m_fiffInfo.projs;
119  computedCov.bads = m_fiffInfo.bads;
120  computedCov.nfree = m_iSamples;
121 
122  // regularize noise covariance
123  computedCov = computedCov.regularize(m_fiffInfo, 0.05, 0.05, 0.1, doProj, exclude);
124 
125 // qint32 samples = rawSegment.cols();
126 // VectorXf mu = rawSegment.rowwise().sum().array() / (float)samples;
127 
128 // MatrixXf noise_covariance = rawSegment * rawSegment.transpose();// noise_covariance == raw_covariance
129 // noise_covariance.array() -= samples * (mu * mu.transpose()).array();
130 // noise_covariance.array() /= (samples - 1);
131 
132 // std::cout << "Noise Covariance:\n" << noise_covariance.block(0,0,10,10) << std::endl;
133 
134 // printf("%d raw buffer (%d x %d) generated\r\n", count, tmp.rows(), tmp.cols());
135 
136  m_lData.clear();
137  m_iSamples = 0;
138 
139  return computedCov;
140  } else {
141  qWarning() << "[RtCov::estimateCovariance] Number of samples equals zero. Regularization not possible. Returning empty covariance estimation.";
142  return FiffCov();
143  }
144 
145 }
146 
147 //=============================================================================================================
148 
149 RtCovComputeResult RtCov::compute(const MatrixXd &matData)
150 {
151  RtCovComputeResult result;
152  result.mu = matData.rowwise().sum();
153  result.matData = matData * matData.transpose();
154  return result;
155 }
156 
157 //=============================================================================================================
158 
159 void RtCov::reduce(RtCovComputeResult& finalResult, const RtCovComputeResult &tempResult)
160 {
161  if(finalResult.matData.size() == 0 || finalResult.mu.size() == 0) {
162  finalResult.mu = tempResult.mu;
163  finalResult.matData = tempResult.matData;
164  } else {
165  finalResult.mu += tempResult.mu;
166  finalResult.matData += tempResult.matData;
167  }
168 }
#define FIFFV_MNE_NOISE_COV
FIFFLIB::FiffCov estimateCovariance(const Eigen::MatrixXd &matData, int iNewMaxSamples)
Definition: rtcov.cpp:71
RtCov class declaration.
FiffCov regularize(const FiffInfo &p_info, double p_fMag=0.1, double p_fGrad=0.1, double p_fEeg=0.1, bool p_bProj=true, QStringList p_exclude=defaultQStringList) const
Definition: fiff_cov.cpp:323
QStringList bads
Definition: fiff_cov.h:197
fiff_int_t dim
Definition: fiff_cov.h:193
Eigen::MatrixXd data
Definition: fiff_cov.h:195
fiff_int_t nfree
Definition: fiff_cov.h:198
QList< FiffProj > projs
Definition: fiff_cov.h:196
fiff_int_t kind
Definition: fiff_cov.h:190
static void reduce(RtCovComputeResult &finalResult, const RtCovComputeResult &tempResult)
Definition: rtcov.cpp:159
static RtCovComputeResult compute(const Eigen::MatrixXd &matData)
Definition: rtcov.cpp:149
QStringList names
Definition: fiff_cov.h:194
covariance data
Definition: fiff_cov.h:77