GMP-selected dual and lensed AGNs: selection function and classification based on near-IR colors and resolved spectra from VLT/ERIS, KECK/OSIRIS, and LBT/LUCI

The Gaia-Multi-Peak (GMP) technique can be used to identify large numbers of dual or lensed AGN candidates at sub-arcsec separation, allowing us to study both multiple SMBHs in the same galaxy and rare, compact lensed systems. The observed samples can be used to test the predictions of the models of SMBH merging once 1) the selection function of the GMP technique is known, and 2) each system has been classified as dual AGN, lensed AGN, or AGN/star alignment. Here we show that the GMP selection is very efficient for separations above 0.15'' when the secondary (fainter) object has magnitude G<20.5. We present the spectroscopic classification of five GMP candidates using VLT/ERIS and Keck/OSIRIS, and compare them with the classifications obtained from: a) the near-IR colors of 7 systems obtained with LBT/LUCI, and b) the analysis of the total, spatially-unresolved spectra. We conclude that colors and integrated spectra can already provide reliable classifications of many systems. Finally, we summarize the confirmed dual AGNs at z>0.5 selected by the GMP technique, and compare this sample with other such systems from the literature, concluding that GMP can provide a large number of confirmed dual AGNs at separations below 7 kpc.


Introduction
The existence of dual AGNs, i.e., pairs of active supermassive black holes (SMBH) at kpc separations inspiralling in the same galaxy, is one of the fundamental predictions of current models of galaxy and SMBH formation and evolution (e.g.Tremmel et al. 2017).The study of these systems allows us to test the models of formation and evolution of binary SMBHs.These systems are particular important because mergers of SMBHs in the mass range 10 4 − 10 9 M ⊙ will dominate the gravitational wave ⋆ e-mail: filippo.mannucci@inaf.it(GW) events detected by the future Laser Interferometer Space Antenna (LISA) space mission and by the Pulsar Timing Arrays (e.g., Arzoumanian et al. 2018;Amaro-Seoane et al. 2023).Dual AGNs allow us to investigate the starting conditions of the merging process, and interpret future GW results (e.g.DeGraf et al. 2023;Dong-Páez et al. 2023;Chen et al. 2022a).
Models of galaxy and SMBH co-evolution predict that dual AGNs at small separations (less than 10 kpc) are quite common, reaching a few percent of all existing AGNs at z > 0.5, where merging activity is expected to be larger than in the local Universe (Volonteri et al. 2022) and therefore must exist in a significant fraction of galaxies.Despite these expectations and a decade of active searches, only a handful of confirmed dual AGNs at small separations are currently known (see Chen 2021, for a recent compilation of the existing data).
These objects can only be detected through imaging with a spatial resolution considerably better than one arcsec, typically requiring observation from space (see De Rosa et al. 2019, for review).Being rare systems, they also require very wide field surveys.The ESA mission Gaia (Gaia Collaboration et al. 2016) is currently the only mission providing high-resolution data over large areas of the sky, therefore can be exploited to detect dual and lensed AGNs.The presence of multiple Gaia sources associated, within a few arcsec, to known AGNs (the "multiplicity" method) have been used by several authors (e.g.Lemon et al. 2018;Chen et al. 2022b) to detect both dual and lensed systems.The extra astrometric jitter induced by AGN variation in unresolved pairs was also exploited to select several multiple AGN candidates (e.g.Shen et al. 2019;Hwang et al. 2020;Chen et al. 2022b).
Recently, in Mannucci et al. (2022), we developed the GMP method which exploits the Gaia catalog to select a large number of dual AGN candidates at sub-arcsec separations.This method is based on the detection of multiple peaks in the light profiles of otherwise unresolved AGNs.Using spatially-resolved imaging and spectroscopy, we showed that the GMP is a very effective method to select dual systems, while in Ciurlo et al. (2023) and Scialpi et al. (2023) we built the first sample of spectroscopically confirmed, GMP-selected, dual AGNs at sub-arcsec separations.
The GMP technique is expected to provide a large number of confirmed dual AGNs, especially at z > 0.3 where the host galaxy does not contribute significantly to the Gaia light profile.These results will be compared with the expectations of the models of galaxy/SMBH formation and co-evolution, which predict various quantities such as the distribution in separation, luminosity and mass ratio, redshift, etc (Steinborn et al. 2016;Capelo et al. 2017;Rosas-Guevara et al. 2019;Volonteri et al. 2022;Chen et al. 2022a).The samples of GMP-selected dual AGNs will provide strong constraints on these models once the efficiency of this method in selecting multiple systems is known as a function of separation and luminosity of the components.
QSOs showing multiple components need to be classified as either a pair of AGNs (either a dual system or a single AGN lensed by a foreground galaxy) or an alignment between an AGN and a foreground star.Different possibilities exist: 1.The classification is usually obtained by spatially resolved spectroscopy from space (e.g Junkkarinen et al. 2001;Mannucci et al. 2022) or by adaptive-optics (AO) assisted spectroscopy from the ground (e.g Ciurlo et al. 2023;Scialpi et al. 2023).The result classification can be very reliable, but obtaining resolved spectroscopy is very telescope-expensive and limits the number of confirmed systems.2. Spatially unresolved spectroscopy, usually from the ground in seeing-limited conditions, can be analyzed to detect stellar features at zero velocity, revealing the presence of a star projected close to the AGN (Shen et al. 2023;Scialpi et al. 2023).These data cannot distinguish dual and lensed systems, since in both cases the spectral differences are not expected to be large enough to be seen in the total spectrum.3. Since stars and AGNs have, in most cases, different colors, multi-band spatially-resolved imaging can be used to study the nature of the two components (e.g.Chen et al. 2022b;Chen 2021).Intrinsic luminosities can also be used to iden-tify lensed systems, in many cases much brighter than the typical AGN at the same redshift.
By applying these three different classification methods to the same targets, in this paper we show that spatially-resolved imaging and spatially-unresolved spectroscopy can complement resolved spectroscopy in classifying at least part of the systems, distinguishing between dual AGNs, lensed systems, and star/AGN alignments.
In Sect. 2 we discuss the statistical properties of the GMPselected systems in terms of luminosity ratios and projected separation.In Sect. 3 we present and compare the different classification techniques.In Sect.3.1 we describe the Large Binocular Telescope (LBT) observations of seven targets aimed at measuring the near-IR colors of each component.Sect.3.2 describes how the analysis of the integrated spectra can reveal the presence of stars projected close to the AGN.In.Sect.3.3 we present spatially-resolved ESO Very Large Telescope (VLT) and Keck spectra of five GMP-selected targets, which include two confirmed dual AGNs and one compact lensed system.In Sect.3.4 we compare the results of the three classification.Finally, in Sect. 4 we summarize the state of GMP target observations and discuss the results in the broader context of observed dual AGNs at small separations.

GMP selection efficiency and contamination
The GMP technique is based on the detection of multiple peaks in the light profiles of known AGNs in the multiple Gaia scans.The Gaia DR3 catalog (Gaia Collaboration et al. 2022) provides the parameter ipd_ f rac_multi_peak (FMP) that contains the percent fraction of scans showing a multiple peak for a given target.The value of this parameter, which ranges from 0 (multiple peak never detected) to 100 (always detected), is determined by several factors.Separation and luminosity ratios are important parameters to consider, as pairs with too small separations and/or too large luminosity ratios will not produce a measurable double peak.Since Gaia obtains sky scans in multiple directions, an object pair will not show a double peak when the scanning direction is perpendicular to the separation line.In addition, when scanning objects with separations δ larger or of the order of the Gaia photometric window (0.7 ′′ , Gaia Collaboration et al. 2016), the components may be separated into different catalog entries, each one with a single peak and a low value of FMP.In contrast, since Gaia projects on the same array two different fields 106.5 • apart (Gaia Collaboration et al. 2016), double peaks can be occasionally detected in light profile of isolated objects due to the alignment with a source in the other field.Thus, the selection requires a statistical approach, and needs to tackle the following complementary challenges: 1. Define the minimum value of FMP to keep the level of contamination by false positives (single systems with values of FMP above threshold) low.2. Measure the efficiency of the GMP method of detecting multiple objects as a function of magnitudes and angular separations δ of the components.et al. (2018).These images have depth and spatial resolution better than Gaia and, as a consequence, can be used to obtain a sample of projected pairs and test the performance of the GMP method.Even if the central parts of the GCs are avoided, some of these fields are very crowded.We have limited our analysis to the 84 fields having densities such that the average distance of a source to the closest object is above 5 ′′ .Higher densities tend to have too many objects at close projected distances, a regime that is not present in the high galactic latitude fields used for the GMP selection, and that could affect Gaia detection properties.
Restricting to lower densities would result in a low number of projected pairs at low separation.Regions within 5 ′′ of bright stars (G<15.5)were excluded because the rectangular Gaia primary mirrors produce long diffraction spikes that could produce the detection of secondary peaks in nearby sources.

Contamination
The  1: Main properties of the targets used to compare the various techniques, and source of the observational data.For each object we report the total Gaia G-band magnitude, with the exception for J0812+2007, which is separated into two targets in the Gaia catalog, and or which we report the magnitude of both objects.Integrated spectra are from SDSS DR16 (Lyke et al. 2020), LAMOST (Yao et al. 2019), and 2QZ (Croom et al. 2004).

Efficiency
The efficiency of pair detection can be estimated by the ratio of the number of recovered pairs vs. the total number of pairs present in the HST data, as a function of magnitude and separation δ.Our catalog contains about 1300 projected pairs with δ < 1.2 ′′ and with primary G<21 mag.The fraction of recovered objects (FMP>8) with G<20.5 mag as a function of separation is shown in the right panel of Fig. 1.The pairs that are not recovered (red histogram) dominate the counts at δ < 0.15 ′′ .Above this limit, virtually all the projected pairs are recovered within our magnitude limit.GMP objects that are isolated in the Gaia catalog, i.e., with no other object within 1.5 ′′ of separation, are shown in green.System corresponding to two nearby entries in the catalog, and that therefore can also be identified without the use of the GMP method, are shown in blue.The "isolated" population dominates separations between 0.15 ′′ and 0.6 ′′ , while above this limit most of the objects are separated in the catalog.
Of course, the "separated" population can also be identified by simply looking for pairs of objects at close separations in the Gaia archive, without considering the FMP parameter.In particular, this population becomes dominant at δ > 0.6 ′′ , while below this limit most of the targets are "isolated" and can only be recovered by the GMP method.
For G<20.5 mag and δ > 0.15 ′′ , the selection efficiency does not depend critically on magnitude and separation.This is shown in Fig. 2, where we plot the efficiency as a function of both primary and secondary magnitudes at different separations.The estimate of the G-band magnitude of the secondary stars undetected by Gaia is derived from the photometry in the HST F814W band by applying a color correction derived from the detected sources.Provided that the fainter object has G<20.5 mag and the separation is δ > 0.15 ′′ , the efficiency is high and not strongly dependent on these parameters.This is a critical point that will allow us to compare the observational results with the models.Fig. 2 also shows that, if the primary component has G>17 mag, the GMP method can efficiently select systems with a secondary as faint as G=20.5 mag.The range of brighter magnitudes, G<17 mag, is not sampled in our fields but it is not very relevant because the majority of AGNs at z > 0.5 are fainter than this limit.As a consequence, the range of luminosity ratio sampled by the GMP method depends on the luminosity of the primary source.For G=18.5 mag, a typical value for bright AGN at z > 1, the magnitude difference with the secondary member is limited to ∼2 mag, a factor of ∼ 6 in luminosity.This limits the number of detected sources but guarantees that both components give some contribution to the integrated spectra, as discussed in the next section.Larger luminosity ratios can be covered with the ESA satellite EUCLID (Euclid Collaboration et al. 2022).Multiple AGN detections with Euclid will be discussed in a future paper (Ulivi et al., in preparation).

Classifying multiple systems in different ways
Since the GMP technique only provides the evidence of the presence of multiple components in an AGN, the classification of these components requires additional observations.As shown in Mannucci et al. (2022), Ciurlo et al. (2023), andScialpi et al. (2023), resolved spectroscopy is the ultimate technique but is observationally very demanding and requires high spatial resolution.For this reason, in this section we investigate how spatiallyresolved colors (Sect.3.1) and total spectra (Sect.3.2) can provide a classification of the components.To do this, we have considered a sample of nine GMP-selected systems (see 1) with archival total spectra.We have obtained spatially-resolved colors of seven of them with LBT/LUCI, and spatially resolved spectra of five of these systems with Keck/OSIRIS and VLT/ERIS.Target selection is described in Sect.3.1.1 and 3.3 for the LBT and Keck/VLT samples, respectively.The classifications based on these three data sets are compared in Sect 3.4.

LBT imaging
We obtained Adaptive-optics (AO) assisted imaging of seven GMP-selected systems at the LBT, including the five multiple systems presented in Mannucci et al. (2022).We extracted the targets from the Milliquas catalog v7.2 (Flesch 2021), selecting objects with secure AGN classification and redshifts, and within ∼ 30 ′′ of suitable natural-guide stars (NGS) to drive the AO system.After cross-correlating with the Gaia archive, we selected systems with FMP≥ 10, a more conservative limit that the value FMP≥ 8 defined in Sect. 2. Six out of seven targets are 'isolated', i.e., have no companion in the Gaia catalog, while J0812+2007 has a known companion 0.66 ′′ apart.
Observations were performed on March 9th 2022 for objects J0732+3533, J0812−0040, J0812+2007, J0927+3512, and J0950+432, and on June 3rd 2022 for J1318−0136 and J1510+5959.We used the LUCI1 imager and spectrograph (Mandel et al. 2007) with the SOUL AO module (Pinna et al. 2016(Pinna et al. , 2021)) In all cases, we used an image sampling of 15 mas pixel −1 .Data reduction was performed with our custom made python reduction routine (PySNAP).The seeing was rather variable during the June 2022 run, and we automatically rejected single images showing a point-spread-function (PSF) significantly worse than the average.All the systems show two components at sub-arcsec separation.The K s -band imaging of these sources was presented in Mannucci et al. (2022).In Table 2 we summarize the observed magnitudes and the main properties of the observed objects.All the sources are well resolved in both H and K s bands, except for J0812+2007 in the H-band (see below).This image was obtained under bad seeing conditions, and the resulting PSF (FWHM∼ 0.65 ′′ ) is actually larger than the separation (0.54 ′′ ).As a consequence, the relative photometry of the two components (see below) in this band is unreliable.
Photometry was performed by PSF fitting.Since the two sources in each system are at very small distance, they are expected to have the same PSF, thus we fitted two elliptical Gaussian functions with the same shape (FWHM and position angle), centers defined on the image with the best resolution, and free normalization.This procedure provides reliable values of the parameters even when the PSF FWHM is not significantly smaller than the separation, and the sources are not well resolved.Since two systems are bright enough to be in the 2MASS catalog (Skrutskie et al. 2006), we used their total magnitudes to obtain an absolute flux calibration.In all cases where the sources are well resolved, the uncertainties are dominated by flux calibrations, estimated to be ∆(mag) ∼0.05.J1318-0136, and J1510+5959 are the two most compact systems, with separations of 0.27 ′′ and 0.28 ′′ , similar to the PSF FWHM.In these cases we estimate larger photometric errors of 0.10 and 0.15 mag.For the bad-seeing H-band photometry of J0812+2007 we estimate errors on the relative photometry of ∆(mag) ∼0.30.

Spatially resolved colors
Fig. 3 shows the color-magnitude diagram for each component of the systems of the LBT sample, compared with QSOs and stars from the literature.For QSOs we used the 2MASS magnitudes of the SDSS DR16Q quasar catalog 2 (Lyke et al. 2020), and the much fainter QSOs in COSMOS by Civano et al. (2016) 3 .These two catalogs probe the entire luminosity range covered by the current sample (15 < H < 21).For stars we used observed magnitudes of a high galactic latitude star sample from the SDSS, and the model magnitudes of high galactic latitudes 2 https://www.sdss.org/dr16/algorithms/qso_catalog/ 3 https://cdsarc.cds.unistra.fr/viz-bin/cat/J/ApJ/819/62stars from the TRILEGAL Galaxy model (Girardi et al. 2005;Vanhollebeke et al. 2009).Stars show relatively blue H − K s colors of about 0.0-0.2, while QSOs span a larger range and mostly have redder colors 0.4 < H − K s < 1.2 with little evolution with magnitude (and therefore redshift).
The observed colors may be affected by the presence of dust in the AGN or in its host galaxy, whose effect depends not only on its column density, but also on the redshift of the QSO, since our images sample progressively bluer wavelengths at higher redshift.As a reference, Fig. 3 shows the extinction arrows for A V = 1 for three different values of redshift.
The 14 components of the 7 LBT systems span a large part of this diagram.Information on the nature of each component can be derived from its position in this color-magnitude diagram, i.e, by comparing it with the distribution of stars and QSO.A very small fraction of high-galactic latitude stars have colors H − K s ≳ 0.35, therefore objects with significantly redder colors can be identified as QSO, while object with bluer colors are most likely to be stars.For all systems, except for J0812+2007, we do not know which one dominates the optical spectrum and is thus responsible for the classification as a QSO.
One component of J0812-0040, J1318-0136, and J1510+5959 (systems 2, 6, and 7 in Fig. 3) show combinations of very blue colors and H-band magnitudes that are consistent with stars, while the other component has red H − K s colors typical of AGNs.These systems are therefore better described by AGN/star alignments, occurrences expected to be present in our sample at the 30% level for separations larger than 0.5 ′′ (Mannucci et al. 2022).It is very unlikely, although not totally impossible, that these systems are constituted by two AGNs.
J0732+3533 and J0812+2007 (systems 1 and 3) have one component compatible with being an AGN, and the other component consistent with both AGNs and stars.As a consequence, colors cannot provide a classification of these two objects.
Both components of J0950+4329 (system 5) have colors that are only compatible with QSOs, and therefore this is either a dual or a lensed AGN.This is a very bright object at z = 1.77, at the upper luminosity envelope of the bright QSOs at this distance in Croom et al. (2009), and several magnitudes brighter than the typical COSMOS QSOs at that redshift (Civano et al. 2016).For this reason, J0950+4329 is most likely a lensed AGN, whose luminosity is boosted by lensing magnification.The slightly different observed H − K s colors of the two components could be due to different amounts of extinction in the lensing galaxy.Given the separation and the redshift of the system, the lensing galaxy is expected to have H ∼ 19−20 and therefore would be outshined by the lensed QSOs.
Finally, both components of J0927+3512 (system 4) show very blue H − K s colors; the brighter component is only compatible with being a star, while fainter one could be either a star or  (Girardi et al. 2005) and from the SDSS, respectively, while orange and red contours are AGNs from SDSS DR16Q (Lyke et al. 2020) and COSMOS (Civano et al. 2016), respectively.a very bright (for its tabulated redshift of 1.149) and very blue AGN.Magnitude and colors of this system are more compatible with being a double star, despite its classification as an AGN based on the LAMOST spectrum.This point is further discussed in Sec.3.2.
Objects fully resolved by Gaia, i.e., having a value of G magnitude for each component, can also be classified in a color-color plane involving both optical and near-IR magnitudes.This would offer the advantage of sampling a larger wavelength range and would provide a better classification of any companion star.In our sample, only J0812+2007 is among these resolved targets, but the large error of its H-band magnitude prevents us from deriving any additional useful information.
In Sect.3.4 these classification based on resolved colors will be compared with those obtained from both total and spatially resolved spectroscopy derived in the next sections.

Classification using unresolved spectroscopy
Ground-based spectra are primarily used to discover AGNs and measure their redshifts.Most of the available QSO spectra come from the SDSS DR16Q (Lyke et al. 2020), 2QZ (Croom et al. 2004) and LAMOST (Yao et al. 2019) surveys.As shown in Sec. 2, in GMP-selected systems the luminosity ratio between the primary and the secondary object is usually not large, and it is below a factor of six for all objects with G>18.5.For this reason, ground-based spectra, which are usually blend of all the components, can also be used to identify the systems containing a star if they can be reproduced with the superposition of an AGN and a star.
Following the procedure described in Scialpi et al. (2023), we have performed the deconvolution of the spectra of all the targets in Tab. 1.The observed spectra are fitted with a combination of an AGN (Vanden Berk et al. 2001;Temple et al. 2021) and a stellar template (Covey et al. 2007) .The fitting parameters, determined with a χ 2 minimization, are the normalization of the two spectra, the redshift of the AGN, the radial velocity of the star (limited to ±300km/sec), and the level of dust extinction of the AGN template.LBT observations were obtained without performing any spectroscopic analysis, thus the LBT sample is unbiased towards systems constituted by 2 AGNs (either dual or lensed) and AGN/star associations.In contrast, Keck observations were targeted towards systems with no clear stellar features in the blended ground-based spectra, except J0812+2007 that was observed because already included in the LBT sample.Also the ERIS target J1318-0136 was observed because part of the LBT sample.
The results of the deconvolution for the objects in Table 1 are reported in Table 5 and plotted in the appendix.They show that: -J0732+3533, J0950+4329, J1048+4541, and J1103+2348 can be fitted well by a single AGN spectrum and do not show any stellar feature.-In contrast, J0812+2007, J0812-0040 and J1510+5959 show clear stellar features, revealing the presence of a projected star.-J1318-0136 has a low-S/N, uncalibrated spectrum from 2QZ survey (Croom et al. 2004), with problematic features.Our estimate of the redshift is based on the presence of a line at ∼ 3850Å, identified as CIV]λ1549 based on the presence of two low-significance features at 4800 Åand 6950 Å, possibly identified as CIII and MgII.Even if our estimate is in agreement with that originally provided by 2QZ, the redshift determination is not certain.Also, the quality of the spectrum does not allow us to look for stellar features, therefore this total spectrum does not provide any classification.-J0927+3512 is a peculiar object.It was observed by LAM-OST and classified as a QSO at z=1.149 based on a broad feature at 4050 Å, nevertheless this feature has low S/N and no other convincing emission line is present.Also, the MgII line expected at ∼ 6000Åfor this redshift is not present, and the continuum shape shows a peculiar shape, with a minimum in F λ around 7500Å.This spectrum does not convincingly show that an AGN is at all present in this system, therefore J0927+3512 remains unclassified.4. To optimize the visualization, some of the spectra have been multiplied by arbitrary factors, which are reported in the top-left corner of the corresponding panels.Vertical dotted lines show the position of the main expected emission lines.The appendix reports the spectra and their spectral deconvolution.The best fit is shown in red, the QSO in blue, the stellar emission in orange, with pink vertical stripes of width 50Å covering the main stellar absorption features of the best-fitting stellar type.

Classification using resolved spectroscopy
We observed four GMP-selected systems with the integral-field spectrograph (IFS) OSIRIS on the Keck telescope, and one target with ERIS/SPIFFIER on the VLT.On both cases the use of laser-assisted adaptive optics allowed us to fully resolve the sys-tems and obtain independent spectra for the two components.These sources were selected to have FMP≥ 10, a more conservative value then the minimum FMP=8 derived in Sec. 2. We extracted the targets from the Milliquas v7.9 catalog (Flesch 2021), considering known AGNs with a spectroscopic redshift such as to put either Hα or Hβ inside one of the available near-IR bands (0.85 < z < 1.11, 1.28 < z < 1.85, and 2.03 < z < 2.65 for OSIRIS; 0.70 < z < 1.02, 1.27 < z < 1.73, and 1.98 < z < 2.73 for ERIS).We selected targets at high galactic latitude (b > 20 • ) and near a star bright enough to drive the AO system.We selected five targets, three of which (J0812+2007, J0950+4329, and J1318-0136) are in the LBT sample described in Sec.3.1.1.Table 3 summarizes the obtained observations, which are described in the next sections.

Keck/OSIRIS IFU spectra
Keck spectra were obtained on on February 8th, 2023, using a pixel scale of 50 mas.In all cases the AO module was driven by a laser guide star (LGS) and a nearby natural tip-tilt (TT) star to correct the low orders.The positions of the two components in J0812+2007 and J0950+4329 were known before the observations from the Gaia catalog and the LBT observations.The knowledge of the position angle in these two cases allowed us to use a smaller field-of-view (0.8 ′′ × 3.2 ′′ ) with a broader wavelength coverage (OSIRIS filters Hbb and Jbb, respectively).In contrast, J1048+4541 and J1103+2348 are isolated sources in the Gaia database and the position of the companion was not known before the observations.As a consequence, we used a narrower wavelength coverage (Hn3 on both cases) to have a larger field-of-view (1.6 ′′ × 3.2 ′′ ), well matched to the expected separations below 0.7 ′′ .In addition to the science targets, each night we also observed a standard star of spectral type A for telluric calibration.All data cubes were assembled and reduced using the standard OSIRIS pipeline (Lockhart et al. 2019), see Ciurlo et al. (2023) for more details.Fig. 4 shows the images of the observed systems (left panels) and the spectra (right panels) obtained via a weighted sum in the shown apertures.
All the systems show the presence of two point-like, unresolved sources (see Fig. 4).J0812+2007 is constituted by an AGN, revealed by the broad Hα line, and a second object characterized by a featureless continuum.The most probable interpretation for the latter is a star, which is expected in about 30% of the systems (Mannucci et al. 2022).Both components of J0950+4320 show QSO spectra, characterized by broad Hβ and Hγ and bright [OIII]4959,5007 lines.The two spectra become virtually identical by applying a small correction for dust extinction to component A for A V = 0.4, and a normalization factor of 2.65.The similarity of the spectra and the luminosity of the object (see Fig. 3) suggest that this is a lensed QSO, with the different dust extinction produced by the different light paths for the two images.The lensing galaxy is not detected but this is compatible with the relatively low sensitivity of our spectra to extended and continuum-dominated objects close to two bright point sources.Both J1048+4541 and J1103+2348 are constituted by two components each, and all the spectra show the presence of bright, very broad Hα lines.The line profiles in each pair are very different, both in center and shapes.To quantify the difference, we fit a Gaussian line and a constant continuum to each broad emission line.The formal errors on both line center and FWHM are small, a few km/sec, but the real uncertainties, due to the non-Gaussian profiles of the line and the partial coverage of some of the lines, are larger.For this reason we estimate typical uncertainties of ∼100-200 km/sec on the centers and 200-500 km/sec on the FWHM.The Hα central wavelengths of the two components of J1048+4541 differ by 500 ± 300 km/s, and the two lines have very different line widths, i.e.FWHM = 2500±200 km/s for component A and FWHM=5300±500 km/s for component B. The significance of the velocity offset is hampered by the lack of most of the blue part of the line of component B, while the difference in line width is highly significant.The two components of J1103+2348 have similar line widths (FWHM=2950 ± 400 km/s in both cases) but the central wavelengths are offset by 430 ± 200 km/s.In addition, the line shapes are very different, with opposite asymmetries.These differences are not compatible with two lensed images of one variable source but with different delays: indeed, in both cases the expected time delay between the two components is a few days at most, while the luminosities of the objects imply BLR sizes of hundreds of days (Bentz et al. 2013).

VLT/ERIS IFU spectra
J1318-0136 was observed on April 14th 2023 with the new Enhanced Resolution Imager and Spectrograph (ERIS, Davies et al. 2023) on the ESO/VLT, during the first run of the INAF GTO time on this instrument, using a LGS and a nearby TT star.We used the SPIFFIER IFU spectrograph (George et al. 2016), upgrade of the SPIFFI spectrograph Eisenhauer et al. (2003), covering a field of view of 3.2 ′′ × 3.2 ′′ with 50 mas sampling.We integrated for 30min using the H-band grating, which provides a spectral resolution of R∼5200 in the wavelength range between 1.45 µm and 1.85 µm.The data were reduced with the automatic pipeline4 .The PSF has FWHM=120 mas, limited by the 50 mas spaxel scale.
In Fig. 5, we show an image of the system, showing the presence of two, well-separated components at 0.27 ′′ of separation, and the spectra of the two sources extracted using circular apertures of 200 mas-diameter, after removing a second-pass background in each spectral plane of the cube.The bright Hα emission line centered at 1.6375 µm corresponds to redshift z = 1.495 and is in good agreement with the tabulated value of z = 1.486 obtained from rest-frame UV lines, clearly revealing the AGN nature of the brightest component A. On the contrary, component B shows no emission lines and a few absorption features, revealing its nature as a foreground star.
In conclusion, Keck and VLT spectroscopy has allowed us to reliably classify five GMP-selected QSOs: two are AGN/star projection, one is a lensed systems, and two are dual AGNs, with the properties reported in Table 4.In the next section these classifications will be compared with those obtained from the near-IR colors (Sect.3.1) and the total spectra (Sect.3.2).

Comparison between different classifications
The three independent classification based on near-IR colors, integrated spectra, and spatially resolved spectroscopy are compared in Table 5.It appears that there is an excellent agreement between all available classifications: -J0812-0040 and J1510+5959 are classified as AGN/star associations by the available methods, i.e, colors and integrated spectra; -J0812+2007 is an AGN/star systems based on the Keck spectrum, in agreement with the integrated spectrum.-The spectrum of J0732+3533 does not show any stellar contribution, while the colors are compatible with both AGN/AGN and AGN/star systems.-J1318-0136 is classified as an AGN/star alignment based on its near-IR colors, while the spectrum is of low quality and cannot be used for a reliable classification.The ERIS spectra confirm this classification.-Finally, J0927+3512 is a peculiar system: while the LAM-OST survey classifies this object as a secure QSO at z = 1.149, our spectrum does not show any clear emission line, and the MgIIλ2800 line expected at 6017 Å is not detected.Thus, the presence of an AGN in this system is not secure, let alone its redshift.Both components show very blue H − K s colors that, combined with the high apparent magnitude (H<15.5)point toward a classification as a star pair.In contrast, its spectrum is not a clear superposition of two stars.The nature of this system remains uncertain.
Summarizing, while spatially-resolved spectra are the ultimate tool to secure a classification, combining integrated spectra, magnitudes, and colors is a useful method to classify the GMP systems with relatively high confidence.

Discussion and conclusions
In Sec. 2 we have shown that the GMP method has a simple and well defined selection function in terms of the separation of the Projected Separation (Kpc) 0.15" (GMP limit)  2023), and this work.This data compilation can be found in computer-readable format in https://tinyurl.com/dualagnssources and their luminosity ratios.This is a key feature to compare the properties of the discovered dual and lensed systems with the theoretical predictions, and to test, for the first time, some of the aspects of the models of hierarchical formation of galaxies and SMBH.
Based on integrated spectra and AO-assisted observations at VLT, Keck and LBT we have shown that combining integrated spectra and spatially-resolved near-IR luminosity and colors can provide indications on the nature of these systems as dual AGNs, lensed QSOs, or AGN/star alignment.Mannucci et al. (2022), Ciurlo et al. (2023), Scialpi et al. (2023), and this work present samples of dual AGNs selected by the GMP method and classified by resolved spectroscopy.In total, 22 GMP systems were classified, obtaining 10 dual AGNs, 5 lensed systems (including the system classified as a lens by Li et al. 2023), and 7 AGN/star associations.One of these systems and three more GMP candidates have been previously identified using other techniques: Junkkarinen et al. (2001) reported the serendipitous discovery with HST of a very compact dual system at z ∼ 0.84, with two components separated by ∼ 0.3 ′′ .Schechter et al. (2017) identified a dual system at z = 1.34 and δ = 0.71 ′′ by spatial deconvolution of ground-based images of AGNs selected to have red WISE W1-W2 colors.Chen et al. (2022b) identified several multiple AGN candidates by looking for SDSS systems associated to more than one Gaia object within a few arcsec.One of these systems was spectroscopically confirmed to be a dual AGN at z = 2.95 and separation 0.46 ′′ (Mannucci et al. 2022).As shown in Fig. 1, this technique becomes very inefficient at δ < 0.5 ′′ , but sporadically can provide targets down to δ ∼ 0.3 ′′ .Finally, Chen et al. (2023) present a detail analysis of one target at z = 2.17 and separation δ = 0.46 ′′ selected using "varstrometry", i.e" the Gaia astrometric extra jitter due to the unrelated variation of the two AGNs.
Fig. 6 shows the distribution in redshift and projected separations of all the known dual AGNs with separations below 30 kpc.It is evident as all the 14 systems with separation below 7 kpc are GMP-selected systems, and that 9 of these are only detected by this method.The only exceptions are the complex system discovered by Lemon et al. (2022) at z=2.36, which is a lensed system with 6 components with large separations (up to ∼ 4 ′′ ), possible attributed to a dual systems with separation of ∼ 0.7 ′′ , a very interesting yet uncertain and peculiar systems, and the dual system at z=1.889 and δ = 2.2kpc serendipitously discovered by Glikman (2023), whose luminosity ratio between the components is too large to be selected by GMP (see Fig. 2).
In conclusion, the GMP method appears to be the best suited technique to provide a large sample of high-redshift dual AGNs residing in the same host galaxies.We are currently in the process of confirming more systems with AO spectroscopic observations, in particular with VLT/MUSE, VLT/ERIS, and Keck/OSIRIS, to obtain a sample large enough to test the prediction of models of galaxy formation and evolution. Summarizing: -We have tested the performances of the GMP selection using dense stellar fields with numerous projected pairs at low separations.This allowed us to estimate the level of contamination (false positives) as a function of the selecting parameter FMP (Fig. 1, left) and the selection efficiency as a function of separation δ (Fig. 1, left) and source magnitudes (Fig. 2).The results show that low levels of contamination and high selection efficiency can be obtained for δ > 0.15 ′′ and G<20.5; -AO-assisted imaging in the H and K s band with SOUL/LUCI at LBT of 7 GMP-selected systems has allowed us to resolve all of them into multiple components, measure the H − K s of all the sources, and to develop a new color-based method to identify the most probable nature of the GMP-selected systems (Fig. 3); we reproduced the integrated, ground-based spectra of the 9 systems of this sample with a combination of an AGN template with a foreground star.Four of these 9 systems are well reproduced by a simple AGN spectrum (Fig. A.1), three show clear evidence for the presence of an AGN and a star, A.2) and two remain unclassified , A.3. we present the spatially-resolved spectra of 5 systems with Keck/OSIRIS (Fig. 4) and VLT/ERIS (Fig. 5).These spectra has allowed us to discover two dual systems with separations of 0.30 ′′ and 0.52 ′′ , a very compact lensed system (separation 0.25 ′′ ), and two AGN/star associations; comparing the various classification methods, we have shown that the combined analysis of resolved colors and integrated spectra, can provide reliable classifications of the systems in most cases.-Fig.6 shows all the confirmed dual AGNs at z > 0.5 and separations below 30 kpc.Fourteen of the fifteen systems with separations below ∼ 7 kpc are GMP systems, and 10 of them are originally selected by this method.In other words, since its introduction (Mannucci et al. 2022), this techniques has already allow us to triple the confirmed systems at small separations.This result shows that the GMP technique is an important step forward to sample the dual AGNs hosted by the same galaxy and test the predictions of the models at δ < 7 kpc.

Fig. 1 :
Fig.1: Left: Level of contamination of false positives (objects with detection of multi peaks that are not projected pairs) as a function of the FMP values in a sample of dense stellar fields, for G< 20.5.No false positives are found for FMP≥ 8, where we plot 90% upper limits, i.e., the levels that would correspond to 2.3 false positive.Right: fraction of recovered projected pairs as a function of separation (δ) in arcsec.The red histogram shows the pairs that are not recovered, which dominate the counts at δ < 0.15 ′′ and are very rare at δ > 0.20 ′′ .The green histogram reports the fraction of recovered 'isolated' objects that correspond to single entries in the Gaia source catalog, and the blue one shows the GMP-selected objects that have nearby companions in the catalog.The combination of these two GMP populations is show by the black line.

Fig. 2 :
Fig.2: Efficiency of the GMP method as a function of G mag of the primary and primary and secondary components, in three separation ranges.The efficiency does not show strong dependencies on the three parameters up to secondary G magnitudes of 20.5, shown by the dashed line.
Fig.3: Near-IR colors and magnitudes of the observed targets.Blue and purple contour plots show the distribution of Galactic stars from the TRILEGAL model(Girardi et al. 2005) and from the SDSS, respectively, while orange and red contours are AGNs from SDSS DR16Q(Lyke et al. 2020) and COSMOS(Civano et al. 2016), respectively.

Fig. 4 :
Fig. 4: Hα emission line maps (left) and spectra (right) of the systems observed with OSIRIS (target name and redshift reported in the right panels).The line maps are oriented with North to the top and West to the right.The spectra shown on the right panels have been extracted over the squared apertures marked on the left panes (with the same color-coding).Each component of the systems is labelled as in Table4.To optimize the visualization, some of the spectra have been multiplied by arbitrary factors, which are reported in the top-left corner of the corresponding panels.Vertical dotted lines show the position of the main expected emission lines.

Fig. 5 :
Fig. 5: Image (left) and uncalibrated H-band spectra (right) of J1318-0136, obtained with VLT/ERIS.The image has 1 ′′ × 1 ′′ size and is oriented with North to the top and West to the right.The spectra on the right are extracted using circular apertures of 200 mas diameter shown on the image.The bright emission line is the Hα of AGN, while the other component only shows absorption features.

Fig. A. 2 :
Fig. A.2: Spectral deconvolution of the systems best reproduced by a combination of an AGN and a star.Observed spectra are shown in black, the best fit in red, the QSO in blue, the stellar component in orange, with red vertical stripes around the main stellar absorption features of the best-fitting stellar type.The four small panels above the main ones are enlargements of the main absorptions features of the stellar companion.

Fig. A. 3 :
Fig. A.3: Spectral deconvolution of the unclassified objects, i.e., the systems whose spectra are not well reproduced neither from an AGN, nor from an AGN plus a star.

Table 2 :
. Each system was observed in the near-IR H and K s Main properties of the targets observed with LBT.filters with integration times of 6 min in H and 24 min in K s in the March run, and 15 min each in H and K s for the June run.

Table 3 :
Main properties of the five targets observed with OSIRIS/Keck and ERIS/VLT, together with the observational setup.Redshift are as reported in the Milliquas catalog.Bands are reported as in the instrument user manuals The table also reports the FWHM of the PSF, calculated on isolated sources.

Table 4 :
Summary of the results from VLT?ERIS and Keck/OSIRIS observations: most probable classification, projected angular and linear distances from the brightest object, and center of the observed lines.

Table 5 :
Classification of all systems presented in this paper, based on different methods: deconvolution of integrated archival spectra, near-IR colors (from LBT images), and IFU spectroscopy (with Keck or ERIS).