<s> from copy import deepcopy <EOL> from math import sqrt <EOL> import numpy as np <EOL> from scipy import linalg <EOL> from . . io . constants import FIFF <EOL> from . . io . open import fiff_open <EOL> from . . io . tag import find_tag <EOL> from . . io . matrix import ( _read_named_matrix , _transpose_named_matrix , <EOL> write_named_matrix ) <EOL> from . . io . proj import _read_proj , make_projector , _write_proj <EOL> from . . io . proj import _needs_eeg_average_ref_proj <EOL> from . . io . tree import dir_tree_find <EOL> from . . io . write import ( write_int , write_float_matrix , start_file , <EOL> start_block , end_block , end_file , write_float , <EOL> write_coord_trans , write_string ) <EOL> from . . io . pick import channel_type , pick_info , pick_types <EOL> from . . cov import prepare_noise_cov , _read_cov , _write_cov , Covariance <EOL> from . . forward import ( compute_depth_prior , _read_forward_meas_info , <EOL> write_forward_meas_info , is_fixed_orient , <EOL> compute_orient_prior , convert_forward_solution ) <EOL> from . . source_space import ( _read_source_spaces_from_tree , <EOL> find_source_space_hemi , _get_vertno , <EOL> _write_source_spaces_to_fid , label_src_vertno_sel ) <EOL> from . . transforms import _ensure_trans , transform_surface_to <EOL> from . . source_estimate import _make_stc <EOL> from . . utils import check_fname , logger , verbose , warn <EOL> from functools import reduce <EOL> class InverseOperator ( dict ) : <EOL> """<STR_LIT>""" <EOL> def copy ( self ) : <EOL> """<STR_LIT>""" <EOL> return InverseOperator ( deepcopy ( self ) ) <EOL> def __repr__ ( self ) : <EOL> """<STR_LIT>""" <EOL> entr = '<STR_LIT>' <EOL> nchan = len ( pick_types ( self [ '<STR_LIT:info>' ] , meg = True , eeg = False ) ) <EOL> entr += '<STR_LIT>' + '<STR_LIT>' % nchan <EOL> nchan = len ( pick_types ( self [ '<STR_LIT:info>' ] , meg = False , eeg = True ) ) <EOL> entr += '<STR_LIT>' + '<STR_LIT>' % nchan <EOL> if self [ '<STR_LIT:src>' ] [ <NUM_LIT:0> ] [ '<STR_LIT:type>' ] == '<STR_LIT>' : <EOL> entr += ( '<STR_LIT>' <EOL> % self [ '<STR_LIT>' ] ) <EOL> elif self [ '<STR_LIT:src>' ] [ <NUM_LIT:0> ] [ '<STR_LIT:type>' ] == '<STR_LIT>' : <EOL> entr += ( '<STR_LIT>' <EOL> % self [ '<STR_LIT>' ] ) <EOL> elif self [ '<STR_LIT:src>' ] [ <NUM_LIT:0> ] [ '<STR_LIT:type>' ] == '<STR_LIT>' : <EOL> entr += ( '<STR_LIT>' <EOL> % self [ '<STR_LIT>' ] ) <EOL> source_ori = { FIFF . FIFFV_MNE_UNKNOWN_ORI : '<STR_LIT>' , <EOL> FIFF . FIFFV_MNE_FIXED_ORI : '<STR_LIT>' , <EOL> FIFF . FIFFV_MNE_FREE_ORI : '<STR_LIT>' } <EOL> entr += '<STR_LIT>' % source_ori [ self [ '<STR_LIT>' ] ] <EOL> entr += '<STR_LIT:>>' <EOL> return entr <EOL> def _pick_channels_inverse_operator ( ch_names , inv ) : <EOL> """<STR_LIT>""" <EOL> sel = [ ] <EOL> for name in inv [ '<STR_LIT>' ] . ch_names : <EOL> if name in ch_names : <EOL> sel . append ( ch_names . index ( name ) ) <EOL> else : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' % name ) <EOL> return sel <EOL> @ verbose <EOL> def read_inverse_operator ( fname , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> check_fname ( fname , '<STR_LIT>' , ( '<STR_LIT>' , '<STR_LIT>' ) ) <EOL> logger . info ( '<STR_LIT>' <EOL> % fname ) <EOL> f , tree , _ = fiff_open ( fname , preload = True ) <EOL> with f as fid : <EOL> invs = dir_tree_find ( tree , FIFF . FIFFB_MNE_INVERSE_SOLUTION ) <EOL> if invs is None or len ( invs ) < <NUM_LIT:1> : <EOL> raise Exception ( '<STR_LIT>' % fname ) <EOL> invs = invs [ <NUM_LIT:0> ] <EOL> parent_mri = dir_tree_find ( tree , FIFF . FIFFB_MNE_PARENT_MRI_FILE ) <EOL> if len ( parent_mri ) == <NUM_LIT:0> : <EOL> raise Exception ( '<STR_LIT>' % fname ) <EOL> parent_mri = parent_mri [ <NUM_LIT:0> ] <EOL> logger . info ( '<STR_LIT>' ) <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_INCLUDED_METHODS ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv = dict ( ) <EOL> inv [ '<STR_LIT>' ] = int ( tag . data ) <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_SOURCE_ORIENTATION ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = int ( tag . data ) <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_SOURCE_SPACE_NPOINTS ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = int ( tag . data ) <EOL> inv [ '<STR_LIT>' ] = <NUM_LIT:0> <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_COORD_FRAME ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = tag . data <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_INVERSE_SOURCE_UNIT ) <EOL> unit_dict = { FIFF . FIFF_UNIT_AM : '<STR_LIT>' , <EOL> FIFF . FIFF_UNIT_AM_M2 : '<STR_LIT>' , <EOL> FIFF . FIFF_UNIT_AM_M3 : '<STR_LIT>' } <EOL> inv [ '<STR_LIT>' ] = unit_dict . get ( int ( getattr ( tag , '<STR_LIT:data>' , - <NUM_LIT:1> ) ) , None ) <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_INVERSE_SOURCE_ORIENTATIONS ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = tag . data <EOL> logger . info ( '<STR_LIT>' ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> tag = find_tag ( fid , invs , FIFF . FIFF_MNE_INVERSE_SING ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = tag . data <EOL> inv [ '<STR_LIT>' ] = len ( inv [ '<STR_LIT>' ] ) <EOL> inv [ '<STR_LIT>' ] = False <EOL> inv [ '<STR_LIT>' ] = _read_named_matrix ( <EOL> fid , invs , FIFF . FIFF_MNE_INVERSE_LEADS , transpose = True ) <EOL> if inv [ '<STR_LIT>' ] is None : <EOL> inv [ '<STR_LIT>' ] = True <EOL> inv [ '<STR_LIT>' ] = _read_named_matrix ( <EOL> fid , invs , FIFF . FIFF_MNE_INVERSE_LEADS_WEIGHTED , <EOL> transpose = True ) <EOL> if inv [ '<STR_LIT>' ] is None : <EOL> raise ValueError ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = _read_named_matrix ( fid , invs , <EOL> FIFF . FIFF_MNE_INVERSE_FIELDS ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = Covariance ( <EOL> ** _read_cov ( fid , invs , FIFF . FIFFV_MNE_NOISE_COV , limited = True ) ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = _read_cov ( fid , invs , FIFF . FIFFV_MNE_SOURCE_COV ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = _read_cov ( fid , invs , <EOL> FIFF . FIFFV_MNE_ORIENT_PRIOR_COV ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = _read_cov ( fid , invs , <EOL> FIFF . FIFFV_MNE_DEPTH_PRIOR_COV ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = _read_cov ( fid , invs , FIFF . FIFFV_MNE_FMRI_PRIOR_COV ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT:src>' ] = _read_source_spaces_from_tree ( fid , tree , <EOL> patch_stats = False ) <EOL> for s in inv [ '<STR_LIT:src>' ] : <EOL> s [ '<STR_LIT:id>' ] = find_source_space_hemi ( s ) <EOL> tag = find_tag ( fid , parent_mri , FIFF . FIFF_COORD_TRANS ) <EOL> if tag is None : <EOL> raise Exception ( '<STR_LIT>' ) <EOL> mri_head_t = _ensure_trans ( tag . data , '<STR_LIT>' , '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = mri_head_t <EOL> inv [ '<STR_LIT:info>' ] = _read_forward_meas_info ( tree , fid ) <EOL> if inv [ '<STR_LIT>' ] not in ( FIFF . FIFFV_COORD_MRI , <EOL> FIFF . FIFFV_COORD_HEAD ) : <EOL> raise Exception ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] = <NUM_LIT:1> <EOL> inv [ '<STR_LIT>' ] = _read_proj ( fid , tree ) <EOL> inv [ '<STR_LIT>' ] = [ ] <EOL> inv [ '<STR_LIT>' ] = [ ] <EOL> inv [ '<STR_LIT>' ] = [ ] <EOL> inv [ '<STR_LIT>' ] = [ ] <EOL> nuse = <NUM_LIT:0> <EOL> for k in range ( len ( inv [ '<STR_LIT:src>' ] ) ) : <EOL> try : <EOL> inv [ '<STR_LIT:src>' ] [ k ] = transform_surface_to ( inv [ '<STR_LIT:src>' ] [ k ] , <EOL> inv [ '<STR_LIT>' ] , <EOL> mri_head_t ) <EOL> except Exception as inst : <EOL> raise Exception ( '<STR_LIT>' % inst ) <EOL> nuse += inv [ '<STR_LIT:src>' ] [ k ] [ '<STR_LIT>' ] <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> return InverseOperator ( inv ) <EOL> @ verbose <EOL> def write_inverse_operator ( fname , inv , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> check_fname ( fname , '<STR_LIT>' , ( '<STR_LIT>' , '<STR_LIT>' ) ) <EOL> logger . info ( '<STR_LIT>' % fname ) <EOL> fid = start_file ( fname ) <EOL> start_block ( fid , FIFF . FIFFB_MNE ) <EOL> write_forward_meas_info ( fid , inv [ '<STR_LIT:info>' ] ) <EOL> start_block ( fid , FIFF . FIFFB_MNE_PARENT_MRI_FILE ) <EOL> write_string ( fid , FIFF . FIFF_MNE_FILE_NAME , inv [ '<STR_LIT:info>' ] [ '<STR_LIT>' ] ) <EOL> write_coord_trans ( fid , inv [ '<STR_LIT>' ] ) <EOL> end_block ( fid , FIFF . FIFFB_MNE_PARENT_MRI_FILE ) <EOL> _write_proj ( fid , inv [ '<STR_LIT>' ] ) <EOL> if '<STR_LIT:src>' in inv : <EOL> _write_source_spaces_to_fid ( fid , inv [ '<STR_LIT:src>' ] ) <EOL> start_block ( fid , FIFF . FIFFB_MNE_INVERSE_SOLUTION ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> write_int ( fid , FIFF . FIFF_MNE_INCLUDED_METHODS , inv [ '<STR_LIT>' ] ) <EOL> write_int ( fid , FIFF . FIFF_MNE_COORD_FRAME , inv [ '<STR_LIT>' ] ) <EOL> udict = { '<STR_LIT>' : FIFF . FIFF_UNIT_AM , <EOL> '<STR_LIT>' : FIFF . FIFF_UNIT_AM_M2 , <EOL> '<STR_LIT>' : FIFF . FIFF_UNIT_AM_M3 } <EOL> if '<STR_LIT>' in inv and inv [ '<STR_LIT>' ] is not None : <EOL> write_int ( fid , FIFF . FIFF_MNE_INVERSE_SOURCE_UNIT , udict [ inv [ '<STR_LIT>' ] ] ) <EOL> write_int ( fid , FIFF . FIFF_MNE_SOURCE_ORIENTATION , inv [ '<STR_LIT>' ] ) <EOL> write_int ( fid , FIFF . FIFF_MNE_SOURCE_SPACE_NPOINTS , inv [ '<STR_LIT>' ] ) <EOL> if '<STR_LIT>' in inv : <EOL> write_int ( fid , FIFF . FIFF_NCHAN , inv [ '<STR_LIT>' ] ) <EOL> elif '<STR_LIT>' in inv [ '<STR_LIT:info>' ] : <EOL> write_int ( fid , FIFF . FIFF_NCHAN , inv [ '<STR_LIT:info>' ] [ '<STR_LIT>' ] ) <EOL> write_float_matrix ( fid , FIFF . FIFF_MNE_INVERSE_SOURCE_ORIENTATIONS , <EOL> inv [ '<STR_LIT>' ] ) <EOL> write_float ( fid , FIFF . FIFF_MNE_INVERSE_SING , inv [ '<STR_LIT>' ] ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> _write_cov ( fid , inv [ '<STR_LIT>' ] ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> _write_cov ( fid , inv [ '<STR_LIT>' ] ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> _write_cov ( fid , inv [ '<STR_LIT>' ] ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> _write_cov ( fid , inv [ '<STR_LIT>' ] ) <EOL> if inv [ '<STR_LIT>' ] is not None : <EOL> _write_cov ( fid , inv [ '<STR_LIT>' ] ) <EOL> write_named_matrix ( fid , FIFF . FIFF_MNE_INVERSE_FIELDS , inv [ '<STR_LIT>' ] ) <EOL> if inv [ '<STR_LIT>' ] : <EOL> kind = FIFF . FIFF_MNE_INVERSE_LEADS_WEIGHTED <EOL> else : <EOL> kind = FIFF . FIFF_MNE_INVERSE_LEADS <EOL> _transpose_named_matrix ( inv [ '<STR_LIT>' ] ) <EOL> write_named_matrix ( fid , kind , inv [ '<STR_LIT>' ] ) <EOL> _transpose_named_matrix ( inv [ '<STR_LIT>' ] ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> end_block ( fid , FIFF . FIFFB_MNE_INVERSE_SOLUTION ) <EOL> end_block ( fid , FIFF . FIFFB_MNE ) <EOL> end_file ( fid ) <EOL> fid . close ( ) <EOL> def combine_xyz ( vec , square = False ) : <EOL> """<STR_LIT>""" <EOL> if vec . ndim != <NUM_LIT:2> : <EOL> raise ValueError ( '<STR_LIT>' ) <EOL> if ( vec . shape [ <NUM_LIT:0> ] % <NUM_LIT:3> ) != <NUM_LIT:0> : <EOL> raise ValueError ( '<STR_LIT>' ) <EOL> n , p = vec . shape <EOL> if np . iscomplexobj ( vec ) : <EOL> vec = np . abs ( vec ) <EOL> comb = vec [ <NUM_LIT:0> : : <NUM_LIT:3> ] ** <NUM_LIT:2> <EOL> comb += vec [ <NUM_LIT:1> : : <NUM_LIT:3> ] ** <NUM_LIT:2> <EOL> comb += vec [ <NUM_LIT:2> : : <NUM_LIT:3> ] ** <NUM_LIT:2> <EOL> if not square : <EOL> comb = np . sqrt ( comb ) <EOL> return comb <EOL> def _check_ch_names ( inv , info ) : <EOL> """<STR_LIT>""" <EOL> inv_ch_names = inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] <EOL> if inv [ '<STR_LIT>' ] . ch_names != inv_ch_names : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> data_ch_names = info [ '<STR_LIT>' ] <EOL> missing_ch_names = list ( ) <EOL> for ch_name in inv_ch_names : <EOL> if ch_name not in data_ch_names : <EOL> missing_ch_names . append ( ch_name ) <EOL> n_missing = len ( missing_ch_names ) <EOL> if n_missing > <NUM_LIT:0> : <EOL> raise ValueError ( '<STR_LIT>' % n_missing + <EOL> '<STR_LIT>' % missing_ch_names ) <EOL> @ verbose <EOL> def prepare_inverse_operator ( orig , nave , lambda2 , method , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> if nave <= <NUM_LIT:0> : <EOL> raise ValueError ( '<STR_LIT>' ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv = orig . copy ( ) <EOL> scale = float ( inv [ '<STR_LIT>' ] ) / nave <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] = scale * inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> if inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] . ndim == <NUM_LIT:1> : <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] = inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] = np . eye ( len ( inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] ) ) <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] = scale * inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] = scale * inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> if inv [ '<STR_LIT>' ] : <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] = sqrt ( scale ) * inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' % ( inv [ '<STR_LIT>' ] , nave ) ) <EOL> inv [ '<STR_LIT>' ] = nave <EOL> sing = np . array ( inv [ '<STR_LIT>' ] , dtype = np . float64 ) <EOL> inv [ '<STR_LIT>' ] = sing / ( sing ** <NUM_LIT:2> + lambda2 ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> inv [ '<STR_LIT>' ] , ncomp , _ = make_projector ( inv [ '<STR_LIT>' ] , <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) <EOL> if ncomp > <NUM_LIT:0> : <EOL> logger . info ( '<STR_LIT>' <EOL> % ncomp ) <EOL> else : <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if not inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] : <EOL> inv [ '<STR_LIT>' ] = np . zeros ( ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] , <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) ) <EOL> eig = inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] <EOL> nzero = ( eig > <NUM_LIT:0> ) <EOL> inv [ '<STR_LIT>' ] [ nzero , nzero ] = <NUM_LIT:1.0> / np . sqrt ( eig [ nzero ] ) <EOL> inv [ '<STR_LIT>' ] = np . dot ( inv [ '<STR_LIT>' ] , inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> % ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] - np . sum ( nzero ) ) ) <EOL> else : <EOL> inv [ '<STR_LIT>' ] = np . diag ( <NUM_LIT:1.0> / <EOL> np . sqrt ( inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] . ravel ( ) ) ) <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> % ncomp ) <EOL> if method in [ "<STR_LIT>" , '<STR_LIT>' ] : <EOL> if method == "<STR_LIT>" : <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> noise_weight = inv [ '<STR_LIT>' ] <EOL> else : <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> noise_weight = ( inv [ '<STR_LIT>' ] * <EOL> np . sqrt ( ( <NUM_LIT:1.> + inv [ '<STR_LIT>' ] ** <NUM_LIT:2> / lambda2 ) ) ) <EOL> noise_norm = np . zeros ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) <EOL> nrm2 , = linalg . get_blas_funcs ( ( '<STR_LIT>' , ) , ( noise_norm , ) ) <EOL> if inv [ '<STR_LIT>' ] : <EOL> for k in range ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) : <EOL> one = inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] [ k , : ] * noise_weight <EOL> noise_norm [ k ] = nrm2 ( one ) <EOL> else : <EOL> for k in range ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) : <EOL> one = ( sqrt ( inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] [ k ] ) * <EOL> inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] [ k , : ] * noise_weight ) <EOL> noise_norm [ k ] = nrm2 ( one ) <EOL> if inv [ '<STR_LIT>' ] == FIFF . FIFFV_MNE_FREE_ORI : <EOL> noise_norm = combine_xyz ( noise_norm [ : , None ] ) . ravel ( ) <EOL> inv [ '<STR_LIT>' ] = <NUM_LIT:1.0> / np . abs ( noise_norm ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> else : <EOL> inv [ '<STR_LIT>' ] = [ ] <EOL> return InverseOperator ( inv ) <EOL> @ verbose <EOL> def _assemble_kernel ( inv , label , method , pick_ori , verbose = None ) : <EOL> eigen_leads = inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> source_cov = inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] [ : , None ] <EOL> if method != "<STR_LIT>" : <EOL> noise_norm = inv [ '<STR_LIT>' ] [ : , None ] <EOL> src = inv [ '<STR_LIT:src>' ] <EOL> vertno = _get_vertno ( src ) <EOL> if label is not None : <EOL> vertno , src_sel = label_src_vertno_sel ( label , inv [ '<STR_LIT:src>' ] ) <EOL> if method != "<STR_LIT>" : <EOL> noise_norm = noise_norm [ src_sel ] <EOL> if inv [ '<STR_LIT>' ] == FIFF . FIFFV_MNE_FREE_ORI : <EOL> src_sel = <NUM_LIT:3> * src_sel <EOL> src_sel = np . c_ [ src_sel , src_sel + <NUM_LIT:1> , src_sel + <NUM_LIT:2> ] <EOL> src_sel = src_sel . ravel ( ) <EOL> eigen_leads = eigen_leads [ src_sel ] <EOL> source_cov = source_cov [ src_sel ] <EOL> if pick_ori == "<STR_LIT>" : <EOL> if not inv [ '<STR_LIT>' ] == FIFF . FIFFV_MNE_FREE_ORI : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> is_loose = <NUM_LIT:0> < inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] [ <NUM_LIT:0> ] < <NUM_LIT:1> <EOL> if not is_loose : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> eigen_leads = eigen_leads [ <NUM_LIT:2> : : <NUM_LIT:3> ] <EOL> source_cov = source_cov [ <NUM_LIT:2> : : <NUM_LIT:3> ] <EOL> trans = inv [ '<STR_LIT>' ] [ : , None ] * reduce ( np . dot , <EOL> [ inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] , <EOL> inv [ '<STR_LIT>' ] , <EOL> inv [ '<STR_LIT>' ] ] ) <EOL> if inv [ '<STR_LIT>' ] : <EOL> logger . info ( '<STR_LIT>' ) <EOL> K = np . dot ( eigen_leads , trans ) <EOL> else : <EOL> logger . info ( '<STR_LIT>' ) <EOL> K = np . sqrt ( source_cov ) * np . dot ( eigen_leads , trans ) <EOL> if method == "<STR_LIT>" : <EOL> noise_norm = None <EOL> return K , noise_norm , vertno <EOL> def _check_method ( method ) : <EOL> if method not in [ "<STR_LIT>" , "<STR_LIT>" , "<STR_LIT>" ] : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> return method <EOL> def _check_ori ( pick_ori ) : <EOL> if pick_ori is not None and pick_ori != '<STR_LIT>' : <EOL> raise RuntimeError ( '<STR_LIT>' <EOL> % pick_ori ) <EOL> return pick_ori <EOL> def _check_reference ( inst ) : <EOL> """<STR_LIT>""" <EOL> if _needs_eeg_average_ref_proj ( inst . info ) : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if inst . info [ '<STR_LIT>' ] : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> def _subject_from_inverse ( inverse_operator ) : <EOL> """<STR_LIT>""" <EOL> return inverse_operator [ '<STR_LIT:src>' ] [ <NUM_LIT:0> ] . get ( '<STR_LIT>' , None ) <EOL> @ verbose <EOL> def apply_inverse ( evoked , inverse_operator , lambda2 = <NUM_LIT:1.> / <NUM_LIT> , <EOL> method = "<STR_LIT>" , pick_ori = None , <EOL> prepared = False , label = None , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> _check_reference ( evoked ) <EOL> method = _check_method ( method ) <EOL> pick_ori = _check_ori ( pick_ori ) <EOL> nave = evoked . nave <EOL> _check_ch_names ( inverse_operator , evoked . info ) <EOL> if not prepared : <EOL> inv = prepare_inverse_operator ( inverse_operator , nave , lambda2 , method ) <EOL> else : <EOL> inv = inverse_operator <EOL> sel = _pick_channels_inverse_operator ( evoked . ch_names , inv ) <EOL> logger . info ( '<STR_LIT>' % len ( sel ) ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> K , noise_norm , vertno = _assemble_kernel ( inv , label , method , pick_ori ) <EOL> sol = np . dot ( K , evoked . data [ sel ] ) <EOL> is_free_ori = ( inverse_operator [ '<STR_LIT>' ] == <EOL> FIFF . FIFFV_MNE_FREE_ORI and pick_ori is None ) <EOL> if is_free_ori : <EOL> logger . info ( '<STR_LIT>' ) <EOL> sol = combine_xyz ( sol ) <EOL> if noise_norm is not None : <EOL> logger . info ( '<STR_LIT>' ) <EOL> sol *= noise_norm <EOL> tstep = <NUM_LIT:1.0> / evoked . info [ '<STR_LIT>' ] <EOL> tmin = float ( evoked . times [ <NUM_LIT:0> ] ) <EOL> subject = _subject_from_inverse ( inverse_operator ) <EOL> stc = _make_stc ( sol , vertices = vertno , tmin = tmin , tstep = tstep , <EOL> subject = subject ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> return stc <EOL> @ verbose <EOL> def apply_inverse_raw ( raw , inverse_operator , lambda2 , method = "<STR_LIT>" , <EOL> label = None , start = None , stop = None , nave = <NUM_LIT:1> , <EOL> time_func = None , pick_ori = None , buffer_size = None , <EOL> prepared = False , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> _check_reference ( raw ) <EOL> method = _check_method ( method ) <EOL> pick_ori = _check_ori ( pick_ori ) <EOL> _check_ch_names ( inverse_operator , raw . info ) <EOL> if not prepared : <EOL> inv = prepare_inverse_operator ( inverse_operator , nave , lambda2 , method ) <EOL> else : <EOL> inv = inverse_operator <EOL> sel = _pick_channels_inverse_operator ( raw . ch_names , inv ) <EOL> logger . info ( '<STR_LIT>' % len ( sel ) ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> data , times = raw [ sel , start : stop ] <EOL> if time_func is not None : <EOL> data = time_func ( data ) <EOL> K , noise_norm , vertno = _assemble_kernel ( inv , label , method , pick_ori ) <EOL> is_free_ori = ( inverse_operator [ '<STR_LIT>' ] == <EOL> FIFF . FIFFV_MNE_FREE_ORI and pick_ori is None ) <EOL> if buffer_size is not None and is_free_ori : <EOL> n_seg = int ( np . ceil ( data . shape [ <NUM_LIT:1> ] / float ( buffer_size ) ) ) <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' % ( n_seg ) ) <EOL> n_times = data . shape [ <NUM_LIT:1> ] <EOL> sol = np . empty ( ( K . shape [ <NUM_LIT:0> ] // <NUM_LIT:3> , n_times ) , <EOL> dtype = ( K [ <NUM_LIT:0> , <NUM_LIT:0> ] * data [ <NUM_LIT:0> , <NUM_LIT:0> ] ) . dtype ) <EOL> for pos in range ( <NUM_LIT:0> , n_times , buffer_size ) : <EOL> sol [ : , pos : pos + buffer_size ] = combine_xyz ( np . dot ( K , data [ : , pos : pos + buffer_size ] ) ) <EOL> logger . info ( '<STR_LIT>' <EOL> % ( pos / buffer_size + <NUM_LIT:1> , n_seg ) ) <EOL> else : <EOL> sol = np . dot ( K , data ) <EOL> if is_free_ori : <EOL> logger . info ( '<STR_LIT>' ) <EOL> sol = combine_xyz ( sol ) <EOL> if noise_norm is not None : <EOL> sol *= noise_norm <EOL> tmin = float ( times [ <NUM_LIT:0> ] ) <EOL> tstep = <NUM_LIT:1.0> / raw . info [ '<STR_LIT>' ] <EOL> subject = _subject_from_inverse ( inverse_operator ) <EOL> stc = _make_stc ( sol , vertices = vertno , tmin = tmin , tstep = tstep , <EOL> subject = subject ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> return stc <EOL> def _apply_inverse_epochs_gen ( epochs , inverse_operator , lambda2 , method = '<STR_LIT>' , <EOL> label = None , nave = <NUM_LIT:1> , pick_ori = None , <EOL> prepared = False , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> method = _check_method ( method ) <EOL> pick_ori = _check_ori ( pick_ori ) <EOL> _check_ch_names ( inverse_operator , epochs . info ) <EOL> if not prepared : <EOL> inv = prepare_inverse_operator ( inverse_operator , nave , lambda2 , method ) <EOL> else : <EOL> inv = inverse_operator <EOL> sel = _pick_channels_inverse_operator ( epochs . ch_names , inv ) <EOL> logger . info ( '<STR_LIT>' % len ( sel ) ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> K , noise_norm , vertno = _assemble_kernel ( inv , label , method , pick_ori ) <EOL> tstep = <NUM_LIT:1.0> / epochs . info [ '<STR_LIT>' ] <EOL> tmin = epochs . times [ <NUM_LIT:0> ] <EOL> is_free_ori = ( inverse_operator [ '<STR_LIT>' ] == <EOL> FIFF . FIFFV_MNE_FREE_ORI and pick_ori is None ) <EOL> if not is_free_ori and noise_norm is not None : <EOL> K *= noise_norm <EOL> subject = _subject_from_inverse ( inverse_operator ) <EOL> for k , e in enumerate ( epochs ) : <EOL> logger . info ( '<STR_LIT>' % ( k + <NUM_LIT:1> ) ) <EOL> if is_free_ori : <EOL> sol = np . dot ( K , e [ sel ] ) <EOL> if is_free_ori : <EOL> logger . info ( '<STR_LIT>' ) <EOL> sol = combine_xyz ( sol ) <EOL> if noise_norm is not None : <EOL> sol *= noise_norm <EOL> else : <EOL> if len ( sel ) < K . shape [ <NUM_LIT:0> ] : <EOL> sol = ( K , e [ sel ] ) <EOL> else : <EOL> sol = np . dot ( K , e [ sel ] ) <EOL> stc = _make_stc ( sol , vertices = vertno , tmin = tmin , tstep = tstep , <EOL> subject = subject ) <EOL> yield stc <EOL> logger . info ( '<STR_LIT>' ) <EOL> @ verbose <EOL> def apply_inverse_epochs ( epochs , inverse_operator , lambda2 , method = "<STR_LIT>" , <EOL> label = None , nave = <NUM_LIT:1> , pick_ori = None , <EOL> return_generator = False , <EOL> prepared = False , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> _check_reference ( epochs ) <EOL> stcs = _apply_inverse_epochs_gen ( epochs , inverse_operator , lambda2 , <EOL> method = method , label = label , nave = nave , <EOL> pick_ori = pick_ori , verbose = verbose , <EOL> prepared = prepared ) <EOL> if not return_generator : <EOL> stcs = [ stc for stc in stcs ] <EOL> return stcs <EOL> '''<STR_LIT>''' <EOL> @ verbose <EOL> def _prepare_forward ( forward , info , noise_cov , pca = False , rank = None , <EOL> verbose = None ) : <EOL> """<STR_LIT>""" <EOL> fwd_sol_ch_names = forward [ '<STR_LIT>' ] [ '<STR_LIT>' ] <EOL> ch_names = [ c [ '<STR_LIT>' ] for c in info [ '<STR_LIT>' ] <EOL> if ( ( c [ '<STR_LIT>' ] not in info [ '<STR_LIT>' ] and <EOL> c [ '<STR_LIT>' ] not in noise_cov [ '<STR_LIT>' ] ) and <EOL> ( c [ '<STR_LIT>' ] in fwd_sol_ch_names and <EOL> c [ '<STR_LIT>' ] in noise_cov . ch_names ) ) ] <EOL> if not len ( info [ '<STR_LIT>' ] ) == len ( noise_cov [ '<STR_LIT>' ] ) or not all ( b in noise_cov [ '<STR_LIT>' ] for b in info [ '<STR_LIT>' ] ) : <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> n_chan = len ( ch_names ) <EOL> logger . info ( "<STR_LIT>" % n_chan ) <EOL> noise_cov = prepare_noise_cov ( noise_cov , info , ch_names , rank ) <EOL> eig = noise_cov [ '<STR_LIT>' ] <EOL> nzero = ( eig > <NUM_LIT:0> ) <EOL> n_nzero = sum ( nzero ) <EOL> if pca : <EOL> whitener = noise_cov [ '<STR_LIT>' ] [ nzero ] / np . sqrt ( eig [ nzero ] ) [ : , None ] <EOL> logger . info ( '<STR_LIT>' % n_nzero ) <EOL> else : <EOL> whitener = np . zeros ( ( n_chan , n_chan ) , dtype = np . float ) <EOL> whitener [ nzero , nzero ] = <NUM_LIT:1.0> / np . sqrt ( eig [ nzero ] ) <EOL> whitener = np . dot ( whitener , noise_cov [ '<STR_LIT>' ] ) <EOL> gain = forward [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] <EOL> fwd_idx = [ fwd_sol_ch_names . index ( name ) for name in ch_names ] <EOL> gain = gain [ fwd_idx ] <EOL> info_idx = [ info [ '<STR_LIT>' ] . index ( name ) for name in ch_names ] <EOL> fwd_info = pick_info ( info , info_idx ) <EOL> logger . info ( '<STR_LIT>' % n_nzero ) <EOL> return fwd_info , gain , noise_cov , whitener , n_nzero <EOL> @ verbose <EOL> def make_inverse_operator ( info , forward , noise_cov , loose = <NUM_LIT> , depth = <NUM_LIT> , <EOL> fixed = False , limit_depth_chs = True , rank = None , <EOL> verbose = None ) : <EOL> """<STR_LIT>""" <EOL> is_fixed_ori = is_fixed_orient ( forward ) <EOL> if fixed and loose is not None : <EOL> warn ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> loose = None <EOL> if is_fixed_ori and not fixed : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if fixed : <EOL> if depth is not None : <EOL> if is_fixed_ori or not forward [ '<STR_LIT>' ] : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> elif forward [ '<STR_LIT>' ] is False : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if depth is not None : <EOL> if not ( <NUM_LIT:0> < depth <= <NUM_LIT:1> ) : <EOL> raise ValueError ( '<STR_LIT>' ) <EOL> if is_fixed_ori : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if not forward [ '<STR_LIT>' ] : <EOL> forward = convert_forward_solution ( forward , surf_ori = True ) <EOL> assert forward [ '<STR_LIT>' ] <EOL> if loose is not None : <EOL> if not ( <NUM_LIT:0> <= loose <= <NUM_LIT:1> ) : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if loose < <NUM_LIT:1> and not forward [ '<STR_LIT>' ] : <EOL> raise ValueError ( '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> gain_info , gain , noise_cov , whitener , n_nzero = _prepare_forward ( forward , info , noise_cov , rank = rank ) <EOL> forward [ '<STR_LIT:info>' ] . _check_consistency ( ) <EOL> if depth is not None : <EOL> patch_areas = forward . get ( '<STR_LIT>' , None ) <EOL> depth_prior = compute_depth_prior ( gain , gain_info , is_fixed_ori , <EOL> exp = depth , patch_areas = patch_areas , <EOL> limit_depth_chs = limit_depth_chs ) <EOL> else : <EOL> depth_prior = np . ones ( gain . shape [ <NUM_LIT:1> ] , dtype = gain . dtype ) <EOL> if fixed : <EOL> if depth is not None : <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> if not is_fixed_ori : <EOL> depth_prior = depth_prior [ <NUM_LIT:2> : : <NUM_LIT:3> ] <EOL> forward = convert_forward_solution ( <EOL> forward , surf_ori = forward [ '<STR_LIT>' ] , force_fixed = True ) <EOL> is_fixed_ori = is_fixed_orient ( forward ) <EOL> gain_info , gain , noise_cov , whitener , n_nzero = _prepare_forward ( forward , info , noise_cov , verbose = False ) <EOL> logger . info ( "<STR_LIT>" <EOL> % len ( gain_info [ '<STR_LIT>' ] ) ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> source_cov = depth_prior . copy ( ) <EOL> depth_prior = dict ( data = depth_prior , kind = FIFF . FIFFV_MNE_DEPTH_PRIOR_COV , <EOL> bads = [ ] , diag = True , names = [ ] , eig = None , <EOL> eigvec = None , dim = depth_prior . size , nfree = <NUM_LIT:1> , <EOL> projs = [ ] ) <EOL> if not is_fixed_ori : <EOL> orient_prior = compute_orient_prior ( forward , loose = loose ) <EOL> source_cov *= orient_prior <EOL> orient_prior = dict ( data = orient_prior , <EOL> kind = FIFF . FIFFV_MNE_ORIENT_PRIOR_COV , <EOL> bads = [ ] , diag = True , names = [ ] , eig = None , <EOL> eigvec = None , dim = orient_prior . size , nfree = <NUM_LIT:1> , <EOL> projs = [ ] ) <EOL> else : <EOL> orient_prior = None <EOL> logger . info ( '<STR_LIT>' ) <EOL> gain = np . dot ( whitener , gain ) <EOL> logger . info ( '<STR_LIT>' ) <EOL> source_std = np . sqrt ( source_cov ) <EOL> gain *= source_std [ None , : ] <EOL> trace_GRGT = linalg . norm ( gain , ord = '<STR_LIT>' ) ** <NUM_LIT:2> <EOL> scaling_source_cov = n_nzero / trace_GRGT <EOL> source_cov *= scaling_source_cov <EOL> gain *= sqrt ( scaling_source_cov ) <EOL> source_cov = dict ( data = source_cov , dim = source_cov . size , <EOL> kind = FIFF . FIFFV_MNE_SOURCE_COV , diag = True , <EOL> names = [ ] , projs = [ ] , eig = None , eigvec = None , <EOL> nfree = <NUM_LIT:1> , bads = [ ] ) <EOL> logger . info ( '<STR_LIT>' <EOL> '<STR_LIT>' ) <EOL> eigen_fields , sing , eigen_leads = linalg . svd ( gain , full_matrices = False ) <EOL> logger . info ( '<STR_LIT>' % np . max ( sing ) ) <EOL> logger . info ( '<STR_LIT>' % trace_GRGT ) <EOL> eigen_fields = dict ( data = eigen_fields . T , col_names = gain_info [ '<STR_LIT>' ] , <EOL> row_names = [ ] , nrow = eigen_fields . shape [ <NUM_LIT:1> ] , <EOL> ncol = eigen_fields . shape [ <NUM_LIT:0> ] ) <EOL> eigen_leads = dict ( data = eigen_leads . T , nrow = eigen_leads . shape [ <NUM_LIT:1> ] , <EOL> ncol = eigen_leads . shape [ <NUM_LIT:0> ] , row_names = [ ] , <EOL> col_names = [ ] ) <EOL> nave = <NUM_LIT:1.0> <EOL> has_meg = False <EOL> has_eeg = False <EOL> ch_idx = [ k for k , c in enumerate ( info [ '<STR_LIT>' ] ) <EOL> if c [ '<STR_LIT>' ] in gain_info [ '<STR_LIT>' ] ] <EOL> for idx in ch_idx : <EOL> ch_type = channel_type ( info , idx ) <EOL> if ch_type == '<STR_LIT>' : <EOL> has_eeg = True <EOL> if ( ch_type == '<STR_LIT>' ) or ( ch_type == '<STR_LIT>' ) : <EOL> has_meg = True <EOL> if has_eeg and has_meg : <EOL> methods = FIFF . FIFFV_MNE_MEG_EEG <EOL> elif has_meg : <EOL> methods = FIFF . FIFFV_MNE_MEG <EOL> else : <EOL> methods = FIFF . FIFFV_MNE_EEG <EOL> if depth is None : <EOL> depth_prior = None <EOL> inv_op = dict ( eigen_fields = eigen_fields , eigen_leads = eigen_leads , <EOL> sing = sing , nave = nave , depth_prior = depth_prior , <EOL> source_cov = source_cov , noise_cov = noise_cov , <EOL> orient_prior = orient_prior , projs = deepcopy ( info [ '<STR_LIT>' ] ) , <EOL> eigen_leads_weighted = False , source_ori = forward [ '<STR_LIT>' ] , <EOL> mri_head_t = deepcopy ( forward [ '<STR_LIT>' ] ) , <EOL> methods = methods , nsource = forward [ '<STR_LIT>' ] , <EOL> coord_frame = forward [ '<STR_LIT>' ] , <EOL> source_nn = forward [ '<STR_LIT>' ] . copy ( ) , <EOL> src = deepcopy ( forward [ '<STR_LIT:src>' ] ) , fmri_prior = None ) <EOL> inv_info = deepcopy ( forward [ '<STR_LIT:info>' ] ) <EOL> inv_info [ '<STR_LIT>' ] = [ bad for bad in info [ '<STR_LIT>' ] <EOL> if bad in inv_info [ '<STR_LIT>' ] ] <EOL> inv_info . _check_consistency ( ) <EOL> inv_op [ '<STR_LIT>' ] = '<STR_LIT>' <EOL> inv_op [ '<STR_LIT:info>' ] = inv_info <EOL> return InverseOperator ( inv_op ) <EOL> def compute_rank_inverse ( inv ) : <EOL> """<STR_LIT>""" <EOL> eig = inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] <EOL> if not inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] : <EOL> rank = np . sum ( eig > <NUM_LIT:0> ) <EOL> else : <EOL> ncomp = make_projector ( inv [ '<STR_LIT>' ] , inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] ) [ <NUM_LIT:1> ] <EOL> rank = inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] - ncomp <EOL> return rank <EOL> @ verbose <EOL> def estimate_snr ( evoked , inv , verbose = None ) : <EOL> """<STR_LIT>""" <EOL> from scipy . stats import chi2 <EOL> _check_reference ( evoked ) <EOL> _check_ch_names ( inv , evoked . info ) <EOL> inv = prepare_inverse_operator ( inv , evoked . nave , <NUM_LIT:1.> / <NUM_LIT> , '<STR_LIT>' ) <EOL> sel = _pick_channels_inverse_operator ( evoked . ch_names , inv ) <EOL> logger . info ( '<STR_LIT>' % len ( sel ) ) <EOL> data_white = np . dot ( inv [ '<STR_LIT>' ] , np . dot ( inv [ '<STR_LIT>' ] , evoked . data [ sel ] ) ) <EOL> data_white_ef = np . dot ( inv [ '<STR_LIT>' ] [ '<STR_LIT:data>' ] , data_white ) <EOL> n_ch , n_times = data_white . shape <EOL> n_zero = ( inv [ '<STR_LIT>' ] [ '<STR_LIT>' ] <= <NUM_LIT:0> ) . sum ( ) <EOL> logger . info ( '<STR_LIT>' <EOL> % ( n_ch , n_zero , n_ch - n_zero ) ) <EOL> signal = np . sum ( data_white ** <NUM_LIT:2> , axis = <NUM_LIT:0> ) <EOL> noise = n_ch - n_zero <EOL> snr = signal / noise <EOL> lambda2_est = np . empty ( n_times ) <EOL> lambda2_est . fill ( <NUM_LIT> ) <EOL> remaining = np . ones ( n_times , bool ) <EOL> bad = ( snr <= <NUM_LIT:1> ) <EOL> lambda2_est [ bad ] = <NUM_LIT> <EOL> remaining [ bad ] = False <EOL> lambda_mult = <NUM_LIT> <EOL> sing2 = ( inv [ '<STR_LIT>' ] * inv [ '<STR_LIT>' ] ) [ : , np . newaxis ] <EOL> val = chi2 . isf ( <NUM_LIT> , n_ch - <NUM_LIT:1> ) <EOL> for n_iter in range ( <NUM_LIT:1000> ) : <EOL> f = sing2 / ( sing2 + lambda2_est [ np . newaxis ] [ : , remaining ] ) <EOL> f [ inv [ '<STR_LIT>' ] == <NUM_LIT:0> ] = <NUM_LIT:0> <EOL> ew = data_white_ef [ : , remaining ] * ( <NUM_LIT:1.0> - f ) <EOL> err = np . sum ( ew * ew , axis = <NUM_LIT:0> ) <EOL> remaining [ np . where ( remaining ) [ <NUM_LIT:0> ] [ err < val ] ] = False <EOL> if not remaining . any ( ) : <EOL> break <EOL> lambda2_est [ remaining ] *= lambda_mult <EOL> else : <EOL> warn ( '<STR_LIT>' ) <EOL> snr_est = <NUM_LIT:1.0> / np . sqrt ( lambda2_est ) <EOL> snr = np . sqrt ( snr ) <EOL> return snr , snr_est </s>
