Keywords: Interferometry, radio, LOFAR, algorithms, extragalactic

Mythopoeia, that "He sees no stars who does not see them first, of living silver made that sudden burst, to flame like flowers beneath the ancient song" 1 . In his defense of myth-making, he formulates the argument that the attribution of meaning is an act of creation -that "trees are not 'trees' until so named and seen"and that this capacity for creation defines the human creature. The scientific endeavour, in this context, can be understood as a social expression of a fundamental feature of humanity, and from this endeavour flows much understanding. This thesis, one thread among many, focuses on the study of astronomical objects as seen by the radio waves they emit.

Il faut imaginer Sisyphe heureux.

-Albert Camus

When you make the finding yourself -even if you're the last person on Earth to see the light -you'll never forget it.

-Carl Sagan

What are radio waves? Electromagnetic waves were theorised by James Clerk Maxwell (Maxwell 1864) in his great theoretical contribution to modern physics, their speed matching the speed of light as measured by Ole Christensen Rømer and, later, James Bradley. It was not until Heinrich Rudolf Hertz's 1887 experiment that these waves were measured in a laboratory, leading to the dawn of radio communications -and, later, radio astronomy. The link between radio waves and light was one of association: light is known to behave as a wave (Young double-slit experiment), with the same propagation speed as electromagnetic radiation. Light "proper" is also known to exist beyond the optical regime: Herschel's experiment shows that when diffracted through a prism, sunlight warms even those parts of a desk which are not observed to be lit (first evidence of infrared light). The link between optical light and unseen electromagnetic radiation is then an easy step to make, and one confirmed through countless technological applications (e.g. optical fiber to name but one). And as soon as this link is established, a question immediately comes to the mind of the astronomer: what does the sky, our Universe, look like to the radio "eye"? Radio astronomy has a short but storied history: from Karl Jansky's serendipitous observation of the centre of the Milky Way, which outshines our Sun in the radio regime, in 1933, to Grote Reber's hand-built back-yard radio antenna in 1937, which successfully detected radio emission from the Milky Way itself 2 , to such monumental projects as the Square Kilometer Array and its multiple pathfinders, it has led to countless discoveries and the opening of a truly new window on the Universe. The work presented in this thesis is a contribution to this discipline -the culmination of three years of study, which is a rather short time to get a firm grasp of radio interferometry both in theory and in practice. The need for robust, automated methods -which are improving daily, thanks CONTENTS to the tireless labour of the scientists in the field -is becoming ever stronger as the SKA approaches, looming large on the horizon; but even today, in the precursor era of LOFAR, MeerKAT and other pathfinders, it is keenly felt. When I started my doctorate, the sheer scale of the task at hand felt overwhelming -to actually be able to contribute to its resolution seemed daunting indeed! Thankfully, as the saying goes, no society sets for itself material goals which it cannot achieve. This thesis took place at an exciting time for radio interferometry: at the start of my doctorate, the LOFAR international stations were -to my knowledge -only beginning to be used, and even then, only tentatively; MeerKAT had not yet shown its first light; the techniques used throughout my work were still being developed. At the time of writing, great strides have been made. One of the greatest technical challenges of LOFAR -imaging using the international stations -is starting to become reality. This technical challenge is the key problem that this thesis set out to address. While we only achieved partial success so far, it is a testament to the difficulty of the task that it is not yet truly resolved.

One of the major results of this thesis is a model of a bright resolved source near a famous extragalactic field: properly modeling this source not only allows the use of international LOFAR stations, but also grants deeper access to the extragalactic field itself, which is otherwise polluted by the 3C source's sidelobes. This result was only achieved thanks to the other major result of this thesis: the development of a theoretical framework with which to better understand the effect of calibration errors on images made from interferometric data, and an algorithm to strongly mitigate them.

The structure of this manuscript is as follows: we begin with an introduction to radio interferometry, LOFAR, and the emission mechanisms which dominate for our field of interest. These introductions are primarily intended to give a brief overview of the technical aspects of the data reduced in this thesis. We follow with an overview of the Measurement Equation formalism, which underpins our theoretical work. This is the keystone of this thesis.

We then show the theoretical work that was developed as part of the research work done during the doctorate -which was published in Astronomy & Astrophysics. Its practical application -a quality-based weighting scheme -is used throughout our data reduction. This data reduction is the next topic of this thesis: we contextualise the scientific interest of the data we reduce, and explain both the methods and the results we achieve.

Astrophysical Interest

Radio astronomy is a young but extremely dynamic field with a large range of science cases. From the study of the epoch of reionisation, which involves looking for polarised, redshifted 21cm line emission to study the emergence of the first luminous sources1 , to that of variable radio sources such as pulsars and (theoretically) exoplanets, going through the study of ultra-high-energy cosmic rays, solar science, space weather science, the study of cosmic magnetism, and the creation of deep extragalactic surveys, radioastronomy opens up an extremely fertile window into our Universe.

In Fig. 1.1, for example, much structure is immediately visible: supernova remnants in the galactic plane (the "bubbles"), gas following galactic magnetic field lines (the filamentary structure), and emission from diffuse gas structure around the core of the Milky Way. Other sources are also visible, however: some unresolved sources, and very 6 CHAPTER 1. THESIS OVERVIEW: ASTROPHYSICS & TECHNIQUES many AGN sources. This is one particularly striking example of what radio telescopes allow us to see of astrophysical sources, but is only one among many. The work presented here falls under the auspices of the deep extragalactic surveys science project. The underlying scientific motivation behind this Key Science Project or KSP is the study of, among other things, radio-active galactic nuclei, or radio AGNs. A radio AGN is a type of radio galaxy whose emission is dominated by its nucleus. This nucleus emission is non-thermal, and assumed to be emitted from both collimated jets springing from the accretion of matter onto a galaxy's supermassive black hole and from a hot accretion disk surrounding said black hole. AGNs are then commonly split into different types, based on spectral line behaviour (as in Seyfert 1943) or morphological characteristics (as in Fanaroff & Riley 1974).

So-called unified AGN models, which aim to explain the underlying physical object which is assumed to lie behind all observed examples of AGN, were developed to explain the diversity in observed AGN types (e.g. Antonucci & Miller 1985;Tasse et al. 2008;Netzer 2015;Morabito et al. 2017). These models would include a central engine, usually a supermassive black hole with hot accretion (innermost accretion regions estimated to be around 10 7 K by Pringle & Rees 1972), emitting thermally; an optically thick disk torus, typically of size 10 to 100pc, blocking light emitted from the central engine along the equatorial plane, and the cavity inside which would be the region from which broad emission lines would emit; and finally, ionisation cones in the direction of the poles, where the emission from the central engine would collimate into large-scale jets (which would also, at the scale of a few tens of parsecs, create the region of narrow spectral line emission).

The interaction between AGN jets and the surrounding intergalactic medium can result in structure such as that shown in Fig. 1.2 below: the bright nucleus of a radio galaxy in the centre of two extremely radio-bright jets, with lots of diffuse emission from their lobes. These lobes show where the accretion-accelerated particles begin to slow down -and therefore emit -due to contact with the intergalactic medium. Hence, their sometimes peculiar shape: this medium can be complex and inhomogeneous (if the radio galaxies lie in a cluster, for example) which can lead to the sort of jet structure shown in this image. Of course, not all galaxies have active nuclei, and those which do also have other forms of radio emission: for example, free-free radiation from the interstellar medium, ionised by light from newborn stars, will often contribute to radio-galactic flux. Some galaxies with active nuclei do not have strong jet emission: these are called radio-quiet AGNs. We can thus identify two components to a radio AGN: the host galaxy and its jets. Peak star formation is thought to have occurred around z ∼ 2, and so the study of AGN activity around this period is of strong scientific interest. Could the jets from active AGNs have contributed to the decline in star formation? Are there large populations of older AGNs, and if so, how old? These are the astrophysical questions which motivate the study of AGN sources (Best 2007;Xu et al. 2018;[START_REF] Ueda | [END_REF], which is the driving scientific concern which motivates the present doctorate.

The Low-ν Sky: Emission Mechanisms

This section aims to bring a brief introduction to the emission mechanisms which dominate at low frequencies, and thus determine the physics accessible to extragalactic astronomers working in this band. Specifically, it will briefly describe free-free radiation (which is emitted from ionised plasma outflows, and thus a tracer of various physical objects e.g. young stellar objects -cf. Coughlan et al. (2017) and references thereinand starburst regions -cf. Varenius et al. (2015)) and synchrotron radiation, which is a tracer of more violent and energetic processes (and thus associated with supernova remnants, AGN, or radio halos -cf. Cassano et al. (2010)).

Of course, each individual galaxy will have differing contributions from different mechanisms -in practice, when studying galactic populations, the overall flux at a given frequency will simply be summed up for each galaxy and become one point in a Spectral Energy Distribution, or SED. However, when studying individual objects, it can be critical to understand which emission mechanism dominates at what frequencies. For example, the radio and far-infrared spectrum for nearby M82 are shown in Fig. 1.3 2 : https://www.cv.nrao.edu/course/astr534/FreeFreeEmission.html. The flat curve corresponds to free-free emission, while synchrotron radiation (negative slope) and thermal dust emission (postive slope) dominate at low and high frequencies respectively. 2 Figure and work taken from the NRAO website. For further information, see here

Synchrotron Radiation

Synchrotron radiation is the dominant mode of emission for extragalactic sources at low wavelengths. As such, it will be covered in greater detail than the other mechanisms above, which are important in other bands (and therefore relevant for multi-spectral analysis) but much less relevant to LOFAR observations. This section draws heavily from Garret Cotter's high-energy astrophysics lectures3 , as well as Alan Loh's doctoral thesis (Loh 2016) and the NRAO online course4 . Synchrotron radiation (or "magnetobremsstrahlung") occurs when the trajectory of a relativistic charged particle is guided by a magnetic field. As the German name suggests, it is the magnetic equivalent of free-free radiation. It is the tracer of very energetic processes, such as the interaction between AGN jets and the intergalactic medium. In the low-energy regime, it is also known as cyclotron radiation, after the device in which it was first measured in a laboratory.

The synchrotron power spectrum of a population of electrons with exactly the same energy is a function of its Lorentz factor and its emission angle5 [START_REF] Longair | High energy astrophysics[END_REF]). The typical radiation spectrum for a single electron is shown in Fig. 1.4 as a function of the critical frequency, which is itself a function of the emission angle θ, the Lorentz factor γ and the magnetic field strength B: ν c = 3 2 γ 2 eB 2πmec sin(θ). e and m e are the electron charge and mass, respectively.

Of course, in reality, we know that electron populations in astrophysical sources do not all share exactly the same energy when emitting. The effect of the electron energy distribution must therefore be taken into account. Because the population in question is energised through non-thermal processes, its energy distribution is not Maxwellian, and is generally written as a power law: S(E) ∝ E α → S ν ∝ ν -α , where α is the object's spectral index. Note that the index of the energy distribution (which relates the number of electrons at a given energy) is different from the synchrotron index, with the two related by the following equation: α = -(γ -1) 2. The question now becomes: what α is appropriate for what kind of object? Does every source have the same α at all frequencies? Do all synchrotron sources have the same spectral index at a given frequency, and if not, what does this tell us?

To answer these questions, we must recall that so far, we have assumed that every photon emitted through synchrotron radiation reaches us. This is not necessarily the case in practice: it is possible for a given photon to be scattered (i.e. absorbed & re-emitted) by surrounding high-energy electrons. The likelihood of this occurring at a given frequency is a function of the absorption cross-section, and is a complex quantity to calculate (see Longair (1994) pp258-260). Suffice to say that at longer wavelengths (lower frequencies) absorption is more effective, and at shorter wavelengths (higher frequencies)

As we can see, this spectrum falls off sharply with decreasing frequency. As such, while it is expected to dominate over thermal emission in the absence of synchrotron radiation, it is synchrotron emission that is expected to dominate overall -if present -at the low frequencies of our observations, which are done with LOFAR. Free-free emission acts as a tracer for star formation and other "gentler" physical processes detected at low radio frequencies.

LOFAR: The LOw-Frequency Array

Our observations have all been made with the LOw Frequency Array LOFAR (van Haarlem et al. 2013). In this section, we describe its technical properties and its current state of the art. In particular, the distinction between "Dutch" LOFAR and "international" LOFAR -and the technical problems associated with each -will be made explicit in this section.

LOFAR is a SKA pathfinder instrument, which means that it serves not only as a cutting-edge instrument in its own right, but does so with the explicit aim of serving as testing grounds for technologies & techniques which could be usefully implemented in the SKA. It is in this context, for example, that trailblazers such as NenuFAR [START_REF] Zarka | Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics[END_REF], a low-frequency extension of LOFAR, are tested. LOFAR is an interferometric array, with each of its elements made up of antennas combined into a phased array, themselves forming an interferometer. It thus consists of antennas which are combined to form stations, which are themselves distributed throughout the Netherlands and Europe.

There are two bands to LOFAR, which are known as LOFAR-HBA (High-Band Antennas) and LOFAR LBA (Low-Band Antennas). Figs. 1.7a and 1.8 show the layout of the 6 innermost core stations and the French international LOFAR station, respectively.

We see the presence, in both cases, of two very different antenna types. One of these antenna types is not a single antenna, but rather a phased array: 16 antenna dipoles distributed in a 4 × 4 array. International LOFAR stations include 96 such tiles, distributed in a "filled-disk" pattern, where the disk fill ratio is determined by the size of the tile. Dutch "remote" stations have 48 such tiles spread into 2 disks (HBA INNER and HBA OUTER), while "core" station HBA tiles are split into two 24-tile crosses. This pattern is shown in Fig. 1.9a. The antennas from each tile are combined into a phased array with a single "tile beam", and all tile beams are themselves combined as a phased array into a station beam 6 .

LBA antennas, use a very simple -and cost-effective -design. A typical LBA antenna is shown in Fig. 1.10. 96 such antennas are spread in a semi-random pattern in each LOFAR station 7 . In Dutch stations, 48 of these are used at all times. They are combined as a single phased array. LOFAR-HBA is sensitive to higher frequencies, from 120 MHz to 240 MHz 8 . It has a smaller field of view than LOFAR-LBA, but a better angular resolution.

(a) LOFAR core, known as the Superterp.

(b) Distribution of the so-called LOFAR "core" stations, which include the Superterp.

(c) LOFAR with both "core" and "remote" stations. The dipole design frees observers from the need to physically point antennas at all: the final station pointing is achieved by digitally introducing delays in the recorded phase before averaging the data. In this sense, the pointing is achieved in a similar way as individual HBA tile pointing, albeit . The LBA antennas are receptive to signals emitted in the 30-80 MHz frequency range.

At the time of writing, use of the LBA data is still relatively new, as its calibration is a very tricky problem. For similar reasons, international LOFAR (i.e. the full LOFAR array) has not been used, at the time of writing, to create wide-field survey images. A large part of the work described in this manuscript consists of reaching a point where full use can be made of international LOFAR, in a streamlined and repeatable way.

This Thesis

The work done in this thesis can be summarised as the development of a new algorithm for mitigating calibration artefacts in radio-interferometric images, and its application to a peculiar famous extragalactic field. Because our data volumes put radio interferometry in the domain of "Big Data" (there are ∼245TB of LOFAR data centred on the Extended Groth Strip, for example), there is a very real need for fast, reliable methods with which to improve interferometric images. This algorithm follows from an improved understanding of the noise properties of interferometric images, and can be understood as a way to better constrain the flux distribution in images by minimising the presence of artefacts. This, in turn, helps self-calibration convergence by improving the image resulting from each imaging step of the self-calibration process.

We then outline the application of this work, along with other modern tools, to the creation of a large, deep and high-resolution image of a rich patch of sky. While this approach is potentially extremely rewarding in terms of the science one can do with its resulting maps, it is also uniquely complicated. There exist many LOFAR HBA (∼120-240 MHz) sky surveys (cf. Williams et al. 2016;Shimwell et al. 2017), but they do not yet make use of international LOFAR, and are thus limited in resolution (∼ 5 ′′ ). And while there are surveys made at LOFAR frequencies (∼150MHz) with other instruments (e.g. Intema et al. 2017), they do not match the resolution attainable with international LOFAR. The aim of our image is then to show that the international LOFAR stations can be used to create a radio image that matches HST/ACS resolution -which would allow for significant cross-matching and source analysis, since there is HST data available for the field being imaged, the Extended Groth Strip (EGS). This would not be a true survey, but more of a proof of concept: the application of new tools and the development of an appropriate strategy to image a large patch of the sky with the full LOFAR array. This image would have better resolution than FIRST (see Helfand et al. 2015, and references therein), which has a 4 ′′ resolution, relatively low sensitivity, and no short baselines (thereby missing diffuse emission and extended structure) at ∼1.4 GHz, while still being sensitive to diffuse emission due to the use of the LOFAR core stations.

Making this image is possible only because LOFAR observations of the Extended Groth Strip includes an extremely bright source, 3C295, in their field of view. This is both a blessing and a curse: unless this source is properly modeled and subtracted, it will pollute images of the EGS. However, because it is so bright, it can serve as a good calibrator source for the LOFAR international stations -if a good enough model of 3C295 can be made and used to calibrate these stations, then high-resolution images of the EGS becomes attainable. This would give access to a large population of radio AGNs at high resolution, observed at a large range of frequencies, and many of which lie at redshifts near the peak of star formation; this improvement of our low-frequency understanding of these objects could, in turn, lead to an improved understanding of these objects more generally.

To really be able to perform wide-field images using the full LOFAR array, we cannot rely on the happy circumstances we find ourselves in for our applications project: a much more thorough and reliable approach is necessary. The work of the LOFAR long-baseline team (e.g. Jackson et al. 2016) shows that kind of approach: in our project, we test specific strategies and approaches which are complimentary to their work, but cannot be applied more generally without their results.

Algorithmic Developments

A significant portion of the work done as part of this thesis was on investigating the relationship between calibration gain errors and calibration artefacts in the image-plane. What started as an attempt to properly write a weighting scheme to down-weigh poorlycalibrated visibilities9 eventually became a formal relationship between these quantities, showing both the behaviour of image-plane variance as a function of direction but also of pixel covariance over the full image. This work is shown in much greater detail in Chapter 3.

The key finding of our algorithmic work is a fundamental relationship between the covariance matrices of residual visibilities and the covariance of pixels in the image-plane: the "Cov-Cov relationship" between visibility covariance and image-plane covariance. The pixel variance and covariance in the image-plane are determined by a "noise-PSF", convolved with every pixel in the sky. This also means that thermal noise in interferometric images is correlated between pixels, and that this correlation is given exactly by the PSF: in the limit of an empty, non-deconvolved sky, we therefore properly model the thermal noise properties of its interferometric images. This noise-PSF is the product of the Fourier transforms of the calibration gain covariance matrix with each cell mapped not from uv space to lm coordinates10 but rather between their respective covariance spaces -from a new differential Fourier plane (henceforth "(δuδv)"-plane) to the imageplane covariance space δlδm. This image-plane covariance space describes the variance in each pixel and the covariance between pixels11 . It describes the expected calibration artefacts and noise level around each source, does not vary as a function of direction, and is convolved with each source in the field to yield the final error map. Because all unwanted (in our case, unphysical) signal can be thought of as noise, we will refer to the pixel variance map as the "noise-map".

Based on this theoretical finding, we develop a weighting scheme algorithm with which to improve the noise properties of interferometric images by down-weighting either visibilities with high variance in their calibration residuals, or visibilities which are strongly correlated to others in their calibration residuals. We implement a robust version of the first weighting scheme, which relies on an antenna-based estimate of the visibility variance.

Imaging the Full Primary Beam

In this section, we briefly describe previous observations of the Extended Groth Strip (EGS), an extragalactic field with a rich multi-wavelength coverage described in Table 1.1. It has been long observed as part of the All-Wavelength Extended Groth Strip International Survey collaboration (AEGIS, Davis et al. (2007)), which later became part of the CANDELS collaboration (Grogin et al. (2011)). The field is centred at α = 14 h 17 m , δ = +52 o 30 ′ , placing it between the tail of Ursa Major and Draco. Its size is 0.7 ′ × 0.1 ′ . It has notably been the subject of very deep Hubble Space Telescope observations, which the international LOFAR could attempt to match. Our data consist of an observation centred on the EGS. The LOFAR beam (i.e. how much of the sky LOFAR is sensitive to at a given time) is much larger than the EGS. As such, our measured visibilities includes contributions not just of sources in the EGS, but also from those all around it -including a very bright 3C 12 source, which will be discussed at length below. Due to technical limitations (i.e. available memory), it is not feasible to image the full LOFAR primary beam at once at international LOFAR's 0.1 ′′ resolution. As such, we begin by imaging at a lower resolution, and thus do not use the international stations at first. This allows us to create a low-resolution model of all the sources in the primary beam before proceeding to international station calibration, giving us a sanity check for core and remote station gains. It also allows us to subtract out-of-field sources when imaging the EGS at high resolution. Initial direction-independent calibration is done using Prefactor, the same pipeline as for the Lofar Two-metre Sky Survey, followed by direction-dependent calibration with the killMS-DDF facet-based pipeline (LOTSS, Shimwell et al. 2017). This work, and its results, are shown in Chapter 4. Note that the AEGIS observations summarised in Table 1.1 do not cover the full LOFAR primary beam. The sources observed at low resolution are therefore matched with Sloan Digital Sky Survey (SDSS, York et al. (2000)) images. Overlays are shown to highlight the wealth of interesting sources, and their optical counterparts, which can be picked up as a byproduct of highresolution maps of the radio sky. We select 12 particularly interesting sources and show the results of source matching for these candidates. Overlays of the same SDSS images with NVSS (NRAO VLA Sky Survey, Condon et al. (1998)) postage stamps are also shown for comparison. Many of the chosen sources do not lie within the Extended Groth Strip, and are therefore only imaged at low resolution, but high-resolution followup images and observations could be made should the need arise.

3C295 and the Extended Groth Strip

The brightest source in the primary beam of our observation is 3C295. It is one of the brighter sources in the Northern radio sky. It lies less than a degree away from the centre of the Extended Groth Strip, which has two important consequences: first, that it can potentially be used as a calibrator source for LOFAR; second, that unless well-modeled and subtracted at the observing instrument's maximum resolution, its sidelobes will pollute the EGS such that very little scientifically useful information might be recovered from a given observation. If we wish to image the EGS with international LOFAR, we therefore cannot rely on a low-resolution model: small errors in the exact location of 3C295's flux could corrupt most of the very long-baseline visibilities, thus making them useless. A high-resolution model of 3C295 is necessary for us to proceed.

The first step towards acquiring a deep, high-resolution image of the Extended Groth Strip is therefore to create a good, high-resolution model of 3C295 at LOFAR frequencies and for long baselines. Creating this model is one of the key results of this thesis. Starting from a high-resolution VLA model of 3C295, subbands selected across the LOFAR bandwidth are self-calibrated in iteratively greater numbers until a given noise threshold is reached. The flux scale is then boot-strapped based on the results of Scaife & Heald (2012). This ensures that 3C295 will be adequately subtracted from all corrected visibilities -including international baselines -and thus that images made using these visibilities will not be polluted by sidelobes from 3C295. The work done towards this goal is shown in Chapter 5.

Before launching the EGS imaging run using all calibrated datasets, we measure the impact of direction-dependent effects on known sources in the field. Two directiondependent effects are expected to dominate: decorrelation and differential gains. The first is modeled by the algorithms used in this thesis, whereas the second is a function of the gain variability as a function of distance from the calibrator. Eight calibrator sources, lying at various distances and directions from the calibrator source, are examined. They are summarised in Table 5.2, which is reproduced below: note that number 6 is actually 3C295. Their distribution in the sky is shown in Fig. 1.11.

If directional gains are measured to have little or no impact, then direction-dependent calibration of the international stations will not be required. This is the ideal scenario, as it does not require the use of fringe fitting (Pearson & Readhead 1984). If they are measured to have a significant impact, then direction-dependent calibration will be required for international LOFAR. We expect the angular scale of differential gain evolution to be of the order of a degree, and so most of the EGS should be relatively unaffected.

Chapter 2

Techniques of Interferometry Contents

In this section, we contextualise the work done as part of this PhD by discussing the difficulties introduced by the combination of sparse uv-coverage and weak a priori constraints on the sky brightness distribution. The conceptual framework which we rely on throughout this manuscript, the Radio Interferometer's Measurement Equation, is described in greater detail in Section 2.7. We will begin by discussing radio interferometry itself. We will then discuss the so-called imaging problem, and end by discussing the calibration problem. Note that both of these are artifical subsets of the general inverse problem of radio interferometry. While in practice calibration is done before imaging, it is conceptually helpful to begin with imaging. Until we reach our introduction to radio interferometric calibration, we will therefore assume that calibration has been successfully carried out, and that we are working on the corrected visibilities, i.e. gain-corrected visibilities.

We start by linking interferometry to more concrete concepts: specifically, we will give a (very!) brief introduction to radio antennas and their characteristics. We will then use this introduction to show the advantages and drawbacks of radio interferometry, along with the concrete technical problems that the method introduces.

A Brief Introduction to Radio Astronomy

Radio astronomy consists of observing the electromagnetic emission of astrophysical1 sources at very long wavelengths by means of radio antennas. These antennas measure a voltage proportional to variations in the electromagnetic field in all the directions they are sensitive to. Since astrophysical signals are weak, a good sensitivity is necessary. Achieving good sensitivity with radio antennas means having a very large collecting area -in this respect, they behave exactly the same way as optical telescopes. Similarly, when well-designed, they are diffraction-limited -this means that, once again like optical telescopes, their resolution is limited by their diameter.

However, because radio frequencies are so much lower, achieving a resolution comparable to e.g. the HST requires extremely large dishes. While there exist telescopes, both old and new, which work on this principle (from Arecibo Observatory, shown in Fig. 2.1a, to the upcoming FAST telescope in China, shown in Fig. 2.1b), the associated technical difficulties (pointing is complex business, as are the associated optics requirements; maintenance costs are high, etc...) have made the single-dish approach prohibitive. One way to move past this technical limitation is to resort to interferometry: this is the topic of this section.

Interferometry: Bypassing the Diffraction Limit

There are two quantities of interest to all astronomers: sensitivity and angular resolution 2 . An instrument's sensitivity is a function of its collecting area 3 . Resolution, for well-designed instruments 4 , is limited by diffraction in the absence of atmospheric effects. This introduces specific issues in the radio domain. Radio waves have very long in the sky, one must introduce a phase delay in each antenna (or, for fundamental interferometer elements such as the LOFAR tiles or NenuFAR, by playing with the cable length between each dish and the correlator). This is done by three means: delay tracking (which is exact only for the delay tracking centre6 , and corresponds to playing with the cable length between antennas in the interferometric array), antenna pointing7 , and fringe stopping8 . This means that there are actually three different "field centres", which are usually (but not always) set to be the same: the pointing (direction of maximum antenna sensitivity), the delay (which can be introduced in cable length) and the phase tracking (which consists of applying a known signal to the measured signal and recording the resulting heterodyne wave so as to keep a reasonable amount of post-correlation data without losing information on sources outside of the phase delay centre). The position of sources are then given in terms of the cardinal9 angle σ: this is the difference between the position of a source and the phase centre, which is where the array is pointing. For example, in Fig. 2.4, if we point antennas and introduce phase delays such that we are pointing the array towards Source 2, the position of Source 1 will be σ.

The uv-plane

For the sake of brevity, radio astronomers tend to talk of a uv-plane rather than uvw-space to describe visibility space. The set of uvw-values for all the baselines of an interferometric array is known as its uv-coverage, and defines the array's properties entirely. For the VLA, for example, the instantaneous uv-coverage when observing the zenith will be as shown in Fig. 2.5.

Individual antennas of an array can be pointed mechanically, and so the loss of antenna sensitivity in the direction of interest (and therefore the loss of net interferometric array sensitivity) can be minimised. But what happens to the array itself? It is useful here to go back to the illustration of Fig. 2.2b. Think of each dish in the array representing a "filled" segment of a massive but empty dish. By pointing our observation in a given direction, this dish is transformed -the absolute positions of the "filled" segments remain constant, and so their relative position after transformation is changed. Our "sparse dish" goes from circular to elliptical. This is shown in Fig. 2.6.

The Point-Spread Function

We have seen that the purpose of an interferometric array is to overcome the diffraction limit of single-dish antennas. We have described visibilities, which are the quantities measured by an interferometric array. What remains is to describe how these measurements are related to the sky brightness distribution.

From Dirty to Clean: Deconvolving the PSF

Raw images made from radio interferometric data consist of the underlying flux distribution convolved to the array's PSF. The quantity of interest to astronomers is the underlying flux distribution. To recover that information from images, deconvolution is needed.

In practice, performing a full deconvolution on large images (easily over 25 million pixels in total) is not computationally viable. Scientists therefore resort to algorithms and techniques to accelerate imaging and deconvolution. For example, the use of Direct Fourier Transforms (DFT), which are extremely slow, is avoided when going from visibility space to image space -Fast Fourier Transforms (FFTs) are preferred. Similarly, one seeks to avoid having to perform full deconvolution on the dirty image, since many pixels in the image do not actually contain astrophysical signal. In this chapter, we will briefly cover the most common method of deconvolution used in radio astronomy.

We use dynamic range (DR) as a quality metric in radio images. DR is the ratio of the flux of the brightest source in the field to the flux of the faintest source in the field detectable above the noise. It is therefore a function of both source brightness and noise, and so we define dynamic range as:

DR = S max max(S min , σ) (2.10)
where S max is the flux of the brightest source in the image, S min the flux of the faintest source detectable above the noise, and σ is the noise in the image. The higher the DR, the deeper we image the sky, and thus the better the image. There are two big constraints to reckon with: firstly, we do not want to deconvolve the PSF from all pixels in the image, but only from those within which large amounts of underlying flux fall. Secondly, if we only deconvolve some pixels rather than performing a full deconvolution, then care must be taken not to treat artefacts in the field (which can be created from overlapping PSF sidelobes from neighbouring sources, even with perfect calibration) as true flux. Doing so would lead to an attempt to deconvolve a PSF from a pixel while assuming an incorrect underlying flux value for this pixel and thus result in the introduction of further artefacts in the image.

The dominant family of deconvolution algorithms used in radio astronomy at the time of writing is CLEAN10 (see Högbom 1974;Pearson & Readhead 1984), which is understood in compressed sensing theory as a matching-pursuit algorithm optimising L2 (i.e. a least-squares fit) with the addition of an L1 regularisation (i.e. the total flux in the model is minimised). The various CLEAN algorithms seek out the brightest pixel(s) in the image, potentially applying a mask to the image first in cases where some a priori knowledge on flux distribution is available. It then deconvolves those pixels up to a predefined threshold (either a fraction of the initial pixel value or some factor of the estimated noise value). It does this some predefined number of times. Then, it collects the aggregate model (the flux estimate for each deconvolved pixel up to this point), Because the imaging problem (i.e. the problem of deconvolution) is linear, it can only be convex -though it is not necessarily well-conditioned. As for the calibration problem, experience strongly indicates that it is nearly convex and ill-conditioned. For calibration convexity, Cyril Tasse (pers. comm.) has found that the Hessian matrix of the gain solutions is positive semidefinite, i.e. that all its eigenvalues are positive and non-zero except for one. This behaviour likely explains why self-calibration is, in the practical experience of radio-interferometric practicioners, neither well-conditioned nor convex: if the problem is poorly regularised, then the same initial model may give different end results, and to truly verify the convergence property of self-calibration would require both a full Monte Carlo simulation and enough time to allow the very slowly-converging procedures to finish converging. It is not a practical thing to test rigorously.

Because the net calibration & imaging problem has so many degrees of freedom, it is convenient to break it into two parts: calibration and imaging. By doing this, instead of exploring a massive set of hyperparameters (e.g. relative flux distribution at a given frequency, overall spectral index, the 8 values of a complex 2 × 2 matrix per antenna per time and frequency interval, ...), a smaller subset of these hyperparameters is explored while all else is kept constant. By ensuring that these subsets each give reasonable results individually, it becomes possible to solve the full problem. This practice of imposing specific priors to each subset of the problem is a form of regularisation: it can turn a poorly-conditioned problem into a better-conditioned one, provided that the priors are well-chosen.

In effect, this means going from exploring the sky-gains parameter space all at once to exploring it along one "axis 12 " at a time. This is illustrated in Fig. 2.12.

In this context, a good deconvolution algorithm is one that is able to minimise how much it is affected by calibration errors while maximising how much true underlying flux it can recover into a model, while a good calibration algorithm is one that is able to minimise how much it is affected by sky model errors while maximising how much true underlying gain structure it can recover.

Weighting Schemes

We will finish our introduction to the inverse problem of interferometry by discussing one method by which the conditioning of deconvolution can be improved. Consider two extreme cases: a field containing a few point sources(Fig. 2.13a), and a field containing a diffuse and turbulent source, with very complex flux distribution at all scales (Fig. 2.13b). Deconvolving the PSF from the first image using CLEAN algorithms will be child's play, but doing the same with the second will be extremely complex. The first case is wellconditioned, and the second to ill-conditioned.

The better the conditioning, the easier the deconvolution and the closer to the ground truth its result. As such, when the true underlying structure of a source is poorly-known, We close this section by noting that data flagging can be considered a form of weighting scheme. This consists of assigning a null weight to those data points considered too corrupted by various instrumental or atmospheric effects (typically, Radio Frequency Interference -RFI -cf. Offringa (cf. 2010) for an example) to be scientifically useful. Those weights which are deemed too corrupted are assigned a weight of zero. This corresponds to simply dropping those data points which would corrupt the image reconstruction by introducing unphysically strong fringes in the image.

Flagging is also useful to remove corrected visibilities for which near-zero gains are found: the associated corrected visibilities will consist of some number divided by a nearzero number, for which numerical errors can introduce very large errors. "Clipping" these visibilities after calibration helps improve the final images.

An Observer's Perspective of Radio Interferometry

Astronomers who request observation time on interferometric arrays are generally not experts in the theory and operation of these arrays. This is a simple consequence of scientific division of labour. This thesis is written from the perspective of an "expert" user: one who specialises in understanding interferometric arrays and reducing their data, rather than their subsequent analysis. However, it can be very helpful to contextualise this perspective in the experience of other scientists, which is what this section aims to do. This is simply to highlight the conceptual structure of a radio interferometric observation, and where the concepts outlined in this chapter fall in the overall pattern. In other words, it is best understood as a guide to see where a given section's topic slots into the overall production of an image made from visibilities.

LOFAR data can be acquired through observation or via the Long-Term Archive. The worst of this data will already be flagged. Users can thus generally proceed directly to calibration. Because there is no guarantee of sufficient signal-to-noise in the "target field" (the object that the astronomer is actually interested in), it is standard practice to observe a calibrator source for 5 minutes before and after each observation. Such sources typically need to be bright and unresolved for the interferometer being used. Because they are bright and at phase centre for the 5-minute observation, the SNR when solving for calibration solutions for these 5 minutes will tend to be quite good. These solutions will then be interpolated between the 5 minutes onto the 8-hour observation. Calibration will be described in Section 2.6. So far, the astronomer will have worked entirely with visibilities. However, generally speaking, the aim will be to get information on the object's spatial brightness distribution at the observing frequency. This will require imaging: going from visibilities to the image-plane. This can be done through a variety of tools, but will generally involve deconvolution. All that was covered so far was imaging: specific examples of deconvolution itself were given in Section 2.2.5.

After obtaining the initial image, it may be desirable to perform a few rounds of selfcalibration (cf. Brogan et al. 2018) on the target field. This consists of extracting a model from the initial image and using it as a new calibration model, then re-imaging. Doing this will usually dramatically improve the final image. This now requires an understanding of interferometric calibration, which is the topic of the next chapter.

Calibration Methods in Radio Interferometry

In this section, we will discuss the implementation of interferometric array calibration. Our analysis is based on the RIME formalism, described in Section 2.7. One key metric of calibration quality is the dynamic range, mentioned in Sec. 2.2.5. High dynamic ranges mean that a high contrast has been obtained, and fainter sources can be reached.

Here, because the difference between thermal noise and artefacts becomes relevant, we redefine dynamic range as follows

DR = S max max(S min , σ thermal , σ artef acts ) (2.11)
where S max is the flux of the brightest source, S min the flux of the faintest source detectable above the noise, σ thermal the thermal noise in the image, and σ artef acts the noise associated with calibration artefacts. This gives a definition for dynamic range useful even in fields with a single source. We can thus use it as a metric for calibration quality specifically.

The distinction between these two noise sources is crucial; one can never go 'below noise' for a given observation, no matter the quality of calibration. Astronomers typically observe for longer periods of time in order to reduce σ thermal in their images, but this will not reduce the artefacts caused by poor calibration solutions. Uncorrected Direction-Dependent Effects13 will not go away on their own, no matter how long the integration time. Similarly, there is a limit to how much improving calibration will improve the final image: eventually, more data is required to drive noise down.

We can distinguish three 'generations' of calibration methods, of increasing complexity (Noordam & Smirnov 2010). We will describe them in terms of the RIME, showing how each generation increases in generality to account for more exotic effects. They are referred to interchangeably as 'nth-generation calibration' or 'nGC' methods.

Generational Analysis

First-Generation: Open-Loop Calibration First-generation calibration methods (1GC methods) consist of open-loop calibration. This relies entirely on instrument stability, and thus imposes significant design constraints on radio telescopes. It consists of briefly observing an external calibrator before and after each observation run to find calibration solutions for those calibrator observations, and then interpolate between them over the "target" observation. In RIME terms (see Section 2.7), this consists of solving for a very basic form of G p :

G p = a p I (2.12)
where a p is a complex constant solved for during open-loop calibration and I is the unit matrix. By interpolating between the two calibrator observations, it is possible to have a linear time-dependent gain estimate.

While values for a p can in theory be found for both autocorrelation and both crosscorrelations (e.g. XY and YX for linear receptors), low signal-to-noise means that in practice, a single set of values is solved for per antenna14 . With these techniques, one can achieve dynamic ranges of about 100:1 ((Noordam & Smirnov 2010)).

Second-Generation: Self-Calibration

Second-generation calibration methods (2GC methods) are defined by self-calibration15 , commonly referred to as self-cal (Pearson & Readhead (1984)). As described in Section 2.11.2, this method can only be deployed if the brightness matrix of the sky is the same for all baselines16 (i.e. that we are not affected by direction-dependent effects).

The first instance of self-calibration (and adaptive optics in radio interferometry) is a paper published in the era of 1GC by Jennison [START_REF] Jennison | [END_REF], expanding on his PhD work, which was published in 1951). With sufficient signal-to-noise, he showed that phase closure could be calculated and errors due to the atmosphere thus mitigated.

Self-calibration and adaptive optics are essentially the same thing: an attempt to correct for phase errors in the signal measured by an instrument, which have been introduced by the atmosphere (i.e. correct time-varying and direction-dependent processes). This is shown in Fig. 2.15. Self-calibration and adaptive optics nevertheless meet different constraints: most notably, because interferometric data is digital, radio astronomers can perform their corrections after the observation.

If amplitude and phase gains can be written as antenna-dependent, then each antennabased error is estimated N times (once per each baseline which includes this antenna). By estimating, and correcting for, these errors, one can use a simple source model to infer an improved one. This is why the method is referred to as self-cal; by calibrating 'on' a good calibrator source (bright, compact and unresolved), one can drastically improve one's source model, along with one's calibration solutions.

By extending this idea to VLBI (see e.g. (Rogers et al. 1974)), amplitude closure was introduced to the field along with phase closure (as described in Rogers et al. 1983). These quantities are immune to antenna-based effects.

The great advantage of radio interferometry over optics, however, is that we can iterate over progressively improved source models. This is because we have the ability to record phase information.

The RIME Formalism

In this section, we introduce the mathematical framework which forms the basis for our algorithmic work. While radio interferometry has historically been a complex business [START_REF] Thompson | Interferometry and Synthesis in Radio Astronomy[END_REF], modern radio interferometry is based on the Radio Interferometer's Measurement Equation, a powerful and elegant formalism which underpins modern calibration and imaging algorithms. In this section, I will give a simplified account of this formalism as described in Smirnov (2011a) and its companion papers (Smirnov 2011b,c,d). This will form the theoretical basis on which further sections will expand.

For a mathematically rigorous description of the RIME, the references used in this work (particularly Smirnov (2011c) and Smirnov (2011d)) are obviously the first place to look. We frame the RIME in the continuity of previous theoretical frameworks of radio interferometry, while also describing it in the context of contemporary software packages and algorithmic tools.

Setting up the RIME: a single point source

The Radio Interferometer's Measurement Equation, or RIME, is a formalism which allows for an elegant and efficient formulation of the physical processes which affect signal propagation, from astrophysical effects to instrumental effects. Its fundamental underlying hypothesis is linearity: that transformations along the signal are linear with respect to the basis chosen to represent our signal. In other words, this means assuming that propagation effects are separable by antenna.

Consider a single quasi-monochromatic point source in the sky. Its signal at a point in time and space can be described by a complex 2-vector e -this assumes plane waves. We can then represent physical processes which affect this signal's propagation using Jones matrices. In other words:

v = J e (2.13) e = ⎛ ⎜ ⎝ e x e y e z ⎞ ⎟ ⎠ (2.14)
where v is the voltage measured by our antenna, and J the Jones matrix describing the net propagation effects -atmospheric, astrophysical, and instrumental -which affect our source's signal. e is the electromagnetic wave emitted by a source. J can be written as a series of matrix products, each individual matrix describing a physical phenomenon. They are then called a Jones chain, and the resulting J is often called the total Jones matrix.

For e to be a 2-vector, rather than a 3-vector, we must choose to use xyz as our coordinates, with z the direction of propagation of our signal. Hamaker et al. (1996) then show that the following relation holds:

2ee H = 2 < e x e ⋆ x > < e x e ⋆ y > < e y e ⋆ x > < e y e ⋆ y > = I + Q U + iV U -iV I -Q = B (2.15)
where B is the coherency matrix of the source's emission.

A visibility is the correlation of voltage from two antennas. In other words,

V pq = 2⟨v p (v q ) H ⟩ (2.16)
V pq = 2⟨J p e(J q e) H ⟩ (2.17)

V pq = 2⟨J p ee H J H q ⟩ (2.18)
where J p corresponds to the signal propagation chain of antenna p, and J q for antenna q. (...) H corresponds to taking the Hermitian conjugate of a matrix, and ⟨...⟩ an average over some interval. Note that a factor of 2 has been introduced here, as a matter of convention. This is due to the definition of the Stokes parameter in relation to our e outer product.

We can thus write:

V pq = J p BJ H q (2.19)
This gives us an elegant formulation of the relationship between the signal emitted by an astrophysical source and its corresponding measured visibility. The Stokes parameters18 of our source can be directly calculated from our measured visibility, provided J p and J q are known. Of course, in practice, they are not, and must therefore be modelled.

To do this, we must expand each Jones matrix into a corresponding Jones chain.

Expanding our Jones Chain

Having written our RIME as a matrix multiplication problem, we can now begin to differentiate between different types of Jones matrices. J p describes the aggregate effects which occur over the course of the signal's propagation to our antenna p. It can be split into a Jones chain of n separate matrices, corresponding to different effects:

J H sp = n 1 J H nsp = J H 1sp J H 2sp ...J H nsp (2.20) J sp = J nsp ...J 2sp J 1sp (2.21)
where s is a particular direction in the sky: these Jones matrices are therefore directiondependent. Going from Eq. ( 2.20) to Eq. ( 2.21) simply consists of applying the Hermitian operator to both sides of Eq. ( 2.20) to recover J nsp . Since the Jones terms are applied sequentially to the brightness matrix, starting with J 1sp , the terms to the right of our decomposition of J sp can be said to occur 'at the source', i.e. before those to the left (said to occur 'at the antenna').

Our aim is thus to find the Jones chain which most concisely describes the different types of effects affecting our signal. In order of increasing complexity, there are three: so-called geometric effects, direction-independent effects, and direction-dependent effects.

Geometric Effects: the Phase Delay Jones Matrix

The physical quantity measured by an interferometric array, the visibility, is not a function of the absolute phase of our signal, but rather a function of the phase difference between our measured voltages v p and v q . In other words, it is a function of the phase difference on baseline pq, which consists of antennas p and q. The direction which minimises this difference is called the phase centre19 .

We shall use the conventional coordinate system and notations (Smirnov 2011a;Thompson et al. 2001, etc). We thus have a z axis pointing towards the phase centre, and an antenna p at coordinates u p = (u p , v p , w p ). The phase difference between p and phase centre (where, by definition, u = 0) for a signal arriving from direction σ σ σ is then, assuming a small field of view:

k p = 2πi((u p l + v p m + w p (n -1))) (2.22) where n = √ 1 -l 2 -m 2 .
Here, l, m, n are the direction cosines20 of σ σ σ, and u is defined in units of wavelength21 .

Having done this, we can now introduce a scalar K-Jones matrix to represent the effect of this phase delay. Since it is scalar, it will commute with all other Jones matrices (see Smirnov 2011a). It can be be written as:

K p = e -ikp = e -2πi(upl+vpm+wp(n-1))
(2.23)

The Coherency Matrix

The commutation property of K p with other Jones matrices allow us to always place it at the rightmost end of our Jones chain. We can thus write our Jones chain as follows:

J p = G p K p (2.24) (2.25)
Here, we have decoupled the nominal phase delay term (K p ) from all 'corrupting' effects (G p ). We can take this one step further by defining the source coherency, X pq , as follows:

X pq = K p BK H p (2.26) V pq = G p X pq G H q (2.27)
This is useful on a conceptual level: the process of calibration is now to find a best estimate for G p , which contain all non-'geometric' physical effects on the signal's propagation from an astrophysical source to our voltage measurement.

On a practical level, using the coherency matrix X allows us to write our RIME more concisely, by absorbing the geometric effects into the core of our RIME. While not strictly necessary, it can be a useful thing to do.

Multiple Point Sources

So far, B represents the signal from a single point source. However, the formalism can be extended to multiple point sources. This is because electric fields -and therefore radio signals -are additive, and that signal from different sources is not coherent (see Eq. (2.5)), which means that visibilities are additive. If our sky contains s point sources, then our RIME becomes:

V pq = s J sp B s J H sq (2.28)
We now expand our Jones chain as before. Some, but not all, elements in the chain will only depend on the antenna. These will be the same for all sources. We can thus write our Jones chain as:

J sp = G p E sp K sp (2.29)
where E sp are direction-dependent effects, which vary for each source. Note that we have kept the source-dependent effects (i.e. those which vary depending on the line of sight) are to the right of our chain, and the antenna-dependent (G p ) effects to the left. K sp is a scalar and can thus be moved at will; it is therefore placed at the far right to recover the coherency matrix if necessary.

Our multiple-source RIME can thus be written in the following form:

V pq = G p s E sp K sp B s K H sq E H sq G H q (2.30)
This 'onion' form of the RIME -so-called because, like an onion, it is formed in layers -is usually sufficient to describe signal propagation in radio-interferometric arrays. We find the three types of Jones matrices promised earlier: the geometric effects accounted for by K , the direction-dependent effects modeled into E , and the direction-independent effects (antenna-based effects) modeled into G.

There are still limitations to this RIME; notably, it assumes that the sky consists entirely of point sources (absence of diffuse emission), and that Jones matrices are timeindependent. These are of course unphysical hypotheses, though they are a very good first approximation. We shall now address these two points.

Diffuse Emission: the Full-sky RIME

Deriving the Fully-Sky RIME

If we want to be able to take diffuse emission properly into account, we must go from our set of discrete point sources B s and integrate the underlying continuous brightness distribution B(σ σ σ) over all possible directions. This gives us the following expression:

V pq = 4π J p (σ σ σ)B(σ σ σ)J H q (σ σ σ)dΩ (2.31)
Here, effects such as our array's beam lie in the Jones matrices. By performing the analysis shown in Sect. 3.1 of [START_REF] Thompson | Interferometry and Synthesis in Radio Astronomy[END_REF], we can rewrite this expression in terms of a sine projection of the sphere onto the (l, m) plane tangential at field centre. The integral then becomes:

V pq = l,m J p (σ σ σ)B(σ σ σ)J H q (σ σ σ) dldm n (2.32)
Having written this, and using l l l as shorthand (l, m), we can once again decompose our Jones matrix into a Jones chain containing a direction-independent term G, a directiondependent term E , and the phase term K . Substituting this into our integral, we get the following expression:

V pq = G p l l l 1 n E p (l l l)K p (l l l)BK H q (l l l)E q (l l l)dldm G H q (2.33) V pq = G p l l l 1 n E p (l l l)BE q (l l l)K p (l l l)K H q (l l l)dldm G H q (2.34) = G p l l l 1 n E p (l l l)BE q (l l l)e -2πi upql+vpqm+wpq(n-1) dldm G H q (2.35)
where

u pq = u p -u q v pq = v p -v q w pq = w p -w q (2.36)
This is a three-dimensional Fourier transform. It can be reduced to a two-dimensional Fourier transform by writing the n-dependence in the exponential (and normalisation) as a Jones matrix, thus removing the non-coplanarity terms from the explicit RIME. We can thus define the following:

W p = 1 √ n e -2πiwp(n-1)
(2.37)

B pq = E p W p BW H q E h q (2.38)
and write our RIME in its final form:

V pq = G p l l l B pq (l l l)e -2πi upql+vpqm+wpq(n-1) dldm G H q (2.39)
Note that B was the actual sky brightness distribution as a function of direction, but now B pq is the sky brightness as seen by baseline pq -in the presence of directiondependent effects, this is not the same for all baselines.

TIME-DEPENDENT JONES MATRICES, SMEARING, AND DECORRELATION49

Recovering the van Cittert-Zernike Theorem

Earlier calibration models and algorithms rely on the hypothesis that all baselines see the same sky. This is the fundamental premise of self-calibration. The topic of self-cal is covered in Section 2.6.1. This premise only holds when all Direction-Dependent Effects (henceforth referred to as DDEs) are identical for all baselines (and therefore for all antennas). We then effectively have E s,p = 1 for all p. This is a necessary condition for the apparent sky B pq to be the same for all baselines22 , as it allows us to write:

E p (l l l)W p (l l l) = E (l l l)W (l l l) (2.40) B pq (l l l) = B app (l l l) = E (l l l)W (l l l)B(l l l)W H (l l l)E H (l l l) (2.41)
If this condition is met, we can then rewrite the full-sky RIME as:

V pq = G p X pq G H q (2.42)
where X pq = X(u pq , v pq ) = X(u u u). Our coherency matrix includes the beam attenuation implicitly, as B app is now the coherency of the apparent sky brightness distribution. Comparing Eqs. (2.26) and ( 2.39) directly shows that the coherency matrix is simply the two-dimensional Fourier transform of B app (l l l). We can thus call X(u u u) the sky coherency, in keeping with our nomenclature for the RIME of a point source.

Note that deriving the van Cittert-Zernike theorem ((van Cittert 1934), covered in [START_REF] Thompson | Interferometry and Synthesis in Radio Astronomy[END_REF])) from the RIME was simply of matter of treating phase as a Jones matrix in its own right 23 . A significant consequence of this is that DDEs can be incorporated within the RIME formalism.

2.12 Time-dependent Jones matrices, Smearing, and Decorrelation

Time-dependence of Jones matrices

Until this point, we have not considered any kind of sky variability in time. In effect, Eq. ( 2.42) describes the snapshot taken by an interferometer for a single measurement. Limiting ourselves to this scenario, however, strongly limits the sampling of the uv-plane which can be performed. Indeed, interferometers generally rely on the Earth's rotation to rotate their baselines pq, thus increasing uv-coverage at no additional cost.This technique is known as super-synthesis (Joardar et al. (2010) 24 ) For Eq. (2.42) to hold throughout an observation, we must therefore assume that B app remains constant in time over the course of the observation. In effect, this means assuming that

E p (t, l l l) = E p (l l l) = E (l l l) ∀(t, p) (2.43)
The above means that direction-dependent effects must be time-independent for Eq. ( 2.42) to hold. We can thus split direction-dependent effects into two categories: trivial DDEs, which satisfy this condition and are a simple multiplicative correction, and non-trivial DDEs, which are not. An example of a trivial DDE is the primary gain; and example of a non-trivial (i.e. time-variable) DDE is ionospheric turbulence.

Smearing and Decorrelation

Of course, our Jones matrices are not the only time-dependent element of the RIME formalism. As mentioned above, to image with interferometric arrays, radio astronomers rely on the Earth's rotation to better sample uv-space. This means that u u u becomes a time-dependent variable.

Furthermore, if we choose to define u u u in terms of wavelength, u u u p becomes dependent in frequency. Our time-independent, frequency-independent RIME thus becomes only valid for infinitesimal increments in time and frequency.

In practice, measurements made by any interferometric array are actually averaged over the integration time of the measurement, as well as the frequency channel size. A better formulation is therefore to explicitly write the RIME as an integration over time and frequency bins:

⟨V pw ⟩ = 1 ∆t∆ν t 1 t 0 ν 1 ν 0 V pq (t, ν)dvdt (2.44) = 1 ∆t∆ν t 1 t 0 ν 1 ν 0 J p (t, ν)BJ H q (t, ν) (2.45)
J (t, ν) describes the propagation of an electromagnetic signal through various physical processes. This signal has a complex phase, which is variable in both time and frequency. Since it is a rotating vector, integrating over any time or frequency band will always result in a net loss in amplitude in the measured ⟨V⟩. This mechanism is commonly referred to as time/bandwidth decorrelation or smearing (Smirnov 2011a).

The general effect is referred to as decorrelation25 , and the specific case of decoherence caused by the K term is referred to as smearing.

Smearing increases with baseline length (u u u pq ) and distance from phase centre (l, m). The noise amplitude, however, does not decrease. This means that smearing results in a decrease in sensitivity, which becomes significant when attempting to do wide-field, high-resolution images (see Section 5.6).

In the context of the RIME, smearing is an element-by-element operation, i.e. different elements do not affect each other. Treating smearing within the RIME is thus a trivial expansion of the scalar equations.

Assuming that ∆t and ∆ν are small enough that the amplitude of V pq remains constant over the integration band (i.e. that V pq is well-sampled by our instrument), then a useful first-order approximation26 is:

⟨V pw ⟩ ≃ sinc ∆Ψ Ψ Ψ 2 sinc ∆Φ Φ Φ 2 V pq (t mid , ν mid ) (2.46)
where

t mid = (t 0 + t 1 ) 2 (2.47) ν mid = (ν 0 + ν 1 ) 2 (2.48) ∆Ψ Ψ Ψ = arg V pq (t 1 , ν mid ) -arg V pq (t 0 , ν mid ) (2.49) ∆Φ Φ Φ = arg V pq (t mid , ν 1 ) -arg V pq (t mid , ν 0 ) (2.50)
This equation is convenient to determine the impact of decorrelation for a source at a given distance from phase centre. Note, however, that it can only -by definition -be valid for a single source at a time. If using a RIME such as Eq. ( 2.30), then it can be calculated for each source individually.

Noise

The RIME presented thus far has been for an ideal case, free of noise in our measured visibility. In practice, each complex visibility will be affected by uncorrelated Gaussian noise 27 , which means that a RIME of the same form as Eq. ( 2.30) would have to be formuled at follows:

V pq = G p s E sp K sp B s K H sq E H sq G H q + σ pq (2.51)
where σ pq is a (2, 2) matrix with real and complex parts in each component. A similar addition can be made for any of the RIMEs presented here. This noise sets a hard limit on the sensitivity achieved for a given observation, and 'reaching the noise' becomes the gold standard of calibration (Smirnov 2011b). One must then take σ pq into account properly when solving for the Jones matrices.

Introduction

Interferometers sample Fourier modes of the sky brightness distribution corrupted by instrumental and atmospheric effects rather than measuring the sky brightness directly. This introduces two problems2 for astronomers to invert: calibration and imaging. Both of these problems are ill-conditioned.

The problem of imaging consists of correcting for the incomplete uv-coverage of any given interferometer by deconvolving the instrument's Point-Spread Function (PSF) from images. Its poor conditioning comes from our limited a priori knowledge of the sky brightness distribution, combined with large gaps in our uv-coverage, which prevents us from placing strong constraints on image deconvolution. It can be better-conditioned in different ways, including through the use of weighting schemes (see Briggs 1995; Yatawatta 2014, and references therein) to improve image fidelity at the start of deconvolution. When inverting the imaging problem, we often assume that the sky is stable within the domain (i.e. is constant in time and frequency). There are exceptions, such as wide-band deconvolution algorithms (e.g. Rau & Cornwell 2011) that explicitly take into account the sky's frequency-dependence, but still assume that the sky brightness distribution does not vary with time.

The problem of calibration is what concerns us in this section. It consists of estimating and correcting for instrumental errors (which includes effects such as antenna pointing errors, but also the phase-delays caused by ionospheric activity, troposphere, etc). Calibration consists of solving for gain estimates, where a gain models the relationship between the electromagnetic field of an astrophysical source and the voltage that an antenna measures for this source. Because measurements are noisy, calibration often involves some fine-tuning of solution intervals, to ensure that the solutions are well-constrained while the solution intervals stay as small as signal-to-noise allows. The calibration inverse problem involves three competing statistical effects: thermal noise in the measurements, true gain variability, and sky model incompleteness. If gain solutions are sought individually for each measurement, then (unless the observed field is exceptionally bright) calibration estimates will be dominated by thermal noise, and will not adequately describe the actual gains. Similarly, if a single gain estimate is fitted to too many measurements, the intrinsic gain variability will be "averaged out"; for example, a choice of time and frequency interval that is too large will cause the solver to estimate a constant gain while the underlying function varies quickly, thereby missing much of the gain structure. This will introduce error which will be correlated in time and frequency. This occurs, for example, when solving for ionospheric phase delays: in the most extreme case, where the solution interval is significantly larger than the scale of ionospheric fluctuations, its varying phase can average out to zero over the interval in time and frequency. Finally, if the model being fitted is incomplete, unmodeled physical flux will likely be absorbed unpredictably into both the gain solutions and the residual visibilities: this absorption of physical flux into gain solutions is known as source suppression (see Grobler et al. 2014;Kazemi & Yatawatta 2013, and references therein).

In practice, it is reasonable to assume that there exist a scale at which gain variation is acceptably low: we can then reduce the noise of our gain estimates by finding a single gain solution for a small number of measurements, assuming that the underlying gain variation is very small and stable over short intervals. This is generally a valid hypothesis, but the specific value for the variation scale can be contentious. Indeed, while the noise level can often be treated as constant throughout an observation, the gain variability itself is generally not constant: there will be time periods where the gains will tend to remain constant for longer, and others where variability will be very quick. This means that, for any choice of calibration interval, some gain estimates will be better than others, and almost all could have been improved (at a cost to others) by a different choice of time (and frequency) intervals.

Since we have measurements which are better-calibrated than others (in that better estimates for their gains were obtained through chance alone), we could, in principle, take inspiration from "lucky imaging" (an optical-domain method for making good images: for more details, see Fried 1978, and references therein) to weigh our visibilities according to their calibration quality. Those weights would in effect be an improvement of currently existing methods such as clipping noisy residual visibilities: in the extreme case where all visibilities are equally-well calibrated except a few which are extremely noisy, it should be equivalent to clipping. Otherwise, the weights should show at least a slight improvement over clipping.

The key finding of our algorithmic work is a fundamental relationship between the covariance matrices of residual visibilities and the map of the covariance in the image-plane: the "Cov-Cov relationship" between visibility covariance to image-plane covariance. We show that the pixel statistics in the image-plane are determined by a "noise-PSF", convolved with each source in the sky (modeled or not). This noise-PSF is the product of the Fourier transforms of the gain covariance matrix with each cell mapped not from uv space to lm coordinates but rather between their respective covariance spaces -from a new differential Fourier plane (henceforth "(δuδv)"-plane) to the image-plane covariance space δlδm. This image-plane covariance space describes the variance in each pixel and the covariance between pixels3 . It describes the expected calibration artefacts and thermal noise around each source, does not vary as a function of direction, and is convolved with each source in the field to yield the final error map. Because all unwanted (in our case, unphysical) signal can be thought of as noise, we will refer to the pixel variance map as the "noise-map".

The notion of a (δuδv)-plane arises organically from the framework of radio interferometry: we are associating a correlation between visibilities to coordinates in covariance space, just as we associate the visibilities themselves to the uv-domain. The (δuδv)-plane is the natural domain of these correlations. As previously stated, even if all sources in the field are perfectly known and modeled, a poor choice of calibration interval can introduce correlated noise in the residuals, which would then introduce larger variance near sources in the noise-map. Conversely, if calibration is perfect, the noise-map should be completely flat (i.e. same variance for all pixels), as there would be no noise-correlation between pixels.

The main result of this section is a new adaptive, quality-based weighting scheme based on this insight. Using the Cov-Cov relationship, we can create a new weighting scheme by estimating the residual visibility covariance matrix in a given observation. By weighting visibilities so as to change their covariance matrix, one can change the shape of the noise-PSF and thus improve the final image: this manifests as either decreased noise or decreased calibration artefacts. Note that this weighting is applied after calibration, but before image deconvolution: applying it will therefore not only improve the residual noise in the image, and thus the sensitivity achievable with a given pipeline, but will also improve deconvolution by minimising calibration artefacts in the field: it should thus effectively remove spurious, unphysical emission from final data products. Estimating the covariance matrix is the main difficulty of our framework: we do not know the underlying covariance matrix, and the conditioning of our estimation thereof is limited by the number of measurements within each solution interval. As such, we have no guarantee that our estimate of the corrected visibility covariance matrix is accurate. This problem can be alleviated, for example by estimating the covariance matrix for the antenna gains themselves, and use it to build the visibility covariance matrix: this effectively improves conditioning (cf Sec. 3.4).

The remainder of this chapter is split into four main sections. In Section 3.2, we derive the Cov-Cov relationship. With its newfound insights, we propose quality-based weighting schemes with which to improve radio interferometric images in Section 3.3. We follow in Section 3.4 by showing how to estimate, from real data, the covariance matrix from which the quality-based weights are derived. Our approach seems to give good results. Finally, we close the chapter on a discussion of the applicability and limitations of the quality-based weighting scheme.

Building the Noise Map

In this section, we derive our first fundamental result: the Cov-Cov relationship, Eq. 3.25, which describes how the statistics of residual visibilities (and thus the antenna calibration solutions, henceforth "gains") relate to the statistics of the image plane, i.e. of images made using the associated visibilities. The dimensions of the matrices (denoted by boldface capital letters) and vectors (denoted by boldface lowercase letters) used in this work are given in Table 3.1, along with the scalar numbers used to denote specific dimensions. All other variables are scalars.

The Cov-Cov Relationship in the δuδv plane

Let us begin by defining visibility gains. Using the Radio Interferometry Measurement Equation formalism for a sky consisting of a single point source (Hamaker et al. (1996), Smirnov (2011a), and companion papers), we can write the following relation between the sky and the signal as measured by a single baseline at time t and frequency ν:

V tν pq = d K d p,tν J d p,tν B d ν J d q,tν H K d q,tν H + N (3.1)
All the quantities above are 2 × 2 matrices. Eq. 3.1, implies a linear relationship between the coherency matrix B d ν and the visibilities recorded by a given baseline (V tν pq ), with the addition of a thermal noise matrix N , which is also of shape 2 × 2 and contains different complex-valued realisations of the noise in each cell. Since electric fields are additive, the sky coherency matrix can be described as the sum of the contributions from individual sources in directions d in the sky. We also assume that the sky does not vary over time, i.e. that B d ν is not a function of time. ourselves to the scalar case, which corresponds to assuming that emission is unpolarised. We assume that B d = s d I, where s is the flux of our single point source and I is the 2 × 2 identity matrix. We also assume that J p,tν = g p,tν I, where g tν p is the complex-valued gains of antenna p at time t and frequency ν. This means that we assume that the gains are direction-independent, and so J d ...,tν becomes J ...,tν . Similarly, K d p,tν = k d p,tν I, the Fourier kernel in the direction of the source, d. N has 1 realisation of ǫ in each cell, where:

ǫ ∼ N (0, σ) + iN (0, σ) (3.2)
where σ is the variance of the thermal noise. Let us denote each (t, ν) pair by τ , and ignore the sky's frequency-dependence. The following scalar formulation is then equivalent to Eq. 3.1:

V τ pq = d s d k d p,τ k d q,τ g τ p g τ q + ǫ (3.3) k d p,τ = exp (2πi (u p,τ l d + v p,τ m d + w p,τ (n p -1))) (3.4)
Calibration is the process of finding an accurate estimate of g τ p for all antennas p, at all times t and frequencies ν. Since we are in a direction-independent regime, the quality of our calibration then determines the statistical properties of the residual visibilities (and the image-plane equivalent, the residual image). The residual visibilities associated with calibration solutions are defined as our measured visibilities minus the gain-corrupted model visibilities. ĝτ p then denotes our calibration estimate for g τ p . We now begin to limit the generality of our framework by assuming that all sufficiently bright sources have been modeled and subtracted: unmodeled flux is then negligible. We can then write the residual visibilities as:

rτ pq = d s d k d p,τ k d q,τ g τ p g τ q -ĝτ p ĝτ q + ǫ (3.5)
The flux values in the image-plane pixels4 are the Fourier transform of the visibility values mapped onto each pixel. This can be written as follows:

ỹ = ⎛ ⎜ ⎝ ⋮ ∑ pq I r pq,lm ⋮ ⎞ ⎟ ⎠ (3.6) I r pq,lm = τ ω pq,τ rτ pq k lm pq,τ (3.7)
where lm are the directional cosine positions of a given pixel, and k lm pq,τ = k d p,τ k d q,τ the Fourier coefficient mapping a point in Fourier space to a point on the image-plane. ω pq,τ is the weight associated to a given visibility.

Let us now write this using a matrix formalism. The contribution of a single visibility b = (pq, τ ) to the image-plane residuals can be written as:

ỹb =F H S b ω b (κ b γ + ǫ) (3.8) ỹ = b ỹb (3.9)
where ǫ is a vector of the ǫ of Eq. 3.2 and ỹ is a vector of size n pix , with

κ b = d s d κ d b (3.10) κ d b = k d p,τ k d q,τ = k lm pq,τ (3.11) γb = g τ p g τ q -ĝτ p ĝτ p (3.12) γ = ⎛ ⎜ ⎝ ⋮ γb ⋮ ⎞ ⎟ ⎠ (3.13)
and w b is the scalar weight associated with each visibility. By default, w b = 1 n b : all visibilities then have the same weight, and ỹ then becomes the average of all ỹb . γ is a vector of all γb , and thus of size n b . F is the Fourier kernel, of size n pix × n pix . S b is a matrix of size n pix × n b : its purpose is to encode the uv-coverage. Each S b contains only a single non-zero cell, different for different S b . The height (number of rows) of S b is determined the size of the uv-grid, and its length (number of columns) by the number of visibilities.

The order of operations is thus: each residual visibility (κ b γ + n) is assigned some weight w b and its uv-coordinates are set by S b . The inverse Fourier transform (F H ) is then applied to this grid, and so we recover its image-plane fringe. By averaging over all fringes, we recover the dirty image.

The residual image will thus depend on three quantities: the residual gains, the flux in the image, and the weighting scheme. Let us consider the relationship between the statistics of residual visibilities and the variance at a given point in the corresponding residual image.

Statistical Analysis

In the following analysis, we treat our gain solutions and thermal noise as random variables in order to compute the covariance matrix of our residual image, Cov{ỹ}. The diagonal of this matrix gives the variance for each pixel, while the wings give the covariance between pixels. Using the property that Cov{Ax} = ACov{x}A H , we can apply the Cov{} operator to Eq. 3.9 to write:

Cov{ỹ} = bb ′ F H w b w b ′ κ b κ b ′ def = φ bb ′ S b Cov{γ}S T b ′ F (3.14) + b w 2 b σ 2 F H S b IS T b F def = C b (3.15) = bb ′ φ bb ′ F H S b Cov{γ}S T b ′ F + b w 2 b σ 2 C b (3.16)
So far, we have only applied definitions. The net effect of S b Cov{γ}S T b ′ (dimensions of n pix × n pix ) is to encode where a given baseline samples the uv-plane, and map one cell at matrix coordinates (b, b ′ ) from the correlation matrix Cov{γ} onto the visibility grid. S b is not the gridding kernel, but rather the sampling matrix, which determines where we have measurements and where we do not. We can thus write that

S b Cov{γ}S T b ′ = [Cov{γ}] bb ′ S b 11 T S T b ′
, where [Cov{γ}] bb ′ is the value from the appropriate cell and 1 is the vector-of-ones of appropriate length (here, n b ). This allows us to write:

Cov{ỹ} = bb ′ φ bb ′ [Cov{γ}] bb ′ F bb ′ + b w 2 b σ 2 C b (3.17) with F bb ′ = (F b ) H F b ′ def = 1 T S T b ′ F (3.18)
Here, C b is a Toeplitz matrix, i.e. a convolution matrix, associated with baseline b.

The set of all C b defines the convolution kernel which characterises the Point-Spread Function (henceforth PSF) associated with a given uv-coverage, of size n pix × n pix . F bb ′ , meanwhile, is not generally Toeplitz. Its cells can be written as:

F bb ′ [d, d ′ ] = e 2iπ(u b l d -u b ′ l d ′ +v b m d -v b ′ m d ′ +(n d -1)w b -(n d ′ -1)w b ′ ) (3.19)
Let us investigate how the sky brightness distribution (i.e. d-dependence) affects the noise-map. We can write the sum over bb ′ as two sums: one over b, b ′ = b and one over b, b ′ ≠ b. Thus:

Cov{ỹ} = b φ bb [Cov{γ}] bb + w 2 b σ 2 C b + b,b ′ ≠b φ bb ′ [Cov{γ}] bb ′ F bb ′ (3.
φ bb ′ = w b w b ′ κ b κ b ′ (3.21) = w b w b ′ d s d κ d b d ′ s d ′ κ d ′ b ′ (3.22) ≈ d w b w b ′ s 2 d κ d b κ d b ′ (3.23) φ bb ′ ≈ d φ d bb ′ (3.24)
Note that F bb = C b , since those are the coordinates along the diagonal: for these values, the matrix-of-ones at the centre of F bb ′ becomes the identity matrix. Note also that

φ d bb = w 2 b s 2 d , since κ d b κ d b = 1.
We can then write Eq. 3.20 as:

Cov{ỹ} = d ⎛ ⎝ b φ d bb [Cov{γ}] bb + w 2 b σ 2 φ bb C b + b,b ′ ≠b φ d bb ′ [Cov{γ}] bb ′ F bb ′ ⎞ ⎠ (3.25)
where we have now limited our formalism to the case where the sky is dominated by sufficiently-separated point-like sources. This is our fundamental result: assuming unpolarised emission coming from wellseparated point sources and normally-distributed thermal noise, it gives a direct relationship between the covariance of the residual visibilities and the covariance of the residual image-pixel values. We thus call it the Cov-Cov relationship. It describes the statistical properties of the image-plane as the result of a convolution process changing an average noise level at different points in the image-plane, allowing us to describe the behaviour of variance and covariance in the image. Let us focus on the first By applying the Diag {} operator (which returns the diagonal of an input matrix as a vector) to both sides of Eq.3.25, we can find an expression for the variance map in the image-plane:

Var{ỹ} =Diag {Cov{ỹ}} (3.26) = d b φ d bb [Cov{γ}] bb + w 2 b σ 2 Diag {C b } =1 + b,b ′ ≠b φ d bb ′ [Cov{γ}] bb ′ Diag {F bb ′ } (3.27)
where

l d = (l d , m d , (n d -1)) (3.28) δu bb ′ = (δu bb ′ , δv bb ′ , δw bb ′ ) (3.29) = (u b -u b ′ , v b -v b ′ , w b -w b ′ ) (3.30)
In Eq 3.27, we have:

Diag {F bb ′ } [d] = e 2iπl d ⋅δu bb ′ (3.31)
For b ≠ b ′ , the diagonals of F bb ′ are the Fourier kernels mapping δuδv space to δlδm. F bb ′ can then be thought of as a Fourier transform. It is not a diagonal matrix. It behaves as a covariance fringe, allowing us to extend standard interferometric ideas to covariance space: each fringe can be thought of as a single "spatial filter" applied to the pixel covariance matrix. Just as a given baseline has coordinates in uv-space, a given correlation between baseline residual errors has coordinates in uv correlation space, which we will henceforth refer to as δuδv-space.

This δuδv space warrants further discussion: Fig. 3.1 shows, for a given uv-track (Fig. 3.1b), both the corresponding δuδv domain (Fig. 3.1c) and point-spread function (Fig. 3.1a). The symmetric, negative uv-track is treated as a separate track, and thus ignored. This means that we do not fully constrain the noise-PSF (since the covariance matrix of the symmetric track is simply the Hermitian of the first), but we do not seek to constrain it in this section, but rather to show that our results hold. We can see that the δuδv-tracks are symmetrical about the origin. The δuδv space corresponding to a given uv-track can thus be most concisely described as a "filled uv-track", with its boundaries defined by the ends of the uv-track. The set of F bb ′ , each of which maps one value of the covariance matrix to a fringe in the image-plane, would then characterise a PSF equivalent for the noise distribution, which we refer to as the noise-PSF. In our formalism, the only sources of statistical effects in the field are calibration errors and thermal noise. The average variance in all pixels will be given by the diagonal of by each baseline will determine the noise-map, and the true image-plane artefacts will then be one set of realisations of this underlying distribution.

Adaptive Quality-based Weighting Schemes

As discussed in Section 3.1, some intervals of an observation will have lower gain variability. These will show up in the gain covariance matrix as intervals with lower variance. Similarly, those with larger intrinsic gain variability will have greater error in their gain estimate. By giving greater weights to the former, and lower weights to the latter, we expect to be able to improve image reconstruction. We thus talk of adaptive quality-based weighting, as the weights will adapt based on the calibration quality.

The pixel variance is determined by the visibility covariance matrix, as shown in Eq. 3.27. The diagonal of the visibility covariance matrix will add a flat noise to all pixels, while its wings will determine the calibration artefact distribution, which will be convolved to the sky brightness distribution. We thus have two sources of variance in the image-plane. Minimising the far-field noise (i.e. the variance far from sources) in an image would involve down-weighting noisier calibration intervals while up-weighting the more quiescent ones, without taking noise-correlation between visibilities into account. This is because the far-field noise will be dominated by the diagonal component of the covariance matrix (cf. Eq. 3.27). By the same token, minimising calibration artefacts would involve down-weighting measurements with strongly-correlated noise, and up-weighting the less-correlated. This would not, however, minimise the diagonal component: in fact, it will likely exaggerate its up-weighting and down-weighting. As such, it will increase the constant level of the noise-map, but flatten the noise-PSF's contribution. There are thus two competing types of noise that we seek to minimise: uncorrelated noise, which corresponds to δuδv = 0 (i.e. the diagonal components of the gain covariance matrix), and correlated noise, which corresponds to δuδv ≠ 0 (i.e. its wings). Minimising the first will minimise far-field noise without optimally reducing artefacts, while minimising the last will minimise noise near sources at a cost to far-field noise. In the following sub-sections, we will discuss weighting schemes used to accomplish this.

Optimising sensitivity

The Cov-Cov relationship (Eq. 3.25) tells us that, far from any sources, the variance map (Eq. 3.27) is dominated by a constant term: the contribution from thermal noise and the diagonal of the residual visibility covariance matrix. Maximising sensitivity far from sources therefore implies minimising Diag {Cov{γ}}. This is equivalent to assigning visibilities weights inversely proportional to their variance:

w b = 1 Var{γ b } (3.37)
For each baseline, those times with larger variance in the residuals will be downweighted, and those with smaller variance will be up-weighted; this scheme does not require information about the underlying gains, only the error on our solutions. Since we are treating σ 2 n as a constant for all antennas and all times, those times where our gains estimate is closer to the true gains will be up-weighted, and those moments where they are farther from the actual gains will be down-weighted: hence the term "adaptive quality-based weighting". Note that the diagonal of the weighted residuals' covariance matrix should therefore become constant: this weighting scheme explicitly brings the residuals closer to what is expected in the case of perfect calibration, assuming uncorrelated noise. For the remainder of this paper, we will refer to these weights as sensitivity-optimal weighting.

Minimising Calibration Artefacts

Minimising calibration artefacts -i.e. optimising the sensitivity near bright sourcesmeans flattening the noise-map. Since the noise-map can be understood as a noise-PSF convolved with all the modeled sources in the sky modulating the background variance level, it will be flattest when its peak is minimised. From the Cov-Cov relationship (Eq. 3.25), we can see that, at the peak of the noise-PSF (which would be the variance at the pixel where l = m = 0), the Fourier kernel is unity: the variance for that pixel is thus the sum of all the cells in the covariance matrix. By accounting for normalisation, we can write the variance at the centre of the noise-PSF as:

V (w) = w T Cov{γ}w w T 11 T w (3.38)
Our optimality condition is then, after some algebra:

0 = ∂ ∂w (V ) (3.39) ↔ Cov{γ}w = 11 T w w T 11 T w -1 w T Cov{γ}w (3.40)
We find that one w which satisfies the above is:

w = Cov{γ} -1 1 (3.41)
where 1 is a vector of ones. Provided that the hypotheses made so far hold6 , these weights depend only on calibration quality: badly-calibrated cells will include spurious time-correlated signal introduced by trying to fit the noise n on visibilities. Downweighting these cells helps suppress artefacts in the field, at the cost of far-field sensitivity. This weighting scheme is thus only a function of the relative quality of calibration at different times, boosting better-calibrated visibilities and suppressing poorly-calibrated visibilities. For the remainder of this paper, we will refer to these weights as artefactoptimal weighting.

Estimating the Covariance Matrix

In our simulations, we have worked from a known covariance matrix and shown that our predictions for the residual image's behaviour hold. With real data, however, we do not have access to this underlying covariance matrix. Since our weights are extracted from said matrix, estimating it as accurately as possible remains a challenge: this is in turn limited by the number of samples which can be used for each cell.

Each cell in the covariance matrix is built by averaging a number of measurements, or samples. The more samples are available, the better our estimate becomes: once we have more samples than degrees of freedom, we say that our estimation is well-conditioned. Otherwise, it is poorly-conditioned. In this section, we will discuss ways in which we can improve the conditioning of the covariance matrix estimation.

Baseline-based Estimation

One way to improve the conditioning of our covariance matrix estimation is to make the same hypothesis as the calibration algorithm: we can treat the underlying gains as constant within each calibration interval. Provided this interval is known, this allows us to find a single estimate for each interval block of the covariance matrix, turning a n b × n b matrix into a smaller n intervals × n intervals equivalent, where n intervals is the number of solution intervals used for to find the gain solutions. We then improve our conditioning by a factor of n int , which is the number of samples in a calibration interval. The estimate Cov{γ}of the covariance matrix Cov{γ} is built by applying the covariance operator:

Cov{γ} def = Ĉγ = 1 n int i∈n int (γ i -⟨γ⟩) (γ i -⟨γ⟩) H (3.42)
where the ⟨⋯⟩ operator denotes taking the average over the full vector. If the calibration solver's gain estimates are unbiased (i.e. E{ĝ} = g) and the model of the sky is sufficiently complete, this quantity should be zero. Having created Ĉγ , which will be of size n b × n b , its cells can now be averaged over blocks of n int × n int . This allows us to estimate the weights for each baseline and each time.

Mathematically, we retrace the steps of Section 3.2. In the absence of directiondependent effects, we define the residual visibilities as before, and use them to define the normalised residual visibilities ρ b :

r b =w b κ b γb + ǫ (3.43) ρ b = r b k b (3.44)
We then organise the residuals in cells:

r C = ⎛ ⎜ ⎝ ⋮ ρ b∈C ⋮ ⎞ ⎟ ⎠ (3.45) R = ⋯ r C ⋯ (3.46)
R corresponds to a matrix containing all the residual visibilities within one calibration cell C, i.e. for b ∈ C where g C = const. It is therefore of size n intervals × n C , where n intervals is the number of calibration intervals in the observation. Normalising the residual visibilities by k b allows us to recover the underlying covariance matrix by multiplying the residual visibility matrix R with its Hermitian conjugate:

Ĉγ [b ∈ C, b ′ ∈ C ′ ] = R H R [C, C ′ ] (3.47)
Note that we have divided the noise term by the flux model κ b , which can be very small in some cells. As such, care must be taken not to cause the relative thermal noise contribution to explode: those cells where this would occur are dominated by thermal noise, and information on the covariance matrix cannot be recovered from them.

In this framework, we simply treat the index C as containing all the times and frequencies, for individual baselines, corresponding to a single calibration interval. Ĉγ is then an estimate of the residual visibility covariance matrix.

Antenna-based Estimation

In the subsection above, we assumed that finding one solution per interval will give us strong enough constraints to make the problem of estimating the covariance matrix well-conditioned: this may not be true in all cases. Conditioning may then need to be improved further: in this subsection, we show one way in which this can be done. There are others, e.g. using the rank of the matrix itself to find better-conditioned estimates of the covariance matrix at a lower resolution (i.e. a single estimate for a greater number of cells). They will not be presented in this paper, but are a possible avenue future work on this topic.

In estimating the covariance matrix for each baseline and each calibration cell, we are severely limited by the small number of samples in each cell. One way to overcome this problem is to find estimates for the variance of antenna gains, and use these to return to the baseline variances. In this formalism, we extend C -the set of visibilities in a calibration cell -to include all visibilities pointing at a single antenna at a given time. Let us begin by writing an expression for the gain vector, which contains the gains for all antennas and all calibration cells:

ĝc = ⎛ ⎜ ⎝ ⋮ ĝτ∈c p ⋮ ⎞ ⎟ ⎠ (3.48) Ĝ = ⋯ ĝc ⋯ (3.49)
and the variance on each antenna in each calibration cell is then:

=E{ĝ c ĝH c } -E{ĝ c }E{ĝ c } H =g c g H c (3.50)
As we can see, Eq. 3.50 is simply a vector form of Eq. 3.12. The residual gains of Eq. 3.12 can now be understood as random samples of the covariance between the gains for antennas p and q at a given time, assuming complete skymodel subtraction. We can thus define the variance sample matrix as an estimate of the variance matrix :

V C = Var{g c } (3.51) = τ ∈C ĝτ ĝH τ -g τ g H τ (3.52)
We define the residual matrix as:

r τ = d s d K d,τ ĝτ ĝH τ -g τ g H τ K H d,τ + ǫ (3.53)
where we explicitly place ourselves in the limits of our formalism, i.e. that we do not have direction-dependent gains. We now see that at the core of Eq. 3.53 lies Vĝ τ , where ∑ τ ∈C Vĝ τ = V C . The K-matrix is defined as follows:

K d,τ = ⎛ ⎜ ⎝ k d p,τ 0 k d q,τ 0 ⋱ ⎞ ⎟ ⎠ (3.54)
Since the residual matrix depends on the gains, we define the residual visibility vectors as:

r C = ⎛ ⎜ ⎝ ⋮ r τ ∈C ⋮ ⎞ ⎟ ⎠ (3.55) R = ⋯ r c ⋯ (3.56)
r C corresponds to a matrix containing all the residual visibilities within one calibration cell C, i.e. for τ ∈ C where g C = const. Let us define n C as the number of elements in each calibration cell. The residual variance sample matrix can now be built by multiplying the residual visibility matrix with its Hermitian conjugate:

⩔ =R H R (3.57)
Note that we do this because it allows us to turn a single noise realisation ǫ into a statistical quantity σ. We can relate ⩔ to the variance of individual antenna gains:

⩔ = τ ∈C ⎛ ⎝ d,d ′ s d K d,τ V C H K H d,τ s d ′ K d ′ ,τ V C K H d ′ ,τ + Iσ 2 ⎞ ⎠ (3.58)
To reach this point, in Eq. 3.23, we made the hypothesis that the sky brightness distribution is dominated by emission from well-separated point sources. Applying this hypothesis again here, we can make the approximation that the cross-terms in the sum over d, d ′ average to zero: ∑ d,d ′ ≠d ≈ 0. We then have:

⩔ ≈ τ d s 2 d K d,τ V C H V C K H d,τ + Iσ 2 (3.59) = V C 2 ○ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ τ d s 2 d k d,τ k H d,τ def = S ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ + n C Iσ 2 (3.60) V C = S ○-1 (⩔ -n C Iσ 2 ) (3.61)
where ○ denotes the Hadamard or entrywise product and k = Diag {K}. Thus, ⩔ allows us to estimate the variance of each antenna and for each calibration cell by using all the visibilities pointing to that antenna within that calibration cell. With this information, we can then rebuild the baseline-dependent matrix, having improved our sampling by a factor of n ant .

constant covariance matrix if this matrix is zero everywhere outside of the diagonal. In effect, the noise-map becomes flatter, but much broader.

Applying the Correction to Real Data

In this section, we show the effect of adaptive quality-based weighting on real data. The dataset used in this section is a single sub-band from an 8-hour LOFAR observation centred on the Extended Groth Strip (α=14:19:17.84,δ=52:49:26.49). The observation was performed on August 28th, 2014. The subband includes 8 channels of width 24.414 kHz each, for a total bandwidth ranging from 150.2 to 150.4 MHz. The data have been averaged in time to 1 data point per second. The data was calibrated using Wirtinger calibration (see Tasse 2014;Smirnov & Tasse 2015, and references therein) and a sky model consisting only of a nearby calibrator source, 3C295. A reference image (a cutout of which is shown in Fig. 3.7a) was made by calibrating the data according to best practice for LOFAR survey data: 1 calibration solution per 8 seconds and per 4 channels. The residual data was then corrected by the gain solutions and imaged using Briggs weighting (robust=0), pixel size of 1.5 ′′ , and deconvolved using the default devoncolution algorithm in DDFacet (Tasse et al. 2017).

The data was then time-averaged to create a new, 2.4 GB dataset with 1 data point per 8 seconds. Deliberately poor calibration was then performed on this dataset, solving for 1 calibration solution per 2 minutes (caeteris paribus). The resulting corrected residual data was imaged using the same imaging parameters as the reference image, and a cutout of the result is shown in Fig. 3.7b. As expected, the very long calibration intervals introduce calibration artefacts in the image. The brightest sources are still visible, but much of the fainter emission is buried under these artefacts. We are then in a case where our residual visibilities are dominated by calibration error rather than sky model incompleteness.

Weights were then calculated based on the badly-calibrated residual visibilities. Fig. 3.7c was made using the same visibilities as Fig. 3.7b and applying baseline-based, sensitivity-optimal weight. Similarly, Fig. 3.7d used the poorly-calibrated residual visibilities with the application of baseline-based, artefact-optimal weighting. These weights are likely to be poorly-conditioned. In both cases, all other parameters were conserved.

Note that applying antenna-based sensitivity-optimal weighting to the badly-calibrated data (shown in Fig. 3.8b) allows us to recover the reference image with only a very small increase in rms (increased by a factor of 1.14). Further testing on complex field simulations will be required to ascertain the usefulness of artefact-optimal weighting: it is likely that it fails to correct the image fully due to the poor conditioning of the covariance matrix used here.

The pixel histograms show us that the weights do not completely mitigate the poor calibration interval choice, but certainly give a dramatic improvement over the unweighted, poorly-calibrated residuals. This is compatible with our statement that the weights give similar residuals in the image with a dramatic improvement in time at some cost in

Discussion

This chapter began by investigating the use of an algorithm inspired by "lucky imaging" to improve images made using radio interferometric data. By investigating the statistics of residual visibilities, we have made the following findings:

• A relationship between the statistics of residual visibilities and residual pixel values (the "Cov-Cov relationship").

• A description of the noise-map in the image plane as a constant variance level modulated by a noise-PSF convolved with the sources in the field. This gives the variance in the flux of the image as a function of distance from the sources in the sky for a given calibration.

• An Adaptive Quality-Based weighting scheme, which reduces the noise in the image (and the presence of calibration artefacts) by minimising either the constant noise term or the noise-PSF.

While our results are not a panacea for poor calibration, they show that we can not only improve images made with well-calibrated data, but also mitigate the worst effects of poorly-calibrated visibilities in otherwise well-calibrated datasets. Provided that the gain variability timescale is long enough at certain points of the observation, we can effectively get images of similar quality using both the "standard" best-practice calibration interval for LOFAR survey data (calibration solution interval of 8 seconds) and a significantly larger solution interval of 2 minutes (frequency interval unchanged). Of course, if no such stable interval exists, there will be no good intervals to upweigh, and we will be left only with equally-poor data chunks. This means that, in the right conditions, net pipeline time can be sped up by a factor of nearly three (by increasing the size of the calibration solution intervals), at a slight cost in sensitivity. This increase will be greater than what could be achieved with existing comparable methods such as "clipping". However, it does not aim to replace such methods, but rather to act as a complement to them: clipping will be better at flagging a handful of bad visibilities, while this weighting scheme will help push the general visibility statistics (and therefore the image-plane statistics) towards a better distribution.

We emphasize that the adaptive quality-based weighting schemes work because the noise-map describes the underlying noise distribution, of which calibration artefacts are one single realisation. To fully characterise the artefacts, the correlation between different pixels (i.e. off-diagonal elements of Cov{ỹ}) must be computed; this has not been done in this thesis. Nevertheless, lesser constraints on the spatial distribution of artefacts can be found using only the diagonal elements of Cov{ỹ}. The weighting schemes merely seek to minimise this spatial distribution as much as possible: the end result is fewer artefacts, which can be distributed across a much larger area. This is the source of the dramatic improvement from 3.7b to Fig. 3.7c. We have simply down-weighted those visibilities where spurious signal was introduced by the calibration solutions, and up-weighted those visibilities where such signal was lesser. Since this spurious signal is the source of calibration artefacts, downweighting the associated visibilities reduces it dramatically. The poor improvement from Fig. 3.7b to 3.7d is likely due to limits in the conditioning of our estimation of the covariance matrix.

The work presented here can be improved upon, notably by working on improving the conditioning of our covariance matrix estimate: for real observations, it is impossible to have more than one realization of each gain value for all antennas. By treating each visibility within a calibration interval as a realization of the true distribution, we can better estimate the covariance matrix per baseline, and thus reach a better estimate of the variance in the image-plane. Of course, in practice, we can never access to the true, underlying time-covariance matrix for each baseline. Significant hurdles remain:

• The impact of sky model incompleteness (since calibration requires a sky model) is ignored in this paper; we start by assuming that we have a complete sky model.

In practice, of course, acquiring a complete sky model is often a key science goal in and of itself. The impact of this hypothesis therefore ought to be investigated in future work.

• The conditioning of our covariance matrix estimation remains a concern. By using an antenna-based approach, we can improve conditioning by a factor of n ant , but this is only one approach among many. Further work is needed to investigate which method, if any, proves optimal.

Once the conditioning issues of the artefact-optimal weights are resolved, a Briggs-type implementation of the overall adaptive quality-based weighting scheme with a sliding parameter could conceivably be implemented. As things stand, however, there is no good reason to use the currently-existing artefact-optimal weighting scheme over the well-conditioned sensitivity-optimal weighting scheme.

Chapter 4

Imaging the full LOFAR field of view

Aims & Methodology

In this section, we describe the work done to make a wide-field image of the entire LOFAR primary beam, centred on the Extended Groth Strip (EGS). A large-field, image of this field has never yet been made at such low frequencies. Due to technical limitations, the maximum attainable image size (in terms of number of pixels) will limit resolution to about 5 ′′ . As such, we do not use the international stations for this part of the project, as their contribution would be averaged out within a single pixel. They are only used in later parts of the project, when imaging 3C295 and calibrator sources in the field at high resolutions. A more complete sky model can be created from this ∼ 5 ′′ -resolution image, and bright out-of-field sources can be subtracted from high-resolution images of the EGS so that their sidelobe emission do not contaminate the final images. To achieve this, we require direction-dependent calibration. Because we do not use the international LOFAR stations at this point, we do not yet need a high-resolution model of 3C295.

We aim to achieve a sensitivity of a few hundreds of µJy using visibilities recorded between 110 and 180 MHz. This corresponds to a sensitivity of ∼ 20 µJy at 1.4 GHz for typical synchrotron radio sources, which have spectral indices of α ∼ -0.7 (where S ν ∝ ν α ). Reaching this sensitivity over the full LOFAR primary beam requires direction-dependent calibration. Our methodology is therefore to calibrate the full LO-FAR bandwidth using the LOFAR Surveys KSP pipeline, which performs a facet-based direction-dependent calibration. We then perform source analysis on the most interesting objects in the field, using ancillary data from NVSS (NRAO VLA Sky Survey, Condon et al. 1998), SDSS (Sloan Digital Sky Survey, York et al. 2000) and WISE (Wright et al. 2010) surveys and cross-referenced using SIMBAD (Wenger et al. 2000). Overlays for those sources will be shown in Section 4.4 in a "radio-source zoo" intended to show the reader what sources a typical LOFAR extragalactic field contains.

Data Reduction

Data & Observation Properties

The dataset used for this PhD project is part of the LOFAR Surveys KSP, and consists of a number of 8-hour pointings covering as much of the sky visible to LOFAR as possible. We analyse a subset of this dataset, a single one of these pointings, an observation performed on the 28th of August 2014 and centred on the EGS. It is an 8-hour-long observation; there exist another four 8-hour pointings in the same direction, but these data are not reduced in this project. We limit ourselves to analysing only the HBA observation, meaning that we use 365 sub-bands which sample a total bandwidth of 110-182 MHz. We use all Core and Remote LOFAR stations. When imaging 3C295 and other sources in the field with the international stations during Chapter 5 and Section 5.6, we only use those which were online at the time (specifically the German stations 1-5 and 7, along with the Swedish and British stations).

The observation is centred on the EGS, and not on 3C295. This means that our calibrator source is not at phase centre, which can potentially introduce problems. Thankfully, new developments in imaging (Tasse et al. 2017) ensure that the direction-dependent PSF (one of the problems introduced by an off-centre calibrator) is properly modeled and corrected for. This means that direction-independent calibration will be exactly correct in the direction of 3C295, and a good model of this source can be made.

The data was retrieved from the LOFAR Long-Term Archive, and is thus pre-processed and flagged for RFI using the standard tools (NDPPP and AOFlagger; see [START_REF] Van Diepen | DPPP: Default Pre-Processing Pipeline[END_REF]Offringa 2010, respectively). The next step is therefore to calibrate the data.

Reducing the Data

Reducing LOFAR data is a complex business for a number of reasons. When seeking to make high-fidelity, wide-field images, the need to minimise decorrelation (cf. Section 2.12.2) means that data can only be averaged with care. As visibilities are averaged in time and frequency 1 , the relative contribution of signal picked up on longer baselines from sources far from the observation centre is decreased. This translates to a blurring of sources as distance from phase centre increases.

What's more, Jones matrices can vary over short durations or show peculiar spectral behaviour. Even if interferometric data were averaged massively using baselinedependent window functions, information on this underlying gain behaviour would be lost -this means that averaging data to reduce its size (and associated reduction wall time) is a risky proposition.

As such, we seek to reduce an eight-hour LOFAR HBA observation, averaged down to 1 measurement per second and per 24 kHz2 . This corresponds to a total of 7 terabytes of raw data: 2920 frequency measurements each second, for 28800 second, of 4 correlations, for each baseline. The observation mode used was HBA DUAL INNER, in which all Dutch stations operate with 24 of their 48 tiles, and these substations are then correlated separately with the rest of the array. This was chosen because this configuration does not negatively impact the array density of short baselines (allows for better sensitivity to, and thus recovery of, diffuse emission), nor does it introduce the calibration difficulties that non-uniform beams would cause.

The calibration strategy consists of a first round of direction-independent calibration followed by direction-dependent self-calibration. In effect, we use the following RIME:

V pq = G p s E sp K sp B s K H sq E H sq G H q (4.1)
which is the same as Eq. (2.30). Here, K sp = k sp I. K sp therefore commutes with all other Jones matrices, and we can write: K sp K H sq = k s,pq To account for the faceting explicitly, we can write:

V pq = G p Ω E Ωp s∈Ω B s k s,pq E H Ωq G H q (4.2)
where Ω indices denote different facets, and s indices denote individual sources. Because we make the approximation that the E term is constant within a facet, we directly write E Ωp rather than E sp . Each source in a given facet is assumed to have similar directiondependent gains, and their cumulative signal-to-noise allows for better estimations of this direction-dependent gain than would be possible if solving for each source individually. Thus:

V pq = G p Ω s∈Ω E Ωp B s k s,pq E H Ωq G H q (4.3)
We begin with direction-independent calibration to find an estimate of G p ∀p. This first step is done because this part of the Jones chain is the same throughout the sky, and its inverse solution can therefore be applied directly to the visibilities as follows:

V corr pq = Ĝ-1 p V pq ĜH q -1 (4.4) = Ĝ-1 p G p Ω s∈Ω E Ωp B s k Ω,pq E H sq G H q ĜH q -1 (4.5)
In so doing, we assume that E sp = 1 ∀s, p when solving for Ĝ. In general, this is of course invalid. For the calibrator source, however, we have E sp = 1 ∀p by construction.

Because we can choose the calibrator source, the presence of a sufficiently bright source (in our case, an unresolved 3C295) will give extremely good estimations of G. We then have Ĝ-1 p G p ∼ I, and

V corr pq ∼ Ω s∈Ω E Ωp B s k s,pq E H Ωq (4.6)
and direction-dependent calibration then consists of solving for Ê Ωp = Ĝ-1 p G p E Ωp ∼ E Ωp for all facets -this means that, if there are errors in the direction-independent gain solutions, they will bias all direction-dependent gain solutions.

Direction-independent calibration is done with Prefactor (van Weeren et al. 2016), using an initial sky model generated using the VLSSr (VLA Low-frequency Sky Survey -redux, Lane et al. 2012), WENSS (WEsterbork Northern Sky Survey Rengelink et al. 1997) and the NVSS (NRAO VLA Sky Survey, Condon et al. 1998), with 3C295 dominating the model. This model also serves as the initial model for the direction-dependent facet-based self-calibration on the corrected visibilities, performed using the killMS-DDF pipeline (Shimwell et al. 2017). The quality-based weighting scheme discussed in Chapter 3 was used to improve self-calibration convergence and the final image.

Direction-dependent calibration is required, as the field of view is large enough that differential ionospheric errors are introduced as a function of direction: the signal from different sources cross different sections of the ionosphere, and the actual gains for these sources are thus different from the gains for the calibrator (see Eq. (4.6)). They must thus be solved for individually. Assuming that the scale of ionospheric change is of the order of a degree (which is what is observed empirically), then a facet-based approach will give a "good enough" result for most of the field. The direction-dependent calibration consisted of a physics-based approach, solving for XX and YY correlations along with a rotation matrix simulating differential Faraday rotation. These gain solutions were then smoothed in time and frequency to provide a time-independent, frequency-dependent amplitude solution for each station. The full LOFAR bandwidth is large enough to allow for Clock/TEC separation, which consists of explicitly modeling the two directionindependent phenomena which dominate in the gain's spectral structure. The clock difference between stations3 introduces a phase error proportional to ν, while differences in Total Electron Content (TEC) along the line of sight of different stations introduce a phase error inversely proportional to ν. With sufficient spectral coverage, these two effects can therefore be separated and modeled individually, thereby improving the phase solution conditioning.

Direction-dependent self-calibration dramatically improves the noise level in the final image, as shown in Fig. 4.2. Note that not only do artefacts centred on 3C295 (which dominate the noise level in Fig. 4.2a) decrease dramatically, but many of the "diffuse" sources in the field also seem to disappear: this is because their associated artefacts also decrease dramatically, and only physical emission remains.

These sources nevertheless remain in the field, as shown in Fig. 4.3. Note that the noise level is not constant throughout the image: this is because the array is not uniformly sensitive to the entire sky, but has a specific beam response. Correcting for this beam response in the image-plane gives integrated flux rather than apparent flux.

Effect of Quality-based Weighting Scheme

Finally, the quality-based weighting scheme developed in Chapter 3 was applied throughout this procedure. Noise maps made with and without the weighting scheme are shown in Fig. 4.4. A few things are important to note in this comparison. Firstly, the weighting schemes decrease the noise by a factor of ∼ 10: for the noise-map to be visible without weighting schemes, the contrast of the weighted image needed to be pushed way low. For reference, Fig. 4.4b is the same image (albeit without the facets shown explicitly) as Fig. 4.1 -though its levels have been changed to better illustrate the effect of the quality-based weighting scheme. This is why we can't see much at all in that image: it is the contrast between Fig. 4.4b and Fig. 4.4a which is relevant here. The noise distribution is more clearly visible in that figure. Relatedly, by comparing that image and Fig. 4.4a, we note that the unweighted image has a different noise per facet, unlike the weighted image. This is due to the fact that a single visibility will be calibrated better or worse depending on the signal-to-noise of individual facets, constrained by the direction-independent calibration which uses the highest signal-to-noise available. It is no coincidence that the facet containing 3C295 is the only one with a noise even remotely comparable to the weighted image's: it is the facet with the best signal-to-noise. This becomes significant when taking into account the fact that the last 2 hours of this observation are much noisier than the first 6: the quality-based weights will decrease the overall contribution of these two hours of data throughout the image, whereas not using them will lead to different facets having stronger or lesser contributions from these noisy visibilities -and therefore different noise levels. This is also why fringes can clearly be seen in some facets: a single visibility (or a set of a few visibilities) are particularly badlycalibrated, and the associated fringe dominates in the facet. Linking back to Chapter 3, because we have a gain solution per visibility and per facet, we have one noise-PSF per facet; this suggests that the quality-based weighting scheme could likely be extended to a direction-dependent formalism. Table 4.1. Table recapitulating the names, positions and likely associations of all 12 chosen sources in the primary beam. These sources were primarily chosen because they had peculiar or interesting diffuse emission. The associated LOFAR image thumbnails are also given for reference.

For each source, we show four images: a cutout of the LOFAR image around the source, an SDSS cutout of the same area with the LOFAR image superimposed as an overlay, the same SDSS cutout with an NVSS image overlayed, and finally the same SDSS image again with a WISE image overlay. The first image aims to show why the source was selected. The second shows optical counterparts to LOFAR emission. The third serves to check whether we pick up similar structure and emission in NVSS images as the LOFAR wide-field image. Note that the NVSS images are at a much lower resolution, but can be a useful way to ensure that we are not choosing artefacts or spurious emission as sources to investigate. Finally, the WISE overlay attempts to see whether infrared emission can be associated to these sources. Note that the sensitivity of SDSS images are not homogeneous, as not all patches of the sky have received the same coverage.

Images & Cross-matching

see that the emission can be successfully matched to no less than two individual NVSS sources. In other words, if those NVSS sources are linked to the LOFAR source -and given their spatial distribution, they almost certainly are -then they are two components of a single, deeper structure. The infrared image tells us very little.

Note that an optical source lies along the line linking the two diffuse poles of EGS-3: this source is SDSS J142933.44+544335.2, the brightest galaxy in a cluster. If EGS-3 is a radio galaxy with jet emission, which it almost certainly is (the rough symmetry of the two "cotton balls" on either side of a "rod"-like component strongly indicate that this is a very warped jet), then its optical emission must lie somewhere within the "rod". That this emission would be associated with violent interactions between jetaccelerated particles and ambiant inter-cluster and intra-cluster media seems likely: it would explain the complex, turbulent distribution of emission and the strong "warping" of the jet structure into the whirls and eddies seen in Fig. 4.8a. EGS-3 can therefore be fairly unambiguously associated to a galaxy cluster.

Discussion

The work done in this section had the primary goal of creating a decent model of the sky seen by LOFAR's primary beam. To do this, a wide-field, 6 ′′ -resolution image of the fully primary beam was made: this image makes no use of the LOFAR international stations due to technical limitations. Using PyBDSM [START_REF] Mohan | PyBDSF: Python Blob Detection and Source Finder[END_REF], a sky model of the full field of view was extracted from the aforementioned image made after a few rounds of direction-dependent self-calibration.

However, the scientific value of this image is not limited to its use in imaging the Extended Groth Strip using international LOFAR. Much of the image will not be recreated using the international stations, in part due to computational constraints (long run time for such a large task) and in part due to technical limitations (performing direction-dependent calibration using international LOFAR would likely require the use of techniques such as fringe-fitting, which are as-of-yet unimplemented in the DDF-kMS pipeline).

Further work on this image is of course both possible and highly desirable: the work done here was artisanal, hand-picking sources and associating them individually to equivalents in other bands. Modern survey work, however, requires automation. Applying Bayesian fitting methods to find most likely counterparts to all the sources extracted in the field could yield useful information, as we would then have access to population statistics. However, doing this would require running yet more self-calibration to remove remaining artefacts, whereas our primary aim is to produce a high-resolution map of the Extended Groth Strip, which has the richest multi-frequency coverage and is therefore the most likely to give effective matches to all the sources in the field. We therefore do not proceed to do this work on the low-resolution, wide-field image of the Extended Groth Strip and its neighbourhood, but are likely to return to the data to do this at a later date.

Aims & Methodology

Our aim in this section is to create a high-resolution (0.1 ′′ ) model of 3C295, which did not exist at LOFAR HBA frequencies when this project began. This means creating a model that will allow us to find good phase-calibration solutions for LOFAR international baselines, and thus use the international stations to create a high-resolution image of the EGS. Because we only solve for gains in the direction of 3C295, we must then estimate the impact of directional gain errors over the rest of the field.

One major constraint in this project is that 3C295 is fully resolved at these resolutions, which makes deconvolution difficult. As such, we consistently use a uv-cut of 100km throughout our imaging in this section, to limit the impact of diffuse emission (which is very hard to deconvolve properly -poor deconvolution introduces artefacts in the restored image) on our images. This choice of uv-cut is a compromise between removing diffuse emission in the images and providing sufficiently good conditioning to the calibration inverse problem1 . Provided that the phases of short-baseline gains do not evolve over the self-calibration procedure, then the information on spatial brightness distribution associated to these spatial fringes is not lost.

Another constraint is that the model needs to be applicable for the full LOFAR bandwidth. We therefore select 6 sub-bands out of the total LOFAR HBA bandwidth, evenly spread throughout the bandwidth as shown in Fig. 5.1. This approach allows us to benefit from the improved conditioning of high-frequency data (since the diffraction limit on angular resolution is proportional to frequency) throughout the LOFAR bandwidth. By starting calibration with a model extracted from a VLA observation of 3C295 with a similar angular resolution, self-calibrating over the initial 6 subbands, then self-calibrating once with 60 subbands (centred on the initial 6) and then the full bandwidth, we expect to achieve good signal-to-noise. Then, we bootstrap the model such that the spectral index of its integrated flux is compatible with Scaife & Heald (2012).

The resulting model is expected to be reliable enough to calibrate our Groth Strip data for the full LOFAR array and the entire HBA bandwidth. It would theoretically be possible to use the direction-dependent gain solutions from Chapter 4 to improve our images, but these are only available for core and remote stations, and not international stations -direction-dependent effects are therefore expected to still arise in the final images due to this. As such, we focus on direction-independent calibration from this point on. Quantifying the impact of the resulting directional gain errors is the subject of Section 5.6.

Data Reduction

Our RIME

As outlined in Section 2.7, calibration is the process of solving for Jones matrices, which are associated to individual antennas. Here, the only source in our model (the source that dominates, by far, over other sources in the field) is 3C295; we can thus limit ourselves to direction-independent calibration, which is equivalent to performing direction-dependent calibration in the direction of 3C295 only. The associated RIME can be written as:

V ν pq = G ν p s=3C295 E ν sp K ν sp B K ν sq H E ν sq H G ν q H (5.1) = G ν p E ν 3C295,p B ν 3C295 k s,pq E ν 3C295,q G ν q H (5.2) = J ν 3C295,p B ν 3C295 k s,pq J ν 3C295,q H (5.3)
where we mute the time-dependence of all terms, and explicitly model the LOFAR beam during both calibration and imaging, so as to ensure that the true sensitivity of each baseline to this out-of-phase-centre calibrator is taken into account appropriately.

calibration (after the successive self-calibration passes). Similarly, there is little apparent change in phase structure (note that the phase wraps from -π to π). It is important to note, however, that there is a difference in the level of the amplitude curve: in the first pass of calibration, the average amplitude value is generally higher than in the last pass of calibration. This is expected behaviour, as there is less diffuse emission in the uv-cut imaging model of 3C295 than in the initial calibration model. This is equivalent to applying a choice of weighting scheme which down-weighs diffuse emission (i.e. a form of uniform weighting).

As for the international station gains, they are shown below for three of the seven international stations used for our observation5 , in Figs. 5.8 to 5.10. We show only one of the four German stations. As we can see, while structure very clearly remains in some of these amplitude curves (notable in the last two hours, e.g. Fig. 5.8b), some of the structure disappears from these amplitude curves. It should be noted that the last two hours of the observation are generally much noisier than the first 6; the last quarter of these plots show this. Quality-based weighting schemes were found to compensate for this effect quite well, and minimise the associated image deterioration6 . Note that the phase wraps extremely fast in the international station gain plots: this is because phase is calculated relative to CS001HBA0, the first antenna in our array. Because the international stations are much farther away from this core station than other core or remote stations are, we expect the associated phase solutions to wrap faster.

Spectral behaviour of gain solutions

So far, we have shown the gain solutions for one subband, out of a total of 324 used subbands. The spectral behaviour of the amplitude and phase solutions can give us insight into whether we are absorbing true flux or calibration artefacts into our model: if our calibration model becomes more accurate to the true underlying flux, we expect e.g. the phase of core station calibration solutions to remain the same (as we are improving high-resolution details of our model) while our international station calibration solutions would become flatter in amplitude. Similarly, we expect to see large-scale frequency-dependent behaviour in the phase solutions of international stations, with smaller fluctuations around this behaviour: the effect of differential TEC, ionospheric turbulence, etc. Indeed, while the core stations (and, to a lesser extent, the remote stations) can be approximated as having the "same ionosphere" above them, this is not the case for international stations: this is a major source of the difficulties in calibrating these stations. As such, for international stations, we expect to see complex structure in the phase gain solutions, characteristic of the differences between the local ionosphere above the core and the international stations. For core stations, we only expect to see structure characteristic of scintillation (regular stripes with no frequency-dependent structure) in the phase gain solutions7 .

Discussion

The plots of pre-self-calibration gain amplitude and phases as a function of frequency are a strong indicator that the gain structure due to unmodeled flux is dominated by a strong, bright source (ring-like structure in e.g. Fig. 5.16). However, the remaining presence of diffuse ring-like structure, even after the self-calibration process, for the international station gain plots is also a strong indicator that we do not recover all of the true sky brightness distribution in the model, as some of the flux is still clearly absorbed into the gains. The most straightforward way to improve on this would be to repeat the procedure while including a full low-resolution sky model: while 3C295 is by far the brightest source in the field at a recorded ∼86.9 Jy in the catalog of the sources in our wide-field image (Table 4.1), it is far from the only contributor (the net flux in the field, according to that same catalog, is 531.1 Jy). Since calibration artefacts are introduced by both incomplete sky model and poor calibration inverse problem conditioning, improving one of these two conditions would likely result in a net improvement of the model.

Because of time/computational constraints, we used calibration solution intervals of 8 channels and 1 minute in time and frequency respectively. Given the brightness of 3C295, this was clearly not justified on the basis of signal-to-noise; it did, however, mean that a single pass of self-calibration on one subband took 10 minutes rather than the hour and a half with solution intervals of 4 channels and 8 seconds in time and frequency respectively (factor of 9 improvement in time). With the advent of new tools like cubical (Kenyon et al. 2018a,b), a much faster direction-independent implementation of the algorithm used for killMS8 which was in development throughout my PhD, it is likely possible to dramatically speed the calibration part of the self-calibration loop (the imaging time constraint would dominate as more subbands are added). Because of the use of quality-based weighting schemes, it is unlikely that a dramatic improvement would be seen when self-calibrating over a few subbands, as the relative contribution of poorly-calibrated visibilities would be down-weighted by the quality-based weighting scheme: as more and more data is added, however, the relative contribution of the artefacts introduced by the choice of calibration intervals would likely begin to overcome the purely thermal noise in the images, which is expected to decrease as more data is added.

All the same, this project served as a useful exercise in both developing a good working methodology and in testing calibration strategies for international LOFAR datasets. A lot of the work which went into this project led to dead ends, in the sense that they gave no directly useful results -but they still proved to be an extremely helpful context in which to understand how to improve the LOFAR pipeline. In particular, gain smoothing 9 has shown promising result on LOTSS images, but was not used in this project because of the full-Jones solutions: properly implementing this smoothing in a framework capable of accounting for polarisation could provide the kind of regularisation that our gain solutions need.

TESTING DIRECTIONAL GAIN ERRORS -AIMS & METHODOLOGY 129

Testing Directional Gain Errors -Aims & Methodology

Having made both a wide-field image of the EGS without international stations and a high-resolution model of the brightest source in the LOFAR primary beam, 3C295, in order to calibrate the international LOFAR station, we can now begin to proceed to imaging the EGS with the international LOFAR stations. Before we do so, however, we must quantify whether the direction-independent calibration on the international baselines is sufficient when imaging the EGS. We expect to see some direction-dependent effects arising in the neighbourhood of 3C295, and also expect that the SNR in the field will be insufficient to perform direction-dependent calibration on the international stations: are these direction-dependent effect so severe that it is not worth imaging the EGS?

In other words, we seek to estimate whether, in practice, the EGS can, in part or in full, be placed on the same facet as 3C295 -and with what loss in accuracy to directional gain errors. The first direction-dependent effect we must take into account is decorrelation, also known as smearing (see Section 2.12). This effect causes signal loss on the longest baselines for sources far from the phase centre. To do this, we make use of the LOFAR Long-Baseline Calibrator Survey (LBCS, Jackson et al. 2016). This is a survey of suitable compact sources to use as calibrators for international LOFAR baselines: these sources are expected to consist mainly of self-absorbed synchrotron radio cores. The survey's source selection procedure is as follows: it selects all sources identified as single sources in WENSS (Rengelink et al. 1997). Then, based on the work of Moldón et al. (2015), which found that the main predictors for compact structure were a high total WENSS flux density and a flat low-frequency spectrum (in line with the expectation above), a selection was made on the basis of spectral index and WENSS flux density for all these sources.

# RA [hms]

Dec A significant preliminary test is therefore to check that we do in fact see any of the LBCS sources. However, this is not the end of the problem: decorrelation causes flux from a source to be "smeared" (hence its other name) around the source. The DDFacet suite can account for this effect by modelling the direction-dependent PSF, i.e. a PSF as a function of direction is called directional gain errors, as it is caused by directiondependent effects.

Assuming that direction-dependent effects (beam, ionosphere) are relatively smooth and dependent only on distance from 3C295, we attempt to characterise the impact of directional gain errors and their associated image-plane DDEs on various calibrator sources around 3C295 and the EGS.

We select 8 LBCS sources, which are all the "good" sources within 5 degrees of the centre of the EGS. Fig. 5.20 shows the position of these calibrator sources with respect to the EGS. If we compare Fig. 5.20 and Fig. 4.5, which are the same image at different levels of contrast, we can immediately see that the EGS lies in two separate facets of the direction-dependent calibration we performed without the international baselines: these are the eastern neighours of the 3C295 facet. Neither of these EGS facets contains a LBCS source (indeed, if they did, our work would have perhaps been easier); however, nearby facets do contain one or multiple LBCS calibrators. If the effect of directional gain errors for these sources is comparable to those for the EGS facets, then we can estimate their impact on an eventual high-resolution image of the EGS.

Testing Directional Gain Errors -LBCS Sources in the Primary Beam

The LBCS sources were chosen automatically using a routine written by Leah Morabitoh and shared in a personal communication. This was done in early 2017. Only the brightest reliable sources were kept, provided they had a good reliability factor (the LBCS catalog includes a reliability factor for each source; we have not verified that these reliability factors have not changed since). We thus had 8 LBCS sources around the field, described in Table 5.2. Their position on the sky is shown in Fig. 5.20. They will henceforth be referred to as LBCS1, LBCS2, and so on -note that LBCS6 is actually 3C295.

For each LBCS source, we create an image made using the same data and imaging parameters as our final 3C295 self-calibration pass. This means that we simply use our final restored image of 3C295 for LBCS6, but make a new image (with the same parameters except the direction) for each other LBCS source. We then show overlays of the results onto our wide-field image of the full primary beam. Note that both the international LOFAR images (∼ .4 ′′ resolution) and the widefield images (∼ 20 ′′ ) are oversampled: their pixel sizes are 0.1 ′′ /pixel and 5 ′′ /pixel respectively.

Because we want to show both the impact of DDEs and the positions of the sources, we show two different overlays for each source. One starts at 5σ and the other starts at 15σ; both increase exponentially up to the maximum value in the image. The first thus gives more information on artefact structure around the sources, while the second gives more information on source structure itself. The σ-level cutoffs were chosen arbitrarily. Note that the σ is calculated using random pixels in the image: it is therefore not the "local RMS" used in e.g. PyBDSM extraction (which will be higher, as it includes

The most encouraging result shown above, however, is for LBCS4: this source appear to be only very slightly affected by directional gain errors. Although LBCS3 seems affected by directional gain errors, it is not so affected that it is smeared out completely: we can thus expect there to still be signal at that distance from the calibrator and phase centre. Since the EGS lies between these two calibrator sources and 3C295, this is extremely positive. Similarly, the fact that LBCS4 resolves into two separate sources (as does LBCS1, for that matter) would indicate that even neighbouring bright sources do resolve properly -this is not a rigorous test, but it is a fortunate serendipitous observation to be able to make. We conclude this thesis manuscript with an overview of the results achieved over the course of this doctorate and a discussion of individual items in this larger context. We follow the structure of the manuscript. These conclusions & discussions are intended to supplement the conclusions & discussions at the end of each chapter with a synthetic overview of the results of this doctoral thesis.

The projects undergone as part of this doctoral degree have strong opportunities for future work, both in terms of the algorithmic work & its application. These will be outlined in detail in the relevant section. It is worth immediately pointing out that only preliminary results were obtained for the application to data by the end of the doctorate's three-year period. These preliminary results are shown in this manuscript, and while they are not of science-grade quality (no astrometric corrections, proper flux bootstrapping, strategy likely needs revisiting for best results) they are nevertheless useful markers of the wealth of information available in the EGS.

Algorithmic Work

One of the aims of this PhD was to develop new interferometric tools & techniques and to apply them to LOFAR data. This was arguably accomplished successfully:

• We found a relationship between the statistics of residual visibilities and residual image-plane pixel values (the "Cov-Cov relationship").

• From this relationship, we see that calibration artefacts are directly caused by correlated calibration gain errors. In the absence of correlated calibration errors, no artefacts are present.

• This allowed for a description of the noise-map in the image plane as a constant variance level modulated by a noise-PSF convolved with the sources in the field.

• This led to the development of an adaptive, quality-based weighting scheme, which reduces the noise in the image (and the presence of calibration artefacts) by minimising either the constant noise term or the noise-PSF.

• One version of this adaptive, quality-based weighting scheme was successfully applied throughout our data reduction & wide-field imaging. Other versions are possible, but at the time of writing suffer from issues due to poor conditioning: the development & deployment of robust versions of these will be the subject of future work.

The Cov-Cov Relationship

The theoretical framework of the cov-cov relationship, which links the variance and covariance in the visibilities to the variance and covariance between pixels in images made with these visibilities, was shown to hold on real data by improving images made with both direction-independent and direction-dependent calibration. This relationship tells us that the variance in images made from interferometric data is the result of two contributions:

Cov{ỹ} = d ⎛ ⎝ b φ d bb [Cov{γ}] bb + w 2 b σ 2 φ bb C b + b,b ′ ≠b φ d bb ′ [Cov{γ}] bb ′ F bb ′ ⎞ ⎠ (6.1)
which is a repeat of Eq. (3.25). The first term above, which corresponds to the sum over b, can be thought of as "thermal" (i.e. uncorrelated) noise in the visibilities, giving rise to a constant variance throughout the image. The second, which corresponds to sky brightness distribution absorbed into the gain solutions -i.e. correlated noise in the calibration residual visibilities -gives rise to what we call a noise-PSF: a distinctive shape which is convolved with every source in the field (in the case where the true gains are direction-independent), and which gives the distribution function of which calibration artefacts are a single realisation. It corresponds to the sum over b, b ′ ≠ b above -i.e. to those cells of the visibility covariance matrix which are off-diagonal, or correlations between residual visibilities.

In Eq. (6.1), Cov{ỹ} is the covariance matrix of the residual image pixel values, and Cov{γ} is the covariance matrix of the residual visibilities, which is of size N b × N b .

[Cov{γ}] bb ′ is a scalar quantity equivalent to the value of the covariance matrix at coordinates b, b ′ . φ d bb ′ corresponds to the product of the weighted 1 visibilities, associated to a source in direction d 2 , as seen by baselines b and b ′ . φ d bb is therefore the square of the magnitude 3 of this quantity for a single baseline and direction. σ 2 is the thermal noise in the visibilities, which is assumed to be constant for all of them. Finally, C b is the image-space fringe corresponding to a single baseline, while F bb ′ is slightly too complex to summarise quickly, but its diagonal maps a single visibility covariance value to a δuδv fringe, thus determing the artefact distribution around sources in the field.

Taking the diagonal the matrices shown above gives the "var-var" relationship, which gives the variance map of the image as a function of the variance in the visibilities. If the visibility covariance matrix is diagonal (i.e. we have no correlated noise in the visibilities), then the second term will be null, and the variance of each pixel in the image-plane will be a constant value (since Diag {C b } = I∀b). This is an extremely strong prediction, but one that was validated through simulations: calibration artefacts in interferometric images are due only to correlated noise in the visibilities. In the absence of such correlated noise, there will be no calibration artefacts in an image.

Note that this relationship also confirms, among other things, that "thermal" noise in interferometric images is actually correlated between pixels, and that this correlation is given exactly by the PSF. This is an expected result, and can only be corrected by fully-deconvolving the PSF from an interferometric image 4 . If this is done, then an image made using visibilities containing only uncorrelated noise will itself have noise uncorrelated between pixels, but will not otherwise.

Mapping the variance of images made with interferometric data

This relationship was strongly tested through simulations, verifying the noise-PSF behaviour of calibration artefacts in Fig. 3.5, reproduced in Fig. 6.1. The visibilities for three point sources were simulated for a single uv-track and frequency, and these were multiplied by multiple realisations of time-correlated residual gains. By taking the DFT of each realisation to create dirty maps of this simulated residual map, and calculating the variance for each pixel across our realisations to find the variance map, we see that we do in fact see a PSF-like behaviour in the noise-map, in agreement with our predictions.

1 e.g. Briggs-weighted, quality-weighted, etc 2 i.e. the visibility which the interferometer would recore if this was the only source in the sky. 3 i.e. the product of the visibility and its own complex conjugate. 4 Here, when talking about "fully deconvolving the PSF from an interferometric image", we mean to systematically deconvolve the PSF from each pixel in the image simultaneously. This would be an impractically computationally expensive task, and is therefore not a realistic proposal. Regardless, the point being made here is that the noise of images made using visibilities with purely thermal noise is correlated between pixels -even though the noise is not correlated between visibilities.

calibration with a strong imaging uv-cut;

• Preliminary investigation of the impact of directional gain errors from directionindependent calibration using the model of 3C295 obtained from the step above;

• Creating a high-resolution image of the EGS after subtracting all sources from the widefield image & the high-resolution model of 3C295.

Of these four points, only the first three were successfully carried out. We therefore limit our discussion to the first three points, treating the last as the subject for future work.

Imaging the LOFAR primary beam

We begin by creating a decent model of the sky seen by LOFAR's primary beam. To do this, a wide-field, 6 ′′ -resolution image of the fully primary beam was made: this image makes no use of the LOFAR international stations due to technical limitations (memory use). Using PyBDSF [START_REF] Mohan | PyBDSF: Python Blob Detection and Source Finder[END_REF], this model was extracted from an image made after a few rounds of direction-dependent self-calibration. This directiondependent self-calibration was done using the killMS/DDF pipeline, which fully models direction-dependent errors due to decorrelation: its application to data is one of the novel contributions of this work. This model is still in a very preliminary stage: astrometric errors have not been corrected for, nor has integrated flux bootstrapping been performed. As such, this image (and the model extracted from it) is not of science-grade quality. Even so, and preliminary as it is, it is sensitive and wide-field enough to pick up a good number of interesting diffuse sources. 12 of these were selected to form a so-called "primary beam galactic zoo" of sources, with the aim of showing the reader the contents of a typical LOFAR deep extragalactic field. They were then matched to NVSS, NVSS & WISE images to find counterparts: the results are shown in Table 6.1, which is a repeat of Table 4.1.

The true interest of modeling these sources at the wide-field resolution achieved is that will allow us to subtract (most of) the emission associated with them when imaging the EGS at high resolution. This, in turns, will allow us to image the EGS itself without pollution from the sidelobes of nearby bright sources. Note that there is nothing inherently special about this field compared to LOTSS fields (Shimwell et al. 2017): its interest is in what it can allow us to achieve with the EGS.

Imaging the EGS with LOFAR-VLBI

Using 3C295 to self-calibrate the international stations Over the course of this doctoral thesis, I was able to develop an extensive familiarity with the problem of interferometry. Specifically, the work done on self-calibrating 3C295 Radio galaxy, z = .292 shows that, while we seem to converge towards the same model regardless of initial conditions (2 point sources, full VLA model, 1 point source5 ), this apparent convergence is very slow. As discussed in Section 2.3, the problem of imaging is convex (if not necessarily well-conditioned) but the problem of calibration is not fully convex: this introduces the possibility of a net non-convexity in the problem of interferometry. Regardless of whether my experience was that of non-convexity or of poor conditioning leading to impractically slow convergence time, the solution is to regularise.

One obvious avenue of inquiry likely to yield results would therefore be in the application of amplitude & phase smoothing on our gain solutions: this corresponds to a regularisation of the calibration inverse problem, improving its conditioning by applying a prior which encourages smooth gain solutions. This prior is physically-motivated, in the sense that we expect the true, underlying gain solutions to be smooth, and so jagged structure in the gains can generally be assimilated to source flux being absorbed into the gain solutions. The tools to perform this smoothing are already developed and being tested on LOTTS DR2 by Cyril Tasse, Tim Shimwell, and Martin Hardcastle. They were not applied here because we used full-Jones calibration solutions (which were found to give better gain solutions with our available SNR early on), whereas the DR2 pipeline solves for scalar gain solutions. Once the smoothing framework is generalised to polarisation, however, its application could result in dramatic improvements in international station gain solutions. This work is complementary to that of the LOFAR long-baseline team, and relies heavily on their research (e.g. Jackson et al. 2016). In essence, we rely on the fact that 3C295 is a bright, dominant source in the EGS primary beam: assuming that the overall contribution from fainter sources in the field (which peak at ∼ 4Jy, to 3C295's ∼ 99Jy at that frequency) is negligible, then self-calibration on 3C295 ought to yield both an improved model of that source and good gain solutions for the international LOFAR stations. One point of originality of this work is that we consistently take into account direction-dependent errors due to decorrelation (but not DDEs) -this has only become possible recently, with the appropriate modeling of decorrelation both in imaging (through the direction-dependent PSF) and calibration (by modeling the decorrelation of model sources properly when creating model visibilities). This removes one of the sources of bias in solving the overlal problem of interferometry, and allows calibrator sources outside phase centre to be used reliably.

As it turns out, this was a bold assumption to make: the model resulting from our self-calibration is still rife with calibration artefacts. Based on our understanding of the cov-cov relationship, many of these can be understood as being due to sky brightness being absorbed into the gain solutions: gain inspection supports this theory, showing unexpected and correlated structure in our international station gain amplitude solutions.

Other avenues remain open to go beyond this issue: using the full low-resolution primary beam sky model from Chapter 4 could be one way to improve our images. However, we would still find ourselves with much unmodeled emission on the baselines including an international station. Another possible avenue of improvement would be through the use the EHT-imager tools: these rely on fitting model data to closure quantities (amplitude & phase), which are unaffected by antenna-based errors. This is currently being tested by the long baseline working group.

All these avenues of research are, of course, only possible because we find ourselves in very specific circumstances: the presence of a very bright resolved source in a famous (and therefore relatively well-known) extragalactic field. As such, we have the freedom to explore methods which supplement those of the LBCS team, who are more strongly constrained in that their aim is to provide a solid, reliable catalog of sources suitable for calibrating the international LOFAR stations. This, in turns, requires a solid, reliable strategy, in ways that our relatively straightforward application to a single field does not.

Quantifying the impact of international station directional gain errors

This work is strongly linked to the self-calibration on 3C295: depending on the quality of both the model and the gain solutions attained (i.e. of how the quality of our solution to the overall inverse problem of radio interferometry), directional gain errors will be more or less present across the primary beam.

We choose 8 nearby sources from the LBCS, a catalogue of compact radio sources in the Northern sky, spread around both the EGS and 3C295. These are reprinted below, in copies of Table 5. 2 and Fig. 5.20. Two direction-dependent effects are expected to arise, one as a function of distance from phase centre (i.e. from the EGS) and the other as a function of distance from the direction-independent calibration source (i.e. from 3C295).

The direction-dependent effect which is function of distance from the EGS, decorrelation, is modeled in the DDF/kMS pipeline. As such, we do not expect this effect to dominate. The other direction-dependent effect, function of distance from 3C295, however, is not modeled: this effect corresponds to directional gain errors. Directionindependent calibration is only exact in the direction of the calibrator source when a single source is present in the field of view; this is because, in reality, the astrophysical signal from each source in the field of view reaches us through a slightly different ionosphere, are in areas of the antenna beam which are of varying sensitivity, etc. Some of these models can be modeled and corrected to some degree (e.g. LOFAR beam), while others cannot at the time of writing (ionospheric effects).

With our final model of 3C295 and gain solutions, flawed as they are, we find that there is little direction-dependent variation in calibration artefact intensity for the various LBCS sources in the field. This is a strong encouragement to continue the project and move on to imaging the EGS using international LOFAR stations.

Future Work

As mentioned in Section 6.1.3, there is much work ahead of us with regards to the algorithmic work developed in this doctoral thesis. Not only can the theoretical framework be expanded to include more complex cases, but those results which it already provides are limited by the conditioning of our estimation of the visibility covariance matrix. There is therefore room for both expanding on and improving existing results.

As for the application to data, it is very encouraging to see that not only are multiple LBCS sources resolved to show physically-sensible morphologies (FRI, FRII, headtail radio galaxies, etc), but in some instances are even resolved to show not just one, but multiple sources with physically-sensible morphology. Many of the LBCS candidates appear to be double or even triple sources. This would indicate that the directional gain errors for many of the LBCS sources are not so strong that they deteriorate the image beyond the recovery of such morphological information; since calibration artifacts, DDEs included, are more present near bright sources, the fact that they are not strong enough to pollute LBCS sources is a promising indicator that even fainter sources in the EGS have recoverable morphological information, given our data and calibration strategy. However, we were not able to create images of this field before the end of this doctoral degree: creating this image (and performing the associated quality-control and source extraction) is therefore very high indeed on the priority for future work. The preliminary results shown in this manuscript are not of science-grade quality (no astrometric corrections, proper flux bootstrapping, strategy likely needs revisiting for best results), but they are nevertheless useful markers of the wealth of information available in the EGS as seen with the full LOFAR array.

Chloe, Prudence, Elise, Sacha, and many more besides; I am every day thankful to count myself your friend. I am also extremely grateful to those teachers who first helped me develop an interest in science and mathematics: without Nasreddine Nouar, Katrina Kemper, and Judith Byrd, I would not be what I am today. I hope to be as good a teacher to my own students as you were to me. Similarly, I would like to thank those who encouraged me to pursue this interest in higher education: Catherine Heymans, Andreas Zech and Jacques le Bourlot. Yours is an example I strive to emulate.

Finally, to close this manuscript on a light note, I want to extend my heartfelt gratitude to the cats who have graced my life, at home and abroad. May they purr blissfully for many years to come.

Gr âce à une analyse statistique de l'Equation de la Mesure des Interf érom ètres Radio, un sch éma de pond ération adaptatif est d ériv é, bas é sur la qualit é de calibration des donn ées d'un instrument interf érom ètrique.

Ce sch éma est utilis é sur une observation d'un champ extragalactique, l'Extended Groth Strip, observation qui contient une source radio-vive (3C295) dans son champ de vue. Cette source est r ésolue avec LOFAR-VLBI; un mod èle de source est cr é é afin de calibrer les stations LOFAR internationales. Cela permettra d'imager le champ à une r ésolution comparable à celle du Hubble Space Telescope, dont des donn ées sont disponibles pour ce champ extragalactique.

Mots Cl és

Interf érom étrie, radio, LOFAR, algorithmie, extragalactique Abstract By performing a statistical analysis of the Radio Interferometer's Measurement Equation, we derive adaptive quality-based weighting schemes. These are deployed on an observation of the Extended Groth Strip, which includes a bright 3C source in the field of view. This source, which is resolved for LOFAR-VLBI, is modeled and used as a calibrator source for the Extended Groth Strip. This will allow the field to be imaged with a resolution matching the Hubble Space Telescope's, of which data are available for this field.

  Interest . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.2 The Low-ν Sky: Emission Mechanisms . . . . . . . . . . . . . . . . . . 8 1.2.1 Synchrotron Radiation . . . . . . . . . . . . . . . . . . . . . . . . 9 1.2.2 Free-Free Radiation . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.3 LOFAR: The LOw-Frequency Array . . . . . . . . . . . . . . . . . . . . 13 1.4 This Thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 1.4.1 Algorithmic Developments . . . . . . . . . . . . . . . . . . . . . . 17 1.4.2 Imaging the Full Primary Beam . . . . . . . . . . . . . . . . . . 18 1.4.3 3C295 and the Extended Groth Strip . . . . . . . . . . . . . . . 19

Figure 1 . 1 .

 11 Figure 1.1. MeerKAT image of the galactic centre, made by SARAO. This is a low-resolution version of the true image.

Figure 1 . 2 .

 12 Figure 1.2. Radio optical montage of 3C31, by Alan Bridle: https: // www. cv. nrao. edu/ ~abridle/ 3c31xopt. htm .

Figure 1 . 3 .

 13 Figure 1.3. Radio and far-infrared spectrum for galaxy M82, as estimated by the NRAO online course (https: // www. cv. nrao. edu/ course/ astr534/ FreeFreeEmission. html ).The flat curve corresponds to free-free emission, while synchrotron radiation (negative slope) and thermal dust emission (postive slope) dominate at low and high frequencies respectively.

  (d) International LOFAR. Newer stations (one in Ireland, three in Poland) are not shown here.

Figure 1 . 7 .

 17 Figure 1.7. Geographic location and distribution of LOFAR stations, explicitly showing what is meant by Superterp, core, remote and international stations. All images from the official ASTRON website: https: // www. astron. nl/ radio-observatory/ astronomers/ users/ technical-information/ lofar-array-configuration/ lofar-array-conf

  Figure 1.9. HBA tile layout and image. Note that both cases show a 24-tile HBA array: each core station includes two such arrays. International stations include a single larger tile array.

Figure 1 .

 1 Figure 1.10. Picture of a single LBA antenna. Source: LOFAR technology website.

( a )

 a Arecibo telescope, in Puerto Rico. (b) FAST telescope, in China

Figure 2 . 1 .

 21 Figure 2.1. Examples of large single-dish radio telescopes.

  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Building the Noise Map . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.1 The Cov-Cov Relationship in the δuδv plane . . . . . . . . . . . 3.2.2 Statistical Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.3 Noise Map Simulations . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Adaptive Quality-based Weighting Schemes . . . . . . . . . . . . . . . 3.3.1 Optimising sensitivity . . . . . . . . . . . . . . . . . . . . . . . . . 3.3.2 Minimising Calibration Artefacts . . . . . . . . . . . . . . . . . . 3.4 Estimating the Covariance Matrix . . . . . . . . . . . . . . . . . . . . . 3.4.1 Baseline-based Estimation . . . . . . . . . . . . . . . . . . . . . . 3.4.2 Antenna-based Estimation . . . . . . . . . . . . . . . . . . . . . . 3.5 Applying the Correction to Simulated Data . . . . . . . . . . . . . . . 3.6 Applying the Correction to Real Data . . . . . . . . . . . . . . . . . . . 3.7 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

  20) Note that the only direction-dependent terms in the above are s d and κ d b , which are both inside φ bb ′ (for both b = b ′ and b ≠ b ′ ). By making the approximation that the Fourier kernels of different sources are orthogonal (i.e. that κ d b κ d ′ b ′ = (κ d b ) 2 δ bb ′ ) 5 we can write:

  Contents 4.1 Aims & Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.2 Data Reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 4.2.1 Data & Observation Properties . . . . . . . . . . . . . . . . . . . 80 4.2.2 Reducing the Data . . . . . . . . . . . . . . . . . . . . . . . . . . 80 4.3 Effect of Quality-based Weighting Scheme . . . . . . . . . . . . . . . . 84 4.4 Zoo of Radio Sources . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 4.4.1 Selected Sources . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 4.4.2 Images & Cross-matching . . . . . . . . . . . . . . . . . . . . . . 89 4.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103

Fig

  & Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 5.2 Data Reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.2.1 Our RIME . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.2.2 Calibration Strategy . . . . . . . . . . . . . . . . . . . . . . . . . 107 5.2.3 Testing the impact of initial model choice . . . . . . . . . . . . . 108 5.2.4 Frequency-dependent Scale Error . . . . . . . . . . . . . . . . . . 110 5.2.5 Examining gain solutions . . . . . . . . . . . . . . . . . . . . . . . 112 5.3 Overlays & Spectral Analysis . . . . . . . . . . . . . . . . . . . . . . . . 124 5.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 5.5 Testing Directional Gain Errors -Aims & Methodology . . . . . . . . 129 5.6 Testing Directional Gain Errors -LBCS Sources in the Primary Beam 132

Fig

  The Jones matrices (J d ...,tν ) contain the antenna gain information in matrix form, while K d ...,tν is the Fourier kernel. Let us limit Fourier kernel for direction d, antenna p and one (t, ν) pair. Size 2 × 2. Baseline selection matrix, which picks out 1 visibility out of the full set. Size n pix × n b C b n pix × n pix convolution kernel that defines the PSF. F bb ′ Convolution matrix mapping one δuδv to δlδm. The set of all F bb ′ determines the noise-PSF. Size n pix × n pix .

	Scalars
	n pix	Total number of pixels in image-plane & number of cells in uv-grid
	n ant	Total number of antennas in the array
	n b	Number of visibilities
	b τ	Index for a single visibility. Equivalent to (pq, tν) Equivalent to (t, ν)
	Vectors
	ỹ	Residual image vector, size n pix
	ǫ γ 1	Vector of ǫ, size n pix Contains gain products, size n b . See Eq. 3.12 Vector containing 1 in every cell, size n b
	δu bb ′	Vector of coordinates in differential Fourier
		plane, of length 3.
	l d	Vector of sky coordinates, of length 3.
	Matrices
	V tν
	K d p,tν

pq Visibility seen by a baseline pq at time and frequency t, ν. Size 2 × 2. J p,tν Jones matrix for antenna p for one (t, ν) pair. Contains the gains information. Size 2 × 2 B Sky brightness distribution matrix, of size 2 × 2. N Noise matrix, of size 2 × 2. Contains a single realisation n of the thermal noise in each cell. F Fourier transform matrix, of size n pix × n pix . S b

Table 3 .

 3 1. Table recapitulating the meaning and dimensions of vectors and matrices used in Sec. 3.2.1. Only scalars which give matrix dimensions or indices are given here.

Table 5 . 2 .

 52 Table giving the positions of our 8 LBCS sources and their distance from both the centre of the EGS and 3C295.

	[dms] Dist. from EGS [deg.] Dist. from 3C295 [deg.]
	1 14:30:18.72 52:17:29.80	2.041	2.904
	2 14:19:44.44 54:23:04.58	1.928	2.517
	3 14:21:20.05 53:03:46.00	0.864	1.743
	4 14:21:09.41 51:22:32.46	1.294	1.728
	5 14:11:50.32 52:49:02.66	0.844	0.619
	6 14:11:20.23 52:12:04.30	0.915	0.000
	7 14:08:07.00 52:55:11.36	1.409	0.869
	8 14:08:09.76 52:44:46.56	1.354	0.680

  Algorithmic Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 6.1.1 The Cov-Cov Relationship . . . . . . . . . . . . . . . . . . . . . . 142 6.1.2 Mapping the variance of images made with interferometric data143 6.1.3 Quality-based Weighting Schemes . . . . . . . . . . . . . . . . . 144 6.2 Imaging the Extended Groth strip . . . . . . . . . . . . . . . . . . . . . 145 6.2.1 Imaging the LOFAR primary beam . .

	Chapter 6
	Conclusion & Future Work
	Contents
	6.1

. . . . . . . . . . . . . . 146 6.2.2 Imaging the EGS with LOFAR-VLBI . . . . . . . . . . . . . . . 146 6.2.3 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149

Table 6 .

 6 1. Table recapitulating the names, positions and likely associations of all 12 chosen sources in the primary beam. These sources were primarily chosen because they had peculiar or interesting diffuse emission. The associated LOFAR image thumbnails are also given for reference.

See here for the full text: http://home.agh.edu.pl/ ~evermind/jrrtolkien/mythopoeia.htm.

https://www.nrao.edu/whatisra/hist_reber.shtml

See https://www.haystack.mit.edu/ast/science/epoch/ for more information -the ultimate goal of this topic is to learn more about the process of structure formation in the universe: what exactly are the evolutionary links between the very smooth matter distribution revealed by CMB studies and the highly structured universe of galaxies and clusters thereof which we know today?

http://www-astro.physics.ox.ac.uk/ ~garret/teaching/

https://www.cv.nrao.edu/course/astr534/SelfAbsorption.html

The emission angle is the angle between the electron's total velocity (including z-axis) and its radial velocity.

i.e. to justify weighting visibilities by a function inversely proportional to the variance of the associated residual visibilities.

10 uv space is where interferometric measurement coordinates lie, and lm coordinates are the cardinal sines on which the inverse Fourier transform of the visibilities are

mapped. 11 The noise-PSF also relates δw to δn, as shown in our matrix formalism, but we do not explicitly reference this in the text for the sake of brevity, as visibility space is usually referred to as "the UV-plane" in literature, rather than "the UVW-space". For more information, see Chapter 3

And, famously, attempting to observe "intelligent" sources (cf Enriquez et al.

2017)2 Other kinds of resolution -in time and frequency, for example -are also extremely important to many astronomers, but are not what concern us here.3 It is also a function of technological factors, of course, but ceteris paribus, a more sensitive telescope means a telescope with a wider collecting area4 By this, we mean that we assume that an instrument is also designed to optimise resolution.

Source: http://www.astron.nl/eris2013/Documents/2_laing_fundamentals_interferometry. pdf, https://science.nrao.edu/science/meetings/2016/15th-synthesis-imaging-workshop/ SISS15Advanced.pdf

This maximises antenna sensitivity in a given direction, but is not done with LOFAR

Fringe stopping consists of down-converting the signal received by antennas by mixing with a local oscillator. The geometric delay compensation is usually applied on this down-converted signal. Source: http://www.gmrt.ncra.tifr.res.in/gmrt_hpage/Users/doc/WEBLF/LFRA/node74.html

So-called because, in the case where these angles are small, they correspond to positions along the cardinal directions.

The other main family of deconvolution algorithms, Maximum Entropy Methods (MEM), are not very widely-used at the time of writing, but remain an active area of research (cf.Chael et al. (2018)).

See Section 2.10

This reduces calibration to solving only for the intensity gains: the data can then only be used for intensity mapping (e.g.[START_REF] Jennison | IAU Symposium[END_REF]). This practice therefore precludes polarimetry.

On the discovery of self-calibration and its evolution in parallel to adaptive optics, see the chapter titled "The Almost Serendipitous Discovery of Self-Calibration" in "[START_REF] Kellerman | Serendipitous discoveries in radio astronomy : proceedings of a workshop held at the National Radio Astronomy Observatory[END_REF].

This is not to say that the brightness matrix must contain only point sources, but rather that moving our interferometer 500m to the East should not change its measured visibilities

Stokes parameters determine the polarisation of light; we will not get into more details in this manuscript, as it is hardly relevant to the topic at hand. Suffice to say that the 2-vector formalism helps account for the polarised nature of light.

Bear in mind the limitations of this, described in Section 2.2.2.

The n -1 term in the exponential is due to the fact that, at phase centre, kp = 0 for all u and n = √ 1 -0 -0 = 1.

This means that, for a given baseline pq, u will vary as a function of frequency.

The "apparent" sky is the 'true' sky attenuated by the "power beam", which corresponds to our instrument's directional sensitivity in the traditional view.

Noordam, J.E. 1996, AIPS++, Note 185: The Measurement Equation of a Generic Radio Telescope

http://www.jpier.org/PIERB/pierb22/05.10032105.pdf

Decorrelation can equivalently be referred to as "decoherence".

see Smirnov (2011b) and companion papers.

Detailed treatment given in[START_REF] Thompson | Interferometry and Synthesis in Radio Astronomy[END_REF], section 6.2, which shows how to obtain the ideal radiometer equation

This section is based onBonnassieux et al. (2017), which was published as part of this PhD.

Which are simply an artificial split of a single original problem; see Section 2.3 for more details.

The noise-PSF also relates δw to δn, as shown in the matrix formalism, but this is not explicitly referenced in the text since visibility space is usually referred to as "the UV-plane" in literature, rather than "the UVW-space".

As opposed to the Fourier-plane pixels, which are the elements of the grid onto which the measured visibilities are mapped for imaging.

This hypothesis is equivalent to assuming that the sky is dominated by well-separated point sources, where we take "well-separated" to mean that the sources are multiple PSF Full-Width Half-Maximum apart

Specifically, that the residual visibilities do not contain any true signal, i.e. that all sources have been properly subtracted from the visibilities

Provided that simple averaging is performed; the use of baseline-dependent window functions can alleviate this effect to some extent

This is the bandwidth for each channel of our observation.

Each station has a relatively-accurate clock associated with it, save for the Superterp stations which all share a single clock. Like all clocks, they have finite accuracy, and desynchronisation errors thus creep in over time.

The more uncorrupted (i.e. flagged for RFI) visibilities are available, the better the estimation of gain solutions -seeSmirnov & Tasse (2015) for the demonstration.

Namely DE602, SE607 and UK608 -these are a German, Swedish and British station respectively.

The evidence for this is in the previous chapter, as the test data is the same as used here; nevertheless, an extra test was carried out in the context of these calibration and imaging tests for the sake of completeness, and found that quality-based weighting schemes resulted in better results than simply flagging the last 2 hours.

This approach was taken due to time constraints: verifying that our model was compatible with the closure phases and amplitudes would have been a better test to perform, but would also likely have

In the context of this work, it would literally be an equivalent but faster implementation of exactly the same algorithm.

i.e. the application of a smoothing operation in gain phase and amplitude, justified by the fact that they are expected to be smooth and so jagged behaviour is likely erroneous.

This initial model was tested, but not followed through to the same extent as the 2-point-source initial model once it became clear that it was both biased away from the ground truth compared to the other two (as we could expect), but that it seemed to converge -albeit very slowly -towards the same ground truth as them.

Acknowledgements

A few paragraphs can never be room enough to express the debt of gratitude that I owe, let alone to the veritable crowd to whom I owe it. This section is thus by necessity incomplete -but I will nevertheless do my best to account, as much as possible, for those who made this thesis possible. I'll begin with thanks to Cyril Tasse, Philippe Zarka and Oleg Smirnov. Though each will claim to have only done a small part of the work in guiding and teaching me over the course of this PhD, they have all contributed monumentally -it was not always easy, but it was certainly worthwhile. I would also like to thank all the members of my thesis defense, many of whom crossed borders to attend. I would particularly like to thank my examiners -Neal Jackson, Patrick Charlot and Ciara Ferrari -whose thorough comments and invaluable insights greatly improved the quality of this work.

Similarly, I owe an infinite debt to the post-doctoral fellows and colleagues I had the pleasure to work with and befriend. I would particularly like to thank Marcel Atemkeng, Julien Girard, Roger Deane, Kshitij Thorat, Trienko Grobler, Arun Aniyan and Alan Loh -both for putting up with my ridiculous imprecations to my computer terminal, and for helping me refine my understanding of my own work. Of course, I am just as grateful to have known and worked with many other fine people: Corentin Louis, Lucas Grosset, and Lisa Bardou in France; Theophilius Ansah-Narh, Kelachukwu Iheanetu, Benjamin Hugo, Ulrich Mbou Sob, Kwazi Mthembu and all the RATT Masters students during my time in South Africa.

I would also like to take the time to thank the people who made this thesis possible on a logistical level: nothing could have been done without Emmanuel Thetas and the nancep cluster. But on a more fundamental level, I would like to thank secretaries in both France and South Africa: without Jacqueline Plancy, Géraldine Gaillant, and Ronel Groenwald, neither the Observatoire de Paris nor the RATT would be what they are now.

Outside the lab, I would like to extend my warmest thanks to my wonderful partner Christine, for also putting up with my ridiculous imprecations to my computer at home and for holding on while I was in South Africa. I also thank my family, whose support was often invaluable -my parents Beatrice and Marc, and my sisters Violette and Coline. To my friends scattered across the world: Laith, Bob, Erith, Alex, Cameron, Craig, Ru, Alasdair, Rachel, Sammie, Michael, Jois, Poppy, Paul, Celine, Philippe, Helene, Sarah,

https://www.astron.nl/radio-observatory/ astronomers/technical-information/http://www.lofar.org/about-lofar/system/lofar-numbers/ LOFAR.org website: http://www.lofar.org/about-lofar/system/lofar-numbers/