v2.0.0
Loading...
Searching...
No Matches
MNELIB::MNEInverseOperator Class Reference

MNE-style inverse operator. More...

#include <mne_inverse_operator.h>

Public Types

typedef QSharedPointer< MNEInverseOperatorSPtr
typedef QSharedPointer< const MNEInverseOperatorConstSPtr

Public Member Functions

 MNEInverseOperator ()
 Constructs an empty inverse operator with invalid sentinel values.
 MNEInverseOperator (QIODevice &p_IODevice)
 Constructs an inverse operator by reading from a FIFF IO device.
 MNEInverseOperator (const FIFFLIB::FiffInfo &info, const MNEForwardSolution &forward, const FIFFLIB::FiffCov &noiseCov, float loose=0.2f, float depth=0.8f, bool fixed=false, bool limit_depth_chs=true)
 Constructs an inverse operator from measurement info, forward solution, and noise covariance.
 MNEInverseOperator (const MNEInverseOperator &other)
 Copy constructor.
 ~MNEInverseOperator ()
 Destructor.
bool assemble_kernel (const FSLIB::FsLabel &label, const QString &method, bool pick_normal, Eigen::MatrixXd &K, Eigen::SparseMatrix< double > &noise_norm, QList< Eigen::VectorXi > &vertno)
 Assemble the inverse kernel matrix.
bool check_ch_names (const FIFFLIB::FiffInfo &info) const
 Verify that inverse-operator channels are present in the measurement info.
Eigen::MatrixXd cluster_kernel (const FSLIB::FsAnnotationSet &annotationSet, qint32 clusterSize, Eigen::MatrixXd &D, const QString &method=QStringLiteral("cityblock")) const
 Cluster the inverse kernel by cortical parcellation.
Eigen::MatrixXd & getKernel ()
 Access the most recently assembled kernel (mutable).
Eigen::MatrixXd getKernel () const
 Access the most recently assembled kernel (const).
bool isFixedOrient () const
 Check whether the inverse operator uses fixed source orientations.
MNEInverseOperator prepare_inverse_operator (qint32 nave, float lambda2, bool dSPM, bool sLORETA=false) const
 Prepare the inverse operator for source estimation.
void write (QIODevice &p_IODevice)
 Write the inverse operator to a FIFF file.
void writeToStream (FIFFLIB::FiffStream *p_pStream)
 Write the inverse operator into an already-open FIFF stream.

Static Public Member Functions

static MNEInverseOperator make_inverse_operator (const FIFFLIB::FiffInfo &info, MNEForwardSolution forward, const FIFFLIB::FiffCov &noiseCov, float loose=0.2f, float depth=0.8f, bool fixed=false, bool limit_depth_chs=true)
 Assemble an inverse operator from a forward solution and noise covariance.
static bool read_inverse_operator (QIODevice &p_IODevice, MNEInverseOperator &inv)
 Read an inverse operator from a FIFF file.

Public Attributes

FIFFLIB::FiffInfoBase info
FIFFLIB::fiff_int_t methods
FIFFLIB::fiff_int_t source_ori
FIFFLIB::fiff_int_t nsource
FIFFLIB::fiff_int_t nchan
FIFFLIB::fiff_int_t coord_frame
Eigen::MatrixXf source_nn
Eigen::VectorXd sing
bool eigen_leads_weighted
FIFFLIB::FiffNamedMatrix::SDPtr eigen_leads
FIFFLIB::FiffNamedMatrix::SDPtr eigen_fields
FIFFLIB::FiffCov::SDPtr noise_cov
FIFFLIB::FiffCov::SDPtr source_cov
FIFFLIB::FiffCov::SDPtr orient_prior
FIFFLIB::FiffCov::SDPtr depth_prior
FIFFLIB::FiffCov::SDPtr fmri_prior
MNESourceSpaces src
FIFFLIB::FiffCoordTrans mri_head_t
FIFFLIB::fiff_int_t nave
QList< FIFFLIB::FiffProjprojs
Eigen::MatrixXd proj
Eigen::MatrixXd whitener
Eigen::VectorXd reginv
Eigen::SparseMatrix< double > noisenorm

Friends

std::ostream & operator<< (std::ostream &out, const MNEInverseOperator &inv)
 Stream-output operator for diagnostic printing.

Detailed Description

MNE-style inverse operator.

Encapsulates the SVD decomposition of the whitened and weighted lead-field matrix together with noise and source covariance, priors, and the source space. Supports reading/writing FIFF files, preparation for application (regularisation, whitening, noise normalisation), and gain-matrix clustering.

Definition at line 137 of file mne_inverse_operator.h.

Member Typedef Documentation

◆ ConstSPtr

Const shared pointer type for MNEInverseOperator.

Definition at line 141 of file mne_inverse_operator.h.

◆ SPtr

Shared pointer type for MNEInverseOperator.

Definition at line 140 of file mne_inverse_operator.h.

Constructor & Destructor Documentation

◆ MNEInverseOperator() [1/4]

MNEInverseOperator::MNEInverseOperator ( )

Constructs an empty inverse operator with invalid sentinel values.

Definition at line 76 of file mne_inverse_operator.cpp.

◆ MNEInverseOperator() [2/4]

MNEInverseOperator::MNEInverseOperator ( QIODevice & p_IODevice)
explicit

Constructs an inverse operator by reading from a FIFF IO device.

Parameters
[in]p_IODeviceIO device (e.g. QFile) containing the inverse operator.

Definition at line 98 of file mne_inverse_operator.cpp.

◆ MNEInverseOperator() [3/4]

MNEInverseOperator::MNEInverseOperator ( const FIFFLIB::FiffInfo & info,
const MNEForwardSolution & forward,
const FIFFLIB::FiffCov & noiseCov,
float loose = 0.2f,
float depth = 0.8f,
bool fixed = false,
bool limit_depth_chs = true )

Constructs an inverse operator from measurement info, forward solution, and noise covariance.

Parameters
[in]infoMeasurement info specifying channels; bad channels in info.bads are excluded.
[in]forwardForward operator.
[in]noiseCovNoise covariance matrix.
[in]looseTangential-component weight in [0, 1] (ignored when fixed is true).
[in]depthDepth-weighting exponent in [0, 1]; 0 disables depth weighting.
[in]fixedIf true, use fixed source orientations normal to the cortical mantle.
[in]limit_depth_chsIf true, restrict depth weighting to gradiometers (or magnetometers/EEG as fallback).

Definition at line 107 of file mne_inverse_operator.cpp.

◆ MNEInverseOperator() [4/4]

MNEInverseOperator::MNEInverseOperator ( const MNEInverseOperator & other)

Copy constructor.

Parameters
[in]otherInverse operator to copy.

Definition at line 122 of file mne_inverse_operator.cpp.

◆ ~MNEInverseOperator()

MNEInverseOperator::~MNEInverseOperator ( )
default

Destructor.

Member Function Documentation

◆ assemble_kernel()

bool MNEInverseOperator::assemble_kernel ( const FSLIB::FsLabel & label,
const QString & method,
bool pick_normal,
Eigen::MatrixXd & K,
Eigen::SparseMatrix< double > & noise_norm,
QList< Eigen::VectorXi > & vertno )

Assemble the inverse kernel matrix.

Applies eigenfield/eigenlead decomposition, regularisation, whitening, projection, and optional noise normalisation to produce the kernel matrix that maps sensor data to source estimates.

Parameters
[in]labelCortical label restricting the source space (empty for all sources).
[in]methodInverse method: "MNE", "dSPM", or "sLORETA".
[in]pick_normalIf true, keep only the surface-normal component (requires free orientation).
[out]KAssembled kernel matrix (nSources x nChannels).
[out]noise_normDiagonal noise-normalisation matrix (empty for MNE).
[out]vertnoVertex numbers per hemisphere.
Returns
true on success.

Definition at line 158 of file mne_inverse_operator.cpp.

◆ check_ch_names()

bool MNEInverseOperator::check_ch_names ( const FIFFLIB::FiffInfo & info) const

Verify that inverse-operator channels are present in the measurement info.

Parameters
[in]infoMeasurement info whose channel names are checked.
Returns
true if all channels match, false otherwise.

Definition at line 311 of file mne_inverse_operator.cpp.

◆ cluster_kernel()

MatrixXd MNEInverseOperator::cluster_kernel ( const FSLIB::FsAnnotationSet & annotationSet,
qint32 clusterSize,
Eigen::MatrixXd & D,
const QString & method = QStringLiteral("cityblock") ) const

Cluster the inverse kernel by cortical parcellation.

Groups source-space points within each parcellation label using KMeans and returns a reduced kernel together with the cluster operator.

Parameters
[in]annotationSetAnnotation set for left and right hemisphere.
[in]clusterSizeMaximum number of sources per cluster.
[out]DCluster operator matrix (nSources x nClusters).
[in]methodDistance metric: "cityblock" or "sqeuclidean".
Returns
Clustered kernel matrix.

Definition at line 356 of file mne_inverse_operator.cpp.

◆ getKernel() [1/2]

Eigen::MatrixXd & MNELIB::MNEInverseOperator::getKernel ( )
inline

Access the most recently assembled kernel (mutable).

Returns
Reference to the kernel matrix.

Definition at line 382 of file mne_inverse_operator.h.

◆ getKernel() [2/2]

Eigen::MatrixXd MNELIB::MNEInverseOperator::getKernel ( ) const
inline

Access the most recently assembled kernel (const).

Returns
Copy of the kernel matrix.

Definition at line 389 of file mne_inverse_operator.h.

◆ isFixedOrient()

bool MNELIB::MNEInverseOperator::isFixedOrient ( ) const
inline

Check whether the inverse operator uses fixed source orientations.

Returns
true if orientation is fixed-normal.

Definition at line 396 of file mne_inverse_operator.h.

◆ make_inverse_operator()

MNEInverseOperator MNEInverseOperator::make_inverse_operator ( const FIFFLIB::FiffInfo & info,
MNEForwardSolution forward,
const FIFFLIB::FiffCov & noiseCov,
float loose = 0.2f,
float depth = 0.8f,
bool fixed = false,
bool limit_depth_chs = true )
static

Assemble an inverse operator from a forward solution and noise covariance.

Performs depth weighting, source-covariance construction, whitening, and SVD decomposition of the weighted lead-field matrix.

Parameters
[in]infoMeasurement info; bad channels in info.bads are excluded.
[in]forwardForward operator (modified internally for fixed orientation if needed).
[in]noiseCovNoise covariance matrix.
[in]looseTangential-component weight in [0, 1].
[in]depthDepth-weighting exponent in [0, 1]; 0 disables.
[in]fixedUse fixed surface-normal orientations.
[in]limit_depth_chsRestrict depth weighting to gradiometers (or fallback).
Returns
Assembled inverse operator.

Definition at line 635 of file mne_inverse_operator.cpp.

◆ prepare_inverse_operator()

MNEInverseOperator MNEInverseOperator::prepare_inverse_operator ( qint32 nave,
float lambda2,
bool dSPM,
bool sLORETA = false ) const

Prepare the inverse operator for source estimation.

Scales covariances by the number of averages, builds the regularised pseudo-inverse, the whitener, the SSP projector, and (optionally) the noise-normalisation factors for dSPM or sLORETA.

Parameters
[in]naveNumber of averages (scales the noise covariance).
[in]lambda2Regularisation parameter (SNR^{-2}).
[in]dSPMCompute dSPM noise-normalisation factors.
[in]sLORETACompute sLORETA noise-normalisation factors.
Returns
A prepared copy of the inverse operator ready for apply_inverse.

Definition at line 884 of file mne_inverse_operator.cpp.

◆ read_inverse_operator()

bool MNEInverseOperator::read_inverse_operator ( QIODevice & p_IODevice,
MNEInverseOperator & inv )
static

Read an inverse operator from a FIFF file.

Parameters
[in]p_IODeviceFIFF IO device (e.g. QFile or QTcpSocket).
[in,out]invReceives the loaded inverse operator.
Returns
true on success, false on failure.

Definition at line 1030 of file mne_inverse_operator.cpp.

◆ write()

void MNEInverseOperator::write ( QIODevice & p_IODevice)

Write the inverse operator to a FIFF file.

Parameters
[in]p_IODeviceIO device to write to.

Definition at line 1263 of file mne_inverse_operator.cpp.

◆ writeToStream()

void MNEInverseOperator::writeToStream ( FIFFLIB::FiffStream * p_pStream)

Write the inverse operator into an already-open FIFF stream.

Parameters
[in]p_pStreamOpen FIFF stream.

Definition at line 1278 of file mne_inverse_operator.cpp.

◆ operator<<

std::ostream & operator<< ( std::ostream & out,
const MNEInverseOperator & inv )
friend

Stream-output operator for diagnostic printing.

Parameters
[in]outOutput stream.
[in]invInverse operator to print.
Returns
The output stream.

Definition at line 403 of file mne_inverse_operator.h.

Member Data Documentation

◆ coord_frame

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::coord_frame

Coordinate frame of the inverse (MRI or head).

Definition at line 354 of file mne_inverse_operator.h.

◆ depth_prior

FIFFLIB::FiffCov::SDPtr MNELIB::MNEInverseOperator::depth_prior

Depth prior (depth-weighting coefficients).

Definition at line 363 of file mne_inverse_operator.h.

◆ eigen_fields

FIFFLIB::FiffNamedMatrix::SDPtr MNELIB::MNEInverseOperator::eigen_fields

Left singular vectors (eigenfields), channels x components.

Definition at line 359 of file mne_inverse_operator.h.

◆ eigen_leads

FIFFLIB::FiffNamedMatrix::SDPtr MNELIB::MNEInverseOperator::eigen_leads

Right singular vectors (eigenleads), sources x components.

Definition at line 358 of file mne_inverse_operator.h.

◆ eigen_leads_weighted

bool MNELIB::MNEInverseOperator::eigen_leads_weighted

True if eigenleads already include R^{0.5} weighting.

Definition at line 357 of file mne_inverse_operator.h.

◆ fmri_prior

FIFFLIB::FiffCov::SDPtr MNELIB::MNEInverseOperator::fmri_prior

fMRI prior (if available).

Definition at line 364 of file mne_inverse_operator.h.

◆ info

FIFFLIB::FiffInfoBase MNELIB::MNEInverseOperator::info

Lightweight measurement info (channel names, types, etc.).

Definition at line 349 of file mne_inverse_operator.h.

◆ methods

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::methods

Modality flag: FIFFV_MNE_MEG, _EEG, or _MEG_EEG.

Definition at line 350 of file mne_inverse_operator.h.

◆ mri_head_t

FIFFLIB::FiffCoordTrans MNELIB::MNEInverseOperator::mri_head_t

MRI-to-head coordinate transformation.

Definition at line 366 of file mne_inverse_operator.h.

◆ nave

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::nave

Number of averages (default 1).

Definition at line 367 of file mne_inverse_operator.h.

◆ nchan

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::nchan

Number of channels (= number of singular values).

Definition at line 353 of file mne_inverse_operator.h.

◆ noise_cov

FIFFLIB::FiffCov::SDPtr MNELIB::MNEInverseOperator::noise_cov

Noise covariance matrix.

Definition at line 360 of file mne_inverse_operator.h.

◆ noisenorm

Eigen::SparseMatrix<double> MNELIB::MNEInverseOperator::noisenorm

Diagonal noise-normalisation matrix (dSPM / sLORETA).

Definition at line 372 of file mne_inverse_operator.h.

◆ nsource

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::nsource

Number of source-space points.

Definition at line 352 of file mne_inverse_operator.h.

◆ orient_prior

FIFFLIB::FiffCov::SDPtr MNELIB::MNEInverseOperator::orient_prior

Orientation prior (loose-constraint weighting).

Definition at line 362 of file mne_inverse_operator.h.

◆ proj

Eigen::MatrixXd MNELIB::MNEInverseOperator::proj

SSP projector matrix applied to data.

Definition at line 369 of file mne_inverse_operator.h.

◆ projs

QList<FIFFLIB::FiffProj> MNELIB::MNEInverseOperator::projs

SSP projectors read from the file.

Definition at line 368 of file mne_inverse_operator.h.

◆ reginv

Eigen::VectorXd MNELIB::MNEInverseOperator::reginv

Diagonal regularised inverse of the singular values.

Definition at line 371 of file mne_inverse_operator.h.

◆ sing

Eigen::VectorXd MNELIB::MNEInverseOperator::sing

Singular values of the whitened lead field.

Definition at line 356 of file mne_inverse_operator.h.

◆ source_cov

FIFFLIB::FiffCov::SDPtr MNELIB::MNEInverseOperator::source_cov

Source covariance matrix (depth + orientation weighting).

Definition at line 361 of file mne_inverse_operator.h.

◆ source_nn

Eigen::MatrixXf MNELIB::MNEInverseOperator::source_nn

Source-normal vectors (nsource x 3, or nsource*3 x 3 for free).

Definition at line 355 of file mne_inverse_operator.h.

◆ source_ori

FIFFLIB::fiff_int_t MNELIB::MNEInverseOperator::source_ori

Source orientation constraint (FIFFV_MNE_FREE_ORI / FIXED_ORI).

Definition at line 351 of file mne_inverse_operator.h.

◆ src

MNESourceSpaces MNELIB::MNEInverseOperator::src

Source space (surfaces and/or volume).

Definition at line 365 of file mne_inverse_operator.h.

◆ whitener

Eigen::MatrixXd MNELIB::MNEInverseOperator::whitener

Whitening matrix (derived from noise covariance).

Definition at line 370 of file mne_inverse_operator.h.


The documentation for this class was generated from the following files: