ABSTRACT

Imaging Atmospheric Cherenkov Telescopes (IACTs) currently in operation feature large mirrors and order of 1 ns time response to signals of a few photo-electrons produced by optical photons. This means that they are ideally suited for optical interferometry observations. Thanks to their sensitivity to visible wavelengths and long baselines optical intensity interferometry with IACTs allows reaching angular resolutions of tens to microarcseconds. We have installed a simple optical setup on top of the cameras of the two 17 m diameter Major Atmospheric Gamma-Ray Imaging Cherenkov (MAGIC) IACTs and observed coherent fluctuations in the photon intensity measured at the two telescopes for three different stars. The sensitivity is roughly 10 times better than that achieved in the 1970s with the Narrabri interferometer.

1 HISTORICAL INTRODUCTION

Hanbury-Brown & Twiss (1954, 1956) demonstrated that at sufficiently short baselines, the spatial correlation measurements of thermal light sources, such as stars, should exhibit a photon bunching signal for measurements within a defined area on the ground. This area depends on the angular size of the star so sampling the intensity of the correlated signal allows measuring the size of the star.

After initial tests on Sirius a two 6.5 m diameter telescope interferometer with adjustable baseline was built in Narrabri, Australia (Hanbury-Brown 1974). This interferometer was used to measure the diameters of the 32 brightest stars in the Southern hemisphere. Because the interferometer operated at visible wavelengths and the distance between the two telescopes was as large as 188 m Narrabri’s interferometer could measure diameters as small as ∼400 μas.

Such ‘intensity interferometry’ has hardly developed in the following decades probably because, contrary to ‘amplitude interferometry’ (Glindemann et al. 2000; ten Brummelaar et al. 2015), it calls for large photon collection areas, sensitivity to a few photons, and fast time response in the order of a few nanoseconds down to picoseconds. The situation may have changed recently due to the availability of telescopes with similar requirements: the so-called Imaging Atmospheric Cherenkov Telescopes (IACTs). IACTs (Hillas 2013) detect faint (few photons) and fast (ns) pulses of near UV to blue Cherenkov light emitted by particle showers produced by cosmic or gamma rays above tens of GeV. Current and future IACTs may be modified to work as intensity inferferometers (Dravins et al. 2012, 2013; Matthews et al. 2019).

Besides, conventional photomultipliers (PMTs) have time resolutions exceeding 1 ns but other photodetectors are even faster, potentially reaching <1 ns resolutions. Avalanche photodiodes have been recently applied to intensity interferometry (Guerin et al. 2018).

2 SHORT THEORETICAL REVIEW

The reader can find a detailed discussion on the theory of the intensity interferometer and the methods to analyse the correlation signal in the book by Hanbury-Brown (1974). Here, we will only review a few definitions and expressions that are essential for our analysis.

Consider two telescopes pointing at a star, as sketched in Fig. 1. Let I1(t) and I2(t) be the signal intensities registered at each of the telescopes as a function of time.

Sketch of the MAGIC telescopes and the new components that have been installed to enable it as an optical intensity interferometer.
Figure 1.

Sketch of the MAGIC telescopes and the new components that have been installed to enable it as an optical intensity interferometer.

The second-order coherence function g(2) is defined as:
(1)
where 〈〉 represents averaging over time and τ is a time delay introduced by the different light and signal paths in the two telescopes.
For a randomly polarized light source and a light detector whose electronic bandwidth Δf is much smaller than the optical bandwidth Δν it can be shown that:
(2)
where the complex visibility V12 is the Fourier transform of the source brightness distribution.
In the remaining of this work we will rather use the normalized contrast c, defined as:
(3)
It is straightforward to show that:
(4)
For a star that is described as a disc of uniform intensity up to an angular diameter θ it can be proven that:
(5)
where B1 is the Bessel function of the first order, λ is the light wavelength, and d is the distance between the telescopes (baseline).
Dividing c for a given baseline d and c at baseline 0 (or a small enough baseline) returns |V12(d)|2:
(6)
c(0) is in effect a constant of the setup.

In conclusion, c decreases with increasing baseline following equation (5). A measurement of the curve c(d)/c(0) versus d allows us to determine the angular size of the star.

3 DESCRIPTION OF MAGIC AND THE TECHNICAL MODIFICATIONS TO IMPLEMENT THE INTERFEROMETER

MAGIC is a system of two IACTs located at the Roque de los Muchachos Observatory (Observatorio del Roque de los Muchachos, ORM) on the island of La Palma in Spain (Aleksic et al. 2016a). Equipped with 17 m diameter mirror dishes and fast PMT cameras, the telescopes record images of extensive air showers in stereoscopic mode, enabling the observation of Very High Energy (VHE) γ-ray sources at energies ≳50 GeV (Aleksic et al. 2016b).

The telescope reflector follows a parabolic shape to minimize the time spread at the focal plane. The focal ratio is 1. The reflector is formed by ∼250 1 m2 spherical mirror tiles. Two actuators behind each tile can be used to shift its orientation and correct the reflector for deformations in its parabolic shape as the telescope changes zenith. This so-called Active Mirror Control (AMC) adjusts the mirror in a few seconds and runs typically every 20 min during standard observations. Approximately 70 per cent of the light of a point source is focused on a pixel. It must be noted that individual mirror tiles are slightly staggered thus introducing an extra time spread in the order of 200 ps.

Each telescope is equipped with a 1039-pixel PMT camera at the primary focus. The PMTs are 25.4 mm diameter. A hexagonal shape Winston cone is mounted on top on each PMT. The distance between PMT centres is 30 mm. Each pixel covers a 0.1° Field Of View (FOV).

The PMTs are Hamamatsu, type R10408, with a hemispherical photocathode and 6 dynodes. The PMT bias voltages for the cathode and dynodes are generated by a Cockroft–Walton DC–DC converter, which can provide up to 1250 V peak voltage. The electrical signals are amplified (AC coupled, ∼25 dB amplification) and then transmitted via independent optical fibers by means of vertical cavity surface emitting lasers (VCSELs). The average pulse width signal is measured to be 2.5 ns (full width at half-maximum, FWHM).

The PMTs operate at a rather low gain of typically 3–4 × 104 in order to also allow observations under moonlight without damaging the dynodes.

Fig. 1 shows a sketch of the two MAGIC telescopes along with a few of the elements that were installed for the interferometric observations.

For interferometry observations we only use the signal from the central PMT of each of the two telescopes. We operate the detector in so-called ‘analogue mode’, that is, we continuously sample and record the signal at the PMTs in a regime where each ∼2–4 ns sample corresponds to tens of photoelectrons. The sensitivity of the detector is proportional to the photon flux detected by the PMT and in first approximation does not depend on the optical bandpass of the detected light (see equation 16 below).

However we have added a filter in front of the central pixels to make the anode currents of the PMTs (DCs) low enough that PMTs are not damaged when pointing at bright stars or during bright Moon observations. The filter is held 15 mm in front of the Winston cone using a mechanical frame previously developed for bright Moon VHE observations (Ahnen et al. 2017). This frame must be mounted and dismounted manually, which generally prevents standard VHE observations for the whole night.

We selected the interference filter Semrock 432/36 nm BrightLine of 32 mm diameter and 3.5 mm thickness. The filter transmission is >90 per cent for 414–450 nm. However light enters the filter directly from the mirror so it follows a wide range of incident angles up to 26.6°. As a result the effective filter transmission curves shifts to shorter wavelength. The filter has a mild dependence on incident angle (effective refractive index 2.02 for both polarizations). Fig. 2 shows the effective transmission range as obtained with a simulation that takes into account the shape of the telescope reflector, the geometry of both filter and pixel, the filter transmission for normal incidence and the dependence on transmission wavelength with incident angle. The centre wavelength shifts to ∼427 nm and the FWHM of the transmission curve is roughly the same as the original one but the transmission curve develops two tails at short and long wavelengths. This actually reduces the sensitivity of the interferometer by about 15 per cent (see Section 5).

Effective transmission of the filter installed in front of the MAGIC central pixels as a function of wavelength, as calculated using a Monte Carlo simulation of the optical system.
Figure 2.

Effective transmission of the filter installed in front of the MAGIC central pixels as a function of wavelength, as calculated using a Monte Carlo simulation of the optical system.

During normal VHE observations, the reflector is focused at a distance of ∼10 km so that the images of atmospheric showers produced by γ-rays are as sharp as possible. For our observations, we would rather focus the telescope at the star i.e. at infinity. AMC allows to move the focus in all three directions. We shifted it along the optical axis in small steps of a few mm and measured DC in the central pixels until we reached a maximum DC when focusing about 33 mm beyond the nominal focus for VHE observations. The filter has essentially the same physical size of a Winston cone and is situated some 20 mm ahead of it so a fraction of the light is shadowed. A simulation shows that this fraction is 18 per cent. This fraction was confirmed comparing DC measurements with and without the filter frame.

The MAGIC cameras transmit the analog signal of each pixel as an optical signal via fiber optics to a separate readout location. The fibres are ∼162 m long. The analogue optical signal from each of the two telescopes is generated by means of a multimode VCSEL present in each pixel of the camera. The 780-nm multimode optical signal is converted to an electrical signal at the readout end of the fiber using a biased photodiode and a 40 dB wide-band amplifier located at the end of the optical fiber. The signal is then digitized using two channels of a Rohde+Schwarz RTO2044 oscilloscope. The 2.5 ns single photo-electron pulse width response of the PMTs corresponds to an input signal bandwidth of approximately 110 MHz. The oscilloscope was used in high definition mode, which provides an accuracy of 13 bits with an input bandwidth of 200 MHz.

Sampling frequencies between 250 Msamples per second (250 MSps) and 1 GSps were considered for signal acquisition. The data were saved in records consisting of 10 times 100 MS records, for a total recording time of 1 s @ 1 GSps to 4 s @ 250 MSps. The 4 GB records were transferred via USB3 from the oscilloscope to an external solid state disc (SSD), resulting in a duty cycle of approximately 10 per cent at 500 MSps. Since the rate of data transfer to disc remains more or less constant independent of sampling rate, the duty cycle and thus the observation time in the sensitivity calculations scales with the sampling rate. The achievable bandwidth also depends on the sampling rate. As a result, the sensitivity is theoretically independent of the sampling rate chosen, since an increase in duty cycle is offset by a decrease in bandwidth. In practice, lower bandwidths (sampling rates) are less sensitive to delay differences due to staggering of mirror tiles, delay calculations, etc.

Before committing to observing on site, different sampling rates and the implementation of offline analysis software were investigated using an optical signal model in the laboratory. A test signal was generated using a Nichia laser diode which normally is used to produce coherent light at 405 nm. The diode was used as an LED at currents well below the laser threshold of the diode, but capable of fast modulation. A DC current was used to produce thermal emission, and an additional RF noise source with a 450 MHz high-pass filter was used to add an additional optical component with a short random temporal correlation. This mimics the expected optical signal whose temporal coherence is also much shorter than the response time of the PMTs.

This test signal was used to verify the data analysis software and compare different sampling rates by means of the operation of two optical modules as used in the MAGIC cameras in a dark box with connection to the acquisition setup described (see Fig. 3).

Sketch of the laboratory setup that was used to test the readout.
Figure 3.

Sketch of the laboratory setup that was used to test the readout.

A low pass flat time delay filter with a cut-off frequency of 117 MHz (Mini circuits SBLP-117+) was added to each channel between the amplifier and the oscilloscope to match the phase and frequency response of both channels and as a precaution against external interference signals.

4 ANALYSIS METHOD

4.1 Contrast and Pearson’s correlation function

As mentioned above the electronic chain of the MAGIC telescopes is AC coupled. Frequencies below 10 kHz are greatly attenuated. This means that we do not know the DC value of I1 or I2. Therefore, we cannot calculate an absolute value for g(2) or c.

However our goal is not to measure c but rather to determine the curve c(d)/c(0) versus d and fit it to equation (5) in order to extract the star’s angular diameter. As we shall see we can extract this curve even if we do not know 〈Ii〉.

To begin with we measure the Pearson’s correlation coefficient ρ:
(7)
because Ii(t) − 〈Ii〉 are fast variations with respect to the pedestal i.e. they correspond to our actual measurement. We can use ρ and its uncertainty to establish the strength of the correlation signal.
On the other hand 〈I1〉 and 〈I2〉 are respectively proportional to the anode current of the photodetectors (conventionally referred to as DC1 in MAGIC-1 and DC2 in MAGIC-2) with a fixed proportionality factor as long as the gains of the PMTs remain constant:
and the fluctuations of I1 and I2 are Poissonian, that is, proportional to the square root of the anode currents:
As a consequence:
(8)
where K is constant if the gain of the PMTs (i.e. the high voltage) remain constant. During our observations, the gain actually changed so we must write:
(9)
where G1 and G2 are factors that are obtained from dedicated sets of DC measurements at different HV settings: see Table 3.

DC1 and DC2 hardly change in time scales of a few seconds but depend strongly on the star and observation conditions on a longer time-scale. In MAGIC DC1 and DC2 are measured roughly every second. We will calculate ρ on the shortest possible time-scales (subruns). To combine the data on longer time-scales we will use the DC1 and DC2 measurements to determine c/K. The same procedure was used in the Narrabri observations (see Hanbury-Brown 1974, chap. 10).

The final goal is to determine c(d)/K for a given star over a range of baselines. Using equation (6) would then allow us to determine the star’s angular diameter.

4.2 Computational method

The computation time is driven by the calculation of the Pearson’s coefficient ρ. Each data subrun is split into N non-overlapping time windows of a fixed number of samples, S. Within each time window, the time is discretized, and runs from i = 1 to S, therefore the signals are given by I1(i) and I2(i). For each one of these time windows, ρ is estimated by applying equation (7) to its samples. Subsequently, the final estimate of ρ for a single data run is obtained by averaging the result obtained for each time window:
(10)
where ρi is the Pearson’s coefficient computed for time window i. In addition, the statistical error of ρ(τ) is estimated as:
(11)
where Δ(ρi(τ)) is the rms of the sample of the ρ estimates obtained from the N time windows.
The only non-trivial averaging that has to be computed in the definition of ρ (equation 7) to obtain equation (10) is the cross term 〈I1(t)I2(t + τ)〉. For non-overlapping windows, it is obtained as
(12)
where S is the number of samples, τ is given in discrete sample units, and W(t) is a windowing function that is identical to one if and only if t is within the time window and zero otherwise. The windowing function must be introduced to keep the computation of the estimate independent between different time windows. As a drawback, the statistics available for each value of τ is different, as it is reflected in the denominator of equation (12), which can be trivially computed to give
(13)
where |x| denotes the absolute value of x.
The direct computation of the numerator of equation (12) for all non-trivial values of τ within a time window is time consuming, therefore we make use of the Convolution Theorem for discrete Fourier transforms (Oppenheim et al. 1999), that can be formulated as that given two periodic discrete functions f and g, with period P, then
(14)
where |$\mathcal {F}(x)$| denotes discrete Fourier’s transform of x, x* denotes complex conjugate of x, and |$\mathcal {F}^{-1}$| denotes the inverse discrete Fourier’s transform. Moreover, the discrete Fourier’s transforms can be efficiently computed using any of the fast Fourier Transform (FFT) family of algorithms.
In order to make use of equation (14), it is necessary to adapt equation (12). For a given time window this is achieved by defining two periodic functions |$\tilde{I}_1(i)$| and |$\tilde{I}_2(i)$| of period 2S, identically to I1(i) and I2(i) respectively for imod 2S ≤ S and 0 for imod 2S > S. With this definition, if τ is restricted to |τ| < S − 1, the equation (12) can be written as
(15)
which is the final procedure used to compute the cross-term in the ρi estimates in equation (10).

The size of the time windows used in the computation of ρ obtained by optimizing the time required by the analysis. With the computational architecture and the FFT algorithm we have used, the optimal has been found to be S = 16 384 samples.

4.3 Correction for the Moon

The registered DCs are in fact produced by the combination of diffuse Moonlight and light coming from the star. During our actual observations (described in Section 6), the DC produced by the Moon are up to ∼15 per cent of the DC of the star.

We will subtract the DCs produced by the Moon in each individual telescope from the measured DCs to calculate the DCs expected only from the star. We will apply the latter to equation (9) to determine the normalized contrast due to the star alone.

In fact, except for a single measurement on the first night we did not measure systematically the DCs produced only by the Moon. Instead we used a model to calculate the Moon DCs that has been validated with previous VHE Moon observations (Britzger 2009).

5 SENSITIVITY ESTIMATE AND SELECTION OF TARGET STARS

An observational campaign was planned for the period around the Full Moon of 2019 April, that is, during nights when the sky brightness is too high to schedule regular VHE observations.

We can calculate the significance of the measured signal and the sensitivity of our setup using exclusively ρ. Hanbury-Brown (1974, equation 5.17) allows to calculate the significance (signal over noise) of the correlation for a given experimental setup and unpolarized light. The equation can be written as:
(16)
where A is mirror area, α(λ0) is the quantum efficiency (QE) of the PMTs for the filter’s central wavelength λ0, q(λ0) is the QE of the remaining optics, n0) is the star’s differential photon flux, bv is the effective cross-correlation electrical bandwidth, F is the excess noise factor of the PMT, and T is the observation time. Finally, σ is the normalized spectral distribution of the light after the filter as defined in formula (5.6) of Hanbury-Brown (1974). σ would be equal to 1 if the filter transmission curve is a boxcar function and the spectrum of the light is flat, whereas it reduces to 0.87 for the transmission curve in Fig. 2.

We have assumed that both telescopes are identical and we have neglected the effect of the Moon and additional noise in the readout chain. As ρ is proportional to c for short time periods we may use the same expression to calculate the significance of ρ for a subset of data of a few seconds.

We used equation (16) to calculate the expected detection time of candidate stars. Table 1 shows the parameters of our setup.

Table 1.

Setup parameters used to estimate the expected detection times: mirror area, QE of PMTs, QE of remaining optics, effective cross-correlation electrical bandwidth, normalized spectral distribution of the light, and inverse of excess noise factor of the PMT. We assume that bv is similar to the electrical bandwidth of the individual channels.

A250 m2
α29.5 percent
q0.236
bv110 MHz
σ0.87
F−10.87
A250 m2
α29.5 percent
q0.236
bv110 MHz
σ0.87
F−10.87
Table 1.

Setup parameters used to estimate the expected detection times: mirror area, QE of PMTs, QE of remaining optics, effective cross-correlation electrical bandwidth, normalized spectral distribution of the light, and inverse of excess noise factor of the PMT. We assume that bv is similar to the electrical bandwidth of the individual channels.

A250 m2
α29.5 percent
q0.236
bv110 MHz
σ0.87
F−10.87
A250 m2
α29.5 percent
q0.236
bv110 MHz
σ0.87
F−10.87

We used the Jean–Marie Mariotti Center Stellar Diameter Catalogue (Cheli 2016) to identify the brightest stars in B observable from La Palma during that period and having angular diameters between 0.3 and 1.0 mas so as to maximize the expected correlation signal. Table 2 shows the time that is necessary to reach a signal of 5σ significance following equation (16) for the three best candidate stars. The differential flux after extinction n was calculated from the magnitude in Ducati (2002) and the atmospheric extinction calculated for ORM and the zenith angle of the star. We have assumed an extinction coefficient k = 0.23 mag per airmass for the filter’s wavelength based on King (1985). d has been calculated from the positions of the telescopes and the star’s zenith and azimuth angles.

Table 2.

Brightness and angular diameter of the three best candidates for interferometric observations, together with the zenith angle and baselines of the observation campaign, the expected visibility, and the expected and measured time for a 5σ detection signal.

StarmBθTimezddV2ExpectedMeasured
(mas)(UTC)(°)(m)T (s)T (s)
Adhara, ϵ CMa1.290.77 ± 0.0521:006734.10.803638 / 35
22:007625.90.885457 / 64
23:008628.00.882160
Benetnasch, η UMa1.670.82 ± 0.0321:005355.30.50126
22:004266.00.36210
23:003373.20.27336354
00:002578.60.22510
Mirzam, β CMa1.730.50 ± 0.0321:006441.30.876073
22:007439.80.88102146
23:008644.50.859780
StarmBθTimezddV2ExpectedMeasured
(mas)(UTC)(°)(m)T (s)T (s)
Adhara, ϵ CMa1.290.77 ± 0.0521:006734.10.803638 / 35
22:007625.90.885457 / 64
23:008628.00.882160
Benetnasch, η UMa1.670.82 ± 0.0321:005355.30.50126
22:004266.00.36210
23:003373.20.27336354
00:002578.60.22510
Mirzam, β CMa1.730.50 ± 0.0321:006441.30.876073
22:007439.80.88102146
23:008644.50.859780
Table 2.

Brightness and angular diameter of the three best candidates for interferometric observations, together with the zenith angle and baselines of the observation campaign, the expected visibility, and the expected and measured time for a 5σ detection signal.

StarmBθTimezddV2ExpectedMeasured
(mas)(UTC)(°)(m)T (s)T (s)
Adhara, ϵ CMa1.290.77 ± 0.0521:006734.10.803638 / 35
22:007625.90.885457 / 64
23:008628.00.882160
Benetnasch, η UMa1.670.82 ± 0.0321:005355.30.50126
22:004266.00.36210
23:003373.20.27336354
00:002578.60.22510
Mirzam, β CMa1.730.50 ± 0.0321:006441.30.876073
22:007439.80.88102146
23:008644.50.859780
StarmBθTimezddV2ExpectedMeasured
(mas)(UTC)(°)(m)T (s)T (s)
Adhara, ϵ CMa1.290.77 ± 0.0521:006734.10.803638 / 35
22:007625.90.885457 / 64
23:008628.00.882160
Benetnasch, η UMa1.670.82 ± 0.0321:005355.30.50126
22:004266.00.36210
23:003373.20.27336354
00:002578.60.22510
Mirzam, β CMa1.730.50 ± 0.0321:006441.30.876073
22:007439.80.88102146
23:008644.50.859780

The best candidates for a detection are Adhara (ϵ CMa), Benetnasch (η UMa), and Mirzam (β CMa). The next three stars are η CMa, δ Sco, and γ Crv but the expected detection times exceeded 5 min. Since we were uncertain about a few of the setup parameters (bv or F−1 for the specific central PMTs) we decided to take data only on Adhara, Benetnasch, and Mirzam.

Unfortunately, there were no stars of similar brightness with significantly smaller θ. As we shall see this prevented us from making an absolute measurement of θ for the observed stars.

6 OBSERVATIONS

We installed the interferometry setup in MAGIC and took data on the three selected stars for five nights with bright and Full Moon in 2019 April. A log of the observations can be found on Table 3.

Table 3.

Log of the interferometric observations made with MAGIC in 2019 April. Tabulated are the start day of the observations, the star name, time in UTC, general weather conditions, sampling rate in MSamples per second, anode current, and high voltage in the central pixels of MAGIC-1 (M1) and MAGIC-2 (M2).

CampaignDateStarTimeSamplingM1 DCM2 DCM1 HVM2 HV
dayMJD(UTC)rate(μA)(μA)(V)(V)
12019/04/15Adhara21:59–22:34500 MSps37–2133–1510211063
MJD 58588Benetnasch22:47–23:58...42–4441–44961983
22019/04/16Adhara21:11–21:41...UnstableUnstable10211063
MJD 58589Benetnasch21:45–23:53...UnstableUnstable961983
32019/04/17Adhara21:02–22:51...39–635–7911954
MJD 58590Benetnasch22:55–23:58...40–4439–42961983
42019/04/18Adhara20:55–22:46250 MSps43–743–8931983
MJD 58591Benetnasch22:50–23:58...46–4342–30961983
52019/04/19Mirzam20:56–22:38...44–940–1010211063
MJD 58592Benetnasch22:43–23:59...45–4442–30961983
CampaignDateStarTimeSamplingM1 DCM2 DCM1 HVM2 HV
dayMJD(UTC)rate(μA)(μA)(V)(V)
12019/04/15Adhara21:59–22:34500 MSps37–2133–1510211063
MJD 58588Benetnasch22:47–23:58...42–4441–44961983
22019/04/16Adhara21:11–21:41...UnstableUnstable10211063
MJD 58589Benetnasch21:45–23:53...UnstableUnstable961983
32019/04/17Adhara21:02–22:51...39–635–7911954
MJD 58590Benetnasch22:55–23:58...40–4439–42961983
42019/04/18Adhara20:55–22:46250 MSps43–743–8931983
MJD 58591Benetnasch22:50–23:58...46–4342–30961983
52019/04/19Mirzam20:56–22:38...44–940–1010211063
MJD 58592Benetnasch22:43–23:59...45–4442–30961983
Table 3.

Log of the interferometric observations made with MAGIC in 2019 April. Tabulated are the start day of the observations, the star name, time in UTC, general weather conditions, sampling rate in MSamples per second, anode current, and high voltage in the central pixels of MAGIC-1 (M1) and MAGIC-2 (M2).

CampaignDateStarTimeSamplingM1 DCM2 DCM1 HVM2 HV
dayMJD(UTC)rate(μA)(μA)(V)(V)
12019/04/15Adhara21:59–22:34500 MSps37–2133–1510211063
MJD 58588Benetnasch22:47–23:58...42–4441–44961983
22019/04/16Adhara21:11–21:41...UnstableUnstable10211063
MJD 58589Benetnasch21:45–23:53...UnstableUnstable961983
32019/04/17Adhara21:02–22:51...39–635–7911954
MJD 58590Benetnasch22:55–23:58...40–4439–42961983
42019/04/18Adhara20:55–22:46250 MSps43–743–8931983
MJD 58591Benetnasch22:50–23:58...46–4342–30961983
52019/04/19Mirzam20:56–22:38...44–940–1010211063
MJD 58592Benetnasch22:43–23:59...45–4442–30961983
CampaignDateStarTimeSamplingM1 DCM2 DCM1 HVM2 HV
dayMJD(UTC)rate(μA)(μA)(V)(V)
12019/04/15Adhara21:59–22:34500 MSps37–2133–1510211063
MJD 58588Benetnasch22:47–23:58...42–4441–44961983
22019/04/16Adhara21:11–21:41...UnstableUnstable10211063
MJD 58589Benetnasch21:45–23:53...UnstableUnstable961983
32019/04/17Adhara21:02–22:51...39–635–7911954
MJD 58590Benetnasch22:55–23:58...40–4439–42961983
42019/04/18Adhara20:55–22:46250 MSps43–743–8931983
MJD 58591Benetnasch22:50–23:58...46–4342–30961983
52019/04/19Mirzam20:56–22:38...44–940–1010211063
MJD 58592Benetnasch22:43–23:59...45–4442–30961983

The weather conditions were poor during the second night of observations so in what follows we will ignore the two data samples taken on that night.

Data were taken with two sampling rates: 250 and 500 MSps to investigate the effect of the sampling rate in the sensitivity. We recorded a total of 8.2 TB of data over the whole campaign.

Different sets of PMT high voltages HV were applied with the goal to keep the DCs below a safe value of ∼40 μA but as high as possible to improve the PMT single photo-electron response and to reduce the impact of electronic noise.

In segmented acquisition mode, the oscilloscope takes data of both channels simultaneously for 100 million samples and then stops digitizing while it saves the data to internal memory. 10 such ‘subruns’ of 100 Msamples are successively recorded and these 10 subruns are saved together as a ‘run’ to the SSD adding even more dead time. The global dead time is about 90.7 per cent for 500 MSps sampling date and 81.4 per cent for 250 MSps. In other words, the duty cycle is correspondingly about 10 per cent and 20 per cent.

Adhara and Mirzam are observed long after transit so the DC strongly decreases along the observations whereas Benetnasch is observed during transit and the DC is roughly constant.

7 ANALYSIS RESULTS

7.1 Significance of the correlation signals

We have calculated Pearson’s correlation coefficient ρ for all data subruns as a function of time delay τ.

Once the variable term of τ is corrected for, a clear correlation signal is detected in all samples at a fixed τ ∼4 ns. As an illustration, Fig. 4 shows the mean ρ for one of Adhara’s data samples. The significance can be estimated either with respect to the rms of ρ far from the signal region or with respect to the rms of ρ calculated by dividing the total sample in smaller subsamples. Both methods are consistent. For this specific data sample the significance is 15.4σ.

ρ as a function of τ between the two telescopes for Adhara’s observation on the third night. The variable delay term due to telescope orientations and positions has been subtracted. A signal with a significance of around 15σ is observed near zero delay for an effective observation time of ∼600 s.
Figure 4.

ρ as a function of τ between the two telescopes for Adhara’s observation on the third night. The variable delay term due to telescope orientations and positions has been subtracted. A signal with a significance of around 15σ is observed near zero delay for an effective observation time of ∼600 s.

Table 4 tabulates the significance of the correlation signal measured for all data samples.

Table 4.

Significance of the correlation signal measured for the four nights with good sky conditions.

StarCampaignEffective timeSignificance
day(s)(σ)
Adhara12005.3
Adhara360015.4
Adhara4104012.6
Benetnasch14025.8
Benetnasch33504.9
Benetnasch47607.3
Benetnasch58805.1
Mirzam510729.3
StarCampaignEffective timeSignificance
day(s)(σ)
Adhara12005.3
Adhara360015.4
Adhara4104012.6
Benetnasch14025.8
Benetnasch33504.9
Benetnasch47607.3
Benetnasch58805.1
Mirzam510729.3
Table 4.

Significance of the correlation signal measured for the four nights with good sky conditions.

StarCampaignEffective timeSignificance
day(s)(σ)
Adhara12005.3
Adhara360015.4
Adhara4104012.6
Benetnasch14025.8
Benetnasch33504.9
Benetnasch47607.3
Benetnasch58805.1
Mirzam510729.3
StarCampaignEffective timeSignificance
day(s)(σ)
Adhara12005.3
Adhara360015.4
Adhara4104012.6
Benetnasch14025.8
Benetnasch33504.9
Benetnasch47607.3
Benetnasch58805.1
Mirzam510729.3

Both Adhara and Mirzam were setting. The zenith angle increases quickly but the baseline hardly changes. As a consequence |V|2 only changes by ∼10 per cent for Adhara and even less for Mirzam. In the contrary, the baseline and consequently |V|2 of Benetnasch change by almost a factor 2 during the first 2 h of the night. We only targeted this star after 22:30 UTC because we prioritized Adhara and Mirzam for which we expected a stronger signal. In conclusion, the observations cover only a small range in |V|2 for each individual star and differences in detection time for each of them are dominated by changes in extinction, i.e. depend on zenith angle.

For Adhara and Mirzam, the predicted times are in the order of 1 min at the beginning of the night and degrade to tens of minutes in a matter of two hours. For Benetnasch, the times are longer but the dependence with zenith angle is milder.

Let us compare the strength of the signal with the expectations in Table 2. Since the detection times changes rather quickly along the night we make the comparison for short samples around the times in the table and we restrict ourselves to the data samples with more statistics, that is, Adhara nights 3 and 4, Benetnasch night 4 and Mirzam. The times that were needed to measure a 5σ signal are tabulated in the last column of Table 2. We have generally underestimated the detection times, but the expected and measured times match within 20 per cent. The only exception, a measured time 40 per cent longer than the expected one, is found for Mirzam around 22:00 UTC. Considering that we have not measured some of the parameters affecting the sensitivity (bandwidth, noise factor, atmospheric extinction, etc.) the match is rather satisfactory.

Hanbury-Brown (1974, see equation 11.12 and fig. 11.10) concluded that the significance of an unresolved star at 45° elevation measured with the Narrabri interferometer scales with observation time and magnitude in B-band mB as:
(17)
where T0 is the time of observation in seconds. A comparison with our results shows that our current interferometer is about 10 times more sensitive than Narrabri assuming that our duty cycle would be 100 per cent.

As a general indication, with this sensitivity we would detect a 5σ signal from a mB = 5 mag unresolved star in about 3 h. Resolving the star however demands to sample the correlation function down to, say, |V|2 = 0.3. For an mB = 4 mag star this would be possible in about 5 h.

It must be said that our longest observation was 15 min and we only targeted a limited number of stars so we have little control over our systematic errors and this sensitivity remains somewhat speculative.

7.2 Contrast and correlation function

We have used expression (8) to calculate c in arbitrary units for all runs.

Fig. 5 shows c (in arbitrary units) as a function of run number for Adhara and the third night. A fit to a straight line shows that there is no significant dependence on run number even if DC decreased by a factor ∼5 during the observations. This is a good crosscheck that equation (8) can be used to estimate c. The decrease in DC and correspondingly photon flux from the star also results in an increase in the uncertainty of c.

Contrast c as a function of run number for Adhara’s observation on the third night. The units in the Y-axis are arbitrary. p0 and p1 are the free parameters in the linear fit c = p0 + p1 · Run.
Figure 5.

Contrast c as a function of run number for Adhara’s observation on the third night. The units in the Y-axis are arbitrary. p0 and p1 are the free parameters in the linear fit c = p0 + p1 · Run.

A linear fit was applied to all samples and used to estimate c for each star and night. All fits were successful although a small dependence with run number was still left in a few samples. This probably indicates that the effect of the Moon has not been fully removed.

Even if we do not know c in absolute terms we can crosscheck if it correlates linearly with |V|2 for all stars and data samples under study, as expected from equation (6).

Fig. 6 shows c as a function of |V|2 for all stars and observation nights. We have made a linear fit to the points. The results of the fit are shown on the same figure.

Contrast c calculated from the data with equation (8) as a function of |V|2 determined from the baselines during the observations and star diameters in Table 2. p0 and p1 are the free parameters in the linear fit c = p0 + p1 · |V|2.
Figure 6.

Contrast c calculated from the data with equation (8) as a function of |V|2 determined from the baselines during the observations and star diameters in Table 2. p0 and p1 are the free parameters in the linear fit c = p0 + p1 · |V|2.

For each star c is constant within the statistical error. This is consistent with the fact that the baseline hardly changes during the observations.

The linear fit is successful and the fit goes through the origin, as expected from equation (6).

This plot may actually be used to calibrate our setup and measure the diameter of other stars. All in all these results are encouraging and indicate that MAGIC can already be used as an intensity interferometer.

However it must be stressed that we have only taken a few minutes of data on a few nights of observations so we have hardly any control over our systematics. Only to mention a few effects we may suffer from sporadic electronic noise, telescope tracking may be unstable for some positions of the sky or the time delay between the two telescopes may slightly change over time-scales of days.

We must quantify these effects and in general take longer data samples on a number of stars up to at least magnitude 4, sampling the whole sky and different observation conditions, before producing scientific results.

8 CONCLUSIONS AND OUTLOOK

We have searched for time correlation in the fluctuations of the signals produced by three stars of θ ∼0.5–1 mas in the central PMTs of the two MAGIC telescopes.

We have observed a clear correlation signal at different telescope baselines and four different observation nights. The significance of the signal is consistent with the parameters of our instrumental setup.

We have also estimated the time correlation function for all data samples and compared with estimates based on the star diameter and telescope baselines. They are both consistent within the systematics of our instrument.

These results prove that MAGIC can be used as an optical intensity interferometer. Even with this simple instrumental setup we expect to detect a significant correlation signal for an unresolved star of mB ∼5 in 3 h of observations.

Our setup is currently limited by the duty cycle of the acquisition system so our plan is to upgrade it to a system with a duty cycle in excess of 90 per cent using commercial digitizing cards that allow to transfer the data continuously to a computer. A higher duty cycle will generate a significant amount of data though. We aim at reducing the data online either at the digitizing cards or at the computer.

We also plan to improve our mechanical setup. For safety reasons, the current filter frame must be mounted and dismounted during the day so we are essentially limited to observing during Full Moon nights when VHE observations are impractical. We will build a new frame that can be deployed at any time. This allows not only to take additional data interspersed with regular VHE observations but also to test the response of the interferometer under very different light conditions.

On a longer time-scale, the interferometer can be improved in a number of ways. Since we digitize the signal using the standard photosensors already available in the VHE camera and the standard electronic path down to a central electronics room we can easily increase the number of channels. We are considering adding other spectral channels and two polarizations.

PMTs can also be easily upgraded because the cameras of the MAGIC telescopes are modular. Each individual ‘cluster’ of 7 PMTs is designed to be easily swapped with a new cluster for maintenance. Clusters equipped with new photosensors may be built as long as they are mechanically and electrically compatible: in fact clusters based on SiPMs have been already tested (Hahn et al. 2018; Guberman et al. 2017; Rando et al. 2015). SiPMs show higher photodetection efficiency. New photosensors can also be selected to improve the time resolution or to match the spectral response of new spectral channels.

ACKNOWLEDGEMENTS

This work would have been impossible without the support of our colleagues of the MAGIC collaboration. We would also like to thank the Instituto de Astrofísica de Canarias for the excellent working conditions at the ORM. The financial support of the German Bundesministerium fuer Bildung und Forschung and Max Planck Gesellschaft, the Italian Istituto Nazionale di Fisica Nucleare and Istituto Nazionale di Astrofisica, the Spanish Ministerio de Economia y Empresapartly using European Research Development Funds (FPA2017-82729-C6-6-R and C6-3-R) is gratefully acknowledged. This work was also supported by the Spanish Unidad de Excelencia ‘María de Maeztu’ (MDM-2014-0369) and by the 2018 Leonardo Grant for Researchers and Cultural Creators from Banco Bilbao Vizcaya Argentaria Foundation.

REFERENCES

Ahnen
M.
et al. .,
2017
,
Astropart. Phys.
,
94
,
29

Aleksić
J.
et al. .,
2016a
,
Astropart. Phys.
,
72
,
76

Aleksić
J.
et al. .,
2016b
,
Astropart. Phys.
,
72
,
61

Britzger
D.
,
2009
,
Studies of the Influence of Moonlight on Observations with the MAGIC Telescope
.
Diploma thesis, Technische Universitat München

Cheli
A.
et al. .,
2016
,
A&A
,
589
,
A112

Dravins
D.
,
LeBohec
S.
,
Jensen
H.
,
Nũnez
P. D.
,
2012
,
New Astron. Rev.
,
56
,
143

Dravins
D.
,
LeBohec
S.
,
Jensen
H.
,
Núñez
P. D.
,
CTA
Consortium
,
2013
,
Astropart. Phys.
,
43
,
331

Ducati
J. R.
,
2002
,
Catalogue of Stellar Photometry in Johnson’s 11-colour system, CDS/ADC Collection of Electronic Catalogues
.

Glindemann
A.
et al. .,
2000
, in
Léna
P.
,
Quirrenbach
A.
, eds,
Proc. SPIE Conf. Ser. Vol. 4006, Interferometry in Optical Astronomy
.
SPIE
,
Bellingham
, p.
2

Guberman
D.
, et al. .,
2017
, in
Proc. 35th International Cosmic Ray Conference, Busan, Korea
.
Proceedings of Science
Trieste (Italy)
,
p. 816

Guerin
W.
,
Rivet
J.-P.
,
Fouche
M.
,
Labeyrie
G.
,
Vernet
D.
,
Vakili
F.
,
Kaiser
R.
,
2018
,
MNRAS
,
480
,
245

Hahn
A.
,
Dettlaff
A.
,
Fink
D.
,
Mazin
D.
,
Mirzoyan
R.
,
Teshima
M.
,
2018
,
Nucl. Instr. Methods A
,
912
,
259

Hanbury-Brown
R.
,
1974
,
The Intensity Interferometer: Its Application to Astronomy
.
Taylor & Francis
,
London

Hanbury-Brown
R.
,
Twiss
R. Q.
,
1954
,
Phil. Mag.
,
45
,
663

Hanbury-Brown
R.
,
Twiss
R. Q.
,
1956
,
Nature
,
178
,
1046

Hillas
A. M.
,
2013
,
Astropart. Phys.
,
43
,
19

King
D. L.
,
1985
,
Atmospheric Extinction at the Roque de los Muchachos Observatory, La Palma
.
RGO/La Palma, Technical Note 31, Royal Greenwich Observatory
,
London (UK)

Matthews
N.
et al. .,
2019
, in
Proc. 36th International Cosmic Ray Conference, Madison, Wisconsin
.
Proceedings of Science
Trieste (Italy)
, p.
740

Oppenheim
A. V.
,
Schafer
R. W.
,
Buck
J. R.
,
1999
,
Discrete-Time Signal Processing
.
Prentice Hall
,
Upper Saddle River, NJ

Rando
R.
et al. .,
2015
, in
Proc. 34th International Cosmic Ray Conference, The Hague, Netherlands
.
Proceedings of Science
,
Trieste (Italy)
,
p. 940

ten Brummelaar
T. A.
et al. .,
2015
,
ApJ
,
628
,
453

This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)