The Gaia Ultracool Dwarf Sample – IV. GTC/OSIRIS optical spectra of Gaia late-M and L dwarfs

As part of our comprehensive, ongoing characterisation of the low-mass end of the main sequence in the Solar neighbourhood, we used the OSIRIS instrument at the 10.4 m Gran Tele-scopio Canarias to acquire low-and mid-resolution (R ≈ 300 and R ≈ 2500) optical spectroscopy of 53 late-M and L ultracool dwarfs. Most of these objects are known but poorly investigated and lacking complete kinematics. We measured spectral indices, determined spectral types (six of which are new) and inferred effective temperature and surface gravity from BT-Settl synthetic spectra fits for all objects. We were able to measure radial velocities via line centre fitting and cross correlation for 46 objects, 29 of which lacked previous radial velocity measurements. Using these radial velocities in combination with the latest Gaia DR3 data, we also calculated Galactocentric space velocities. From their kinematics, we identified two candidates outside of the thin disc and four in young stellar kinematic groups. Two further ultracool dwarfs are apparently young field objects: 2MASSW J1246467 + 402715 (L4 𝛽 ), which has a potential, weak lithium absorption line, and G 196–3B (L3 𝛽 ), which was already known as young due to its well-studied primary companion.


INTRODUCTION
Ultracool dwarfs (UCDs) are objects with effective temperatures  eff ⪅ 2700 K (spectral type ≳ M7 V, Kirkpatrick et al. 1999) con-tinuing on from the low-mass tail of the main sequence, that consist of spectral types late-M, L, T and Y dwarfs.These UCDs consist of a combination of low-mass stars and brown dwarfs.Brown dwarfs are sub-stellar objects incapable of hydrogen fusion and are defined by mass, between the deuterium minimum mass burning limit, ∼13 Jupiter masses (Saumon et al. 1996;Chabrier et al. 2000) and the hydrogen minimum mass burning limit, ∼72 Jupiter masses (Chabrier & Baraffe 1997; Baraffe et al. 1997).The majority of known UCDs are within the Solar neighbourhood (e.g.Smart et al. 2019;Kirkpatrick et al. 2021;Sarro et al. 2023) with typically dim apparent optical magnitudes (Gaia  ⪆ 17 mag).The closest stars to the Sun have been catalogued throughout the history of astronomy.For example, the Catalogue of Nearby Stars (CNS) from Gliese (1957) has been updated with every all-sky photometric and astrometric survey, including the most recent release using Gaia DR3 data (CNS5, Golovin et al. 2023).This Solar neighbourhood has been further described in the 'The Solar Neighborhood' series by the Research Consortium on Nearby Stars (RE-CONS 1 ) team with publications from Henry et al. (1994) to Vrĳmoet et al. (2022).Specifically, M dwarfs within 30 pc were covered in another series of articles from Delfosse et al. (1999) to Crifo et al. (2005).Volume limited samples such as the recent Gaia Collaboration et al. (100 pc, 2021b), Kirkpatrick et al. (20 pc, 2021) and Reylé et al. (10 pc, 2021) works provide important constraints on the initial mass function (Salpeter 1955;Scalo 1986;Kroupa 2001;Chabrier 2003), which underpins all aspects of astrophysics from stars to galaxies to cosmology.
Spectral features of low mass stars, M, L and T dwarfs, and their definitions were initially described by Tinney & Reid (1998), Kirkpatrick et al. (1999), Martín et al. (1999b), Burgasser et al. (2002), Geballe et al. (2002) and Kirkpatrick (2005).The bulk of the flux emitted by L dwarfs lies in the near infrared (NIR) and continues strongly towards the mid-infrared spectral regions for later spectral type UCDs.However, several features of youth, e.g. a weak sodium doublet, 8183,8195 Å (Schiavon et al. 1997a), are apparent in mid-to high-resolution optical spectra.Additionally, in the optical regime features such as the 9850-10200 Å FeH Wing-Ford band (Schiavon et al. 1997b) can be seen, which can be indicative of low or high metallicity.Optical spectra have an advantage in that there are fewer and weaker telluric absorption bands than in ground-based infrared spectra, where water and oxygen bands can dominate (Reiners et al. 2007;Smette et al. 2015).However, only the closest and brightest UCDs can be observed with optical spectroscopy due to the low relative flux; further and fainter UCDs require large aperture telescopes and long exposure times.
The photometry of UCDs is important because the change in colour across the optical and NIR regime (Leggett et al. 2002) correlates with physical and atmospheric properties.These changing processes, such as dust, condensate cloud formation and subsequent clearing as an atmosphere cools, are well covered in the literature (e.g.Marley et al. 2002;Dahn et al. 2002;Saumon & Marley 2008).Understanding a changing atmosphere for different ages with a range of masses has allowed the computing of 'cooling tracks' (Burrows et al. 1997; Baraffe et al. 2015).Accounting for theoretical atmospheric physics has been used in model grids such as BT-Settl (Allard et al. 2011), or Sonora (Marley et al. 2021;Karalidi et al. 2021), and when interpreting the results of retrieval techniques (e.g.Burningham et al. 2017;Calamari et al. 2022).Being able to constrain the mass and/or age has underpinned modern observational UCD astronomy, but is challenging due to the mass/age degeneracy (Burrows et al. 1997).For example, benchmark systems (e.g.Pinfield et al. 2006;Dupuy et al. 2009) allow us to constrain the age of a brown dwarf with the coeval main sequence primary.The metallicity and surface gravity of an object of a given spectral type are the major variables affecting the photometric colour (Stephens et al. 2009), see references to 'blue' and 'red' L dwarfs (e.g.Faherty et al. 2009;Schmidt et al. 2010).Any works that infer spectral type, surface gravity and effective temperature must take into account the atmospheric physics, as these directly correlate with observable features.
Gaia is a European Space Agency mission, launched in 2013 to make high-precision measurements of positions, parallaxes, and proper motions of well over a billion sources and photometry in three different photometric filters ( BP , ,  RP ).The third Gaia data release (EDR3 and DR3 -Gaia Collaboration et al. 2021aCollaboration et al. , 2023, respectively) , respectively) containing astrometric and photometric measurements, was in December 2021, with the remaining measurements and inferred parameters, including spectra, in June 20222 .
Obtaining the full 6D (right ascension, declination, proper motions, parallax, radial velocity: , ,   cos ,   , ,   ) positional and kinematic information is fundamental to fully characterise the populations of UCDs within a volume limited sample (e.g.Best et al. 2021).Precise measurements of radial velocities (RVs) are obtained from high signal-to-noise observations taken with high resolution spectrographs with resolving powers of R∼100 000, leading to uncertainties ∼1-5 m s −1 .This has only been achievable for the nearest, brightest UCDs (e.g.Zechmeister et al. 2019).Blake et al. (2010) achieved   ≈ 50-200 m s −1 with the Keck Near-Infrared Spectrometer (NIRSPEC), which had a resolution of R≈25 000.The 'Brown Dwarf Kinematics Project' has gathered further UCD RVs (Burgasser et al. 2015;Hsu et al. 2021) with both the NIRSPEC and the Magellan Echellette (MagE, R∼4100,   ≈ 2-3 km s −1 ) spectrographs.By comparison, the lower-resolution spectroscopy such as those discussed in this work (R≈2500) is only capable of theoretical minimum uncertainties of ≳5 km s −1 ; this is still useful when constraining the kinematics of the Solar neighbourhood.Parallaxes and proper motions of UCDs were historically gathered from ground based time-domain campaigns (e.g.PARSEC: Andrei et al. 2011;Marocco et al. 2013;Smart et al. 2018) that have been generally superseded by Gaia for the brightest objects,  ⪅ 20 mag.In the case of most late-L and T dwarfs, ground-based astrometry is still the predominant source (e.g.Vrba et al. 2004;Dupuy & Liu 2012;Liu et al. 2016;Best et al. 2018).For dimmer objects, beyond mid-L dwarfs, parallaxes and proper motions are gathered by space-based infrared surveys and are analysed in-depth by Kirkpatrick et al. (2021).Young moving groups are constrained using these complete kinematics.See the BANYAN Σ series and references therein for detail on nearby young moving groups and clusters (Gagné et al. 2014, to Gagné & Faherty 2018) or similarly, the LACEwING code (Riedel et al. 2017), designed around young objects in the Solar neighbourhood.Subdwarfs, meanwhile, are characterised by their statistically higher space velocities indicative of the older population (e.g.Lodieu et al. 2005;Burgasser et al. 2007;Lodieu et al. 2017;Zhang et al. 2017).This is the fourth item in the Gaia UltraCool Dwarf Sample series (GUCDS, Smart et al. 2017Smart et al. , 2019;;Marocco et al. 2020) and is an ongoing, international, multi-year programme aimed at characterising all of the UCDs visible to Gaia.Gaia DR3 produced astrophysical parameters for ≈470 million sources (Creevey et al. 2023), including effective temperatures,  eff .The ≈94 000 Gaia DR3  eff values relating to UCDs by Creevey et al. (2023) were provided under the teff_espucd keyword.The full sample of UCDs detected by Gaia with Gaia DR3  eff values were documented and analysed by Sarro et al. (2023).In our analysis, we will use the values from these Gaia DR3 derivative works to compare with the equivalent values directly measured by us.There is significant overlap between the Sarro et al. (2023) sample and the GUCDS, although the majority of UCD sources as seen by Gaia are as yet not characterised through spectroscopic follow-up.A subset of this Sarro et al. (2023) sample has public Gaia RP spectra (see the Gaia xp_summary table 3 ), which covers the  RP passband (Δ ≈ 6200-10420 Å, Riello et al. 2021).This subset from Sarro et al. (2023) was further analysed for spectroscopic outliers by Cooper et al. (2024).The internally calibrated Gaia RP spectra and processing were discussed thoroughly by Carrasco et al. (2021), De Angeli et al. (2023) and Montegriffo et al. (2023).
The aim of this work is to complement the literature population with measurements and inferences from low-and mid-resolution optical spectroscopy.In Section §2 we explain the target selection ( §2.1) and observation strategy ( §2.2).Different reduction techniques with a test case are discussed in Section §3.Section §4 explains our techniques for determining spectral types ( §4.1), astrophysical parameters ( §4.2), and kinematics ( §4.3) including membership in moving groups ( §4.4).Section §5 follows a discussion of our results for spectral types ( §5.1), kinematics ( §5.2) and astrophysical parameters ( §5.3).We also discuss individual objects ( §5.3.1)before summarising the overall conclusions in Section §6.

DATA COLLECTION
We obtained optical spectroscopy of 53 unique UCDs using the OSIRIS (Optical System for Imaging and low-intermediate Resolution Integrated Spectroscopy -Cepa 1998) instrument on the 10.4 m Gran Telescopio Canarias (GTC) at El Roque de los Muchachos in the island of La Palma, Spain, under proposal IDs GTC54-15A and GTC8-15ITP (PIs Caballero and Marocco, respectively).The objects were observed in semesters 2015A, 2015B and 2016A.
The observed data from the GTC were complemented with Gaia DR3.Gaia also carries a radial velocity spectrometer, although this was unsuitable for our purposes as all of our targets were fainter than the Gaia selection limit (Katz et al. 2023, 𝐺 < 14 mag,).We acquired 63 spectra in which we observed 53 unique objects, shown in Table 1.These 63 observations are shown in Table A1, including the air mass and humidity of the observation.Of the 63 spectra, 46 were observed with the R2500I volume phased holographic grating (hereafter VPHG), whilst 17 were observed with the R300R grism.Ten of the 53 objects were observed with both dispersive elements.
Twenty of the 53 objects already had full 6D positional and kinematic information in the literature.Fifty-one had proper motions, 43 had parallaxes, and two had only  and .All values along with their provenance are given in Table 1.In the next sub-sections we discuss the target list selection and observations.

Target selection
Our targets were drawn from a combination of two samples: benchmark systems (system with a star and a UCD, Pinfield et al. 2006) and known L dwarfs with poor or no available spectroscopy.The targets were selected by Marocco et al. (2017) and Marocco et al. (2020), and here we briefly summarise their selection criteria.Both samples were chosen with the aim of gathering low-and midresolution spectra, mostly to achieve radial velocities and to confirm their status as L dwarfs.Benchmark system selection used the procedure of Marocco et al. (2017, their section 4).To summarise, primary systems consisting of possibly metal-rich or metal-poor stars were selected with metallicity cuts of [Fe/H] < −0.3 and [Fe/H] > 0.2 dex from a number of catalogues (Marocco et al. 2017, their table 2).If more than one value of [Fe/H] was available, the one with the smallest uncertainty was used; Marocco et al. (2017) did not investigate if there were any systematic offsets between different catalogues, as this was beyond the scope of that work.The companions to these systems were filtered by a series of colour, absolute magnitude and photometric quality cuts from 2MASS, SDSS (the Sloan Digital Sky Survey, York et al. 2000) and ULAS (United Kingdom Infrared Telescope Deep Sky Survey, Large Area Survey, Lawrence et al. 2007) photometry in equation (1).These colour cuts in equation ( 1) are taken directly from Marocco et al. (2017) as that work created part of the target list used in this work.Magnitudes from 2MASS were converted into UKIRT/WFCAM magnitudes via the equations of Stephens & Leggett (2004).
These companions were determined as being candidate benchmark systems with a maximum matching radius of 3 arcmin, i.e. the maximum separation to the primary object.The remaining targets, known L dwarfs, were already spectroscopically confirmed bright L dwarfs that were predicted to be visible to the astrometry and photometry in (at the time, upcoming) Gaia data releases.These known L dwarfs should be single systems.They would, however, not be bright enough for the Gaia radial velocity spectrometer (Katz et al. 2023), and thus were chosen to determine radial velocities for, as a complement to the 30 pc volume-limited sample.This list was complemented with additional targets too dim for Gaia photometry and astrometry, which were detected in UKIDSS, and by a few wellknown L dwarfs, such as G 196-3B, which could serve as template standards.

Cross-matching
All observed targets (Table 1) were cross-matched with Gaia, 2MASS, and AllWISE.These surveys were chosen because they are all-sky and we were aiming for completeness in this process.The targets were also cross-matched with Pan-STARRS (50/53 successful matches), for the additional optical components for those sources within the Pan-STARRS footprint.This sample of 53 objects was then also cross-matched against the astrophysical parameter and xp_summary tables from Gaia DR34 .Thirty-eight of these objects had a teff_espucd value, and 28 had a public RP spectrum.Internally calibrated Gaia RP spectra were then extracted from the Gaia archive with a linearly dispersed grid from 6000 Å to 10500 Å using the gaiaxpy.convert(Ruz-Mieres 2022) and gaiaxpy-batch (Cooper 2022a) codes.We also searched for common proper motion systems within Simbad (Wenger et al. 2000) with the selection criteria given in the GUCDS, specifically equation (1) of Marocco et al. (2020): In equation (2),  is the separation in arcseconds,  is the proper motion position angle in degrees, whilst  (milli-arcseconds) and  (milli-arcseconds per year) are our target list's Gaia DR3 parallax and proper motion, respectively.Like with the photometric selection, equation (1), the common proper motion selection was taken directly from Marocco et al. (2020).This is because the target list in this work is drawn from the same wider target list used in the GUCDS.In effect, this selection is creating a widest possible physical separation of 100 000 AU (see the discussion on binding energies by Caballero 2009).

Observations
The OSIRIS instrument used a 2 × 1 mosaic of 2048 × 4096 pixel (photosensitive area) red-optimised CCDs (Marconi MAT-44-82 type) with a 7.8 × 7.8 arcmin 2 unvignetted field of view.We used the standard operational mode of 2 × 2 binning, which has a physical pixel size of 0.254 arcsec pixel −1 .For our purposes, we used the 7.4 arcmin long slit with a width of 1.2 arcsec.We had variable seeing between 0.6 and 2.5 arcsec, with the vast majority having seeing <1.2-1.5.The undersampling of the Full Width at Half-Maximum (FWHM) when the seeing is significantly less than the slit width would cause uncertainty in the wavelength calibration.In the worst cases, this can approach the resolution element.This was then included in the systematic uncertainty estimate on the radial velocities.We used the R300R and R2500I grisms and purely read off CCD 2 due to the instrument calibration module having a strong gradient from CCD 1 to 2 in the flat fields.The R300R grism has a wavelength range of ≈4800-10 000 Å with a dispersion of ≈7.74 Å pix −1 for a resolution of ≈350 whilst the R2500I VPHG has a wavelength range of ≈7330-10 000 Å with a dispersion of ≈1.36 Å pix −1 for a resolution of ≈2500, as per the online documentation 5 .Both dispersive elements experience an increase in fringing at wavelengths ≳9200 Å to ≥5 per cent.The R300R grism however, had second order light from 4800 to 4900 Å contaminating the 9600 to 9800 Å region.This is because standards, but not UCDs, have flux in the blue regime, hence affecting the flux calibration in the red regime.As a result, the R300R spectra were conservatively truncated to 9000 Å.Our standards were a selection of white dwarfs plus two well-studied bright main sequence dwarf stars, all with literature flux calibrated spectra and spectral types: Ross 640 (DZA6, Oke 1974;McCleery et al. 2020); Hilt 600 (B1, Hamuy et al. 1992Hamuy et al. , 1994)); GD 153 (DA1, Bohlin et al. 1995Bohlin et al. , 2014)); G191-B2B (DA1, Oke 1990;Bohlin et al. 1995Bohlin et al. , 2014)); GD 248 (DC5, Tremblay et al. 2011;McCleery et al. 2020), GD 140 (DA2, Tremblay et al. 2011;McCleery et al. 2020) andG 158-100 (dG-K, Oke 1990).We took a series of short exposures for the brightest objects to avoid saturation and non-linearity.The majority of observations had a bright moon whilst the sky condition varied from photometric to clear with humidity typically ≲50 per cent.All calibration frames were taken at the start and end of each night, the arc lamps being used to solve the wavelength solution were: Hg-Ar, Ne and Xe.The full observing log is given in Table A1.

DATA REDUCTION
We aimed to determine spectral types, spectral indices and radial velocities from directly measuring the GTC spectra.Furthermore, we inferred astrophysical parameters (effective temperature,  eff [K]; surface gravity, log  [dex]; and metallicity, [Fe/H] [dex]) from comparisons with atmospheric models.Our adopted PypeIt6 (Prochaska et al. 2020a;Prochaska et al. 2020b) reduction procedure applied to every object was as follows: master calibration files were created by median stacking the relevant flat, bias and arc frames.Basic image processing was performed including bias subtraction, flat fielding, spatial flexure correction and cosmic ray masking via the L.A. Cosmic Rejection algorithm (van Dokkum 2001).We then manually identified the arc lines using the median stacked master arc.These arc lines were used to manually create a wavelength solution through pypeit_identify with typical RMS values of ≈0.0804 Å for the R2500I VPHG and ≈0.1394 Å for the R300R grism.The R2500I wavelength calibration solution was a 6 th order polynomial, whilst the R300R solution was only 3 rd .The information inside the object headers (observation date, object sky position, longitude and latitude of the observatory) was used to heliocentric correct the wavelength solution.The PypeIt wavelength solution was defined in vacuum.
The standard frames were median stacked before the global sky was subtracted and corrected for spectral flexure (to account for fringing).Both the stacked standard and object were then extracted using both boxcar (5 pixel) and optimal (Horne 1986) extraction methods, with the latter being the presented spectra.
We then fitted a function to account for the sensitivity, CCD quantum efficiency and zeropoint.The telluric regions listed by Reiners et al. (2007) and Smette et al. (2015) were masked out.We divided each standard by its corresponding flux calibrated spectrum from the literature, as listed above.This sensitivity function was then applied to the reduced standard and object to flux calibrate the extracted spectra.If an observation had more than one science frame, those were co-added after wavelength and flux calibration.
The standards observed under the R2500I VPHG were used to create a telluric model from a high resolution atmospheric grid derived at Las Campanas.This grid was interpolated through to find the best match across airmass and precipitable water vapour.The telluric model was applied back to the flux calibrated standard and object.This telluric corrected standard was visually checked to confirm that the telluric model was behaving appropriately.The configuration files used in our reduction procedure are given in Appendix A6.
It is important to mention here that we made a comparison between this PypeIt reduction and that of a customised reduction (both the full basic image and spectral reductions) using standard IRAF tasks.This was done with the aim of validating the quality of the PypeIt data against that from a well proven reference source.In Appendices A2 and A3 we describe this procedure in detail for one suitably chosen test object from our selection sample, and which is common to both independent reductions: J1745−1640.
A comparison between the PypeIt reduction, and that which used standard IRAF routines, is shown in the normalised spectra of J1745−1640 in Figure 1.We show good agreement in the flux profile up to ∼8900 Å.The IRAF reduced spectra is brighter in the broad H 2 O region, due to the differing telluric correction methods.The MagE spectrum was not telluric corrected whilst the IRAF spectrum was telluric corrected using a blackbody, instead of Ross 640 (the corresponding white dwarf standard).This difference does not affect the model fitting of the spectra, as this is done in localised, small, chunks.All spectra then agree at wavelengths ⪆9800 Å.

ANALYSIS
Here, we discuss the analysis of the reduced spectra, in order to produce spectral types, astrophysical parameters and kinematics.We discuss our measurements of astrophysical parameters first because the cross-correlation technique used to measure RV requires the use of a best-fitting model derived template, obtained from the best fit of astrophysical parameters.The code used for both estimating astrophysical parameters and calculating RV is rvfitter (Cooper 2022b).This program was developed to effectively recreate in python older codes (e.g.IRAF.Fxcorr, IRAF.Splot, IDL.gaussfit) designed for allowing a user to manually cross-correlate spectra and fit line centres with different profiles.All wavelengths discussed in this Section are in standard air, hence we converted our PypeIt spectra from vacuum to air.This was performed via the specutils package, using the corrections by Edlen (1953).

Spectral typing
We spectral typed both the R300R and R2500I spectra using the classifyTemplate method of the kastredux (Burgasser 7500 8000 8500 9000 9500 10000 R2500I spectra for J1745−1640, normalised at 8100-8200 Å, comparing two independent reduction procedures: PypeIt in black and IRAF in orange.In blue, the heliocentric corrected MagE spectra (Burgasser et al. 2015) for the same object is shown (which is not telluric corrected).The Earth symbol indicates the telluric bands present in the spectra.
The spectra had all been smoothed in wavelength with a Gaussian 5 kernel, and we only compared the regions from 8000 to 8500 Å for R2500I and 7000 to 8000 Å for R300R.This was decided through experimentation, which deliberately excluded regions with telluric features, as those features can cause poorer solutions.Each object was also visually checked against known standards (Kirkpatrick et al. 1999), the spectral sub-types by which we refer to as 'by eye'.Any spectra with indicators of youth are given optical gravity classes as defined by Cruz et al. (2009), from , ,  in order of decreasing surface gravity.The kastredux spectral types were our adopted spectral types.

GTC spectral sequence
The 46 spectra from the R2500I VPHG, ordered by our adopted spectral type, are shown in Figures 2 and 3.All spectra are heliocentric corrected, such that the relative motion of the Earth has been removed.Each spectrum shown had an outlier masking routine applied such that points within a rolling ≈15 Å (ten data points) chunk are removed if they had a difference greater than the standard deviation from the median.Additionally, some objects had problematic O 2 A-band tellurics.In those cases, we interpolated over the region 7540-7630 Å from the maximum of the first ≈7.5 Å to minimum of the last ≈7.5 Å.Where appropriate, spectra were co-added.All spectra appear noisy in the primary H 2 O band of ≈9200-9600 Å.The 17 heliocentric corrected, reduced spectra from the R300R grism are shown in Figure 4.The R300R spectra were trimmed from 6500 <  < 9000 Å due to (a) the lack of signal in the blue regime and (b) to constrain to purely the first order light.Unlike the R2500I spectra, the R300R spectra were not telluric corrected.

Fundamental astrophysical parameters
We used the rvfitter.crosscorrelatecode on our R300R and R2500I spectra with BT-Settl CIFIST model grids from 1200 ≤  eff ≤ 4000 K and 4.5 ≤ log  ≤ 5.5 dex (Allard et al. 2011).Lower surface gravity grids were available but not routinely used as the focus was on RV measurement with an a priori expectation of field surface gravity, ≈5 dex.These models assume a solar metallicity with no variation and are linearly dispersed in steps of 100 K and 0.5 dex.This code allowed us to visually select the best fitting model from the array of model grids and for each spectral line from Table 2.
We used these chosen lines rather than correlating against the entire model because the models do not exactly match the flux profile of ground based spectra.It was also known that the BT-Settl grids were generated using a different line list to our selected alkali lines, taken from the NIST database (Kramida et al. 2021).For efficiency purposes, each model when being loaded into the code, was interpolated onto the wavelength array of the object being compared against.The models could optionally be Gaussian smoothed, The first 24 of the R2500I VPHG spectra with a linear offset applied, sorted by spectral sub-type.We show the short names and the spectral sub-types from this work, attached to each spectrum.At the top of the figure are grey lines denoting a selection of spectral features typical to L dwarfs, plus the two main telluric bands.Table 2.The list of atomic alkali metal lines used when estimating astrophysical parameters and calculating radial velocities.Wavelengths are as measured by Kramida et al. (2021) and are defined in standard air.
which was helpful for fitting any 'messy' regions of models (e.g. telluric bands in models with  eff ⪆ 2000 K).We normalised the model and data by their respective medians in a given variably sized 'chunk' around each spectral line.We noted that around certain lines, particular models appeared almost identical to each other, e.g.around 7000-8000 Å, the 1900 and 2000 K models are not visually distinct.This means there is a higher uncertainty for effective temperatures within the 1900-2000 K region.Not every spectral line was used for each object as some have poorly resolved features or low signal-to-noise.Our selected  eff was the mean  eff from each line measurement, as was log .To determine the error on each  eff and log  final value, we chose to use the standard deviation from each independent line fit divided by square root of the number of lines used.This error was added in quadrature with half of the separation between each grid, i.e.50 K for  eff and 0.25 dex for log .
Additionally, we created an 'expected' effective temperature,  eff , using the Filippazzo, sixth order field  eff relation (Filippazzo et al. 2015) and our adopted spectral types.The errors on  eff correspond with the mean difference in  eff across ±0.5 spectral subtypes (our spectral sub-type uncertainty), plus the quoted relation RMS of 113 K.

Calculating the radial velocities
Only our R2500I spectra were used to determine RVs as the features in R300R spectra are mostly blended/unresolved.We used two methods by which to measure an adopted RV: line centre fitting and cross correlation.We note that our seeing (Table A1, corrected for airmass) was almost always smaller than the slit width, which affects the RV offset as the slit is not fully illuminated.The full width at half-maximum was typically 3-4 pixels, corresponding to ≈0.75-1 arcseconds.Most observations were seeing-limited, whilst a few, taken in poorer conditions, were slit-limited.The following methods were performed only on heliocentric corrected spectra, hence any quoted RV values are heliocentric corrected.

Line centre fitting
Using the same atomic absorption lines listed in Table 2, we applied the rvfitter.linecenteringcode to interactively fit Gaussian, Lorentzian and Voigt profiles with the minimum possible width.This minimum possible width is equal to the number of free parameters plus one (although this does not guarantee a successful fit).We used these different profiles to obtain the best fit for a particular line given its underlying absorption characteristics and the available signal-to-noise of the spectral region.The fitting technique used was least-mean-square7 minimisation.For each spectral line, we subtracted a linear continuum from the data.The continuum corresponds to the medians of selected regions to the blueward and redward sides of the spectral line.Each continuum region is chosen to follow the shape of the spectra with a minimum width of ≈50 Å within 100-200 Å of the spectral line.Also shown during the fitting routine is a fifth order spline, as a visual aid; the minima of the spline does not necessarily correspond to the line position.A example of this routine is given for J1745−1640 in Figure 5.The fits were only accepted if they appeared to accurately represent the spectral lines profile upon visual inspection.In general, the most consistently reliable lines were the rubidium lines, sodium doublet and first caesium line.The potassium doublet often was affected by rotational broadening whilst the second caesium line was often affected by neighbouring features.The uncertainty for each line, was the value in the diagonal of the covariance matrix corresponding to centroid position from the least-squares fit, plus the wavelength calibration RMS for that object, Doppler shifted into RV space.
After measuring every line, we then calculated the overall weighted mean ( LC ) and weighted standard deviation ( LC ), the weights were the inverse of the uncertainties of each line used, squared.The uncertainty from the vacuum to air conversion was negligible (≪0.1 km s −1 ) compared to the fitting uncertainties calculated from the eight (or less, if rejected) aforementioned lines.The final line centre RV standard error was the weighted standard deviation divided by the square root of the number of lines fit.

Cross-correlation
In addition to estimating the astrophysical parameters with rvfitter.crosscorrelate in Section §4.2, we also used the same package to measure RV by manually shifting the best fitting BT-Settl model as a template.No smoothing was applied to the model template to match the spectral resolution of the object spectrum.This was because smoothing could confuse where the centroid of a line was, when looking by-eye.Likewise, there was no continuum subtraction applied to the object spectrum.The RV shift was in steps of 5, 10, 100 km s −1 , which in turn defined the RV uncertainty on each line (2.5, 5, 50 km s −1 , i.e. the margin of error).These RV errors are added to the wavelength calibration RMS for the given object (Doppler shifted into an RV error).Not all atomic lines were always used, only in the cases where the model appeared to closely match the apparent line profile.The typical technique was to select a broad region (Δ = 100-200 Å) around each spectral line, find the best fitting template in terms of  eff and log , then narrow that region (Δ ≈ 50 Å) to then find an RV.This was a predominantly by-eye technique, although root-mean-square deviation divided by interquartile range (RMSDIQR) values were computed as a numerical guide when comparing models.We also show a fifth order spline, as with the line centering method, as a visual aid.This initial broad region is shown for J1745−1640 in Figure 6.
As in Section §4.3.1, the overall cross-correlated weighted mean RV value ( XC ) and weighted standard deviation ( XC ) was calculated using all of the manually selected lines.We used the same method to estimate the uncertainty in final cross-correlation derived RVs as for the line centre results, by finding the standard error of the mean.

Adopted RV
We created an adopted RV by constructing a weighted mean, using the deviation in each method as the weighting.The different RV values for each line, method and the corresponding probability distribution functions (PDFs) are shown in Figure 7, for J1745−1640.We also note that our final adopted RV for J1745−1640 obtained from combining the results of the two measurement techniques (32.7 ± 6.5 km s −1 ) is in agreement with the values obtained from both the customised IRAF reduced data and the value reported by Burgasser et al. (2015), within their respective uncertainties.See Appendix A3 for a full description.
The adopted RV was the mean ( RV ) whilst the standard error ( RV ) was equal to the standard deviation ( RV ) divided by √ 2. The mean and standard deviation was calculated through the inverse variance weighting equations (5 and 6).Typically, we found that the cross-correlation technique was more precise (being more controlled by-eye) and robust.The line centre fitting was often more accurate, however, and performed best on the higher quality spectra.In the bottom panel, the orange curve is the cross-correlated PDF; the blue curve is the line centre PDF; and the black curve is the adopted PDF.The dotted vertical lines are the mean RV values as associated with each PDF.

Kinematics
Galactic UVW velocities were calculated using our adopted RVs plus Gaia astrometric measurements, using the equations from astrolibpy.We corrected for the Local Standard of Rest (LSR) using the values from Coşkunoǧlu et al. (2011) where U, V, W = (−8.50,+13.38, +6.49) km s −1 .These equations follow the work by Johnson & Soderblom (1987), except that U is orientated towards the Galactic anti-centre.We also used BANYAN Σ (Gagné et al. 2015a(Gagné et al. , 2018)), which provided moving group classification with associated probability.When using BANYAN Σ, we checked the resultant probabilities both with and without RV.This was because RV has by far the lowest precision, thus could reduce a likely membership candidate into a field object in error.Our final values are the ones which include RV.Notably, when using velocities in the Galactic reference frame, one can select a Galactic component with  total (where  total is the total space velocity,  total = √  2 +  2 +  2 ).We followed the work by Nissen & Schuster (2010) and define thick disc and halo objects as having  total > 70 km s −1 and  total > 180 km s −1 respectively.This definition, especially for separating the thin and thick disc, is indicative of metallicity; see the Besançon Galaxy models (Czekaj et al. 2014;Lagarde et al. 2021).

RESULTS
In this Section, we present the spectral types, radial velocities and astrophysical parameters.In Table A2, we provide photometry from the Gaia, 2MASS and ALLWISE catalogues.We discuss individually interesting objects and objects where our measured results differ significantly from published values.

Spectral types
In Table 3 we list published spectral types based on optical spectra, near-infrared spectra and the 'by eye' and kastredux methods discussed in Section §4.1.This work has produced the first spectral type estimates for six of the 53 objects.
The 47 objects with known spectral types have a standard deviation of 0.5 sub-types between the published values and the 'by eye'/kastredux results, which we adopt as the error on the new spectral sub-types.When the literature values for a given object differ we adopted the optical spectral type.Our spectral types across the two methods are displayed against the adopted literature spectral types in Figure 8.All objects except J1004−1318 have sub-type differences between the spectral type derived in this work and the adopted literature spectral type of less than two sub-types.J1004−1318, has an optical (Opt) spectral sub-type of L0 (Martín et al. 2010) whilst Marocco et al. (2013) found a sub-type of L1 using nearinfrared (NIR) spectra; we find a sub-type of L3.However, a more recent study, Robert et al. (2016), found a sub-type of L4 (NIR), which is more consistent with our result.The fit statistic from kastredux is about twice larger for L1 than for L3.In Figure 2, J1004−1318 does not seem dissimilar to the neighbouring objects, whereas the L0/L1 spectra appear different (e.g.weaker alkali lines).The different spectral typing of J1004−1318 may be due to lower signal-tonoise (S/N) ratios of some observations.For example, Martín et al. ( 2010) exposed for 2400 s at the 2.56 m Nordic Optical Telescope, while we exposed for 1500 s, and with moderately good seeing and low airmass, with a telescope with an aperture over 16 times larger.

Radial velocity analysis
We have derived RVs for 46 of the observed 53 objects, the seven objects that we did not measure RVs were only observed with the R300R grism.For 20 of the 53 objects, there are published RVs and for 17 of these we have measured RVs.The objects J1004+5022, J1441−0945 and J1617+7733B are candidate members of benchmark systems (Section §2.1), and we adopt the RVs of their primary stars as a comparison with our measured values for the secondary, for a total of 20 comparison RVs.In Figure 9, we plot histograms of the 20 published and the 46 measured values.We also show   2008).The ':' after a spectral type indicates uncertainty of ±1 whilst 'p' indicates peculiarity.The surface gravity flag  is given when appropriate, and is discussed in Section §5.3.1.The adopted spectral type is the kastredux method, only overwritten where there are gravity flags in the 'by eye' method.In addition, J1246+4027 has been typed as having a potential Li i detection (6708 Å), which can only be seen in the R300R spectra.
the difference between the published and measured values of the 20 overlapping objects.If there is more than one literature value, we take the weighted mean RV and standard error on the mean, to compare against the adopted RV from this work.We show literature measurements with respect to their resolutions and define these as: low,  < 2 500; mid, 2 500 ≤  ≤ 25 000; high,  > 25 000.The error used to define  are the quadrature summed errors from the literature and our adopted RV.
Our 46 RVs in the heliocentric reference frame are presented in Table 4.This reference frame has been experimented with, in that the heliocentric/barycentric correction via pypeit has been compared with a manual barycentric correction using barycorrpy (Kanodia & Wright 2018).Resultant RV differences from the manual barycentric correction to the pipeline barycentric correction differ by ≈0.1 km s −1 .The difference between heliocentric and barycentric correction is 0.5 km s −1 in the case of J1745−1640.
The median difference between our adopted RVs and the literature RVs was 7.8 km s −1 .This 7.8 km s −1 was then added in quadrature to our adopted RV error.We used this value to account for systematic uncertainties such as night-to-night instrumental drift and any FWHM undersampling.A S/N ratio of 20-30 also correlates with an RV uncertainty of ≈8 km s −1 , which was the typical S/N ratio seen around our alkali lines.Some lines, such as the potassium doublet, had lower S/N ratios and lower local resolutions due to a combination of wider features and lower flux values.All objects except J0940+2946 and J1221+0257 have an adopted and literature RV difference less than twice the sum of the respective errors in quadrature.J0940+2946 was 2.69 from the weighted mean literature value.Of the two literature values constructing this weighted mean, our value is <2 from the value from Kiman et al. (2019), which is notably larger than the value from Schmidt et al. (2010).J1221+0257 was 2.08 from the weighted mean literature value.Our RV value was closest to the value from Kiman et al. (2019), with less agreement shown with the value from Schmidt et al. (2010), which itself was most similar to the values from Burgasser et al. (2015) and Hsu et al. (2021).We note for both of these objects that the RV values from Schmidt et al. (2010) utilised considerably lower resolution spectra, hence a worse agreement being shown.Any objects in Table 4 which have known primary stars with literature RVs are discussed below: J1004+5022: G 196-3B is the binary companion to G 196-3A (Kirkpatrick et al. 2008).G 196-33A has a mean RV of −1.6 ± 0.4 km s −1 (Shkolnik et al. 2012;Schlieder et al. 2012b;Binks & Jeffries 2016;Gaia Collaboration et al. 2018a).This mean RV of the primary is 0.1 away from the RV of the secondary companion from this work.
J1441−0945: DENISJ144137.2−094558 is the binary companion to G 124-62 (Bouy et al. 2003;Seifahrt et al. 2005).G 124-62 has an RV of −41.65 ± 5.91 km s −1 (Gaia Collaboration et al.The RV values from the literature on the  axis with our adopted RV values, on the .We show a one-to-one relation, over which our 20 comparison RVs are plotted.Squares are from low-resolution literature measurements, whereas circles and diamonds are mid-and high-resolution literature measurements respectively.Orange points are like-for-like comparisons and blue points are for the three benchmark systems, i.e., comparisons between our measured secondary RV against the literature RV of the primary.2018a), which is within 1.4 of the companion (which had large uncertainties).

Moving groups
Our results for UVW Galactic kinematic components are presented in Table 5 with each object's moving group classification and associated probability from BANYAN Σ.When accounting for RV in BANYAN Σ, the resultant probability was often lower than the calculation without RV.This was due to the Bayesian probabilities being designed for a higher recovery rate (moving from 82 per cent to 90 per cent) when accounting for the RV (see the BANYAN Σ cautionary note 8 , Gagné et al. 2018).In addition, the RV uncertainties from this work are much higher than proper motion or parallax uncertainties from Gaia.

Galactic components
Thin disc objects were differentiated from thick disc and halo objects using the LSR corrected UVW Galactic velocities; the thick disc and halo objects were those with  total > 70 km s −1 and  total > 180 km s −1 respectively (Nissen & Schuster 2010). total is the total space velocity.We calculated upper and lower bounds for UVW Galactic velocities using the propagated parallax, proper motion, and RV errors; these UVW velocities with associated uncertainties are shown in Figure 10.The objects J1109−1606 ( total = 103 ± 5 km s −1 ) and J1539−0520 ( total = 69 ± 4 km s −1 ) are found using the above criteria to be most likely thick disc objects, and are highlighted in Figure 10.J1539−0520, is a borderline thick disc object, within 1 of the thick disc cut-off.Considering that a nearby object is most likely within the thin disc (Holmberg et al. 2009), J1539−0520 is a reasonable thick disc candidate, hence the inclusion here.It was also assigned a 64.6 per cent probability of being in the thick disc by Cooper et al. (2024), although it did not pass the conservative subdwarf candidate selection criteria in that work.Without metallicity information, an object being in the thick disc is not a direct inference on age.These objects are worth visiting with higher resolution spectroscopy to gain metallicity information, to confirm any potential subdwarf candidacy.This future work would also involved gathering NIR spectra, as in work by Zhang (2018); Zhang et al. (2018, and references therein).

Astrophysical parameters
We present the  eff and log  values from the model fitting (Section §4.2) in Table 6 along with  eff , assuming our adopted spectral type and equation ( 4) by Filippazzo et al. (2015) and teff_espucd values from Gaia DR3.In Figure 11, we plot the difference be-  2 used, 0 otherwise.Quoted RVs are already heliocentric corrected.A ' †' symbol next to an RV means the RV is that of the primary star in the common proper motion system a given object is part of.
tween our value and the expected value.In the cases of objects with both R2500I and R300R spectra available, we default to the higher resolution result.
Although the best-fitting surface gravity values can be indicative of youth, they are quite degenerate and without corresponding metallicity values, therefore they are not relied upon in our discussion below.The best fitting spectral sub-types and BT-Settl models are shown in a spectral sequence for R2500I VPH spectra in Fig- ures A1 and A2.
Figure 12 shows a set of colour-absolute magnitude diagrams (CAMD), 2MASS  −   , 2MASS   , AllWISE 1 − 2 and AllWISE  1 .Parallaxes from Gaia were used to generate the absolute magnitudes.Highlighted here are the objects with spectral features that are indicative of youth.These are compared to known young UCDs from Faherty et al. (2016) and Liu et al. (2016, 'VL-G' or 'Young'), as well as the full sample from the GUCDS.These young objects tend to be over-bright, although the effect varies across filters and is further complicated by intrinsic scatter plus variability.Burgasser et al. (2015), 3. Best et al. (2020), 4. Weinberger et al. (2016), 5. Blake et al. (2010).U is in the direction of the Galactic anti-centre.Derived using this work's adopted radial velocity in combination with Gaia DR3 kinematics unless otherwise indicated.We also show the predicted Galaxy component, taken from the UVW velocities and  total cuts in Nissen & Schuster (2010).†: J1213−0432 had an additional probability (26 per cent) of being a member of Argus, for a total non-field probability of 98 per cent.

Individual objects
We further discuss here objects we have indicated as being nontypical, with interesting features or results.We check for any age classifications, based on the moving group membership from BANYAN Σ and location on the CAMD in Figure 12.There are additional objects which exist in the same colour space as our highlighted objects in Figure 12 which are not discussed below.This is because there can be large implicit colour scatter due to unresolved binarity, metallicity and dust.Hence, only objects which are interesting either spectrally or kinematically are discussed.The following four objects were found to be members of the moving groups listed above, in Section §5.2.1.
J1109-1606 J1539-0520 Figure 10.Toomre diagram, as done by Bensby et al. (2005), using Gaia DR3 astrometry in combination with our calculated RVs.V is on the  axis, against the velocity dispersion ( √ U 2 + W 2 ) on the  axis.Black circles are UVW velocities calculated with the RVs from this work, with associated error-bars given.We show the respective thick disc and halo selection lines at  total > 70 km s −1 and  total > 180 km s −1 respectively.
The change of from Gaia DR2 to Gaia DR3 in isolation did not alter the confidence (99.2 per cent), whereas the inclusion of our adopted RV value dropped this to 98.9 per cent.Our adopted RV was 15.0 ± 8.3 km s −1 , which is within 1 of the 'optimal' RV from BANYAN Σ, 21.5 ± 1.5 km s −1 .From Figure 12, we see J0453−1751 (a) is photometrically similar to known young objects.Its  eff of 1850 ± 70 K is in good agreement with  eff and teff_espucd, although is cooler than the 2100 K from Gagné et al. (2015a).We can conclude that this object is an L3 within  Pictoris.
J0502+1442: 2MASS J05021345+1442367, an L0, we find as a member of the Hyades cluster with a 99 per cent probability.This improves the membership confidence by Gagné & Faherty (2018, 75 per cent) and concurs with the classifications by Gaia Collaboration et al. (100 per cent confidence, 2018b); Cantat-Gaudin et al. (100 per cent confidence, 2020, using the Melotte 25 name).Works by Oh & Evans (2020) and Spina et al. (2021) also placed this object in Melotte 25 with 96 per cent and 99 per cent confidences, respectively.It also agrees with the classification by Lodieu et al. (2019), which had a 'c parameter' of 5.88, well within their Hyades membership limit, c < 25.9. Figure 12, places J0502+1442 (b) also as photometrically similar to known young objects, being somewhat over-bright, although there is considerable overlap with standard M-L sequence.With a  eff of 2212 ± 126 K, J0502+1442 is an L0 object in the Hyades cluster.
J1058−1548: Another L3 object, SIPS J1058−1548, is classified with 93 per cent confidence as a member of Argus.Gagné et al. (2015b) had the same classification with a much lower probability (35 per cent).Gaia DR2 astrometry in isolation gave a confidence of 96.3 per cent, whilst Gaia DR3 reduced this to 94.8 per cent, the inclusion of our adopted RV value further dropped this to 93.1 per cent.Our adopted RV was −0.9 ± 11.1 km s −1 , which is within 1 of the 'optimal' RV from BANYAN Σ, 8.5 ± 1.4 km s −1 .Specifically, in Gaia DR2, this was  = 54.6 ± 0.5 mas,   cos  = −258.1±0.8 mas yr −1 and   = +31.1±0.7 mas yr −1 .In Gaia DR3,  = 55.1 ± 0.3 mas,   cos  = −258.6± 0.3 mas yr −1 and   = +30.8± 0.3 mas yr −1 .These values are in broad agreement with non-Gaia works, where  ranges from 49.2 mas-66.5 mas, ), but is not as convincingly over-bright as neighbouring known young objects, see (c) in Figure 12.Sanghi et al. (2023) conclude that for J1058−1548, "it is probable that the YMG assignment [Argus] is incorrect", because their spectrum well matched L-dwarf FLD-G standards, although the log  value of 4.27 dex was an outlier and more typical of a VL-G object (their figure ( 21)).The log  value in this work was 5.0 ± 0.3 dex, although this less robust than that from Sanghi et al. (2023), who also had a much lower  eff =1570 K, which itself is more akin to a cooler object, ≈L5.We would argue that this a probable L3 member of Argus but more high resolution spectra and modelling is required to ascertain youth.J1213−0432: 2MASS J12130336−0432437 (L4) we classify as a member of Carina-Near or Argus (98 per cent), which is an update on the 75 per cent classification of being in Carina-Near by Gagné & Faherty (2018).Just using Gaia DR2 astrometry gave a confidence of 68.5 per cent (with a 30.6 per cent likelihood of being in Argus), whilst Gaia DR3 increased this to 74.3 per cent (24.7 per cent for Argus), the inclusion of our adopted RV value (with large uncertainty) updated this to 72.0 per cent, with a 26.0 per cent likelihood of being in Argus.Our adopted RV was −25.3 ± 22.4 km s −1 , which is within 1.5 of the 'optimal' RV from BANYAN Σ, 2.4 ± 0.8 km s −1 .In Gaia DR2, it was  = 59.5 ± 1.0 mas,   cos  = −368.1 ± 2.2 mas yr −1 and   = −34.6±1.4 mas yr −1 .In Gaia DR3,  = 59.1±0.6 mas,   cos  = −367.9± 0.7 mas yr −1 and   = −34.0± 0.5 mas yr −1 .The work by Best et al. (2020) is also in good agreement,  = 56.3± 3.8 mas,   cos  = −380.9± 2.7 mas yr −1 and   = −33.4± 2.4 mas yr −1 .Figure 12 (d) shows this object as being under-bright compared with known young objects, with a  eff of 1783 ± 143 K. Being the age of Carina-Near could explain this relative under-brightness, as it should be tending towards field-like behaviour.This object can be classified then as an L4 member of Carina-Near.
There are two further field objects that we have highlighted as interesting due to their spectral features: J1246+4027: The L4 dwarf, 2MASSW J1246467+402715, observed at the two resolutions, is of interest due to the potential Li i detection at ≈6708 Å.As this feature is only in the wavelength regime of the R300R spectra, this is not definitive enough a detection to confirm lithium (see discussion by Martín et al. 2018, using the equation from Cayrel (1988)).Higher resolution (R ⪆2000) spectra would be required for confirmation (Gálvez-Ortiz et al. 2014).Assuming a true detection, employing the lithium test (Rebolo et al. 1992) alongside our fitted effective temperature of  eff = 1750 ± 91 K would identify this object as being substellar.This  eff is in good agreement with the expected temperature of  eff = 1717 ± 116 K and the Gaia DR3  eff of 1780 ± 162 K.This substellar argument is in line with discussion by Basri (1998), Martín et al. (1999a) and Kirkpatrick et al. (1999), because our  eff is in the range 2670 >  eff > 1400 K. Figure 12 suggests J1246+4027 (e) neighbours some known young objects.The best fitting model had a surface gravity of log  = 4.6 ± 0.3 dex, although we have no complementary metallicity information.BANYAN Σ finds no correlation with any known young moving groups.J1246+4027 could be classed as an L4 field object.
J1004+5022: G 196-3B is known to be a low gravity brown dwarf (Rebolo et al. 1998;Kirkpatrick et al. 2008;Allers & Liu 2013), to which we concur, with a spectral sub-type of L3.Our log  value is 4.5 ± 0.2 dex ( eff = 1740 ± 113 K), as would be expected from the already known young nature.This object sits extremely red and over-bright in Figure 12 (f), even more extremely than most known young objects.It is a companion to the well known G 196-3A M3 star, to which we compared our kinematics in Section §5.2, finding a 0.1 difference.There is much deeper discussion on this benchmark system by Zapatero Osorio et al. (2010), which measures an angular separation of  = 15.99 ± 0.06 ′′ .Combined with a Gaia DR3 parallax of  = 46.1952± 0.5452 mas (in agreement with the 49.0 ± 2.3 mas and 41.0 ± 4.1 mas measurements by Liu et al. 2016;Zapatero Osorio et al. 2014, respectively), this implies a projected separation of  = 739 ± 1 AU.This is slightly more than the projected physical separation range calculated by Zapatero Osorio et al. (2010), 285-640 AU.We found a probability of the secondary being a field object of 99.9 per cent, which is an increase on the 32 per cent probability of being a member of AB Doradus by Gagné et al. (2014).Liu et al. (2016) kinematically confirmed that G 196-3B is a young field object.This is also in agreement with the 50 per cent classification of the primary being a member of AB Doradus by Schlieder et al. (2012a), which was later downgraded to 0 per cent by Binks & Jeffries (2016); however, the primary was also classified as being a member of the controvertible Castor moving group (Barrado y Navascues 1998) with 75 per cent confidence (Klutsch et al. 2014).The Castor moving group was not included in BANYAN Σ, hence not being included in our analysis.We classify this object as an L3 object.

SUMMARY AND CONCLUSIONS
We have presented the low and mid resolution optical GTC/OSIRIS spectra of 53 objects observed between 2015 and 2016.Our data reduction was non-standard, using a pipeline package, PypeIt; this reduction was validated with an independent IRAF spectral extraction and calibration for one of the objects.We used kastredux to create 53 automated spectral types, six of which are for objects not yet spectrally typed, alongside the established technique of comparing against spectral standard template spectra.We found that our chosen spectral reduction package, PypeIt, introduced some non-optimal artefacts during reduction.One example is a spike appearing near the O 2 A band from the telluric correction procedure, which required interpolating over for visualisation purposes (it does not affect wavelength solutions).
In addition to using new data reduction software, we also used novel analysis software, rvfitter, that we developed to perform manual line centering and cross-correlation (against BT-Settl CIFIST models).The rvfitter code also used an uncertaintyweighted mean to create an adopted RV.This produced 46 RVs, 29 of which are new, which we have validated against standard IRAF and IDL software techniques.There were 17 RVs which were compared against literature values, showing good agreement with a median difference of 7.8 km s −1 , adopted as our systematic uncertainty.Our median RV uncertainty was 11.2 km s −1 , indicating that further high-resolution spectroscopy would be necessary to validate our RV values and conclusions.The cross-correlation also produced mean  eff and log  values for all 53 objects.In this work, we performed further analysis on our spectral types, RVs and  eff values by making comparisons to the literature where appropriate and ensuring all results were within two spectral sub-types, ΔRV < 2 and Δ eff < 2 (against  eff and Gaia DR3 teff_espucd).We then discussed any measurements which did not conform with these standards, such as J0940+2946, which had a ΔRV = 2.69.There were four objects that we classified through BANYAN Σ as being a member of a young moving group: SIPS J1058-1548 (J1058−1548), 2MASS J04532647-1751543 (J0453−1751), 2MASS J12130336-0432437 (J1213−0432), and 2MASS J05021345+1442367 (J0502+1442).There were two objects we placed as members of the thick disc: SIPS J1109-1606 (J1109−1606) and 2MASS J15394189-0520428 (J1539−0520).
Finally, by relating to gravity sensitive alkali lines and the aforementioned young moving group members, we discuss the interesting young candidates J1246+4027 and J1004+5022.2MASSW J1246467+402715 (J1246+4027) has a potential lithium indication and is otherwise an L4 field object.G 196-3B (J1004+5022) is confirmed as a young object, as was known from its primary companion.
In conclusion, this work was part of the GUCDS series of papers.A search of the GUCDS yields 145 known L dwarfs with measured RVs, excluding those from the SDSS.The 29 new L dwarf RVs presented in this work are therefore an ≈20 per cent increase to the number of 6-D complete L dwarfs.A number of interesting objects were identified or confirmed, either into young moving groups or young field objects.We used novel open-source techniques at all stages of our procedure, which we make available to the astronomical community.These techniques have been compared with established and accepted techniques in order to generate a baseline of trust.The observation campaign to complete the 30 pc sample is ongoing, with predominantly NIR spectrographs.This campaign will continue to produce work discussing, expanding and exploring this 30 pc sample.

A2 Comparison with standard routines
In the reduction we use two procedures based on IRAF and Python packages with a comparison target (J1745−1640, DENIS J174534.6−164053,Phan-Bao et al. 2008) as a sanity check.A full image and spectral reduction was carried out using standard tasks within the IRAF package on one of our target objects (J1745-1640) plus complimentary flux standard (Ross 640).This was done to assess both the quality of the data and to ascertain the necessary required reduction steps to maximise data quality.The results from this bespoke reduction method served as a reliable reference by which to measure the performance of a python pipeline (with support for the GTC/OSIRIS instrument recently added), which was later applied to all objects within our sample.

A2.1 Bespoke IRAF Reduction
Our IRAF reduction was applied to the science and calibration frames of J1745−1640 (L1-1.5)and Ross 640 (DZA6) as appropriate using the following tasks, beginning with basic image reduction: CCDPROC: Pre-scan bias level and bias structure removal; flatfielding; illumination correction; data section trimming.
Illumination and CCDPROC: Correction for spatial axis illumination gradients, made from the extensive sky lines of a well exposed object frame.
IDENTIFY , FITCOORDS and TRANSFORM: Correction for geometric image distortion (curvature) along the spatial axis sky background.
For the spectral reduction: APALL: Trace and extraction using both optimal and fixed-width aperture summing using image distortion corrected arc frames.
IDENTIFY and DISPCOR: Wavelength calibration to a linear wavelength dispersion using image distortion corrected arc frames.
STANDARD, SENSFUNC and CALIBRATE: Flux calibration from the flux standard Ross 640 taken on same night as the target.
In addition to the IRAF tasks mentioned above, two extra reduction software tools were utilised during the reduction process: DeFringFlat: An IDL routine aquired from the NASA IDL Astronomy library (Landsman 1993) was used to provide capabilities in de-fringing the flat field frames (DeFringFlat.pro;Rojo & Harrington 2006).
SKYCALC : ESO Sky Model Calculator provides additional telluric correction during flux calibration.A telluric sky model was queried using meteorological (e.g.moon phase, precipitable water vapour) and astrometric parameters (e.g.altitude, angular separation) appropriate for the object in question.
During the bias subtraction we discovered that the pre-scan region of the second CCD containing the spectrum displayed a gradient across it in ADU.A carefully chosen restricted section of the pre-scan region was used (∼3 pixels wide), which was found to be reliable for row-by-row bias level subtraction, before the 2D image bias structure was removed.
To correct for illumination gradients evident along the spatial axis of the 2D image introduced by the slit illumination function, we utilised the extensive sky lines of the well exposed object frames as a pseudo twilight sky flat (no sky flats were available).The IRAF Illumination task provided this functionality for correction, and we estimate that, after the correction was applied, the error introduced by the slit illumination gradient was reduced to a maximum of ∼1.5 per cent in the flat-field frames.
The latter, longer wavelength half of the flat-field frames showed evidence of fringing between wavelengths of approximately 8500 Å to 10, 000 Å, coincident with the area of the CCD containing the spectra of interest.We used the IDL routine DeFringFlat as mentioned above to attempt to remove as much of the fringing as possible and found the best fit using the Morlet 'wavemother' model, and near default parameters.We estimate from measuring the cleaned flat-fielded image that the amplitude of the fringing was reduced from an original level of approximately 7 per cent, to a maximum of about 1.7 per cent.
A combined arc frame was made from the three arcs available from the night of observation to cover the entire wavelength region of the spectrum.An initial wavelength solution was created and applied as part of the geometric image distortion correction, which resulted in a wavelength solution with an RMS error of 0.016 Å.A second wavelength calibration was subsequently made after additional reduction steps to ensure no systematic errors had been introduced, resulting in a more reasonable final RMS to the fitted wavelength solution of 0.025 Å.The final wavelength corrected spectrum had a linear dispersion 1.396 Å pixel −1 over the entire extracted range of 7339 Å -10, 155 Å.
Two separate flux calibrations were then made: one used a blackbody to represent the DZ white dwarf flux standard with an effective temperature 8070 K (Blouin et al. 2018) and with an -band magnitude of 13.66 mag (Bergeron et al. 2001); the second used the low resolution calibrated flux standard spectrum of Ross 640 contained in the IRAF database.In both cases, the sensitivity functions were created by interpolating over the affected telluric regions, and regions of intrinsic absorption features.Both of these sensitivity functions provided flux calibrations with almost identical results.A correction for atmospheric extinction and telluric features to the target was included during the flux calibration.An initial extinction correction was made from using a file containing tabulated extinction magnitudes as a function of wavelength applicable to the observatory site, that was provided on the GTC instrument website.However, an improved extinction correction was obtained from the much higher spectral resolution telluric sky model mentioned above (via the ESO Sky Model Calculator).The improvement is particularly evident over the wavelength regions containing the potassium K i 7665,7699 Å doublet and the H 2 O band at about 9500 Å.

A3 Radial velocity method validation
In keeping with our strategy outlined in Section A2 we again invoked an independent check, this time to validate our methods by helping to identify any problems with our RV measurements relating to the PypeIt reduced data set.The techniques used to measure RVs via the centres of atomic neutral alkali lines and through crosscorrelation of spectra were employed by Burgasser et al. (2015), and we adopt a similar twin measurement approach to derive our final RVs.We achieved this through the use of both IRAF and custom prepared routines within IDL to measure the RV via the Fourier cross-correlation and the line centre fitting methods.This analysis was conducted on the bespoke IRAF reduced data of our test object J1745−1640.We then used our validated RVs to classify any objects into young moving groups and stellar associations.

A3.1 Line centres
Two interactive methods were employed here: the first using routines in IDL to measure the 1D centroids of fitted Gaussian profiles to the atomic lines of J1745−1640, while the second used the IRAF task Splot to again measure the same lines but via fitting Voigt profiles.
In the first case, sub-sections of the spectrum surrounding the line features to be measured were extracted and interpolated onto a ten times finer wavelength grid, to facilitate the manual fitting of Gaussian profiles with a different number of terms via the Gaussfit.proroutine.Best fitting model profiles to spectral features were initially determined by eye, and determined by how closely the profile matched the feature with more emphasis being given around the line centre region.The reported RMS error and FWHM of fitted profiles were also taken into account for when the different Gaussian profiles produced similar results, such that the number of terms which fitted with the least error and narrowest FWHM were chosen.The measured wavelength shifts from labo-ratory rest-frame line centres (in standard air: Kramida et al. 2021) were then converted to Doppler RVs.
Secondly, and by using Splot, Voigt profiles were fitted to the same line features of appropriately pseudo-continuum subtracted sub-sections of the spectrum, and Doppler RVs were then found in the same manner as previously from the reported line centres.We obtained results for all eight line features from both measurement sets.However, it was apparent that four of the measurements gave the least error and particularly consistent results between both sets, these being Rb i-a, Rb i-b, Na i-a, Cs i-a with mean values for RV found from these four selected for each measurement set.The RV derived from the Gaussian fitted profiles (IDL) was found to be 35.1 km s −1 , and via Voigt profiles (Splot) 29.0 km s −1 (all test results are Heliocentric: barycentric correction calculated using baryvel.pro).Typically, we found that Gaussian profiles were more reliable to fit but Voigt profiles were best for lines which could be successfully fit.From the spread among the individually measured line shifts we place more confidence in the latter result, and assign uncertainties based on the 1- standard deviation of the respective RV measurements of 4.3 km s −1 and 3.8 km s −1 .
The RV as measured by our line centering method using the PypeIt reduced data for J1745−1640 is 36.2 ± 4.4 km s −1 (see Table 4) which is in broad agreement with those from this independent measurement test.The RV measured via line centre fitting as reported by Burgasser et al. (2015) is 28 ± 9 km s −1 .Thus, we have confidence in our RV results derived from our chosen method, which contribute to the final adopted values.

A3.2 Cross-correlation
To validate this second technique of measuring RVs as part of our adopted method, and to ascertain the best way forward in its application, we used the Fourier cross-correlation task Fxcor within IRAF to conduct tests.Our choice of RV rest-frame models were a BT-Settl model spectrum and custom-made synthetic atomic absorption spectra.Our object was again the bespoke IRAF reduced J1745−1640 spectrum.
The BT-Settl spectrum used was the best fitting model with the physical parameters of  eff = 2000 K, log  = 5 dex and Fe/H= 0 dex, corresponding to ≃L1 in spectral type.We smooth the spectrum using a Gaussian kernel to match the dispersion and resolution of the J1745−1640, and appropriate FITS header keywords added for the Fxcor task to recognise the template spectrum as rest-frame.
To help highlight any potential systematic wavelength shifts introduced by the use of the BT-Settl model, and therefore to help assess its suitability as an RV template, we measured the line centre locations of the most reliable Rb i-b and Cs i-a lines by fitting Voigt profiles in Splot.BT-Settl is known to generate models using a different line list to those selected in this work, where we used the NIST database.A maximum difference compared to laboratory rest-frame line centres of 0.13 Å was found, corresponding to 4.5 km s −1 .This shift is similar to the uncertainty found earlier from the fitted line profiles suggesting that the BT-Settl model is reliable for use as a template, however, we add this uncertainty in velocity units in quadrature to the subsequent Fxcor individual RV region measurements.
To facilitate the most accurate RV measurements we extracted sections of both object and template spectra into discrete spectral regions, each covering the main atomic absorption features as well as the FeH Wing-Ford band at ∼9900 Å, then each region was appropriately pseudo-continuum subtracted and normalised.
During the RV measurements, we interactively adjusted the sample test wavelength range around the features of interest to reduce noise in Fourier space domain.Next, the width of the crosscorrelation function (CCF) fit was changed to facilitate a best-fit (Gaussian fit to the CCF was used).The results of these changes to the CCF height, the goodness-of-fit 'R-value' and fit error were noted, until the best RV estimate was obtained.The shape of the CCF profile was also informative to this end, it tended to be broad, with no apparent double peaks seen.No Fourier filtering was applied as it was not found to be beneficial.
For this test, three regions gave consistent results covering both of the rubidium lines, the first caesium line (≈ 8500 Å) and the FeH Wing-Ford band.The average of these individual results gave an RV of 21.2 ± 5.2 km s −1 .
For our second test, we created a noise-free synthetic absorption spectrum of unity continuum with line widths and depths as measured by Voigt profiles of the neutral atomic lines in of J1745−1640, with no attempt to include the FeH band.The line centres were fixed to the laboratory rest-frame wavelength values.Results from all four regions were averaged which covered both of the rubidium lines, the sodium doublet and both caesium lines.Including the potassium doublet gave a similar result for that region but gave a very large increase in uncertainty, so was not included.We find a resulting RV of 24.6 ± 1.7 km s −1 .
Our final test was conducted to ascertain the intrinsic level of uncertainty in RV from the application of this method through the use of Fxcor on a representation of our spectral data.This involved making a cross-correlation between two noise-free synthetic absorption spectra: the same RV rest-frame template as used above in the second test, and with the object being a wavelength shifted version of the same synthetic spectrum, with the FITS header updated accordingly.The shift in wavelength was set at a value corresponding to the adopted RV presented in Burgasser et al. (2015), of 26.2 ± 2.3 km s −1 .We found the average combined RV of the four measured regions used to be 26.7 ± 1.2 km s −1 , indicating that 1.2 km s −1 is our base level uncertainty in using this method.This is, however, in addition to any uncertainty introduced from a real object spectrum (i.e.J1745−1640).
Both of these cross-correlation RV test results for J1745−1640 are in agreement with the equivalent value presented in Burgasser et al. (2015), within their respective uncertainties.The measured RV for J1745−1640 using the cross-correlation package we adopted and apply to our data set (see Section §4.3) has a value of 28.8 ± 4.7 km s −1 .Again, the results of this cross-correlation test validate our method and provide us with confidence in the separately derived RVs as well as in our final adopted values combined from both methods (see Section §4.3.3).

Figure 3 .Figure 4 .
Figure 3. Same as Figure 2 but for the second half of the R2500I VPHG sample.

Figure 5 .
Figure5.J1745−1640 RV calculation via different line profiles (orange: solid -Gaussian; dash-dot -Voigt) against the data (black squares) and fifth order spline fit (blue) in the regime around the eight listed line centres.The flux uncertainty is smaller than the height of each square.The shift from the laboratory line position (vertical dashed grey line) is shown as the vertical solid black line.The horizonal black line (solid or dash-dot, depending on the fitted line profile as above) is the continuum, as is subtracted from the data.A grey band is given, corresponding to the region of data the line profiles are fitted to.The shown region is between the inner edges of the continuum regions.

Figure 6 .
Figure 6.J1745−1640 RV calculation via the manually shifted BT-Settl model (orange) against the data (black squares) and fifth order spline fit (blue).The flux uncertainty is smaller than the height of each square.The laboratory line position (vertical dashed grey line) has been manually shifted by the RV given on the sub-plot title (vertical solid black line).Effective temperature, gravity and metallicity are also indicated on each features title.

PDFFigure 7 .
Figure 7. J1745−1640 RV values for each given line.In the top panel, orange squares are cross-correlated RVs, blue diamonds are line centre RVs; each spectral feature has been indicated on the  axis.In the bottom panel, the orange curve is the cross-correlated PDF; the blue curve is the line centre PDF; and the black curve is the adopted PDF.The dotted vertical lines are the mean RV values as associated with each PDF.

Figure 9 .
Figure 9. [Left Panel]: Histograms of the RVs calculated in this work (orange) and from the literature (blue) to show the relevant population densities.The dashed vertical lines indicate the means of the associated distributions.[Right Panel]:The RV values from the literature on the  axis with our adopted RV values, on the .We show a one-to-one relation, over which our 20 comparison RVs are plotted.Squares are from low-resolution literature measurements, whereas circles and diamonds are mid-and high-resolution literature measurements respectively.Orange points are like-for-like comparisons and blue points are for the three benchmark systems, i.e., comparisons between our measured secondary RV against the literature RV of the primary.

Figure 11 .
Figure11.The expected  eff (calculated via spectral type through a Filippazzo relation,Filippazzo et al. 2015) on the  axis and the best-fitting BT-Settl model mean  eff on the  axis.Blue crosses are for objects with a fit to the R300R spectra whilst black crosses are objects with a fit to the R2500I spectra.

Figure 12 .
Figure 12.CAMDs of 2MASS and AllWISE photometry, focused on the majority of this sample (an inset of the full sequence is shown in the upper right).The 2MASS  −   colour is on the  axis for the first column, with the AllWISE 1 − 2 colour on the  axis on the second column.Absolute 2MASS  magnitude is on the  axis for the first row whilst AllWISE  1 is the  axis of the second row.Underlying the plots as grey points is the full UCD sequence from the GUCDS.Known young objects fromFaherty et al. (2016) andLiu et al. (2016) are displayed as black diamonds.Each object is coloured by our adopted spectral type, with absolute magnitude error shown.Coloured diamonds are the young candidates discussed in Section §5.3.1.Key: a-J0453−1751, b-J0502+1442, c-J1058−1548, d-J1213−0432, e-J1246+4027, f-J1004+5022.
04. Cs i-b -Gaussian Profile with 2.0 Å; 11.1 ± 11.4 km s −1 ; RMSDIQR=0.25.RV Line Centre =36.2 ± 4.4 km s −1 .J1745−1640, Cross Correlation: K i-a -2200 K, log  = 5.0 dex ; 30.0 ± 5.0 km s −1 ; RMSDIQR=0.48.K i-b -2200 K, log  = 5.0 dex ; 20.0 ± 5.0 km s −1 ; RMSDIQR=0.20.Rb i-a -2200 K, log  = 5u x c a l i b ] e x t i n c t _ c o r r e c t = F a l s e f l u x r e a d . ./ S c i e n c e / < s p e c 1 d −s t a n d a r d .f i t s > s e n s f u n c .f i t s . ./ S c i e n c e / < s p e c 1 d −o b j e c t .f i t s > s e n s f u n c .f i t s f l u x end A6.4 Coadding [ c o a d d 1 d ] c o a d d f i l e = . ./ S c i e n c e / < s t a n d a r d .f i t s > c o a d d 1 d r e a d . ./ S c i e n c e / < s p e c 1 d −s t a n d a r d .f i t s > SPAT0240−SLIT0457−DET02 c o a d d 1 d end [ c o a d d 1 d ] c o a d d f i l e = . ./ S c i e n c e / < o b j e c t .f i t s > c o a d d 1 d r e a d . ./ S c i e n c e / < s p e c 1 d −o b j e c t .f i t s _ s h i f t _ b o u n d s = −10., 1 0 .

Table 1 :
The 53 targets observed at the GTC with OSIRIS and presented in this work.

Table 3 .
Our spectral types compared with the literature optical and near-infrared types for each object.

Table 4 .
RVs measured in this work and compared to the literature.

Table 5 .
The UVW velocities and BANYAN Σ classification (with associated probability) from this work.

Table 6 .
Effective temperatures and surface gravities from this work.