# Removing non-stationary, non-harmonic external interference from gravitational wave interferometer data

## I Introduction

Gravitational wave (gw) research started in the early 1960’s thanks to the pioneering work of Weber . Since that time, there has been an ongoing research effort to develop detectors of sufficient sensitivity to allow the detection of these waves from astrophysical sources.

The effect of a gw of amplitude $`h`$ is to produce a strain in space given by $`\mathrm{\Delta }L/L=h/2`$. The magnitude of the problem facing researchers in this area can be appreciate if from the fact that theory predicts that for a reasonable ‘event rate’ one should aim for a strain sensitivity of $`10^{21}`$ to $`10^{22}`$. This means that if we were monitoring the separation of two free test masses of one meter apart, the change in their separation would be $`10^{21}`$ m. Such figures show the size of the experimental challenge facing those developing gravitational wave detectors, and make it clear that the analysis of the data must realize as much of the detector sensitivity as possible.

The different types of detectors can be classified into two major categories: those using laser interferometers with very long arms  and those using resonant solid masses that may be cooled to ultra-low temperatures . The first have the ability to measure the gravitational wave induced strain in a broad frequency band (expected to range from 50 Hz up to perhaps 5 kHz), while the latter measure the gravitational wave Fourier components around the resonant frequency (usually near 1 kHz), with a bandwidth currently of order a few Hz. For resonant mass antennas, the fundamental limitation to their sensitivity comes from the thermal motion of the atoms that can be reduced by cooling them to temperatures of order 50 mK; existing antennas can then achieve a sensitivity $`10^{19}`$ to $`10^{20}`$.

In the early 1970’s, the idea emerged that laser interferometers might have a better chance of detecting gravitational waves. Detailed studies were carried out by Forward and his group  and by Weiss . Since then, several groups have develop prototype interferometric detectors at Glasgow (10 meter Fabry-Perot), Garching (30 m delay line), MIT, Caltech (40 m Fabry-Perot) and Tokyo. Projects to build long arm laser interferometers have also been funded. These are LIGO project  to build two 4 km detectors in USA, VIRGO project  to build a 3 km detector in Italy, GEO-600 project  to build a 600 m detector in Germany and the TAMA-300 project to build a 300 m detector in Japan. The construction of the above detectors has already started. First observations may come as early as 2000. These observations have a potential to see the full range of gravitational wave sources: periodic, burst, quasiperiodic and stochastic.

In the development of data analysis techniques, it is useful to examine data already available from prototypes. Today, prototype interferometers are routinely operating at a displacement noise level of a few times $`10^{19}`$ m$`/\sqrt{\mathrm{Hz}}`$ over a frequency range from 200 Hz to 1000 Hz corresponding to an rms gravitational-wave amplitude noise level of $`h_{rms}2\times 10^{19}`$ .

There are different noise sources that limit the sensitivity of the laser-interferometer detectors. The stochastic noise can be modeled as a sum of six main contributions : photon shot noise, seismic noise, quantum noise (which follows from the indeterminacy of the position of the test masses due to the Heisenberg uncertain principle), the vibration of the suspension wires (‘violin modes’), and thermal noise from the vibration of the test masses and from the low frequency oscillations of the pendulum suspensions. The above noise sources may be considered among the most important but other sources of noise cannot be ignored.

In the measured noise spectrum of the different prototypes, in addition to the stochastic noise, we observe peaks due external interference, where the amplitudes are not stochastic. The most numerous are powerline frequency harmonics. We have shown how to model and remove these very effectively using a technique we call coherent line removal (clr) . Other lines are clearly related to the power supply but appear at non-harmonic frequencies. In this paper, we will describe their characteristics and we will present a procedure to remove them as well. Our goal is to remove as many interference features as possible, so that the interferometer sensitivity is limited by the genuinely stochastic noise.

clr is an algorithm able to remove interference present in the data while preserving the stochastic detector noise. clr works when the interference is present in many harmonics and they remain coherent with each other. In , we applied clr to some interferometric data and the entire series of wide lines corresponding to the electricity supply frequency and its harmonics were completely removed even when the frequency of the supply was not independently known.

In addition to those lines appearing at multiples of the electricity supply frequency (50 or 60 Hz), there are other interference lines whose frequencies change in step with the supply frequency, but not at the harmonic frequencies. From a data-analysis point of view, we try to develop a technique able to remove this interference while producing a minimum disturbance to the underlying noise background. In order to remove these lines, other methods are also available  but these methods will remove the noise and any underlying real signals as well.

We are not sure if this type of interference will be absent in large-scale interferometers, since the physical process that generates all these lines is so far unknown. But it is present in prototype data and therefore it is important to be prepared to remove it from full data.

This paper is based on a study of Glasgow interferometer data taken in March 1996. The method we propose makes use of a reference wave-form signal corresponding to the fundamental harmonic of the electrical interference. We can obtain it directly from the supply voltage or we can construct it using the clr algorithm from the true harmonics present in the data. The method we propose is an adaptive procedure that is tuned in such a way that the electrical interference can be removed and ‘single-line’ signals masked by them can be recovered to at least the $`75\%`$ level.

The rest of the paper is organized as follows: In section II, we describe the electrical interference present in the data. In section III, we present different models of the interference. In section IV, we summarize the principle of the coherent line removal algorithm and we explain how to construct a reference wave-form of the incoming electricity signal from the data. In section V, we present an algorithm to remove the electrical interference, and not only the harmonics of the reference wave-form. The algorithm is recursively applied for small stretches of data. Thus, it allows the parameters to change and to adapt themselves in order to be able to remove the interference with a minimum disturbance of the noise background. Finally, in section VI, we discuss the results obtained.

## II The electrical interference in the prototype data

In this paper, we will focus our attention to the data produced by the Glasgow laser interferometer in March 1996. The data set consists of 19857408 points, sampled at 4 kHz and quantized with a 12 bit analog-to-digital converter with a dynamical range from -10 to 10 Volts. The data are divided into 4848 blocks of 4096 points each. The first 18 minutes of data were rendered useless due to a failure of the autolocking. Thus, in our analysis we ignore the first 1153 blocks. (For preliminary studies of these data see .)

In the study of the prototype data, we observe in the power spectrum many instrumental lines. Some of them are rather broad and appear at multiples of 50 Hz. But there are other ones appearing at different frequencies. In the data, the lines at 1 kHz have a width of 5 Hz. Therefore, we can ignore these sections of the power spectrum or we can try to understand the interference and, if possible, remove it in order to be able to detect any possible gravitational wave signal previously masked by it.

We have already shown how to remove lines at integer multiples of the supply frequency . In long-term Fourier transforms, these lines are broad, and the structure of different lines is similar apart from an overall scaling proportional to the frequency. In smaller length Fourier transforms, the lines are narrow, with central frequencies that change with time, again in proportion to one another. It thus appears that all these lines are harmonics of a single source (e.g., the electricity supply) and that their broad shape is due to the wandering of the incoming electricity frequency. These lines have been observed in different interferometer prototypes .

But this is not the end of the story. Further analysis of the prototype data reveals the presence of many other features that are related to the incoming electricity frequency. These other lines are not as powerful as the harmonics. In many cases, they are just slightly above the stochastic noise level.

The easiest way to detect their presence is by studying in detail the spectrogram (i.e., the magnitude of the time-dependent Fourier transform versus time). This is a time-dependent frequency analysis in which the whole data set is split into small segments and for each of them a discrete Fourier transform is computed. By examining the spectrogram, we identify a large number of lines that present similar time-frequency evolution to the harmonics of 50 Hz. Therefore, all these line are related in some way with the incoming electricity supply. An example is shown in figure 1.

These lines are spread over the whole spectrum. A big population of them lies below 140 Hz. There are also some isolated ones around 222, 238, 444, 575, 887, 1105 and 1275 Hz and, after 1750 Hz, there is again a large population of them. Our guess is that this kind of line might be present through the whole spectrum, but many are buried in the stochastic noise.

An interesting feature is that some of these lines appear to fall in to their own harmonic families. For example, we have found the families: (222, 443, 886, 1772 Hz), (48, 96 Hz), (66, 132 Hz), (72.5, 145 Hz), (116, 232 Hz). In all cases, the width of the lines (which is the interval of the frequency wandering) seems to scale as the frequency. The amplitude of the lines is time-dependent as well.

Although the nature of the process that generates these lines remains unknown, many of them may be intrinsic to the incoming electricity signal. We have analyzed some a-posteriori electricity supply voltage (from Glasgow University) and we have observed several features, and not just harmonics of 50 Hz. External switches or the running of other electrical devices (as computers) can generate some of these lines. This could explain their non-linear and non-stationary character. Another possible sources are ground motions due to mechanical motors running at different frequencies than 50 Hz (e.g., trains), but whose frequency wanders in a way related somehow to the power supply.

Remembering the extremely small amplitude of the disturbances these lines represent, it follows that it may be difficult for experimenters to exclude these lines from the data of detectors now under construction. It is therefore important to understand how to remove them if they are there.

## III Modeling the lines

A first model for these non-harmonic lines is that they could be beats between stationary frequencies and the supply harmonics. We have examined this in detail and we believe it is unlikely. In particular, the width of these lines seems to be proportional to their frequency, which could not be the case for a beat, and also we would expect the lines to appear in couples, that it is not the case.

Using this model,

$$h(t)=\alpha M(t)^n\mathrm{exp}(i2\pi ft),$$
(1)   
where $`\alpha `$ is a complex amplitude, $`M(t)^n`$ is a supply harmonic, and $`f`$ is the beat frequency, we have analyzed several lines. For example, the line at 99.7 Hz, we find that the best match corresponds to the second harmonic $`n=2`$ with a beat frequency of $`f=0.7391`$ Hz. Using these values, we try to remove the interference by applying a least square method, but as we show in figure 2 the interference is not cancelled. The same occurs with other lines.

After rejecting the possibility of the beats, the next simplest model is to assume a non-integer harmonic of the supply:

$$h(t)=\alpha (t)M(t)^q,$$
(2)   
where $`\alpha (t)`$ is a slowly varying complex amplitude, $`M(t)`$ is a reference wave form corresponding to the fundamental harmonic of the electrical interference, but now $`q`$ is a real number, not only an integer as in the case of the harmonics. The reference signal, $`M(t)`$, can be obtained directly from the supply voltage, or it can be obtained directly using clr as we describe below. The justification for this model is simply the success we have in removing the interference.

The model is not perfect, and our line-removal is sometimes incomplete. It appears that the index $`q`$ may also be a slow function of time. But we will not pursue this refinement here.

## IV Coherent Line Removal

For any of the previous models, we need to know the reference wave-form, $`M(t)`$, corresponding to the fundamental harmonic of the interference. In the case of electrical interference, this wave-form can be obtained from the data applying the clr algorithm. In this section, we summarize the principle of clr. For further details we refer the reader to .

clr works when the interference is present in many harmonics, and assumes that the interference has the form

$$y(t)=\underset{n}{}a_nm(t)^n+\left(a_nm(t)^n\right)^{},$$
(3)   
where $`a_n`$ are complex amplitudes, $`m(t)`$ is a nearly monochromatic function near a frequency $`f_0`$ and $``$ is the complex conjugate. The idea is to use the information in the harmonics of the interference to construct the reference function $`M(t)`$ that is as close a replica as possible to $`m(t)`$.

Assuming additive noise, the data produced by the system is just

$$x(t)=y(t)+n(t),$$
(4)   
where $`y(t)`$ is the interference given by Eq. (3) and the noise $`n(t)`$ in the detector is a zero-mean stationary stochastic process. The algorithm consists in defining a set of functions $`\stackrel{~}{z}_k(\nu )`$ in the frequency domain as

$$\stackrel{~}{z}_k(\nu )\{\begin{array}{cc}\stackrel{~}{x}(\nu )& \nu _{ik}<\nu <\nu _{fk}\\ 0& \text{elsewhere},\end{array}$$
(5)   
where $`\stackrel{~}{}`$ denotes the Fourier transform, $`(\nu _{ik},\nu _{fk})`$ correspond to the upper and lower frequency limits of the harmonics of the interference and $`k`$ denotes the harmonic considered. These functions are equivalent to

$$\stackrel{~}{z}_k(\nu )=a_k\stackrel{~}{m^k}(\nu )+\stackrel{~}{n}_k(\nu ),$$
(6)   
where $`\stackrel{~}{n}_k(\nu )`$ is the noise in the frequency band of the harmonic considered. Their inverse Fourier transforms yield

$$z_k(t)=a_km(t)^k+n_k(t).$$
(7)   
Since $`m(t)`$ is supposed to be a narrow-band function near a frequency $`f_0`$, each $`z_k(t)`$ is a narrow-band function near $`kf_0`$. Then, we define

$$B_k(t)\left[z_k(t)\right]^{1/k},$$
(8)   
that can be rewritten as

$$B_k(t)=(a_k)^{1/k}m(t)\beta _k(t),$$
(9)   
where

$$\beta _k(t)=\left[1+\frac{n_k(t)}{a_km(t)^k}\right]^{1/k}.$$
(10)   
All these functions, $`\{B_k(t)\}`$, are almost monochromatic around the fundamental frequency, $`f_0`$, but they differ basically by a certain complex amplitude. These factors, $`\mathrm{\Gamma }_k`$, can easily be calculated, and we can construct a set of functions $`\{b_k(t)\}`$

$$b_k(t)=\mathrm{\Gamma }_kB_k(t),$$
(11)   
such that, they all have the same mean value. Then, $`M(t)`$ can be constructed as a function of all $`\{b_k(t)\}`$ in such a way that it has the same mean and minimum variance. If $`M(t)`$ is linear with $`\{b_k(t)\}`$, then statistically the best choice for $`M(t)`$ is

$$M(t)=\left(\underset{k}{}\frac{b_k(t)}{\mathrm{Var}[\beta _k(t)]}\right)/\left(\underset{k}{}\frac{1}{\mathrm{Var}[\beta _k(t)]}\right),$$
(12)   
where

$$\mathrm{Var}[\beta _k(t)]=\frac{n_k(t)n_k(t)^{}}{k^2|a_km(t)^k|^2}+\text{corrections}.$$
(13)   
In practice, we approximate

$$|a_km(t)^k|^2|z_k(t)|^2,$$
(14)   
and we assume stationary noise. Therefore,

$$n_k(t)n_k(t)^{}=_{\nu _{ik}}^{\nu _{fk}}S(\nu )𝑑\nu ,$$
(15)   
where $`S(\nu )`$ is the power spectral density of the noise.

The amplitude of the different harmonics of the interference can be obtained then applying a least square method.

## V Removal of non-harmonic external interference

clr can remove the integer harmonics of a reference signal but it does not remove other interference present in the data. We can use clr to construct a reference waveform $`M(t)`$, but we have to design another technique able to get read off the electrical interference present at non-harmonics of the 50 Hz line frequency.

The models of signals we propose are buried in noise. Therefore, we face the problem of detecting signals and estimating their parameters.

### A Maximum likelihood detection

A standard method is the maximum likelihood detection which consists of maximizing the likelihood function $`\mathrm{\Lambda }`$ with respect to the parameters of the signal. If the maximum of $`\mathrm{\Lambda }`$ exceeds a certain threshold, we say that the signal is present. (See  for signal analysis theory in the context of gravitational wave broadband detectors.)

We assume that the noise $`n(t)`$ in the detector is an additive, zero-mean, Gaussian and stationary random process. Then the data $`x(t)`$ (if the expected signal model $`h(t)`$ is present) can be written as

$$x(t)=n(t)+h(t).$$
(16)   
The logarithm of the likelihood function has the form

$$\mathrm{ln}\mathrm{\Lambda }=\underset{k=0}{\overset{N1}{}}\frac{\stackrel{~}{x}_k\stackrel{~}{h}_k^{}}{S_k}\frac{1}{2}\underset{k=0}{\overset{N1}{}}\frac{|\stackrel{~}{h}_k|^2}{S_k},$$
(17)   
where $`S_k`$ is the power spectral density of the noise, $`k`$ is the frequency index running from $`0`$ to $`N1`$ and $`N`$ is the number of sampled points. The likelihood ratio $`\mathrm{\Lambda }`$ depends on the particular set of data $`x(t)`$ only through the sum

$$G=\underset{k}{}\frac{\stackrel{~}{x}_k\stackrel{~}{h}_k^{}}{S_k}.$$
(18)   
This sum is called the detection statistic for the signal $`h`$. Its variance is

$$d^2=\underset{k}{}\frac{|\stackrel{~}{h}_k|^2}{S_k}.$$
(19)   
If there is no signal present, the mean of $`G`$ is zero, but if the signal $`h`$ is present the mean of $`G`$ will be equal to its variance.

The thresholds on the likelihood function $`\mathrm{\Lambda }`$ must be set with regard to the false alarm probability. For a detection statistics $`G`$ and a variance $`d^2`$, the false alarm probability is

$$P_F=\frac{1}{2}\mathrm{erfc}\left(\frac{G}{\sqrt{2}d}\right).$$
(20)   
This is equivalent to study the output signal-to-noise ratio (snr), i.e., the value of the detection statistics divided by the standard deviation,

$$\mathrm{snr}\frac{G}{d}=\underset{k}{}\frac{\stackrel{~}{x}_k\stackrel{~}{h}_k^{}}{S_k}/\sqrt{\underset{k}{}\frac{|\stackrel{~}{h}_k|^2}{S_k}}.$$
(21)   
Notice that snr (and therefore the false alarm probability) is independent of the value of the amplitude of the model signal $`h(t)`$ used for this pattern-matching procedure. Of course, the snr is proportional to the amplitude of whatever multiple of $`h(t)`$ is contained in $`x(t)`$ . This is thus a linear detector.

### B The parameter space

Assuming the form for the electrical interference, $`h(t)=M(t)^q`$, we have to construct as many filters as different values of $`q`$ need to be considered and, for each of them, calculate their snr.

In order to set the parameter space, we can consider $`M(t)`$ as a monochromatic signal at a frequency $`f_0`$. Hence, $`M(t)^q`$ will be a monochromatic signal at $`qf_0`$. Thus, the maximum value of $`q`$ to be considered corresponds to

$$q_f=\frac{f_{Nyquist}}{f_0},$$
(22)   
where

$$f_{Nyquist}=\frac{f_s}{2},$$
(23)   
and $`f_s`$ is the sampling frequency.

The frequency resolution is $`\mathrm{\Delta }\nu =1/T`$. Therefore, we can resolve two signals if the separation in $`q`$ is of the order

$$\mathrm{\Delta }q=\frac{\mathrm{\Delta }\nu }{f_0}.$$
(24)   
This is the maximum separation in $`q`$ we can allow. Note that the size of the parameter space

$$N_q=\frac{q_f}{\mathrm{\Delta }q}=\frac{f_s}{2}T,$$
(25)   
increases with the observation time $`T`$.

For 128 blocks of the Glasgow data (this corresponds to the usual stretch of data we work at once), we obtain

$$\mathrm{\Delta }q0.00015.$$
(26)   
We can calculate the minimum number of filters and the number of floating points operations needed to calculate the snr for all of them.

### C An approximately matched filter

The calculation of the exact matched filter for all the parameter space of $`q`$ is computationally expensive. Given the reference function $`M(t)`$, for each value of $`q`$, we need to calculate $`h(t)=M(t)^q`$, perform its Fourier Transform, and then, compute the output filter via the formula of the snr given by Eq. (21).

Since we assume $`M(t)`$ is a nearly monochromatic function, all the functions $`h(t)=M(t)^q`$ are also going to be nearly monochromatic but at different frequencies. This means that the values of their Fourier transforms are just relevant in small frequency bands. Hence, the first approximation we can perform is to reduce the index summation in Eq. (21) to a small interval. This is equivalent to consider that the Fourier transform of the template is zero outside this interval.

But this is not enough. What is really computationally expensive is the construction of the templates. Since the reference signal $`M(t)`$ is almost monochromatic and changes frequency smoothly in time, we can approximate the values of $`\stackrel{~}{M^q}(\nu )`$ in a certain interval, $`\nu _{iq}<\nu <\nu _{fq}`$, by those of $`\stackrel{~}{H^q}(\nu )`$ defined by

$$\stackrel{~}{H^q}(\nu )\stackrel{~}{M^n}\left(\nu \frac{q}{n}\right),$$
(27)   
where $`n`$ can be chosen to be an integer, i.e., we can build $`\stackrel{~}{H^q}`$ as a similar function to $`\stackrel{~}{M^n}`$, where $`M^n(t)`$ corresponds to an harmonic of the reference signal.

If we construct $`\stackrel{~}{H^q}(\nu )`$ via Eq. (27) using the nearest harmonic, we will expect $`\stackrel{~}{H^q}(\nu )`$ to be close to $`\stackrel{~}{M^q}(\nu )`$. Therefore, we just need to calculate the Fourier transforms of all the harmonics (which is a small number in comparison to all the possible values of $`q`$), calculate $`\stackrel{~}{H^q}(\nu )`$ from the nearest harmonic and substitute their values in equation (21).

The power spectral density, $`S_k`$, can be estimated from the data using Welch Method  averaging over shorter periodograms.

With all these simplifications, we apply the approximate matched filter to 128 blocks (approximately two minutes) of prototype data. By setting a threshold of snr=4.5, we have found several values of $`q`$ for which the electrical interference might be present. See table I.

If we assume that the interference follows the model

$$h(t)=\alpha M(t)^q,$$
(28)   
with $`\alpha `$ being constant, for each $`q`$ we can find the value of the complex amplitude $`\alpha `$ by applying a least square method in the time domain, i.e., $`\alpha `$ is the value of $`\beta `$ that minimizes the quantity $`|x(t)\beta M(t)^q|^2`$.

With this method and using the same 128 blocks of data, we try to remove the interference. The result is that the interference is attenuated but not removed. Therefore, we think that this model is too simple and that to remove real interference we need a more complicated model.

### D Amplitude modulation

The extra complication is to allow the amplitude to vary slowly with time. Therefore, we assume that the interference takes the form

$$h(t)=\alpha (t)M(t)^q,$$
(29)   
where $`q`$ are real constants, $`\alpha (t)`$ is a slow changing function of time, and the values of $`q`$ do not differ much from those calculated with the matched filter as described in the previous subsection.

We have to find a procedure to determine $`\alpha (t)`$, the amplitude modulation. Since the signal is buried in noise, we cannot find its exact value. What we do is to find an approximate function $`\beta (t)`$ to $`\alpha (t)`$ such that allows us to remove the interference and, at the same time, keep the intrinsic level of noise present in the interferometer.

The method we use consists of splitting the data into small pieces (i.e., small number of blocks), and for each piece calculating a value $`\beta `$ as if it was a constant. Then, we construct $`\beta (t)`$ as the succession of those values obtained.

In practice, we separate the data into sets of $`2n_b+1`$ blocks with overlaps of $`n_b`$ blocks. For each of them we calculate the $`\beta `$ value and we associate it to the block $`n_b+1`$. In this way, we construct $`\beta (t)`$ as a set of discrete values that change smoothly in time. The number of blocks, $`n_b`$, used must be tuned according to the data. If the value of $`n_b`$ is too big, it does not allow enough amplitude modulation and the interference is not cancelled. By contrast, if $`n_b`$ is too small, the function $`M(t)^q`$ for those small number of blocks will be almost monochromatic, in the sense that it will affect to a few frequency bins. Then, this method will behave as an adaptive multitaper method  and hence, it will remove whatever is in those frequency bins (any signal plus noise).

For the prototype data, we have found that these values, $`n_b`$, need to depend on $`q`$. Thus, we build the function $`n_b(q)`$ satisfying the following requirements:

* it must be able to remove the interference,
* it should leave a noise level comparable to the intrinsic noise background of the interferometer,
* an artificial ‘single-line’ signal with an amplitude of the order of the electrical interference should not be attenuated more than 25 %.

To satisfy all these requirements, we propose the following function

$$n_b(q)=\{\begin{array}{cc}7,\hfill & q<2\\ 5,\hfill & 2q<4\\ 4,\hfill & 4q<6\\ 3,\hfill & 6q<8\\ 2,\hfill & 8q<16\\ 1,\hfill & q16.\end{array}$$
(30)   
Using this function, $`n_b`$, and the values of $`q`$ displayed in table I, we have succeeded in removing the interference present in 128 blocks. The result are shown in figures 3 and 4.

### E Frequency drift

We proceed now to remove the interference from the whole data stream. We split the data into fragments of 128 blocks and we apply the previous method using the amplitude modulation.

As a first attempt, we assume that the values of $`q`$ remain constant during all time. As a result, some of the lines are removed but some others are not after a certain time.

This is not surprising. Visual inspection of the spectra shows that these lines drift,

$$q(t)=q_0+\delta q(t),$$
(31)   
where $`\delta q(t)`$ is small in comparison with $`q_0`$.

Therefore, for the longer data set, we allow small changes in $`q`$. We apply to each fragment of 128 blocks the approximate matched filter, as described before, but using a much reduced parameter space, i.e., just allowing a maximum variation of $`10\mathrm{\Delta }q`$ about the $`q`$ values of the previous fragment of 128 blocks.

In order to construct the templates, we make use again of Eq. (27), but we use the old values of $`q`$ (those obtained in the previous data fragment). Then we choose the new values of $`q`$ by maximizing the snr.

Using this procedure, we have removed all the lines corresponding to those initial values of $`q`$ listed in table I. See figure 5.

From the evolution of $`q`$, we observe that some values of $`q`$ remain almost constant, but many others change in time, obtaining a maximum variation of 0.005 for 3695 blocks of data.

Assuming the incoming signal is monochromatic, we can estimate the best timescale to perform the matched filter, i.e., the length of the stretch of data for which the frequency resolution is greater than the maximum expected frequency variation due to the drift of $`q`$. This yields 120 blocks, while we were using 128.

## VI Discussion

We have described an algorithm able to remove any kind of interference related with the incoming electricity main supply. The study is based on the data produced by the Glasgow interferometer prototype. In the data, we have observed many interference lines that are highly non-linear and non-stationary.

The form of the interference can be modeled by $`h(t)=\alpha (t)M(t)^{q(t)}`$, where $`M(t)`$ corresponds to the fundamental harmonic of the incoming electrical signal, $`\alpha (t)`$ is a slow varying function of time, and $`q(t)`$ is almost constant, but it can drift in time, i.e., $`q(t)=q_0+\delta q(t)`$, where $`\delta q(t)`$ is small in comparison to $`q_0`$.

In order to detect those signals and estimate their parameters, we use the method of maximum likelihood detection to determine the values of $`q`$ (at a certain instance). Then, we apply an adaptive procedure to determine the amplitude modulation and we repeat it recursively for the whole data stream.

The result is that all those lines which were initially detected (i.e., those with enough snr) have been tracked and completely removed. This method is able to recover monochromatic signals that are buried by the interference. The signal distortion is less than $`25\%`$. Thus, this procedure can assist in the search for continuous waves and also clean the statistics of the noise in the time domain. As we pointed out in an earlier paper , the removal of lines like these can reduce the level of non-Gaussian noise. Therefore, line removal is important since it can raise the sensitivity and duty cycle of the detectors to short bursts of gravitational waves, as well.

###### Acknowledgements.

We would like to thank J. Hough and the gravitational waves group at Glasgow University for providing their gravitational wave interferometer data for analysis. This work was partially supported by the European Union, TMR Contract No. ERBFMBICT972771.