MNE-CPP  0.1.9
A Framework for Electrophysiology
mne.cpp
Go to the documentation of this file.
1 //=============================================================================================================
38 //=============================================================================================================
39 // INCLUDES
40 //=============================================================================================================
41 
42 #include "mne.h"
43 #include <fiff/fiff.h>
44 #include <iostream>
45 
46 //=============================================================================================================
47 // QT INCLUDES
48 //=============================================================================================================
49 
50 #include <QFile>
51 #include <QDebug>
52 
53 //=============================================================================================================
54 // USED NAMESPACES
55 //=============================================================================================================
56 
57 using namespace MNELIB;
58 using namespace Eigen;
59 using namespace FIFFLIB;
60 
61 //=============================================================================================================
62 // DEFINE MEMBER METHODS
63 //=============================================================================================================
64 
65 bool MNE::read_events(QString t_sEventName,
66  QString t_fileRawName,
67  MatrixXi& events)
68 {
69  QFile t_EventFile;
70  qint32 p;
71  bool status = false;
72 
73  if (t_sEventName.isEmpty()) {
74  p = t_fileRawName.indexOf(".fif");
75  if (p > 0) {
76  t_sEventName = t_fileRawName.replace(p, 4, "-eve.fif");
77  } else {
78  printf("Raw file name does not end properly\n");
79  return 0;
80  }
81 
82  t_EventFile.setFileName(t_sEventName);
83  if(!MNE::read_events_from_fif(t_EventFile, events)) {
84  printf("Error while read events.\n");
85  return false;
86  }
87  printf("Events read from %s\n",t_sEventName.toUtf8().constData());
88  } else {
89  // Binary file
90  if (t_sEventName.contains(".fif")) {
91  t_EventFile.setFileName(t_sEventName);
92  if(!MNE::read_events_from_fif(t_EventFile, events)) {
93  printf("Error while read events.\n");
94  return false;
95  }
96  printf("Binary event file %s read\n",t_sEventName.toUtf8().constData());
97  } else if(t_sEventName.contains(".eve")){
98 
99  } else {
100  // Text file
101  printf("Text file %s is not supported jet.\n",t_sEventName.toUtf8().constData());
102 // try
103 // events = load(eventname);
104 // catch
105 // error(me,mne_omit_first_line(lasterr));
106 // end
107 // if size(events,1) < 1
108 // error(me,'No data in the event file');
109 // end
110 // //
111 // // Convert time to samples if sample number is negative
112 // //
113 // for p = 1:size(events,1)
114 // if events(p,1) < 0
115 // events(p,1) = events(p,2)*raw.info.sfreq;
116 // end
117 // end
118 // //
119 // // Select the columns of interest (convert to integers)
120 // //
121 // events = int32(events(:,[1 3 4]));
122 // //
123 // // New format?
124 // //
125 // if events(1,2) == 0 && events(1,3) == 0
126 // fprintf(1,'The text event file %s is in the new format\n',eventname);
127 // if events(1,1) ~= raw.first_samp
128 // error(me,'This new format event file is not compatible with the raw data');
129 // end
130 // else
131 // fprintf(1,'The text event file %s is in the old format\n',eventname);
132 // //
133 // // Offset with first sample
134 // //
135 // events(:,1) = events(:,1) + raw.first_samp;
136 // end
137  }
138  }
139 
140  return true;
141 }
142 
143 //=============================================================================================================
144 
145 bool MNE::read_events_from_fif(QIODevice &p_IODevice,
146  MatrixXi& eventlist)
147 {
148  //
149  // Open file
150  //
151  FiffStream::SPtr t_pStream(new FiffStream(&p_IODevice));
152 
153  if(!t_pStream->open()) {
154  return false;
155  }
156 
157  //
158  // Find the desired block
159  //
160  QList<FiffDirNode::SPtr> events = t_pStream->dirtree()->dir_tree_find(FIFFB_MNE_EVENTS);
161 
162  if (events.size() == 0)
163  {
164  printf("Could not find event data\n");
165  return false;
166  }
167 
168  qint32 k, nelem;
169  fiff_int_t kind, pos;
170  FiffTag::SPtr t_pTag;
171  quint32* serial_eventlist_uint = NULL;
172  qint32* serial_eventlist_int = NULL;
173 
174  for(k = 0; k < events[0]->nent(); ++k)
175  {
176  kind = events[0]->dir[k]->kind;
177  pos = events[0]->dir[k]->pos;
178  if (kind == FIFF_MNE_EVENT_LIST)
179  {
180  t_pStream->read_tag(t_pTag,pos);
181  if(t_pTag->type == FIFFT_UINT)
182  {
183  serial_eventlist_uint = t_pTag->toUnsignedInt();
184  nelem = t_pTag->size()/4;
185  }
186 
187  if(t_pTag->type == FIFFT_INT)
188  {
189  serial_eventlist_int = t_pTag->toInt();
190  nelem = t_pTag->size()/4;
191  }
192 
193  break;
194  }
195  }
196 
197  if(serial_eventlist_uint == NULL && serial_eventlist_int == NULL)
198  {
199  printf("Could not find any events\n");
200  return false;
201  }
202  else
203  {
204  eventlist.resize(nelem/3,3);
205  if(serial_eventlist_uint != NULL)
206  {
207  for(k = 0; k < nelem/3; ++k)
208  {
209  eventlist(k,0) = serial_eventlist_uint[k*3];
210  eventlist(k,1) = serial_eventlist_uint[k*3+1];
211  eventlist(k,2) = serial_eventlist_uint[k*3+2];
212  }
213  }
214 
215  if(serial_eventlist_int != NULL)
216  {
217  for(k = 0; k < nelem/3; ++k)
218  {
219  eventlist(k,0) = serial_eventlist_int[k*3];
220  eventlist(k,1) = serial_eventlist_int[k*3+1];
221  eventlist(k,2) = serial_eventlist_int[k*3+2];
222  }
223  }
224  }
225 
226  return true;
227 }
228 
229 //=============================================================================================================
230 
231 void MNE::setup_compensators(FiffRawData& raw,
232  fiff_int_t dest_comp,
233  bool keep_comp)
234 {
235  // Set up projection
236  if (raw.info.projs.size() == 0) {
237  printf("No projector specified for these data\n");
238  } else {
239  // Activate the projection items
240  for (qint32 k = 0; k < raw.info.projs.size(); ++k) {
241  raw.info.projs[k].active = true;
242  }
243 
244  printf("%d projection items activated\n",raw.info.projs.size());
245  // Create the projector
246 // fiff_int_t nproj = MNE::make_projector_info(raw.info, raw.proj); Using the member function instead
247  fiff_int_t nproj = raw.info.make_projector(raw.proj);
248 
249  if (nproj == 0) {
250  printf("The projection vectors do not apply to these channels\n");
251  } else {
252  printf("Created an SSP operator (subspace dimension = %d)\n",nproj);
253  }
254  }
255 
256  // Set up the CTF compensator
257 // qint32 current_comp = MNE::get_current_comp(raw.info);
258  qint32 current_comp = raw.info.get_current_comp();
259  if (current_comp > 0)
260  printf("Current compensation grade : %d\n",current_comp);
261 
262  if (keep_comp)
263  dest_comp = current_comp;
264 
265  if (current_comp != dest_comp)
266  {
267  qDebug() << "This part needs to be debugged";
268  if(MNE::make_compensator(raw.info, current_comp, dest_comp, raw.comp))
269  {
270 // raw.info.chs = MNE::set_current_comp(raw.info.chs,dest_comp);
271  raw.info.set_current_comp(dest_comp);
272  printf("Appropriate compensator added to change to grade %d.\n",dest_comp);
273  }
274  else
275  {
276  printf("Could not make the compensator\n");
277  return;
278  }
279  }
280 }
281 
282 //=============================================================================================================
283 
284 bool MNE::read_events_from_ascii(QIODevice &p_IODevice,
285  Eigen::MatrixXi& eventlist)
286 {
287  if (!p_IODevice.open(QIODevice::ReadOnly | QIODevice::Text)){
288  return false;
289  }
290  QTextStream textStream(&p_IODevice);
291 
292  QList<int> simpleList;
293 
294  while(!textStream.atEnd()){
295  int iSample;
296  textStream >> iSample;
297  simpleList.append(iSample);
298  textStream.readLine();
299  qDebug() << "Added event:" << iSample;
300  }
301 
302  eventlist.resize(simpleList.size(), 1);
303 
304  for(int i = 0; i < simpleList.size(); i++){
305  eventlist(i,0) = simpleList[i];
306  }
307  return true;
308 }
FIFF raw measurement data.
Definition: fiff_raw_data.h:78
static bool make_compensator(const FIFFLIB::FiffInfo &info, FIFFLIB::fiff_int_t from, FIFFLIB::fiff_int_t to, FIFFLIB::FiffCtfComp &ctf_comp, bool exclude_comp_chs=false)
Definition: mne.h:207
#define FIFF_MNE_EVENT_LIST
QList< FiffProj > projs
Definition: fiff_info.h:268
QSharedPointer< FiffTag > SPtr
Definition: fiff_tag.h:152
static bool read_events_from_ascii(QIODevice &p_IODevice, Eigen::MatrixXi &eventlist)
Definition: mne.cpp:284
FIFF File I/O routines.
Definition: fiff_stream.h:104
qint32 get_current_comp()
Definition: fiff_info.cpp:137
qint32 make_projector(Eigen::MatrixXd &proj) const
Definition: fiff_info.h:278
Eigen::MatrixXd proj
QSharedPointer< FiffStream > SPtr
Definition: fiff_stream.h:107
MNE class declaration, which provides static wrapper functions to stay consistent with mne matlab too...
FIFF class declaration, which provides static wrapper functions to stay consistent with mne matlab to...
static bool read_events_from_fif(QIODevice &p_IODevice, Eigen::MatrixXi &eventlist)
Definition: mne.cpp:145
void set_current_comp(fiff_int_t value)
Definition: fiff_info.h:292