Skip to main content

InvHpiFit

Namespace: INVLIB  ·  Library: Inverse Library

Python equivalent

mne.chpi in MNE-Python.

#include <inv/inv_hpi_fit.h>

class INVLIB::InvHpiFit

Drives one HPI fit: per-coil magnetic-dipole localisation, frequency matching against the digitised coil layout and Procrustes solve for the dewar-to-head transform.

The class is reusable across consecutive measurement blocks and caches the sensor geometry to avoid recomputing it on every call.

Drives one HPI fit (per-coil dipole localisation, coil ordering, dewar-to-head transform).

src/examples/ex_inv_hpi/main.cpp
InvHpiDataUpdater updater(info); // good MEG channels, their sensors and the digitised coil positions
updater.prepareDataAndProjectors(segment, MatrixXd::Identity(info->nchan, info->nchan));
// The coil frequencies are known but not their order; the first fit with bOrderFrequencies finds it.
const InvHpiModelParameters model(QVector<int>({154, 158, 161, 166}), info->sfreq, info->linefreq, true);
InvHpiFit hpi(updater.getSensors());
HpiFitResult fit;
hpi.fit(updater.getProjectedData(), updater.getProjectors(), model, updater.getHpiDigitizer(), true, fit);

Public Methods​

InvHpiFit()​

Default constructor.

InvHpiFit(sensorSet)​

Constructs the HPI from a InvSensorSet.

Parameters:

  • sensorSet : const InvSensorSet & The MEG sensorSet used for the hpi fitting.

checkForUpdate(sensorSet)​

Checks if InvSensorSet has changed and updates member InvSensorSet.

Parameters:

  • sensorSet : const InvSensorSet & The MEG sensorSet used for the hpi fitting.

fit(matProjectedData, matProjectors, hpiModelParameters, matCoilsHead, hpiFitResult)​

Perform one single HPI fit.

Parameters:

  • matProjectedData : const Eigen::MatrixXd & Data to estimate the HPI positions from. Projectars should be already applied.

  • matProjectors : const Eigen::MatrixXd & The projectors to apply.

  • hpiModelParameters : const InvHpiModelParameters & The model parameters to use for the Hpi Fitting, especially to compute the coil amplitudes.

  • matCoilsHead : const Eigen::MatrixXd & The hpi coil locations in head space.

  • hpiFitResult : HpiFitResult & The fitting results.

fit(matProjectedData, matProjectors, hpiModelParameters, matCoilsHead, bOrderFrequencies, hpiFitResult)​

Static Methods​

storeHeadPosition(fTime, matTransDevHead, matPosition, vecGoF, vecError)​

Store results from dev_Head_t as quaternions in position matrix.

The format is the same as you get from Neuromag's MaxFilter.

Parameters:

  • fTime : float The corresponding time in the measurement for the fit.

  • matTransDevHead : const Eigen::MatrixXf & The device->head transformation matrix.

  • matPosition : Eigen::MatrixXd & The matrix to store the results.

  • vecGoF : const Eigen::VectorXd & The goodness of fit per coil.

  • vecError : const QVector< double > & The Hpi estimation Error per coil.

compareTransformation(mDevHeadT, mDevHeadDest, fThreshRot, fThreshTrans)​

Compares two 4x4 transformation matrices by evaluating their rotation and translation differences against given thresholds.

Parameters:

  • mDevHeadT : const Eigen::MatrixX4f & The first transformation matrix.

  • mDevHeadDest : const Eigen::MatrixX4f & The second (destination) transformation matrix.

  • fThreshRot : const float & The rotation threshold in degrees.

  • fThreshTrans : const float & The translation threshold in meters.

Returns:

  • bool — True if the movement exceeds either threshold.

Example​

Used in examples​

ex_hpiFit​

Example for cHPI fitting on raw data with SSP. The result is written to a .txt file for comparison with MaxFilter's .pos file.

src/examples/ex_hpiFit/main.cpp, line 182
InvHpiDataUpdater hpiDataUpdater = InvHpiDataUpdater(pFiffInfo);
InvHpiFit HPI = InvHpiFit(hpiDataUpdater.getSensors());

ex_dsp_rt​

DSP real-time workers: averaging, covariance, noise spectrum, connectivity, HPI fitting and inverse operators.

src/examples/ex_dsp_rt/main.cpp, line 297
hpiWorker.doWork(updater.getProjectedData(), updater.getProjectors(), model, updater.getHpiDigitizer());
HpiFitResult direct;
InvHpiFit(updater.getSensors()).fit(updater.getProjectedData(), updater.getProjectors(), model, updater.getHpiDigitizer(), direct);
ok &= expect(workerFit.GoF.size() == 4 && workerFit.GoF.minCoeff() > 0.98 && workerFit.devHeadTrans.trans.isApprox(direct.devHeadTrans.trans, 1e-6f),
             QString("RtHpiWorker fits 4 coils (GoF >= %1) like InvHpiFit").arg(workerFit.GoF.size() ? workerFit.GoF.minCoeff() : 0.0, 0, 'f', 4));

ex_inv_hpi​

Inverse library: continuous head position (HPI) fitting.

src/examples/ex_inv_hpi/main.cpp, line 107
updater.prepareDataAndProjectors(segment, MatrixXd::Identity(info->nchan, info->nchan));
// The coil frequencies are known but not their order; the first fit with bOrderFrequencies finds it.
const InvHpiModelParameters model(QVector<int>({154, 158, 161, 166}), info->sfreq, info->linefreq, true);
InvHpiFit hpi(updater.getSensors());
HpiFitResult fit;
hpi.fit(updater.getProjectedData(), updater.getProjectors(), model, updater.getHpiDigitizer(), true, fit);
// MaxFilter head position at t = 0: translation (-0.814, -11.046, 53.572) mm. One 200 ms window, so allow 1 mm
// (test_hpiFit_integration holds the mean over the whole recording to 0.3 mm).
MatrixXd maxfilter;
IOUtils::read_eigen_matrix(maxfilter, data + "/Result/ref_hpiFit_pos.txt");
const Vector3d translation = fit.devHeadTrans.trans.block<3, 1>(0, 3).cast<double>();
const Vector3d reference = maxfilter.block<1, 3>(0, 4).transpose();
ok &= expect(fit.hpiFreqs == QVector<int>({166, 154, 161, 158}) && fit.GoF.minCoeff() > 0.98 && (translation - reference).norm() < 0.001,

Authors of this file​