{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n\n# Annotate movement artifacts and reestimate dev_head_t\n\nPeriods, where the participant moved considerably, are contaminated by low\namplitude artifacts. When averaging the magnetic fields, the more spread the\nhead position, the bigger the cancellation due to different locations.\nSimilarly, the covariance will also be affected by severe head movement,\nand source estimation will suffer low/smeared coregistration accuracy.\n\nThis example uses the continuous head position indicators (cHPI) times series\nto annotate periods of head movement, then the device to head transformation\nmatrix is estimated from the artifact-free segments. The new head position will\nbe more representative of the actual head position during the recording.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Authors: Adonay Nunes <adonay.s.nunes@gmail.com>\n#          Luke Bloy <luke.bloy@gmail.com>\n# License: BSD-3-Clause\n# Copyright the MNE-Python contributors."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import mne\nfrom mne.datasets.brainstorm import bst_auditory\nfrom mne.io import read_raw_ctf\nfrom mne.preprocessing import annotate_movement, compute_average_dev_head_t\n\n# Load data\ndata_path = bst_auditory.data_path()\ndata_path_MEG = data_path / \"MEG\"\nsubject = \"bst_auditory\"\nsubjects_dir = data_path / \"subjects\"\ntrans_fname = data_path / \"MEG\" / \"bst_auditory\" / \"bst_auditory-trans.fif\"\nraw_fname1 = data_path_MEG / \"bst_auditory\" / \"S01_AEF_20131218_01.ds\"\nraw_fname2 = data_path_MEG / \"bst_auditory\" / \"S01_AEF_20131218_02.ds\"\n# read and concatenate two files, ignoring device<->head mismatch\nraw = read_raw_ctf(raw_fname1, preload=False)\nmne.io.concatenate_raws(\n    [raw, read_raw_ctf(raw_fname2, preload=False)], on_mismatch=\"ignore\"\n)\nraw.crop(350, 410).load_data()\nraw.resample(100, npad=\"auto\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Plot continuous head position with respect to the mean recording position\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Get cHPI time series and compute average\nchpi_locs = mne.chpi.extract_chpi_locs_ctf(raw)\nhead_pos = mne.chpi.compute_head_pos(raw.info, chpi_locs)\noriginal_head_dev_t = mne.transforms.invert_transform(raw.info[\"dev_head_t\"])\naverage_head_dev_t = mne.transforms.invert_transform(\n    compute_average_dev_head_t(raw, head_pos)\n)\nfig = mne.viz.plot_head_positions(head_pos)\nfor ax, val, val_ori in zip(\n    fig.axes[::2],\n    average_head_dev_t[\"trans\"][:3, 3],\n    original_head_dev_t[\"trans\"][:3, 3],\n):\n    ax.axhline(1000 * val, color=\"r\")\n    ax.axhline(1000 * val_ori, color=\"g\")\n\n# The green horizontal lines represent the original head position, whereas the\n# red lines are the new head position averaged over all the time points."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Plot raw data with annotated movement\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "mean_distance_limit = 0.0015  # in meters\nannotation_movement, hpi_disp = annotate_movement(\n    raw, head_pos, mean_distance_limit=mean_distance_limit\n)\nraw.set_annotations(annotation_movement)\nraw.plot(n_channels=100, duration=20)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "After checking the annotated movement artifacts, calculate the new transform\nand plot it:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "new_dev_head_t = compute_average_dev_head_t(raw, head_pos)\nraw.info[\"dev_head_t\"] = new_dev_head_t\nfig = mne.viz.plot_alignment(\n    raw.info,\n    show_axes=True,\n    subject=subject,\n    trans=trans_fname,\n    subjects_dir=subjects_dir,\n)\nmne.viz.set_3d_view(fig, azimuth=90, elevation=60)"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.12.2"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}