Mathieu Kohli
email: mathieu.kohli@protonmail.ch

Price structure of traded honey, pollen and wax linked to winter survival of a honey bee colony

la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

The interaction between Apis mellifera honeybees and beekeepers raises the issue of how they can mutually benecially cooperate? While bees struggle to be resilient through tough weather conditions, beekeepers need to sell honey, wax and pollen at a reasonable price. Here we investigate from a mathematical perspective how a resilient bee colony would price its own honey and other commodities depending on its internal state and on environmental factors.

We start by assuming the environment is a deterministic function of time. Then we allow randomness on the environment to understand how that changes the bees' perception of the price of honey.

If a beekeeper has this price information in mind, he only needs to add the price of his own labor to compute the price at which he can set his honey. Moreover, since this model focuses on an ideal resilient bee colony, it can be compared to real colonies to discriminate which ones have a high resilience behavior in the climate in which they live. Another potential application of the model is the prediction of swarming behavior.

Value of honey resources

The starting point of our study is the knowledge that a specic ower eld or bush or tree in the vicinity of a hive has a limited potential to feed the hive. Inasmuch as the owers only provide nectar to bees (as opposed to pollen), this potential can be dened as the net energy intake gained by the colony per time unit as it sends one forager bee back and forth between the hive and that resource. We call it the potential eciency of a forager bee on this nectar resource. It depends on parameters that are identied in [Kohli 2023, Kohli and Anstett]:

the greater the distance d separating the hive and the oral resource, the more energy is spent to transport nectar [Ausseuil et al. 2018] and the longer it takes to travel between these two points. To this, we must add the eect of wind, which is equivalent to a distortion of distances on the geographic map [Kohli and Anstett]. Distances can be compared to the maximal distance from the hive d max at which bees can be found foraging [START_REF] Beekman | Long-range foraging by the honey-bee, Apis mellifera L[END_REF], Greenleaf et al. 2007, Ratnieks and Shackelton 2015]. The owers can be highly or scarcely rewarding for bees. We measure this quality through a coecient q f that is dened as the energy contained in the stomach of a forager that is full with nectar from ower f , to which we substract the energy needed to transform this quantity of nectar into honey through evaporation of superuous water [Mitchell 2019]. Quantum q f can be compared to q 0 that corresponds to a reference ower which is the most rewarding ower possible. The ower species also inuences the average number of owers m f a forager must visit to ll its stomach with nectar, as well as the average time β f a bee spends on a ower to forage. The higher the density ρ f of owers, the less energy and time is lost collecting nectar within a small range. Density ρ f should be compared to a critical density ρ f,crit below which owers are too scattered to form a resource. Some parameters do not depend on the oral ressource but must be deduced for experiments. Namely V , the speed at which foragers travel between the hive and the owers [START_REF] Nachtigall | Flight of the honey bee VII: metabolic power versus ight speed relation[END_REF], as well as v hop , the average speed at which bees hop from one ower to another. The last relevant parameter that is included in the model we consider is the average time t hive a forager spends inside the hive between two trips to a oral resource [START_REF] Fewell | Individual foraging eort as a function of colony population in the honey bee, Apis mellifera L[END_REF], Colin et al. 2022]. Once the parameters that aect the potential eciency η f of a forager bee on a nectar resource have been stated, we can express it through the following formula [Kohli 2023, Kohli and Anstett]:

η f = q 0 q f q0 -d dmax - q f q0 ρ crit,f ρ f t hive + 2 d V + m f β f + 1 2v hop √ ρ f . (1.1)
If we assume the bee colony optimizes the use of its foragers on nectar resources, it should send as much foragers as possible to the resource where the foragers have the maximal forager eciency. If some foragers remain unemployed, they should be aected to the second best resource in terms of eciency of the bees and so on until no forager is left idle.

To implement the previous process, the knowledge of how many foragers ϕ f a oral resource can constantly feed is required. It can be derived from the same parameters that were used to express the eciency of foragers [Kohli 2023, Kohli and Anstett]:

ϕ f = Sρ f λ f q f t h + 2d V + m f β f + 1 2v hop √ ρ f , (1.2)
where S is the surface of the oral resource and λ f represents the power oered to foragers by one individual ower [START_REF] Baude | Historical nectar assessment reveals the fall and rise of oral resources in Britain[END_REF]].

Relative price of pollen compared to honey

To include pollen as a resource in the allocation process previously presented, there is one main problem. Unlike honey, pollen does not provide the energy the hive needs. It provides the matter (amino acids, essential fatty acids, vitamins,...) necessary to form bees [Di Pasquale et al. 2013].

Therefore, the net gain for the bee colony brought back by a forager bee is a quantity of matter minus an energy, which can not be directly reduced to one unity. Even if we are able to compare how eciently bees can forage on two distinct pollen resources, another diculty follows: how can we compare a nectar resource to a pollen resource? To simultaneously overcome both these challenges, we admit that each day during the foraging season, bees are able to assign a value to both nectar and pollen based on their current needs. The ratio between the value of high quality pollen and the value of honey is the relative price of pollen in honey.

To implement this in practice, we replace q f in formula 1.1 by

price × q0 qf q0
where q0 is the maximal mass of pollen that a forager can transport on its legs, price is the energetic equivalent of the maximal mass of pollen q0 provided this pollen is of the highest quality possible pollinic ower, in the sense that bees that eat this pollen live the longest (see [Di Pasquale et al. 2013]) qf q0 is the ratio between the average lifespan of a bee fed with pollen of ower f and the average lifespan of a bee fed with the reference ower [Di Pasquale et al. 2013]. To summarize, a relative price of pollen in honey corresponds to the distribution of foragers on the local oral resources, which in turn corresponds to a certain inow of nectar and pollen in the hive. Now the nectar and pollen inow can be linked to the variation of the reserves of both these commodities by substracting the losses to the inows. We write the balance sheet for energy as it is stated in [Kohli 2023, Kohli and Anstett]:

dE dt = µ δM f orag dt -α δN death dt -πN -h (T, N ) .
(2.1) Each term in the previous equation has a specic meaning:

δM f orag dt stands for the net quantity of honey produced with the nectar inow entering the hive through foraging. This term can be computed as

µ δM f orag dt = f exploited nectar resource η f ϕ f , (2.2)
where ϕ f is computed through 1.2 except for the resource where the foragers are the least ecient, where ϕ f simply is the number of remaining foragers that are not working on a better resource. δN death dt is the death rate in the bee colony. See [Di Pasquale et al. 2013] for a study of this death rate related to the pollen the colony feeds on. πN represents the "basis metabolism" inside the hive, where π stands for the average metabolic rate of a bee that is neither ying nor taking part in the thermoregulation of the hive. [START_REF] Kovac | Respiration of resting honeybees[END_REF] h represents the energy used to thermoregulate the hive [Stabentheiner 2021]. It is a function of the ouside temperature and the number of worker bees, but also of the number of males, the number of nurses and the size of the brood. A similar equation holds for pollen and the number of bees this pollen is transformed into. It has terms equivalent to those in Equation 2.1 except for the term h in Equation 2.1 that has no pollinic equivalent.

Winter survival rate and quantitative state of the colony

At the begining of winter, there are two major factors that depend on the work of the bees during the previous month and that account for the bees' ability to successfully reach the next spring alive.

Those factors are the quantity of honey enclosed in the hive and the number of bees surounding the queen.

Either through a thermodynamic simulation of energy consumption as a function of probable winter conditions outside the hive [Omholt and Lonvik 1986, Myerscough 1993, Ocko and Mahadevan 2014, Kohli 2023, Kohli and Anstett] or by collecting experimental data [Döcke et al. 2019], one can establish the probabilty P (t = f all winter, N, H) of surviving through winter if at the time of transition between fall and winter, there are N bees and H kilograms of honey inside the hive.

At times during the foraging season, in order to describe the quantitative state of the bee colony, it can be relevant to decompose N into several components that list the number of bees of dierent age groups and H into an actual honey component and a wax component (since bees process honey to obtain wax) as well as a component of heat accumulated in the form of hive temperature.

We now endeavour to compute, at any time t between the begining of spring and the end of fall, P (t, N, H, Envir). The previous quantity P (t, N, H, Envir) is dened as the probability to survive the following winter if on day t, there are H kilograms of honey and N bees in the hive, if all the environmental parameters from spring to fall which include ower ressources and weather are known in advance.

The function

Envir : [ spring ; fall ] -→ { weather } × { oral maps }
contains all environmental information during the foraging season including future times after t.

The complete knowledge of the environment between time t and the begining of winter implies determinism in the optimal choices of the colony between time t and the begining of winter, so that all the randomness comes from winter. And at any time before winter, we can write d dt P (t, N (t), H(t), Envir) = 0.

(3.1)

This translates into

N ′ (t) ∂P ∂N + H ′ (t) ∂P ∂H = - ∂P ∂t . (3.2)
In the previous equation, term ∂P ∂t corresponds to the loss of probability of survival in winter that would arise if the hive was somehow "put onto pause mode" between time t and time t + dt. The optimization of the survival rate in winter requires to maximize -∂P ∂t through an optimal choice of foraged resource, which is equivalent to an optimal choice of price:

∂ ∂price ∂P ∂t = 0, (3.3) so by Equation 3.2, ∂ ∂price H ′ (t) N ′ (t) , --→ grad (P) = 0, (3.4)
and

∂ ∂price H ′ (t) N ′ (t) , --→ grad (P) = 0. (3.5) Although v Envir (t) def = H ′ (t) N ′ (t)
depends on the foraging process and also on the internal expenditures of the hive, ∂v Envir ∂price (t) only depends on foraging choices.

The problem of computing P for any time t during the foraging season can now be solved. Since we assumed that P is known at the begining of winter, we only need to deduce P(t -dt) from P(t). We start by nding out which price satises 3.5, which allows us to discover H ′ (t) and N ′ (t). Then we use 3.2 that leads us to P(t -dt).

Randomness in the environment

Of course, the future is not xed and there are dierent scenarios Envir that can play out with dierent probabilities. We denote by S the probability on the set of environmental conditions.

To understand what can happen when there are dierent paths leading from time t to the begining of winter, let us assume that there are just two dierent functions Envir1 and Envir2 that describe the same local conditions up to time t r but diverge at time t r . The hive must bet on Envir1 or Envir2 and follow the optimal path in the state space (H, N ) up to time t r + dt when it will be able to reassess if it did the correct choice or not. But if there is not only one but several hives in the same state (H, N ), some will choose Envir1 and others will choose Envir2, which means that the state (H, N ) will split into two states that increasingly diverge until time t r +dt. To know in which proportions both choices are made, we refer to Appendix A. Now at time t r + dt, all hives realize which choice was right, for example Envir1, and both splited states follow the optimal path in the state space (H, N ) knowing they are in Envir1 from there on. This leads to a splited state at the begining of winter.

By knowing the proportion of hives in each state and the probability to survive winter in each case, we deduce the global probability to survive winter along the environmental conditions Envir1, ramied at time t r if at time t we are in the state (H, N ) that we denote by P ramif (t, N, H, Envir1).

The previous notation stands for any ramication structure, even if it is more complex than just one ramication at time t r . To replace a simple state in the (H, N )-space by any splited state, we introduce probability densities on the (H, N )-space. A simple state then becomes a Dirac density. Now at time t, the actual probability to survive winter can be written as

T(t, N, H) = P ramif (t, N, H, Envir)dS(Envir),
where it is understood that dS is conditioned to the knowledge of what happend in the season before time t.

Absolute prices

Now that the meaning of a relative price of pollen in honey has been claried, we are ready to dene an absolute price of both honey and pollen in a currency such as dollars, euros,... Indeed, let us assume a beekeeper has only one hive and the market price of a bee swarm is for example 200 euros. This is the price a beekeeper would have to pay to maintain his activity if he is certain that his bee colony will not survive winter. Conversely, if the beekeeper is certain that his hive will live through winter, that saves him 200 euros. This means that the value of T = 1 is 200 euros.

This means that the price of a kilogram of honey is given by α ∂T ∂H , where we took for example α = 200 euros. This price in turn translates into an equivalence between currency and pollen through the relative price of pollen in honey.

Let us notice that wax prices can also be computed. The production of wax does not compete with the production of honey or pollen. Indeed, it is when they are between 5 and 20 days old that bees produce wax and build wax structures, while they usually start foraging after they are three weeks old [START_REF] Clement | Le traité rustica de l'apiculture, chapitre 2 : la vie sociale de la colonie[END_REF]]. The production of 1g of wax requires around 8.4g of honey according to [START_REF] Rodney | Dietary requirements of individual nectar foragers, and colony-level pollen and nectar consumption: a review to support pesticide exposure assessment for honey bees[END_REF] and since wax production does not interfere with any other production, the price of wax is simply 8.4 as much as the price of honey.

The best interest of both the beekeeper and the bees is that the honey, wax or pollen harvests happen at the moment when these commodities are at their lowest price. The beekeeper will add the cost of his own labour to this price when he sells his honey.

Swarming conditions

For a hive to form a swarm, the probability that both halfs of the initial colony die before the begining of the next spring should be inferior to the probability that the initial colony dies before that same moment. This translates into

1 -T t, N 2 , H 2 2 ⩽ 1 -T(t, N, H). (6.1)
Here, we did not take into account the dissymmetry between the old queen that leaves with the danger of wandering and the new queen that stays inside the hive with all the wax structure already built.

Let us notice that the beekeeper has one major advantage over the bee colony when it comes to computing T: he has access to meteorological predictions for the next few days. At a pivotal moment of the year such as swarming period, that means that he can anticipate swarming behavior.

Hornets

Bees are to hornets what pollen is to bees: food to feed their larva. From that perspective, we can apply the pricing results that apply to honey bees to hornets similarily. This can oer insight to nd out if hornet traps near beehives makes hunting bees an unsustainable activity or not. The value of lost hornets must be considered pro rata temporis of their remaining life expectancy and should be compared to the number of hornet larvae that can be fed with the bees that get captured around the hive.

Appendix A Betting on two opposite teams simultaneously

Let us assume a player is given the choice to bet between two outcomes.

Gain if

Gain if Team 1 wins Team 2 wins Player bets on Team . This strategy is the most resilient one, as it brings at its maximum the amount that is won with certainty.

If it so happens that the amount that should be bet on one team is negative, which is the case if and only if a < c or d < b, this leads to betting entirely on the opposite team. This corresponds to the "preparing for the worst, hoping for the best" strategy.

In the case we are studying, we translate this using the dictionary a P ramif (t, N, H, Envir1 dt (v Envir2 -v Envir1 ) P ramif (t, N, H, Envir2)

Here, dt denotes the time it takes for the hive to realize which of the two environmental scenarios actually is taking place. In the case where c < a and b < d, if dt tends to zero, the gain that can be assured at the ramication point is the translation of ad-bc (a-c)+(d-b) : P r (2) (v 1 -v 2 ) P r (1) -P r (1) (v 1 -v 2 ) P r (2) (v 1 -v 2 ) (P r (1) -P r (2)) .

(A1)

On the other hand, if a < c or d < b, randomness must still be taken into account and the hive bets on the worst case scenario.

  , we impose b < a, c < d. If instead of betting all on one team, the player decides to bet d-b (a-c)+(d-b) on Team 1 and a-c (a-c)+(d-b) on Team 2 then the player knows he will always win ad-bc (a-c)+(d-b)

  Envir1 -v Envir2 ) P ramif (t, N, H, Envir1) d-c


Jean-Pierre Rivet
email: jean-pierre.rivet@oca.eu

Farrokh Vakili

Olivier Lai

David Vernet

Mathilde Fouché

William Guerin

Guillaume Labeyrie

Robin Kaiser

Optical long baseline intensity interferometry : prospects for stellar physics

Keywords: Temporal and spatial photon bunching, micro-arc-second interferometry, in the optical wavelengths

More than sixty years after the first intensity correlation experiments by Hanbury Brown and Twiss, there is renewed interest for intensity interferometry techniques for high angular resolution studies of celestial sources. We report on a successful attempt to measure the bunching peak in the intensity correlation function for bright stellar sources with 1 meter telescopes (I2C project). We propose further improvements of our preliminary experiments of spatial interferometry between two 1 m telescopes, and discuss the possibility to export our method to existing large arrays of telescopes.

1 Introduction: historical background and revival Like classical amplitude interferometry, intensity interferometry is a highresolution observing technique used in astronomy to extend the angular resolution far beyond the classical resolution limit of single aperture instruments (optical or radio telescopes).

The principle underlying intensity interferometry is to monitor the space and/or time correlations of the electromagnetic wave intensity fluctuations in a given, limited bandwidth.

This idea was first proposed by Robert Hanbury Brown and theorized by Richard Twiss in the early 1950', in order to estimate the sizes of radio-sources beyond the resolution capabilities of classical amplitude radio-interferometry.

This method was first tested on-sky by Hanbury Brown et al (1952) to determine the size of "Cygnus A" and "Cassiopea A" radio-sources. Independently, similar diameter values were obtained for these objects by classical amplitude interferometry by [START_REF] Mills | Apparent angular sizes of discrete radio sources: Observations at Sydney[END_REF] and [START_REF] Smith | Apparent angular sizes of discrete radio sources: Observations at Cambridge[END_REF].

In order to extend their concept to the optical domain, Hanbury Brown and Twiss (HBT hereafter) first performed successful laboratory experiments on an artificial light source (Mercury arc lamp) with two photomultiplier tubes (Hanbury Brown and Twiss, 1956a). They demonstrated further the efficiency of their method in real, less than ideal atmospheric conditions, by measuring the diameter of the star Sirius with a prototype intensity interferometer involving two 1.56 m searchlight projectors used as photons collectors (Hanbury Brown and Twiss, 1956b). This first preliminary experiment paved the way for a more ambitious experiment : the Narrabri stellar intensity interferometer (NSII). The experimental setup involved two 6.5 m reflectors on a circular track with diameter 188 m. This configuration allowed adjusting the interferometric baselines up to 188 m, whilst maintaining the baseline projection constant (Hanbury Brown et al, 1967a,b;Hanbury Brown, 1968). From 1964 to 1972, the NSII has been used to perform ambitious experiments and led to the first consistent catalog of diameter measurements for 32 stars from spectral type O5 to F8 (Hanbury Brown et al, 1974a).

Classical amplitude interferometry requires to control the time-dependent sidereal optical path difference (OPD) between the arms of the interferometer with accuracies better than the wavelength of the radiation under consideration, on sub-millisecond time scales. This constraint is not considered as a challenge anymore, at least for long wavelengths (radio-interferometers) or in the optical domain for baselines of at most a few hundred meters, like for VLTI, CHARA and NPOI astronomical interferometers (Gomes et al, 2017;[START_REF] Mourard | Spectral and spatial imaging of the Be+sdO binary φ Persei[END_REF]Garcia et al, 2016). In the latter case, high precision optical delay lines are required to transport coherently the light signal from individual telescopes to the focal plane of the interferometer. For amplitude interferometry with very long baselines at short wavelengths (optical domain), the OPD compensation is still a challenging issue.

On the contrary, intensity interferometry only requires OPD control with accuracies given by the bandwidth of the bandpass filter (equivalently, by the time resolution of the acquisition chain), rather than by the wavelength of the radiation itself. For example, an intensity acquisition chain with 100 ps timeresolution requires an OPD control with millimetric accuracy only. Moreover, the link between individual telescopes only involves signal cables and no optical delay lines. In addition, the role of the telescopes is only to collect photons onto the sensitive surfaces of the sensors. This requires correct guiding and if possible fast tip-tilt correction, but the optical quality of the telescopes can be less-than-ideal. Besides, as emphasized by [START_REF] Tan | Optical intensity interferometry through atmospheric turbulence[END_REF], operation is still possible in real conditions, that is, through atmospheric turbulence and even with some sky background contamination (full Moon).

Intensity interferometry is thus much easier to implement, especially with long baselines and is less sensitive to optical defects and atmospheric turbulence. However, it suffers from a lack of sensitivity [START_REF] Brown | An attempt to detect a corona around beta Orionis with an intensity interferometer using linearly polarized light[END_REF], requiring large collecting areas. For this reason, intensity interferometry has not been developed further since the early 1970', until the last decade where a significant jump was performed by the technology of light detectors and signal digital correlators. Single photon avalanche diodes (SPAD) are now commercially available in the visible domain, with unprecedented sensitivity, short dead time and high time resolution. Moreover, time to digital converters (TDC) which can process the electric pulses delivered by SPADs with count rates as high as several Mcps (millions of counts per second) and time resolution of a few picoseconds are now also available. A review of commonly used SPADs and TDCs can be found in table 1 from [START_REF] Pilyavsky | Single-photon intensity interferometry (SPIIFy): utilizing available telescopes[END_REF]. These recent technological advances shed a new light on intensity interferometry for astronomy. For this reason, and despite its lack of sensitivity, intensity interferometry is now considered as a possible alternative to standard amplitude interferometry for very long baselines (kilometric and more) in the visible and near infrared domain.

This solution has been proposed as an alternative observing mode for the future Cherenkov Telescope Array [START_REF] Lebohec | Optical intensity interferometry with atomospheric Cherenkov telescope array[END_REF][START_REF] Nuñez | High angular resolution imaging with stellar intensity interferometry using air Cherenkov telescope array[END_REF]Dravins et al, 2013). With several tens of 4 m to 23 m telescopes spread over kilometric-sized areas, this facility will offer an unprecedented photon collecting area with sub-milliarcsecond angular resolutions and dense (u, v) plane coverage. However, these telescopes are not designed to be diffraction-limited and have a very large aperture ratios. This is likely to introduce technical difficulties to couple the telescopes with the photon detection devices.

Another approach to revive intensity interferometry methods is to implement it on existing diffraction-limited optical telescopes. This is the strategy followed for example by the "Asiago group" (Capraro et al, 2010;[START_REF] Zampieri | Intensity interferometry with Aqueye+ and Iqueye in Asiago[END_REF]. They have implemented intensity interferometry on the Copernicus and Galileo telescopes in Asiago (1.82 m and 1.22 m respectively), coupled with the cutting edge high velocity and high efficiency photon-counting photometers "Aqueye+" and "Iqueye" [START_REF] Naletto | Iqueye, a single photon-counting photometer applied to the ESO new technology telescope[END_REF].

In 2017, Pilyavsky et al. conducted extensive numerical simulations of intensity interferometry on seven 0.9 m to 4.0 m telescopes at Kitt Peak observatory [START_REF] Pilyavsky | Single-photon intensity interferometry (SPIIFy): utilizing available telescopes[END_REF]. Their simulations suggest that valuable science throughput on stellar diameters of hot stars can be achieved with commercially available components and without excessive technical complexity.

This renewal of interest around intensity interferometry has led our research team to conduct new experiments with state-of-the-art avalanche photodiodes and time-to-digital converters. In 2016, we had performed laboratory experiments with artificial light sources (Dussaux et al, 2016), to assess the capability of our setup to detect the photon bunching phenomenon. In this article, we describe further experiments we have conducted in 2017 on bright stars with one, then with two 1 m telescopes. We also discuss further more ambitious experiments that we plan to perform in a near future. We finally present possible implementations of these techniques on existing arrays of telescopes.

2 From amplitude to intensity interferometry on the French Riviera Following the pioneer optical long baseline interferometers on the French Riviera [START_REF] Labeyrie | Interference fringes obtained on Vega with two optical telescopes[END_REF][START_REF] Labeyrie | Fringes obtained with the large 'boules' interferometer at CERGA[END_REF]) and given our currently operating focal instruments (VEGA-CHARA at Mount Wilson [START_REF] Mourard | Spectral and spatial imaging of the Be+sdO binary φ Persei[END_REF] and AMBER at Paranal ESO-VLTI [START_REF] Petrov | AMBER, the near-infrared spectro-interferometric three-telescope VLTI instrument[END_REF] or MATISSE [START_REF] Lopez | MATISSE status report and science forecast[END_REF], recently commisioned), we have initiated pilot experiments in 2015, aimed at attempting intensity interferometry in the visible and the near infrared using the two 1 meter optical telescopes of the C2PU facility (https://www.oca.eu/index.php/fr/accueil-c2pu) in view of implementing this technique on future arrays of optical telescopes. The C2PU facility is located at the Calern observing site of the Observatoire de la Côte d'Azur (hereafter OCA). We call this project "I2C" for Intensity Interferometry at Calern. I2C follows a stepwise approach to validate the focal instrument concept and components, including analysis and interpretation methods, and to reliably estimate the achievable sensitivity.

Like other groups, we have considered application to Cerenkov arrays which can be used for intensity interferometry during bright time, especially CTA for its photon-collecting surface which is almost equivalent to a 100 m dish with a maximum baseline of 2500 m, once the array is completed. However the angular spread of a Cerenkov telescope focal image is large, so that the light can only poorly be coupled to a single multimode fiber. Also, the sky background photon noise, integrated over the entire Cerenkov telescope angular spread, will affect the limiting magnitude. It is possible to introduce optical components to correct for the Cerenkov telescope primary mirror aberrations, but this dramatically complicates the design and complexity of these arrays primarily dedicated to high energy astrophysics. Another possible approach would be to use large surface detectors such as photomultipliers [START_REF] Matthews | Development of a digital astronomical intensity interferometer : laboratory results with thermal light[END_REF], at the expense of a lower electronic bandpass, therefore lower intrinsic setup visibility.

Our technical choice is based on optical fibers connected to SPADs rather than photomultipliers. So, there are advantages to using existing optical telescopes to implement intensity interferometry techniques. Optical telescopes, whatever their size, generally provide seeing limited focal images, so that stellar photons can be transported by optical fibers to photon-counting detectors in order to estimate the intensity correlation function g 2 (τ, ρ) defined by :

g (2) (τ, ρ) = < I(t, r)I(t + τ, r + ρ) > < I(t, r) > 2 . ( 1
)
With seeing-limited telescopes, at least 50% coupling efficiencies are achievable with multimode optical fibers, which are cheap and commercially available.

However, the very narrow bandpass required to preserve coherence of the photons implies that a large fraction of the light is rejected and not used in the correlation measurements. Therefore as already noted by Hanbury [START_REF] Brown | An attempt to detect a corona around beta Orionis with an intensity interferometer using linearly polarized light[END_REF], the most promising way to significantly improve the instrumental sensitivity of intensity interferometry is to observe spectrally dispersed light from the source with multiple photodetectors simultaneously at different wavelengths i.e., multiple spectral channels. Besides accessing valuable spectral lines to probe stellar surfaces and circumstellar structures, triggered by magnetic fields, mass loss or colliding winds, there is also gain on the intrinsic instrumental contrast, and hence the signal-to-noise ratio (SNR), which improves as the square root of the spectral channels simultaneously correlated [START_REF] Trippe | Optical multi-channel intensity interferometry -or : How to resolve O-stars in the Magellanic clouds[END_REF].

To keep the volume and design of a high dispersion spectrograph reasonable, it is much preferable to use single mode (SM) fibers, which allow the use of fibered Fabry-Pérot etalons, fibered Bragg gratings or integrated (on-chip) spectrograph technology, to work on a large number of simultaneous narrow band spectral channels. To couple the light into SM fibers, the telescopes must be equipped with adaptive optics; however, the optical throughput and the decreasing Strehl ratio at short wavelengths can negatively affect the sensitivity, and a trade study is required to compare different SM coupling technologies, such as photonics lanterns. These can be used to split the light into a number of SM fibers, which can then each be dispersed spectrally. The output of all the SM fibers can then be co-added before or after detection of each spectral channel.

Results from temporal bunching at C2PU

On the nights of February 20, 21, and 22, 2017, we have measured (Guerin et al, 2017) what we believe to be the first intensity correlation measured with the light of a star other than the sun [START_REF] Tan | Measuring temporal photon bunching in backbody radiation[END_REF] since HBT historical experiments. The temporal photon bunching g (2) (τ, r = 0), obtained in the photon counting regime, was measured for 3 bright stars, α Boo, α CMi, and β Gem. The light was collected at the focal plane of the West telescope of the C2PU facility (IAU observatory code : 010, latitude : 43 • 45 ′ 13 ′′ N, longitude : 06 • 55 ′ 22 ′′ E, altitude : 1270 m). The instrument is a F/12.5 pure Cassegrain combination, but to improve the fiber coupling efficiency, we have reduced the focal ratio to F/5.6 by inserting a combination of commercially available focal reducers. In the resulting focal plane, we have placed our compact optical setup. It consists of a dichroic beam splitter, a CCD guiding camera, a combination of filters and a multi-mode fiber (see Fig. 1 in Guerin et al ( 2017)). The light reflected by the dichroic (λ < 650 nm) is sent to the guiding camera. The transmitted light goes through a linear polarizer, then through two cascaded narrow band dielectric filters centered on λ 0 = 780 nm. The resulting bandwidth is 1 nm with a throughput estimated as 61% for the central wavelength.

Table 1 Main data and results of temporal bunching on the three stars observed in February 2017 at C2PU. T is the total integration time, in hours and minutes. The flux F is the number of detected photon counts per second and per detector, averaged over the total integration time. The contrast C = g (2) (0) -1 is the value of the correlation at zero delay given by the amplitude of the Gaussian fit, its uncertainty is the 1σ confidence interval of the fit. The noise is the rms noise on the data evaluated in the wings of the correlation function, where The light is then transported by a 20 m-long multimode graded-index fiber (MMF) with a core diameter of 100 µm, which is connected to a 50:50 fibered splitter such that two photodetectors can be used to compute the autocorrelation function without being limited by the dead time of the detectors, which is much larger than the coherence times. The two output ports of the splitter are connected to SPADs via two other MMFs of 1 m for the first detection channel and 2 m for the second. In order to avoid any spurious correlation induced by electronic cross-talk, we have introduced an electronic delay between the two channels by using a 10 m shielded BNC cable for the second channel and 1 m for the first channel. The total delay (optical and electronic) between the two channels is t 0 ≃ 45 ns and is subtracted in the data processing. For each detected photon, the SPADs produce a 10 ns pulse, whose rising edge is detected and processed by a TDC with a time bin of 162 ps. The TDC is operated in the "Delay-Histogram" mode, which yields a good approximation of the intensity correlation function after proper normalization.

g (2) (τ ) ≃ 1. Star Sp. type R mag I mag T F (×10 6 ) C (×10 -3 ) noise (×10 -3 ) α Boo K0-III C -1.03 -1
For total exposure times of a few hours, we obtained contrast values around 2×10 -3 (see Table 1), in agreement with the theoretical expectation for chaotic sources, given the optical and electronic bandwidths of our setup. Figure 1 displays a typical correlation graph (g (2) (τ ) vs τ ) obtained on the star α Boo (Arcturus).

It is worth comparing the performances of our measurements with the Narrabri intensity interferometry results published by Hanbury Brown et al (1967b). For this sake, we can use the star α CMi, which has been observed by both instruments. For the Narrabri instrument, the contrast of the correlation measured with the shortest baseline is reported with a SNR of 11.5 for a 13.9 hours observing time. We obtained a SNR of 6.8 with only 4.6 hours observing time and with a collecting area approximately 85 times smaller (the Narrabri interferometer used two 6.5 m collectors, whilst we used a single 1 m telescope with 9.7% central obscuration). The improvement in performance is due to the larger electronic bandwidth, the better quantum efficiency of the detectors, the fact that there are no correlation losses due to the separation of the telescopes (the minimal baseline of the Narrabri interferometer was 9.5 m) and our working wavelength of 780 nm, which is more adapted to the spectral type (F5IV) of α CMi than Narrabri's wavelength (443 nm). Moreover, the τ (ns)

g (2) (τ) -4 -3 -2 -1 0 1 2 3 4 0.9990 0.9995 1 1.0005 1.0010 1.0015 1.0020 1.0025
Fig. 1 Temporal intensity correlation function g (2) (τ, r = 0) measured on α Boo (Arcturus) in February 2017. The Gaussian fit (dashed line) is used to estimate the contrast (reported in Table 1) and the FWHM of the bunching peak. The fit window is [-10, 10] ns. performances of the system used by Guerin et al ( 2017) can be improved in two ways. First, a faster TDC and/or the choice of the "photon time tagging" mode, associated to an appropriate data post-processing system, allow for improved measurements of the correlation function. Second, we can increase the efficiency of our telescope-to-fiber coupling strategy by introducing a fast tip-tilt correction stage or, preferably, a low-order adaptive optics.

Preliminary results from spatial bunching

The transition from the temporal bunching experiments of February 2017 (Guerin et al, 2017) to spatial bunching observations has been quite fast and the first successful observations were performed in October 2017. The presence in the same facility (C2PU) of two independent but nearly identical 1 m telescopes, separated by 15 m has made this evolution relatively easy. For this, the fibered splitter used to feed both SPADs from the flux of a single telescope was removed. Instead, two optical fibers coming from both telescopes were connected to the two SPADSs. So, our telescope-to-fiber coupling system has been duplicated (same guiding camera, same set of filters) and both replicas were set in the focal planes of the two telescopes.

Three bright stars (α Lyr, β Ori and α Aur) were observed during 4 nights, from 10 th to 13 th October 2017. The correlation peak emerged easily, provided a suitable delay be taken into account (connection wire lengths difference and time-dependent astronomical optical path difference). Fig. 2 depicts the correlation peak detected on α Lyr for a total of 11.1 hours integration time τ (ns) spread over 3 nights. The on-sky projected baseline ranged from 9.6 m to 14.3 m, with a time-averaged value of 11.9 m. The contrast was measured as (0.89 ± .09) × 10 -3 (1-σ uncertainty, i.e. SNR= 10). Taking into account the temporal resolution of the detection chain, this value is quite consistant with the known value of the uniform disk angular diameter of α Lyr (3.14 mas), taking into account a 11.9 m projected baseline. A more detailed study on other bright stars is under progress and will be presented in a forthcoming detailed paper.

g (2) (τ)

The I2C project : future developments

Following the straightforward success of intensity interferometry techniques on the C2PU facility at OCA, we have initiated a short term plan to overcome its limitations and to improve its robustness for a further deployment on larger facilities in other observatories (see Sect. 3). The immediate plan is to improve our prototype experiment sub-system components to gain on the sensitivity (limiting coherent magnitude, i.e. source apparent brightness and effective visibility) whilst comparing its performances to operating optical direct interferometers like CHARA and/or NPOI through stellar sources with accurate angular diameters.

A first step is based on introducing polarization beam-splitters to observe stars with extended atmospheres and/or wind envelopes simultaneously in two polarizations, parallel and perpendicular to the interferometric baseline. HBT's early observations at Narrabri on Rigel (Hanbury Brown et al, 1974b) Fig. 3 Optical configuration of the I2C interferometer at OCA for stellar polarization and/or differential angular diameter measurements in the continuum versus emission lines such as Hα. T W and T E are respectively the West and East 1 m optical telescopes of the C2PU facility. M 1 and m 2 are the primary and secondary mirrors (both telescopes have nearly identical optical combinations). Seeing limited images at the Cassegrain secondary foci are split in their two linearly polarized components and injected into two separate multimode optical fibers (MMOF), which transport the corresponding lights to SPADs. All four SPADs are connected to four input ports of the correlator X. Here, each polarization from both telescope is correlated to its counterpart from the other telescope so that g (2) (τ, r) estimates are obtained simultaneously in polarized light parallel and perpendicular to the baseline. For C2PU the maximum East-West baseline reaches 15 m, corresponding to an angular resolution of 9.6 mas at 700 nm. Note that F are 1 nm bandwidth filters either centered on Hα line or at the 780 nm continuum to achieve differential visibility measurements.

have proven that substantial gains in accuracy on the visibility are required to measure any deviation from the uniform disk model, according to realistic stellar atmosphere models of early type stars (Cassinelli and Hoffman, 1975). Since then, direct interferometers have been unsuccessful to achieve such measures [START_REF] Vakili | Study of stellar polarization with the CERGA interferometer[END_REF] due to instrumental polarization effects (e.g. time-varying oblique reflections and cross-talk) introduced by the optical complexity of amplitude interferometers. Intensity interferometry is immune by construction to such effects since one can directly separate two polarizations and correlate them independently without any loss in SNR. Moreover, a gain of √ 2 is expected for unpolarized sources, by co-adding the correlation functions of both polarization states.

Fig. 3 depicts the generic configuration of the I2C interferometer project. This layout also allows to observe differential effects between two spectral channels (one in the continuum and one on some emission lines such as Hα). This is interesting for hot stars, for instance Be or LBV like γ Cas or P Cyg, which have been formerly observed in amplitude interferometry at OCA [START_REF] Vakili | Subtle structures in the wind of P Cygni[END_REF].

A second step consists in increasing the angular resolution by performing intensity correlation experiments involving one or two telescopes from C2PU and one telescope of the "MeO" laser-ranging facility [START_REF] Samain | Clock comparison based on laser ranging technologies[END_REF], located at 150 m of the C2PU building. Synchronization between these distant telescopes can be achieved via cable or fiber connection with adequate control or measurement of transmission delays.

Most modern TDCs can be used in two modes : cross-correlation computation and photon time-tagging mode. Both modes can be used for astronomical intensity interferometry. If the electric signals of two (or more) SPADs linked to two (or more) telescopes, are fed into two (or more) channels of the same TDC via coaxial cables, then the cross-correlation computation mode will deliver directly the photon arrival time correlation function with a time delay which needs to be computed and specified a priori. This operating mode requires a physical connection between telescopes through cables. On the contrary, if each telescope has its own SPAD, its own TDC (set in time-tagging mode), which delivers its time-tagging strings to its own computer, then no physical link is required between the telescopes. Correlations can be computed numerically afterwards. Of course, all TDCs need to be accurately synchronized to a common time reference (e.g. UTC). This configuration is more expensive but more straightforward to extrapolate to very long baselines, since no physical connection is required. Another advantage is that the correlation functions can be computed numerically a posteriori from the sets of photon arrival times, recorded independently by different TDCs and stored on independent mass storage devices. Contrary to correlations computed "on the fly" during the observation run, a posteriori correlations computation allows for more elaborate algorithms, for example with time-varying and/or adjustable parameters. This configuration with no physical link between telescopes could be tested as a third step of our experiment. In our case, the necessary time reference for the TDCs would be provided by the local time standard of the "MeO" instrument. This time standard offers a temporal resolution and stability on the order of a few picoseconds, which is significantly better than the coherence length of our intensity correlation measurements (of the order of 3 cm, because of the time-resolution of our SPADs).

Finally, we are also considering 3-telescope intensity correlation experiments between the two C2PU 1 m telescopes and the MeO laser-ranging 1.5 m telescope. The West-North-West direction of MeO with respect to C2PU offers much longer projected baselines than C2PU itself. They range between a few tens of meters and 150 m depending on the target's position. This leads to submas resolutions for visible wavelengths (e.g. in the range of 600-800 nm). Indeed, triple correlations between 3 telescopes would provide the closure phase quantity [START_REF] Malvimat | Intensity interferometry with more than two detectors ?[END_REF] but at the expense of a much higher exposure time for a given SNR. Several hundreds hours of integration would be necessary with the aforementioned OCA telescopes to attain the SNR value required to access any measure of phase information from the triple correlation quantity. However the power of intensity interferometry, despite its intrinsic low SNR and sensitivity, lies in the fact that the correlation signal can be integrated over many nights since no fringe detection or tracking is necessary as in amplitude interferometry.

Table 2 Potential arrays of optical telescopes with the highest photon-collecting surfaces and kilometric baselines suitable for modern intensity interferometry. N T is the number of telescopes with diameters larger than 1.5 m, S is the total collecting surface of the array in m 2 (assuming 10% to 15% losses due to the central obstruction), N b is the number of possible array baselines, Bm is the array theoretical maximum baseline in meters and ρ the potential angular resolution in mas at 700 nm wavelength. Last column gives the potential photon collecting area and/or the longest baseline becoming available if extremely large telescopes are built at or close to the corresponding observing sites.

Observatory

Prospects and discussion

The revival of intensity interferometry half a century after the pioneering work by Hanbury Brown and Twiss suggests that these techniques are promised to open new opportunities for very high angular resolution observations in the future, specially in the visible wavelengths, down to the violet if the site transparency permits. Two approaches have been proposed until now : the first one is based upon the use of existing or foreseen Cerenkov telescope arrays (Dravins and LeBohec, 2008;[START_REF] Matthews | Development of a digital astronomical intensity interferometer : laboratory results with thermal light[END_REF], and the second one is to use existing seeing limited optical telescopes (Zampieri et al, Pilyavsky et al, 2017;Guerin et al, 2017). Simulations and laboratory experiments confirm the imaging potential of HBT arrays whatever Cerenkov or optical. We started and effectively demonstrated that photon bunching, either temporal or spatial can be straightforwardly achieved with two modest 1 m telescopes at OCA using almost "off-the-shelf" quantum optics photonic components, at least if instrumental effects are correctly controlled and calibrated. As argued by [START_REF] Trippe | Optical multi-channel intensity interferometry -or : How to resolve O-stars in the Magellanic clouds[END_REF], multichannel correlations, with hundreds of simultaneous narrow spectral channels, is the most promising way to progress over the intrinsic limitations of intensity interferometry for a given limited number and photon collecting surface. Much experience is to be gained however by using existing 1-2 m class telescopes, coupled later to 8-10 m class telescopes for on-sky observations in order to control the detailed behavior of photonic subsystems and components of an intensity interferometer on various astrophysical sources before speculating on the observation of exotic or extragalactic sources with microarcsecond angular resolutions. Considering this, we propose to consider intensity interferometry on a few major observatories (see Table 2) and a first prototype HBT array based on the four 1.8 m ATs at Paranal/Chile (A. Mérand, private communication) that would be available for such experiments and observations when the VLTI focal instruments MA-TISSE and GRAVITY occupy the Paranal optical delay lines.

Fig. 2

 2 Fig.2Spatial intensity correlation function g (2) (τ, r) measured for α Lyr in October 2017. The Gaussian fit (dark solid line) is used to extract the contrast and the FWHM of the bunching peak. The fit window is[-10, 10] ns.

Acknowledgements The I2C pilot experiment is supported by INPHYNI and Lagrange laboratories, Döblin Federation and grants from OCA and the Excellence Initiative UCA-JEDI from University Côte d'Azur. We are grateful to A. Dussaux for his valuable contri-

bution to this project. We also thank E. Samain, C. Courde and J. Chabé (GeoAzur lab., OCA) for fruitful discussions about space and time metrology. Ph. Bério from Lagrange Laboratory (OCA) is also kindly acknowledged.