MNE-CPP  beta 1.0
rtcov.cpp
Go to the documentation of this file.
1 //=============================================================================================================
38 //*************************************************************************************************************
39 //=============================================================================================================
40 // INCLUDES
41 //=============================================================================================================
42 
43 #include "rtcov.h"
44 
45 #include <iostream>
46 #include <fiff/fiff_cov.h>
47 
48 
49 //*************************************************************************************************************
50 //=============================================================================================================
51 // QT INCLUDES
52 //=============================================================================================================
53 
54 #include <QDebug>
55 
56 
57 //*************************************************************************************************************
58 //=============================================================================================================
59 // USED NAMESPACES
60 //=============================================================================================================
61 
62 using namespace RTINVLIB;
63 using namespace FIFFLIB;
64 
65 
66 //*************************************************************************************************************
67 //=============================================================================================================
68 // DEFINE MEMBER METHODS
69 //=============================================================================================================
70 
71 RtCov::RtCov(qint32 p_iMaxSamples, FiffInfo::SPtr p_pFiffInfo, QObject *parent)
72 : QThread(parent)
73 , m_iMaxSamples(p_iMaxSamples)
74 , m_pFiffInfo(p_pFiffInfo)
75 , m_bIsRunning(false)
76 {
77  qRegisterMetaType<FiffCov::SPtr>("FiffCov::SPtr");
78 }
79 
80 
81 //*************************************************************************************************************
82 
84 {
85  if(this->isRunning())
86  stop();
87 }
88 
89 
90 //*************************************************************************************************************
91 
92 void RtCov::append(const MatrixXd &p_DataSegment)
93 {
94 // if(m_pRawMatrixBuffer) // ToDo handle change buffersize
95 
96  if(!m_pRawMatrixBuffer)
97  m_pRawMatrixBuffer = CircularMatrixBuffer<double>::SPtr(new CircularMatrixBuffer<double>(32, p_DataSegment.rows(), p_DataSegment.cols()));
98 
99  m_pRawMatrixBuffer->push(&p_DataSegment);
100 }
101 
102 
103 //*************************************************************************************************************
104 
105 void RtCov::setSamples(qint32 samples)
106 {
107  m_iNewMaxSamples = samples;
108 }
109 
110 
111 //*************************************************************************************************************
112 
114 {
115  //Check if the thread is already or still running. This can happen if the start button is pressed immediately after the stop button was pressed. In this case the stopping process is not finished yet but the start process is initiated.
116  if(this->isRunning())
117  QThread::wait();
118 
119  m_bIsRunning = true;
120  QThread::start();
121 
122  return true;
123 }
124 
125 
126 //*************************************************************************************************************
127 
129 {
130  m_bIsRunning = false;
131 
132  m_pRawMatrixBuffer->releaseFromPop();
133 
134  m_pRawMatrixBuffer->clear();
135 
136  return true;
137 }
138 
139 
140 //*************************************************************************************************************
141 
143 {
144  quint32 n_samples = 0;
145 
146  FiffCov::SPtr cov(new FiffCov());
147  VectorXd mu;
148 
149  while(m_bIsRunning)
150  {
151  if(m_pRawMatrixBuffer)
152  {
153  MatrixXd rawSegment = m_pRawMatrixBuffer->pop();
154 
155  if(n_samples == 0)
156  {
157  mu = rawSegment.rowwise().sum();
158  cov->data = rawSegment * rawSegment.transpose();
159  }
160  else
161  {
162  mu.array() += rawSegment.rowwise().sum().array();
163  cov->data += rawSegment * rawSegment.transpose();
164  }
165  n_samples += rawSegment.cols();
166 
167  if(n_samples > m_iMaxSamples)
168  {
169  mu /= (float)n_samples;
170  cov->data.array() -= n_samples * (mu * mu.transpose()).array();
171  cov->data.array() /= (n_samples - 1);
172 
173  cov->kind = FIFFV_MNE_NOISE_COV;
174  cov->diag = false;
175  cov->dim = cov->data.rows();
176 
177  //ToDo do picks
178  cov->names = m_pFiffInfo->ch_names;
179  cov->projs = m_pFiffInfo->projs;
180  cov->bads = m_pFiffInfo->bads;
181  cov->nfree = n_samples;
182 
183  // regularize noise covariance
184  *cov.data() = cov->regularize(*m_pFiffInfo, 0.05, 0.05, 0.1, true);
185 
186  emit covCalculated(cov);
187 
188  cov = FiffCov::SPtr(new FiffCov());
189  n_samples = 0;
190  }
191 
192 
193 // qint32 samples = rawSegment.cols();
194 // VectorXf mu = rawSegment.rowwise().sum().array() / (float)samples;
195 
196 // MatrixXf noise_covariance = rawSegment * rawSegment.transpose();// noise_covariance == raw_covariance
197 // noise_covariance.array() -= samples * (mu * mu.transpose()).array();
198 // noise_covariance.array() /= (samples - 1);
199 
200 // std::cout << "Noise Covariance:\n" << noise_covariance.block(0,0,10,10) << std::endl;
201 
202 // printf("%d raw buffer (%d x %d) generated\r\n", count, tmp.rows(), tmp.cols());
203 
204  }
205  }
206 }
#define FIFFV_MNE_NOISE_COV
virtual void run()
Definition: rtcov.cpp:142
RtCov(qint32 p_iMaxSamples, FiffInfo::SPtr p_pFiffInfo, QObject *parent=0)
Definition: rtcov.cpp:71
QSharedPointer< FiffCov > SPtr
Definition: fiff_cov.h:97
RtCov class declaration.
QSharedPointer< CircularMatrixBuffer > SPtr
QSharedPointer< FiffInfo > SPtr
Definition: fiff_info.h:99
void covCalculated(FIFFLIB::FiffCov::SPtr p_pCov)
void push(const Matrix< _Tp, Dynamic, Dynamic > *pMatrix)
void setSamples(qint32 samples)
Definition: rtcov.cpp:105
virtual bool stop()
Definition: rtcov.cpp:128
Matrix< _Tp, Dynamic, Dynamic > pop()
covariance data
Definition: fiff_cov.h:94
Definition: fiff.h:98
virtual bool start()
Definition: rtcov.cpp:113
void append(const MatrixXd &p_DataSegment)
Definition: rtcov.cpp:92
Definition: rtave.h:89
bool isRunning()
Definition: rtcov.h:207
FiffCov class declaration.