Analysis of BMR Tilt from AutoTAB Catalog: Hinting toward the Thin Flux Tube Model?

One of the intriguing mechanisms of the Sun is the formation of bipolar magnetic regions (BMRs) in the solar convection zone (CZ), which are observed as regions of concentrated magnetic fields of opposite polarity on the photosphere. These BMRs are tilted with respect to the equatorial line, which statistically increases with latitude. The thin flux tube model, employing the rise of magnetically buoyant flux loops and their twist by Coriolis force, is a popular paradigm for explaining the formation of tilted BMRs. In this study, we assess the validity of the thin flux tube model by analyzing the tracked BMR data obtained through the Automatic Tracking Algorithm for BMRs. Our observations reveal that the tracked BMRs exhibit the expected collective behaviors. We find that the polarity separation of BMRs increases over their normalized lifetime, supporting the assumption of a rising flux tube from the CZ. Moreover, we observe an increasing trend of the tilt with the flux of the BMR, suggesting that rising flux tubes associated with lower flux regions are primarily influenced by drag force and Coriolis force, while in higher flux regions, magnetic buoyancy dominates. Furthermore, we observe Joy’s law dependence for emerging BMRs from their first detection, indicating that at least a portion of the tilt observed in BMRs can be attributed to the Coriolis force. Notably, lower flux regions exhibit a higher amount of fluctuations associated with their tilt measurement compared to stronger flux regions, suggesting that lower flux regions are more susceptible to turbulent convection.


INTRODUCTION
The Bipolar Magnetic Regions (BMRs), commonly observed in line-of-sight (LOS) magnetograms, emerge on the solar surface in an East-West orientation with a finite tilt having the leading polarity closer to the equator (Hale et al. 1919).Statistically, BMR tilt increases with the latitude of emergence, which is commonly referred to as Joy's law.This behavior was validated by many authors, including, Wang & Sheeley (1991); Howard (1991); Sivaraman et al. (1999) using white-light observations prior to the availability of magnetograms.
anubsreedevi.rs.phy20@itbhu.ac.in karak.phy@iitbhu.ac.in The tilt in the BMRs plays an essential role in the reversal of the existing poloidal field through the dispersal and cancellation of fluxes on the solar surface (Babcock 1961;Leighton 1964).This process is popularly known as the Babcock-Leighton mechanism and is an essential component of solar dynamo (e.g., Cameron & Schüssler 2023;Karak 2023).
Our present understanding attributes the formation of BMRs to the magnetically buoyant, large-scale toroidal flux tubes generated by the solar dynamo process at the base of the convection zone (BCZ; Parker 1955).Numerical simulation studies, assuming the thin flux tube approximation, have focused on understanding the dynamics of Ω-shaped flux tubes and provide constraints to the magnetic field strength in the CZ (Caligari et al. 1995;Fan 2009).The works of D 'Silva & Choudhuri (1993) and Fan et al. (1994) demonstrate that the tilt in the BMR is due to the action of Coriolis force on the diverging flows at the apex of rising flux tubes.Consequently, stronger BMRs, with higher magnetic field strength at the BCZ, are expected to ascend rapidly through the CZ and experience Coriolis force for less time, leading to a reduced amount of tilt in them.Their works also predicted that with the increase of the flux in the tube, the BMR tilt increases due to the effect of the drag.Observational studies, such as as by Tian et al. (2003) using magnetic field observations from Huairou Solar Observatory Station, validate these findings, noting an initial increase in tilt angle with flux followed by a decrease for higher flux BMRs.Similar variations in tilt angle with magnetic field strength are observed in simulations by Weber et al. (2011) and in observational studies by Jha et al. (2020).Interestingly, Stenflo & Kosovichev (2012) did not find any systematic dependence of tilt on flux from LOS magnetic field observations of the Michelson Doppler Imager (MDI).
The theory of Coriolis force as the reason behind tilted BMR is very promising and has been studied extensively.If this theory is true, then BMRs should emerge on the photosphere with a definite tilt.Analyzing 715 BMRs from MDI magnetograms during 1996-2008, Kosovichev & Stenflo (2008) found Joy's law behavior of BMR tilts at their mid-emergence phase, concluding, 'The observations indeed show the predicted latitudinal dependence (Joy's law) and indicate that the tilt is formed below the surface.'Later, Schunker et al. (2020) noted that BMRs emerge with zero tilt and develop tilt in accordance with Joy's law at a later stage in their lifespan, challenging the thin flux tube model's prediction.They suggested that the observed Joy's law behavior is because of the inherent north-south separation of BMRs as they reach the surface.
Joy's law is statistical, and thus, the latitudinal dependency of tilt is evident after averaging over large data samples.Persistent, significant scatter around Joy's law is consistently observed in various observational studies across different data sets (e.g.,Howard 1996;Wang & Sheeley 1991;Dasi-Espuig et al. 2010;McClintock & Norton 2013).Fan et al. (1994) and Longcope & Fisher (1996) hint towards the role of turbulent convection on rising flux tubes as the possible cause of the scatter in tilt angle.This observed behaviour was further supported by the simulation of Weber et al. (2013).Recently, a thorough analysis on the inconsistencies in Joy's law has been explored by Will et al. (2023).
Currently, observational support for the thin flux tube model is very limited and not adequate enough to establish their existence.Lack of dependence of BMR's tilt on flux has been used to strongly rule out the thin flux model based on studies presented in Kosovichev & Stenflo (2008); Schunker et al. (2019).While the latter study was based on a limited data set, the study of Stenflo & Kosovichev (2012) did not include the tracked information of BMRs and hence they counted every detection of BMR as a new one, giving a higher weightage to longliving BMRs (bigger ones).On the basis of all these results, we are still not in a position to completely rule out the thin flux tube model.
This manuscript aims to understand the origin mechanism of BMRs using the tracked BMR information from the Automated Tracking Algorithm for BMRs (Au-toTAB; Jha et al. 2021;Sreedevi et al. 2023) catalog, which contains the tracked information of 12,610 unique BMRs for the period 1996 -2023.We study the general behavior of BMR properties over their lifetime and explore the validity of the thin flux tube model.Before we present our study, in Section 2, we list what kind of observational signatures we expect based on this theory.Then in Section 4 we present our result followed by conclusions in Section 5.

OBSERVATIONAL EXPECTATIONS OF THE THIN FLUX TUBE MODEL
The simplest explanation for the formation of the BMR is given by the thin flux tube model (Parker 1955).The model assumes a thin, untwisted flux tube with a diameter negligible compared to the length scale of perturbation in the tube, anchored in the deep CZ.Upon becoming magnetically buoyant, tubes rise through the CZ to the photosphere.If the flux tubes remain anchored in the CZ, they rise in the form of an Ω-shaped loop structure.As the flux tubes rise, the draining plasma from the apex of the tube is subjected to Coriolis force, which leads to the observed tilt in the BMRs (D 'Silva & Choudhuri 1993).
Therefore, if Ω-shaped loops do exist and they are responsible for the formation of the BMRs, one can verify this idea by tracking the separation of BMR's opposite polarities during the initial phases of their evolutions.We expect the polarities to separate as BMRs mature, indicating the rise of the Ω loop when observed in LOS magnetograms.
Second, if the Coriolis force acts on the rising flux tubes, then the leading polarity of the BMR is expected to be closer to the equator, and BMR's tilt will have a latitudinal dependency (e.g., see Section 2 of D 'Silva & Choudhuri (1993) for rough calculations).In such a scenario, we anticipate the BMRs will emerge with a definite tilt, and they are expected to follow Joy's law due to the pronounced effect of the Coriolis force in the CZ.
Finally, BMRs with higher flux are expected to have higher tilt as predicted by the empirical relation given in Fan et al. (1994) for the tilt angle (γ), (1) Here, B 0 is the initial magnetic field in the flux tube (forming BMR) at the BCZ, Φ is the magnetic flux inside the rising flux tubes, and λ is the emerging latitude of the flux tube.We note that, in the above model (Equation 1), B 0 and the Φ are made independent of each other, while in the Sun this may not be true.Nevertheless, based on the above relation, we expect that the tilt (γ) decreases with the increase of magnetic field strength (B 0 ).Although the magnetic flux that is measured in the BMR (Φ m ) is not exactly the same as Φ, they are expected to be related.However, for the magnetic field it is not obvious; field that is observed in the BMR on the solar surface is quite different from the initial field B 0 of the flux tube.Despite this, we shall also check the dependence of tilt with the measured magnetic field in addition to the dependency on the flux in the BMR.To the best of our knowledge, no definitive observational evidence supporting the thin flux tube model was conclusively confirmed.In Section 4, through the catalog of AutoTAB, we evaluate if the BMR properties align with the mentioned findings from the thin flux tube model.

DATA AND METHOD
We start with a brief description of the AutoTAB catalog (Sreedevi & Jha 2023;Sreedevi et al. 2023) analyzed in this study.The catalog encompasses tracked information of 12,173 BMRs during 1996 -2023, generated using AutoTAB.These tracked BMRs are observed in the LOS magnetograms from MDI (1996MDI ( -2011;;Scherrer et al. 1995, full cadence) and Helioseismic and Magnetic Imager (HMI: 2010-present;Scherrer et al. 2012, 96 minutes cadence) onboard the Solar and Heliospheric Observatory (SOHO) and Solar Dynamic Observatory (SDO), during the period of September 1996 to December 2023 which includes complete Cycles 23, 24, and early 25.The working of AutoTAB is summarized below while the detail is published in Sreedevi et al. (2023) 1 .
To automatically detect the BMRs from the LOS magnetogram, AutoTAB uses a similar method prescribed by Stenflo & Kosovichev (2012), which was also adopted by Jha et al. (2020) with slight modifications.All the detected regions satisfy the flux balance condition, similar to Stenflo & Kosovichev (2012) and Jha et al. (2020) to ensure that the amount of positive and negative flux are balanced in each BMR.This flux balance condition only checks for the amount of flux (positive and negative), rather than the distribution of flux in BMRs.Hence, the AutoTAB catalogue (this work) also includes the multi-polar regions which are very frequent and thus important in contributing to the polar field in Sun (Yeates 2020).The detected BMRs are saved as binary files, and AutoTAB uses these files to track the detected regions.These regions undergo pre-processing steps before the tracking for improved tracking efficiency, followed by the technique of feature association in successive binary files to track the BMR in future instances; see Sections 2.2 and 2.3 of Sreedevi et al. (2023) for details of AutoTAB.It has to be noted that the detection and tracking algorithms operate independently.Therefore, AutoTAB gives users freedom to choose a different methods of detection for BMR or any other features, which can be efficiently tracked by AutoTAB.
AutoTAB tracks BMRs through their evolution on the nearside of the Sun, and the tracked BMRs have a wide range of lifetime2 , which includes those tracked from emergence to decay, those tracked in their evolutionary stages only, and a small category of BMRs that live for less than 8 Hrs (mostly the ephemeral regions).These three groups were respectively classified as "Lifetime (LT)", "Diskpassage (DP)", and "Short-lived (SL)" by Sreedevi et al. (2023).The BMRs falling in the latter class are excluded from this study and will be explored in detail in the forthcoming manuscript.Lifetime BMRs are mostly composed of moderately small BMRs within the flux range of 10 20 -10 22 Mx; some of them may not produce enough contrast in white-light images to be visible as sunspots.Meanwhile, the DP class constitutes BMRs, which have been only tracked during their appearance on the near side of the Sun.Therefore, this class includes (i) BMRs that are first detected near the east limb (< 45 • E) and disperse on the nearside; (ii) BMRs that emerge on the nearside of the Sun but cross the west limb (> 45 • W); and (iii) BMRs that first detected near the east limb (< 45 • E) and cross the west limb (> 45 • W).Hence, this class mainly comprises larger and stronger BMRs, exhibiting consistent evolution of magnetic properties throughout their lifetime.The majority of the tracked BMRs fall into this category (Sreedevi et al. 2023).
AutoTAB catalog includes the tracked information of all BMRs and their physical parameters, such as maximum magnetic field (B max ), average magnetic field (B mean ), total unsigned flux (Φ m ), tilt (γ), as well as positional information (latitude λ and longitude ϕ), at each instance during their lifetime.The tilt of a BMR is calculated similarly to Stenflo & Kosovichev (2012) and Jha et al. (2020), which varies between ±90 • .The convention followed here is that the BMRs, which strictly follow Joy's law, have a positive tilt in the Northern Hemisphere, whereas it is negative in the Southern Hemisphere.Furthermore, the tilt of the BMRs emerging in the southern hemisphere is multiplied by a factor of −1 to match with the tilts of the northern hemisphere assuming hemispheric symmetry in tilt distribution.This convention is followed throughout the analysis.

Assigning Physical Parameters to BMRs
In most of the previous studies, every observation of a BMR is treated as a new BMR.However, since the AutoTAB catalog provides the tracked information, we need to find a way to study the statistical properties of BMRs by assigning a unique number for BMR parameters.Thus, to get the single representative parameter for LT class BMRs, we picked the time period during which the measured Φ m exceeds 80% of its maximum value during tracking.Following that, we calculate the average properties only during that period.To avoid projection effects, the maximum flux Φ m of a BMR is considered only when it resides within 55 • E-W.In Figure 1(a), we illustrate the evolution of flux of a typical BMR from the AutoTAB catalog.Here, shaded region represents the time window when Φ m is more than 80% of maximum Φ m , and thus, the physical parameters (B mean , B max , Φ m , γ, and λ) are averaged over this time window.
Using this approach, we can easily assign a single value of a physical parameter for the LT class BMRs.However, assigning a single value to a BMR in the DP class, may not be appropriate as they are tracked at different evolutionary stages and we may not have the maximum Φ m during this period.For the tilt angle, it may be more questionable to assign a single value as the fluctuations in the tilt is more pronounced compared to other parameters (refer to Figure 9 for case studies).
To gauge these fluctuations in tilt angles over BMR's tracking period, we calculate the standard deviation σ t (γ) in tilt for a BMR over its lifetime for all the DP class BMRs.The variation of the mean of σ t (γ) in each flux bin of 10 22 Mx is presented as a function of Φ m in Figure 1  While the fluctuation in the measured tilt over the tracked lifetime (σ t (γ)) is even larger than the mean tilt, the trend in this figure suggests a relatively stable variation with respect to the flux.We however observe a potential decrease in σ t (γ) for higher flux BMRs.The larger error bars at higher flux ranges result from a smaller number of BMRs in these bins.Hence, in this study, DP class BMRs can be uniformly analyzed and the same method for assigning the parameter values for the LT class can be followed for the DP class as well.
Figure 1(b) also suggests that BMRs having higher flux are less affected by convective buffeting, possibly due to quenched convection around them and a larger amount of magnetic tension; more on this will be discussed in Sections 4.2.4 and 4.2.5.
After outlining the method for assigning representative values to each tracked BMR, we proceed to assess the statistical behavior of BMRs.

Statistical Properties of tracked BMRs
We start validating the AutoTAB catalog by studying its statistical properties and comparing them with the anticipated collective behaviors of BMRs based on previous studies (e.g.,Howard 1991Howard , 1996;;Sivaraman et al. 2007).This serves as a preliminary step before investigating the observational signature for the thin flux tube model.

Distribution of tilt
One of the best-known properties of tilt is its distribution.Hence, we start by evaluating whether the assigned tilts of BMRs shows the familiar distribution, as observed in earlier studies by Wang & Sheeley (1991); Dasi-Espuig et al. ( 2010) and many others.In Figure 2(a), we show the distribution of tilts for all the BMRs from all latitudes.From this figure, we note that the tilt shows the well-known Gaussian distribution as reported in earlier works.The least-square fit to the distribution with Gaussian function (shown by blue solid line in Figure 2a) gives mean (µ) of 7.78 • and a standard deviation (σ) of 16.46 • .These values closely align with those obtained by Jha et al. (2020), where they have not tracked BMRs and treated each detection as an independent measurement.We also compared our distribution with Stenflo & Kosovichev (2012) by only considering BMRs in the latitude range of 15 • -20 • in both hemispheres and we find µ = 9.27 • and σ = 14.90 • , which show an excellent agreement.
We further analyzed the tilt distribution for individual cycles (Cycles 23 & 24) and separately for each hemisphere and the fitting parameters are presented in Table 1.Obtained µ and σ values do not show a significant variation.Therefore, for our analyses, we combined the data from two hemispheres.

Joy's law
Another well-known property of BMRs is that they show a systematic dependence on latitude, i.e., Joy's law.Therefore, we plotted the mean tilts of BMRs in each 5 • latitude bin by folding both hemispheres; see Figure 2(b).Here, we calculate the mean by fitting a Gaussian function in the distribution of tilt in each latitude bin.
We fit these mean values with the standard Joy's law function, i.e., γ = γ 0 sin λ + b, (blue solid line), which yields Joy's law slope γ 0 of 29.56 • ± 3.11 • , in agreement with previous reports (Stenflo & Kosovichev 2012;Jha et al. 2020).If we force the line to pass through the origin (i.e., b = 0), the γ 0 gets modified to 28.98 • ±3.10 • .Sometimes, instead of sinusoidal dependence, a linear dependence (γ = m Joy λ + b) is also used for the tiltlatitude relation (e.g., Wang & Sheeley 1991;Sivaraman et al. 1999).In Figure 2(b), we also fit the mean γ with this linear fit, which is shown by the dashed blue line (almost on top of red).The slope of the fitted line is 0.49 when b ̸ = 0 and 0.50 when b = 0.These values are significantly higher than the previously reported values of 0.26 and 0.28, based on white-light observations at Mount Wilson and Kodaikanal Observatories, respectively (Dasi-Espuig et al. 2010).One reason behind this discrepancy is the difference in data type, as magnetogram tends to give a higher slope of Joy's law than the white light (Wang et al. 2015).Since there is no significant difference between the fitted lines in Figure 2b, we use the sinusoidal dependence as the standard Joy's law.

Flux vs magnetic field
Assigned B max and Φ m values of the tracking data from AutoTAB reveal a correlation between them, depicted in a scatter plot of Figure 3.The (Pearson) correlation coefficient between the quantities is 0.72, suggesting a good correlation between B max and Φ m .A similar trend (not shown) is observed between B mean and Φ m as well.It is important to note that the measured B max values are affected by the saturation limit of MDI and HMI, but the saturation effects will not significantly influence the measurement of Φ m of the regions (Hoeksema et al. 2014).Therefore, the observed stabilization of B max with Φ m beyond B max of 3 kG in Figure 3 may reflect this phenomenon and measurements in high field regions will be significantly constrained by instrument saturation.The quadratic fit is slightly better and the relation becomes B max below 3 kG is given by B max = (0.83 ± 0.02)Φ 2 m + (0.27 ± 0.04)Φ m .(If we fit with cubic relation, then the equation reads B max = (0.08 ± 0.03)Φ 3 m + (0.56 ± 0.1)Φ 2 m + (0.48 ± 0.09)Φ m .)The cut-off of 3 kG is set to address the instrument saturation.Given this limitation, the relation in Figure 3 suggests stronger BMRs are associated with higher magnetic flux.

Validation of thin flux tube model
After discussing the collective behaviors of the tracked regions, in this section, we delve into an examination of whether the tracked information from AutoTAB aligns with the expected outcomes proposed by the thin flux tube model as discussed in Section 2.

Evolution of footpoint separation of BMRs
To assess the validity of the rising flux tube model, i.e., the rise of Ω loop, we calculate foot-point separation (D) at each time step of observation throughout the lifetime of the LT class BMRs.Here, D is defined as the angular distance between the two polarities of the BMR at each instance of the tracking of a BMR calculated using, (2) The heliographic location of the polarities is calculated based on the flux density-weighted mean location of BMR's polarity.In Figure 4(a), we show the footpoint separation (D) of BMR polarities as a function of time, normalized by their lifetimes to bring them on the same footing.
In Figure 4(a), we see a steady, rapid increase in D during the initial phase of the evolution of the BMRs, followed by a slow increase in the later stage.However, the evolution of D of a typical BMR, shown in Figure 1 of Kosovichev & Stenflo (2008), depicts a downward trend in the footpoint separation towards the disintegration phase of the BMR, which is not evident in Figure 4(a).To further explore the behavior observed by Kosovichev & Stenflo (2008), we look at a few individual BMRs (Figure 9 in Appendix).We observe that in some cases, there is indeed a downward trend in D towards the later phase of BMR life (see Figure 10(a), (f), (g), (i), (k)).However, we find a statistically consistent increase of D along with saturation in the later part of BMR's lifetime, which aligns with the expectation of the rise of Ω shaped flux tube.
A closer inspection of Figure 4(a) shows a rapid increase in the footpoint separation occurs during the initial phase (5% to 30% of lifetime).This accelerated growth phase in BMR evolution was previously reported as Phase-1 (acceleration) by Schunker et al. (2019).The  4(a) is primarily contributed by Bin-II, and the trend also suggests that higher flux regions are associated with higher D, which is consistent with previous findings (Schunker et al. 2019).The evolution of ∆ϕ in panel (c) reflects the rapid increase seen in D in higher flux ranges (Bin-II and Bin-III); however, in Figure 4(d), such a trend is not evident in the evolution of ∆λ.We also would like to point out that the trend depicted in Figure 4(a) persists across the latitudinal bins, with their patterns overlapping.
It has been suggested in previous simulation studies by D 'Silva & Choudhuri (1993); Schüssler & Rempel (2005) that the flux tubes get tethered from the CZ as they rise up.If the tethering happens, the reflection of the same can be seen as a rapid change in longitude separation and a steady change in latitude separation due to imparted circular motion from the twisted flux tubes.Now, the question arises, when does the tethering happen in a flux tube?One can assume that flux tubes with lower strength and lower flux content will tether from the CZ at an earlier time compared to stronger flux tubes.
The rapid change in ∆ϕ in higher flux bins of Figure 4(c) may be attributed to this effect.However, we fail to find any statistically significant variation in ∆λ in Figure 4(d).An important consideration here is that the BMRs in LOS magnetogram data are detected only once the flux tubes emerge from the photosphere radially.Additionally, AutoTAB tracks them exclusively when both polarities have emerged and they hold a strong flux balance condition.Hence, the dataset lacks information about the onset emergence phase of the BMR.The absence of a rapid increase in lower flux regions could be because tethering had already occurred before AutoTAB started tracking BMRs.This is in contrast to higher flux regions, where tethering occurs at a later stage, allowing AutoTAB to track them during the tethering period ef-fectively and reflected as the rapid increase in ∆ϕ, seen in Figure 4(c).

Tilt angle at the first detection
The evolution of D lends support to the assumption of ascending flux tubes associated with the formation of BMRs.One debated point has been whether the BMRs emerge with definite tilt or they acquire the observed tilt after emergence.According to the thin flux tube model, BMRs should have acquired tilt as the flux tubes rise through the CZ.Hence, we anticipate a definite tilt angle in BMRs at the emergence.To evaluate this, we collected all those BMRs that emerged on the near side of the Sun and looked for their tilt and Joy's law behavior at the very first detection.We emphasize that the latitude and tilt values considered here correspond to the BMR's first detection.To avoid the projection effect, we restrict the BMR emergence between 45 • E-W longitudes.AutoTAB tracks 5635 such BMRs, which lie within the flux range of 10 19 -10 23 Mx with the median at 4.05 × 10 21 Mx.
In Figure 5, first, we plot the Gaussian mean in each latitude bin of 5 • and then fit Joy's law function using the least square fitting method.Here we note that the BMRs show a definite tilt at their first detection, which increases with latitude, as we expect from Joy's law.However, Joy's law amplitude of γ 0 = 27.17• implies a somewhat weaker latitude dependence compared to Figure 2(b).Comparing this γ 0 with that obtained during the matured phases, we find that the tilt increases in the later phase of BMR's life, which is consistent with the previous findings (Kosovichev & Stenflo 2008;Schunker et al. 2019).This is also visible in the individual case studies in Figure 10.Nevertheless, based on the behavior observed in Figure 5, we can say that, statistically, BMRs emerge with a definite tilt in accordance with Joy's law.Schunker et al. (2019Schunker et al. ( , 2020) ) suggested that BMRs emerge with zero tilt and develop tilt in accordance with Joy's law in the later part of their life from the analysis of the Solar Dynamic Observatory Helioseismic Emerging Active Region (SDO/HEAR) survey, thereby arguably ruling out the Coriolis force as a possible cause for the tilt in the BMRs.The active regions chosen for their study mainly had two criteria: (i) The regions should appear in the continuum observation, and (ii) Regions should emerge in the quiet Sun region and not in a region close to existing AR.We argue that these criteria may lead to selection bias, and the result may be affected by this.To validate if the BMR emerges with zero tilt or not, we conducted some case studies of BMR.We remind that AutoTAB tracks the BMR only for the time period in which the flux balance condition holds, and the initial developing signatures of BMRs might have been missed in tracking.Therefore, we select some BMRs at various stages of solar cycles and flux ranges and go back in time (time before they are detected by AutoTAB) to observe their tilt at the emergence phase.While making our selection, we make sure that the latitude of the emergence of these BMRs is between ±30 • and with no prominent flux emergence nearby.Thereafter, we look for the signature of emergence for these BMRs in the magnetogram by going back in time before they are detected by AutoTAB.Snapshots of 12 such BMRs and their evolutions are shown in Figure 9 along with the evolution of Φ, γ, and D in Figure 10.Our observations based on these selected BMRs reveal that 8 of them emerge with definite tilt, while the remaining emerge with nearly zero tilt.Furthermore, we also note that the BMR signatures emerging with or without tilt do not depend on the latitude of their emergence, the time it takes to be detected by AutoTAB, or their flux strength.The factors determining a BMR emerge with or without tilt remain unclear and require further exploration.Nevertheless, it is noted that nearby flux emergence can influence the tilt of the emerging BMRs.
In summary, as observed in case studies, and based on Figure 5, we confirm that Joy's law trend is clearly evident at the first detection of tracking.This suggests that some tilt is imparted to the BMRs below the surface before they are observed in the magnetograms, and the major cause of the BMR tilt can be the Coriolis force.As the BMR flux is observed to show a wide variation in magnitude (more than three orders), we expect to detect a variation of tilt with the magnetic flux, which was one of the predictions of the thin flux tube model as discussed in Section 2. To assess the tilt dependence on flux, we compute γ 0 = ⟨γ⟩/⟨sin λ⟩ in each 10 22 Mx flux bin and plot them against the Φ m in Figure 6(a).We note that in computing γ 0 , a normalization factor ⟨sin λ⟩ is used to get rid of the latitudinal dependency in Joy's law.However, this quantity, γ 0 , is strictly not the slope of Joy's law (γ 0 ).Taking the average on both sides of Equation 1, we can regard this factor γ 0 as the mean slope of Joy's law.Instead of computing it in this way, if we compute γ 0 in the traditional method (i.e., by fitting tilt v/s latitude variation) in each flux bin, then we get a statistically insignificant value of γ 0 due to limited data in some bins (however see the next section).Despite this, we observe a consistent increase in γ 0 with Φ m until 7×10 22 Mx, and the best fit is found to be linear (γ 0 = aΦ m + b).

Flux dependence of tilt
This dependence is much stronger than the one predicted by the thin flux tube model of Fan et al. (1994); see Equation 1.However, we note that the measured total unsigned flux of the BMR might be an averaged representation of flux inside the tube, and a precise simulation-like behavior cannot be anticipated from observational data.Nevertheless, any observed dependence of γ 0 on Φ m indicates a significant dependence of tilt on the flux of the BMR and thus supports the thin flux tube model.After 7 × 10 22 Mx, we observe an indication of a decrease of γ 0 with Φ m (assuming that there is no saturation limit in the measurements).The reliability of these data points in the graph is however compromised due to the limited number of high-flux BMRs.Nevertheless, this decrease could be due to the dominance effect of the magnetic field at a large flux (as B max increases with Φ m ; Figure 3).While the tilt increases with the magnetic flux, it decreases with the field in the flux tube; see Equation 1.In high flux regime, the effect of magnetic field dominates over the flux and it causes to decrease the tilt.This decrease of tilt is referred to as the tilt quenching by Jha et al. (2020) (also see the inclusion of tilt quenching in dynamo models Karak & Miesch 2017, 2018), who observed a decrease in γ 0 with B max (see Figure 4(a) in Jha et al. 2020).In our tracked BMR data also we find a similar trend between B max and γ 0 as observed by Jha et al. (2020)-an initial increase of γ 0 with B max followed by a mild descend beyond 1.75 kG.However, the reliability of the data points for higher B max bins is limited for the same reasons, as highlighted in Section 4.1.3(saturation limit in measurements).
The trends in Figure 6 indicate a dynamic interplay of forces on rising flux tubes.Small and moderate field and flux might represent a region where the drag and Coriolis forces dominate over the magnetic buoyancy, while the regime with B max > 1.75 kG (and Φ m ≳ 7 × 10 22 ) represents the magnetic buoyancy dominated regime, resulting in reduced tilt.

Flux dependence of Joy's law
To further explore the magnetic flux dependence of tilt, here we revisit Joy's law.We have seen in Figure 2 that the tracked BMRs from the AutoTAB catalog collectively obey Joy's law.However, we have already noted that the tilt data is extremely noisy (e.g., Figure 1(b)), and Joy's law trend is evident only after averaging.A possible contribution to the scatter in the tilt is due to the turbulent convection affecting the flux tube as it rises through the CZ (e.g, Longcope & Choudhuri 2002;Weber et al. 2013), and thus, the effect of scatter is more prominent in lower flux bins.Already, we have seen some evidence of it in Figure 1(b).We, therefore, segregate the BMRs into three bins with an equal number of data points in each bin.To keep the same number of data points in each bin, the flux ranges in these bins become: 1.72 × 10 19 -3.13 × 10 21 Mx (Bin-I), 3.13 × 10 21 -2.38 × 10 22 Mx (Bin-II), and 2.38 × 10 22 -9.92 × 10 24 Mx (Bin-III).We note that although the maximum flux of BMR in Bin III goes to a very large value, there are only a few BMRs above 2 × 10 23 Mx; see Figure 3 (in fact, the median flux values in each bin are 6.64 × 10 20 Mx, 9.81 × 10 21 Mx, and 5.28 × 10 22 Mx, respectively.)Figure 7 displays tilt as a function of latitude, along with Joy's law fit, for the BMRs in these bins.We observe that as we move to the higher flux regimes (Bins-II and III), the scatter decreases, which agrees with the theoretical expectation that the stronger BMRs are less buffeted by convection.
Also, seeing a large scatter around Joy's law and a large error in the fitted parameters in Bin-I, having BMR flux < 3.13 × 10 21 Mx, we raise a doubt whether Joy's law behavior is valid only for stronger BMRs.Interestingly, if we discard this bin, then from panels (b) and (c), we observe that the slope of Joy's law γ 0 increases as we move from Bin-II (having mean flux 1.11 × 10 22 Mx) to Bin-III (mean flux 4.93 × 10 23 Mx).This γ 0 increase is in agreement with the theoretical prediction and, in particular, with the result presented in the previous section (Section 4.2.3) that γ 0 = ⟨γ⟩/⟨sin λ⟩ increases with flux.
We note that in Figure 7, the result also remains consistent if we include only the BMRs within the flux range of 10 20 -10 23 Mx (i.e., if we cut the two tails of the flux distributions) and bin the data in equal numbers.Interestingly, if we bin the data in equal flux bin such that the flux ranges become 10 20 -10 21 Mx (Bin-I), 10 21 -10 22 Mx (Bin-II), and 10 22 -10 23 Mx (Bin-III), then the value of γ 0 in Bin-II and Bin-III are comparable.However, the decrease in scatter from Bin I to Bin II remains consistent in this case as well.This suggests that Joy's law fitting is sensitive to how we bin the data points.This could be why previous authors could not find a systematic increase of γ 0 with the flux (Kosovichev & Stenflo 2008;Stenflo & Kosovichev 2012).Therefore, ⟨γ⟩/⟨sin λ⟩ may be a better quantity when measuring the flux dependent of BMR tilt as done in Section 4.2.3.

Variation of tilt fluctuations with footpoint separation
Based on the thin flux tube model, it is expected that as the tubes rise through CZ, they are buffeted by convection.Including turbulence in the numerical model of thin flux tube, Longcope & Fisher (1996) predicted dependence of the tilt angle fluctuations i.e., root-meansquare value of the tilt RMS(γ) with the footpoint separation D. They showed that as D increases, (i) the averaging effect of small-scale wiggles over the flux tube is greater, leading to straight rise of the tube and (ii) the flux increases, which consequently decreases the rise time or the time for the interaction with convection.They showed that RMS(γ) scales as 1/D.We investigate this behavior by showing the variation of RMS(γ) with D from the AutoTAB catalog in Figure 8.We observe that RMS(γ) decreases with increasing footpoint separation exactly reported by Longcope & Fisher (1996, RMS(γ) = a/D + b).We note that from whitelight observations, Fisher et al. (1995) also found a similar dependence: RMS(γ) ≃ 10 • (D/100) −3/4 .

CONCLUSIONS
The thin-flux tube model provides the simplest explanation for the formation of BMRs, postulating magnetically buoyant tubes anchored at the convection zone.According to this model, the Coriolis force imparts tilt to the emerging BMRs.Despite its popularity, the theory lacks observational backing.This study undertakes an investigation into the validity of the model's assumptions using the AutoTAB catalog.
AutoTAB is an in-house developed algorithm to automatically detect and track the BMRs from MDI (1996-2011) and HMI (2010-2023) LOS magnetogram data.The resulting comprehensive catalog documents the evolution of 12,173 BMRs on the nearside of the Sun.As a sanity test, we first assign a single representative values for BMR properties (B max , B mean , Φ m , γ) by averaging their values during the times when their flux exceed 80% of the peak values.We then show that BMRs follow (i) Joy's law trend (γ 0 = 29.56• ), (ii) Gaussian-like tilt distribution (µ = 7.78 ± 0.43), and (iii) magnetic flux vs field dependence, all are consistent with previous studies.With these basic tests, we proceed with detailed analyses of our tracked BMRs to validate the thin flux tube model as a theory behind the BMR formation.
The buoyant rise of flux tube assumption is validated by examining footpoint separation (D) over the BMR's tracked lifetime.Our findings indicate a rapid increase in D during the initial phase of BMR evolutions, followed by a gradual rise and eventual saturation towards the end of their lifetime.The rapid increase is primarily attributed to the longitude separation (∆ϕ), particularly pronounced in higher flux regions.This may imply a connection between rapid D growth and footpoint tethering from the convection zone, where the immediate effect of the same manifests as a sudden increase in ∆ϕ.
In line with the thin-flux tube model, BMRs are expected to appear with tilt at the onset of the emergent phase due to the effect of Coriolis force during the rise of flux tubes.Our analysis of Joy's law trend in tracked BMRs of LT class supports this expectation by demonstrating a clear Joy's law trend in their first detection.This was further evaluated through case studies of individual BMRs selected from different phases of cycles and strengths.Our findings reveal a nuanced scenario, where the signatures of some BMRs emerge with zero tilt and develop it at later phases, while the rest exhibit a significant tilt from the beginning.This tilt behavior is independent of the emergence latitude or cycle phase, indicating a potential contribution of Coriolis force to a part of observed tilt in BMRs.
Further, based on the thin flux tube model, we expect the tilt to increase with the increase of magnetic flux (D 'Silva & Choudhuri 1993;Fan et al. 1994); also see Equation 1.We explored this flux dependence of tilt by computing γ 0 = ⟨γ⟩/⟨sin λ⟩.Our results reveal a linear increase in γ 0 with Φ m , signifying a pronounced tilt dependence on magnetic flux.A similar trend is observed concerning B max , with an initial increase in γ 0 for lower B max regions, followed by a decrease after 1.75 kG.We further observed that Joy's law variation is flux dependent.The Joy's law is less significant (and has a large scatter around the mean trend) at the small flux bin (below about 10 21 Mx).The slope of Joy's law increases as we move the BMR flux bin from 10 21 Mx to 10 23 Mx.This result is again in agreement with the prediction of the thin flux tube model.
Finally, based on the thin flux tube model, we expect that the tilt fluctuation (rms value) to depend on the footpoint separation (e.g, Longcope & Fisher 1996).BMRs having large footpoint separation have less chance to be buffeted by the convection.From our data, we find that the tilt fluctuation decreases inversely with the footpoint separation as predicted by the model of Longcope & Fisher (1996).
Although our analysis provides some support to the thin flux tube model and hints at Coriolis force as the reason behind a part of the tilt observed in the BMR, further study using richer data is needed to strengthen our conclusion.Notably, we have to carefully identify the emergent phase of a BMR and automatically compute the tilt of a large number of BMRs to check the statistical reliability of tilt at the very early phase of BMR.We also need to explore the flux dependence of the tilt angle using a longer dataset.
A.S. sincerely expresses her gratitude to ARIES, Nainital, for the warm hospitality extended during the project's initial phase.B.B.K. acknowledges the financial support from the Department of Science and Technology (SERB/DST), India, through the Ramanujan Fellowship (project No. SB/S2/RJN-017/2018). A.S. and B.B.K. thank Sami Solanki, Natalie Krivova, and Robert Cameron for their valuable comments and discussions during their visit to MPS, which helped in improving the manuscript.The observational data, the LOS magnetograms from MDI and HMI, used in this article were obtained via JSOC, courtesy of the NASA/SOHO and NASA/SDO science teams.

APPENDIX
A. EXAMPLES OF EVOLUTIONS OF BMRS Panels (a)-(f) and (g)-(l) in Figure 9 present snapshots to illustrate the evolution of 12 distinct BMRs tracked by AutoTAB from Solar Cycle 23, starting from the point when the BMR signatures first emerged.The snapshots with red rectangles are produced from the times when AutoTAB could track the BMRs, while the snapshots with blue rectangles show their evolutions back in time.The criteria for region selection were as follows: (i) Individual BMRs' emerging signatures should fall within ±55 • longitude and ±30 • latitude, and (ii) There should be no significant flux emergence in the vicinity.The "T" value at the bottom indicates the evolutionary stage of the BMR."T:0" signifies the initial detection by AutoTAB.Pre-and post-detection times are indicated in hours, denoted by negative and positive signs.The time and date of first detection is mentioned at the top of the panel of "T:0".Bracketed numbers denote the mean heliographic latitude and longitude for each evolutionary stage.
Among the selected regions, we observed that eight BMRs emerged with a significant tilt at the initial phase itself, while others emerged with nearly zero tilt.The evolution of γ, Φ m , and D during their tracked lifetime for all the selected BMRs are illustrated in panels (a)-(l) in Figure 10. (b).

Figure 1 .
Figure 1.(a) Evolution of magnetic flux (Φm) of a representative BMR tracked by AutoTAB.The shaded area indicates the duration during which the measured Φm is greater than 80% of the maximum Φm recorded for the BMR.(b) Mean of the standard deviation of BMR tilt, σt(γ) in each 5 × 10 21 Mx flux bin plotted against Φm of the DP class BMRs.

Figure 2 .
Figure 2. (a) Tilt Distribution: Number of BMRs in 5 • tilt bins are shown as bars.The blue solid line represents the Gaussian fitted curve (with an offset) with fitting parameters mentioned in the panel.The vertical solid blue line represents the 0 • tilt, and the dashed blue line represents the Gaussian fitted mean at 7.78 • .(b) (Gaussian) Mean tilt in each 5 • latitude bin as a function of the latitude.Blue solid and red dashed lines represent Joy's law (γ = γ0 sin λ + b) and straight line fits (γ = mJoyλ + b) with fitting parameters mentioned in the panel along with χ 2 value for the fit at the bottom right.The numbers appearing below the points display the total number of BMRs in the associated bins.

Figure 3 .
Figure 3. Scatter plot between the measured Bmax and Φm.The dashed blue line represents the quadratic fit (see text).Flux values are capped at 3×10 23 Mx, and the outliers beyond this range attribute to the effects of defective pixels and faulty detection.

Figure 4 .
Figure 4. (a) General evolutions of footpoint separation of BMR (D) over its lifetime.(b) Shows the same but in three different flux ranges: Bin-I: 10 20 -10 21 Mx (blue), Bin-II: 10 21 -10 22 Mx (red), and Bin-III: 10 22 -10 23 Mx (green).(c) and (d) show the same as in (b) but for the longitude separation (∆ϕ) and latitude separation (∆λ), respectively.Note that the time is normalized to their lifetimes, and the quantities (D, ∆ϕ, and ∆λ) are averaged over all BMRs in each 96 minute cadence.rate of rise in footpoint separation slows down in the later phase, suggesting that this can be Phase-2 (deceleration) as suggested by Schunker et al. (2019).To further investigate this behavior, we segregate the BMRs into three different flux ranges based on the assigned flux, Bin-I: 10 20 -10 21 Mx, Bin-II: 10 21 -10 22 Mx, and Bin-III: 10 22 -10 23 Mx and we study the evolution of D and the major contributors to D, i.e., ∆ϕ, and ∆λ in each of the flux bins in Figure 4(b), (c), and (d), respectively.Figure 4(b) represents the evolution of D in various flux bins.The increase in D observed in Figure4(a) is primarily contributed by Bin-II, and the trend also suggests that higher flux regions are associated with higher D, which is consistent with previous findings(Schunker et al. 2019).The evolution of ∆ϕ in panel (c) reflects the rapid increase seen in D in higher flux ranges (Bin-II and Bin-III); however, in Figure4(d), such a trend is not evident in the evolution of ∆λ.We also would like to point out that the trend depicted in Figure4(a) persists across the latitudinal bins, with their patterns overlapping.It has been suggested in previous simulation studies by D'Silva & Choudhuri (1993);Schüssler & Rempel (2005) that the flux tubes get tethered from the CZ as Figure 5. Joy's law plot at the first detection of BMR emerging between 45 • E-W.(Gaussian) Mean tilt in each 5 • latitude bin as a function of the latitude.Blue dashed lines represent Joy's law fit (γ = γ0 sin λ + b) with parameters and χ 2 value for the fit mentioned in the panel.Numbers appearing below the points mark the total number of BMRs in the associated bins.

Figure 6 .
Figure 6.(a) Mean γ0 is calculated in each flux bin (Φm) of length 1 × 10 22 Mx and is plotted as a function of Φm.The dashed blue line represents the straight line fit (γ0 = aΦm + b), excluding the last two points with parameters mentioned on top of the panel.Here, the error bars represent the standard error in each flux bin, which is bigger in the high flux bins because of less number of BMRs.(b) Mean γ0 is calculated in each Bmax bin of length 0.5 kG and is plotted against Bmax.Dashed black vertical line represents Bmax at 1.75 kG.

Figure 8 .
Figure 8. Mean root-mean-square of the tilt (RMS(γ)) calculated in each footpoint separation (D) bin of length 25 Mm and is plotted as a function of D.

Figure 9 .
Figure 9. Snapshots of the evolution of 12 BMRs tracked by AutoTAB with 1 day cadence.See text for details. in the brackets at the bottom of panels represent the mean latitude and longitude of polarities of the BMR.

Figure 10 .
Figure 10.Time evolutions of flux, tilt, and footpoint separation of 12 BMRs whose evolution are shown in Figure 9.

Table 1 .
Fitting parameters of the tilt distribution for different cycles and hemispheres.