A Study of the Radio Spectrum of Mrk 421

We present the results of a spectral analysis using simultaneous multifrequency (22, 43, 86, and 129 GHz) very long baseline interferometry (VLBI) observations of the Korean VLBI Network on BL Lac object, Markarian 421. The data we used were obtained from 2013 January to 2018 June. The light curves showed several flux enhancements with global decreases. To separate the variable and quiescent components in the multifrequency light curves for milliarcsecond-scale emission regions, we assumed that the quiescent radiation comes from the emission regions radiating constant optically thin synchrotron emissions (i.e., a minimum flux density with an optically thin spectral index). The quiescent spectrum determined from the multifrequency light curves was subtracted from the total CLEAN flux density, yielding a variable component in the flux that produces the time-dependent spectrum. We found that the observed spectra were flat at 22–43 GHz, and relatively steep at 43–86 GHz, whereas the quiescent-corrected spectra are sometimes quite different from the observed spectra (e.g., sometimes inverted at 22–43 GHz). The quiescent-corrected spectral indices were much more variable than the observed spectral indices. This spectral investigation implies that the quiescent-spectrum correction can significantly affect the multifrequency spectral index of variable compact radio sources such as blazars. Therefore, the synchrotron self-absorption B-field strength (B SSA) can be significantly affected because B SSA is proportional to the fifth power of turnover frequency.


INTRODUCTION
Active Galactic Nuclei (AGN) are the central regions of galaxies that exhibit unusually bright luminosity across multiwavelengths, from radio to gamma rays.It is thought that the mechanism responsible for this luminosity is a supermassive black hole (10 6 − 10 10 M ⊙ ) at the center of the AGN, which accretes enormous amounts of matter from its surrounding environment.A relativistic jet is formed by the accreted matter (Ulrich et al. 1997).
A blazar is a subclass of AGN that exhibits extremely energetic and highly variable observational behavior.The source is among the most luminous and powerful objects, and features a relativistic jet that points directly toward us with a very small viewing angle.This alignment causes a blazar to appear as a point-like source in the sky that exhibits rapid and dramatic variability in its brightness (Urry & Padovani 1995).Blazars encompass both BL Lac objects and flat-spectrum radio quasars (FSRQs).While FSRQs and BL Lac objects share similar characteristics, such as being highly variable and bright, their spectra differ.The spectrum of a BL Lac object has weak or no broad emission, unlike that of an FSRQ.
Markarian 421 (Mrk 421), is one of the most extensively studied high-synchrotron-peaked blazars and one of the closest sources, located at a redshift z=0.031.Its high luminosity enables us to investigate its emission (Katarzyński et al. 2003;Charlot et al. 2006;Horan et al. 2009;Acciari et al. 2011;Aleksić et al. 2015;Arbet-Engels et al. 2021).Its source Mrk 421 was the first extragalactic source detected in the TeV energy bands, as reported by Punch et al. (1992) using the Whipple 10-m Cherenkov telescope.Mrk 421 is classified as a BL Lac object.Acciari et al. (2011) reported no correlation between TeV γ-ray and optical/radio light curves obtained from observations performed from December 2005 to May 2008 with VERITAS/Whipple gamma-ray data and Metsahovi and UMRAO radio data.On the other hand, when focusing on the very high energy (VHE) E > 100 GeV, they found a strong correlation between TeV gamma-ray and X-rays and suggested that this can be explained by the one-zone synchrotron self-Compton model (Zhu et al. 2016).In a multiwavelength study conducted on Mrk 421 over 5.5-yr period, from December 2012 to April 2018, Arbet-Engels et al. (2021) performed several cross-correlation analyses.They found that the TeV and X-ray light curves were very well correlated, with a lag of less than 0.6 days.The GeV and 15 GHz radio light curves also showed correlation, with a lag of over 30-100 days.However, no correlation was observed between the TeV and the radio light curves.They suggested that the variability at high energy may be due to the leptonic mechanism, given the short time scale of variability.
Although extensive multiwavelength studies on Mrk 421 have been extensively conducted to date, research on its radio spectrum has been limited.In particular, simultaneous multifrequency very long baseline interferometry (VLBI) imaging observations with high time resolution have not been performed yet.In this study, we present various analyses of the radio spectrum of Mrk 421 to explain the correlation between the spectrum and the flux density variations.
In this paper, we report results from simultaneous multifrequency VLBI observations of Mrk 421 carried out using the Korea VLBI Network (KVN) at 22, 43, and 86, from January 2013 to June 2018.The paper is organized as follows.
In Section 2, we present descriptions of the observations for both the KVN and all archival data used in this paper.Section 3 and 4 describe the results and various analyses of flux density variability at multiple frequencies, respectively.Finally, Section 5 includes a discussion on the results and analysis.

KVN observations
We conducted simultaneous multifrequency VLBI imaging monitoring observations of Mrk 421 using the Korean VLBI Network (KVN) at 22, 43, 86, and 129 GHz as part of the iMOGABA(Interferometric Monitoring of Gamma-ray Bright active galactic Nuclei) project (Lee et al. 2016(Lee et al. , 2017b;;Algaba et al. 2018;Lee et al. 2020;Kang et al. 2021).The KVN is a VLBI network system consisting of three 21-meter radio telescopes located in Seoul (KVN Yonsei), Ulsan (KVN Ulsan), and Jeju (KVN Tamna), Korea.The monitoring observations were primarily carried out every month.However, we excluded the 129 GHz data from the analysis because there were no valid measurements, higher than 3-σ of errors.The observation period was from January 2013 to June 2018, except for maintenance periods from June to August.Measurements were made in 41 epochs over a period of 5.5 years.At 22, 43, and 86 GHz, 40 epochs, 33 epochs, and 21 epochs were measured, respectively.The observing frequencies were 21.700-21.764, 43.400-43.464, 86.800-86.864, and 129.300-129.364GHz, with a total bandwidth of 256 MHz.The angular resolutions of the KVN were 6, 3, and 1.5 mas at 22, 43, and 86 GHz, respectively.We used a left circular polarization observation mode at a recording rate of 1 Gbps.For more details of the KVN observations, see Lee et al. (2016Lee et al. ( , 2017a)).

Data calibration
The observed data were processed using the DiFX software correlator in Daejeon, Korea (Deller et al. 2007;Lee et al. 2015).The DiFX correlator generated a cross-correlation function spectrum with a resolution of 0.125 MHz and an accumulation period of one second.After the correlation, the data was calibrated with the AIPS (Astronomy Image Processing System) software package, provided by the National Radio Astronomy Observatory (NRAO), following standard procedures, including phase and amplitude calibration, fringe fitting, and bandpass calibration.We used the KVN pipeline developed by Hodgson et al. (2016) to perform the data calibration.Sensitivity at high frequencies was improved by transferring phase solutions from lower to higher frequencies (i.e., FPT) as described by Rioja & Dodson (2011) and Algaba et al. (2015).Amplitude calibration was conducted using the system temperatures measured during the observations.The system temperatures were adjusted for atmospheric opacity based on sky-tipping curve measurements conducted every hour at each radio telescope.Renormalization of the fringe amplitudes was conducted to correct the amplitude distortion due to quantization, and losses from quantization and re-quantization were corrected  (April 23, 2014).The left, middle, and right panels indicate source structures at 22, 43, and 86 GHz, respectively.The source shows a compact core-dominated structure at all frequencies.The contour level starts from 3 times the noise and increases by a factor of √ 2. The beam size is indicated by the filled ellipse in the bottom left corner.(Lee et al. 2015).After the amplitude calibration and the re-quantization loss calibration, the uncertainty of the amplitude calibration was within 5 % at 22 and 43 GHz and 10-30 % at 86 and 129 GHz.

Imaging and modelfitting
After the phase and amplitude calibration, we used DIFMAP to make a CLEAN contour map.To make a CLEAN map with the phase self-calibration, for the first time we performed the startmod command.This uses a point source model of 1 Jy to conduct the initial phase self-calibration.Then, the visibility data at 22, 43, and 86 GHz were averaged at 30-second intervals.This averaging time aligned with the typical coherence timescales of KVN observations.We flagged bad amplitude and phase data points and outliers in the vplot.Then, the repetitive CLEANing and phase self-calibration within the central emitting regions were performed until the noise root-mean-square (rms) level was no longer significantly lowered.The standard CLEANing and self-calibration within the central emitting regions were conducted until no significant flux density was added compared to the image rms level.
Since Mrk 421 is a compact source observed on mas-scales in the KVN observations, finding the best model is straightforward.Figure 1 displays the CLEAN contour maps of Mrk 421 at 22, 43, and 86 GHz for a representative epoch (well observed at all frequencies) showing compact core-dominated components.After making the CLEAN contour map, we computed the image quality factor ξ r .This is the ratio of the image rms noise and its mathematical expectation, ξ r = S r /S r,exp , where S r is the maximum absolute flux density in the residual image, and S r,exp is the expectation of S r .For a more comprehensive understanding of this image quality evaluation scheme, please refer to Lobanov et al. (2006) and Lee et al. (2016Lee et al. ( , 2017a)).The ξ r values derived in this paper range from 0.65 -0.97, suggesting that the images sufficiently represent the structure identified in the visibility data.
To establish the model of the source, we fit the model in the uv -data.In this process, we used a single circular two-dimensional Gaussian model to obtain the model parameters, the information of the model, such as flux, size, and degree of the major axis.The model fitting was repeated until reduced χ 2 approached one, and model fitting stops when the reduced χ 2 does not decrease anymore.
Table 1 lists the fitted parameters, including the observing frequencies, restoring beam size B maj,min , position angle B PA , total CLEAN flux density S KVN , peak flux density S p in units of Jy per beam, rms in the image σ, the dynamic range (ratio of peak to rms) of the image D, and image quality factor ξ r .

Multifrequency light curves
We present the results of the multifrequency radio observations of Mrk 421. Figure 2 displays the multifrequency light curves of Mrk 421 at 22, 43, and 86 GHz observed between January 16, 2013 (MJD 56308) and June 1, 2018 (MJD 58270).During the 5.5-year observation period, the flux density gradually decreases and several small flux enhancements are shown (which we call flares).The light curves show a similar trend at all frequencies.Although we can't show the 15 GHz light curve of the OVRO in this study, the 15 GHz light curve exhibits flux enhancements peaking at approximately MJD 56400, 56900, and 57200 (see Figure 1 in Arbet-Engels et al. 2021).In the KVN light curves, only the flare at MJD 56900 is seen due to the relatively sparse data points.Because of that, the flare at MJD 56900 looks like a large flare.But the flare is very sharp indeed.Q -43 GHz band; W -86 GHz band; D -129 GHz band; 4-6 -restoring beam: 4 -major axis; 5 -minor axis; 6 -position angle of the major axis; 7 -total CLEAN KVN flux density; 8 -peak flux density; 9 -off-source RMS in the image; 10 -dynamic range of the image; 11 -quality of the residual noise in the image (i.e., the ratio of the image root-mean-square noise to its mathematical expectation).

Radio spectrum
The simultaneous multifrequency monitoring KVN observations enable us to investigate source spectra accurately without the time uncertainty.The multifrequency spectral analyses of variable blazars using arcsecond-scale observations by Kim et al. (2022) and Jeong et al. (2023) show that variable spectral properties can be revealed by subtracting off the quiescent spectrum in the multifrequency light curves.We calculated the quiescent spectrum by fitting the four local minima obtained on MJD 56716, 57107, 57400, and 57840 -58093 shown in Figure 3.The fourth local minimum we used is obtained on MJD 58093 at 22 GHz, MJD 58079 at 43 GHz, and MJD 57840 at 86 GHz.The quiescent spectrum is obtained as 0.263±0.045Jy, 0.271±0.021Jy, and 0.156±0.027Jy at 22, 43, and 86 GHz, respectively, by extrapolating each of the three frequencies to a common time (MJD 58269), as shown in Figure 4.
The mas-scale spectra of Mrk 421 are shown in Figure 5. Out of the 41 epochs, data existed for at least two frequencies for 33 epochs.The observed spectrum is shown in black, and the quiescent-corrected spectrum (the observed spectrum by subtracting the quiescent spectrum) is shown in red.The spectra generally appear to be flat between 22 and 43 GHz, and seem to be relatively steep between 43 and 86 GHz.The quiescent-corrected spectra are sometimes consistent with the observed spectra and sometimes quite different.Unlike the observed spectra, the quiescent-corrected spectra generally show inverted spectra between 22 and 43 GHz.Although we were not able to precisely calculate the turnover frequency, ν c , of Mrk 421 in this analysis, we clearly find that the quiescent-corrected spectrum becomes steeper between 22 and 43 GHz and flatter between 43 and 86 GHz than the observed spectrum.This implies that the quiescent emission correction for the source spectrum in the mas-scale may be important to study the intrinsic spectral properties of the variable emission regions.3 × 10 1 4 × 10 1 6 × 10 1
In Figure 6, we see that the 22 and 43 GHz quiescent-corrected spectral indices are consistently lower (steeper) than the observed spectral indices, although the errors are large.We noticed that the values of α 22/43 were steeper when the flares were around the peaks on MJD 56901 and MJD 57289, whereas those were flatter at local minima (MJD 56716 and MJD 57400).
The 43 and 86 GHz quiescent-corrected spectral indices showed more variability with spectral indices that were sometimes consistent, sometimes lower, and sometimes higher than the observed spectral indices, although the errors were large.To further investigate the effect of correcting for the quiescent spectrum, we plotted the spectral indices versus flux density.In Figures 7 -8 indices were less correlated with the 86 GHz flux densities than before.

DISCUSSION
Figures 7 and 8 show that correcting for the quiescent spectrum can significantly affect the measurements of the spectral indices.This is of particular importance for calculations that require measurements of the synchrotron selfabsorption frequency, for example, the magnetic field strength, B SSA ∝ ν 5 c .This means that even small errors in measuring ν c can lead to significant errors when measuring the magnetic field strength.This suggests that magnetic field strength estimations based on multifrequency data should be carefully obtained with a proper determination of the optically thin quiescent spectrum of compact radio sources even on mas-scales.
To test the previous work, which did not consider correcting the quiescent spectrum, we re-analyzed the multifrequency spectral properties of OJ 287 (Lee et al. 2020).The minimum flux densities of OJ 287 were 2. 58,2.27,1.83,and 1.36 at 22,43,86,and 129 GHz,respectively.The spectral index between 86 and 129 GHz was −0.73.We used this optically thin spectral index to correct the observed spectrum.We subtracted the minimum flux density obtained from the optically thin spectral index at each frequency for 8 epochs analyzed in Lee et al. (2020).Figure 9 shows the spectra of OJ 287.The black and grey symbols indicate the observed and the quiescent-corrected spectra, respectively.After subtracting the optically thin spectrum, the spectrum was fitted with a curved power-law function (Lee et al. 2016(Lee et al. , 2017b(Lee et al. , 2020)).Then, we obtained the turnover frequency and the peak flux density.The turnover frequencies ob- tained from the quiescent-corrected spectra were in the range of 28-80 GHz, and the peak flux densities were 0.7-3.3Jy.The turnover frequencies were larger by a factor of approximately 2.5 larger than the ones obtained from the observed spectra.Because of the B SSA ∝ ν 5 c , the synchrotron self-absorption magnetic field strength might be underestimated in Lee et al. (2020).This could also potentially change the comparison result between the equipartition magnetic field strength and synchrotron self-absorption magnetic field strength.

SUMMARY
We have presented the results of simultaneous multifrequency observations at 22, 43, and 86 GHz with the KVN.The light curves showed a global decaying and several small radio flares.To separate the variable component and the quiescent component in the light curves, we subtracted the quiescent component, assuming the minimum flux density.We found that the observed spectra were flat at 22-43 GHz, and relatively steep at 43-86 GHz, whereas the quiescent-corrected spectra are sometimes quite different from the observed spectra (e.g., sometimes inverted at 22-43 GHz ).This means that the turnover frequency, ν c is shifted to a higher frequency.Based on the assumption of the synchrotron self-absorption magnetic field strength estimation, B SSA ∝ ν 5 c .Therefore, to properly estimate the magnetic field strength, subtracting the quiescent spectrum in the KVN light curves is important.

Figure 1 .
Figure 1.CLEAN maps of Mrk 421 obtained in Epoch 14 (April 23, 2014).The left, middle, and right panels indicate source structures at 22, 43, and 86 GHz, respectively.The source shows a compact core-dominated structure at all frequencies.The contour level starts from 3 times the noise and increases by a factor of √ 2. The beam size is indicated by the filled ellipse in the bottom left corner.

Figure 4 .
Figure 4.The minimum flux density at each frequency, by fitting the four local minima obtained on MJD 56716, 57107, 57400, and 57840 -58093, and extrapolating each of the three frequencies to a common time (MJD 58269).

Figure 6 .
Figure 6.Spectral index α of Mrk 421.The top panel is α 22/43 and the bottom panel is α 43/86 .The red symbols indicate the quiescent-corrected spectral index.A horizontal line is a guideline of zero.

Figure 7 .
Figure 7. Spectral index as a function of the flux density.The left and right panels indicate 22 and 43 GHz, respectively.a) and b) indicate α vs. flux density.c) and d) indicate the αcor vs. flux density.ρ and p values indicate coefficient and probability, respectively.

Figure 8 .
Figure 8. Spectral index as a function of the flux density.The left and right panels indicate 43 and 86 GHz, respectively.a) and b) indicate α vs. flux density.c) and d) indicate the αcor vs. flux density.ρ and p values indicate coefficient and probability, respectively.