Numerical and analytical considerations
The detection scheme consists of a quadrant photodiode (QPD) in the image plane of a Newtonian telescope. The QPD is focused at a known distance, r_{foc}, which coincides with the position of a black termination cavity. As such, the QPD is unfocussed at the aperture but is gradually focused along the FOV toward r_{foc}. Expressed differently, the FOV of the QPD segments overlap entirely at the telescope aperture, and are gradually separated with distance until imaged sharply at r_{foc}. In this setup, the FOV is assumed to be evenly illuminated by the sunlight and is therefore equivalent to the probe volume. The properties of the imaging system are illustrated in Fig. 1.
We further make the following assumptions:

1.
The observed insects are small compared to the probe volume (corresponding to point sources when illuminated by homogeneous sunlight).

2.
Throughout the duration of an observation, the insect velocity vector is constant.

3.
The coaxial movement of the observed insects is limited compared to the length of the probe volume.
When an insect transits the probe volume, sunlight is scattered into the telescope and an intensity increase is recorded. The signal includes a lowfrequency envelope from the insect body and an oscillatory component from the insect wings [30, 31]. As observed in Fig. 1d), the signal recorded in two adjacent detector elements is identical at close range, differs but blends together at medium range, and is sharply resolved at far range. The overlap can be quantified through timelag correlation of the detector segments, which is the crosscorrelation between the two time vectors obtained by the detector elements as a function of sliding delay. The delay time, τ, then corresponds to the maximum correlation. τ increases linearly with range, as seen in Fig. 1d), and is inversely proportional to the flight speed of insects. The flight speed is eliminated by dividing τ with the transit time Δt, yielding a quotient that is unique for each range. The distance \(\hat{r}\) to an insect transiting the probe volume is estimated according to Eq. 1.
$$\hat{r} = \frac{{\tau {\o}_{\text{tel}} f}}{{\tau \left( {\frac{{{\o}_{\text{tel}} f}}{{r_{\text{foc}} }}  d_{\text{s}} } \right) + \Delta t\frac{{d_{\text{s}} }}{2}}}.$$
(1)
In Eq. 1, ø_{tel} is the telescope aperture, f is the telescope focal length and d_{s} is the quadrant sensor width. The validity of Eq. 1 can be evaluated with the simulated insect signals in the raytracing model, where the distance is known. Equation 1.1 is obtained by rearranging Eq. 1. d_{foc} is the width of the sensor image in the focal plane as calculated with the lens equation (see Fig. 1c). α = ø_{tel}/d_{foc} is the quotient of the width of the FOV at the telescope location and the width of the FOV at termination, i.e., the inverse of the linear scaling coefficient of the FOV.
$$\frac{{\hat{r}}}{{r_{\text{foc}} }} = \frac{{\alpha \frac{2\tau }{\Delta t}}}{{\left( {\alpha  1} \right)\frac{2\tau }{\Delta t} + 1}}$$
(1.1)
The FOV converges if it is terminated in the nearfield and diverges if it is terminated in the farfield. If the FOV is terminated at the limit between the near and farfield, r_{lim}, it maintains a constant width along the monitored transect. In other words, if r_{foc} = r_{lim}, the quotient α is equal to 1. The distance to termination therefore has implications on the ranging properties; see Fig. 2. The setup parameters given in Fig. 1 yield a limit between near and farfield r_{lim} = 121.2 m.
Aside from r_{foc}, there are other setup and observation parameters that affect the ranging accuracy. Equation 1 includes the quotient τ/Δt and is therefore invariant to the flight speed of insects. τ is equal to 0 at the aperture and Δt/2 at termination (see Fig. 1d). Δt depends on the cruise altitude of insects at close range due to the round aperture of the telescope, but is invariant to cruise altitude at far range due to the quadratic shape of the sensor. The quotient τ/Δt is equal to 0 at the aperture and increases linearly with range to 1/2 at termination, indicating that the ranging equation as a whole is invariant to cruise altitude. The climb and heading angles with which insects transit the probe volume can affect the results of Eq. 1. However, this effect is minute since the insect transit distance along the optical axis is negligible in comparison to the length of the probe volume (assumption 3).
Field measurements and insect observation properties
Experiments were carried out on 20130720 between 12:00 and 15:20 at Stensoffa in Southern Sweden (N55.70, E13.45). A quadrant photodiode (Hamamatsu S4349 Si PIN) was used to monitor an air volume through a Newtonian telescope with a focal length of 120 cm and a diameter of 30 cm. The detector has a bandwidth of 2.9 kHz and was sampled at 20 kHz with a data acquisition board (NI USB6211). The telescope was aimed roughly north in order to collect insect data in backscatter mode. A black termination cavity was placed at a distance of 142.3 m from the telescope, where the detector was focused. The setup, measurement principle, and site are shown in Fig. 3. More details on the experiments are found in [32].
The sunlight propagates through the atmosphere before impinging on the FOV of the detector. As such, turbulence and other atmospheric phenomena can cause rapid variations in the optical background. By filtering the signal with a sliding median, the optical background was obtained. The noise level was obtained as the difference between the optical background and a sliding minimum filtered signal. A detection threshold with signaltonoise ratio SNR = 2 was set.
A set of harmonic basis functions ψ_{
n
}(t) was constructed based on test frequencies f_{0test}. This is described by Eq. 2. The backscattered signal intensity from an insect contains light scattered in the insect body and in the insect wings, I_{tot}(t) = I_{body}(t) + I_{wing}(t). These two signal components were separated using a sliding minimum filter with a window width equivalent to the period of f_{0test}. The estimated wingbeat signal \(\hat{I}_{\text{wing}}\)(t) was expressed as a linear combination of basis functions with coefficients a_{
n
} as expressed in Eq. 3. \(\hat{I}_{\text{wing}}\)(t) was fitted against the original wingbeat data I_{wing}(t). By scanning f_{0test} and minimizing the residuals of the fit, a final value for the wingbeat frequency f_{0} was obtained. This process is further detailed in [31].
$$\psi_{n} \left( t \right) = \sin \left( {2\pi nf_{{0{\text{init}}}} t} \right) + \cos \left( {2\pi nf_{{0{\text{init}}}} t} \right),\quad n = 0,1,2, \ldots$$
(2)
$$\hat{I}_{\text{wing}} (t) = \mathop \sum \limits_{n} a_{n} \psi_{n} ,\quad a_{n} \left( t \right) = \left( {\psi_{n}^{T} \psi_{n} } \right)^{  1} \psi_{n}^{T} I_{\text{wing}}$$
(3)
I_{tot} vectors can be obtained separately in the signals from the four detector elements. The east–west and up–down time delays, τ_{ew} and τ_{ud}, were obtained through timelag correlation of the corresponding I_{tot} vectors. \(\hat{r}\) was calculated according to Eq. 1 with τ obtained as a weighted average of τ_{ew} and τ_{ud}, see Eq. 4.
$$\tau = \frac{{\tau_{\text{ew}} \cdot I_{{{\text{body}}\_{\text{e}}}} \cdot I_{{{\text{body}}\_{\text{w}}}} + \tau_{\text{ud}} \cdot I_{{{\text{body}}\_{\text{u}}}} \cdot I_{{{\text{body}}\_{\text{d}}}} }}{{I_{{{\text{body}}\_{\text{e}}}} \cdot I_{{{\text{body}}\_{\text{w}}}} + I_{{{\text{body}}\_{\text{u}}}} \cdot I_{{{\text{body}}\_{\text{d}}}} }} ,$$
(4)