1. Introduction

The first-order correlation function, g ⁽¹⁾, is widely used in optical applications because a simple experimental arrangement makes it easily accessible, e.g. Young’s double-slit interferometer [1]. In Young’s experiment, two points on the wavefront separated by a distance x ₁₂ produce a fringe pattern with resultant intensity $I_1+I_2+2\sqrt{I_1{I}_2}\mid g^{\left(1\right)}({\mathbf{x}}_{12},\tau )\mid \cos \left(\Delta {\phi }_{12}\right)$ at a screen location where the phase difference is Δφ ₁₂ and the propagation time-difference is τ. The fringe visibility 𝓥 measures the magnitude of the first-order correlation function $\mid g^{\left(1\right)}({\mathbf{x}}_{12},\tau )\mid =\genfrac{}{}{0.1ex}{}{I_1+I_2}{2\sqrt{I_1{I}_2}}𝓥$ between field points at the double pinhole. When the light-state distribution, |g ⁽¹⁾(x₁₂,τ)|, is known then the phase difference is all that is required to define relationship between the field points. But in the case of unknown field distributions, measuring this quantity alone is ambiguous because it cannot distinguish between light states like: entangled photons, incoherent light, or coherent and thermal light (See Table 1). The ability to resolve such potential ambiguities is highly significant, for example, for the current debate over the existence of microcavity-polariton Bose-Einstein condensation (BEC).[2]–[10]

Table 1 shows that to distinguish coherent [11] and thermal (chaotic) light states [12] or incoherent light state and entangled photons (through this paper we referee to entangled two-photon state produced as a result of spontaneous parametric down conversion), the second-order correlation function $g^{\left(2\right)}({\mathbf{x}}_{12},\tau )=\genfrac{}{}{0.1ex}{}{〈I_1\left(t\right)I_2\left(t+\tau \right)〉}{〈I_1\left(t\right)〉〈I_2\left(t+\tau \right)〉},$, associated with Hanbury Brown and Twiss (HBT) correlated intensity fluctuations [13], must also be measured. Here, I ₁,₂(t) is the light intensity at a point $\pm \genfrac{}{}{0.1ex}{}{1}{2}{\mathbf{x}}_{12}$ and time t. However, recalling some of the early controversies surrounding photon-correlation measurements [14, 15, 16, 17], reminds us that second-order correlation measurements are technically more difficult to accomplish in the optical range and requires dedicated single-photon detectors as well as a coincidence-counting apparatus. An additional difficulty arises due to the impossibility of taking single, simultaneous image of the g ⁽²⁾(x, t,x, t) distribution. Instead, a pair of detectors and a beam splitter must be used to sample the distribution over a range of detector-pair positions. More recent correlation measurements, such as those for BEC cavity polaritons, often require the ability to resolve the spatial dependence of g ⁽²⁾ in a single run. Taken together, all these requirements demand an integrated, monolithic, photon detector (as flexible as a camera), that is capable of imaging intensity noise-correlations.

Elsewhere [18], we have presented such a device and applied it in a table-top implementation of a stellar HBT interferometer [19], which in their original paper was used to measure the apparent angular diameter of the star Sirius from the current-noise correlations of two detectors separated along a variable baseline. As seen from the Table 1, row for a thermal light source, the angular width of Sirius star θ=0.0063″ assumes that detector separation baseline x ₁₂ can be as large as a 10⁶×λ, yet enabling to observe second-order correlations.

However for a table-top implementation, for which the angular width of a typical source is θ~1000″ (e.g., a 200 µm spot seen from a distance of 2 cm), the detectors have to be placed at a distance x ₁₂ no more than a few tens of wavelength λ in order to observe correlations. It was thus not possible to measure g ⁽²⁾ in the optical spectral range with stand alone detectors located in one plane. A conventional table-top implementation necessitates the introduction of a beam splitter together with two detectors at equivalent planes related by reflection, as in the original HBT implementation [13]. At the same time for such delicate experiments demanding a clear-cut answer whether the intensity correlations exist or not, as in the polariton BEC experiments, a strong sensitivity to variations in the relative detector position might render the measurement procedure impractical and the results biased by an error in the relative detector position of the few tens of wavelength scale. Therefore, we proposed [18] the use of single photon

Table 1. Values of first and second order correlation functions for incoherent, coherent and thermal light states as well as entangled two-photon state produced as a result of spontaneous parametric down conversion [12]. Single-mode states are considered, θ is the angular width of the source, λ is the wavelength, τ_c is the coherence length. Integration effects due to limited detector response times and resolution of coincidence counter are not indicated.

View Table

avalanche detectors (SPADs) integrated into a monolithic array with a simple data-treatment algorithm for realizing table-top (laboratory) measurements of an HBT interferometer in the original “stellar”-like configuration.

 

Fig. 1. Micrograph of the 4×4 SPAD array (a), schematics of the SPAD pixel structure (b) and electronic readout circuit (c). For a typical operation conditions V _OP=-21V, which by 4 V exceeds the breakdown voltage, V _DD=3.3V and VBIAS=0V.

Download Full Size | PDF

Here, we report on the peculiarities of our SPAD detector-array for probing local g ⁽²⁾ correlations. The structure of this paper is as follows. In Sec. 2, we describe our detector system, which comprises a 4×4 array of Si SPADs implemented in CMOS technology (Fig. 1) and a simple algorithm to treat the multiphoton time-of-arrival distributions from different SPAD pairs, implemented off-chip. The device also incorporates on-chip high-bandwidth I/O circuitry which facilitates the external data processing. We discuss the data treatment procedure to correct for spurious correlations in measurements (Sec. 3) and we demonstrate the operation of our device for cases of very different light statistics (Sec. 4): measuring the beat note of a multimode coherent source, measuring the depth of intensity oscillations in a superposition of incoherent and thermal source, and measuring the instrument response function of a table top interferometer as a function of the source angular-width.

2. Structure of g⁽²⁾ imager

Our g ⁽²⁾ imager is based on a single photon avalanche diode (SPAD) array [20], whose photomicrograph is shown in Fig. 1(a). The array pitches are 30 and 43µm in the horizontal and vertical directions, respectively. The chip is fabricated in a 0.35 µm CMOS technology. Each pixel of 4×4 SPAD array is independently accessible and has a structure as indicated in Fig. 1(b). The SPADs p+-n-well junction is reverse biased above breakdown by a voltage known as excess bias voltage. In this mode of operation, called Geiger mode, single photons may be detected and counted using auxiliary electronics. The anode is enclosed in the p+-diffusion obtained from the standard implant used to define a PMOS channel. To prevent premature edge breakdown in the periphery of the multiplication region, a 3.5 µm guard ring of lightly doped p-type implant is designed to surround the anode.

The p⁺ anode of each SPAD pixel is kept at a negative voltage V _OP of 21 V, which by 4 V exceeds the breakdown voltage [Fig. 1(c)]. When an avalanche is triggered, e.g. as a result of photoabsorption in the depleted n-well region or thermally induced carrier release from a trap defect, the avalanche is subsequently quenched via a resistive current path in the PMOS transistor structure with the transistor gate kept at a ground potential (V _BIAS=0). This ballast resistivity brings the bias voltage across the p⁺-n junction below the breakdown voltage and is used to read out the photodetection events. After avalanche quenching followed by a recharge of the depletion region, the detection cycle is completed and can be initiated again by arriving photons or thermally excited carriers. Each SPAD pixel comprises a high-bandwidth inverter so as the Geiger pulses are converted into digital signals and transmitted for off-chip data processing. This approach considerably improves the signal-to-noise ratio, however as discussed in section 3, the data should be corrected for spurious correlations due to a possible crosstalk between SPAD pixel outputs.

During the time required to complete a detection cycle, the detector is insensitive to arriving photons. In our structures, the dead time is 12–15 ns, so as a theoretical limit on the count rate set by detector saturation (pile-up) effects is about 30 MHz.

Tunneling of minority carriers through the shallow p-doped guard ring and/or trap-assisted carrier release may trigger an avalanche and thus create spurious Geiger pulses. The lowest detectable photon flux is set by the dark count rate (DCR) of SPADs, which at room temperature conditions, is in the 5–10 Hz range. Such low DCR is achieved by implementing n-wells of small diameter d=3.5 µm. The lowest detectable photon flux density in our experiments is thus 10^-8 photons/s/cm² and the dynamic response range is 64 dB.

At 4V excess bias voltage, the measured photon detection probability (PDP) of such 0.35µm-CMOS SPAD pixels is 40 % at a wavelength 450 nm and exceeds 25 % in the wavelength range at 400–550 nm. PDP is proportional to excess bias. The spectral sensitivity curve and pixel statistics on dark count rate can be found in [20]. For our experiments, the chip was mounted on a standard ceramic case and was wire-bounded with 20 µm-thick aluminum wires of 6–7 mm length spaced by 100 µm. The bonding wires can be seen on the peripheries of the chip in Fig. 1(a). For what follows it is important to know that the bonding wires of individual detector data lines are ordered with the SPAD pixel number i indicated on the microphotographic image.

Since all 16 detectors in the array have separate parallel outputs, $\left(\genfrac{}{}{0ex}{}{16}{2}\right)=120\mathrm{HBT}$ HBT interferometer measurements are possible simultaneously between various detector pairs at temporal resolution set by the SPAD timing jitter characteristics (80ps). Measurements of pairwise second order correlation are based on implementing expression [21]:

with a simple algorithm. Here, integers i, j (i≠j) enumerate detector pixels D_i and D_j in Fig. 1, â_{i, j} is the photon annihilation operator at two detectors, ρ̂ is the density operator for the field and detectors, and trace is taken over the detector and field states [22]. At this stage, the pairwise intensity noise correlations g⁽²⁾(x_ij,τ) are computed from multiphoton arrivals at different detector pairs of the array, using a programmable four-channel 6-GHz bandwidth digital oscilloscope (Wavemaster 8600A, LeCroy) and implementing a simple expression

where X_i and X_j are discrete random variables whose values 0 (no event) or 1 (photon detection) correspond to the binary data stream emanating from any pair of detectors D_i and D_j, respectively. The spatial lag x_ij is set by the separation of the detector pair within the SPAD array. Time-lag increments τ=lT are set by multiples of temporal resolution T (equivalent to the bin width in conventional methods), where NT is the width of the measurement window and M is the overall number of measurements series. Note that both addition and the bitwise AND operator (∧ in the numerator of Eq. (2)) could be implemented easily using on-chip digital electronic circuits.

As can be seen from Eq. (2), the concept of our g ⁽²⁾ imager [18] is drastically different to conventional approaches based on measuring start-stop times and computing a timing histogram of delayed photon arrivals at two detectors, followed by renormalization of the histogram according to a statistical hypothesis for the field (e.g., a hypothesis like g ⁽²⁾(τ=0)=1 or g ⁽²⁾(τ=±NT/2)=1). Our g ⁽²⁾-imager algorithm (Eq. (2)) does not suffer from missing detection events in case of subsequent photon arrivals at the same detector and operates correctly with multiphoton arrivals. Here, we will show that differently from the conventional approach, it operates correctly with modulated intensity signals at any count rates within the dynamical response range of SPADs and any width NT of the temporal window of interest.

An inherent limitation of our current prototype detector-array is the off-chip data treatment implemented with a general purpose digital oscilloscope that enables user-defined Matlab data processing subroutines. As a result of the relatively slow Matlab interpreter, a data trace of 5 µs length (N=5000 points at T=1ns resolution) requires 75 ms of computations to calculate the second-order correlation function in 200 points (in the interval -100ns<τ<100ns). To achieve a low standard variation of the processed g̃⁽²⁾ data, about M=10⁶ such data traces needs to be acquired and treated. The undesirable latency of this numerical processing step could be reducesd dramatically by integrating our robust algorithm into the chip.

3. Multiphoton correlation measurements

The performance of our g ⁽²⁾ imaging device (Fig. 1) was tested by measuring the distributions g ⁽²⁾(x_ij,τ) when the 4×4 SPAD array detector is uniformly illuminated with an incandescent light bulb (incoherent broadband light source). Imperfections in g ⁽²⁾ imaging with the SPAD array chip might be caused by several reasons including the concept of monolithic implementation itself. Figure 2 shows two time-resolved correlation patterns g̃⁽²⁾ _ij(τ) (left panels) and correlation maxima images g̃⁽²⁾ _ij (0) (right panels) acquired with respect to detectors D8 (top panels, j=8) and D9 (bottom panels, j=9) respectively. In the figure, the pairwise correlations g̃⁽²⁾ _i8(τ) (g̃⁽²⁾ _i9 (τ)) are ordered according to the detector index i in Fig. 1. For an ideal g ⁽²⁾-imager, one would expect to measure no correlations (g ⁽²⁾(x_ij,τ)=1) and a uniform correlation maxima map for all pairwise detector combinations. However, Fig. 2 reveals a deviation from the expected distribution.

 

Fig. 2. Second-order correlation functions g ⁽²⁾(x_ij,τ) of SPAD array pixels (a) and the maps of correlation maxima g ⁽²⁾(x_ij,0) [(b) and (c)] measured for incoherent light using detector D8 [(a) and (b)] or D9 (c) as a reference. The time resolution is 1 ns. The panels corresponding to individual SPAD pixels are arranged in the same order in which they appear in the imager pattern and corresponds to SPAD pixel position in Fig. 1. Green scale curves in (a) show autocorrelation curves g ⁽²⁾(0,τ) in linear (olive curve, left axis) and logarithmic (green curve, right axis) scales

Download Full Size | PDF

Measured second-order auto-correlation function of a detector g̃⁽²⁾ _ii(τ) shows a large coincidence peak g ⁽²⁾(0)≫1 at zero time lag [green scale curves in (a)]. In fact, measurements with one and the same detector at τ=0 do not represent statistical correlations between field states projected on two independent detectors in Eq. (1) [21, 22]. Instead, the statistical algorithm in Eq. (2) is applied to one and the same field state projection. The excess correlations at zero time-lag can be understood by noting that the bitwise AND-ing of the binary data-streams in Eq. (2) yields 〈X ^(m) _i(n)∧X ^(m) _i(n)〉=〈X ^(m) _i(n)〉=µT for the average detection rate µ and time window NT. The result is then ${\tilde{g}}_{ii}^{\left(2\right)}=\left(0\right)=\genfrac{}{}{0.1ex}{}{N\left(\mu NT\right)}{{\left(\mu NT\right)}^2}=1⁄\mu T.$. In Fig. 2, the correlation curves are acquired at resolution time T=1 ns and the average rate of detections 6 MHz, yielding an estimate g̃⁽²⁾ _jj(0)=166, in perfect agreement with the acquired data in Fig. 2 (olive curve, left axis).

In the vicinity of the self-correlation peak, the dead-time bands of g̃⁽²⁾ _ii (τ)=0 are due to a recharge of the detector depletion region after each Geiger pulse. Measurements in our SPAD array indicate the dead time τ_D of 12–15 ns (Fig. 2(a), green curve, right axis). At larger time lags |τ|>τ _D, the self-correlation function g̃⁽²⁾ _ii (τ) again shows an excess of delayed coincidences, exponentially decaying to g ⁽²⁾=1. [Note a logarithmic scale of the right axis in (a).] This time, i.e., the excess of delayed coincidences, is due to the thermal release of minority carriers from traps in the depletion region, that were captured during the preceding avalanche current pulse. A casuality relationship between afterpulses and their generating photodetection events explains the measured increase of delayed coincidences for the incoherent light source [23, 24]:

where p_A(t)dt is the probability to detect an afterpulse in the time interval (t, t +dt) and ${\epsilon =\int }_{\tau D}^{\infty }{P}_A\left(t\right)dt$ is the overall probability of afterpulsing. We stress that the impact of afterpulsing effects on auto-correlation function (3) varies with the incident photon flux and is particularly pronounced at small count rates µ. Measuring auto-correlation curve (3) for a larger time lag (not shown in the figures) shows the characteristic decay time of afterpulsing probability 40 ns with an integral effect ε=7% at 4 V excess above threshold. The afterpulsing probability ε is determined by the number of trap defects in the active (n-well) region material, the number of carriers involved in an avalanche and the recharge time. Its impact on (3) can be drastically reduced by increasing the deadtime of the detector, e.g. by using a smaller recharge current for the depletion region. We have also noticed a slight decrease of afterpulsing probability with excess of the biasing voltage above the avalanche breakdown. In any case, the auto-correlation measurements, as in (3), do not allow the second-order correlation statistics of a source to be accessed but can be used to measure the count rate, dead-time, and afterpulsing probability of SPAD detectors.

The results of pairwise measurements in Fig. 2 deviate slightly from the expected cross-correlation function g ⁽²⁾ _ij(τ)=1 as well. For i≠j, the axis scale is limited to the range 0.8<g ⁽²⁾ _ij<1.3. In each of the correlation maps in Figs. 2(b) and (c), one can distinguish discrepancies of two kinds: (i) There are a couple of SPAD pixel-pairs exhibiting a bunching of Geiger pulses at zero time-lag (red curves), while (ii) the rest of detector pairs exhibit anti-correlations at τ=0 (black curves).

Thus in Fig. 2(a), which was acquired using SPAD pixel D8 as a reference, detectors D7 and D9 exhibit correlation maxima g ⁽²⁾(0)~1.3 measured at temporal resolution T=1ns. At higher temporal resolution and shorter integration time (T=100ps), its height increases while the natural FWHM measures 320 ps in width (not shown in the figures). Measuring these spurious correlations at various excess voltages above the avalanche threshold and incident light intensities, we have not seen any impact. Such behavior attests for the electrical cross-talk origin of spurious correlations that occur at the SPAD output data lines and have no impact on the detection process. In particular, the detectors D7 and D9 have data line bonding pads on the CMOS chip adjacent to the wire of SPAD pixel D8. (The bonding wires of SPAD data outputs seen at the chip edges in Fig. 1 are ordered with the SPAD pixel index i.) These wire lines are subjected to high currents when the logical state changes between 0 (0 V) and 1 (3.3 V). The measured excess of correlations can thus be attributed to electrical cross-talk between the SPADS pixels with neighboring bonding wires. This is confirmed in measurements of pairwise correlations with a reference detector D9 [Fig. 2(c)] showing the excess of correlations for detector pairs D8–D9 and D9–D10. The patterns in Fig. 2 indicate that there is no direct optical crosstalk [25] between the adjacent SPADs within the array, which would appear in Figs. 2(b) and (c) as a cross pattern centered at the reference pixel.

All other SPAD pixel pairs in Fig. 2 with larger index difference (|i-j|≥2) have non-nearest neighbor bonding wires and show small anti-correlations g ⁽²⁾ _ij(0)~0.85. The anti-bunching dip of Geiger pulses is of 2 ns width followed by strongly damped oscillations of 4 ns period (they can be seen in Fig. 3(b)). This feature corresponds to reflections in the 90 cm length RF cables transmitting photodetection events data to the processing timing electronics (oscilloscope Wavemaster 8600A, LeCroy). In particular, we attribute this to the impedance mismatch between CMOS output inverters of the SPAD array chip and the 50 Ohm RF cable. Measurements of spurious anti-bunching dips at various illumination power and excess above threshold conditions have not revealed any variations in the dip height or width, attesting the non-detector origin of measured g̃⁽²⁾,^bg_ij(τ) bias.

Spurious correlations and afterpulsing effects have to be taken into account in the photon coincidence measurements utilizing our SPAD arrays. Starting from the results of Ref. [23] and spurious background correlations g̃⁽²⁾,^bg_ij(τ) between X_i and X_j data lines, one can conclude that for incoherent and coherent light sources as well as for thermal, entangled or single photon sources with coherence time smaller than the detector dead time τ_D, the measured correlation function (2) is

where g ⁽²⁾ _ij is the correlation function for photons that would be measured by ideal detectors. The second term takes a reduction of the photon bunching peak (anti-bunching dip) due to the spurious afterpulses of one detector uncorrelated with photon detection events and afterpulses of another detector. Note that spurious uncorrelated dark counts have the same effect (see also Eq. (5)), however for illumination intensity levels used in the experiments discussed here, their impact is negligible. The third term is the correction due to background correlation caused by data transmission lines from the SPAD chip to the timing electronics.

We shall stress that in agreement with Eq. (4), there were no variations observed in the impact of afterpulsing effects on the cross-correlation function with increasing (decreasing) incident light power. This feature allows the detector’s pairwise correlation measurements (2) to be used as an efficient tool for accessing the field statistics g ⁽²⁾(x_ij,τ). Unlike conventional detection methods based on start-stop timing histograms of delayed single photon arrivals, which are not capable to treat multiple photon arrivals at one detector during the start-stop interval measurements, our approach implements properly normalized multiphoton distribution. As such it is robust against missing detection events, the impact of Poisson-like distribution decay ∝exp(-µτ) (see Ref [18] for details) and intensity modulation. As opposed to standard methods, Eq. (2) permits arbitrary count rates and temporal window of interest and does not require a statistical hypothesis to normalize g ⁽²⁾.

4. Experimental results and discussion

The imager was tested by measuring statistical properties of several light states including multimode coherent state, near field of an extended intensity-modulated thermal light source and the far field of stationary thermal light source.

4.1. Multimode coherent state

As a model system for a multimode coherent state we use the emission of a He-Ne laser operating on a fundamental Gaussian mode of the cavity at 633 nm wavelength transition. The cavity length of L=21 cm assumes a large separation of the longitudinal modes c/2L=710 MHz, such that a small optical gain and narrow linewidth allows only two or three longitudinal modes to reach the lasing threshold (in function of detuning from the gain line center). No special provision for stabilization of the cavity length or lasing modes has been made. The relative phases of modes vary rapidly in time such that no stationary intensity modulation at the intermodal frequency can be observed in the emitted output radiation.

The output beam of the laser impinging the g ⁽²⁾ imaging detector was attenuated to reduce the count rate at detectors down to 2 MHz. The temporal resolution of the timing electronics was set to T=100 ps so as to observe dynamics of lasing modes making a cavity roundtrip in 1.4 ns.

 

Fig. 3. a) shows the correlation function g̃⁽²⁾ _5,9(τ) [Eq. (2)] measured by two detectors in the middle of array (detectors D5 and D9) at 30 µm baseline. The background component due to a spurious crosstalk at the same detector pair g̃⁽²⁾,^bg _5,9 (τ), which was measured with an incoherent broadband light source, is shown in (b). The corresponding correlation function of the multimode laser beam g(2)(x5,9,τ) calculated from Eq. (4) after corrections for background correlations and afterpulsing effects in the detectors is plotted in (c). It can be seen that this procedure is an efficient tool for removing spurious (anti-) correlations seen as a dip in both the measured correlations g̃⁽²⁾ _ij(τ) and reference background g̃⁽²⁾,^bg_ij(τ).

Download Full Size | PDF

According to the Table 1, for each of the modes of the laser, one measures g ⁽²⁾(τ)=1. However, as evidenced by the measured second order correlation function in Fig. 3(c), when these modes are simultaneously excited in the lasing spectrum, one can see the effective beat signal of intensity and phase noise correlations in the modes [30]. The numerical fit with a cosine function (green curve) reveals the noise beat period of 1.4 ns, that is the time of a roundtrip in the laser cavity. The amplitude of harmonic oscillations 0.2 indicates that the noise of different modes is just partially correlated. We stress that no such regular oscillations are observed in the intensity of the output laser beam.

4.2. Intensity-modulated extended thermal light source

As a model system for an extended quasi-monochromatic thermal (chaotic) light source, the green line of mercury (546 nm) from the emission spectrum of a low-pressure Hg-Ar discharge lamp (CAL2000, Ocean Optics) is used. The lamp is driven with an original power supply incorporating a 30-KHz AC source producing 1.8 KV pulses to initiate the discharge and 120V in a steady state operation conditions. The green line of mercury is filtered with a 10 nm band pass filter (FL543.5- 10, Thorlabs).

Such a quasimonochromatic chaotic light source exhibits the first order correlations |g ⁽¹⁾|>0 in the Young’s interference experiment (Fig. 4(a)). It is characterized by Gaussian distribution of photons in the Glauber P-representation (Plank distribution in n-representation) so as one should expect measuring second order correlations $g^{\left(2\right)}\left(x,\tau \right)=1+\genfrac{}{}{0.1ex}{}{1}{q}\mid g^{\left(1\right)}\left(x,\tau \right){\mid }^2$ with q being the number of photon modes (of equal intensities) reaching the detectors [11]. Here we report on measurements in the near field of the source such that the number of modes is large. (The angular width of the source θ~1 in Table 1). Therefore, for a detector pair with the baseline of the order of λ and more, that is practically in all cases, the second order correlation function should be close to one, indicating no second order coherence.

 

Fig. 4. Multimode coherent state emitted by a He-Ne laser at 632.8 nm wavelength. (a): Acquired second order coherence function g̃⁽²⁾ _5,9(τ) [Eq. (2)]. (b): Spurious correlations background g ⁽²⁾,^bg _5,9 (τ) measured with the help of an incandescent light bulb. (c): Second order correlation function of the field g ⁽²⁾(x_5,9,τ) corrected for spurious correlations.

Download Full Size | PDF

Eq. (4)]. The laser cavity length is L=21 cm with the cavity roundtrip time of 1.4 ns (the period of oscillation seen in (c)). The detectors D5 and D9 separated by x_5,9=30 µm baseline are used at a temporal resolution of timing electronics T=100ps.

 

Fig. 5. Extended quasi-monochromatic thermal light source (546nm line of mercury). (a): Young’s interference fringes indicate phase correlation |g ⁽¹⁾|>0. (b) and (c): Measured second order correlation function g ⁽²⁾ _5,9 in the near field of the source obtained after corrections (4) at temporal resolution 100ps (b) and 100ns (c). The oscillations are due to the AC power supply of the Hg-Ar discharge lamp.

Download Full Size | PDF

To measure the second order correlations in the near field of such chaotic light source, a 3 mm size domain of the lamp discharge region was imaged onto the SPAD array at magnification 0.44. The image created by a 150 mm focal distance lens was overilluminating the entire SPAD array such that we were able to test correlations at two pairs of detectors, one in the middle of the array (SPAD pixels D5 and D9) and one pair with detectors at the opposite array corners (pixels D0 and D15). These detector pairs with baselines of 30 µm and 158 µm measured correlations g ⁽²⁾(x,τ) at the points of the lamp discharge that are separated by 68 and 360 µm distance, respectively. We were expecting to see no correlations (g⁽²⁾≡1), since the correlation measurements were done between different emission points of the extended source. However, in both cases, we have observed excess of coincidence counts. Figure 4 details the correlation function g ⁽²⁾ _5,9(τ) measured at the detector pair in the array center. The correlations measured with SPAD pixels at the array corners exhibit the same correlation curve with excess of coincidences g ⁽²⁾ _0,15(0)=2.7 (not shown in the figure).

The difference between the usual statistical data analysis and our method, can be understood as follows. The distribution g ⁽²⁾(τ) is assumed to be an even function about the origin with boundary condition g ⁽²⁾(τ→±∞)~1. However, when applied to a given measurement window of interest that is at τ=±NT/2, the boundary condition g ⁽²⁾(±NT/2)=1, that is used in practice, can be invalid if the window is not large enough. Consider the high resolution data in Fig. 4(b) sampled at T=100ps. Those data, which may be regarded as a uniformly distributed histogram of two-photon arrival times, can be renormalized to meet the above boundary condition, viz., g ⁽²⁾(τ)=1 at the edges of the sample window. But because of the renormalization, g ⁽²⁾ will not have the required value of 2.7.

By limiting the temporal resolution to T=100 ns and increasing the width of the measurement window to NT=50 µs, the origin of the excessive correlations can be made visible [Fig. 4(c)]. The measured second order correlation function just reveals the fact of long-period intensity oscillations of the lamp at the double frequency of the AC power supply, at about 60 KHz. Because of these oscillations and the large angular width of the source, the photon bunching at τ=0 due to Planck’s statics [13] is not visible. The change in the width of the measurement window and the temporal resolution (in Figs. 4(b) and (c)) has no impact on the measured excess value g ⁽²⁾(0)=2.7, attesting for correct normalization of the correlation function when multi-photon arrivals are treated with the algorithm (2)–(4). Interestingly, after correction for the afterpulsing effects at the detectors [ε in Eq. (4)], the minima of g ⁽²⁾(τ) function are close to zero, in agreement with fluorescence lifetime considerations.

Modulation depth of the correlation function can be tailored in a mixed state. Figure 5 details second-order correlation measurements for a superposition of incoherent light (as in Fig. 2) and a quasimonochromatic chaotic light beam (as in Fig. 4) at various probabilities to detect photons in chaotic light state. The incandescent lamp is used as a source of incoherent emission, which is superimposed on mercury green line fluorescence from the intensity-modulated Hg-Ar lamp. The probability P(Hg) that a detected photon has been emitted by mercury atoms, is given by the intensity ratio I_Hg/I estimated from the count rates at the detectors. For example, P(Hg)=0.38 when I_Hg=6.8 kHz, I =|17.8kHz.

One can easily show that for a superposition of mutually uncorrelated photonic states, characterized by partial intensities I ^(α) and correlations g ⁽²⁾,^α (x_ij,τ), the second order coherence function $g^{\left(2\right)}({\mathbf{x}}_{\mathrm{ij}},\tau )=\genfrac{}{}{0.1ex}{}{〈I_i\left(t\right)I_j\left(t+\tau \right)〉}{〈I_i\left(t\right)〉〈I_j\left(t+\tau \right)〉}$ for the overall photon fluxes I_i,j=∑α I^(α _i,j impinging the detectors assumes the expansion :

This expression can be simplified in the case of equal count rates of the detectors, yielding the partial weights 〈I ^(α)〉²/〈I〉² for state contributions of g ⁽²⁾,^α (τ)-1 to the overall excess (or lack) of coincidences of the correlation function. For example, for the overall rate 〈I〉=(1+e)〈I(^α)〉 of photodetection events caused by a process of intensity 〈I ^(α)〉 and statistics g ^(2),α _ij(τ) accompanied by spurious (incoherent) afterpulses occurring with probability ε, Eq. (5) predicts a reduction of the bunching peak (or antibunching dip). The second order coherence function of the entire process is then 1+(g ^(2),α _ij(ei )-1)/(1+ε)², in agreement with Eq. (4). Along the same lines, the effect of spurious dark counts of detectors can be taken into account at low count rates 〈I _i,j〉.

Fig. 5.Intensity correlations measured for the superposition of incoherent light (from incandescent lamp) and intensity-modulated quasimonochromatic thermal light (green 546nm line of mercury from Hg-Ar lamp CAL2000) at various relative intensities. Temporal resolution is T=100ns. (a): Second order correlation function g ⁽²⁾(τ). (b): Correlation function maximum g ⁽²⁾(0) plotted as squared relative intensity of thermal light incident on detector pair.

For the experimental results reported in Fig. 5, the incoherent light constituent does not contribute to the excess of correlation (the correlation function g ⁽²⁾,^incoh(τ)≡1). The average intensities at the detectors (SPAD pixels D5 and D9) are equal so as the contribution from the chaotic light source of Fig. 4 is reduced as 1+(g ⁽²⁾,^Hg _5,9 (τ)-1)〈I_Hg〉²/〈I〉² [Fig. 5(a)]. The maximum values of the second order correlations g ⁽²⁾ _5,9(0)-1 linearly grows with the squared probability of detecting a thermal light photons 〈I_Hg〉²/〈I〉² [Fig. 5(b)].

4.3. Quasi-monochromatic thermal light source

So far we have not yet observed the HBT photon bunching in the emission from a quasimonochromatic chaotic light source. To make it visible, we replace the power source of the lamp and conduct measurements in the far field. In particular, we drive the Hg-Ar lamp bulb (the same as used in previous experiments) from a DC voltage source. This lamp bulb (CAL-2000-B, Ocean Optics) having a cold cathode and U-folded discharge and designed for operation in the AC regime, also operated well with a DC power supply (160V @ 15 mA). To start the discharge, we used the original AC power supply of the lamp, which was connected in parallel with a DC source via a filter (2H inductance) such that the AC supply was gradually turned off, whereas the DC source was gradually turned on.

 

Fig. 6. (a) Experimental setup of the table-top stellar HBT interferometer. (b) Correlations measured at various detector separation x_{i j} and distance L to the fiber end. The Fresnel parameter $FN=\genfrac{}{}{0.1ex}{}{w{x}_{ij}}{\lambda L}$ is indicated in the panels. (c) Measured (points) and modeled (curves) second-order correlations in function of detector separation for the model Eq. (6) (dashed blue curve) and Eq. (8) (solid red curve). The inset shows the corresponding near field distributions at the fiber end (indicated with same type). (d) shows the corresponding modelled $g^{{(2)}_{m}ax}$ (r) patterns (top line of panels) and far field patterns (bottom line of panels). (e) Imaged second-order correlation maxima along the row of the array. Its position in the g ⁽²⁾ plane is indicated in (d) with green solid lines. It is assumed that g ⁽²⁾=2 along the diagonal. The green line of Hg (546nm) is used. Temporal resolution.

Download Full Size | PDF

As in the experiment of Fig. 4, the light emitted by the lamp was transmitted through a 10nm-width bandpass filter, which keeps only the emission at the green line of mercury (546 nm). The inherent difficulty is related to the fact that the number q of uncorrelated photon modes contributing to the intensity interference at SPAD array detectors $g^{\left(2\right)}\left(x,\tau \right)=1+\genfrac{}{}{0.1ex}{}{1}{q}\mid g^{\left(1\right)}\left(x,\tau \right){\mid }^2$ has to be kept small. (The requirement q~1 can be also seen from Eq. (5) for a superposition of q modes of equal intensities.) For a light beam of transversal size w seen from detectors at a distance L in the ray cone of angle Θ=x ₁₂/L, the number of photon modes is $q=\int \genfrac{}{}{0.1ex}{}{d^2\mathbf{r}{d}^2\mathbf{k}}{{\left(2\pi \right)}^2}\sim \genfrac{}{}{0.1ex}{}{1}{64}{\left(kw\Theta \right)}^2=\genfrac{}{}{0.1ex}{}{1}{16}{\left(\pi \genfrac{}{}{0.1ex}{}{w{x}_{12}}{\lambda L}\right)}^2.$. For a SPAD array of a few 100×100 µm² size, which defines a characteristic baseline x12 of detectors, the requirement q~1 imposes the angular width of such extended source to be θ=w/L~10-2 or less (see also Table 1).

To compromise the brightness of the source and its angular width at the detectors, the light from the lamp bulb was injected into a multimode fiber of 1 m length and core diameter w~200 µm [Fig. 6(a)]. The other end of the fiber was used to illuminate the SPAD array in the far field zone of the fiber end, at a variable distance L (1–3 cm) from the fiber such that the whole 4×4 SPAD array is over illuminated. Such extended thermal light source is of the small angular width w/L=10-2 rad and exhibits first-order correlations g ⁽¹⁾, when the SPAD array is replaced by Young double-pinhole interferometer as in Fig. 4(a), which was taken in the same setup configuration but with the Hg-Ar lamp driven from the AC power supply.

In intensity interference measurements with a quasi monochromatic thermal light source, the coherence width and Doppler-broadened spectral line have the same impact as in the case of the Young double slit interferometer. The classical expression for second-order spatio-temporal correlations for a non-polarized single-mode chaotic light source are determined by the first-order correlations [11, 21, 26] with coherence time ${\tau }_c=2\sqrt{2\pi \ln 2}⁄\Delta \omega$ due to the inhomogeneous broadening Δω (FWHM) of the emission line

where the second term in the right hand side takes into account the decorrelation effects due to unpolarized light (the coefficient 1/2), the zero-delay degree of spatial coherence and the Gaussian profile of the delayed first-order correlation function for an inhomogeneously-broadened line.

Figure 6(b) shows several traces g ⁽²⁾(τ) measured at different distances to the fiber end at several detector pairs (the Fresnel like parameter $FN=\genfrac{}{}{0.1ex}{}{w{x}_{ij}}{\lambda L}$ is indicated for each trace). At small FN<0.5 (large distance to detectors and small baseline), the detector counts show a correlation peak. At FN=1, it disappears and then at FN~1.2, it reappears but with a smaller amplitude of correlation excess. At first sight, such behavior corresponds well to the expected oscillations of the sinc function in Eq. (6). However measuring correlations at various distances from the fiber end with several detector pairs and correcting the data for afterpulsing effect (Fig. 6(c), points), we find large discrepancies with expected dependence (Fig. 6(c), dashed blue curve). The difference is outside the instrumental error and the fast growth of g ⁽²)max at FN~1.2 cannot be compensated for with a numerical fit to the model in Eq. (6).

To interpret the experimental results in Fig. 6, one shall admit that g ⁽²⁾, being a function of two variables, can be tailored not only in the time domain, as exemplified in previous sections, but it can also be tailored via spatial modulation of the near field distribution. This fact is often omitted in experiments on photon bunching, leading to wrong data interpretation. In the same way as the time domain variations of intensity (intensity noise) renders invalid the exponential term of the famous expression (6), the sinc term in that expression is obsolete if the source has a nonuniform intensity distribution in the near field.

Using Glauber’ expression g ⁽²⁾ ₁₂=1+G ⁽¹⁾ ₁₂ G ⁽¹⁾ ₂₁/G ⁽¹⁾ ₁₁ G ⁽¹⁾ ₂₂ [11] with G ⁽¹⁾ _ij=〈E_i E*_j〉 for an extended source with near field intensity distribution I(x,y) and detector pair with baseline x_ij centered at x=0, after integration over the wavelet contributions $E_{1,2}\propto \genfrac{}{}{0.1ex}{}{1}{\lambda L}\mid E_0\left(x,y\right){\mid }^2\exp (-ik{x}_{12}x⁄L)dxdy$ from mutually incoherent elements dxdy of the extended source, one finds

which is nothing else than the Van Cittert-Zernike theorem applied for visibility of the second order correlations. Here, as in Eq. (6), the coefficient 1/2 takes the decorrelation effects due to unpolarized light into account. We shall stress that only the 1-D Fourier transform of the near field intensity in the direction of the detector baseline is involved. For a uniform intensity distribution, integration over the source area from -w/2 to w/2 yields expression (6).

Field distribution in multimode fiber highly depends on excitation conditions and uniformity of the fiber core material [27, 28, 29]. Due to cylindrical symmetry of the fiber, the intensity distribution has an intensity bump or a dip at the core center. To simplify analysis, we consider only the main radial harmonic of frequency Ω, assuming the approximation I(r)∝1+γ cos(Ωr):

where γ is the intensity modulation depth and we have used Fresnel parameter FN=^wx _ij/λL to shorten the expression. The solid red curve in Fig. 6(c) shows the result of the numerical fit with the model (8), reporting γ=-0.3 and spatial frequency Ω=2.8 π/w. The corresponding near field (shown in the inset) has a small intensity dip at the center. Note that the difference to uniform intensity distribution is remarkable only in the correlation maxima distribution g ⁽²⁾ _max(r) (Fig. 6(d), top line panels) while it is not visible at all in the far field distribution impinging the detectors (Fig. 6(d), bottom line panels).

Being limited by the number of acquisition channels, we were able to record simultaneous correlations between four independent detectors form one row of the SPAD array. In Fig. 6(e), g ⁽²⁾(x _ij,0) measured along the array row is plotted as a pairwise correlation map g ⁽²⁾(i, j) for measured (left panel) and calculated with the models (8) and (6) correlation maxima (two right panels). In this image map, the spatial oscillations of the coherence factor are clearly visible. The measured excess of correlations for the corner pixels on the secondary diagonal corresponds better to the model (8), confirming previous conclusion on impact of the spatial modulation of the light source.

5. Conclusion

We have presented a g ⁽²⁾-imager built with conventional CMOS technology, which is capable of measuring second-order spatio-temporal correlated photons. We have discussed several important aspects related to temporal and spatial modulation of the source intensity when conducting HBT correlation measurements with such imager. Such detectors will find various applications. This approach allows the functions g ⁽²ⁿ⁾ of other even orders to be implemented as well.

The application in Sec. 4.3 used a green line emission of mercury in a Hg-Ar discharge lamp as a reliable source of HBT photons pairs. The lamp was emitting around hundred µW of thermal light at the green-line transition, only small fraction of which can be coupled into a multimode fiber. Yet the photon flux µ emanating from the fiber was quite high, about 10¹¹ photons/sec. If not decorrelations due to the angular width of such source, one would need detectors and coincidence electronics with integration time T~10ps to detect coincidences at the highest rate R_c~2µ ² T of about 10¹¹ photon pairs/sec in the near field of the fiber.

In practice, because of the low collection efficiency of detectors in the far field of the fiber end, the count and coincidence rates are very low, about µ=100 kHz and R_c=10 Hz respectively, for integration time T=1 ns. Using detectors of large diameter d will not overcome the problem, yielding a reduction of measured excess of correlations (g ⁽²⁾(τ)-1) by a factor of ~w ² d ²/λ ² L ² [see Sec. 4.3], like long integration time T results in a decoherence factor of τ ² _c/(τ ² _e+T ²)~τ ² _c/T ² [31].

Although the alternative of building a quasi-thermal source by using a He-Ne laser beam impinging a rotating sintered disk offers higher intensity, larger coherence time, and hence much faster acquisition, we scrupulously avoided this type of source so as not to introduce the undesirable possibility of spurious correlations such that to observe a clear impact of detector imperfections, non-stationarity and non-uniformity of the field on measured second-order correlations.

It is interesting to stress the impossibility of detecting a normalized correlation function g ⁽²⁾(τ) for entangled biphotons (Table 1) [12]. A state of the art entangled photon source based on spontaneous parametric down conversion in periodically-poled nonlinear crystal shows conversion efficiency of 10-⁷ and is capable of producing bi-photons at rates of 10⁸ photon pairs/sec per miliwatt of pump power. The principal difference to a thermal light source illuminating detectors is that all generated bi-photon pairs produce coincidence detections, provided perfect collection efficiency and no decoherence due to the angular width of the source. Technical difficulties arise from low overall power and a very short coherence time of bi-photons (typically, τ _c~0.1 ps) so as to measure the peak value g ⁽²⁾(0)=1 for entangled photons (see Table 1) one will need detectors with the response time T~0.1 ps, which do not exist yet. On practice, one just measures non-normalized second order correlations G ⁽²⁾(τ) with the time lag τ defined by the optical path difference to two detectors.

Future work will include the development of larger arrays of SPADs, the integration of on-chip data processing based on equation (1), and the extension to other g ⁽²⁾-imaging applications.