This is an automated email from the git hooks/post-receive script. yoh pushed a commit to tag 0.4 in repository python-mne.
commit fc767b8957a6cfc27180a4006d0b561276777169 Author: Alexandre Gramfort <[email protected]> Date: Thu May 24 19:10:22 2012 +0200 ENH : expose qrs_detector params in mne_compute_ecg_proj.py (add --tstart option) --- bin/mne_compute_proj_ecg.py | 14 +++++++++++++- bin/mne_compute_proj_eog.py | 14 +++++++++++++- mne/artifacts/ecg.py | 36 +++++++++++++++++++++++++++------- mne/artifacts/eog.py | 10 +++++++--- mne/preprocessing/ssp.py | 47 ++++++++++++++++++++++++++++++++++++++------- 5 files changed, 102 insertions(+), 19 deletions(-) diff --git a/bin/mne_compute_proj_ecg.py b/bin/mne_compute_proj_ecg.py index acdc104..b05cc0d 100755 --- a/bin/mne_compute_proj_ecg.py +++ b/bin/mne_compute_proj_ecg.py @@ -42,6 +42,12 @@ if __name__ == '__main__': parser.add_option("--h-freq", dest="h_freq", type="float", help="Filter high cut-off frequency in Hz", default=100) + parser.add_option("--ecg-l-freq", dest="ecg_l_freq", type="float", + help="Filter low cut-off frequency in Hz used for ECG event detection", + default=5) + parser.add_option("--ecg-h-freq", dest="ecg_h_freq", type="float", + help="Filter high cut-off frequency in Hz used for ECG event detection", + default=35) parser.add_option("-p", "--preload", dest="preload", help="Temporary file used during computation (to save memory)", default=True) @@ -85,6 +91,8 @@ if __name__ == '__main__': help="ID to use for events", default=999) parser.add_option("--event-raw", dest="raw_event_fname", help="raw file to use for event detection", default=None) + parser.add_option("--tstart", dest="tstart", type="float", + help="Start artifact detection after tstart seconds", default=0.) options, args = parser.parse_args() @@ -101,6 +109,8 @@ if __name__ == '__main__': n_eeg = options.n_eeg l_freq = options.l_freq h_freq = options.h_freq + ecg_l_freq = options.ecg_l_freq + ecg_h_freq = options.ecg_h_freq average = options.average preload = options.preload filter_length = options.filter_length @@ -116,6 +126,7 @@ if __name__ == '__main__': event_id = options.event_id proj_fname = options.proj raw_event_fname = options.raw_event_fname + tstart = options.tstart if bad_fname is not None: bads = [w.rstrip().split()[0] for w in open(bad_fname).readlines()] @@ -146,7 +157,8 @@ if __name__ == '__main__': tmin, tmax, n_grad, n_mag, n_eeg, l_freq, h_freq, average, filter_length, n_jobs, ch_name, reject, - bads, avg_ref, no_proj, event_id) + bads, avg_ref, no_proj, event_id, + ecg_l_freq, ecg_h_freq, tstart) raw.close() diff --git a/bin/mne_compute_proj_eog.py b/bin/mne_compute_proj_eog.py index 08c4968..770348d 100755 --- a/bin/mne_compute_proj_eog.py +++ b/bin/mne_compute_proj_eog.py @@ -48,6 +48,12 @@ if __name__ == '__main__': parser.add_option("--h-freq", dest="h_freq", type="float", help="Filter high cut-off frequency in Hz", default=35) + parser.add_option("--eog-l-freq", dest="eog_l_freq", type="float", + help="Filter low cut-off frequency in Hz used for EOG event detection", + default=1) + parser.add_option("--eog-h-freq", dest="eog_h_freq", type="float", + help="Filter high cut-off frequency in Hz used for EOG event detection", + default=10) parser.add_option("-p", "--preload", dest="preload", help="Temporary file used during computation (to save memory)", default=True) @@ -88,6 +94,8 @@ if __name__ == '__main__': help="ID to use for events", default=998) parser.add_option("--event-raw", dest="raw_event_fname", help="raw file to use for event detection", default=None) + parser.add_option("--tstart", dest="tstart", type="float", + help="Start artifact detection after tstart seconds", default=0.) options, args = parser.parse_args() @@ -104,6 +112,8 @@ if __name__ == '__main__': n_eeg = options.n_eeg l_freq = options.l_freq h_freq = options.h_freq + eog_l_freq = options.eog_l_freq + eog_h_freq = options.eog_h_freq average = options.average preload = options.preload filter_length = options.filter_length @@ -118,6 +128,7 @@ if __name__ == '__main__': event_id = options.event_id proj_fname = options.proj raw_event_fname = options.raw_event_fname + tstart = options.tstart if bad_fname is not None: bads = [w.rstrip().split()[0] for w in open(bad_fname).readlines()] @@ -148,7 +159,8 @@ if __name__ == '__main__': tmin, tmax, n_grad, n_mag, n_eeg, l_freq, h_freq, average, filter_length, n_jobs, reject, bads, - avg_ref, no_proj, event_id) + avg_ref, no_proj, event_id, + eog_l_freq, eog_h_freq, tstart) raw.close() diff --git a/mne/artifacts/ecg.py b/mne/artifacts/ecg.py index 18b7c43..456097a 100644 --- a/mne/artifacts/ecg.py +++ b/mne/artifacts/ecg.py @@ -5,7 +5,7 @@ from ..filter import band_pass_filter def qrs_detector(sfreq, ecg, thresh_value=0.6, levels=2.5, n_thresh=3, - low_pass=5, high_pass=35): + l_freq=5, h_freq=35, tstart=0): """Detect QRS component in ECG channels. QRS is the main wave on the heart beat. @@ -22,10 +22,12 @@ def qrs_detector(sfreq, ecg, thresh_value=0.6, levels=2.5, n_thresh=3, number of std from mean to include for detection n_thresh: int max number of crossings - low_pass: float + l_freq: float Low pass frequency - high_pass: float + h_freq: float High pass frequency + tstart: float + Start detection after tstart seconds. Returns ------- @@ -34,12 +36,16 @@ def qrs_detector(sfreq, ecg, thresh_value=0.6, levels=2.5, n_thresh=3, """ win_size = round((60.0 * sfreq) / 120.0) - filtecg = band_pass_filter(ecg, sfreq, low_pass, high_pass) - n_points = len(filtecg) + filtecg = band_pass_filter(ecg, sfreq, l_freq, h_freq) absecg = np.abs(filtecg) init = int(sfreq) + n_samples_start = int(init * tstart) + absecg = absecg[n_samples_start:] + + n_points = len(absecg) + maxpt = np.empty(3) maxpt[0] = np.max(absecg[:init]) maxpt[1] = np.max(absecg[init:init * 2]) @@ -73,10 +79,13 @@ def qrs_detector(sfreq, ecg, thresh_value=0.6, levels=2.5, n_thresh=3, a = np.array(numcross)[b] clean_events = time[b[a < n_thresh]] + clean_events += n_samples_start + return clean_events -def find_ecg_events(raw, event_id=999, ch_name=None): +def find_ecg_events(raw, event_id=999, ch_name=None, tstart=0.0, + l_freq=5, h_freq=35, qrs_threshold=0.6): """Find ECG peaks Parameters @@ -89,6 +98,15 @@ def find_ecg_events(raw, event_id=999, ch_name=None): The name of the channel to use for ECG peak detection. The argument is mandatory if the dataset contains no ECG channels. + tstart: float + Start detection after tstart seconds. Useful when beginning + of run is noisy. + l_freq: float + Low pass frequency + h_freq: float + High pass frequency + qrs_threshold: float + Between 0 and 1. qrs detection threshold Returns ------- @@ -122,7 +140,11 @@ def find_ecg_events(raw, event_id=999, ch_name=None): ecg, times = raw[ch_ECG, :] # detecting QRS and generating event file - ecg_events = qrs_detector(info['sfreq'], ecg.ravel()) + ecg_events = qrs_detector(info['sfreq'], ecg.ravel(), tstart=tstart, + thresh_value=qrs_threshold, l_freq=l_freq, + h_freq=h_freq) + import ipdb; ipdb.set_trace() + n_events = len(ecg_events) average_pulse = n_events * 60.0 / (times[-1] - times[0]) print ("Number of ECG events detected : %d (average pulse %d / min.)" diff --git a/mne/artifacts/eog.py b/mne/artifacts/eog.py index 15ff8da..c9d7431 100644 --- a/mne/artifacts/eog.py +++ b/mne/artifacts/eog.py @@ -5,7 +5,7 @@ from .. import fiff from ..filter import band_pass_filter -def find_eog_events(raw, event_id=998): +def find_eog_events(raw, event_id=998, l_freq=1, h_freq=10): """Locate EOG artifacts Parameters @@ -14,6 +14,10 @@ def find_eog_events(raw, event_id=998): The raw data event_id : int The index to assign to found events + low_pass: float + Low pass frequency + high_pass: float + High pass frequency Returns ------- @@ -49,8 +53,8 @@ def find_eog_events(raw, event_id=998): indexmax = np.argmax(temp) - # easy to detect peaks with this filtering. - filteog = band_pass_filter(eog[indexmax], sampRate, 1, 10) + # easier to detect peaks with filtering. + filteog = band_pass_filter(eog[indexmax], sampRate, l_freq, h_freq) # detecting eog blinks and generating event file diff --git a/mne/preprocessing/ssp.py b/mne/preprocessing/ssp.py index 711df3f..60713af 100644 --- a/mne/preprocessing/ssp.py +++ b/mne/preprocessing/ssp.py @@ -14,7 +14,8 @@ from ..artifacts import find_ecg_events, find_eog_events def _compute_exg_proj(mode, raw, raw_event, tmin, tmax, n_grad, n_mag, n_eeg, l_freq, h_freq, average, filter_length, n_jobs, ch_name, - reject, bads, avg_ref, no_proj, event_id): + reject, bads, avg_ref, no_proj, event_id, + exg_l_freq, exg_h_freq, tstart): """Compute SSP/PCA projections for ECG or EOG artifacts Note: raw has to be constructed with preload=True (or string) @@ -79,6 +80,15 @@ def _compute_exg_proj(mode, raw, raw_event, tmin, tmax, event_id: int ID to use for events + exg_l_freq: float + Low pass frequency applied for filtering EXG channel + + exg_h_freq: float + High pass frequency applied for filtering EXG channel + + tstart: float + Start artifact detection after tstart seconds. + Returns ------- proj : list @@ -108,10 +118,12 @@ def _compute_exg_proj(mode, raw, raw_event, tmin, tmax, if mode == 'ECG': print 'Running ECG SSP computation' events, _, _ = find_ecg_events(raw_event, ch_name=ch_name, - event_id=event_id) + event_id=event_id, l_freq=exg_l_freq, + h_freq=exg_h_freq, tstart=tstart) elif mode == 'EOG': print 'Running EOG SSP computation' - events = find_eog_events(raw_event, event_id=event_id) + events = find_eog_events(raw_event, event_id=event_id, + l_freq=exg_l_freq, h_freq=exg_h_freq) else: ValueError("mode must be 'ECG' or 'EOG'") @@ -162,7 +174,8 @@ def compute_proj_ecg(raw, raw_event=None, tmin=-0.2, tmax=0.4, average=False, filter_length=2048, n_jobs=1, ch_name=None, reject=dict(grad=2000e-13, mag=3000e-15, eeg=50e-6, eog=250e-6), bads=[], avg_ref=False, no_proj=False, - event_id=999): + event_id=999, ecg_l_freq=5, ecg_h_freq=35, + tstart=0.): """Compute SSP/PCA projections for ECG artifacts Note: raw has to be constructed with preload=True (or string) @@ -224,6 +237,15 @@ def compute_proj_ecg(raw, raw_event=None, tmin=-0.2, tmax=0.4, event_id: int ID to use for events + ecg_l_freq: float + Low pass frequency applied for filtering ECG channel + + ecg_h_freq: float + High pass frequency applied for filtering ECG channel + + tstart: float + Start artifact detection after tstart seconds. + Returns ------- proj : list @@ -236,7 +258,8 @@ def compute_proj_ecg(raw, raw_event=None, tmin=-0.2, tmax=0.4, projs, ecg_events = _compute_exg_proj('ECG', raw, raw_event, tmin, tmax, n_grad, n_mag, n_eeg, l_freq, h_freq, average, filter_length, n_jobs, ch_name, - reject, bads, avg_ref, no_proj, event_id) + reject, bads, avg_ref, no_proj, event_id, + ecg_l_freq, ecg_h_freq, tstart) return projs, ecg_events @@ -246,7 +269,7 @@ def compute_proj_eog(raw, raw_event=None, tmin=-0.2, tmax=0.2, average=False, filter_length=2048, n_jobs=1, reject=dict(grad=2000e-13, mag=3000e-15, eeg=500e-6, eog=np.inf), bads=[], avg_ref=False, no_proj=False, - event_id=998): + event_id=998, eog_l_freq=1, eog_h_freq=10, tstart=0.): """Compute SSP/PCA projections for EOG artifacts Note: raw has to be constructed with preload=True (or string) @@ -308,6 +331,15 @@ def compute_proj_eog(raw, raw_event=None, tmin=-0.2, tmax=0.2, event_id: int ID to use for events + eog_l_freq: float + Low pass frequency applied for filtering E0G channel + + eog_h_freq: float + High pass frequency applied for filtering E0G channel + + tstart: float + Start artifact detection after tstart seconds. + Returns ------- proj : list @@ -320,6 +352,7 @@ def compute_proj_eog(raw, raw_event=None, tmin=-0.2, tmax=0.2, projs, eog_events = _compute_exg_proj('EOG', raw, raw_event, tmin, tmax, n_grad, n_mag, n_eeg, l_freq, h_freq, average, filter_length, n_jobs, None, - reject, bads, avg_ref, no_proj, event_id) + reject, bads, avg_ref, no_proj, event_id, + eog_l_freq, eog_h_freq, tstart) return projs, eog_events -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-mne.git _______________________________________________ debian-med-commit mailing list [email protected] http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/debian-med-commit
