{'fraction_non_alphanumeric': 0.039540217901266744, 'fraction_numerical': 0.016603242454256373, 'mean_word_length': 4.203716262037163, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 18, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Systematic analysis of SNR in bipartite Ghost Imaging with classical and quantum light

We present a complete and exhaustive theory of signal-to-noise-ratio in bipartite ghost imaging with classical (thermal) and quantum (twin beams) light. The theory is compared with experiment for both twin beams and thermal light in a certain regime of interest.

Among them a great interest has been attracted by the so called ghost imaging [26] (GI), based on the correlation in spatial intensity fluctuations (speckles [27]). In this technique a light beam crosses (or it is reflected) by an object to be imaged. However, the beam that crossed the object is detected by a detector without any spatial resolution (bucket detector). The image of the object is retrieved when the bucket detector signal is correlated with the signal of a spatially resolving detector measuring a light beam whose noise is spatially correlated to the previous beam (reference beam). The first demonstration of this technique was achieved with nonclassical states of light [28], known as twin beams, produced by parametric down conversion (PDC). Then it was shown, both theoretically and experimentally, that this result can be achieved also with beam-split thermal light [29][30][31][32][33], although with a smaller visibility.
In practice this technique can be useful in the presence of phase distortions (for example when the beam crosses a diffusive medium, like fog), where intensity correlations with a second beam allow one to retrieve the spatial information [34].
Due to the conceptual and practical interest attracted to GI, a lot of works have been devoted to clearly theoretically describe [21,[35][36][37] and to improve this protocol [38][39][40][41][42]. The main parameters of GI discussed in these works are typical for any kind of imaging: the signal-to-noise ratio (SNR), which shows how well the image of the object is distinguishable from the background, and the resolution. The latter is related to the number of speckles (spatial modes) contained in the image and is basically the number of elementary details of the object that can be reconstructed in the ghost image. Note that various terms are used in the literature: SNR [35,43,45] is sometimes called contrast-to-noise ratio (CNR) [46]; number of spatial modes [37,43] is sometimes called the number of details of the image [40].
Despite the huge amount of literature devoted to GI, several points are still not clear concerning this technique. For instance, in all existing experimental works, GI is performed by measuring either normalized intensity correlation functions (CFs) [40,43] or intensity CFs with the background subtracted (also called CFs of intensity fluctuations) [30][31][32]. At the same time, up to recently, most theoretical papers considered only the simplest CFs, without normalization or background subtraction. In a recent theoretical paper, Chan et al. [46] for the first time considered both the previous procedures and showed that they lead to almost the same SNR value, while GI through the measurement of usual CFs has much smaller SNR. However, they only considered GI with bright thermal light and presented only the results of numerical simulations. In [35] Erkmen et al. developed a theory of the SNR of GI in the analog detection scheme using the not-normalized intensity correlation, providing results when ac-coupling of the photocurrents is performed by introducing signal frequency filtering.
In this paper, for providing the reader with a general discussion of second-order [51] GI in view of practical applications, we generalize the theory of Ref [46] by including into consideration, beyond thermal light ghost imaging (ThGI), the case of twin beams (TwGI) and assuming an arbitrary brightness of the light source. Moreover, we consider another protocol for GI that exploits the variance of the difference signal from the object and reference detectors [37]. It turns out that this method of GI provides the same SNR as the background subtraction and normalization of the CF for middle and high brightness of the source.
We analyze the influence of important aspects like the brightness of the source, the losses, the number of spatiotemporal modes collected by the detector, and the resolution on the performance of GI in terms of SNR.
Finally, in this work we compare the theory with the experiment, for two particularly interesting cases: mesoscopic twin beams (photon numbers per spatiotemporal mode about 0.2) and bright pseudo-thermal light, both detected by CCD cameras. The former is the best GI option if high intensities should be avoided (for instance, for the sake of not damaging the imaged object), and the latter is the simplest way to achieve high SNR if using bright light is not a problem.

II. THEORY
In general, in order to obtain a ghost image, a spatially incoherent beam interacts with the object and then is collected by a bucket detector without any spatial resolution, while the correlated beam, which does not interact with the object, is registered by a spatially resolving detector, namely an array of pixels. This procedure is repeated K times. The image of the object is described by a parameter S(x j ), where x j represents the position of the pixel j in the reference region. Usually, S has the form , a function f that involves the correlation function E[N p 1 N q 2 (x j )], (p, q ≥ 0), of the observable corresponding to the total number of photons collected at the bucket detector and at the j-th pixel of the reference arm, respectively. Here, represents the average over the set of K realizations. For the sake of simplicity we consider an imaged object defined by two levels of transmission, T = 1 and T = 0.
The signal-to-noise ratio of a ghost imaging protocol can be defined as the ratio of the mean "contrast" to its standard deviation (mean fluctuation): where S in and S out are the intensity values of the reconstructed ghost image, when x j is either inside (T = 1) or outside (T = 0) of the object profile, respectively, and δS ≡ S − S is the fluctuation. The mean value, denoted by . , represents the theoretical expectation value, which can be estimated in practice by performing space averages over the regions "in" and "out" of the ghost image. The ghost image can be reconstructed by exploiting different GI protocols, namely different parameters S. Although all the protocols exploit the correlations between the object beam and the reference one, not all off them give the same results in terms of SNR. In particular, in the following we analyze four different protocols based on using: i) the Glauber intensity correlation function, when S(x) = G (2) Here, the quantum expectation value is obviously S(x) = N 1 N 2 (x) , and the standard deviation can be evaluated directly as δ 2 S(x) = δ 2 [N 1 N 2 (x)] /K. The photon numbers are considered as quantum operators and the quantum mean values are evaluated in the appendix, where we describe in detail our theoretical model.
ii) the normalized intensity CF, when S(x) = g (2) We note that, unlike all the other parameters considered at point i, iii and iv, the standard deviation of g (2) measurement was evaluated in our work by means of the uncertainty propagation of the quantities G (2) iii) the covariance, or CF of intensity fluctuations, In this case we have Cov( ). The coefficient depending on the number of realizations, namely the number of acquired images, is the usual one allowing the unbiased estimation of the theoretical correlation. If we assume that E[N 1 ] and E[N 2 ] can be well represented by their expectation values N 1 and N 2 in the calculation of the fluctuations (it happens if K ≫ 1), we have iv) the variance of the intensity difference, . Whatever the parameter S, it involves in some form the CFs, which are second-order intensity moments. Therefore, the evaluation of its fluctuations requires expressions for the first-to fourth-order moments of the intensity. These are calculated in the Appendix, where the details of our model are extensively described.
We consider several parameters that influence the SNR.
• The brightness of the source, which is related to the average number of photons µ per single spatio-temporal mode.
• The total number M of spatio-temporal modes collected by each element (pixel) of the spatially resolving detector in the reference arm. It is approximately given by the product of the number of spatial modes, M sp = max[A pix /A coh , 1], and the number of temporal modes, M t = max[T det /T coh , 1]. Here, A pix , T det are the pixel size and the detection (integration) time, respectively, while A coh , T coh represent the characteristic coherence area (roughly, the speckle size) and the coherence time of the source. From the viewpoint of the light statistics, the total number of modes is M ≥ 1; it is considered equal to 1 when A pix ≤ A coh and T det ≤ T coh . While through the paper we always consider a usual situation in which the integration time is much larger than the coherence time, the case of A pix ≤ A coh and A pix > A coh are both analyzed in the text. We stress that M should not be confused with the number of realizations, for instance the number of acquired frames K: M is the number of modes collected by a pixel in a single acquired frame, so that I ≡ N 2 = η 2 M µ represents the total number of photons detected in the pixel, further called the "illumination level".
• We take into account the overall transmission-collection-detection efficiency of the two channels, 0 ≤ η j ≤ 1 with j = 1, 2, i.e. the probability to detect an emitted photon. According to our model, developed in the appendix, the situation in which the pixel of the reference detector is smaller than the coherence area (A pix ≤ A coh ) corresponds to a reduction of the collection efficiency. Therefore, η 2 includes also a factor that takes into account for the geometrical collection probability of a photon in a certain spatial mode, • Finally, the number of spatial resolution cells of the reconstructed image is represented by R = A in / max[A pix , A coh ], with A in the area of the ghost image where T = 1. When the pixel size is equal or smaller than the size of a single spatial mode, R is determined by the number of spatial modes, otherwise it is given by the number of pixels in the area A in . The ideal condition for maximizing in principle the resolution and the collection efficiency is A pix ≃ A coh , that is, the pixel size should approximately coincide with the speckle size [2].
As an example of a specific relevant case, in Tab. I we report the expressions for the SNR calculated for the ideal case of unity transmission/quantum efficiency and collection efficiency in both optical paths, η 1 = η 2 = 1. The SNR values are normalized to √ K.
TwGI ThGI TABLE I: Expressions for the signal-to-noise ratio (SNR) of the reconstructed ghost image with thermal light (ThGI) and with twin beams (TwGI) for different protocols described in the text, in the lossless case (η = 1).
First we observe that both for TwGI and ThGI the performances of the scheme based on G 2 are drastically worse than for the other methods because the SNR drops as SN R G 2 → 1/R (R ≫ 1) and as 1/ √ M (M ≫ 1) with respect to the resolution parameter R and with respect to the number of modes detected by each pixel M . At the same time, all the other protocols (exploiting g 2 , Cov and V ar) scale as 1/ √ R (R ≫ 1) and asymptotically ∼ const for M ≫ 1. We report the theoretical dependence of SNRs on R in Fig. 6 and Fig.7 in Sec. III C. This fact makes the not-normalized and not-subtracted correlation function G 2 , albeit largely considered in the literature (see, for instance, Refs. [35,45,46]), inappropriate for ghost imaging when compared with the other protocols. Thus, we do not investigate it further through the paper. Fig. 1 presents the SNR for all the protocols as a function of the illumination level when M = 1 and η = 1. In the typical situation of large R (R > 10) it turns out that SN R g 2 , SN R Cov and SN R V ar behave the same for large values of detected photons I ≫ 1, approaching approximatively the upper value of (2R) −1/2 , although they reach this bound at different values of the brightness. In particular, for all ThGI protocols the bound is reached when I ≫ 1. For TwGI based on the variance and covariance, the condition depends on the resolution, namely I ≫ 1/(2R), while our calculation shows that g (2) reaches the flat region as soon as I ≫ 1/(2R 2 ). Concerning the performance of a very low-brightness source, Fig. 1 shows that the advantage of using twin beams is in general very pronounced but, even with the same source, some protocol seems to be more convenient than others. In particular, for ThGI all the protocols scale as ∝ I with the exception of the variance method, for which SNR approaches zero faster, i.e. ∝ I 3/2 . For TwGI all the methods lead to the same asymptotic behavior ∝ I 1/2 . Therefore, we can conclude that for the ghost imaging of a complex object the three protocols exploiting SN R g 2 , SN R Cov and SN R V ar have the same performance at medium and high intensities (for the same source, quantum or classical), while for very low brightness, preference should be given to the normalized correlation function g 2 . Now we focus on the difference between the signal-to-noise ratio of TwGI, SN R T w , that is obtained in the limit of large R, and the corresponding SN R T h obtained with thermal light in the same limit, by using the covariance as the parameter for the GI reconstruction. In order to perform a fair and probably more useful comparison in a practical situation, we already introduced the "illumination" I = ηM µ, i.e. the total number of photons detected in the pixel in the single image. In Fig. 2(a) we report the quantum and classical SNR, for different values of the transmission/detection efficiency, balanced for the two channels, while fixing M = 1 (a single temporal mode and a single spatial mode are detected). Not surprisingly, the curve for the thermal light is not influenced by the losses, while the TwGI curve, in general showing a better performance than the classical one, for high losses approaches the ThGI. The advantage of twin-beam light is evident for a low-brightness source, i.e. when the number of photons per spacetemporal mode µ ≪ 1. On the contrary, when µ ≫ 1, thermal and twin-beam states of light produce the same results (Fig. 2). The physical reason can be understood with the help of Eqs. (6), (7) and (8) for the two-mode statistics. Both classical and quantum light have the same single-mode thermal fluctuations δ 2 n j = ηµ + η 2 µ 2 . However, the photon-number correlations for thermal light, δn 1 δn 2 T H = n 1 n 2 = η 2 µ 2 , are relevant when the number of photons is large, n j ≡ µ > 1, and do not include the shot-noise component ∝ µ, which remains uncorrelated. On the other hand, the correlations in twin beams are δn 1 δn 2 T B = η 2 µ(1 + µ), which shows that even the shot-noise component of the fluctuations is correlated. Moreover,classical two-mode correlations, generated by a beam splitter obeying the Bernoulli statistics, are not reduced by the losses, which are again described by a beam splitter model (see the Appendix and the discussion after Eq. (8)).
It is also important to present ( Fig.2(b)) the unbalanced case, η 1 > η 2 , which in our model includes the typical situation in many experiments where the pixel is smaller than the coherence area. In this case, η 2,coll ∼ A pix /A coh ≤ 1, and the advantage of the TwGI is reduced. We observe that thermal ghost imaging is weakly influenced by the unbalancing. This indicates that, in order to make the best of the TwGI, in the photon-counting regime one should properly set the pixel dimension.
In Fig.3, the dependence of SNR on the number of modes M collected by a single pixel is shown. We see that the TwGI is insensitive to M while the SNR for ThGI is reduced for large M but can be recovered at higher values of illumination I. A discussion about that point can be found in the Appendix after Eq. (10). This means that ThGI performs best in the single-mode regime, while TwGI, even at low brightness, can be easily brought to a high-illumination level by accumulating a large number of temporal modes in a single frame, without any decrease in the SNR. However, provided that the source has the same parameters, i.e. the same brightness µ, the same M and therefore the same I, the advantages of quantum light can be of several orders of magnitude, even for large illumination I. In conclusion, our discussion demonstrates that single-mode (M = 1) ThGI is the best solution, considering the less demanding experimental resources, when one is not limited in the brightness of light. On the other hand, in situations where a low light level, µ < 1, is needed, for instance because of the photosensitivity of the object, quantum light provides a larger SNR (even for a relatively large illumination level).

III. THE EXPERIMENT
For the experimental comparison of GI obtained with thermal light and twin beams, we used two setups. In the first case, GI was performed with bright pseudo-thermal radiation (more than 10 4 photons per radiation mode). In the second case, nonclassical radiation was used, namely, squeezed vacuum obtained via type-II parametric downconversion.
In both experiments, we used the 'numerical mask' method suggested in Ref. [45]. We selected two regions, one in the object beam, where we virtually inserted the mask, and the other in the reference beam. In the bucket channel, we simulated the presence of a completely opaque mask, which fully transmitted light within a certain area A in and fully absorbed light outside it. This was done by simply setting to zero the signals of all pixels outside of this area. This method, which is only possible with an array photodetector in the bucket channel, presents the advantage of arbitrarily varying the shape and size of the mask without any technical difficulties and provides an ideal limit for this part of the set-up in a realistic situation.

A. Experimental setup with pseudo-thermal light
The experiment with pseudo-thermal light ghost imaging was performed with the setup shown in Fig. 4 [43]. Second-harmonic radiation of a Nd:YAG laser with the wavelength 532 nm, pulse duration 10 ns, and repetition rate 47 Hz was incident on an Arecchi's rotating ground-glass disc generating pseudo-thermal light [44]. The radiation scattered by the disc was split by a non-polarizing beamsplitter and then both output beams were registered by two different parts of the matrix of a commercial digital photographic camera Sigma SD14. The far-field speckle patterns, symmetrical in the two parts, were focused on the matrix by a lens with the focal length 50 mm. As a result, the typical speckle size on the matrix was about 0.230 mm. With the pixel size of the matrix being 7.8 µm, a single speckle occupied, on the average, an area of 30 × 30 pixels. This was confirmed by the measurement of the spatial intensity correlation function.
The rotation rate of the disc was chosen in such a way that the speckle pattern did not change noticeably during a single pulse but changed completely from pulse to pulse. The camera was taking frames with exposure times 1/50 s, which was much larger than the pulse duration but still less than the distance between the pulses. As a result, most of the frames contained the speckle pattern from a single pulse ( [43], Fig. 4) and few were 'empty'. The total number of captured frames was 5000. In each frame, we selected square regions in the two symmetrical parts (Fig.  6 of Ref. [43]), further referred to as the object and reference channels. As before, in the object channel, the 'masks' were introduced numerically, by simply ignoring the readings of the pixels outside them.
The intensity of the laser light was reduced by using a Glan prism; by measuring the signals from separate pixels of the camera versus the input intensity we verified that the camera was operating in a linear regime and was not saturated. It is important that the output data of the camera was written in the RAW format. The mean values of signals from separate pixels were on the order of 300, which indicates that the number of incident photons was at least as much as that. Hence, per single speckle, i.e. for a single radiation mode, we had photon numbers as high as 3 · 10 4 . As the pulse duration was much smaller than the typical time of the speckle pattern fluctuation, and the pixel size was much smaller than the speckle size, the number of modes was M = 1. The transmission/quantum efficiency was quite low (less than 0.01), due to the beam attenuation, but, as discussed in section II, this is not important for ThGI.

B. Experimental setup with type-II PDC
In contrast to most ghost imaging experiments reported in the literature that are realized with faint PDC, here we set an experiment in the mesoscopic regime using a mesoscopic PDC source and a CCD camera as a photon number resolving detector. This represents an important step in view of practical applications of this method. The setup is reported in Fig.5.
For generating the twin beams we used a type II BBO non-linear crystal (l = 7 mm) pumped by the third harmonic (355 nm) of a Q-switched Nd:YAG laser. The pulses had a duration of T p = 5 ns with a repetition rate of 10 Hz and a maximum energy, at the selected wavelength, of about 200 mJ. The spurious non-Gaussian components were eliminated by spatial filtering (a lens with a focal length of 50 cm and a diamond pin-hole, 250 µm of diameter) of the beam that was re-collimated with a diameter of w p = 1.25 mm at the crystal. After the BBO, the pump was stopped by two UV mirrors (with ≃ 98% declared transmission at 710 nm), and by a low-pass frequency filter (≃ 95% transmission at 710 nm). The orthogonally polarized down-converted signal and idler beams were separated by two polarizers (≃ 97% transmission) and finally the far field was registered by a 1340X400 CCD (Princeton Pixis:400BR with pixel size of 20 µm, quantum efficiency of 80% and readout noise 4 electrons/pixel). A mechanical shutter provided the synchronization with the laser so that each acquired frame corresponded to the PDC emission generated by a single shot of the laser. The far field was observed at the focal plane of the lens with 10 cm focus in an 'f − f ' optical configuration. Our setup was developed exploiting the potentiality of the quantum scheme [15], which, according to the theory, should outperform the classical light scheme for the same brightness (number of photon per spatio-temporal mode). For optimizing the detection of quantum correlations, basically two requirements should be fulfilled. First, transmission/detection and collection efficiency must be large. We realized it by reducing the optical losses after the crystal and using a CCD camera with high quantum efficiency, while the collection efficiency can be kept high if the pixel area A pix is larger than the coherence area A coh (approximatively a speckle), i.e. a pixel should collect more than one spatial mode [27]. The second point is that the detection noise must be low compared with the shot noise. For both these reasons, it is very efficient to perform hardware binning of the physical pixels. It consists of grouping the physical pixels in squared blocks, each of them being processed by the CCD electronics as single "superpixel". The number of photons collected by a superpixel is the sum of photon numbers of each pixel, whereas the readout noise is just slightly increased with respect to the one of a single pixel. We chose its dimension comparable to A coh , i.e. the pixel (hereinafter we discard the prefix "super") size was set to 240 × 240(µm) 2 . The "illumination" I, defined as the number of photons detected by a pixel in a single shot image (frame), was about I = 1900. The expected number of temporal modes M t = T p /T coh detected in one frame was 5 · 10 3 , considering the coherence time T coh of PDC around one picosecond. The number of spatial modes detected by the pixel was M sp = A pix /A coh ≃ 4. A measure of the coherence area, basically determined by the pump transverse diameter, is the size of the spatial cross-correlation of the signal and idler intensity patterns that has been evaluated to A coh ∼ 120 × 120(µm) 2 . The total number of modes M = M t × M sp turned out to be compatible with the level of excess noise E ≡ I/M ∼ 0.12 due to the thermal statistics in the single beam. Thus, the number of photons per single space-temporal modes was µ = I/(ηM ) ≃ 0.20. Under these conditions, we recorded 4000 frames.

C. Experimental Results
After having discarded all the frames corrupted by incoming cosmic rays, as well as 'empty' frames in the ThGI setup, we started the analysis by defining two correlated regions in the object and in the reference arm. In the TwGI setup, the regions included 13 × 15 = 195 pixels, while in the ThGI setup, they contained 750 × 750 pixels. These regions were the same for all the frames.
In the twin beams set-up, one of the main problems we had to cope with was a strong instability of the Q-switched laser power from pulse to pulse (about 14%). Since the expectation value of the number of photons per mode µ is proportional to the square of the power P of the pump (µ ∝ sinh 2 (const √ P )), fluctuations of P lead to the fluctuations of µ. As a consequence, the temporal statistics on many pulses is characterized by the meanμ and the variance V (µ). Indeed, using the equations presented in the Appendix, we have obtained the expectation value of Cov inside and outside of the mask taking into account the pump instability: When M and R are large, as in our experiment, the fluctuations of the pump dominate and this introduces a nonzero covariance for pixels in the reference arm that are, in principle, uncorrelated with the mask. The overall result is a strong noise of "background correlation" that hides the PDC spatial correlation, essential in the ghost image reconstruction. As a matter of fact, this effect would lead to an increased noise on the ghost image and a consequent decreasing of SNR.
In order to compensate for the instability of the pump, we normalized each frame by its average value taken over a certain CCD spatial region. This procedure completely removed the problem.
In Fig.6 we present the obtained experimental values of the SNR versus the resolution parameter R of the mask for the case of TwGI. The data set is in agreement with the theory where the experimental parameters M , µ and η have been estimated independently. However, the best fit is obtained when η = 0.42 that is smaller than the value 0.62 that was measured very accurately in [16]. This discrepancy can be explained by the relatively small pixel size (reduced collection efficiency of correlated photons) and the presence of background noise generated by detection and straylight. A detailed analysis of such problems can be found in [16]. We stress that in our experiment we exploited successfully the theoretical results pointing out that the multimode regime is the most appropriate for TwGI. In fact, by collecting M ∼ M sp · M t = 4 · 5000 we reached a level of illumination I = 1900 such that the experimental sources of noise (like electronic read-noise of the CCD and straylight ecc.) were almost negligible. For G 2 we only report the theoretical prediction, which is close to zero. In accordance to this prediction, our experimental data involving 4000 frames was insufficient to produce any ghost image with G 2 protocol.
In summary, we have demonstrated ghost imaging with twin beams in a mesoscopic regime (µ = 0.2, M ≡ M sp · M t = 4 · 5000, I = 1900) by a photon-number resolving detector for such parameters of the source (number of photons per mode and the illumination level) for which thermal ghost imaging provides much worse SNR, according to the values reported in Fig. 3.
On the other hand, in figure 7 shows the results of ThGI experiment. The illumination level was about I ≈ 300, photon number per mode µ ≈ 10 4 , the number of modes M = 1, and the resolution parameter R was varied from R ≈ 4 to R ≈ 200 by changing the size of rectangular 'numerical masks'. One can see that the experimental points agree well with the theoretical dependence for SNR, which, according to Table I, is in this case similar for all three methods (based on variance, covariance, and g (2) ) and given by SN R = √ 2R.

IV. CONCLUSIONS
For providing the reader, in view of practical applications, with a general exhaustive description of ghost imaging, we have discussed both thermal and twin-beam GI for four different protocols, based on the measurement of the following values: the Glauber correlation function G (2) = N 1 N 2 , the normalized Glauber's correlation function g (2) = G (2) / N 1 N 2 , the covariance Cov = δN 1 δN 2 , and the variance of the difference signal V ar = δ 2 (N 1 − N 2 ) .
A first significant result, in accordance with the analysis in [46] for thermal light, is that the G (2) method performs much worse than the others (the SNR falling rapidly with the resolution parameter, i.e. the number of speckles transmitted through the mask in the object channel, and with the number of spatio-temporal modes detected by the reference detector): thus, albeit often considered in the literature for the sake of simplicity, the G (2) method is not worth using in real applications.
On the other hand, the other methods have similar performances. For all of them, twin-beam GI performs largely better than thermal-light GI at low illumination levels (few photons per pixel in a single shot), while they become equivalent at high illumination levels (many photons per pixel in a single shot). Furthermore, TwGI is insensitive to the number of modes collected by a pixel and can be brought to a large illumination level by accumulating a large number of temporal modes in a single shot, even maintaining low brightness, which means low photon number per spatio-temporal mode. This is what is realized in our experimental demonstration, in which about 20000 modes are accumulated in the 5 ns duration of the pump pulse. ThGI performs best in the single-mode regime. In summary, ThGI is preferable whenever one does not have limits on the brightness (e.g. photo-sensible samples), due to its simpler realization.
Finally, we have demonstrated the perfect agreement of the developed general theory of GI with experiment in some particularly interesting (and sometimes unexplored) regimes for both TwGI and ThGI.

V. APPENDIX
In this Appendix we calculate the correlation functions needed for evaluating the SNR of different GI schemes. Both for thermal and PDC light we first consider the single-mode statistics and then we extend it to a more general multi-mode case. Only later we introduce the bucket detector that integrates the light over all the spatial modes in the 'mask' ('object') beam.
PDC [48], in the approximation of a plane and monochromatic pump beam of frequency ω p , can be described by an evolution operator U that is the product of independent operators, each acting on a couple of symmetric modes (q, Ω) and (−q, −Ω) [47]: where ±Ω is the detuning of the signal(1) and idler (2) photons with respect to the degeneracy frequency ω p /2, ±q are the symmetric and correlated transverse wavevectors, a † 1 and a † 2 are the creation operators of the signal and idler photons, respectively, and ξ(q, Ω) is the coupling coefficient. Therefore, it can be seen as a collection of independent processes where couples of symmetric modes are correlated and the output field operators for signal and idler are where U (ξ) = cosh(ξ) and V (ξ) = sinh(ξ) determine the gain of the process by their maximum value and the coherence time by their inverse width in Ω. Here the phases are not considered because they do not play any role in the photon-number statistics.
In our model we assume η j (j = 1, 2) to be the detection probabilities for the photons in the mode j. Here, η j is a coefficient that includes transmission, detection (quantum efficiency, QE, of the detector) and collection efficiency η j,coll of the photon emitted in a certain spatial modes. In particular, in our model the collection efficiency can be less than one when the pixel of the spatially resolving detector is smaller than the spatial mode at the detection plane. In that case, the probability to collect a photon emitted in a certain mode is η j,coll = A pix /A coh . All the losses can be taken into account by a model in which a beam splitter with transmission η j is preceding an ideal detector (QE= 1) [49]. For the calculation of the photon-number correlation function n j = b † j b j , this leads to a simple substitution b j → √ η j b j in the normally ordered expressions for the operators. The single-mode photon-number statistics as well as the two-mode correlation can be obtained through a straightforward calculation by the relations (5), with the input state being the vacuum one: n 2 1 n 2 = µη 1 η 2 1 + 6µ 2 η 1 + µ(2 + 4η 1 ) , n 2 1 n 2 2 = µη 1 η 2 1 + 24µ 3 η 1 η 2 + 6µ 2 (η 1 + η 2 + 4η 1 η 2 ) + 2µ(1 + 2η 1 + 2η 2 + 2η 1 η 2 ) , n 3 1 n 2 = µη 1 η 2 1 + 24µ 3 η 2 1 + 18µ 2 η 1 (1 + η 1 ) + 2µ(1 + 6η 1 ) .
For our purposes we focus on the correlation in the far field obtained as the focal plane of a thin lens of focal length f in an f − f configuration and around degeneracy (Ω = 0). Here, any transverse mode q is associated with a single position x in the detection (focal) plane according to the geometric transformation (2cf /ω p )q → x [4], with c being the speed of light. We observe that the diffraction effects arising from the transverse finite dimension of the optical system (typically the pump profile in the twin beam generation process) lead to a nonzero transverse coherence area, A coh , in the far field, representing the uncertainty in the arrival point of the correlated photons. The fluctuation of the difference photon number is δ 2 (n 1 − n 2 ) T W ≡ δ 2 n 1 + δ 2 n 2 − 2 δn 1 δn 2 = 2ηµ(1 − η) for η = η 1 = η 2 . Thus, in the ideal lossless case (η = 1), a perfect correlation appears in the photon number since δ 2 (n 1 − n 2 ) T W = 0 (sub-shot-noise regime [14,15]).
In the case of thermal radiation we assume a set of independent modes a i (i=1,..) addressed to the input of a 50% beam splitter performing the transformations b 1,i = (a i + a vac,i )/ √ 2 and b 2,i = (a i − a vac,i )/ √ 2, where a vac is the mode at the vacuum input port of the beam splitter. The photon-number correlation of the beams at the two output ports are obtained by exploiting the rule (a † i ) p a p i = p! a † i a i p , for normally ordered moments of thermal light. For the single-beam moments, taking into account the losses, we found again the expressions in Eqs. (6). This is in agreement with the fact that each beam from PDC, considered separately, manifests thermal statistics. On the other hand, the correlations of photon numbers between thermal beams 1 and 2 are n 1 n 2 = 2µ 2 η 1 η 2 , (8) n 2 1 n 2 = 2µ 2 η 1 η 2 (1 + 3µη 1 ) , n 2 1 n 2 2 = 2µ 2 η 1 η 2 1 + 12µ 2 η 1 η 2 + 3µ(η 1 + η 2 ) , n 3 1 n 2 = 2µη 1 η 2 1 + 9η 1 µ + 12µ 2 η 2 1 .
In contrast to the quantum case, the correlation between the photon numbers in the coupled modes is not perfect; hence, for thermal light we have δ 2 (n 1 − n 2 ) T h = 2ηµ, which represents the shot noise limit. At the same time, the excess noise ∼ µ 2 in the second equation in (6) disappears regardless of the losses. According to our experimental realization (but without loss of generality), we start by considering a situation in which both reference and object beams are detected by a spatially resolving array of detectors (the pixels of a CCD camera in our case). Each pixel collects a number M ≥ 1 of independent spatio-temporal modes in the far field and the total number of photons detected by the pixel is N j = Σ M i=1 n j,i , with j = 1, 2. It is easy to show that the moments of the N distribution can be expressed in terms of single-mode moments in Eqs. (6) as Since the modes of signal and idler fields are pairwise correlated, the correlators of the count numbers in two symmetrical pixels (collecting M correlated modes) are The correlators N 1 N 2 2 and N 1 N 3 2 can be obtained by simply exchanging the indexes 1 and 2 in the second and fourth equations of (10). From the first lines of Eqs. (9) and (10) we see that δN 1 δN 2 = M δn 1 δn 2 (η 1 = η 2 = 1). By substituting the two-mode correlators following from Eqs. (6) and (8) we have δN 1 δN 2 T W = N 1 + N 2 1 /M and δN 1 δN 2 T H = N 2 1 /M for quantum and thermal light, respectively. These expressions show that thermal correlations are sensitive to the numbers of modes in both channels, while the shot noise term is not reduced by the increase in M . This provides an intuitive explanation of the different behavior between TwGI and ThGI with respect to the number of modes M . Now we introduce the bucket detector whose output is represented by the operator N 1 = Σ R k=1 N k,1 , i.e. a sum of the photon numbers over R independent sets of spatial cells, each one collecting M ≥ 1 modes. R coincides with the number of pixels of the ghost image if the pixel is sufficiently large to collect more than one spatial mode, otherwise it coincides with the number of spatial modes (speckles) itself.
The correlators between the readings of the bucket detector and the ones of an arbitrary pixel in the reference channel registering light correlated with the one passing through the mask (dubbed with the subscript "in") are The analogous equation for the "out" case is obtained by using the statistical independence of N 1 and N 2,out , i.e. by substituting N p 1 N q 2,out = N p 1 N q 2,out , and then applying Eqs. (11) and (9). By using the definition (1) and Eqs. (11), (12) we arrive to our full analytical expressions for the SNR for the protocols considered in the paper.

VI. ACKNOWLEDGEMENTS
This work has been supported by Regione Piemonte. We thank Ivan Agafonov for providing the set of frames with thermal-light speckle patterns.