Investigation of the correlation between optical and $\gamma$-ray flux variation in the blazar Ton 599

The correlation between optical and $\gamma$-ray flux variations in blazars reveals a complex behaviour. In this study, we present our analysis of the connection between changes in optical and $\gamma$-ray emissions in the blazar Ton 599 over a span of approximately 15 years, from August 2008 to March 2023. Ton 599 reached its highest flux state across the entire electromagnetic spectrum during the second week of January 2023. To investigate the connection between changes in optical and $\gamma$-ray flux, we have designated five specific time periods, labeled as epochs A, B, C, D, and E. During periods B, C, D, and E, the source exhibited optical flares, while it was in its quiescent state during period A. The $\gamma$-ray counterparts to these optical flares are present during periods B, C, and E, however during period D, the $\gamma$-ray counterpart is either weak or absent. We conducted a broadband spectral energy distribution (SED) fitting by employing a one-zone leptonic emission model for these epochs. The SED analysis unveiled that the optical-UV emission primarily emanated from the accretion disk in quiescent period A, whereas synchrotron radiation from the jet dominated during periods B, C, D, and E. Diverse correlated patterns in the variations of optical and $\gamma$-ray emissions, like correlated optical and $\gamma$-ray flares, could be accounted for by changes in factors such as the magnetic field, bulk Lorentz factor, and electron density. On the other hand, an orphan optical flare could result from increased magnetic field and bulk Lorentz factor.


INTRODUCTION
Blazars, the peculiar category of active galactic nuclei (AGN) known for emitting non-thermal variable radiation covering the whole range of the electromagnetic spectrum, are among the most luminous objects in the Universe (luminosity ≈ 10 42 −10 48 erg s −1 ).These objects emit relativistic jets that are oriented in the same direction as the observer's line of sight (Antonucci 1993;Urry & Padovani 1995), resulting in strong Doppler boosting.They reside within elliptical galaxies (Urry et al. 2000) and their energy source is derived from the process of material accreting onto a supermassive black hole.(10 6 -10 10 M⊙) (Shakura & Sunyaev 1973;Shaw et al. 2012).
Blazars are categorized into two groups: flat spectrum radio quasars (FSRQs) and BL Lacs.This classification relies on the rest-frame equivalent width (EW) of the emission lines observed in their optical spectra.FSRQs exhibit strong emis-⋆ E-mail: bhoomikarjpt2@gmail.com sion lines (EW > 5 Å), while BL Lacs exhibit either weak or no emission lines (EW < 5 Å) (Stocke et al. 1991;Scarpa & Falomo 1997).The absence of broad emission lines in the optical spectra of BL Lac objects may be caused by Dopplerboosted non-thermal continuum, which could swamp out the spectral emission lines (Blandford & Rees 1978).Additionally, the presence of EW > 5 Å observed in BL Lac objects might be the outcome of an especially low state of the beamed continuum (Vermeulen et al. 1995).Consequently, the EW by itself is not a reliable indicator of the difference between the two categories of blazars.Ghisellini et al. (2011) and Sbarrato et al. (2012) proposed a categorization scheme for FSRQs and BL Lacs that hinges on the luminosity of the broad-line region (BLR), which is measured in Eddington units.In accordance with this scheme, the luminosity of the broad line region (LBLR) to the Eddington luminosity ratio for FSRQs is greater than 5 × 10 −4 while it is less than 5 × 10 −4 for BL Lacs.
The broadband spectral energy distribution (SED) of blazars has a two-hump structure, where one hump peaks at lower energies in the optical/IR/X-ray bands and the second hump peaks at higher energies in the MeV/GeV bands (Fossati et al. 1998;Mao et al. 2016).The low-energy emission is attributed to the synchrotron emission process of the relativistic electrons in the jets (Urry & Mushotzky 1982;Impey & Neugebauer 1988).Padovani & Giommi (1995) and Abdo et al. (2010c) proposed an alternative method of categorizing blazars by considering the position of the synchrotron peak frequency (ν syn peak ) in their SED.These classes include low-synchrotron peaked blazars (LSP, ν syn peak < 10 14 Hz), intermediate-synchrotron peaked blazars (ISP, 10 14 Hz < ν syn peak < 10 15 Hz), and high-synchrotron peaked blazars (HSP, 10 15 Hz < ν syn peak < 10 17 Hz).Most FSRQs are LSP blazars, while BL Lac sources exhibit all three blazar behaviors.Some BL Lacs are extremely high-synchrotron peaked blazars (EHBL or extreme blazars) with ν syn peak > 10 17 Hz (Foffano et al. 2019).
It is renowned that the synchrotron emission process results in low-energy emission.However, the origin of the higher energy peak remains partially understood.To explain this numerous models have been developed.Leptonic model often provides a satisfactory explanation for the observed broadband SED of the majority of blazars.In the leptonic scenario, the high-energy emission is produced by the inverse compton (IC) process (Abdo et al. 2010c).The seed photons that cause the IC process because of the scattering of the relativistic electron in the jets can originate from a variety of locations.These could originate from within the jets by synchrotron emission process (synchrotron self-Compton or SSC; Konigl 1981;Marscher & Gear 1985;Ghisellini & Maraschi 1989;Tramacere et al. 2009) or from outside the jets (external Compton or EC) mainly from the accretion disk (Boettcher et al. 1997;Dermer & Schlickeiser 2002), from the BLR (Ghisellini & Madau 1996;Dermer et al. 2009) and from the dusty torus (B lażejowski et al. 2000;Sikora et al. 2008).
In addition to the leptonic model, the high-energy peak in the broadband SED of blazars can also be accounted for by the hadronic model (Böttcher et al. 2013).Hadronic model for blazar postulates that protons are accelerated to ultra-relativistic energies and dominate the high-energy emission through proton-synchrotron radiation (Aharonian 2000;Mücke et al. 2003) or photo-pion production via protonproton or proton-photon interaction (Mannheim 1993), which results in the production of ultra-high-energy photons and neutrinos.The hybrid lepto-hadronic model has also been proposed in some cases to explain high-energy emission in blazar (Diltz & Böttcher 2016;Paliya et al. 2016;Cerruti et al. 2019).
In addition to the broadband SED, the study of the temporal correlation between the two major components of the broadband emission i.e. optical and γ-ray emission can give major insights to the high-energy emission in blazars.In the leptonic model of blazar emission, since the relativistic electrons in the jets are responsible for both the optical and γray emission (Böttcher 2007), a close correlation between the optical and γ-ray flux changes is anticipated.The launch of the Fermi Gamma-ray space telescope 1 (referred to as Fermi hereafter; Atwood et al. 2009) in 2008 gave us unusual opportunity to characterise the γ-ray for a large number of sources 1 https://fermi.gsfc.nasa.gov/ssc/data/access/lat/(Abdollahi et al. 2020).To investigate the connection between the optical and γ-ray bands, Fermi data was combined with numerous ground-based optical and IR observations.Several studies have been conducted on this topic, yielding diverse conclusions.For example, correlations between optical and γ-ray changes, with or without lag, have been found in blazars (Bonning et al. 2009;Chatterjee et al. 2012;Liao et al. 2014;Liodakis et al. 2018).These findings firmly support the leptonic model of blazar emission.However, many studies disfavour the one-zone leptonic model due to the absence of a correlation between changes in optical and γ-ray flux.In these studies, either γ-ray flares are present without their corresponding optical flares (orphan γ-ray flares) (Dutka et al. 2013;Cohen et al. 2014;MacDonald et al. 2015;Rajput et al. 2020) or optical flares are present without their γ-ray counterparts (orphan optical flares) (Chatterjee et al. 2013;Rajput et al. 2019, 2020 andreferences therein).
The existing observations suggest that the correlation between variations in optical flux and γ-ray flux is intricate.In this work, we conducted a thorough examination over the long-term duration to understand the diverse behaviour between changes in optical and γ-ray flux through broadband SED modeling in the blazar Ton 599.This blazar was recently observed at its brightest state (Prince 2023;Garrappa & Valverd 2023;Tripathi et al. 2023).The availability of multi-wavelength data over the long-term period renders this source appropriate candidate for conducting study on the optical-GeV connection.Ton 599, also known as 4C +29.45, is named as 4FGL J1159.5+2914 in the 4th catalog of Fermi Large Area Telescope (4FGL; Abdollahi et al. 2022).It is located at a redshift z = 0.725 (Hewett & Wild 2010), with RA = 179 • .882641and Dec = 29 • .245507.This source is Classified as an FSRQ source and it falls under the category of LSP source (Abdo et al. 2010a).It displays high optical variability and exhibits strong polarization (Fan et al. 2006).In the γ-ray band this source was first detected by the Energetic Gamma Ray Experiment Telescope (EGRET) (Thompson et al. 1995) and later by the Fermi (Abdo et al. 2010b).Moreover, this source was detected in the very high-energy band (VHE > 100 GeV) by VERITAS (Mukherjee & VERITAS Collaboration 2017).This source is observed with significant variability in the optical and γ-ray, as well as in other bands of the electromagnetic spectrum (Fan et al. 2006;Patel et al. 2018;Prince 2019;Rajput & Pandey 2021;Bhatta et al. 2023).In 2017, Ton 599 was detected in an unprecedented massive flaring state across the entire electromagnetic spectrum (Cheung et al. 2017;Carrasco et al. 2017;Acosta Pulido et al. 2017).The optical and GeV connection for this source during its 2017 outburst was investigated by Prince (2019), and optical and γ-ray emission regions are found to be co-spatial with the time lag of a few days.This source was modelled using a two-zone leptonic emission model during the same flare by Patel & Chitnis (2020), and it was found that the GeV emission was from EC process, where the seed photons originated from the dusty torus.The largest flare ever recorded by Fermi for this source in a period of around 15 years occurred in January 2023 when Ton 599 was once again found to be flaring, with a flux value of 3.60 ± 0.27 × 10 −6 ph cm −2 s −1 .Here, we present the findings of a multi-wavelength analysis performed on the blazar Ton 599 using data collected over a period of about 15 years, from August 2008 to March 2023.The main objective of this analysis is to determine any po-tential correlation between variations in the optical V-band flux and γ-ray flux and, as a result, putting constraints on the emission process in blazars.We provide the information about the data utilized in this work in section 2. Section 3 provides a detail of the analysis, and section 4 presents the findings and a discussion.The summary of the work is given in the final section.

MULTI-WAVELENGTH DATA REDUCTION
To carry out a multi-wavelength analysis of the source Ton 599, we used publicly accessible data in the γ-ray, X-ray, UV, and optical bands that span around 15 years, from August 2008 to March 2023.

γ-ray data
We analyzed about 15 years of data collected from the Large Area Telescope (LAT) aboard the Fermi.The Fermi-LAT is a pair-conversion telescope and has a very large effective area (> 8000 cm 2 ) (Atwood et al. 2009).It usually operates in the scanning mode and scans the entire sky once every ∼ 3 hours covering an extensive spectral energy range from 20 MeV to TeV energies.For this study, we performed binned likelihood analysis by using the Fermi-LAT data in the energy range of 100 MeV-300 GeV from 2008 August 08 to 2023 March 05, a period of approximately 15 years .To analyze the Fermi data we used the current version of Fermitools 2.2.0.We used PASS 8 data where the photonlike events are categorized as 'evclass=128, evtype=3'.All the events were chosen within a 10 • circular region of interest (ROI) centered on the source and a 90 • zenith angle cut was made to eliminate background γ-rays from the Earth's limb.After that, we generated a good time interval using the suggested criteria "(DATA QUAL > 0)&&(LAT CONFIG==1)".For this analysis, we used the latest instrument response function (IRF) "P8R3 SOURCE V3", Galactic diffuse emission model "gll iem v07" and the extragalactic isotropic diffuse emission model "iso P8R3 SOURCE V3 v1".Our source model includes all the point sources from the 4FGL catalog that fall within 10 • of the source, as well as the Galactic and the extragalactic isotropic diffuse emission components.The spectral shapes of all the sources are adopted from the 4FGL catalog.For the sources lying between 0 • to 5 • , the spectral parameters are left free except the scaling factor and for the source lying between 5 • to 10 • , all the spectral parameters are kept fixed to their 4FGL catalog value.To determine the significance of the detection, we used the maximum likelihood (ML) test statistic, defined as TS = 2∆ log(L), where L represents the likelihood function comparing models with and without a γ-ray point source at the source's location.During the likelihood analysis, the sources with TS < 25, which corresponds to a 5σ detection (Mattox et al. 1996), were removed.This modified model is then employed for the analysis of the light curve and spectra.To generate the one-day binned γ-ray light curve over the time period of interest, we regarded the source as detectable when the value of TS > 9, equivalent to a 3σ detection.

Swift-XRT data
The X-ray telescope (XRT), as detailed in Burrows et al. (2005), mounted on the Swift satellite (Gehrels et al. 2004), operates within the 0.3 to 10 KeV energy range, observed Ton 599 concurrently with Fermi.The X-ray data utilized in this work was taken from the archives at HEASARC 2 .To analyze the data, we adhered the standard procedure given by the instrument pipeline.We used the data obtained only in the photon counting (PC) mode for the light curve and spectrum analysis.We analyzed the data with the xrtpipeline task, utilizing the latest CALDB and response files provided in HEASOFT version 6.29.We employed the standard grade selection 0-12.To generate the energy spectra the calibrated and cleaned events files were added.The source spectra were extracted from a circular area with a radius of 60 ′′ , positioned at the center of the source, whereas, the background spectra were extracted from a region with a radius of 80 ′′ , located away from the source.The exposure maps were combined using the ximage task, and the ancillary response files (ARFs) were created using the xrtmkarf task.The corresponding response matrix files (RMFs) were obtained from CALDB.The RMFs and ARFs files were loaded into the grppha tool together with the source and background spectra in order to combine and group them.To perform the fitting of the final X-ray spectra within XSPEC (Arnaud 1996), we employed an absorbed simple power law model with a Galactic neutral hydrogen column density NH = 1.63×10 20 cm −2 (Dickey & Lockman 1990;Kalberla et al. 2005).

NuSTAR data
NuSTAR, the first focusing hard X-ray telescope, was launched in June 2012.It consists of two co-aligned X-ray detectors pair with their corresponding focal plane modules, called FPMA and FPMB and it operates over a wide energy range from 3 to 79 keV (Harrison et al. 2013).NuS-TAR recorded two observations of Ton 599 in May 2019 and June 2021.For this study, we used one observation (ObsID 60463037004) during the 2021 flare of Ton 599 with an exposure time of 17.5 ks.The data was reduced using the NuS-TAR Data Analysis Software package NuSTARDAS within the HEASOFT3 version 6.29.The calibrated and cleaned level 2 event files were produced using the nupipeline task and the CALDB version 20211202.We extracted the source and background spectra using a circular region of 30 ′′ with the nuproducts script.We grouped the spectra using grppha with at least 20 counts per bin and further processed it into XSPEC.

UV and optical data
We obtained the data in the UV and optical bands from Swift-UV-Optical telescope (Swift-UVOT).We analyzed this data using the online tool4 .To generate the light curves in the UV and optical bands, the magnitudes obtained using the online tool were converted into fluxes using the zero point  method described in Breeveld et al. (2011) without correcting for galactic reddening.However, to examine the spectral characteristics of the source, the magnitudes were corrected for galactic extinction (A λ ), whose values were obtained from NED 5 .In addition to using optical data obtained from Swift-UVOT, we also incorporated optical V-band data from the ground-based survey ASAS-SN 6 and from the Steward ob-0.0 0.2 0.4 55650 55675 55700 55725 55750 MJD(Days) 0.0 0.5 1.0 58007 58032 58057 58082 58107 MJD(Days) 0 5 servatory (SO)7 .The ASAS-SN ground-based survey operates with five stations positioned in both the northern and southern hemispheres.It has the capability to conduct daily observations of the entire visible sky, reaching a depth of g = 18.5 mag (Shappee et al. 2014;Kochanek et al. 2017).The Steward observatory, which is located in Arizona, employs

Multi-wavelength light curves
We have generated multi-wavelength light curves of the source Ton 599 for a period of about 15 years from August 08, 2008 to March 05, 2023 (MJD 54686−60008).The multiwavelength light curves, which include 1-day binned γ-ray, X-ray (PC mode), UV and optical light curves, are shown in Figure 1.From Figure 1 it is conspicuous that Ton 599 has numerous quiescent and flaring stages during this period.We have identified specific time periods within this 15-year duration.These include one quiescent period labeled as epoch A and four flaring periods (when an optical or γ-ray flare is present), marked as epochs B, C, D, and E. These quiescent and flaring epochs have a 100-day time span and flaring epochs are centered on the peak of either the γ-ray or optical light curves.The time periods of these epochs are given in Table 1.These epochs were chosen based on the following criteria: 1) The accessibility of multi-wavelength data 2) optical and γ-ray light-curve data during the quiescent phase should be below the average level 3) optical or γ-ray flare should gradually increase by more than two times from the average level.The values of average optical and γ-ray fluxes are 1.24±0.003× 10 −11 erg cm −2 s −1 and 0.42±0.003× 10 −6 ph cm −2 s −1 , respectively, over the period from August 2008 to March 2023.The average and peak values of optical and γ-ray flux for each epoch are given in Table 1.The details of each of the epochs are given below.

Epoch A (MJD 55650−55750)
The source was in a quiescent state in the γ-ray, X-ray, UV and optical bands during this epoch.The multi-wavelength light curve of this epoch is shown in the upper left panel of Figure 2.

Epoch B (MJD 58007−58107)
During this epoch, the source showed major variation in γray flux.The γ-ray flux increased by about six times from the average level during this epoch.Along with the large changes in the γ-ray flux, changes in the X-ray, UV and optical fluxes were also observed (refer to upper middle panel of Figure 2).The optical flux also increased by nearly 5 times from the average level.During this epoch, two short-term intense γ-ray flares were observed superimposed on the large γ-ray flare at around MJD 58057 and MJD 58077.Visual inspection of the Figure 2 upper middle panel) also gives the hint of corresponding short-term optical flare at around MJD 58057 but there is no hint of an optical flare at around MJD 58077 due to lack of data points.We deduce that during this time period, we noted the same patterns of flux variations in both the optical and GeV energy ranges.

Epoch C (MJD 59355−59455)
The source experienced a significant γ-ray flare during this epoch, this one being around six times greater than the γray flux on average.This substantial flaring event was subsequently accompanied by flares in the X-ray, UV, and optical bands.Notably, the optical flux experienced an increase of approximately five times compared to its average level.Thus, in this epoch too, we found that the variability pattern in the optical and γ-ray bands appeared to be correlated.The multi-wavelength light curve for this epoch is shown in the upper right panel of Figure 2.

Cross-correlation analysis
The initial investigation described in Section 3.1 involved visually examining the correlation between changes in optical and γ-ray emissions over different time periods.Furthermore, we explored this correlation in depth by employing advanced cross-correlation analysis methods.Analyzing the cross-correlation between variations in flux across various energy bands allows for the determination of whether emission in different bands are co-spatial or if they emanate from different regions in the jet.In this study, we conduct a cross-correlation analysis using the discrete cross-correlation function (DCF) as outlined in Edelson & Krolik (1988)) to measure the time lag between unevenly sampled optical (Vband) light curve and γ-ray light curve binned over one day.
The significance of the cross-correlation is evaluated through Monte Carlo simulations, with the assumption of a simple power-law power spectral density (PSD) model.This model is defined as P(ν) ∝ 1/ν β , where ν represents the temporal frequency, and β denotes the slope.A comprehensive explanation of this methodology is provided by Max-Moerbeck et al. (2014b).We have adopted the value of βγ = 1.5 (Max-Moerbeck et al. 2014a) and β optical = 2 (Goyal 2021).The DCF between optical and γ-ray light curves for epochs B, C, D, and E are shown in Figure 3.During these epochs, we detected no significant correlations between variations in optical and γ-ray flux.The results of the DCF analysis are provided in Table 2.

Optical spectral variations
Apart from changes in flux, blazars are also known for exhibiting variations in their spectral characteristics.In general, FSRQs exhibit redder-when-brighter (RWB) trend (Gu et al. 2006;Bonning et al. 2012;Zhang et al. 2015;Rajput et al. 2019;Safna et al. 2020;Negi et al. 2022).However bluer-when-brighter (BWB) pattern is also observed in FSRQ sources (Wu et al. 2011;Rajput et al. 2019;Safna et al. 2020;Rajput et al. 2020).The optical and UV spectral regions (V, B, U, W1, M2, and W2 bands) investigated in this paper include contributions from both non-thermal emission from the relativistic jet and thermal emission from the accretion disk.Hence, by examining spectral variations, one can identify the many factors causing the observed flux variations.To describe the spectral behaviour of the epochs studied in this work we looked for the variations in the B-V colour against the V-band magnitude.The colour-magnitude diagrams of these epochs are shown in Figure 4. Due to lack of data, we have not included the color-magnitude diagram for epoch A.
The trend of spectral behaviour during epochs B, C, D and E was identified by weighted linear least squares fit, which take into account errors in both colour and V-band magnitude.We determined the Spearman rank correlation coefficient (R) and probability of obtaining this correlation under the assumption of the null hypothesis (p-value).If R > 0.5 or R < −0.5 with p-value < 0.05, the source is considered to have displayed any significant colour trend over the epochs.We found no spectral variations during the epochs considered in this work.The results of the weighted linear least squares fitting are given in Table 3.

γ-ray spectral analysis
The inherent distribution of the emitting electrons in blazars may be reflected in the curvature of the γ-ray spectra, which also provides information on the potential acceleration and cooling processes in the jets.The γ-ray spectra of the blazars is either fit by power law (PL) or log parabola (LP) models.
The PL model is defined as: Here dN(E)/dE is the number of photons in cm −2 s −1 MeV −1 , N• is the normalization factor of the energy spectrum, E• is the scaling factor and ΓP is the photon index.
The LP model is defined as Here, α is the photon index at E and β is the curvature index that determines the curvature around the peak.
To evaluate the model that more accurately characterizes the γ-ray spectra (PL against LP), we applied the maximum likelihood estimation technique outlined in Mattox et al. (1996), utilizing the analysis tool gtlike.The gtlike tool maximizes the likelihood by employing different optimizers, which determine the optimal spectral parameter values for the best fit.We also calculated the TScurve value, as specified in Nolan et al. (2012), using the following equation: In the above equation, L is the likelihood function.In order to detect significant curvature in the γ-ray spectra, we employed the threshold value TScurve > 16 (at 4σ level; Mattox et al. 1996).The results of the γ-ray spectral analysis are given in Table 4 and the γ-ray spectra for the selected epochs are shown in Figure 5.The γ-ray spectra are well fit with the LP function for all the epochs, except the quiescent epoch A, where the value of TScurve is less than 16.

Broadband SED modelling
The source has shown varied flux variability behaviour over the last 15 years, during which we have detected instances of correlated optical and GeV variations as well as optical flare without γ-ray counterpart.To understand the diverse behaviour of the source during different epochs we constructed broadband spectral energy distributions and modelled them using one-zone leptonic emission model.To model the broadband SED, we used the publicly available code, JetSet (Tramacere et al. 2009(Tramacere et al. , 2011;;Tramacere 2020).JetSet (version 1.2.2) fits numerical models to the data in order to determine the optimal parameter values that best describe the observed data.
The photometric observations taken in both UV and optical wavelengths during each epoch were averaged filter-wise to generate a single photometric data point for each epoch.On the other hand, the average spectra for X-rays and γrays were created using all of the data that was available for the 100-day span in each of the epochs.Since the EC process primarily accounts for the high-energy γ-ray emission in the FSRQs (Paliya et al. 2015;Shah et al. 2017;Rajput et al. 2019Rajput et al. , 2020)), we used synchrotron, SSC, and EC processes to model our SEDs.In the one-zone leptonic emission model, the spherical blob of radius R, located at a distance of RH form the central black hole.It is postulated that this spherical region is injected with non-thermally accelerated electrons that follow a distribution characterized by a broken power law.The broken power law distribution of electrons is defined as follows: Where, γ is the electron Lorentz factor, N is the electron density in units of cm −3 and K is the normalization constant.p and p1 are the low and high energy spectral slopes, and γmin, γmax and γ break are the Lorentz factors corresponding to the low-energy cut off, high-energy cut off and turn over energy, respectively.The emission region (spherical blob), which is moving down the jet with the bulk Lorentz factor Γ, is permeated with the magnetic field B. The relativistic electrons interact with the magnetic field and produce synchrotron radiation.The SSC emission is produced by the same synchrotron photons generated within the jet (Tramacere et al. 2009).The seed photons for the EC process are external to the jet, which are produced (1) from the accretion disk (Dermer & Schlickeiser 2002;Ghisellini et al. 2009), (2) from broad line region (Dermer et al. 2009), and (3) from dust torus (Sikora et al. 2008).The external photons that directly come from the accretion disk and participate in the EC process have luminosity L disk and temperature T disk .The reprocessed UV photons from the BLR, which is located at a distance of RBLR from the black hole, and IR photons from the dust torus, which is located at a distance of RDT with a temperature of TDT, are involved in the EC process.The fraction of L disk reprocessed by the BLR and dust torus are given as τBLR and τDT.
We have determined the size of the emission region R, which is defined as R ≤ cTminδ/(1+z), where Tmin is the minimum variability timescale and δ denotes the Doppler factor.We have calculated the minimum variability timescale for the γ-ray light curve using the methodology that models flares with an exponential profile, as outlined in Burbidge et al. (1974): Here, ∆t = | t2-t1| and F2 and F1 are the flux values at a time t2 and t1, respectively.The minimum variability timescale is calculated to be ∼ 7.23 hours.We have adopted the Doppler factor value for the source as δ = 18.2, which has been determined using the one-zone leptonic emission model by Paliya et al. (2017).Using the values of Tmin and δ, we have approximated the size of emission region to be 8.24 × 10 15 cm ∼ 8 × 10 15 cm.Furthermore, assuming that this blob covers the whole cross-section of the jets, we have estimated the distance of this blob from the central black hole using the relation RH = R/ϕ, where ϕ represents the opening angle.The intrinsic opening angle for the Ton 599 is estimated as 0.58 • (∼ 0.01 radian) based on observations conducted with the 15.4 GHz Very Long Baseline Array (VLBA) from the 2 cm VLBA MOJAVE program (Pushkarev et al. 2009).Utilizing the values of R and ϕ, we estimate that the location of the emission region is located at roughly 8 × 10 17 cm away from the central black hole.
We have considered a value for L disk equal to 4.5 × 10 45 erg/sec, which was determined through the application of a one-zone leptonic emission model in studies by Ghisellini et al. (2010); Paliya et al. (2017).The inner radius of the spherical shell like BLR (RBLR in ) is considered to be 0.01 pc ∼ 3.1 × 10 16 cm (Donea & Protheroe 2003).The outer radius of the BLR (RBLR out ) is estimated to be around 2.1 × 10 17 cm using the relationship RBLR = 10 17 L 1/2 disk,45 cm, where L disk,45 represents the luminosity of the accretion disk in units of 10 45 erg/sec.This estimation is based on the assumption that the BLR radius scales in relation to the square root of the luminosity of the accretion disk (Ghisellini & Tavecchio 2008).Moreover, it is posited that the spherical shell like torus is situated at a distance RDT = 2.5 × 10 18 L 1/2 disk,45 cm (Ghisellini & Tavecchio 2008), which is approximately equivalent to 5 × 10 18 cm.We have assumed that the fraction of the L disk is reprocessed by the BLR (τBLR) is 0.1 and by the dust torus (τDT)is 0.5 (Ghisellini et al. 2010).Furthermore, the disk temperature T disk and dust torus temperature TDT are assumed to be at their typical values of 1.0 × 10 5 K (Patel & Chitnis 2020) and 1.0 × 10 3 K (Donea & Protheroe 2003), respectively, in our study.And the ratio of cold protons to relativistic electrons is taken as 0.1 in our SED analysis (Ghisellini 2012). .PL and LP fits to the γ-ray spectra for epochs A, B, C, D and E.Among these, the PL fit provides a better description of the γ-ray spectra in epoch A, while the LP fit better explains the γ-ray spectra in epochs B, C, D, and E. Refer to Table 4 for detailed results.Table 4. Spectral analysis results for the γ-ray observations.The γ-ray flux values are expressed in units of 10 −6 ph cm −2 s −1 .The PL and LP fits to the γ-ray spectra during these time periods are shown in Figure 5. Table 5.Values of the frozen parameters for the one-zone leptonic model fits to the observed SEDs during the five epochs.
name of the parameters values 8.0 × 10 17 cm cold proton to relativistic electron ratio 0.1 In one-zone leptonic model, the characteristics of the observed SED are primarily determined by nineteen parameters.Five parameters, namely, p, p1, N, B, and Γ, are chosen to be free, while the remaining fourteen parameters are chosen to be frozen.The purpose of parameter freezing is to compare the changes in certain parameters, such as the magnetic field, bulk Lorentz factor, and particle density, etc., between quiescent and flaring epochs.Details of frozen parameters are given in Table 5.In order to investigate the distinct correlation behavior between optical and γ-ray flux variations during different flaring epochs, we initially fitted the quiescent epoch to acquire the values for the parameters γmin, γmax and γ break .The fitting process using the same frozen parameters and the previously acquired values of γmin, γmax, and γ break from the fitting during the quiescent epoch, was repeated for the other epochs B, C, D, and E with the selection of free parameters.Our SED analysis shows that the leptonic scenario provides a good fit to the observed broadband SED across all the epochs.The resultant best fit parameters for the selected epochs are given in Table 6 and the fitted SEDs are shown in Figures 6.

Connection between optical and γ-ray flux variations
Over a span of around 15 years of study, Ton 599 has shown varied correlated behavior between optical V-band flux and γ-ray flux variations.We found two cases, where, (1) both the optical flux and the γ-ray flux showed flaring behavior, (2) an increase in optical flux without a corresponding increase in γ-ray flux.We chose five epochs A, B, C, D, and E based on the distinct behaviour between optical V-band and γ-ray flux variations.During epoch A, the source was in its quiescent state in each waveband.During epochs B, C and E, there was a structural correlation between the variations in optical and γ-ray flux.However, during epoch D, the source exhibited an optical flare without a corresponding γ-ray response.Our investigation involving linear regression analysis on the logarithmic values of γ-ray flux and optical flux suggests a similar trend.The linear least square fits are shown in Figure 7, depicting the logarithm of γ-ray flux plotted against the logarithm of optical flux.The outcomes of the linear regression analysis are given in Table 7.Moreover, we employed the sophisticated Discrete Cross Correlation Function to crosscorrelate the optical V-band and γ-ray light curves, aiming to establish whether there is a correlation between variations in the optical flux and γ-ray flux.We found that there are weakly significant correlations between the optical and γ-ray flux variations.Given the observed structural similarities between these two energy bands, this low significance may result from the relatively low amount of optical data used.
A time lag of a few days within the period of October-December 2017, falling within the timeframe of epoch B considered here, was previously identified by Prince (2019).In the cross-correlation study of bright blazars, Liodakis et al. (2018Liodakis et al. ( , 2019) ) found that the time lag between optical and γ-ray bands are generally small.Moreover, de Jaeger et al. (2023) found that the time lag distribution between optical and γ-ray bands is consistent with zero lags.These results indicate that the optical and γ-ray flux emission regions are cospatial and favour leptonic emission model of blazars, where the same relativistic electron population is accountable for both the optical and γ-ray emission.Given the scenario that both optical and γ-ray emissions are proportionate to the number of emitting electrons, it is plausible that variations in the electron population of the jet lead to the correlation between these two types of flux variations, as proposed by D' Ammando et al. (2011).

γ-ray spectra
The γ-ray spectral shape of blazars can yield valuable information about the inherent characteristics of the jets such as the distribution of emitting particles.In particular, these γray spectra can help us determine if there is any absorption of γ-rays within the emission region and can also reveal whether γ-ray emissions are attributed to a single component of the jet or multiple components.The high-energy γ-ray spectra of FSRQ sources are better explained by the broken power law or log parabola model (Abdo et al. 2010a;Paliya et al. 2015;Rajput et al. 2019;Sahakyan 2020;Rajput et al. 2020;Rajput & Pandey 2021).The curvature observed in the γ-ray spectra of FSRQs may be attributed to external factors, such as the absorption of γ-rays through photon-photon pair production within the BLR (Poutanen & Stern 2010), and the impact of the Klein-Nishina effect on the inverse Compton scattering of BLR photons by relativistic jet electrons with a curved distribution (Cerruti et al. 2013).However, both the effects account for the IC scattering of the BLR photons, whereas in the large number of FSRQs the γ-ray emission site lies outside the BLR and γ-ray emission is not produced by the IC scattering of the BLR photons (Costamante et al. 2018).Our SED analysis also indicates that the γ-ray emission is produced by the IC scattering of the IR photons from the dusty torus.
In our work, the LP model fits the γ-ray spectra in the flaring epochs (and orphan optical flare epoch) adequately, whereas the PL model fits the γ-ray spectra in the quiescent epoch well.The curvature observed in the γ-ray spectra is probably caused by the inherent characteristics of the electron distribution that produces emission.These electrons follow either a cutoff energy distribution or a log parabola energy distribution.However, during quiescent epoch the low quality of photons data available for analysis makes it challenging for us to identify this curvature.In the LP model fits to γ-ray spectra, the parameters α and β convey crucial information about the properties of the γ-ray spectra.Changes in the α and β values are reflected in the shape of γ-ray spectra.We looked for the variations in the α and β parameters against γ-ray flux during epochs B, C, D and E, which are shown in the Figure 8.We conducted an analysis to find a correlation between the α and β parameters and the γ-ray flux, employing the Spearman rank test.We considered the correlation is significant if the Spearman rank correlation coefficient > 0.5 or < −0.5, along with a p-value < 0.05.We observed that the parameter α showed no significant pattern, and the curvature parameter β exhibited a mild negative correlation with γ-ray flux throughout all epochs.The results of Spearman rank test are given in Table 8.Decreasing trend of β with increasing flux was also noticed in the study of 3LAC (The Third Catalog of Active Galactic Nuclei Detected by the Fermi) sources by Ackermann et al. (2015).According to Coogan et al. (2016), the variations in the location of the γray emission area during various activity levels of the sources may be responsible for these distinct spectral behaviors of the γ-ray spectra.Such varied spectral behavior could also be caused by a change in the IC peak frequency, where the external photon field is responsible for producing the γ-ray emission in FSRQs (Rajput et al. 2020).

Broadband spectral energy distribution
In our study, we found instances of both optical and γ-ray flares as well as optical flare without its corresponding γ-ray counterpart.Liodakis et al. (2019) investigated for a correlation between optical and GeV flux changes in 178 blazars and found that in nearly 50% of the blazars, optical flares occurred without their γ-ray counterparts.A Few blazars with complex behavior between changes in optical and γray flux were also studied by Chatterjee et al. (2013)     could provide a possible explanation for optical flux variations without γ-ray counterparts or vice versa (Mücke & Protheroe 2001).Furthermore, the leptonic scenario of blazars favours a close correlation between changes in the optical and GeV flux.However, results from Rajput et al. (2019Rajput et al. ( , 2020) ) indicate that different correlation behaviour between optical and γ-ray flux changes are well explained by one-zone leptonic model scenario.
In this work, we performed broadband SED modeling using a one-zone leptonic emission model to model the five epochs considered here.The leptonic model provides a good fit during all the epochs with the seed photons for the EC process from the torus.The source was in a quiescent condition during epoch A, therefore we first fitted the broadband SED for epoch A and then looked for changes in parameters in other epochs relative to the quiescent epoch.The accretion disk component is clearly evident in the broadband SED during the quiescent phase, but synchrotron radiation dominates the optical-UV emission during epochs B, C, D, and E. The source showed both optical and γ-ray flares during epochs B, C and E. During epochs B and E, the bulk Lorentz factor was roughly 1.5 times higher than the quiescent epoch A. Furthermore, there was an increase in electron density and a decrease in the magnetic field during period B. In contrast, during period E, there was a decrease in electron density and an increase in the magnetic field.The enhanced optical and γ-ray fluxes in these two time intervals may be attributed to the increase bulk Lorentz factor and changes in electron density and magnetic field.In comparison to the quiescent period, the bulk Lorentz factor changed slightly during epoch C, and there was also decrease in the magnetic field.However, the electron density increased during epoch C, which might be accountable for the optical and γ-ray flares during this epoch.
During Epoch D, when the source showed an orphan optical flare, the bulk Lorentz factor increased almost 1.5 times, while the electron density decreased when compared to Epoch A. Additionally, the magnetic field was at its peak during this epoch.Consequently, the orphan optical flare can be the result of the increased magnetic field and bulk Lorentz factor during epoch D. According to Marscher & Gear (1985) the change in magnetic field can affect the synchrotron emission in the jet, however the γ-ray emission will be unperturbed by the change in the magnetic field.The presence of increased optical polarization measurements and X-ray flare can support this argument.A significant X-ray flare is observed during this time (see the Figure 1), which is dominated by SSC process (the same photons emitted by the synchrotron process are involved in the EC process), however, the optical polarization measurements are not available at this time.From SED analysis, we draw the conclusion that the leptonic model explains the optical-GeV variations during the various epochs in Ton 599.

SUMMARY
We investigated the correlation between the changes in optical V-band and GeV flux in the FSRQ Ton 599 for the time period MJD 54686-60008.Our study includes 1) the identification of the epochs exhibiting optical or γ-ray flare, 2) cross-correlation analysis between optical flux and γ-ray flux variations, 3) analysis of optical spectral variations, 4) analysis of γ-ray spectra 5) broadband SED modeling with a one-zone leptonic emission model for the chosen epochs.We summarize the results below (i) In the period of around 15 years, we selected five epochs A, B, C, D, and E. The source was in its quiescent state during epoch A. Our visual inspections indicated that, in epochs B, C, and E, the source showed variations in both optical and γ-ray emissions, showing structural similarities.However, during epoch, D the source showed an optical flare, without a clear γ-ray counterpart.The source was in its highest flaring state during epoch E.
(ii) Through cross-correlation analysis, we identified a weak linear correlation between optical and γ-ray light curves.
(iii) We found no optical spectral variations during the epochs.
(iv) During the quiescent period, the γ-ray spectra of the source was well described by the PL model, while during the flaring periods, the γ-ray spectra was best described by the LP model.
(v) During all the epochs, from SED modelling, we found the γ-ray emission to be produced by EC process, where the source of seed photons is the torus.Our SED analysis showed that during the quiescent period, the optical-UV emission was caused by the accretion disk component, but during the flaring periods, the synchrotron emission was the cause for the optical-UV emission.The varied behaviour observed in our work between the changes in optical and GeV flux, could be explained by changes in the magnetic field, the bulk Lorentz factor, and the electron density.

Figure 1 .
Figure 1.Multi-wavelength light curves of the source Ton 599 for the period from MJD 54686 to MJD 60008.The top panel refers to the one day binned γ-ray (0.1-300 GeV) light curve, the second panel refers to the X-ray (0.3-10 keV) light curve, the third panel refers to the UV light curve, the fourth and fifth panels refer to the optical light curves.The γ-ray fluxes are measured in the units of 10 −6 ph cm −2 s −1 and the UV and optical fluxes are measured in the units of 10 −11 erg cm −2 s −1 .The two vertical black dashed lines represent the quiescent epoch A. The vertical blue dashed lines refer to the flaring epochs B, C and E, when both optical and γ-ray flares are present.The two vertical green dashed lines represent the flaring epoch D, when optical flare is present without its corresponding γ-ray counterpart.The peaks of optical and γ-ray flares are denoted by the vertical cyan dashed lines.The horizontal black dashed lines in the γ-ray and optical V-band light curves represents the average flux from August 2008 to March 2023.

Figure 2 .
Figure 2. Multi-wavelength light curves for the selected epochs A, B, C, D, and E. The light curves for epochs A, B, and C are displayed in the upper left, middle, and right panels, respectively.The bottom left and right panels represent the light curves for epochs D and E, respectively.The detail of the panels are same as in Figure 1.
During this epoch, the γ-ray flux variations are moderate but the changes in the X-ray, UV and optical flux are significant (see the lower left panel of the Figure2).In comparison to the average flux level, the γ-ray flux increased by less than two times, whereas the optical flux increased by about four times compared to the average flux level.Since there was no clear large amplitude variation in the γ-ray light curve at this epoch, we can infer that in this epoch we observed an orphan optical flare.3.1.5Epoch E (MJD 59905−60005)The source underwent the major γ-ray flare during this epoch, which is the largest γ-ray flare over the whole time period examined in this work.Such large flux variations were also noticed in the X-ray, UV and optical bands.The multiwavelength light curve for this epoch is shown in lower right panel of Figure 2.During this epoch the γ-ray flux increased about nine times from the average flux level and the optical flux increased nearly eight times compared to the average flux.We conclude that in this epoch the optical and γ-ray flux changes are structurally correlated.

Figure 3 .
Figure 3.The outcomes of the discrete cross-correlation analysis conducted on the optical V-band and γ-ray light curves during epochs B, C, D, and E. The DCF is indicated by the black dashed line.The significance levels are represented by the yellow dashed-dotted line (68.27%,1σ), green dotted line (95.45%,2σ), and red dashed line (99.73%,3σ).A negative lag indicates the γ-ray flux leads the optical flux variations.

Figure 4 .
Figure 4. Colour-magnitude diagram for Epochs B, C, D and E. The cyan line depicts the weighted linear least squares fit to the data.
Figure5.PL and LP fits to the γ-ray spectra for epochs A, B, C, D and E.Among these, the PL fit provides a better description of the γ-ray spectra in epoch A, while the LP fit better explains the γ-ray spectra in epochs B, C, D, and E. Refer to Table4for detailed results.

Figure 6 .
Figure 6.Observed broad band SED along with the one zone leptonic emission model fit during epochs A, B, C, D and E (top panel) and the residuals (bottom panel).Here the solid line is the best fit to the SED, orange dashed represents the emission from synchrotron process, green dashed line represents the emission from SSC process, red dashed line represents the EC emission from dusty torus and pink dashed line represents the EC emission from BLR. Purple and brown dashed lines represent the thermal emissions from the dusty torus and accretion disk, respectively.The fit to the quiescent state is applied for determining the values of the parameters γ min , γmax, and γ break in the subsequent epochs modelled in this study.MNRAS 000, 1-16 (2015)

Figure 7
Figure 7. γ-ray v/s Optical V-band fluxes during epochs B, C, D, and E.

Figure 8 .
Figure 8. Correlation between the parameters α and β and the γ-ray flux.The blue line represents the least-square fit applied to the data.

Table 1 .
Details of the epochs studied in this work for optical-GeV correlation and broadband SED modelling.In this table the optical flux values are in units of 10 −11 erg cm −2 s −1 , and the γ-ray flux values are in units of 10 −6 ph cm −2 s −1 .

Table 2 .
The outcomes of the DCF analysis between optical flux and γ-ray flux variations during epochs B, C, D, and E.Here τ represents the time lag, with negative values indicating that the γray is leading the optical flux variations.The DCF values provide estimates of the cross-correlation coefficient, the p-value represents the likelihood of observing the correlation under the assumption of the null hypothesis.

Table 3 .
Results of the weighted linear least squares fit and Spearman rank correlation test applied to the colour-magnitude diagram during epochs B, C, D, and E. None of the correlations are significant for a chosen significance level of 0.05.

Table 6 .
Results of the broadband SED fitting during the epochs A, B, C, D, and E.

Table 7 .
The results obtained by applying linear regression on the logarithm of γ-ray flux and optical flux during epochs B, C, D, and E.

Table 8 .
Results of the Spearman rank test for the correlation of α and β parameters with γ-ray flux.The plots of α and β parameters v/s γ-ray flux are shown in Figure 8.