Compensation System for Biomagnetic Measurements with Optically Pumped Magnetometers inside a Magnetically Shielded Room

Magnetography with superconducting quantum interference device (SQUID) sensor arrays is a well-established technique for measuring subtle magnetic fields generated by physiological phenomena in the human body. Unfortunately, the SQUID-based systems have some limitations related to the need to cool them down with liquid helium. The room-temperature alternatives for SQUIDs are optically pumped magnetometers (OPM) operating in spin exchange relaxation-free (SERF) regime, which require a very low ambient magnetic field. The most common two-layer magnetically shielded rooms (MSR) with residual magnetic field of 50 nT may not be sufficiently magnetically attenuated and additional compensation of external magnetic field is required. A cost-efficient compensation system based on square Helmholtz coils was designed and successfully used for preliminary measurements with commercially available zero-field OPM. The presented setup can reduce the static ambient magnetic field inside a magnetically shielded room, which improves the usability of OPMs by providing a proper environment for them to operate, independent of initial conditions in MSR.


Introduction
Electrophysiological phenomena in the human body lead to the occurrence of magnetic field [1]. While electric biosignals are affected by the presence of insulating layers (fat, bones) and varied conductivity of the body tissues, magnetic field can easily penetrate through. Contactless measurements using magnetometers can avoid difficulties related to the electrode's attachment to the skin.
Magnetic fields generated by living organisms are of extremely low amplitude (fT -pT). Biomagnetic measurements can be thus affected by much higher electromagnetic noise or disturbing

Materials and Methods
Zero-field optically pumped magnetometers require a close to zero residual magnetic field. In commercially available OPMs from QuSpin (QZFM), three-axis integrated field cancellation coils are installed around the vapour cell and provide a field zeroing procedure [36]. The procedure is done automatically with the dedicated software at the beginning of measurements. The first generation (Gen-1) of QZFM can cancel the residual magnetic field of 50 nT, and thus an ambient magnetic field below this value is required for OPMs to operate. The measurement dynamic range of SERF OPMs is defined by the zero-field resonance and may differ depending on a manufacturer. In order to operate the SERF OPM in its linear regime, the dynamic range must be significantly smaller than the magnetic resonance width [48]. Full width at half maximum (FWHM) of the zero-field resonance in QuSpin OPMs (QZFM Gen-1) is~30 nT and electronic control unit output provides a measurement dynamic range of ± 5 nT after zeroing procedure [40]. The conceptual transfer curve of the sensor in the operating mode is presented in Figure 1A.
Sensors 2020, 20, x FOR PEER REVIEW 3 of 15 static field compensation system, based on the Helmholtz square coils can be easily installed inside MSR and its total cost does not exceed 200 Euros. The system is used to bring OPMs to their operational range, when the initial conditions in the shielded room are not sufficient.

Materials and Methods
Zero-field optically pumped magnetometers require a close to zero residual magnetic field. In commercially available OPMs from QuSpin (QZFM), three-axis integrated field cancellation coils are installed around the vapour cell and provide a field zeroing procedure [36]. The procedure is done automatically with the dedicated software at the beginning of measurements. The first generation (Gen-1) of QZFM can cancel the residual magnetic field of 50 nT, and thus an ambient magnetic field below this value is required for OPMs to operate. The measurement dynamic range of SERF OPMs is defined by the zero-field resonance and may differ depending on a manufacturer. In order to operate the SERF OPM in its linear regime, the dynamic range must be significantly smaller than the magnetic resonance width [48]. Full width at half maximum (FWHM) of the zero-field resonance in QuSpin OPMs (QZFM Gen-1) is ~30 nT and electronic control unit output provides a measurement dynamic range of ± 5 nT after zeroing procedure [40]. The conceptual transfer curve of the sensor in the operating mode is presented in Figure 1A. Unfortunately, ideal conditions are not easy to achieve in a man-size two-layer shielding room. There are few reasons for that. Firstly, an ambient magnetic field in the measurement volume can be initially higher than 50 nT. In this situation, the magnetic field must be compensated to make the sensor field zeroing procedure possible. Secondly, even after field zeroing process, the operating point can be close to the limit of compensation capabilities. Finally, the magnetic field can change during the measurement more than ± 5 nT dynamic range of the OPM, especially in long lasting sessions, and the operating point chosen initially may not be adequate anymore. These conditions, Unfortunately, ideal conditions are not easy to achieve in a man-size two-layer shielding room. There are few reasons for that. Firstly, an ambient magnetic field in the measurement volume can be initially higher than 50 nT. In this situation, the magnetic field must be compensated to make the sensor field zeroing procedure possible. Secondly, even after field zeroing process, the operating point can be close to the limit of compensation capabilities. Finally, the magnetic field can change during the measurement more than ±5 nT dynamic range of the OPM, especially in long lasting sessions, and the operating point chosen initially may not be adequate anymore. These conditions, assuming that the OPM can still operate, may lead to sensor overdrive and signal saturation, presented in Figure 1B. To overcome difficulties that may appear during measurements using zero-field magnetometers, we decided to design a three axis compensation coils system. The system was supposed to be installed inside the magnetically shielded room located at Warsaw University of Technology, Faculty of Mechatronics. The MSR shell consists of two layers of µ-metal and one layer of aluminium between them and has inner dimensions of 4.0 × 3.0 × 2.4 [m] (Vacoshield Advanced, Vacuumschmelze GmbH and Co. KG, Hanau, Germany [49]). The MSR in our facility allows for static shielding of the magnetic field only. It does not incorporate built-in compensating coils, so an active shielding is not possible without an additional compensating system. The attenuation level of the static magnetic field equals 54 dB, for the 0.01 Hz it is 32 dB, and for 0.1 Hz its value is at the level of 38 dB.
The Helmholtz coil, consisting of two identical coils, can produce a region of uniform magnetic field and thus can be used for cancelling an ambient magnetic field. A Helmholtz pair can be both circle and square. In both cases, coils with an equal number of turns are placed coaxially. The coils are connected in series, which provides the same electric current flow through both of them. The direction of the electrical current in the coils is set to generate the same direction of magnetic flux. Helmholtz spacing, a distance between the coils providing the homogeneity of the field, is equal to the radius of the coil for circle coils and equal to 0.5445 times the length of a side for square coils [50]. The uniformity of the field at the central volume of the coils is reached at the expense of variation in field strength at a distance from the centre. If slight nonuniformity at the centre is acceptable, different spacings are used to obtain an increase in the length of the quasi-uniform field [50]. Changing the coils spacing to, e.g., 1.15 times the Helmholtz spacing would result in an increase in the side-length of the quasi-uniform field volume from 20% to 60% of the coils' side length, at the expense of the magnetic field uniformity [50]. In this case, the magnetic field at the distance (measured from the centre of the quasi-uniform field volume) of 20% of the coils' side length would be~2% higher than the field measured in the centre. However, the magnetic field at the distance of 30% of the coils side length would be~2% lower than the field measured in the centre. In other words, across the above-mentioned distance equal 10% of the coils' side-length, gradient will be equal to~4% of the magnetic field measured in the centre of the quasi-uniform field volume. For 50 nT (the expected magnitude of the static magnetic field to be cancelled with the compensation system), this would result in a~2 nT (~4% of 50 nT) change across 20 cm (as 20 cm is 10% of designed 2 m square coils' side-length), which is equal~10 pT/mm. The magnetic field gradient in a two-layer MSR is typically lower than 30 pT/mm [39], and OPM requires a gradient not higher than 10 nT/3 mm = 3.3 nT/mm to operate properly. Thus, theoretically, even with the coil spacing providing less uniform magnetic field than in case of the Helmholtz spacing, the gradient produced by compensation coils should not deteriorate OPM performance significantly. Nevertheless, we decided to use the Helmholtz spacing in our system.
Square Helmholtz pair provides better uniformity of magnetic field in comparison to circular coils [51]. Furthermore, their shape is more useful in practice: the construction can stand alone and can be easily assembled and disassembled. The CAD model of the designed 3D Helmholtz coil system is presented in Figure 2. The CAD model was created using open-source Free CAD software [52].
For biomagnetic measurements (e.g., magnetoencephalography, magnetocardiography) a large volume of uniform magnetic field is required. Supporting system for the triaxial compensation coils was made of wood with 2-m-length of a side ( Figure 3). The size was limited by the dimensions of the entrance to the magnetically shielding room, but is still big enough for a human subject and can be easily accessed. The volume with the magnetic field homogeneity better than 0.1% is about 1/5 of the entire length [51], which provides a uniform volume with the side length equal to 0.4 m. Significant dimensions help to avoid the situation when movements of the subject could have led to touching of the coils frame, inducing variations in magnetic field generated by the compensating system. Each coil was wound with 30 and 60 turns of enameled wire independently. Such an approach brings some flexibility in terms of the range of the magnetic field it is possible to generate. The 30-turn windings (with higher wire cross-section area) could be used to work with higher currents or the 60-turn windings (with smaller wire cross-section area) while working with smaller currents. It is also possible to connect windings in series to obtain 90-turn windings. Double winding also gives the possibility to use one set of coils to compensate the residual magnetic field and the other one to introduce external signals. It is planned to use such a configuration in further research. For biomagnetic measurements (e.g., magnetoencephalography, magnetocardiography) a large volume of uniform magnetic field is required. Supporting system for the triaxial compensation coils was made of wood with 2-m-length of a side ( Figure 3). The size was limited by the dimensions of the entrance to the magnetically shielding room, but is still big enough for a human subject and can be easily accessed. The volume with the magnetic field homogeneity better than 0.1% is about 1/5 of the entire length [51], which provides a uniform volume with the side length equal to 0.4 m. Significant dimensions help to avoid the situation when movements of the subject could have led to touching of the coils frame, inducing variations in magnetic field generated by the compensating system. Each coil was wound with 30 and 60 turns of enameled wire independently. Such an approach brings some flexibility in terms of the range of the magnetic field it is possible to generate. The 30-turn windings (with higher wire cross-section area) could be used to work with higher currents or the 60-turn windings (with smaller wire cross-section area) while working with smaller currents. It is also possible to connect windings in series to obtain 90-turn windings. Double winding also gives the possibility to use one set of coils to compensate the residual magnetic field and the other one to introduce external signals. It is planned to use such a configuration in further research.
The 30-turn coils were made with 1-mm 2 cross-section area copper wire, resulting in a maximum resistance equal to 6 Ω per coil. This enables a high current flow of up to 1 A using 12 V power supply. The 60-turn coils were made with 0.1-mm 2 cross-section area wire and gave maximum resistance of 150 Ω per coil.  Compensation coils have to be driven with a stable current source. The dedicated circuit, based on a Howland current source, was designed and is presented in Figure 4. The source is batterypowered and provides filtration of the current supplying the coils to ensure current stability. The battery supply was used to avoid introducing unwanted noise from mains. The whole compensation system requires rather small current (20 mA), hence, while powered with two 9 V batteries, it can last for reasonably long time (~10 h). The current value can be controlled manually via potentiometer or voltage control (e.g., via acquisition card). The output may be set in the range from −5 to 5 mA The 30-turn coils were made with 1-mm 2 cross-section area copper wire, resulting in a maximum resistance equal to 6 Ω per coil. This enables a high current flow of up to 1 A using 12 V power supply. The 60-turn coils were made with 0.1-mm 2 cross-section area wire and gave maximum resistance of 150 Ω per coil.
Compensation coils have to be driven with a stable current source. The dedicated circuit, based on a Howland current source, was designed and is presented in Figure 4. The source is battery-powered and provides filtration of the current supplying the coils to ensure current stability. The battery supply was used to avoid introducing unwanted noise from mains. The whole compensation system requires rather small current (20 mA), hence, while powered with two 9 V batteries, it can last for reasonably long time (~10 h). The current value can be controlled manually via potentiometer or voltage control (e.g., via acquisition card). The output may be set in the range from −5 to 5 mA separately for each of three pairs of coils, resulting in an expected magnetic field strength of~500 nT at the center of Helmholtz coils, when using 90-turn windings.   [55]) at a sample rate equal to 5 S/s.
The system was designed to change the bias magnetic field in all directions in the range of ± 500 nT. The QuSpin OPMs Gen-1 provide the compensation of an ambient magnetic field up to 50 nT. The obtained range of the magnetic field adjustment gave a large margin for measurements, exceeding the field zeroing capabilities of OPMs.
To acquire the data from QuSpin sensors the NI-USB(BNC) 6218 DAQ (National Instruments Corporation, Austin, TX, USA [56]) card was used together with the dedicated application developed in LabView ® software. The raw data, four signals from two QuSpin OPM sensors representing the magnetic field in two orthogonal axes (Z and Y axes measured simultaneously), were sampled at 500 S/s with 16-bit resolution and streamed to the binary file (in big-endian order). The data stored in the file contain only fundamental information: sampling frequency (double, 8 Bytes), number of channels (single, 4 Bytes), data length (single, 4 Bytes) and finally the QuSpin data (double, 8 Bytes per sample, without multiplexing). The application, however, was provided with additional functionalities: real time power spectrum monitor, acquisition path evaluation (via signal power and signal to noise ratio calculation while testing signal applied), simple preview of signals and a timer The system was designed to change the bias magnetic field in all directions in the range of ±500 nT. The QuSpin OPMs Gen-1 provide the compensation of an ambient magnetic field up to 50 nT. The obtained range of the magnetic field adjustment gave a large margin for measurements, exceeding the field zeroing capabilities of OPMs.
To acquire the data from QuSpin sensors the NI-USB(BNC) 6218 DAQ (National Instruments Corporation, Austin, TX, USA [56]) card was used together with the dedicated application developed in LabView ® software. The raw data, four signals from two QuSpin OPM sensors representing the magnetic field in two orthogonal axes (Z and Y axes measured simultaneously), were sampled at 500 S/s with 16-bit resolution and streamed to the binary file (in big-endian order). The data stored in the file contain only fundamental information: sampling frequency (double, 8 Bytes), number of channels Sensors 2020, 20, 4563 7 of 14 (single, 4 Bytes), data length (single, 4 Bytes) and finally the QuSpin data (double, 8 Bytes per sample, without multiplexing). The application, however, was provided with additional functionalities: real time power spectrum monitor, acquisition path evaluation (via signal power and signal to noise ratio calculation while testing signal applied), simple preview of signals and a timer which could set the recording time and to visualize the time elapsed. The power spectrum monitor was also accompanied by a set of filters to choose from, which could remove the mains and its harmonics, the residual hum of the data acquisition card and the DC offset. This was intended to evaluate possibilities of signal quality improvement before the final recording. Only the raw data were recorded to preserve the data format simple and flexible in terms of further processing.

Compensation System Characterisation
The designed and built system was installed inside the two-layer magnetically shielded room at Warsaw University of Technology, Faculty of Mechatronics. To investigate the relationship between the magnetic induction generated by the Helmholtz coils and the current flow in wires, the fluxgate magnetometer Fluxmaster with measurement range of ±2 µT and 0.1 nT resolution was used [54]. The sensor was located in the centre of the setup and directed perpendicularly to each of three coils planes. The coils were supplied with the current in the range from −2.0 to 2.0 mA, which provided bias magnetic field changes in the range of ± 100 nT produced using only one set of coils (60 turns).
The relation of the current value vs. the value of the magnetic field was investigated. Measurements were carried out as follows: the current flowing through the compensation coils was changed in such a way as to set magnetic field inside measurement volume at predefined levels: from −100 to 100 nT in 10 nT steps. First the direction of changes was ascending, from −100 to 100 nT, then descending, from 100 to −100 nT. This procedure was repeated two times, giving four datapoints for each interval in total.
As the sensor may have a slight dependence on direction, the measurements were performed in two directions, positive and negative, with respect to the magnetic field produced by the compensation coils. By positive, we assume the sensor oriented in the direction of the magnetic field vector; negative relates to the reversed setting.
Results for each of three Helmholtz coils are presented in Figure 5. Measurements directions (x, y, z) correspond to the coordinate system in Figure 2.
Linear regression calculated using Matlab TM (The MathWorks, Inc., Natick, MA, USA [57]) polyfit function gave equation B[nT] = p 1 I[mA] + p 2 , where p 1 and p 2 are the polynomial coefficients and I[mA] is the current feeding the coils. Linear regression coefficients in descending powers (p 1 and p 2 ) with a standard error (σp1 and σp2) and determination coefficients R 2 calculated for the measurements done in each direction are presented in Table 1.
Approximate value of the residual magnetic field in the centre of the magnetically shielded room was determined from the linear fit of data. The calculation of the relation between the current and the magnetic field produced gave the value of the magnetic field ordinate for the current value set to zero. This value is approximately equal to the residual magnetic induction within the magnetically shielded room. The residual values for all three directions are as follows: x max~1 8 nT, y max~1 8 nT, z max~2 4 nT. The magnetic field variations were highest in the vertical (z) direction, reaching −30 nT on the other days.
As the sensor may have a slight dependence on direction, the measurements were performed in two directions, positive and negative, with respect to the magnetic field produced by the compensation coils. By positive, we assume the sensor oriented in the direction of the magnetic field vector; negative relates to the reversed setting.
Results for each of three Helmholtz coils are presented in Figure 5. Measurements directions (x, y, z) correspond to the coordinate system in Figure 2.

Preliminary Measurements
Preliminary measurements of the ability to compensate the magnetic field in the volume of interest with QuSpin zero-field magnetometers were made ( Figure 6). The optically pumped magnetometer was located in the centre of the compensation coils setup, perpendicularly to the plane of the green coils (−y direction) in Figure 2 (left and right vertical coils in Figure 3). Initially, measurements of the magnetic field were done without any additional active shielding (compensation coils switched off). Subsequently, the coils were switched on and the magnetic field was minimized manually, based on the The amplitudes of the residual magnetic field changes in the measurements done with and without additional active compensation (compensation coils switched on and off, respectively) cannot be compared, as the signal was saturated. Nevertheless, as no additional signal source was present inside MSR, the amplitude was expected to be similar. The difference is a result of shifted bias magnetic field.
To investigate the distribution of the measured residual magnetic field, histograms for both data series from Figure 6 were done (Figure 7). A bin width was set to 0.02 nT, which resulted in 66 bins for measurement done with the active shielding and 43 bins, when the coils were switched off. The number of bins already suggest a different distribution of the magnetic field. Histogram of measurement taken with the coils switched off (blue) is skewed with extreme dominant number of samples of the magnetic induction in the range from −0.08 to −0.06 nT. This corresponds to the areas of signal saturation. Even if the sensor was calibrated properly with the dedicated QuSpin software at the beginning of each measurement, the magnetic field induction measured without additional compensation shows multiple 'flat' areas (indicated in the picture with arrows, Figure 6A). In these areas, the magnetic field entering the centre of MSR is expected to be outside the compensation capabilities of used OPM. Use of the additional compensation coils changed the bias magnetic field, and signal saturation occurred only once.
The amplitudes of the residual magnetic field changes in the measurements done with and without additional active compensation (compensation coils switched on and off, respectively) cannot be compared, as the signal was saturated. Nevertheless, as no additional signal source was present inside MSR, the amplitude was expected to be similar. The difference is a result of shifted bias magnetic field.
To investigate the distribution of the measured residual magnetic field, histograms for both data series from Figure 6 were done (Figure 7). A bin width was set to 0.02 nT, which resulted in 66 bins for measurement done with the active shielding and 43 bins, when the coils were switched off. Histogram of data obtained with the additional compensation (red) is much more symmetric compared to the previous one (blue). This proves that the magnetic induction was distributed in the full range of OPM used. The sensor overdrive did not significantly influence the measurement and the data could be used for further analysis.
Comparison of datasets obtained with and without the additional compensation of magnetic field, indicate no usability of data measured without compensation system switched on in magnetography of human magnetic fields.
It is expected that the higher the current feeding the coils, the more magnetic field noise originating from the current supply. The noise inside MSR at Faculty of Mechatronics, WUT measured with QuSpin OPM in the band 3-100 Hz was as high as 1.5 pT/Hz for frequencies near DC, down to 50 fT/Hz for frequencies above 70 Hz on the day the preliminary measurements were taken. For frequencies below 10 Hz, the noise observed with the compensation system turned on was on the same level, or even lower, as when the system was turned off. For higher frequencies additional compensation resulted in the additional noise, whose level was significantly higher, reaching −300-400 fT/Hz in the worst case. Changing the ambient magnetic field inside MSR to make OPMs operational and the measurement possible is at the expense of a higher noise level.

Discussion
For proper operation, spin exchange relaxation-free, optically pumped magnetometers require a low magnetic field environment. Commercially available OPMs from QuSpin (QZFM Gen-1) can cancel an ambient magnetic field of magnetic induction below 50 nT, but this value is not easily achieved in the most common two-layer magnetically shielded rooms. The surroundings and Histogram of data obtained with the additional compensation (red) is much more symmetric compared to the previous one (blue). This proves that the magnetic induction was distributed in the full range of OPM used. The sensor overdrive did not significantly influence the measurement and the data could be used for further analysis.
Comparison of datasets obtained with and without the additional compensation of magnetic field, indicate no usability of data measured without compensation system switched on in magnetography of human magnetic fields.
It is expected that the higher the current feeding the coils, the more magnetic field noise originating from the current supply. The noise inside MSR at Faculty of Mechatronics, WUT measured with QuSpin OPM in the band 3-100 Hz was as high as 1.5 pT/

√
Hz for frequencies near DC, down to 50 fT/ √ Hz for frequencies above 70 Hz on the day the preliminary measurements were taken. For frequencies below 10 Hz, the noise observed with the compensation system turned on was on the same level, or even lower, as when the system was turned off. For higher frequencies additional compensation resulted in the additional noise, whose level was significantly higher, reaching −300-400 fT/

√
Hz in the worst case. Changing the ambient magnetic field inside MSR to make OPMs operational and the measurement possible is at the expense of a higher noise level.

Discussion
For proper operation, spin exchange relaxation-free, optically pumped magnetometers require a low magnetic field environment. Commercially available OPMs from QuSpin (QZFM Gen-1) can cancel an ambient magnetic field of magnetic induction below 50 nT, but this value is not easily achieved in the most common two-layer magnetically shielded rooms. The surroundings and maintenance of the MSR have a significant impact on the measurement and may lead to sensor overdrive and signal saturation.
The designed triaxial compensation system, based on square Helmholtz coils, enables generation of the magnetic field in the range from −500 to +500 nT, when supplied from the dedicated current source with the current in the range from −5 to 5 mA. The system was used to change the ambient magnetic field inside the magnetically shielded room and made measurements with OPMs possible, independently from the magnetic field variations. The presented setup improves the usability of zero field optically pumped magnetometers in standard MSR. The proposed solution can be used for decreasing an ambient magnetic field inside a magnetically shielded room to achieve a better environment for OPMs operation. While the presented setup was designed for commercially available QuSpin zero-field magnetometers (Gen-1) it could be extended to use with different OPMs operating in SERF regime to move them to their operational range. When OPMs' built-in cancelling coils are not capable to handle the residual magnetic field within a poorly shielded environment, our compensation system is capable to null up to few hundred nT of the static magnetic field to make OPMs operable. The introduced cost-efficient system can be used for the adaptation of two-layer MSRs to the new technology.
Preliminary results show an improvement in measurements using QZFM optically pumped magnetometers. A compensation of the static residual field of the shielded room with developed setup was possible and successfully used to overcome sensor overdrive. Cancelling an ambient magnetic field inside MSR, in order to make OPMs operational and the measurement possible, was at the expense of an increase in the noise level. The noise in the presented solution is non-negligible and an additional technique should be applied to perform biomagnetic measurements. The technique applied should depend on the application. Averaging, which is often used for evoked magnetoencephalography, is a common solution in recording low energy signals generated by living subjects. In other cases, filters should be applied, but again, this could be insufficient for a simple waveform recording. With two sensors, the measurement could be performed in a gradiometer configuration or an additional sensor could be used as a reference in an adaptive filtering technique. The size of the compensation system developed provides a uniform magnetic field in a volume with a side length equal to 0.4 m. The proposed construction is suitable for measurements of biomagnetic fields of human subjects.
Nevertheless, the manual control of the magnetic field compensation may be insufficient in long measuring sessions (e.g., evoked fields). Two-layer MSR provides a 32-38 dB attenuation of an external magnetic field for a frequency range from 0.01 to 0.1 Hz. A slowly changing magnetic field entering the magnetically shielded room may cause signal saturation, and even the working point could leave the range of the built-in OPM compensation. These problems could be avoided by using an automatic compensation system, which is a field of further development.