11institutetext: Department of Physics and Astronomy, FI-20014 University of
Turku, Finland; 11email<EMAIL_ADDRESS>and<EMAIL_ADDRESS>22institutetext:
Armagh Observatory & Planetarium, College Hill, Armagh BT61 9DH, UK; 22email:
<EMAIL_ADDRESS>33institutetext: Department of Physics &
Astronomy, University of Western Ontario, London, Ontario, N6A 3K7, Canada;
33email<EMAIL_ADDRESS>44institutetext: Leibniz-Institut für Sonnenphysik
(KIS), Schöneckstr 6, Freibirg, Germany 44email: svetlana.berdyugina@leibniz-
kis.de

# Highly sensitive search for magnetic fields in white dwarfs using broad-band
circular polarimetry

Andrei V. Berdyugin 11 Vilppu Piirola 11 Stefano Bagnulo 22 John D. Landstreet
2233 Svetlana V. Berdyugina 44

(Received September 15, 2021; accepted October 29, 2021)

Circular polarisation measurements of white dwarfs of various ages and
spectral types are useful to understand the origin and evolution of the
magnetic field in degenerate stars. In the latest stages of white dwarf
evolution, when stars are so cool that spectral lines are no longer formed in
the normal H- or He-dominated atmospheres, magnetic fields can be probed only
by means of circular polarimetry of the continuum. The study of the fields of
featureless DC white dwarfs may reveal whether Ohmic decay acts on magnetic
white dwarfs, or if magnetic fields continue to be generated even several
billion years after white dwarf formation. Compared to spectropolarimetry,
broad-band circular polarisation measurements have the advantage of reaching a
higher accuracy in the continuum, with the potential of detecting magnetic
fields as weak as a fraction of a MG in DC stars, if the telescope size is
adequate for the star’s magnitude. Here we present the results of a first
(short) observing campaign with the DIPol-UF polarimeter, which we have used
to measure broad-band circular polarisation of white dwarfs. Our observing run
was in part aimed to fully characterise the instrument, and in part to study
the relationship between magnetic field strength (when known from
spectropolarimetry) and circular polarisation of the continuum. We also
observed a small number of previously unexplored DC white dwarfs, and we
present the discovery of two new magnetic white dwarfs of spectral class DC,
probably the first discovery of this kind made with broad-band circular
polarimetric techniques since the late 1970s. We also discuss the
characteristics of our instrument, and predict the level of polarimetric
accuracy that may be reached as a function of stellar magnitude, exposure
time, and telescope size.

###### Key Words.:

White dwarfs – Stars: magnetic fields – polarization

##  1 Introduction

White dwarfs (WDs) are the end point of 90 % of stellar evolution. About 20 %
of such stars possess strong magnetic fields, but if we exclude stars younger
than 0.5 Gyr from the statistics, this frequency rises to almost 30 % (Bagnulo
& Landstreet, 2021). The distribution of the field strength is approximately
uniform per decade, in a range that varies from a few tens of KG to several
hundred MG, and there are no obvious signs of Ohmic decay during WD cooling,
at least until an age of 5 Gyr (Landstreet & Bagnulo, 2019; Bagnulo &
Landstreet, 2021). The origin of these fields is not fully understood. Several
hypotheses have been put forward: for instance, a field could be inherited
from pre-WD evolution, or originate during the merging of a binary system, or
it might be produced by a dynamo acting after crystallisation begins (Isern et
al., 2017). None of the proposed mechanisms is able to explain all observed
magnetic features of degenerate stars, and more than one channel of field
formation may exist.

Magnetic fields are primarily detected through the analysis of the Zeeman
effect on the Stokes profiles of spectral lines produced in the atmospheres of
almost pure H or He found in most WDs. Zeeman splitting may be detected with
high-resolution spectroscopy in DA, DB, and DZ WDs (with prominent H, He, and
sometimes Ca or Fe lines) that have field strengths as low as 50 kG. More
commonly, from low-resolution spectroscopy, the usual detection limit is of
the order of 1 MG. With spectropolarimetry, the detection threshold is 1 – 3
dex lower. Circular polarimetry is sensitive to the component of the magnetic
field along the line of sight averaged over the stellar disk. For this
quantity, usually called mean longitudinal field ($\langle B_{z}\rangle$),
values as small as a few kG may be detected in DA and DB stars (Aznar Cuadrado
et al., 2004; Landstreet & Bagnulo, 2019).

The situation is different in featureless (DC) stars, in which a magnetic
field may be detected only if its strength is sufficiently high to circularly
polarise the continuum. Assuming that a 15 MG field is needed to produce 1 %
of polarisation (Bagnulo & Landstreet, 2020), the weakest fields detected so
far in DC WDs are of the order of several MG. However, it is in principle
possible to firmly detect circular polarisation signals of the order of a few
parts in $10^{-4}$, or even smaller, and thus to detect magnetic fields in DC
WDs as weak as a fraction of a MG (if such fields exist in DC stars). Bagnulo
& Landstreet (2021) have suggested that the field strength may be roughly
uniformly distributed per decade in a range between a few tens of kG and few
hundred MG. Recalling that the magnetic frequency in older stars is of the
order of 30 %, we conclude that with appropriate exposure time and telescope
size, and with an instrument capable of delivering a polarimetric sensitivity
of order $10^{-4}$, one would expect to detect a magnetic field in at least
one out of five DC WDs.

Estimating the frequency of magnetic fields in DC WDs as a function of stellar
age would be very useful to establish whether Ohmic decay operates during DC
WD cooling, and if so, on which timescale; whether the field strength starts
to decrease with time after a certain age; and if there is an age at which
none of the DC stars are found to be magnetic. These data would help to set
strong constraints to the theories that try to explain the origin and
evolution of magnetic fields in degenerate stars.

The circular polarization of the continuum may be measured with a photo-
multiplier tube (PMT) or in imaging mode with a CCD detector, through narrow-
or broad-band filters, or even in white light without any filter at all; or
with spectropolarimetric techniques. Narrow- and broad-band circular
polarisation (BBCP) measurements were extensively used in the 1970s (e.g. Kemp
et al., 1970; Angel & Landstreet, 1970b; Landstreet & Angel, 1971) and the
early 1980s (Liebert & Stockman, 1980; Angel et al., 1981), but have never
been systematically exploited during the past four decades. However,
spectropolarimetry of the continuum for surveys of DC WDs has been used by
Putney (1997) and more recently, by Bagnulo & Landstreet (2020) and Bagnulo &
Landstreet (2021) for a systematic survey of DC WDs of the local 20 pc volume.

Spectropolarimetry potentially has a much higher diagnostic content than BBCP
because the detailed behaviour of the polarised spectrum with wavelength is
dictated by the field strength, atmospheric chemistry, and morphology.
Spectral resolution may be traded with the signal-to-noise ratio, S/N, via
software, allowing one in principle to reach the same S/N as in BBCP after
heavy re-binning, if this is necessary with faint stars. However, Bagnulo &
Landstreet (2020) and Bagnulo & Landstreet (2021) have argued that the
spectropolarimetry of the continuum is generally affected by higher
instrumental and background polarisation than BBCP, making it effectively very
difficult to assess the reliability of any signal smaller than $10^{-3}$, no
matter how high the S/N is. BBCP techniques have been proved to reach a
sensitivity better than $10^{-4}$ and should therefore be regarded as the
technique of choice for the detection of the weakest fields in DCs stars, but
also for larger surveys, to be followed up by spectropolarimetric monitoring
once a magnetic WD has been identified.

We have developed high-precision polarimeters for applications of broad-band
linear polarimetry, ranging from studies of the minute polarisation produced
by interstellar dust and the magnetic field in the solar neighbourhood to
strongly interacting binary stars with a BH component. Our Dipol-2 polarimeter
(Piirola et al., 2014) and the more recent version, DIPol-UF (Piirola et al.,
2021), achieve precision at the $10^{-5}$ level. We expect to reach a similar
performance level in circular polarimetry for sufficiently bright stars. With
this in mind, we initiated a pilot study to explore the capabilities of the
instrument in circular polarimetry of magnetic WDs, and the prospects for
applying the BBCP techniques in general. In this paper we discuss the results
of our study, and report new discoveries of previously unknown magnetic WDs.

##  2 Observations and data reduction

###  2.1 Instrument description

The observations were carried out at the 2.5 m Nordic Optical Telescope (NOT,
located at the Observatorio del Roque de los Muchachos, La Palma, Canary
Islands) on three nights, 2021 July 2-5, using the simultaneous three-colour
($B^{\prime}V^{\prime}R^{\prime}$) polarimeter DIPol-UF. The polarimeter is
described in detail by Piirola et al. (2021).

The light enters DIPol-UF through the modulator, a B. Halle superachromatic
quarter-wave plate (QWP) for circular polarisation measurements. The plate is
rotated to discrete positions with steps of 90$\degr$ for circular
polarimetry, and it remains at a fixed angle while each exposure is in
progress. The optical axis of the QWP is at 45$\degr$ from the polarisation
planes of the analyser. A plane-parallel calcite plate acts as a polarisation
analyser and splits the incident radiation into two orthogonally linearly
polarised rays (ordinary and extraordinary). With the image scale of the NOT,
the o- and e-image separation is about 14.6”. Two dichroic beam splitters
split the light into the three ($B^{\prime}V^{\prime}R^{\prime}$) spectral
passbands, with equivalent wavelengths of 445, 540, and 640 nm, and a full
width at half maximum, FWHM, of 114, 75, and 96 nm, respectively (Fig.1). We
emphasise that our $B^{\prime}V^{\prime}R^{\prime}$ bands are not those of the
Johnson-Cousins system, since broad overlapping passbands cannot be created by
the dichroic beam splitters. Instead, we use sharp-cutoff filters with high
peak transmission. The $V^{\prime}$ band is somewhat narrower than the other
bands due to the particular filter used for that band in addition to the
dichroic beam splitters.

Figure 1: DIPol-UF passbands ($B^{\prime}V^{\prime}R^{\prime}$) determined by
the filters and dichroic beam splitters, in the blue ($B^{\prime}$), green
($V^{\prime}$), and red ($R^{\prime}$) spectral region (continuous lines).
Dashed lines show the effect of the quantum efficiency of the CCD detectors
taken into account.

DIPol-UF uses three ultrafast electron-multiplied (EM) Andor iXon CCD cameras.
The $V^{\prime}$ and $R^{\prime}$ cameras are identical iXon Ultra 897-EX
models, while for the $B^{\prime}$ passband, the iXon Ultra 897-BV model with
the sensitivity optimised for the blue region is used. Faint targets are
recorded in the EM amplifier mode with the EM gain value in the range of 5–20.
In the EM mode, the camera can achieve single-photon sensitivity owing to the
elimination of the readout noise, which significantly increases the S/N. For
bright targets, conventional amplifiers with higher dynamical range are used
in combination with a defocusing technique. The minimum number of the optical
elements, high throughput, negligible systematic effects, and direct optical
elimination of sky polarisation make DiPol-UF a very efficient instrument for
broad-band polarimetry. The sensitivity of DIPol-UF in measuring polarisation
is limited in practice by photon statistics and depends on the number of
registered photons and the achieved S/N.

Each two successive exposures with orthogonal orientations of the QWP provide
one circular polarisation measurement. The uncertainty of circular
polarisation is the inverse of the S/N of the flux measured on both beams
accumulated over all exposures. For a general formulation, see for example
Patat & Romaniello (2006). For targets in the magnitude range 15 $<V<$ 17,
typically 20-40 sec individual exposure times were used. The delay between the
exposures due to camera readout and QWP rotation is $<0.3$ sec and is
therefore negligible. In a typical nightly monitoring interval of 0.5-1.5
hours, the number of individual measurements for one object is in the range 32
- 96. These data provide reliable statistical error estimates.

###  2.2 Target selection

Table 1:  Programme stars and their main physical features.

STAR | $G$ | $d$ | $T_{\rm eff}$ | log $g$ | $M$ | Age | Previous observations
---|---|---|---|---|---|---|---
|  |  | (pc) | (K) | c.g.s. | ($M_{\odot}$) | (Gyr) |
WD 1309$+$853 | G 256-7 | 15.8 | 16.5 | 5300 | 8.12 | 0.66 | 5.46 | DA MWD with $\langle|B|\rangle\simeq 5.4$ MG (7)
WD 1350$-$090 | GJ 3814 | 14.6 | 19.7 | 9580 | 8.13 | 0.68 | 0.81 | DA MWD with $\langle|B|\rangle\simeq 0.45$ MG (5,15)
WD 1503$-$070 | GD 175 | 15.8 | 20.3 | 6640 | 8.00 | 0.59 | 1.80 | DA MWD with $\langle B_{z}\rangle\sim 2$ MG (12)
WD 1556$+$044 | PM J15589+0417 | 16.0 | 22.5 | 6685 | 8.43 | 0.85 | 3.65 | DC (He) never observed
WD 1647$+$591 | DN Dra | 12.3 | 10.9 | 12738 | 8.24 | 0.76 | 0.45 | Non magnetic DA WD (13)
WD 1748$+$708 | G 240-72 | 13.8 | 6.2 | 5570 | 8.34 | 0.79 | 5.86 | Known strongly polarised DQ WD (2,4,9)
WD 1831$+$197 | G 184-12 | 16.3 | 39.3 | 7305 | 8.04 | 0.60 | 1.57 | DQ (He) no polarisation detected (8)
WD 1856$+$534 | LP 141-14 | 17.0 | 24.8 | 4430 | 7.69 | 0.40 | 4.34 | DC (He) never observed
WD 1900$+$705 | Grw +70∘ 8247 | 13.3 | 12.9 | 11880 | 8.54 | 0.93 | 0.92 | Known strongly polarised DA WD (1)
WD 1953$-$011 | LAWD 79 | 13.6 | 11.6 | 7868 | 8.23 | 0.73 | 1.63 | DA6 MWD ($\langle B_{z}\rangle\sim 0.1$ MG with spot) (6,10)
WD 2011$+$063 | G 24-9 | 15.7 | 22.9 | 6455 | 8.10 | 0.63 | 2.34 | DC (He) no polarisation detected (3,8)
WD 2032$+$248 | HD 340611 | 11.6 | 14.8 | 20704 | 8.03 | 0.64 | 0.06 | Non magnetic DA WD (13)
WD 2049$-$253 | UCAC4 325-215293 | 16.0 | 18.0 | 4895 | 7.84 | 0.48 | 4.4 | DC MWD (14)
WD 2049$-$222 | LP 872-48 | 15.0 | 20.3 | 8164 | 8.33 | 0.79 | 1.97 | DC (He) never observed
WD 2117$+$539 | EGGR 378 | 12.4 | 17.3 | 15250 | 7.94 | 0.58 | 0.17 | Non magnetic DA WD (13)
WD 2254$+$076 | G 28-27 | 17.1 | 45.5 | 12847 | 9.26 | 1.29 | $\sim 1.8$ | DA MWD with $\langle|B|\rangle\simeq 16$ MG (11)

111 Key to references: 1: Kemp et al. (1970); 2: Angel et al. (1974a); 3:
Angel et al. (1981); 4: West (1989); 5: Schmidt & Smith (1994); 6: Schmidt &
Smith (1995); 7: Putney (1995); 8: Putney (1997); 9: Berdyugin & Piirola
(1999); 10: Valyavin et al. (2008); 11: Külebi et al. (2009); 12: Gianninas et
al. (2011); 13: Bagnulo & Landstreet (2018), 14: Bagnulo & Landstreet (2020);
15: Bagnulo & Landstreet (2021).

Our target list comprises some WDs that are well studied in circular
polarisation, and some WDs that have never been observed before. Well-studied
targets include (1) WDs for which highly sensitive spectropolarimetric
measurements failed to detect any magnetic field, and which may therefore be
considered as non-magnetic standard WDs, (2) WDs that have relatively weak
fields (from hundreds to several thousand kG), and (3) WDs that are known to
exhibit a signal of BBCP. In addition to this group of 11 stars, we observed 5
DC WDs that were never observed before or in which we wished to confirm a
previous detection. Our programme stars are listed in Table 1, together with a
summary of each star’s main features.

###  2.3 Data reduction

Table 2:  Observing log of bright non-polarised standard stars. STAR | $G$ | DATE | UT | JD – | Exp. |  $V/I$ (%)
---|---|---|---|---|---|---
|  | yyyy-mm-dd | hh:mm | 2400000 | (s) | $B^{\prime}$ | $V^{\prime}$ | $R^{\prime}$
HD 115043 | 6.6 | 2021-07-02 | 21:19 | 59398.389 | 1280 |  $-0.000$$\pm$ | 0.001 |  $-0.001$$\pm$ | 0.001 |  $-0.000$$\pm$ | 0.001
HD 122676 | 6.9 | 2021-07-03 | 21:11 | 59399.383 | 1024 |  $0.001$$\pm$ | 0.001 |  $0.001$$\pm$ | 0.001 |  $-0.000$$\pm$ | 0.001
HD 124694 | 7.0 | 2021-07-04 | 21:19 | 59400.389 | 1152 |  $-0.003$$\pm$ | 0.001 |  $-0.000$$\pm$ | 0.001 |  $0.001$$\pm$ | 0.001
HD 211476 | 6.9 | 2021-07-03 | 05:00 | 59398.708 | 1600 |  $0.002$$\pm$ | 0.001 |  $0.001$$\pm$ | 0.001 |  $-0.000$$\pm$ | 0.001

After the standard CCD image reductions (bias, dark), the fluxes of the o- and
e-images were extracted with an aperture-optimising photometric algorithm that
minimises the standard error of the mean value of the circular polarisation
calculated from the individual measurements. The difference between the o- and
e-beam intensities at the orthogonal orientations of the QWP (Sect. 2.1) yield
one measurement of the Stokes parameter $V/I$, which is usually described as
the net percentage circular polarisation. Depending on the seeing and star
brightness, the optimum aperture radius was in the range 1.2” to 2.8”, and the
sky annulus thickness was 1.5” to 2.2”. To calculate the final value of $V/I$
from individual measurements, we applied the 2$\sigma$ weighting algorithm
(Piirola, V. et al., 2020), which gives a lower weight for points deviating by
$>2\sigma$. Normally, 5-8 % of the points deviated by more than 2$\sigma$ and
were given a weight $W<1$. The remaining points (92-95 %) were equally
weighted ($W=1$). The algorithm helps to suppress the effects of cosmic rays,
transient clouds, or occasional moments of poor seeing.

##  3 Interpretation of the observed circular polarisation

It is important to realise that the presence of a magnetic field in a WD may
lead to detectable broad-band circular polarisation through two quite
different mechanisms. One mechanism is related to the continuum, and the other
to the line opacity.

First, the spectrum of the star may include continuum radiation, outside even
of the wings of possible discrete spectral features, which is circularly (and
perhaps linearly) polarised. These polarised continua are expected to be
produced by the subtle dependence of the sources of continuum opacity in the
WD atmosphere on the magnetic field, so that the continuous opacity at a given
wavelength is greater in one sense of circular polarisation than in the other,
and the emergent continuum radiation has a net polarisation. The continuous
opacity in cool H- or He-dominated WD atmospheres is due to combinations of
opacity produced by free electrons, neutral atoms, negative atomic ions, and
molecular ions. The result is that continuum polarisation is expected to
depend on composition, $T_{\rm eff}$, and on the presence or absence of
metallic impurities. A particularly simple treatment of this basic idea is
described by Shipman (1971).

An example of the presence of true continuum polarisation is shown in Fig. 2
of Putney & Jordan (1995). In this figure, it is clear that in addition to the
circular polarisation spikes coinciding with the $\sigma$ components of
H$\alpha$, there is a baseline of non-zero continuum circular polarisation, in
spite of the absence of any absorption features other than the Zeeman triplet.
This baseline continuum polarisation is produced by the continuous opacity
effects mentioned above.

To provide a first estimate of magnetic field strength from detected non-zero
continuum circular polarisation, we often make use of a simple order-of-
magnitude estimate of the value of the line-of-sight average field $\langle
B_{z}\rangle$ that causes the observed polarisation. This estimate assumes
that the observed continuum value of $V/I$ is very roughly proportional to
$\langle B_{z}\rangle$, with a proportionality constant of about 15 MG per 1 %
circular polarisation. The origin of this estimate and derivation of the
proportionality constant are discussed in detail by Bagnulo & Landstreet
(2020) and references therein.

The resulting estimate of $\langle B_{z}\rangle$ can then be used to provide a
(still less accurate) estimate, or at least a lower limit, of the mean surface
field strength $\langle|B|\rangle$ by noting that $\langle|B|\rangle$ must be
larger than its line-of-sight component ($\langle|B|\rangle>\langle
B_{z}\rangle$), and we generally find $\langle|B|\rangle\geq 2.5\langle
B_{z}\rangle$ from experiments in which these two field strength values were
computed for various dipolar and dipole-like magnetic field structures.

A second mechanism that leads to non-zero circular polarisation in broad-band
filter measurements arises in MWDs that have strong polarised spectral
features, such as polarised Zeeman splitting in the Balmer lines. Two examples
are shown by Putney (1997) for GD 90 (with a field of
$\langle|B|\rangle\approx 9$ MG) in her Fig. 1h and for PG 1658+441 (with a
field of $\langle|B|\rangle\approx 2.3$ MG) in Fig. 1j. For both these stars,
it is clear that the Zeeman pattern of H$\alpha$ is practically symmetric, and
so in a broad-band measurement, the two $\sigma$ components will contribute
approximately equal numbers of polarised photons, which because of the
opposite signs of $\sigma_{+}$ and $\sigma_{-}$ will lead (approximately) to
no net polarised signal in the filter band.

However, at shorter wavelengths, it is no longer obvious that polarised counts
contributed by the blue and red polarised signals from the $\sigma$ components
of each line will cancel: the quadratic Zeeman effect means that the spacing
of various $\sigma$ sub-components is no longer the same in the blue and red
features, and line saturation effects will change the equivalent widths of the
polarised $\sigma$ components so that they no longer cancel out in the filter
measurement. In addition, the filter measurement may be made on a strongly
sloped flux continuum, such as those seen in the intensity spectra of the two
example stars, so that every blue-wing polarised feature will contribute more
photons than the associated red-wing feature. For this reason also, the blue
and red polarisation features will not exactly cancel, and the line spectrum
will contribute to the net polarisation measured in a broad band that includes
lines.

Thus we expect that polarisation in spectral line features may contribute some
net circular polarisation to a broad filter measurement. We note, however,
that while this effect means that a particular broad-band measurement may not
provide an accurate measurement of the continuum polarisation between spectral
lines, the line circular polarisation source still requires the presence of a
magnetic field to operate. A broad-band filter measurement of a non-magnetic
WD should yield a measurement without significant polarisation, and the
detection of a broad-band circular polarisation signal is still a strongly
valid sign of the presence of a magnetic field.

##  4 Results

Instrumental polarisation was determined from observations of four bright
nearby stars, assumed to have zero circular polarisation (Sect. 4.1 and Table
2). The results of our WD observations are given in Table 3 and are presented
in Sects. 4.2 to 4.5.

Table 3:  Observing log of WDs. WD 1900$+$705 and WD 1831$+$197 were observed twice. STAR | DATE | UT | JD – | Exp. |  $V/I$ (%)
---|---|---|---|---|---
| yyyy-mm-dd | hh:mm | 2400000 | (s) | $B^{\prime}$ | $V^{\prime}$ | $R^{\prime}$
WD 1647$+$591 | 2021-07-03 | 00:59 | 59398.541 | 1920 |  $-0.000$$\pm$ | 0.005 |  $0.007$$\pm$ | 0.007 |  $0.004$$\pm$ | 0.006
WD 2032$+$248 | 2021-07-04 | 05:00 | 59399.709 | 1920 |  $0.001$$\pm$ | 0.004 |  $-0.001$$\pm$ | 0.004 |  $-0.002$$\pm$ | 0.005
WD 2117$+$539 | 2021-07-03 | 04:17 | 59398.678 | 1800 |  $0.000$$\pm$ | 0.006 |  $-0.001$$\pm$ | 0.007 |  $-0.007$$\pm$ | 0.008
WD 1748$+$708 | 2021-07-03 | 00:10 | 59398.508 | 1920 |  $0.665$$\pm$ | 0.019 |  $-0.376$$\pm$ | 0.017 |  $-1.126$$\pm$ | 0.014
WD 1900$+$705 | 2021-07-02 | 22:22 | 59398.432 | 640 |  $3.756$$\pm$ | 0.016 |  $3.604$$\pm$ | 0.016 |  $3.827$$\pm$ | 0.019
WD 1309$+$853 | 2021-07-02 | 23:10 | 59398.466 | 3120 |  $0.090$$\pm$ | 0.048 |  $0.109$$\pm$ | 0.048 |  $0.113$$\pm$ | 0.037
WD 1350$-$090 | 2021-07-03 | 22:03 | 59399.419 | 2240 |  $0.048$$\pm$ | 0.015 |  $0.017$$\pm$ | 0.020 |  $0.011$$\pm$ | 0.018
WD 1503$-$070 | 2021-07-03 | 23:17 | 59399.471 | 3360 |  $0.006$$\pm$ | 0.030 |  $0.042$$\pm$ | 0.033 |  $0.061$$\pm$ | 0.025
WD 1953$-$011 | 2021-07-03 | 01:48 | 59398.576 | 2240 |  $0.017$$\pm$ | 0.012 |  $0.015$$\pm$ | 0.011 |  $-0.012$$\pm$ | 0.014
| 2021-07-04 | 01:55 | 59399.580 | 3200 |  $0.016$$\pm$ | 0.007 |  $0.002$$\pm$ | 0.008 |  $0.002$$\pm$ | 0.008
WD 2254$+$076 | 2021-07-05 | 04:04 | 59400.671 | 5760 |  $-0.347$$\pm$ | 0.034 |  $-0.324$$\pm$ | 0.056 |  $-0.529$$\pm$ | 0.055
WD 1556$+$044 | 2021-07-04 | 22:10 | 59400.424 | 3200 |  $-0.218$$\pm$ | 0.033 |  $0.315$$\pm$ | 0.039 |  $-0.131$$\pm$ | 0.030
WD 1831$+$197 | 2021-07-04 | 00:37 | 59399.526 | 3840 |  $-0.003$$\pm$ | 0.038 |  $0.049$$\pm$ | 0.042 |  $0.087$$\pm$ | 0.037
| 2021-07-04 | 23:48 | 59400.493 | 6400 |  $0.051$$\pm$ | 0.025 |  $-0.015$$\pm$ | 0.034 |  $-0.089$$\pm$ | 0.029
WD 1856$+$534 | 2021-07-05 | 01:33 | 59400.565 | 3840 |  $-0.067$$\pm$ | 0.058 |  $-0.005$$\pm$ | 0.060 |  $-0.032$$\pm$ | 0.043
WD 2011$+$063 | 2021-07-03 | 02:48 | 59398.618 | 3120 |  $0.022$$\pm$ | 0.035 |  $0.057$$\pm$ | 0.040 |  $-0.028$$\pm$ | 0.024
WD 2049$-$253 | 2021-07-04 | 03:31 | 59399.647 | 6400 |  $0.442$$\pm$ | 0.044 |  $0.460$$\pm$ | 0.039 |  $0.684$$\pm$ | 0.031
WD 2049$-$222 | 2021-07-05 | 02:39 | 59400.611 | 2240 |  $0.095$$\pm$ | 0.019 |  $0.068$$\pm$ | 0.022 |  $0.142$$\pm$ | 0.020

###  4.1 Observations of non-polarised bright stars

The high S/N measurements of non-polarised stars yield the instrumental
polarisation to a precision better than $10^{-5}$. In the
$B^{\prime}V^{\prime}R^{\prime}$ bands, the values of $V/I$ are 0.0117 $\pm$
0.0006 %, 0.0129 $\pm$ 0.0005 %, and 0.0084 $\pm$ 0.0004 %, respectively. The
instrumental polarisation was subtracted from the observed polarisation of all
targets, including the measurements of the standard stars in Table 2. The WDs
in which non-zero circular polarisation was detected by our measurements are
highlighted in boldface type in Table 3.

###  4.2 Known magnetic WDs with continuum polarisation

Figure 2:  $V/I$ spectropolarimetry of the continuum (empty symbols joined by
solid lines), and new DIPol-UF BBCP (solid circles) of the stars WD 1748$+$708
(top panel) and WD 1900$+$705 (bottom panel). Horizontal error bars represent
the FWHM of the broad-band filters. The polarisation spectra used for
comparison are described in the text.

WD 1748$+$708 and WD 1900$+$705 are well-known high-field MWDs that exhibit
strongly polarised continua. The polarisation of WD 1748$+$708 changes rapidly
with $\lambda$ (Angel et al., 1974b). In Fig. 2 we compare our filter
observations with an unpublished $V/I$ spectrum of WD 1748+708 obtained by
Angel and Landstreet on 1974 September 6.27 UT, using the Oke Multi-Channel
Spectrophotometer on the Mount Palomar 5m telescope, as described for example
by Landstreet & Angel (1975). We note that the sign of $V/I$ reverses within
the $V^{\prime}$ passband. Polarisation values observed in adjacent passbands
can be significantly different with broad-band filters. Accordingly, it is
essential to make observations in more than one passband to search efficiently
for magnetic WDs.

The spectropolarimetric data shown in the lower panel of Fig. 2 were obtained
with the ISIS instrument on the William Herschel Telescope (WHT). The plotted
data are adapted from the detailed search for possible variability of WD
1900+705 reported by Bagnulo & Landstreet (2019). We note that the broad
$B^{\prime}V^{\prime}R^{\prime}$ filters faithfully follow the broad outlines
of the higher-resolution $V/I$ spectra in both stars, clearly revealing the
power of filter polarimetry even for complex polarised spectra.

###  4.3 Non-magnetic WDs

Figure 3:  ISIS spectropolarimetry of the continuum rebinned to 112 Å in the
blue and 63 Å in the red (empty symbols joined by solid lines), and new DIPol-
UF BBCP (solid circles) of three non-magnetic WDs. Star WD 1647$+$591 was
observed three times with grating R600B, and observations obtained in
different epochs are shown with different colours. The other two stars were
observed simultaneously with gratings R600B and R1200R; see Bagnulo &
Landstreet (2018) for more details. Data for different stars have been offset
by 0, $-0.2,$ and $-0.4$ % for display purposes. Our BBCP measurements (zero
polarisation) are consistent with no detectable magnetic field.

Stars WD1647$+$591, WD2032$+$248, and WD2117$+$539 have previously been
observed with highly sensitive spectropolarimetric techniques. All were found
to be non-magnetic with sub-kG accuracy through the analysis of the Zeeman
effect on the Stokes $V$ profiles of their very strong H Balmer lines. All
these stars were observed by Bagnulo & Landstreet (2018) with the ISIS
instrument of the WHT. ISIS has the capability of measuring polarisation in
the continuum, but it has been found that the precision of these measurements
is limited by the difficulty of correcting precisely for spurious background
polarisation, particularly in measurements made in the presence of scattered
moonlight, as discussed in detail by Bagnulo & Landstreet (2020). This effect
means that it can be difficult to establish the reality of continuum
polarisation detected with spectropolarimetry below about the 0.05 % level. It
is therefore very interesting to compare the very weak continuum polarisation
signal observed with ISIS in these three non-magnetic WDs with the results of
our new BBCP measurements. This comparison is shown in Fig. 3. The figure
shows that for all non-magnetic WDs from the sample, DiPol-UF has measured
zero continuum polarisation ($<0.01\%$) in all three passbands. This null
result strongly confirms that our BBCP measurements are free of any systematic
or calibration errors at the level of $<10^{-4}$, and it also confirms the
earlier conclusion about the serious difficulty of measuring continuum
polarisation below about 0.05 % with a spectropolarimeter, especially in the
presence of moonlit background.

###  4.4 Weakly magnetic WDs

Our target list included five WDs that have been discovered to be magnetic
from the Zeeman splitting of their spectral lines. In increasing order of
field strength, these stars are WD 1350$-$090 and WD 1953$-$011, both with
$\langle B_{z}\rangle\simeq 0.15$ MG; WD 1503$-$070, with
$\langle|B|\rangle\sim 2$ MG; WD 1309$+$853, with $\langle|B|\rangle\sim 5.4$
MG; and WD 2254$+$076, with $\langle|B|\rangle\sim 16$ MG. These MWDs have
been included in the study in order to clarify the relation between known
magnetic field properties and observed continuum circular polarisation.

The only star in our sample of weakly magnetic WDs in which we firmly detect a
signal of circular polarisation is WD 2254$+$076, although marginal detections
are obtained in WD 1309$+$853 and WD 1350$-$090\. As WD 2254+076 has the
highest $\langle|B|\rangle$ value of the sample of low-field MWDs, this is not
surprising. The other more weakly magnetic stars in this sample probe the
limits of our ability to identify weakly magnetic WDs.

WD 2254$+$076 was discovered to host a magnetic field based on a flux spectrum
from the Sloan Digital Sky Survey (SDSS) showing Zeeman-split H Balmer lines.
Külebi et al. (2009) estimated a field of $\langle|B|\rangle\sim 16$ MG.
Assuming that the star has a mean longitudinal field between 0 and 10 MG, we
would expect to measure $|V/I|$ up to $\sim 0.7$ %. Earlier observations by
Angel et al. (1981), who had measured a BBCP signal of $-0.16\pm 0.10$ %,
failed to detect a magnetic field in this star, but our measurements
($-0.53\pm 0.06$ % in $R^{\prime}$, $-0.32\pm 0.06$ % in $V^{\prime}$, and
$-0.35\pm 0.03$ % in $B^{\prime}$) are consistent with our prediction and
would serve to definitively identify this WD as magnetic if this were not
already known.

WD 1309$+$853 = G 256–7 was discovered to be magnetic by Putney (1995). From
the spectra shown in her Fig. 2, we have deduced $\langle|B|\rangle\approx
5.4$ MG, while Putney found 4.9 MG. The polarisation spectrum reveals a non-
zero $\langle B_{z}\rangle$. From the observations of Putney (1995), we can
estimate $\langle B_{z}\rangle\sim 1-3$ MG. In this star, we have measured
$0.113\pm 0.037$ in the $R^{\prime}$ filter; in the $B^{\prime}$ and
$V^{\prime}$ filter we have measured similar values of the polarisation, but
only at $\sim 2\,\sigma$ significance level. Because all the Balmer lines in
this very cool DA star are extremely weak, all three filter measurements are
probably dominated by continuum polarisation, which appears to be roughly
constant with wavelength, as is found for the well-known similar DA star WD
0553–053 = G 99–47 (Angel & Landstreet, 1972; Liebert et al., 1975). Both
results are broadly consistent with the estimate that a field for which
$\langle B_{z}\rangle=15$ MG produces approximately 1 % of circular
polarisation in the continuum (Bagnulo & Landstreet, 2020, and references
therein), although more and better calibrations are needed to confirm the
precision and range of usefulness of this approximation.

A field characterised by $\langle|B|\rangle\sim 2$ MG was discovered in WD
1503$-$070 by Gianninas et al. (2011) from Zeeman splitting of the shallow and
weak Balmer lines. Our measurements set an upper limit of about 0.05 % for
circular polarisation. From the estimated $\langle|B|\rangle$ value, we would
expect $|\langle B_{z}\rangle|\la 1$ MG, that is, a signal of circular
polarisation $\la 0.07$ %. This estimate is consistent with our measurements.

In the star WD 1350$-$090, Schmidt & Smith (1995) measured $\langle
B_{z}\rangle=85\pm 9$ kG. Bagnulo & Landstreet (unpublished) have measured the
value of $\langle B_{z}\rangle$ several times, always with values between +105
and +145 kG, and Bagnulo & Landstreet (2021) have suggested that this field is
actually roughly constant with time. We measured $V/I=0.048\pm 0.015$ % in the
$B^{\prime}$ filter, but found values lower than 0.02 % in the other filters.
For such a small $\langle B_{z}\rangle$ value, we would expect a signal of
circular polarisation $\la 0.01$ %

However, this DA star has strong Balmer lines up to about H8 (where the higher
lines merge), and has correspondingly strong line polarisation features
throughout the $B^{\prime}$ filter. Thus we expect that the polarisation
features across individual lines may well have failed to cancel exactly,
either because of saturation effects, and/or because of the systematically
biased weighting given by the continuum slope. Therefore we could reasonably
expect some net polarisation precisely in the filter band that contains five
Balmer lines, but not in the other two filter bands. In addition, the sharp
cutoff of the $B^{\prime}$ filter is located in the blue wing of H8, and
therefore the measurement did not give full weight to the blue $\sigma$
polarisation from this line to cancel the polarised photons from the red
$\sigma$ component. Thus the polarisation detected in this band is almost
certainly not a measurement of true continuum polarisation, but is a reliable
symptom of the presence of magnetic field in this WD. In contrast, the null
polarisation measurements in $V^{\prime}$ and $R^{\prime}$ are consistent with
the prediction of continuum polarisation weaker than about 0.01 %.

WD 1953$-$011 has a constant longitudinal field $\la 0.1$ MG (Valyavin et al.,
2008) and should produce a BBCP signal also $\la 0.01$ %. Our measurements,
which failed to detect circular polarisation, rule out a signal significantly
stronger than this value.

###  4.5 Science targets

We have firmly detected a signal of circular polarisation in three further
WDs: WD 2049$-$253, WD 2049$-$222, and WD 1556$+$044\. The DC star WD
2049$-$253 was discovered to be magnetic by Bagnulo & Landstreet (2020), who
measured $V/I\simeq+0.5$ %, approximately constant between 370 and 510 nm,
with the FORS2 instrument in spectropolarimetric mode. The flux spectrum of
the star shows no obvious spectral line features, but with an effective
temperature slightly below 5000 K, it could have an atmospheric chemistry
dominated by either H or He. Although the spectropolarimetric measurement
showed no obvious problems and was obtained during dark time so that polarised
sky background due to cross-talk from highly polarised scattered moonlight was
not a problem, the first measurement is near the limit of what can reliably be
obtained with FORS. Bagnulo & Landstreet (2020) warned that the discovery
needed to be confirmed by further observations. Our new detection of a signal
at a similar level to that measured by Bagnulo & Landstreet (2020) fully
confirms that the star is indeed magnetic, with a field that is, as estimated
by Bagnulo & Landstreet (2020), probably characterised by $\langle
B_{z}\rangle\sim 7$ MG and so by $\langle|B|\rangle$ of about 20 MG or
slightly higher.

WD 1556$+$044 and WD 2049$-$222 are DC WDs that have never been observed in
polarimetric mode. They represent the first discovery of new MWDs obtained by
means of BBCP since the discoveries in the 1970s that were summarised by Angel
et al. (1981).

WD 2049$-$222 was first identified as a DC star in the survey of hot stars in
the Galactic halo of Bergeron et al. (1992). This classification was confirmed
by Kawka & Vennes (2006), who showed a portion of the spectrum between 6340
and 6780 Å. The completely smooth spectrum segment supports the classification
of this WD as a DC. There is a significant spread in the values of $T_{\rm
eff}$ that are found by different investigators, but they generally range
between 8000 and 9500 K. The absence of the slightest indication of H$\alpha$
in the spectrum segment published by Kawka & Vennes (2006) at this effective
temperature means that the dominant atmospheric chemistry is very probably He.
The polarisation signal that we have measured in WD 2049$-$222 is possibly the
weakest non-zero polarisation signal ever detected in a WD. It is weaker than
that measured in WD 0553$+$053 = G99-47 by Angel & Landstreet (1972), which is
$-0.4$ % at 4000 Å and $-0.2$ % at 8000 Å, and even weaker than that observed
in WD 1116$-$470 by Bagnulo & Landstreet (2021), which was just $\simeq-0.2$ %
in the entire range from 4000 to 9000 Å, but which in fact still awaits
confirmation.

The broadband polarisation measurements of WD 2049–222 all have the same sign,
but the variations appear to be real, with the maximum signal $V/I=0.142\pm
0.020$ % found in the $R^{\prime}$ band. For this star, we estimate $\langle
B_{z}\rangle=1-2$ MG. The fact that the polarisation in all three broad bands
has the same sign is typical of the polarisation spectra of most MWDs showing
non-zero continuum circular polarisation, from WDs with very strong continuum
polarisation, such as WD 1900+705 to the weak continuum polarisation of WD
0553+053. In this star, the available evidence suggests that the observed
polarisation is probably true continuum polarisation.

WD 1556$+$044 was originally thought to be a quasi-stellar object (Véron-Cetty
& Véron, 2010), but was recognised as a DC WD by Limoges et al. (2013), who
showed its blue $I$ spectrum to about 5200 Å. The spectrum is apparently
completely featureless even at high S/N, and with $T_{\rm eff}=6\,685$ K, the
star very probably has an He-rich atmosphere. Like many MWDs, the mass of this
WD, $0.85M_{\odot}$, is considerably higher than those of most of the stars in
Table 1. The polarisation in WD 1556+044 is detected at the $6\sigma$,
$8\sigma,$ and $4\sigma$ levels in the $B^{\prime}$, $V^{\prime},$ and
$R^{\prime}$ filters. In this new MWD, the detection is quite firm, and the
polarisation values, which change sign from $B^{\prime}$ to $V^{\prime}$ to
$R^{\prime}$, indicate strong variations of $V/I$ with wavelength. This is an
unusual characteristic that is found in only a small minority of MWDs for
which the spectrum of $V/I$ is well established. This unusual characteristic
is perhaps most similar to the variation seen in WD 1748+708 = G240-72 (Fig.
2), which is only about 1 000 K cooler. The explanation for the unusual sign
changes in these two MWDs is unknown, but may involve similar chemical
pollution of the nearly pure He atmospheres. This new MWD could have a
magnetic field of $\langle B_{z}\rangle\sim 5$ MG, and a field
$\langle|B|\rangle$ that is at least 12 or 15 MG, and perhaps larger.

We now comment on the remaining three science targets, for which we are able
to set only upper limits for $|V/I|$ and $|\langle B_{z}\rangle|$. Among these
targets, WD 1856$+$534 is of special interest because Vanderburg et al. (2020)
discovered a giant planet orbiting it that has not been tidally disrupted. The
star has never been observed in polarimetric mode, and our measurements
suggest that the star’s circular polarisation is most likely $<0.1$ % (in
absolute value), setting our estimate of the upper limit for $|\langle
B_{z}\rangle\ |$ to 1MG.

We observed WD 1831$+$197 twice. With our first observation, we obtained a
detection of nearly 3 $\sigma$ with the $R^{\prime}$ filter. Taken at a face
value, our second observation represents a 3 $\sigma$ detection (but with the
opposite sign with respect to the previous night). Owing to the lack of
detection in the other filters, we cannot yet accept the star as a MWD, but it
would be worthwhile to re-observe it with a larger telescope, or with a much
longer exposure time. The star could have a field of about 1 MG and might be
variable. At the same time, our observations probably rule out the presence of
a field with a longitudinal component larger than a few MG.

WD 2011$+$063 was one of the first WDs ever observed in polarimetric mode:
Angel & Landstreet (1970a) measured $V/I=0.12\pm 0.09$ %. Our new measurements
are about three times more accurate, but still detect no polarisation signal.
The upper limit for $|\langle B_{z}\rangle\ |$ is again of about 1 MG.

##  5 Polarimetric sensitivity and detection limits

Figure 4:  Relation between magnitude and uncertainty of our measurements of
Table 3 as a function of Gaia $G$ magnitude, normalised to 1 h exposure time.

Our observing run may be used to make quantitative predictions regarding the
level of polarimetric sensitivity that may be achieved for the BBCP
measurements as a function of time, stellar magnitude, and telescope size
using DIPol-UF with our filter set. Fig. 4 shows our measurement uncertainties
in the $B^{\prime}V^{\prime}R^{\prime}$ filters, normalised to 1 h exposure
time, as a function of the Gaia magnitude $G$. Uncertainties depend on the
spectral energy distribution of the star, that is, for hot stars, the
increased flux in the blue reduces errors in the $B^{\prime}$ band, whereas
for cool stars, the best S/N is achieved in the red. For our sample of stars
and their spectral energy dependence, the errors are rather uniformly
distributed (Fig. 4). The high throughput of the $B^{\prime}$ band (Fig.1)
helps to balance the performance of the DIPol-UF instrument in the different
passbands recorded simultaneously.

The relationship obtained for the polarimetric uncertainty as a function of
the input flux indicates that the former is defined by the photon flux with no
detectable systematic errors (both axes in Fig. 4 are logarithmic). This
nicely demonstrates that the DIPol-UF polarimeter provides polarimetric
measurements with precision and accuracy down to 10 parts per million (ppm)
also in circular polarisation, that is, it is comparable to the few ppm
uncertainty achieved in linear polarisation measurements (Piirola et al.,
2021). Thus, using this relation, we can evaluate observational requirements
for a given WD magnetic field with DIPol-UF employed at various telescopes.

From the estimate that a $\langle B_{z}\rangle=15$ MG field produces a
polarisation of approximately 1 % in the continuum (Bagnulo & Landstreet,
2020, and references therein), it follows that a signal of $\sim$ 0.07 % is
expected from a $\langle B_{z}\rangle=1$ MG field. A 3$\sigma$ detection of
this requires an uncertainty of $\sim$ 0.02%. According to the relation in
Fig. 4, this is achieved in one hour for a star with a magnitude of $G\sim
15.5$. For a brighter $G\sim 12$ star, the precision 0.004% in one hour would
yield a 3$\sigma$ detection of $\langle B_{z}\rangle\sim 0.2$ MG. For the
faintest ($G\sim$17) stars in our sample, the precision achieved in one hour
would allow $>3\sigma$ detection of fields $\langle B_{z}\rangle>3$ MG.

The typical seeing values (0.6” - 1.0”) and the sky darkness during our
observations are comparable to the conditions at the best observing sites, and
our detection limits could easily be scaled to larger telescopes. For an 8 m
telescope, the measured polarisation uncertainties would be lower by a factor
of $\sim 3$ for a given exposure time.

##  6 Conclusions

During the past two years, five new MWDs of spectral class DC have been
discovered: WD 0004$+$122, WD 078$-$670, WD 2049$-$253 (discoveries presented
by Bagnulo & Landstreet, 2020, with confirmation for the latter obtained in
this work), WD 1556$+$044, and WD 2049$-$222 (this work), while the discovery
of a $\sim 20$ MG field in WD 1116$-$47 by Bagnulo & Landstreet (2021) still
needs to be confirmed. These are the first discoveries of magnetic fields in
old, featureless WDs since the mid-1990s by Putney (1995). WD 1556$+$044 and
WD 2049$-$222 represent the first discoveries obtained by means of BBCP filter
measurements since the late 1970s. The significance is that we are now able to
explore the magnetic features of stars in the last stages of stellar
evolution. In particular, studying the frequency of the occurrence of magnetic
fields in old WDs will allow us to understand if and when Ohmic decay acts
during the long WD cooling phase, or if magnetic fields are found with high
frequency and high field strength even among the very oldest WDs of our
neighbourhood.

###### Acknowledgements.

Based on observations made with the Nordic Optical Telescope, owned in
collaboration by the University of Turku and Aarhus University, and operated
jointly by Aarhus University, the University of Turku and the University of
Oslo, representing Denmark, Finland and Norway, the University of Iceland and
Stockholm University at the Observatorio del Roque de los Muchachos, La Palma,
Spain, of the Instituto de Astrofisica de Canarias. DIPol-UF is a joint effort
between University of Turku (Finland) and Leibniz Institute for Solar Physics
(Germany). We acknowledge support from the Magnus Ehrnrooth foundation and ERC
Advanced Grant Hot-Mol ERC-2011-AdG-291659. JDL acknowledges the financial
support of the Natural Sciences and Engineering Research Council of Canada
(NSERC), funding reference number 6377-2016.

## References

  * Angel et al. (1981) Angel, J. R. P., Borra, E. F., & Landstreet, J. D. 1981, ApJS, 45, 457
  * Angel et al. (1974a) Angel, J. R. P., Hintzen, P., Strittmatter, P. A., & Martin, P. G. 1974a, ApJ, 190, L71
  * Angel et al. (1974b) Angel, J. R. P., Hintzen, P., Strittmatter, P. A., & Martin, P. G. 1974b, ApJ, 190, L71
  * Angel & Landstreet (1970a) Angel, J. R. P. & Landstreet, J. D. 1970a, ApJ, 162, L61
  * Angel & Landstreet (1970b) Angel, J. R. P. & Landstreet, J. D. 1970b, ApJ, 160, L147
  * Angel & Landstreet (1972) Angel, J. R. P. & Landstreet, J. D. 1972, ApJ, 178, L21
  * Aznar Cuadrado et al. (2004) Aznar Cuadrado, R., Jordan, S., Napiwotzki, R., et al. 2004, A&A, 423, 1081
  * Bagnulo & Landstreet (2018) Bagnulo, S. & Landstreet, J. D. 2018, A&A, 618, A113
  * Bagnulo & Landstreet (2019) Bagnulo, S. & Landstreet, J. D. 2019, MNRAS, 486, 4655
  * Bagnulo & Landstreet (2020) Bagnulo, S. & Landstreet, J. D. 2020, A&A, 643, A134
  * Bagnulo & Landstreet (2021) Bagnulo, S. & Landstreet, J. D. 2021, MNRAS, 507, 5902
  * Berdyugin & Piirola (1999) Berdyugin, A. V. & Piirola, V. 1999, A&A, 352, 619
  * Bergeron et al. (1992) Bergeron, P., Ruiz, M.-T., & Leggett, S. K. 1992, ApJ, 400, 315
  * Gianninas et al. (2011) Gianninas, A., Bergeron, P., & Ruiz, M. T. 2011, ApJ, 743, 138
  * Isern et al. (2017) Isern, J., Garsia-Berro, E., Kulebi, B., & Loren-Aguilar, P. 2017, ApJ, 836, L28
  * Kawka & Vennes (2006) Kawka, A. & Vennes, S. 2006, ApJ, 643, 402
  * Kemp et al. (1970) Kemp, J. C., Swedlund, J. B., Landstreet, J. D., & Angel, J. R. P. 1970, ApJ, 161, L77
  * Külebi et al. (2009) Külebi, B., Jordan, S., Euchner, F., Gänsicke, B. T., & Hirsch, H. 2009, A&A, 506, 1341
  * Landstreet & Angel (1971) Landstreet, J. D. & Angel, J. R. P. 1971, ApJ, 165, L67
  * Landstreet & Angel (1975) Landstreet, J. D. & Angel, J. R. P. 1975, ApJ, 196, 819
  * Landstreet & Bagnulo (2019) Landstreet, J. D. & Bagnulo, S. 2019, A&A, 623, A46
  * Liebert et al. (1975) Liebert, J., Angel, J. R. P., & Landstreet, J. D. 1975, ApJ, 202, L139
  * Liebert & Stockman (1980) Liebert, J. & Stockman, H. S. 1980, PASP, 92, 657
  * Limoges et al. (2013) Limoges, M. M., Lépine, S., & Bergeron, P. 2013, AJ, 145, 136
  * Patat & Romaniello (2006) Patat, F. & Romaniello, M. 2006, PASP, 118, 146
  * Piirola et al. (2014) Piirola, V., Berdyugin, A., & Berdyugina, S. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9147, Ground-based and Airborne Instrumentation for Astronomy V, ed. S. K. Ramsay, I. S. McLean, & H. Takami, 91478I
  * Piirola et al. (2021) Piirola, V., Kosenkov, I. A., Berdyugin, A. V., Berdyugina, S. V., & Poutanen, J. 2021, AJ, 161, 20
  * Piirola, V. et al. (2020) Piirola, V., Berdyugin, A., Frisch, P. C., et al. 2020, A&A, 635, A46
  * Putney (1995) Putney, A. 1995, ApJ, 451, L67
  * Putney (1997) Putney, A. 1997, ApJS, 112, 527
  * Putney & Jordan (1995) Putney, A. & Jordan, S. 1995, ApJ, 449, 863
  * Schmidt & Smith (1994) Schmidt, G. D. & Smith, P. S. 1994, ApJ, 423, L63
  * Schmidt & Smith (1995) Schmidt, G. D. & Smith, P. S. 1995, ApJ, 448, 305
  * Shipman (1971) Shipman, H. L. 1971, ApJ, 167, 165
  * Valyavin et al. (2008) Valyavin, G., Wade, G. A., Bagnulo, S., et al. 2008, ApJ, 683, 466
  * Vanderburg et al. (2020) Vanderburg, A., Rappaport, S. A., Xu, S., et al. 2020, Nature, 585, 363
  * Véron-Cetty & Véron (2010) Véron-Cetty, M. P. & Véron, P. 2010, A&A, 518, A10
  * West (1989) West, S. C. 1989, ApJ, 345, 511


# Comprehensive Study on Denoising of Medical Images Utilizing Neural Network
Based Auto-Encoder

Thoshara Nawarathne1 1Department of Electrical and Electronic Engineering,
University of Peradeniya, Sri Lanaka.  Thanushi Withanage1 1Department of
Electrical and Electronic Engineering, University of Peradeniya, Sri Lanaka.
Samitha Gunarathne1 1Department of Electrical and Electronic Engineering,
University of Peradeniya, Sri Lanaka.
Upekha Delay1 1Department of Electrical and Electronic Engineering, University
of Peradeniya, Sri Lanaka.  Eranda Somathilake2 1Department of Electrical and
Electronic Engineering, University of Peradeniya, Sri Lanaka.  Janith
Senanayake1 1Department of Electrical and Electronic Engineering, University
of Peradeniya, Sri Lanaka.
Roshan Godaliyadda1 1Department of Electrical and Electronic Engineering,
University of Peradeniya, Sri Lanaka.  Parakrama Ekanayake1 1Department of
Electrical and Electronic Engineering, University of Peradeniya, Sri Lanaka.
Janaka Wijayakulasooriya1 1Department of Electrical and Electronic
Engineering, University of Peradeniya, Sri Lanaka.
Chathura Rathnayake3 1Department of Electrical and Electronic Engineering,
University of Peradeniya, Sri Lanaka.

###### Abstract

Fetal motion discernment utilizing spectral images extracted from
accelerometric data incident on pregnant mothers abdomen has gained
substantial attention in the state-of-the-art research. It is an essential
practice to avoid adverse scenarios such as stillbirths and intrauterine
growth restrictions. However, this endeavor of ensuring fetus safety has been
arduous due to the existence of random noise in medical images. This novel
research is an in depth approach to analyze how the interference of different
noise variations affect the retrieval of information in those images. For
that, an algorithm employing auto-encoder-based deep learning was modeled and
the accuracy of reconstruction of the STFT images mitigating the noise has
been measured examining the loss. From the results, it is manifested that even
a substantial addition of the Super-Gaussian noises which have a higher
correlation of the frequencies possessed by the Fetal movement images can be
restored successfully with the slightest error.

###### Index Terms:

Fetal movement, Auto-encoder, Neural network, Denoising, STFT

##  I Introduction

Biomedical image processing is essential to extract the useful information
embedded in most of the biomedical signals. These processed images are useful
for doctors in identifying adverse health conditions [1]. One of the most
current areas of interest in medicine is ensuring the fetus’s health by
discerning the activeness of the baby inside the womb to avoid miscarriage and
intrauterine growth restrictions. Currently, there are several types of
research going on to detect fetal movements from the spectral images
constructed utilizing the signals obtained from the pregnant mother’s
abdomen[2] [3] [4] [5]. However, due to the methods utilized to obtain these
data, there’s a tendency of the original signals being interfered by external
noises such as the addition of mother’s movements, and the effect from the
materials used in the wearable device. This makes the fetal movement
differentiation problem extremely complex resulting in many false-positive
cases. Hence the difficulty in analyzing fetal well-being arises. Therefore,
denoising these raw images is extremely crucial as a preliminary step prior to
information extraction [6].

There are many denoising algorithms in the literature such as the Non-Local
Mean algorithm (NLM) [7] and Total Variation Regularization Algorithm [8], but
the most effective and recent method is Convolution Denoising Auto-Encoders
(CDAE) which is make use of convolution deep auto-encoder method. There are
several research attempts [9] [10] where they have used the auto-encoder
approach to denoise medical images. However those attempts were confined to
analyzing only some specific noise cases and focused on the structure of the
neural network. Hence not much significant research have been conducted to
comprehensively analyze the effect of all noise varieties in image denoising.

For this research first, we selected fetal signals with low noise. Initially
the data were in the time domain. By utilizing short Time Fourier Transform,
images were generated to notice both temporal and frequency domain variations
[11]. Thereafter, in the denoising task, we used an auto-encoder and analyzed
its performance in accurate de-nosing. Subsequently, the noises with different
noise factors were tested. Additionally, how the noise properties such as
colour, kurtosis, and skewness [12] affect the denoising process of the images
were observed. For this purpose, different noise distributions were appended
to the original image, and reconstruction loss was tracked through an auto-
encoder neural network.

##  II Objective

Current fetal movement detection methods include MRI and CT scans. There is a
risk that those may emit harmful radiations to the fetus [13][14]. We found a
data set that comprises accelerometric data obtained in a non-invasive manner
[15]. But it was observed that those signals had the tendency to get
interfered with external noises. In this research, the intention is to employ
an auto-encoder to denoise the images and to utilize this proposed algorithm
as a pre-processing step for the classification of the fetal movements.

##  III Data Set

For this research, a Biomedical data set named ”Fetal Movement Detection
Dataset Recorded Using MPU9250 Tri-Axial Accelerometer” was utilized. In it,
MPU 9250 tri-axial accelerometer data including gyration, external movements,
mothers’ responses have been recorded of 13 singleton mothers with an average
gestation period of 28 to 40 weeks. Each session contains an average of 7
fetal movements.

##  IV Methodology

This section will express the proposed scheme for biomedical image denoising.
It consists of pre-processing steps and auto-encoder based image denoising
methods as shown in Fig. 1.

Figure 1: Denoising Procedure

The data of the previously mentioned data set belong to the time domain. It
can be observed in Fig. 2 the time-domain interpretation of the fetal movement
signals and the impact of the fetal movement is most dominant in the Z-axis of
the accelerometer. Therefore, it was used for further denoising.

Figure 2: Fetal Movement Observation According to Various Axes of the Sensor

Fig. 3 shows the low noise and high noise fetal movement signal represented in
the Time Domain.

Figure 3: (a)Time Domain Low noise Image, (b)Time Domain Noisy image

The selected low noise time-domain signals transformed into frequency-time
variation by utilizing the Short Time Fourier Transform method.

###  IV-A Short Time Fourier Transform

The Short-time Fourier transform (STFT) is a subcategory of Fourier transform.
It is a particularly crucial section of the temporal - frequency distribution
of a signal. STFT can be considered as a method of identifying the complex
amplitude of a signal with respect to its time and frequency, as well as the
phase[16].

From an analytical perspective, the STFT of a signal is computed in two steps.
In the first stage, the total signal is divided into short time sections,
where the length of each section takes as a unique value, which refers to as
the window size. This window length should be selected accordingly to optimize
the output[17]. In the second step, Fourier transform of each shorter sections
were computed separately using equation(1). This explores the Fourier spectrum
of each small segment of the longer signal. Eventually, the visual
representation of the STFT is obtained, which interprets the segment-frequency
content versus time. This output plot is named Spectrogram.

$F(t^{\prime},u)=\sum_{t=-\infty}^{+\infty}f(t).W(t-t^{\prime}).e^{-j2\pi ut}$
(1)

where,

$F(t^{\prime},u)$ \- STFT of $f(t)$ :computed for each window

centered at $t=t^{\prime}$

$\ t^{\prime}$ \- Time parameter

$u$ \- Frequency parameter

$f(t)$ \- Signal to be analyzed

$W(t-t^{\prime})$ \- Windowing function centered at $t=t^{\prime}$

MATLAB was used to generate the STFT images for the above data. For this
purpose, initially, the raw time-domain signal was segmented by selecting the
sample size as 300, considering the average time elapse of the fetal movement
signal. Subsequently, the windowing method, window size, and overlap for the
spectrogram were selected adequately. The transformed sample spectrograms for
both low noisy and high noisy signals in Fig. 3 can be observed in Fig. 4.

Figure 4: (a)Spectrogram Low Noise Image, (b)Spectrogram Noisy Image

Then, the generated STFT images were resized and set as the inputs to the
auto-encoder.

###  IV-B Auto Encoders

An auto-encoder is an unsupervised deep artificial neural network, which can
be effectively utilized for data coding. It is equipped with propagation and
back propagation.

Learning of a representation of a particular set of data is the fundamental
purpose of an auto-encoder. It consists of both encoder and decoder. The
encoder side is particularly known as the reduction side. Specifically,
dimension reduction is achieved by training the network to filter the noise of
the original signal. This can be achieved through convolutional bottleneck
architecture[18]. On the other hand, on the decoding side, the auto-encoder
tries to generate the reduced encoded output, that merely matches with the
original input signal. In practice, there are two basic variations of auto-
encoders: Regularized Auot-encoders and Variational Auto-encoders. Under the
regularized auto-encoder sector it can be further categorized as Sparse auto-
encoders(SAE), Denoising auto-encoders (DAE), and Contractive auto-
encoders(CAE). This paper is based on the denoising auto-encoders.

Denoising auto-encoder focuses on the reduction of noise through identity
mapping. To use the auto-encoder to achieve the above-mentioned goal, it
should first train the network [19].

In this process, initially, the no noise (low noise) images were applied as
inputs ($X$). Then, noises were added to create noisy images ($X^{\prime}$).
After feeding through the denoising auto-encoder, the images can be
reconstructed such that minimum noise was obtained. Encoding algorithm is
represented in equation (2).

$Encoder:Y=A(WX^{\prime}+b)$ (2)

$X^{\prime}$ images were taken by the encoder and it constructed the output
$Y$. Here, $A$ represents the activation function, and $W$ is the weight added
to the hidden layer(s). Besides, a bias $b$ is also added in there.

On the decoder side, the latent representation of $Y$ was mapped again to
reconstruct the images $Z$, which were assumed to be low noisy. The decoder
was modeled as stated in equation (3).

$Decoder:Z=A(W^{\prime}Y+b^{\prime})$ (3)

Where $W^{\prime}$ and $b^{\prime}$ are the weights and the biases added in
the decoder respectively.

Here the auto-encoder parameters such as $W$, $W^{\prime}$, $b$, and
$b^{\prime}$ were tuned to minimize the cost function. It was achieved by
varying the parameters such as the number of layers, number of neurons per
layer, and neuron size.

The performance of denoising auto-encoder was analyzed by varying the fraction
of noise added and characteristics of the noise: whiteness and distribution
etc.

Figure 5: Auto-encoder Architecture of Denoising Images

This encoder-decoder denoiser architecture can be observed in Fig. 5. To feed
the network, initially, 75 low noisy fetal movement signals were selected.
After transforming them into spectrograms those were resized to create images
with size $256\times 256\times 3$. Then the noise addition was done in four
separate sections as listed below.

  1. 1.

Gaussian White Noise Performance

  2. 2.

Gaussian Coloured Noise Performance

  3. 3.

Optimum Colouring Parameter Detection

  4. 4.

Distribution Variation

The used auto-encoder architecture was a sequential model, which was analyzed
using Tensor Flow libraries. This model consisted of two convolutional layers
for the encoder and two for the decoder. In the encoder section, the first
convolutional layer used 32 latent variables with a kernel size of $3\times
3$, with Relu activation. This was followed by the max-pooling layer with a
window size of $2\times 2$. The next hidden layer had 64 neurons, kernel size
of $3\times 3$, with a Relu activation and a max-pooling layer with the same
dimensions as above. Decoder had the parameters with 64 and 32 feature maps,
$3\times 3$ filter size, with $2\times 2$ max-pooling layers activated by Relu
function.

The performance evaluation according to the different noise types and
parameters is discussed in further sections.

####  IV-B1 Gaussian White Noise Performance

To analyze the auto-encoder denoising performance concerning different
characteristics of added noise, initially, a Gaussian-white noise signal was
added. The addition of Gaussian white noise signal represents an ideal
situation when the image gets interfered with noises at every frequency in an
equal power level and the amplitude distribution of those noises represents a
normal distribution [20]. Here, the amount of noise addition was varied by
changing the Noise Factor parameter from 0 to 0.9 by 0.1 intervals, where a
Noise Factor of 0.9 means the images are almost completely noisy. This
relationship is shown in equation (4). The objective of this section was to
identify the optimum noise factor such that it results in the minimum loss of
the algorithm.

$Noise\hskip 2.84526ptFactor=\dfrac{Noise\hskip 2.84526ptPower}{Signal\hskip
2.84526ptPower}$ (4)

####  IV-B2 Gaussian Coloured Noise Performance

The nature of the input fetal movement signals is not white. Therefore, it is
better to explore the denoising algorithm to coloured noise components. To
achieve that goal colouring of the Gaussian-white signal was done using
applying the filter depicted in Fig. 6. By applying that filter transformation
to white Gaussian noise, some noise frequency components are suppressed, while
changing the overall shape of the power spectral density of the Gaussian noise
distribution.

Figure 6: Filter Architecture Used for Colouring Noise

$G(Z)=\dfrac{\sqrt{1-a^{2}}}{1-aZ^{-1}}$ (5)

####  IV-B3 Optimum Colouring Parameter Detection

Next, our approach was to detect the optimum alpha in which the proposed auto-
encoder network presents the highest accuracy in reconstructing the noisy
spectral fetal movement images. So that such types of noises included in the
medical images could be effectively denoised. For that, the coloring
parameter, ‘$a$’ in equation (5) was changed from 0 to 0.9 in steps of 0.1 to
generate Gaussian signals of different color levels. These signals are shown
in Fig. 7. The minimum loss was examined considering both various ’$a$’ values
and noise factors.

Figure 7: Generated Gaussian Noise for Different Colour Levels

####  IV-B4 Distribution variation

Up to now, the above analysis was done considering only normally distributed
signals. However, the added noise for bio-medical signals may be differing
from Gaussian signals. Therefore, it is necessary to explore the behavior of
the auto-encoder denoiser when the input images get added with non-Gaussian
noise. Those non-Gaussian signals are further categorized as Sub-Gaussian and
Super-Gaussian signals, which were generated by changing the features such as
kurtosis and skewness of the noise distributions[21]. Hence, some standard
distributions of varied kurtosis and skewness were used for this purpose.
Uniform, Raised Cosine, and Weiner Semi-Circle distributions were used for
exploring the Sub-Gaussianity nature whereas double exponential, Weibull,
Hyperbolic Secant, Laplace, and Logistic distributions were used to analyze
Super-Gaussian nature. Fig. 8 and Fig. 9 depicts the selected Super-Gaussian
and Sub-Gaussian signal distributions respectively.

Figure 8: Selected Super-Gaussian Distributions

Figure 9: Selected Sub-Gaussian Distributions

For the above distributions, parameters such as the noise factor and the
coloring parameter were tuned. By analyzing the minimum loss in denoising, the
optimum noise value and coloured factor, for each above mentioned
distribution, were obtained.

##  V Results

The minimum loss of the algorithm was varied with the noise factor when the
White Gaussian noise was added to the input. Initially we try to optimize the
minimum loss by varying Noise Factor as well as Number of Epochs. From this
analysis, we obtained that there is a significant drop of loss around Noise
Factor = 0.3 for 10 to 50 Epochs. This variation is depicted in Fig. 10.

Figure 10: Variation of Normalized Loss with Noise Factor and Number of Epochs

Considering the computational cost we selected to use the Number of Epochs as
30. For Number of Epochs = 30 acheivable minimum loss with respective to Noise
Factor is shown in Fig. 11. As it exhibits there, the finest denoising
performance can be obtained when the N.F. is about 0.3 for white Gaussian
noise.

Figure 11: Variation of Minimum Loss vs. Noise Factor for White Gaussian Noise

Then as the analysis extended to checking the performance in colored gaussian
signals the same procedure was carried out by varying the N.F. for different
’$a$’ values. The observed variations are shown in Fig. 12.

Figure 12: Variation of Loss vs. Colouring for Gaussian Noise for N.F. = 0.3

Figure 13: Variation of Loss vs. Noise Factor for Different ’$a$’ Values

In further analysis, it is manifested that, here regardless of the value of
the coloring parameter of the noise signal, the optimum noise factor was also
obtained around 0.3. The obtained minimum loss variation with the noise factor
for sample ‘$a$’ values can be observed in Fig. 13.

The obtained minimum loss was varied when varying the distribution of the
noise. When comparing the Sub-Gaussian and Super-Gaussian white signals,
Super-Gaussian signals perform well in the denoising of the biomedical images.

Even though the added noise proportion was quite high, most of the Super-
Gaussian noise components can be eliminated in a fairly accurate manner than
Sub-Gaussian noise. These results are tabulated in Table I and sample
behaviour for sub and super Gaussian distribution is depicted in Fig. 14.

TABLE I: Noise Factors, Which Obtained Minimum Loss (0.0032) of White Gaussian Noises with Different Probability Distributions Gaussianity | Distribution | Noise | Next Minimum
---|---|---|---
|  | Factor | Loss with N.F.
Sub-Gaussian | Raised Cosine | 0.1 | 0.0034 at N.F. = 0.3
Sub-Gaussian | Uniform | 0.2, 0.3 | 0.0033 at N.F. = 0.5
Sub-Gaussian | Wiener Semi-circle | 0.1, 0.3 | 0.0034 at N.F. = 0.4
Super-Gaussian | Double Exponential | 0.2, 0.3 | 0.0034 at N.F. = 0.8
Super-Gaussian | Hyperbolic Secant | 0.2, 0.4 | 0.0033 at N.F. = 0.6
Super-Gaussian | Laplace | 0.2, 0.3 | 0.0034 at N.F. = 0.8
Super-Gaussian | Logistic | 0.3, 0.4 | 0.0033 at N.F. = 0.5
Super-Gaussian | Weibull | 0.2, 0.5 | 0.0034 at N.F. = 0.7

Figure 14: (a)Minimum Loss for Sub-Gaussan Noise, (b)Minimum Loss for Super-
Gaussan Noise

If the input noise is colored, and its distribution was varied the obtained
minimum loss also varies. Their performances are shown in Table II.

TABLE II: Minimum Average Loss for Different Distributions with Optimum Colouring Parameter Gaussianity | Distribution | Colouring | Average
---|---|---|---
|  | Parameter | Minimum Loss
Sub-Gaussian | Raised Cosine | 0.5 | 0.0036
Sub-Gaussian | Uniform | 0.7 | 0.0038
Sub-Gaussian | Wiener Semi-circle | 0.5 | 0.0037
Super-Gaussian | Double Exponential | 0.4 | 0.0036
Super-Gaussian | Hyperbolic Secant | 0.8 | 0.0034
Super-Gaussian | Laplace | 0.7 | 0.0038
Super-Gaussian | Logistic | 0.7 | 0.0036
Super-Gaussian | Weibull | 0.6 | 0.0036

Also, for a sample Gaussian noise with sub-gaussian distribution the loss
variation is shown in Fig. 15, whereas Fig. 16 depicts the loss varying
behavior of a sample Super-Gaussian distributed Gaussian signal.

Figure 15: Variation of Minimum Loss and N.F. vs. Colouring Parameter for a
Noise with Uniform Distribution

Figure 16: Variation of Minimum Loss and N.F. vs. Colouring Parameter for a
Noise with Double Exponential Distribution

##  VI Conclusion

Image denoising is crucial when analyzing the bio-medical images. Here to
identify the fetal movement signals with considerable accuracy it is better to
use a denoising method as a pre-processing step before feeding into a
classifier. In this research denoising method based on convolutional auto-
encoders was applied to spectral images, and its performance was evaluated.
From the above analysis, it can be observed that when Gaussian-White or
Gaussian-coloured noise was added the optimum denoising was obtained when
Noise Factor is 0.3. Further, the optimum coloring parameter of the denoiser
was alpha = 0.3.

Thereafter, varying the noise distribution from Gaussian to Sub/Super-Gaussian
noise signals with optimum noise factor varied gives promising results.
Specifically for Super-Gaussian signals, even the noise factor more than
0.7-0.8 offers the finest denoising performances displaying high accuracy for
both white and coloured noise components. This may result since most of the
noise signals in nature vary in a smooth manner. So the fetal signals impacted
by these soft variations can be modeled using the Super-Gaussian signals.
Therefore, to denoise the fetal movement biomedical signals it is ideal to use
the above proposed convolutional auto-encoder denoiser.

## References

  * [1] Rådestad, Ingela. (2010). Fetal movements in the third trimester - Important information about wellbeing of the fetus. Sexual & reproductive healthcare : official journal of the Swedish Association of Midwives. 1. 119-21. 10.1016/j.srhc.2010.06.006.
  * [2] W. A. W. S. Wasalaarachchi et al., ”Fetal Movements Identification Based on Non-negative Matrix Factorization and Spectral Clustering,” 2019 14th Conference on Industrial and Information Systems (ICIIS), Kandy, Sri Lanka, 2019, pp. 266-271, doi: 10.1109/ICIIS47346.2019.9063351.
  * [3] U. H. Delay, B. M. T. M. Nawarathne, D. W. S. V. B. Dissanayake, M.P.B. Ekanayake, G.M.R.I Godaliyadda, J. V. Wijayakulasooriya, R. M. C. J. Rathnayake. “Non Invasive Wearable Device for Fetal Movement Detection” 2020 15th International Conference on Industrial and Information Systems (ICIIS), PaperID1570668829 Online Conference held on:11.28.2020 (IN PRESS)
  * [4] W.T. Ruwanga, G.M.R.I. Godaliyadda, M.P.B. Ekanayake, J.V. Wijayakulasooriya, and R.M.C.J.Rathnayake, ”Combining Local Non-negative Matrix Factorization and Spectral Clustering to De-tect Fetal Movements”, Women in Engineering Symposium 2020 (IEEE WIESymp2020),PaperID22, Online Conference held on:11.15.2020.
  * [5] Samitha Gunarathne, Thanushi Withanage, Upekha Delay, Thoshara Nawarathne, Roshan Go-daliyadda, Parakrama Ekanayake, Janaka Wijayakulasooriya, and Chathura Rathnayake,”StandardizedFourth Central Movement Based Three Step Algorithm for Fetal Movements Identification”. On-line International Conference on Advances in Computing and Technology (ICACT 2020 Universityof Kelaniya), PaperID 51, Online Conference held on:11.28.2020 (IN PRESS)
  * [6] Huang HK. Biomedical image processing. Crit Rev Bioeng. 1981;5(3):185-271. PMID: 7023828.
  * [7] K. Leng, ”An improved non-local means algorithm for image denoising,” 2017 IEEE 2nd International Conference on Signal and Image Processing (ICSIP), Singapore, 2017, pp. 149-153, doi: 10.1109/SIPROCESS.2017.8124523.
  * [8] Paul Rodríguez, ”Total Variation Regularization Algorithms for Images Corrupted with Different Noise Models: A Review”, Journal of Electrical and Computer Engineering, vol. 2013, Article ID 217021, 18 pages, 2013. https://doi.org/10.1155/2013/217021.
  * [9] D. Lee, S. Choi, H.-J. Kim, Performance evaluation of image denoising developed using convolutional denoising auto-encoders in chest radiography, Nuclear Inst. and Methods in Physics Research, A (2017), https://doi.org/10.1016/j.nima.2017.12.050.
  * [10] L. Gondara, ”Medical Image Denoising Using Convolutional Denoising auto-encoders,” 2016 IEEE 16th International Conference on Data Mining Workshops (ICDMW), Barcelona, 2016, pp. 241-246, doi: 10.1109/ICDMW.2016.0041.
  * [11] Mateo, Carlos & Talavera, Juan. (2017). Short-Time Fourier Transform with the Window Size Fixed in the Frequency Domain. Digital Signal Processing. 77. 10.1016/j.dsp.2017.11.003.
  * [12] Joanes, Derrick N.; Gill, Christine A. (1998), ”Comparing measures of sample skewness and kurtosis”, Journal of the Royal Statistical Society, Series D, 47 (1): 183–189, doi:10.1111/1467-9884.00122, JSTOR 2988433.
  * [13] Marinac-Dabic, Danica & Krulewitch, Cara & Moore, Roscoe. (2002). The Safety of Prenatal Ultrasound Exposure in Human Studies. Epidemiology (Cambridge, Mass.). 13 Suppl 3. S19-22.10.1097/00001648-200205001-00004.
  * [14] Bulas, Dorothy & Egloff, Alexia. (2013). Benefits and risks of MRI in pregnancy. Seminars in perinatology. 37. 301-4. 10.1053/j.semperi.2013.06.005.
  * [15] Delay, Upekha; Dissanayake, Sajan; Nawarathne, Thoshara; Wasalaarachchi, Wishmi; Hetti Arachchi, Anupama; Abeywardhana, Sachitha ; Withanage, Thanushi; Gunarathne, Samitha; Ekanayake, Mervyn Parakrama; Godaliyadda, G. M. R. I. ; Wijayakulasooriya, J. V. ; Ratnayake, Chathura (2019), “Fetal Movement Detection Dataset Recorded Using MPU9250 Tri-Axial Accelerometer”, Mendeley Data, V2, doi: 10.17632/7svcy4cscy.2
  * [16] Sejdić E.; Djurović I.; Jiang J. (2009). ”Time-frequency feature representation using energy concentration: An overview of recent advances”. Digital Signal Processing. 19 (1): 153–183. doi:10.1016/j.dsp.2007.12.004.
  * [17] Shibli Nisar, Omar Usman Khan, Muhammad Tariq, ”An Efficient Adaptive Window Size Selection Method for Improving Spectrogram Visualization”, Computational Intelligence and Neuroscience, vol. 2016, Article ID 6172453, 13 pages, 2016. https://doi.org/10.1155/2016/6172453.
  * [18] Masci J., Meier U., Cireşan D., Schmidhuber J. (2011) Stacked Convolutional Auto-Encoders for Hierarchical Feature Extraction. In: Honkela T., Duch W., Girolami M., Kaski S. (eds) Artificial Neural Networks and Machine Learning – ICANN 2011. ICANN 2011. Lecture Notes in Computer Science, vol 6791. Springer, Berlin, Heidelberg. https://doi.org/10.1007/978-3-642-21735-7-7.
  * [19] Vincent, Pascal; Larochelle, Hugo (2010). ”Stacked Denoising auto-encoders: Learning Useful Representations in a Deep Network with a Local Denoising Criterion”. Journal of Machine Learning Research. 11: 3371–3408.
  * [20] Howard, Roy. (2015). White noise: A time domain basis. 1-4. 10.1109/ICNF.2015.7288581.
  * [21] Çelikoğlu, Ahmet & Tirnakli, Ugur. (2014). Skewness and kurtosis analysis for non-Gaussian distributions. Physica A: Statistical Mechanics and its Applications. 499. 10.1016/j.physa.2018.02.035.


# THE CONNECTION BETWEEN INTERNETWORK MAGNETIC ELEMENTS AND SUPERGRANULAR
FLOWS

D. Orozco Suárez, Y. Katsukawa<EMAIL_ADDRESS>National Astronomical
Observatory of Japan, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan  L. R.
Bellot Rubio  Instituto de Astrofísica de Andalucía (CSIC), Apdo. de Correos
3004, 18080 Granada, Spain

###### Abstract

The advection of internetwork magnetic elements by supergranular convective
flows is investigated using high spatial resolution, high cadence, and high
signal-to-noise ratio Na I D1 magnetograms obtained with the _Hinode_
satellite. The observations show that magnetic elements appear everywhere
across the quiet Sun surface. We calculate the proper motion of these magnetic
elements with the aid of a feature tracking algorithm. The results indicate
that magnetic elements appearing in the interior of supergranules tend to
drift toward the supergranular boundaries with a non-constant velocity. The
azimuthally averaged radial velocities of the magnetic elements and of the
supergranular flow, calculated from a local correlation tracking technique
applied to Dopplergrams, are very similar. This suggests that, in the long
term, surface magnetic elements are advected by supergranular flows, although
on short time scales their very chaotic motions are driven mostly by granular
flows and other processes.

###### Subject headings:

convection, - polarization, - Sun: photosphere, - Sun: surface magnetism

††slugcomment: Submitted to ApJ

##  1\. Introduction

Figure 1.— _Top-left_ : image from the time sequence of high-resolution
magnetograms. A mpeg animation is available in the electronic edition of the
Astrophysical Journal. It shows the rapid flux appearance and disappearance of
IMEs and their drift toward stronger flux concentrations. Colors outline
magnetic features detected with the tracking algorithm. The magnetograms are
saturated at $\pm$ 30 G. _Top-right_ : Contours of constant distance computed
from the supergranular center (I). _Bottom-left_ : Paths of some of the IMEs
detected and tracked using YAFTA_10 in 13 hours time series of magnetograms.
They are rather lineal in the inner part and more random near the borders.
Color code has no meaning. _Bottom-right_ : Horizontal velocity field (stream-
lines) derived from the motion of granules for a 8 hours window. Axis units
are in arcsecond. Arrows in the top-right panel point to the spatial locations
of centers I, II, and III (see text). The locations of the centers are
represented by red crosses in the rest of the panels.

Convection is probably the most efficient mechanism for energy transport in
the solar surface (Spruit et al., 1990). It works in the same way on all
spatial scales but, due to the large stratification of the solar atmosphere,
its properties change dramatically with height (see Rieutord & Rincon 2010 and
Nordlund et al. 2009 for reviews). For instance, the average diameter of
convective cells varies from the 1 Mm of the granules observed at the surface
to the 20-30 Mm of the supergranular cells at tens of Mm below the surface
that have been inferred from observations, theory, and numerical simulations
(e.g. Rast 2003). The average lifetime also varies from 10 minutes to hours or
days. There is an intermediate scale -mesogranulation- whose existence and
origin are still under debate (see e.g., Yelles Chaouche et al. 2011). First
estimates of the horizontal velocity component of the supergranular flow were
obtained by Simon & Leighton (1964) and Leighton et al. (1962). They reported
velocities of about 0.3 km s-1 with a peak of 0.5 km s-1. Recent
investigations have concentrated on the determination of the weaker vertical
component of the supergranular flow, which is of order 10 m s-1 only (see,
e.g., Duvall & Birch 2010).

It has been known for long that there exists a close relationship between
supergranulation and the fields which make up the quiet Sun magnetic network
(Simon & Leighton, 1964). We also know that most magnetic field elements
pervading the quiet Sun surface are very weak (Orozco Suárez et al., 2007a,
b). These fields can only be advected and transported across the surface by
three different flow patterns: granulation, mesogranulation, and
supergranulation. Indeed, since the discovery of magnetic fields in the
interior of supergranular (network) cells by Livingston & Harvey (1971, 1975)
it is known that internetwork magnetic elements (hereafter IMEs) drift toward
strong magnetic concentrations with a mean horizontal speed111Defined as
$\sqrt{v_{x}^{2}+v_{y}^{2}}$, with $v_{x}$ and $v_{y}$ the two orthogonal
components. of 0.5 km s-1. The proper motion of IMEs can be separated in two
main components222Here we ignore the drifts of IMEs due to mesogranular flows
and the motion of the footpoints of emerging loops: they show a ballistic
velocity of 3 km s-1 during the emergence phase, i.e., when they move from the
granule toward an adjacent intergranular lane (Martínez González & Bellot
Rubio, 2009; Manso Sainz et al., 2011).: one of random nature resulting from
the continual buffeting of the fields by granular convection (Manso Sainz et
al., 2011) and a systematic motion that transports the fields toward strong
flux concentrations delineating the supergranular borders. The horizontal
velocity of IMEs shows an rms fluctuation of about 1.6 km s-1 (de Wijn et al.,
2008). The relative contribution of the random and systematic velocities to
this enormous fluctuation is still under debate. The values reported in the
literature for the net velocity component range from 0.2 km s-1 (de Wijn et
al., 2008) to 0.35 km s-1 (Zirin, 1987). The net velocity may depend on the
location of the IMEs within the supergranular cell and even on the polarity of
the IMEs and the nearby network fields. Zhang et al. (1998) investigated the
spatial dependence of the velocity of IMEs and found a constant radial
velocity (measured from the center of supergranular cells outward) of about
0.4 km s-1 using low resolution (1$\farcs$5) and low cadence (7 minutes)
magnetograms. Recently, de Wijn et al. (2008) have argued that the net
velocity of the IMEs may be larger in the proximity of network fields of
opposite polarity, suggesting the existence of some kind of magnetic
“attraction” between IMEs and the network.

Unambiguous observational proofs that IMEs are advected by convective flows
are still missing. In this Letter we use high spatial resolution magnetograms
acquired with the _Hinode_ satellite (Kosugi et al., 2007) to investigate the
relationship between the net velocity of IMEs and the horizontal supergranular
flow.

##  2\. Observations, magnetic features tracking, and local correlation
tracking

The data used here belong to the Hinode Operation Plan 151, entitled “Flux
replacement in the photospheric network and internetwork”. They consist of
shutterless Stokes I and V filtergrams taken on 2010 April 20 at disk
center333The heliocentric angle varied $\sim$ 0.01 units during the observing
time. with the _Hinode_ Narrowband Filter Imager (NFI; Tsuneta et al. 2008) in
the two wings of Na I D1 589.6 nm, $\pm 16$ pm apart from the line center. The
total integration time for each wavelength and polarization state was 6.4 s.
It took about 38 seconds to complete a 93″$\times$82″ map in the shutterless
mode. The final cadence was set to 80 seconds. The diffraction-limited
resolution of the NFI filtergrams at 589.6 nm is 0$\farcs$24 ($\sim$ 175 Km)
and the pixel size is 0$\farcs$16, thus the data are slightly over-sampled
although the Lyot filtergraph worsen the spatial resolution to about
0$\farcs$3 ($\sim$ 215 Km). We used 13 hours of observations, about half of
the mean lifetime of supergranules (16–23 hr; e.g., De Rosa & Toomre 2004).

Figure 2.— Averaged radial (outward) velocities against radial distance from
the center of the supergranule. Top panels represent the radial velocities
calculated from the tracking of IMEs (left) and from FLCT (right) using three
different centers. The red curves represent fits to a simple supergranular
kinematic model for center I (see the text). The bottom left panel shows a
comparison of the radial velocity calculated from the IMEs (solid line) with
the supergranular flow velocity (dotted and dashed lines, 1 and 8 hr averaged,
respectively) using center I. The fits are put together in the bottom right
panel for comparison.

The data were calibrated as follows: first we obtained magnetograms and
Dopplergrams as $M=0.5\times(V_{b}/I_{b}-V_{r}/I_{r})$ and
$D=(I_{r}-I_{b})/(I_{r}+I_{b})$, with $I$ and $V$ representing the intensity
and the circular polarization signal, and $r$ and $b$ the red and blue wings
of the line. We then removed the (residual) effects of the transmission
profile shift resulting from the satellite orbital variation. Finally, we
applied a subsonic filter (see e.g., Hirzberger et al. 1997) to the data to
remove p-mode oscillations with periods of 5 minutes. The magnetogram signal
$M$ was converted into flux density multiplying it by the calibration constant
$C=9160$ derived from the weak field approximation (Gošić, 2012). No
additional corrections were applied to the Dopplergrams because they are used
only for local correlation tracking. The final rms noise of the magnetograms
is 6.7G, allowing the detection of very weak IMEs. Figure 1 (top left panel)
depicts a single magnetogram of 48″$\times$48″. The stronger circular
polarization signals are concentrated near the borders of the image. These
fields correspond to the network that, in general, outlines supergranular
boundaries. Thus, the image roughly encompasses one supergranular cell. In the
interior, there are small polarization patches that show opposite polarities.
These signals correspond to IMEs.

For detecting and tracking IMEs we used the YAFTA_10 (Yet Another Feature
Tracking Algorithm) software444YAFTA_10, written in IDL, can be downloaded
from the author’s Web site at
http://solarmuri.ssl.berkeley.edu/$\sim$welsch/public/software/YAFTA/
(DeForest et al., 2007). We ignored all magnetic elements with flux densities
below 10 G and smaller than 16 pixel2. Colors in the top-left panel of Figure
1 represent the magnetic elements detected by YAFTA_10 in a single frame. The
animation shows the continuous appearance and disappearance of IMEs of both
polarities in the quiet solar surface and how the tracking software detects
and follows them. Interestingly, the visual impression is one of IMEs moving
radially toward the network. Indeed, as can be seen in the bottom left panel
of Figure 1, the trajectories of the magnetic elements located in the interior
of the supergranule are rather linear and radially aligned with respect to the
center of the field of view (FOV). The paths located in the network look more
random and show little linear trends. Using the tracking results, we derived
the horizontal velocity of the magnetic elements. To compute the velocities,
we neglected IMEs that did not last longer than three frames, leaving us with
a total of 17 322 IMEs, of which only 7 844 are within a 20″ circle from the
center of the FOV.

To calculate the horizontal velocity of the flows we applied a Fourier Local
Correlation Tracking method (FLCT; Fisher & Welsch 2008) to consecutive
Dopplergrams. The algorithm cross-correlates two images by comparing small
subfields. The output is the local displacements of one image with respect to
the other. The subfields were defined by a Gaussian window of full-width at
half-maximum 10 px (about 1$\farcs$5), roughly corresponding to granules that
are the best tracers of large-scale velocity flows (Rieutord et al., 2001).
This process is done for the full time series. The two components of the
horizontal flow velocity are then calculated from running averaged FLCT flow-
field maps. Since the horizontal flow pattern might slightly change depending
on the width of the window and the total number of images used for the
averaging, we have calculated the horizontal velocity for two different
running averaged flow-field maps corresponding to about 1 and 8 hours
duration. The bottom right panel of Fig. 1 displays the stream-lines (flow
map) calculated from the 8 hours averaged FLCT maps. Clearly, the flow goes
from the center of the FOV outward.

Finally, we obtain the radial and transverse components of the horizontal
velocities with respect to the supergranular center, for both IMEs and the
supergranular flow. The radial component is defined to be positive outward. We
locate the center of the supergranule at the position of minimum radial
velocity calculated from the magnetic element tracking, at the position of
minimum horizontal velocity of the 8 hours averaged flow map, and the same for
the 1 hour averaged flow map: C${}_{\rm I}\simeq$ (26$\farcs$6 , 22$\farcs$1),
C${}_{\rm II}\simeq$ (24$\farcs$8 , 21$\farcs$1), and C${}_{\rm III}\simeq$
(25$\farcs$3 , 20$\farcs$2), respectively. Next we calculate azimuthally
averaged values of the radial velocity component as a function of distance
from the center of the supergranular cell (see the top right panel of Figure
1).

##  3\. Results

Figure 2 shows the averaged radial velocity versus radial distance from the
center of the supergranule calculated from the tracking of IMEs and from the
FLCT applied to the Dopplergrams. These values are averages over an extremely
large range of individual radial velocities (the typical rms fluctuation in
individual bins is 0.7 km s-1). The average radial component of the IME
velocity shows a clear dependence with distance to the supergranular cell
center (see top left panel). The magnetic elements near the center of the
supergranule are almost at rest. Their radial velocity increases outwards
until it reaches a maximum of 0.3 km s -1 at 9″ from the supergranular center.
Then, it decreases monotonically. The mean radial velocity, calculated by
averaging the mean radial velocities over radial distance, is about 0.15 km s
-1, smaller than the 0.4 km s -1 estimated555Note that Zhang et al. (1998)
directly averaged the radial velocities of individual IMEs to give a mean
value while de Wijn et al. (2008) use an autocorrelation analysis. by Zhang et
al. (1998) and the 0.2 km s -1 of de Wijn et al. (2008). Remarkably, the
average velocity is always positive, meaning that IMEs tend to move
systematically toward the network. On short time scales, IMEs can move in any
direction with much larger speeds, as revealed by the very large spread of the
individual velocities, which can be both positive and negative.

The top right panel of Fig. 2 shows the radial velocity associated with the
supergranular convection as determined from FLCT with an 8 hour average. The
velocity profile shape and its variation are remarkably similar to those
calculated from the tracking of IMEs. The rms fluctuation in individual bins
is 0.3 km s-1 on average. To stress the similarities between the curves we
represent both radial velocity profiles in the bottom left panel of the same
figure. The radial velocity determined from FLCT is slightly steeper although
the exact variation depends on the adopted time window. Indeed, the velocity
profile becomes gentler as the average window increases from 1 to 8 hours
(dotted and dashed lines, respectively). Overall, Fig. 2 suggests that IMEs
are continuously advected by large-scale convective flows, although on short
time scales the velocity is dominated by granular motions or other effects (as
evidenced by the large rms fluctuations).

To make the radial velocity profiles available to the community we have fitted
them using the supergranular convection kinematic model of Simon & Weiss
(1989) (see also Simon et al. 2001). In their model the radial (outward)
velocity is given by $v(r)=V(r/R)e^{-(r/R)^{2}}$ with $V$ representing the
amplitude of the radial velocity, $R$ the radius of the supergranular cell,
and $r$ the distance from the center of the supergranule. This model
satisfactorily represents the velocity profiles with $V=0.71\pm 0.17$ km s-1
and $R=14.2\pm 2.4$ Mm for the IMEs and $V=0.84\pm 0.14$ km s-1 and $R=12.6\pm
2.6$ Mm for the supergranular radial flow field (see the red curves in Fig.
2). To highlight the similarities and differences in the two fits we represent
them in the bottom right panel of Fig. 2. Thus, according to this model, the
radius of the supergranular cell is about 13 Mm.

##  4\. Conclusions

We have shown that there exists a close correspondence between the average
radial velocity of IMEs and the horizontal velocity of supergranular
convective flows. Such a relation was often taken for granted by the solar
community because weak photospheric fields should be advected by surface
flows, although there was no clear observational proof of it. Our results show
that the dynamic properties of quiet Sun IMEs depend on their location within
supergranular cells. On average, IMEs first accelerate outwards and then
decelerate. This behavior may affect the distribution of magnetic flux within
supergranular cells and also suggests that magnetic diffusion may be more
effective in certain areas of the solar surface.

The authors thank the Japan Society for the Promotion of Science (JSPS) for
the financial support through the postdoctoral fellowship program for foreign
researchers. This work has been partly funded by the Spanish MICINN through
projects AYA2011-29833-C06-04 and PCI2006-A7-0624, and by Junta de Andalucía
through project P07-TEP-2687, including a percentage from European FEDER
funds. _Hinode_ is a Japanese mission developed and launched by ISAS/JAXA,
with NAOJ as domestic partner and NASA and STFC (UK) as international
partners. It is operated by these agencies in co-operation with ESA and NSC
(Norway). Use of NASA’s Astrophysical Data System is gratefully acknowledged.

## References

  * de Wijn et al. (2008) de Wijn, A. G., Lites, B. W., Berger, T. E., Frank, Z. A., Tarbell, T. D., & Ishikawa, R. 2008, ApJ, 684, 1469
  * DeForest et al. (2007) DeForest, C. E., Hagenaar, H. J., Lamb, D. A., Parnell, C. E., & Welsch, B. T. 2007, ApJ, 666, 576
  * De Rosa & Toomre (2004) De Rosa, M. L., & Toomre, J. 2004, ApJ, 616, 1242
  * Duvall & Birch (2010) Duvall, T. L., Jr., & Birch, A. C. 2010, ApJ, 725, L47
  * Fisher & Welsch (2008) Fisher, G. H., & Welsch, B. T. 2008, Subsurface and Atmospheric Influences on Solar Activity, 383, 373
  * Gošić (2012) Gošić, M. 2012, Master Thesis, University of Granada (Spain)
  * Hirzberger et al. (1997) Hirzberger, J., Vazquez, M., Bonet, J. A., Hanslmeier, A., & Sobotka, M. 1997, ApJ, 480, 406
  * Kosugi et al. (2007) Kosugi, T., et al. 2007, Sol. Phys., 243, 3
  * Leighton et al. (1962) Leighton, R. B., Noyes, R. W., & Simon, G. W. 1962, ApJ, 135, 474
  * Livingston & Harvey (1971) Livingston, W., & Harvey, J. 1971, Solar Magnetic Fields, 43, 51
  * Livingston & Harvey (1975) Livingston, W. C., & Harvey, J. 1975, BAAS, 7, 346
  * Martínez González & Bellot Rubio (2009) Martínez González, M. J., & Bellot Rubio, L. R. 2009, ApJ, 700, 1391
  * Manso Sainz et al. (2011) Manso Sainz, R., Martínez González, M. J., & Asensio Ramos, A. 2011, A&A, 531, L9
  * Nordlund et al. (2009) Nordlund, Å., Stein, R. F., & Asplund, M. 2009, Living Reviews in Solar Physics, 6, 2
  * Orozco Suárez et al. (2007a) Orozco Suárez, D., et al. 2007a, ApJ, 670, L61
  * Orozco Suárez et al. (2007b) Orozco Suárez, D., et al. 2007b, PASJ, 59, 837
  * Rast (2003) Rast, M. R. 2003, ApJ, 597, 1200
  * Rieutord & Rincon (2010) Rieutord, M., & Rincon, F. 2010, Living Reviews in Solar Physics, 7, 2
  * Rieutord et al. (2001) Rieutord, M., Roudier, T., Ludwig, H.-G., Nordlund, Å., & Stein, R. 2001, A&A, 377, L14
  * Simon & Leighton (1964) Simon, G. W., & Leighton, R. B. 1964, ApJ, 140, 1120
  * Simon et al. (2001) Simon, G. W., Title, A. M., & Weiss, N. O. 2001, ApJ, 561, 427
  * Simon & Weiss (1989) Simon, G. W., & Weiss, N. O. 1989, ApJ, 345, 1060
  * Spruit et al. (1990) Spruit, H. C., Nordlund, A., & Title, A. M. 1990, ARA&A, 28, 263
  * Tsuneta et al. (2008) Tsuneta, S., et al. 2008, Sol. Phys., 249, 167
  * Yelles Chaouche et al. (2011) Yelles Chaouche, L., Moreno-Insertis, F., Martínez Pillet, V., et al. 2011, ApJ, 727, L30
  * Zhang et al. (1998) Zhang, J., Wang, J., Wang, H., & Zirin, H. 1998, A&A, 335, 341
  * Zirin (1987) Zirin, H. 1987, Sol. Phys., 110, 101