= MEG current phantom (Elekta-Neuromag) =
''Authors: Ken Taylor, John Mosher''

This tutorial explains how to import and process Elekta-Neuromag current  phantom recordings. We decided to release this example for testing and  cross-validation purposes. With these datasets, we can evaluate the  equivalence of various forward models and dipole fitting methods in the  case of simple recordings with single dipoles. The recordings are  available in two file formats (native and FIF) to cross-validate the  file readers available in Brainstorm and MNE. A similar page exists for the [[Tutorials/PhantomCtf|CTF phantom]].

<<TableOfContents(3,2)>>

== License ==
This  tutorial  dataset remains a  property of its authors: Ken Taylor, John Mosher (Epilepsy Center, Cleveland Clinic Neurological Institute, Cleveland, OH USA).

If  you reference this dataset in your publications, please acknowledge them and cite Brainstorm as indicated on the [[http://neuroimage.usc.edu/brainstorm/CiteBrainstorm|website]]. For questions, please contact us through the forum.

== The phantom ==
A current phantom is provided with the Elekta Neuromag for checking the system performance. It contains 32 artificial dipoles and four fixed head-position indicator coils. The phantom is based on the mathematical fact that an equilateral triangular line current produces equivalent magnetic field distribution to that of a tangential current dipole in a spherical conductor, provided that the vertex of the triangle and the origin of the conducting sphere coincide. For a detailed description of how the phantom works, see [[http://neuroimage.usc.edu/paperspdf/1985_IlmoniemHamKnuu_ForwardInversModel.pdf|here]].

The phantom dipoles are energized using an internal signal generator which also feeds the HPI coils. An external multiplexer box is used to connect the signal to the individual dipoles. Only one dipole can be activated at a time. The location of the dipole is recorded relative to the center of the  sphere (0,0,0)m, where X is positive toward the nasion, Y is positive  toward the left ear and Z is positive toward the top of the head (see  the [[CoordinateSystems]] tutorial for more details).

Use of the phantom is shown below. Note that the uncovered version is the phantom that came with the Neuromag-122, which explicitly shows the wiring. The covered version uses the same principle but somewhat different dipole locations. Further details are available in Section 7.2 of the User's Manual.

 . {{attachment:phantom.gif||height="366",width="632"}}

==== References ====
R.J. Ilmoniemi, M.S. Hämäläinen, and J. Knuutila, The Forward and Inverse Problems in the Spherical Model. In: Biomagnetism: Applications and Theory, eds. H. Weinberg, G. Stroink, T. Katil, Pergamon Press, 1985.

Elekta Neuromag System Hardware User's Manual, Revision G, September 2005.

== Description of the experiment ==
Files distributed as part of the CTF phantom tutorial:

 * phantom_ctf/ds/phantom_200uA_20150709_01.ds<<BR>>Current='''200uA''', Moment='''1800nAm''', Frequency='''7Hz''', Location='''(0, -18, 49)mm''' - 03-Jul-2015 <<BR>>This corresponds to a very strong dipole, that could be studied without any averaging.
 * phantom_ctf/ds/phantom_20uA_20150603_03.ds<<BR>>Current='''20uA''', Moment='''180nAm''', Frequency='''23Hz''', Location='''(0, -18, 49)mm''' - 03-Jun-2015<<BR>>This is a weaker dipole, closer to the range of amplitudes we can except from the brain. You will not see the dipole activity emerging from the noise without averaging a few cycles together.
 * phantom_ctf/ds/emptyroom_20150709_01.ds<<BR>>MEG empty room measurements: the phantom is in the MEG helmet but not connected - 09-Jul-2015
 * phantom_ctf/ds/phantom_20160222_01.pos<<BR>>"Head shape" of the phantom digitized with a Polhemus device using the [[Tutorials/TutDigitize|Brainstorm digitizer]].
 * phantom_ctf/fif/phantom_20uA_20150603_03.fif<<BR>>Dataset 20uA converted to FIF format using MNE utility mne_ctf2fiff
 * phantom_ctf/fif/phantom_200uA_20150709_01.fif<<BR>>Dataset 200uA converted to FIF format using MNE utility mne_ctf2fiff
 * All recordings were performed by Elizabeth Bock at the MEG lab at McGill with a CTF MEG system with 275 axial gradiometers.

== Download and installation ==
 * '''Requirements''':  You have followed the introduction tutorials and Brainstorm is installed.
 * Go to the [[http://neuroimage.usc.edu/bst/download.php|Download]] page of this website, and download the file: '''sample_phantom.zip'''
 * Unzip   it in a folder that is not in any of the Brainstorm folders (program   folder or database folder)
 * Start Brainstorm (Matlab scripts or stand-alone version)

Select the menu File > Create new protocol. Name it '''TutorialPhantomElekta''' and select the options:

 * '''No, use individual anatomy''',
 * '''No, use one channel file per acquisition run (MEG)'''.

== Generate anatomy ==
 * In the Matlab command window: type "'''generate_phantom''''''_elekta'''".

This creates a new subject '''Kojak_Sphere''' (so named after the 70's TV show) and generates the "anatomy" for this device: one volume and a few surfaces representing the geometry of the phantom. <<BR>><<BR>>

 . {{attachment:kojak_anat.gif||height="247",width="630"}} <<BR>><<BR>>

 * You can display the MRI and surfaces as presented in the introduction tutorials. <<BR>><<BR>> {{attachment:MRIandSurface.gif||height="427",width="630"}}

== Access the recordings ==
We can now review one of the raw kojak data sets. These have been generated by sequentially activating each of the 32 phantom dipoles in a single raw file.

 * Switch to the '''functional data''' view by clicking the second button under the protocol name.

 * Right-click on the subject folder > '''Review raw file'''.

 * Select the '''200nA''' source.

 * Click '''Event channel''', select '''STI201''', and hit '''OK'''.

Close the figure and double click the link to raw file to open a list of the events. At this point, you may wish to rename your groups so that the ordering remains convenient. To do this, click group 1, click '''Events''' > '''Rename group''', and rename it to '''01'''. Repeat this for 2 - 9.

This step, and others like it, can also be performed by importing the data into '''MATLAB''', running the following code below, and then exporting the data back into Brainstorm:

 .

 . '''          for''' i = 1:32,
  . STR = str2num(link_filename.F.events(i).label);
  . link_filename.F.events(i).label = sprintf('%02.0f',STR);
 '''          end'''

 * Trim away the singleton events by clicking on '''256 (x1)''', shift clicking on '''7936 (x1)''' and pressing the '''delete''' key.
 * Remove the first event from dipole '''01''' by clicking on it, clicking on the '''first''' time instant, and pressing the '''delete''' key.

Events should now appear as per the image below on the right: <<BR>>

 . {{attachment:editing_events.gif||height="444",width="630"}}

__'''Note:'''__''' '''Make sure to save the modifications that you make before proceeding. If you hit the grey X to close all figures and clear memory, you will be prompted to save modifications, and you should click "Yes".

== Importing the dipole events ==
Right click on '''Link to raw file''' and click '''Import in database'''

 * Uncheck '''Apply SSP/ICA projectors'''
 * Check '''Remove DC offset''', set the time range from -100 to 300ms
 * Resample the recordings to 100Hz (change from 1000Hz)
 * Do not create seperate folders
 * Click '''Import'''<<BR>><<BR>> {{attachment:kojak_dipole_import.gif||height="306",width="630"}} <<BR>><<BR>>

With the epochs loaded, we can now calculate the noise covariance from  the prestims. To do this, select all 32 dipoles (click the first then  shift click the last), then right click, '''Noise covariance''' > '''Compute  from recordings'''. Uncheck '''MEG MAG''' and hit '''OK'''. <<BR>><<BR>>

 . {{attachment:kojak_noise_covariance.gif||height="216",width="630"}} <<BR>><<BR>>

Average all of the epochs by selecting all 32 again, dragging them to  the process box, clicking '''RUN'''. Next click the cog, and select '''Average''' > '''Average files'''. Group the files '''By trial group (folder  average)''', and select '''Arithmetic average: mean(x)''', and click "'''Run'''".

== Source estimation ==
Compute a forward model and inverse model for a regular grid inside the phantom volume.

 * In Process1, select the two averages (20uA and 200uA). <<BR>><<BR>> {{attachment:phantom_inverse3.gif||height="102",width="362"}}
 * Select process '''Sources > Compute head model''':  MRI Volume, Regular grid (5mm), Single sphere.<<BR>>Increase the density of the grid for higher spatial accuracy.<<BR>><<BR>> {{attachment:phantom_headmodel1.gif||height="350",width="262"}}  {{attachment:phantom_headmodel2.gif||height="350",width="266"}}
 * Add process '''Sources > Compute sources [2015]''': Kernel only (one per file), Dipole modeling. <<BR>><<BR>> {{attachment:phantom_inverse1.gif||height="250",width="298"}}  {{attachment:phantom_inverse2.gif||height="322",width="227"}}

== Dipole scanning ==
Scan the for the most significant dipole in the grid of computed dipoles estimate previously.

 * In Process1, select the source maps for the two conditions (20uA and 200uA)<<BR>>or leave the averaged recordings selected and click on [Process sources]. <<BR>><<BR>> {{attachment:phantom_dipscan.gif||height="109",width="365"}}
 * Run process '''Sources > Dipole scanning''': Time window='''[0,0]ms'''<<BR>><<BR>> {{attachment:phantom_dipscan2.gif||height="247",width="298"}}  {{attachment:phantom_dipscan3.gif||height="242",width="282"}}
 * Right-click on the dipole file > '''File > View file contents''': <<BR>><<BR>> {{attachment:phantom_dipscan4.gif||height="232",width="471"}}
 * Note that there is no non-linear fitting in this process. This operation selects a dipole in the grid of points available in the grid with a 5mm spacing we used during the computation of the forward model. Therefore the localisation cannot be more precise than the resolution of the grid, all the results have to interpreted with an uncertainty of '''+/-2.5mm'''.
 * For a higher spatial resolution, you just need to recompute another forward model with a denser grid (2mm spacing for a 1mm precision).

== Dipole fitting with FieldTrip ==
Perform a non-linear dipole fit with the function ft_dipolefitting from the FieldTrip toolbox.

 * In Process1, select the average recordings for the two conditions (20uA and 200uA).
 * Run process '''Sources > FieldTrip: ft_dipolefitting''': Time window='''[0,0]''', Sensor type='''MEG'''<<BR>><<BR>> {{attachment:phantom_dipfit.gif||height="350",width="290"}}
 * Right-click on the dipole file > '''File > View file contents''': <<BR>><<BR>> {{attachment:phantom_dipfit2.gif||height="233",width="534"}}
 * The dipoles obtained with this non-linear dipole fitting method generally have a higher goodness of fit and a more accurate location than the dipole scanning previously illustrated. But this precision comes with computation costs that can be significantly higher in the case of more realistic head models.

== Results comparison ==
||<rowstyle="font-weight:bold;">Condition ||Method ||Forward model ||X ||Y ||Z ||GOF ||
|| ||Nominal location || ||0 ||-18 ||49 ||mm ||
||200uA ||Scanning (5mm) ||Single sphere ||-1.00 ||-16.00 ||44.00 ||99.80% ||
|| ||Scanning  (5mm) ||Overlapping spheres ||-1.00 ||-16.00 ||44.00 ||99.81% ||
|| ||Scanning  (5mm) ||OpenMEEG BEM ||-1.00 ||-16.00 ||44.00 ||99.76% ||
|| ||Fitting ||Single sphere ||-1.04 ||-17.00 ||43.98 ||99.95% ||
|| ||Fitting ||Overlapping spheres ||-1.05 ||-16.98 ||44.00 ||99.95% ||
|| ||CTF software || ||-1 ||-17 ||44 ||99.95% ||
|| ||MNE software || ||-0.79 ||-17.00 ||43.98 ||99.9% ||
|| || || || || || || ||
||20uA ||Scanning  (5mm) ||Single sphere ||-1.00 ||-16.00 ||44.00 ||96.94% ||
|| ||Scanning  (5mm) ||Overlapping spheres ||-1.00 ||-16.00 ||44.00 ||96.96% ||
|| ||Scanning  (5mm) ||OpenMEEG BEM ||-1.00 ||-16.00 ||44.00 ||96.90% ||
|| ||Fitting ||Single sphere ||-1.75 ||-17.13 ||44.39 ||98.25% ||
|| ||Fitting ||Overlapping spheres ||-1.78 ||-17.14 ||44.45 ||98.25% ||
|| ||CTF software || ||-1 ||-17 ||44 ||98.38% ||
|| ||MNE software || ||-1.38 ||-16.31 ||44.01 ||99.1% ||


The nominal location indicates where the dipole is supposed to be, relative to the center of the sphere. It is measured with rulers with an overall precision of about 5mm. The range of discrepancy we observe between this nominal location and the position of the dipole estimated from the recordings is acceptable.

<<TAG(Advanced)>>

== Digitized head points ==
The head points collected with the Brainstorm digitizer are usually copied to the .ds folders and imported automatically when loading the recordings. We decided not to include them in this example because in the case of this current phantom, there is no ambiguity in the definition of the anatomical fiducials. As this refined registration with the .pos files is not part of the standard CTF workflow, not including it will make it easier to compare the workflow and results with other programs.

For additional testing purposes, the .pos file for the phantom is included in the sample_phantom.zip package, but you have to add it manually to the recordings. Do not use these points to refine automatically the registration: the fitting algorithm may fail finding the best rotation around the Z axis because the phantom is completely spherical, and the registration is already close to perfection.

 * Right-click on one of the channel files (20uA or 200uA) > '''Digitized head points > Add points'''.<<BR>>Select the file format: "'''EEG: Polhemus'''"<<BR>>Select file: sample_phantom/'''phantom_20160222_01.pos'''<<BR>><<BR>><<BR>> {{attachment:phantom_pos.gif||height="165",width="608"}}
 * Right-click on the channel file > '''MRI registration > Check'''. <<BR>><<BR>> {{attachment:phantom_check.gif||height="215",width="543"}}

== Scripting ==
==== Generate Matlab script ====
Available in the Brainstorm distribution: '''brainstorm3/toolbox/script/tutorial_phantom_ctf.m'''