58 throw std::runtime_error(
"Error during fiff setup raw read");
70 throw std::runtime_error(
"Error during fiff setup raw read");
101 cals = RowVectorXd();
108#include <QElapsedTimer>
114 const RowVectorXi& sel,
117 bool projAvailable =
true;
119 if (this->
proj.size() == 0) {
121 projAvailable =
false;
132 from = this->first_samp;
137 qWarning(
"No data in this range %d ... %d = %9.3f ... %9.3f secs...", from, to, (
static_cast<float>(from)) / this->
info.
sfreq, (
static_cast<float>(to)) / this->info.sfreq);
144 qint32 nchan = this->
info.nchan;
148 using T = Eigen::Triplet<double>;
149 std::vector<T> tripletList;
150 tripletList.reserve(nchan);
151 for (i = 0; i < nchan; ++i)
152 tripletList.push_back(T(i, i, this->
cals[i]));
154 SparseMatrix<double> cal(nchan, nchan);
155 cal.setFromTriplets(tripletList.begin(), tripletList.end());
160 if (sel.size() == 0) {
161 data = MatrixXd(nchan, to - from + 1);
163 if (projAvailable || this->
comp.kind != -1) {
165 mult_full = this->
comp.data->data * cal;
166 else if (this->
comp.kind == -1)
167 mult_full = this->
proj * cal;
169 mult_full = this->
proj * this->
comp.data->data * cal;
172 data = MatrixXd(sel.size(), to - from + 1);
175 MatrixXd selVect(sel.size(), nchan);
179 if (!projAvailable && this->
comp.kind == -1) {
181 tripletList.reserve(sel.size());
182 for (i = 0; i < sel.size(); ++i)
183 tripletList.push_back(T(i, i, this->
cals[sel[i]]));
184 cal = SparseMatrix<double>(sel.size(), sel.size());
185 cal.setFromTriplets(tripletList.begin(), tripletList.end());
187 if (!projAvailable) {
188 for (i = 0; i < sel.size(); ++i)
189 selVect.row(i) = this->
comp.data->data.block(sel[i], 0, 1, nchan);
190 mult_full = selVect * cal;
191 }
else if (this->
comp.kind == -1) {
192 for (i = 0; i < sel.size(); ++i)
193 selVect.row(i) = this->
proj.block(sel[i], 0, 1, nchan);
195 mult_full = selVect * cal;
197 for (i = 0; i < sel.size(); ++i)
198 selVect.row(i) = this->
proj.block(sel[i], 0, 1, nchan);
200 mult_full = selVect * this->
comp.data->data * cal;
209 tripletList.reserve(mult_full.rows() * mult_full.cols());
210 for (i = 0; i < mult_full.rows(); ++i)
211 for (k = 0; k < mult_full.cols(); ++k)
212 if (mult_full(i, k) != 0)
213 tripletList.push_back(T(i, k, mult_full(i, k)));
215 SparseMatrix<double> mult(mult_full.rows(), mult_full.cols());
216 if (tripletList.size() > 0)
217 mult.setFromTriplets(tripletList.begin(), tripletList.end());
221 if (!this->
file->device()->isOpen()) {
222 if (!this->
file->device()->open(QIODevice::ReadOnly)) {
223 qWarning(
"Cannot open file %s", this->
info.filename.toUtf8().constData());
230 MatrixXd one, newData, tmp_data;
231 FiffRawDir thisRawDir;
234 for (k = 0; k < this->
rawdir.size(); ++k) {
235 thisRawDir = this->
rawdir[k];
239 if (thisRawDir.
last > from) {
240 if (thisRawDir.
ent->kind == -1) {
247 one.resize(nchan, thisRawDir.
nsamp);
249 one.resize(sel.cols(), thisRawDir.
nsamp);
253 fid->read_tag(t_pTag, thisRawDir.
ent->pos);
258 if (mult.cols() == 0) {
259 if (sel.cols() == 0) {
261 one = cal * (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
263 one = cal * (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
265 one = cal * (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
267 one = cal * (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
269 one = cal * Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
271 qWarning(
"Data Storage Format not known yet [1]!! Type: %d\n", t_pTag->type);
272 this->
file->device()->close();
277 newData.resize(sel.cols(), thisRawDir.
nsamp);
280 tmp_data = (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
282 for (r = 0; r < sel.size(); ++r)
283 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
285 tmp_data = (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
287 for (r = 0; r < sel.size(); ++r)
288 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
290 tmp_data = (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
292 for (r = 0; r < sel.size(); ++r)
293 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
295 tmp_data = (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
297 for (r = 0; r < sel.size(); ++r)
298 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
300 tmp_data = Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
302 for (r = 0; r < sel.size(); ++r)
303 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
305 qWarning(
"Data Storage Format not known yet [2]!! Type: %d\n", t_pTag->type);
306 this->
file->device()->close();
314 one = mult * (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
316 one = mult * (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
318 one = mult * (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
320 one = mult * (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
322 one = mult * Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
324 qWarning(
"Data Storage Format not known yet [3]!! Type: %d\n", t_pTag->type);
325 this->
file->device()->close();
333 if (to >= thisRawDir.
last && from <= thisRawDir.
first) {
338 last_pick = thisRawDir.
nsamp - 1;
341 }
else if (from > thisRawDir.
first) {
342 first_pick = from - thisRawDir.
first;
343 if (to < thisRawDir.
last) {
348 last_pick = thisRawDir.
nsamp + to - thisRawDir.
last - 1;
355 last_pick = thisRawDir.
nsamp - 1;
364 last_pick = to - thisRawDir.
first;
371 picksamp = last_pick - first_pick + 1;
374 qDebug() <<
"first_pick: " << first_pick;
375 qDebug() <<
"last_pick: " << last_pick;
376 qDebug() <<
"picksamp: " << picksamp;
383 data.block(0, dest, data.rows(), picksamp) = one.block(0, first_pick, data.rows(), picksamp);
391 if (thisRawDir.
last >= to) {
397 if (!this->
file->device()->isOpen()) {
398 this->
file->device()->close();
401 times = MatrixXd(1, to - from + 1);
403 for (i = 0; i < times.cols(); ++i)
404 times(0, i) =
static_cast<float>(from + i) / this->
info.sfreq;
413 SparseMatrix<double>& multSegment,
416 const RowVectorXi& sel,
419 bool projAvailable =
true;
421 if (this->
proj.size() == 0) {
423 projAvailable =
false;
439 qWarning(
"No data in this range\n");
446 qint32 nchan = this->
info.nchan;
450 using T = Eigen::Triplet<double>;
451 std::vector<T> tripletList;
452 tripletList.reserve(nchan);
453 for (i = 0; i < nchan; ++i)
454 tripletList.push_back(T(i, i, this->
cals[i]));
456 SparseMatrix<double> cal(nchan, nchan);
457 cal.setFromTriplets(tripletList.begin(), tripletList.end());
462 if (sel.size() == 0) {
463 data = MatrixXd(nchan, to - from + 1);
465 if (projAvailable || this->
comp.kind != -1) {
467 mult_full = this->
comp.data->data * cal;
468 else if (this->
comp.kind == -1)
469 mult_full = this->
proj * cal;
471 mult_full = this->
proj * this->
comp.data->data * cal;
474 data = MatrixXd(sel.size(), to - from + 1);
477 MatrixXd selVect(sel.size(), nchan);
481 if (!projAvailable && this->
comp.kind == -1) {
483 tripletList.reserve(sel.size());
484 for (i = 0; i < sel.size(); ++i)
485 tripletList.push_back(T(i, i, this->
cals[sel[i]]));
486 cal = SparseMatrix<double>(sel.size(), sel.size());
487 cal.setFromTriplets(tripletList.begin(), tripletList.end());
489 if (!projAvailable) {
490 for (i = 0; i < sel.size(); ++i)
491 selVect.row(i) = this->
comp.data->data.block(sel[i], 0, 1, nchan);
492 mult_full = selVect * cal;
493 }
else if (this->
comp.kind == -1) {
494 for (i = 0; i < sel.size(); ++i)
495 selVect.row(i) = this->
proj.block(sel[i], 0, 1, nchan);
497 mult_full = selVect * cal;
499 for (i = 0; i < sel.size(); ++i)
500 selVect.row(i) = this->
proj.block(sel[i], 0, 1, nchan);
502 mult_full = selVect * this->
comp.data->data * cal;
511 tripletList.reserve(mult_full.rows() * mult_full.cols());
512 for (i = 0; i < mult_full.rows(); ++i)
513 for (k = 0; k < mult_full.cols(); ++k)
514 if (mult_full(i, k) != 0)
515 tripletList.push_back(T(i, k, mult_full(i, k)));
517 SparseMatrix<double> mult(mult_full.rows(), mult_full.cols());
518 if (tripletList.size() > 0)
519 mult.setFromTriplets(tripletList.begin(), tripletList.end());
525 if (!this->
file->device()->isOpen()) {
526 if (!this->
file->device()->open(QIODevice::ReadOnly)) {
527 qWarning(
"Cannot open file %s", this->
info.filename.toUtf8().constData());
536 for (k = 0; k < this->
rawdir.size(); ++k) {
537 FiffRawDir thisRawDir = this->
rawdir[k];
541 if (thisRawDir.
last > from) {
542 if (thisRawDir.
ent->kind == -1) {
549 one.resize(nchan, thisRawDir.
nsamp);
551 one.resize(sel.cols(), thisRawDir.
nsamp);
556 fid->read_tag(t_pTag, thisRawDir.
ent->pos);
561 if (mult.cols() == 0) {
562 if (sel.cols() == 0) {
564 one = cal * (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
566 one = cal * (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
568 one = cal * (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
570 one = cal * (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
572 one = cal * Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
574 qWarning(
"Data Storage Format not known yet [1]!! Type: %d\n", t_pTag->type);
575 this->
file->device()->close();
580 MatrixXd newData(sel.cols(), thisRawDir.
nsamp);
583 MatrixXd tmp_data = (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
585 for (r = 0; r < sel.size(); ++r)
586 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
588 MatrixXd tmp_data = (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
590 for (r = 0; r < sel.size(); ++r)
591 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
593 MatrixXd tmp_data = (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
595 for (r = 0; r < sel.size(); ++r)
596 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
598 MatrixXd tmp_data = (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
600 for (r = 0; r < sel.size(); ++r)
601 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
603 MatrixXd tmp_data = Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
605 for (r = 0; r < sel.size(); ++r)
606 newData.block(r, 0, 1, thisRawDir.
nsamp) = tmp_data.block(sel[r], 0, 1, thisRawDir.
nsamp);
608 qWarning(
"Data Storage Format not known yet [2]!! Type: %d\n", t_pTag->type);
609 this->
file->device()->close();
617 one = mult * (Map<MatrixDau16>(t_pTag->toDauPack16(), nchan, thisRawDir.
nsamp)).cast<double>();
619 one = mult * (Map<MatrixXi>(t_pTag->toInt(), nchan, thisRawDir.
nsamp)).cast<double>();
621 one = mult * (Map<const MatrixXf>(t_pTag->toFloat(), nchan, thisRawDir.
nsamp)).cast<double>();
623 one = mult * (Map<MatrixShort>(t_pTag->toShort(), nchan, thisRawDir.
nsamp)).cast<double>();
625 one = mult * Map<const MatrixXd>(t_pTag->toDouble(), nchan, thisRawDir.
nsamp);
627 qWarning(
"Data Storage Format not known yet [3]!! Type: %d\n", t_pTag->type);
628 this->
file->device()->close();
636 if (to >= thisRawDir.
last && from <= thisRawDir.
first) {
641 last_pick = thisRawDir.
nsamp - 1;
644 }
else if (from > thisRawDir.
first) {
645 first_pick = from - thisRawDir.
first;
646 if (to < thisRawDir.
last) {
651 last_pick = thisRawDir.
nsamp + to - thisRawDir.
last - 1;
658 last_pick = thisRawDir.
nsamp - 1;
667 last_pick = to - thisRawDir.
first;
674 picksamp = last_pick - first_pick + 1;
677 qDebug() <<
"first_pick: " << first_pick;
678 qDebug() <<
"last_pick: " << last_pick;
679 qDebug() <<
"picksamp: " << picksamp;
686 data.block(0, dest, data.rows(), picksamp) = one.block(0, first_pick, data.rows(), picksamp);
694 if (thisRawDir.
last >= to) {
700 if (mult.cols() == 0)
705 if (!this->
file->device()->isOpen()) {
706 this->
file->device()->close();
709 times = MatrixXd(1, to - from + 1);
711 for (i = 0; i < times.cols(); ++i)
712 times(0, i) =
static_cast<float>(from + i) / this->
info.sfreq;
723 const RowVectorXi& sel)
const
728 from = floor(
static_cast<double>(from) * this->
info.sfreq);
729 to = ceil(
static_cast<double>(to) * this->
info.sfreq);
739 const RowVectorXi& picks,
747 int firstSamp = (from >= 0) ? from :
first_samp;
748 int lastSamp = (to >= 0) ? to :
last_samp;
750 if (firstSamp > lastSamp) {
751 qWarning() <<
"[FiffRawData::save] Invalid sample range.";
760 outInfo.
sfreq =
info.sfreq /
static_cast<float>(decim);
767 qWarning() <<
"[FiffRawData::save] Cannot start writing raw file.";
772 int firstOut = firstSamp / decim;
776 const int blockSize = 2000;
777 int blockSamples = decim * blockSize;
779 for (
int samp = firstSamp; samp <= lastSamp; samp += blockSamples) {
780 int nsamp = qMin(blockSamples, lastSamp - samp + 1);
785 qWarning() <<
"[FiffRawData::save] Error reading data at sample" << samp;
786 pStream->finish_writing_raw();
792 int nOut = (nsamp + decim - 1) / decim;
793 MatrixXd decimData(segData.rows(), nOut);
794 for (
int s = 0, idx = 0; s < nsamp && idx < nOut; s += decim, ++idx) {
795 decimData.col(idx) = segData.col(s);
800 pStream->write_raw_buffer(segData, calsOut);
803 pStream->finish_writing_raw();
805 qInfo() <<
"[FiffRawData::save] Saved raw data from sample" << firstSamp
806 <<
"to" << lastSamp <<
"(decim=" << decim <<
")";
FIFF tag: the 16-byte tag header (kind, type, size, next) plus its decoded payload.
FIFF binary tag-stream layer: wraps a QIODevice to read and write FIFF tags, directories,...
#define FIFF_FIRST_SAMPLE
Stim-channel event list (sample, previous value, new value triples) with FIFF read/write helpers.
FIFF continuous raw recording: FiffInfo plus a directory of FIFF_DATA_BUFFER tags for random-access s...
FIFF file I/O, in-memory data structures and high-level readers/writers.
Full FIFF measurement info: per-channel descriptors, sampling and filter setup, projectors,...
bool read_raw_segment(Eigen::MatrixXd &data, Eigen::MatrixXd ×, fiff_int_t from=-1, fiff_int_t to=-1, const Eigen::RowVectorXi &sel=defaultRowVectorXi, bool do_debug=false) const
bool save(QIODevice &p_IODevice, const Eigen::RowVectorXi &picks=Eigen::RowVectorXi(), int decim=1, int from=-1, int to=-1) const
bool read_raw_segment_times(Eigen::MatrixXd &data, Eigen::MatrixXd ×, float from, float to, const Eigen::RowVectorXi &sel=defaultRowVectorXi) const
QList< FiffRawDir > rawdir
QSharedPointer< FiffStream > SPtr
static bool setup_read_raw(QIODevice &p_IODevice, FiffRawData &data, bool allow_maxshield=true, bool is_littleEndian=false)
static FiffStream::SPtr start_writing_raw(QIODevice &p_IODevice, const FiffInfo &info, Eigen::RowVectorXd &cals, Eigen::MatrixXi sel=defaultMatrixXi, bool bResetRange=true)
std::unique_ptr< FiffTag > UPtr