Issue 
A&A
Volume 685, May 2024



Article Number  A16  
Number of page(s)  11  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/202348791  
Published online  30 April 2024 
OGLE2023BLG0836L: The sixth microlensing planet in a binary stellar system
^{1}
Department of Physics, Chungbuk National University,
Cheongju
28644,
Republic of Korea
email: cheongho@astroph.chungbuk.ac.kr
^{2}
Astronomical Observatory, University of Warsaw,
Al. Ujazdowskie 4,
00478
Warszawa,
Poland
^{3}
Korea Astronomy and Space Science Institute,
Daejon
34055,
Republic of Korea
^{4}
Max Planck Institute for Astronomy,
Königstuhl 17,
69117
Heidelberg,
Germany
^{5}
Department of Astronomy, The Ohio State University,
140 W. 18th Ave.,
Columbus,
OH
43210,
USA
^{6}
University of Canterbury, Department of Physics and Astronomy,
Private Bag 4800,
Christchurch
8020,
New Zealand
^{7}
Department of Particle Physics and Astrophysics, Weizmann Institute of Science,
Rehovot
76100,
Israel
^{8}
Center for Astrophysics  Harvard & Smithsonian
60 Garden St.,
Cambridge,
MA
02138,
USA
^{9}
Department of Astronomy and Tsinghua Centre for Astrophysics, Tsinghua University,
Beijing
100084,
PR China
^{10}
School of Space Research, Kyung Hee University,
Yongin,
Kyeonggi
17104,
Republic of Korea
^{11}
Korea University of Science and Technology,
217 Gajeongro,
Yuseonggu,
Daejeon,
34113,
Republic of Korea
^{12}
Department of Physics, University of Warwick,
Gibbet Hill Road,
Coventry
CV4 7AL,
UK
Received:
30
November
2023
Accepted:
12
February
2024
Aims. Light curves of microlensing events occasionally deviate from the smooth and symmetric form of a singlelens singlesource event. While most of these anomalous events can be accounted for by employing a binarylens singlesource (2L 1S) or a singlelens binarysource (1L2S) framework, it is established that a small fraction of events remain unexplained by either of these interpretations. We carried out a project in which data collected by highcadence microlensing surveys were reinvestigated with the aim of uncovering the nature of anomalous lensing events with no proposed 2L 1S or 1L 2S models.
Methods. From the project we found that the anomaly appearing in the lensing event OGLE2023BLG0836 cannot be explained by the usual interpretations, and we conducted a comprehensive analysis of the event. From thorough modeling of the light curve under sophisticated lenssystem configurations, we arrived at the conclusion that a triplemass lens system is imperative to account for the anomalous features observed in the lensing light curve.
Results. From the Bayesian analysis using the measured observables of the event timescale and angular Einstein radius, we determined that the least massive component of the lens has a planetary mass of 4.36_{−2.18}^{+2.35} M_{J}. This planet orbits within a stellar binary system composed of two stars with masses 0.71_{−0.36}^{+0.38} M_{⊙} and 0.56_{−0.28}^{+0.30} M_{⊙}. This lensing event signifies the sixth occurrence of a planetary microlensing system in which a planet belongs to a stellar binary system.
Key words: planets and satellites: detection
© The Authors 2024
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1 Introduction
The light curve of a microlensing event involving a single lens and a single source (1L1S) is represented by (1)
where A(t) is the lensing magnification, F_{s} and F_{b} denote the respective flux values originating from the source and blended light components, and represents the projected lenssource separation normalized to the angular Einstein radius θ_{E}. The lensing magnification varies in time as the lenssource separation changes due to their relative motion as (2)
where 0 and t_{0} represent the minimum lenssource separation (scaled to θ_{E}) and the corresponding time, and t_{E} is the Einstein timescale. The resulting light curve is characterized by a smooth and symmetric form (Paczynski 1986).
Light curves in microlensing events occasionally exhibit deviations from the standard 1L1S form. These deviations are, in most cases, attributed to two primary factors: the potential binarity of the lens, as described by Mao & Paczynski (1991), and the binarity of the source, as noted by Griest & Hu (1992). In the case of a binarylens (2L1S) event, the lensing system creates a complex pattern of caustics. These caustics represent specific positions of the source at which the lensing magnification for a point source diverges to infinity. When a source crosses the caustic, a pair of new images are created or disappear, resulting in a complicated lensing light curve that deviates from the 1L1S form. In the case of a binarysource (1L 2S) event, the lensing magnification corresponds to the mean of the magnifications associated with the individual binary source stars, A_{1} and A_{2}, weighted by the flux contributions of the component source stars, F_{1} and F_{2}: (3)
As a consequence, the light curve of a 1L2S event displays deviations from the standard 1L1S form.
Since 2016, the Korea Microlensing Telescope Network (KMTNet) has been conducting a microlensing survey with frequent observations of stars located in the direction of the Galactic bulge (Kim et al. 2016). Among about 3000 microlensing events that are annually detected from the survey, about 10% of the events exhibit anomalies in the lensing light curves. While the majority of these anomalous lensing events can be explained by applying a 2L1S or a 1L2S framework, it is known that a small fraction of events defy explanation under either of these interpretations. The challenge in precisely characterizing the peculiarities of these events hints at the necessity for more advanced models to interpret the observed anomalies.
We conducted a project that involved revisiting microlensing data collected by the KMTNet survey. The primary goal of this project was to uncover instances of anomalous lensing events for which the conventional 2L1S or 1L2S models had not been previously proposed. Through this investigation, we identified multiple occurrences that required the application of advanced modeling approaches beyond the standard 2L1S or 1L2S frameworks. Han et al. (2019) found that the anomaly appearing in the light curve of the lensing event OGLE2018BLG1011, which corresponds to KMTNet event KMT2018BLG2122, could be explained with a triplelens (3L1S) model in which the lens is composed of two giant planets and their host star. Through a detailed analysis of the light curve for the lensing event OGLE2018BLG1700 (KMT2018BLG2330), Han et al. (2020a) identified the triple nature of the lens by decomposing the anomaly into two parts produced by two binarylens pairs. In one of these binary pairs, the mass ratio of the lens components is approximately q ~ 0.01, while in the other pair the mass ratio is around ~0.3, suggesting that the lens is a planetary system in a binary. Through a careful examination of the central anomaly observed in the lensing curve of the highly magnified event KMT2019BLG1953, Han et al. (2020b) found that the discrepancies from the 2L1S model were significantly reduced when an additional planetary lens companion or a source companion was introduced, although distinguishing these two interpretations was difficult within the precision of the photometric data. In their study, Han et al. (2021a) determined that the anomalous characteristics observed in the lensing light curve of the event KMT2019BLG0797 could be accounted for by a 2L2S model in which both the lens and source are binary systems. By analyzing the event KMT2019BLG1715, for which the lensing light curve exhibited two shortterm deviation features from a causticcrossing 2L1S light curve, Han et al. (2021b) suggested a fivebody (lens+source) model, in which one deviation feature was generated by a planetarymass third body of the lens, and the other feature was generated by a faint source companion, and thus the event is a very complex fivebody system composed of three lens masses (planet + two stars) and two source stars. Han et al. (2021c) found that KMT2018BLG1743 is another planetary lensing event occurring on two source stars. In their analysis of the anomalies observed in the lensing event OGLE2019BLG0304 (KMT2019BLG2583), Han et al. (2021d) put forward two rivaling interpretations between a 3L1S and a 2L2S models. The 3L1S model suggests the presence of a planetarymass third body situated near the primary lens of the binary lens system. The 2L2S model proposes the existence of an additional nearby companion to the source. Zang et al. (2021) found that the central anomaly in the lensing light curve of the event KMT2020BLG0414 was produced by a triplelens system, which consists of an Earthmass planet and its binary host. Through the investigation of the anomalous lensing event KMT2021BLG1077, Han et al. (2022a) identified that the lens of the event is a multiplanetary system in which two giant planets orbit a very lowmass star. It was found by Han et al. (2022b) that the dual bump anomaly feature in the light curve of the lensing event KMT2021BLG1898 could be explained with a 2L2S model, in which the lens contains a giant planet and the source is a binary composed of a turnoff star and a Ktype dwarf. Han et al. (2022c) found that the planetary signal in the lensing light curve of the event KMT2021BLG0240 was deformed either by an extra planetary lens component or by a companion to the source, although the 3L1S and 2L2S interpretations could not be distinguished with the available data. Han et al. (2023a) found that the lensing events OGLE2018BLG0584 and KMT2018BLG2119 were generated by 2L2S lens systems. From the analysis of the lensing event KMT2021BLG1122, Han et al. (2023b) revealed that the anomaly appeared in the light curve was produced by a 3L1S system, which consists of three stars.
In this study, we provide a comprehensive analysis of the microlensing event OGLE2023BLG0836/KMT2023BLG1144, for which no existing model has successfully explained the anomaly observed in the lensing light curve. The anomaly in question presents two distinctive features: a causticcrossing pattern and a strong cuspapproaching peak. Our investigation has led us to the conclusion that the inclusion of a triplemass lens system is imperative to adequately account for all the anomalous features in the lensing light curve of the event.
2 Observations and data
The source of the microlensing event OGLE2023BLG0836/KMT2023BLG1144 is positioned in the direction of the Galactic bulge field, with equatorial coordinates (RA, Dec)_{J2000} = (17:48:44.85, −23:44:29.47), which correspond to the Galactic coordinates (l, b) = (4°.8075, 2°.0912). In this direction, the extinction in the Iband is approximately A_{I} ~ 2.01. The magnification of the source flux caused by lensing was first detected from the Optical Gravitational Lensing Experiment IV (OGLEIV: Udalski et al. 2015) survey on 28 June 2023, corresponding to the reduced Heliocentric Julian Date (HJD′ ≡ HJD − 2460000 = 123). Five days later, the KMTNet group independently found the event and designated it as KMT2023BLG1144. While we initially recognized the anomalous nature of the event through a systematic analysis of the KMTNet data gathered during the 2023 season, we have chosen to label the event as OGLE2023BLG0836, aligning with the reference ID from the OGLE survey that initially detected the event.
The event observations were conducted using the telescopes operated by individual survey groups. The KMTNet group employs three identical telescopes, each featuring a 1.6m aperture and a widefield camera capable of capturing 4 square degrees in a single exposure. In order to ensure continuous coverage of lensing events, the KMTNet telescopes are strategically positioned across the three continents in the southern hemisphere: Siding Spring Observatory in Australia (KMTA), Cerro Tololo Interamerican Observatory in Chile (KMTC), and South African Astronomical Observatory in South Africa (KMTS). Additionally, the OGLE telescope, featuring a 1.3m aperture and a camera yielding a 1.4 square degree field of view, is located at Las Campanas Observatory in Chile. Images from the KMTNet and OGLE surveys were mainly obtained in the Iband with the inclusion of some Vband images taken for source color measurement. Observations of the event by the KMTNet and OGLE surveys were done with cadences of an hour and about two days, respectively. Image reduction and photometry for the lensing event were accomplished using automated pipelines tailored to the individual surveys, with those developed by Albrow (2017) for the KMTNet survey and by Woźniak (2000) for the OGLE survey. For the use of optimal data in the analysis, we conducted rereduction of the KMTNet data using the photometry code developed by Yang (in prep.). Additionally, the error bars of the data were adjusted to ensure consistency with the data scatter, and to set the χ^{2} value per degree of freedom (d.o.f.) for each data set to unity following the procedure outlined in Yee et al. (2012).
Figure 1 presents the lensing light curve of OGLE2023BLG0836, revealing deviations from the typical smooth and symmetric shape observed in a 1L1S event. These deviations are marked by two primary features. The first is the causticcrossing feature, which is characterized by a pair of causticcrossing spikes at HJD′ ~ 115.7 and ~119.6 and a Ushape trough region between the spikes. The second feature is the strong peak, centered at HJD′ ~ 111.3, which is likely to be produced by the source star’s approach to a caustic cusp. The upper panel of Fig. 1 offers a detailed view of this anomalous region. Because caustics in gravitational lensing are formed due to the presence of multiple lensing masses, the causticcrossing feature indicates that the lens is comprised of multiple masses. Furthermore, the approximate symmetry between the ascending and descending segments of the peak anomaly feature suggests that the feature likely originated from the source’s approach to a cusp of a caustic.
Fig. 1 Lensing light curve of OGLE2023BLG0836. The top panel displays an enlarged view focused on the vicinity of the anomaly. 
3 Lensing light curve analysis
3.1 Binarylens analysis
Taking into account the potential involvement of the anomaly features with caustics, our analysis commences by modeling the light curve under the interpretation with a 2L1S lenssystem configuration. This modeling process is conducted to identify a lensing solution, which comprises a set of lensing parameters that best describe the characteristics of the light curve. Under the approximation of rectilinear relative motion between the lens and source, the light curve of a 2L1S event is defined by seven fundamental lensing parameters. The three initial parameters (t_{0}, u_{0}, t_{E}) characterize the approach of the source to the lens. The next three parameters (s, q, α) provide information about the binary lens system: s denotes the projected separation (scaled to θ_{E}) between the lens components M_{1} and M_{2}; q = M_{2}/M_{1} represents the mass ratio of these lens components; and α is the source trajectory angle, defined as the angle between the direction of the relative lenssource proper motion vector µ and the M_{1}−M_{2}, axis. The last parameter, ρ, is defined as the ratio of the angular source radius θ_{*} to the Einstein radius (normalized source radius), that is, ρ = θ_{*}/θ_{E}, and it characterizes how finite source effects contribute to the deformation of the lensing light curve during caustic crossings.
In our pursuit of the lensing solution, we categorized the lensing parameters into two groups. Within the first group, which pertains to the binary parameters (s, q, α), we conducted a gridbased exploration to determine the values of parameters s and q with multiple initial values of α. The other parameters of the second group were determined by minimizing χ^{2} through the use of the Markov chain Monte Carlo (MCMC) method, which employs an adaptive step size Gaussian sampler, as described in Doran & Mueller (2004). To assess the presence of degenerate solutions, we examined the ∆χ^{2} map in the s–q parameter space derived from the grid search. After identifying local solutions, we further refined the lensing parameters of the solutions using a downhill approach. In cases for which the discrepancies in χ^{2} values among these local solutions were marginal, we presented all of them and investigated the causes of degeneracies.
Despite our comprehensive exploration of the parameter space, we were unable to identify a 2L1S solution capable of sufficiently explaining both the causticcrossing and cuspapproaching features within the anomaly. This underscores the necessity of employing a more sophisticated model for a comprehensive understanding of the observed anomaly features.
While the 2L1S model cannot simultaneously accommodate both the causticcrossing and cuspapproaching features, we identified that each anomaly feature can be adequately approximated by a 2L1S model. This is illustrated in Fig. 2, where we present two 2L1S model curves depicting the individual anomaly features. The curve presented in the upper panel represents a 2L1S model derived by fitting the data excluding the observations around the cuspapproaching feature during the time interval 110 ≲ HJD′ ≲ 115, while the curve shown in the lower panel is a model obtained by fitting the data with the exclusion of the data around the caustic feature within the time range 114 ≲ HJD′ ≲ 110. Another example of an anomaly that can be divided into two parts produced by two 2L1S events can be seen in Fig. 2 of Han et al. (2020a). We found a solution with binary parameters (s, q) ~ (0.60, 0.85) from the 2L1S fit to the data excluding the causticcrossing feature. Similarly, we also identified a solution with (s, q) ~ (1.13, 3.4 × 10^{−3}) from the fit to the data excluding those around the cuspapproaching feature. The full lensing parameters of the pair of 2L1S solutions are presented in Table 1. Regarding the lensing parameters in the model describing the causticcrossing feature, we note that the mass ratio of the lens components is on the order of 10^{−3}, indicating that the companion is very likely to be a planetmass object.
Figure 3 displays the lenssystem configurations corresponding to the individual 2L1S models presented in Fig. 2. In each panel, the red figure comprising the concave closed curves represents the caustic, and the arrowed line represents the trajectory of the source. The configuration illustrates that the causticcrossing feature resulted from the source star’s crossing over the planetary caustic induced by a planetary companion. On the other hand, the peak feature of the anomaly was produced by the source star’s close approach to the sharp offaxis cusp of a caustic induced by a binary lens composed of roughly equal masses.
Fig. 2 Division of the anomaly into two parts. The curve presented in the upper panel represents a 2L1S model derived by fitting the data excluding those around the cuspapproaching feature, while the curve in the lower panel is a 2L1S model obtained by fitting the data with the exclusion of the data around causticcrossing feature. 
Lensing parameters of 2L1S solutions.
3.2 Triplelens analysis
Bozza (1999) and Han (2001) pointed out that anomalies produced by a triplelens system, composed of three masses (M_{1}, M_{2}, M_{3}), can often be approximated by combining the anomalies induced by the two binary pairs M_{1}−M_{2} and M_{1}−M_{3} through superposition. Then, the fact that the two anomaly features in the peak region of the OGLE2023BLG0836 light curve are well approximated by two 2L 1S models implies the possibility of the lens system being a triple system. Consequently, we proceeded with a modeling approach based on a 3L 1S configuration for the lens system.
The 3L1S configuration corresponds to the case in which an extra lens component, M_{3}, is present in addition to the 2L 1S configuration. Incorporating this supplementary lens component necessitates the addition of extra lensing parameters in the modeling procedure. These parameters consist of (s_{3}, q_{3}, ψ), which respectively stand for the projected separation and mass ratio between M_{1} and M_{3}, and the orientation angle of M_{3} as measured from M_{1}−M_{2} axis. In order to distinguish these parameters describing M_{3} from those describing M_{2}, we use the notations (s_{2}, q_{2}) to designate the separation and mass ratio between the M_{1}−M_{2} pair.
The 3L1S modeling was carried out using the following procedure. In the first step, we searched for the tertiary lens parameters (s_{3}, q_{3}, ψ) via a grid approach using the parameters of the 2L 1S model as initial values for the other parameters. In our analysis we used the lensing parameters of the 2L 1S solution depicting the cuspapproaching feature. In the second step we identified local solutions appearing in the s_{3}−q_{3}−ψ parameter space, and then refined the individual solutions by gradually minimizing χ^{2} of the fit using a downhill approach.
Through the 3L1S modeling, we identified three distinct solutions resulting from the ambiguity in s_{3}. Figure 4 shows the locations of the individual local solutions in the ∆χ^{2} map on the (log s_{3}, log q_{3}) parameter plane obtained from the grid searches for the tertiary lens parameters. The parameters describing the M_{1}−M_{2} pair, which lie in the ranges of 0.5 ≲ s_{2}≲ 0.6 and 0.59 ≲ q_{2} ≲ 0.88, are similar to those of the 2L1S model describing the cuspapproaching feature. Similarly, the parameters describing the M_{1}−M_{3} pair, which lie in the ranges of 1.11 ≲ s_{3} ≲ 1.16 and 4.9 × 10^{−3} ≲ q_{3} ≲ 5.9 × 10^{−3}, are similar to those of the 2L1S model describing the causticcrossing feature. We label the individual local solutions as “inner,” “intermediate,” and “outer” for the reason discussed below. In the ∆χ^{2} map, it appears that there are two local minima with log q values approximately around 2.5. However, during the refinement of the solutions, we found that these two minima converge into a single solution.
In Table 2, we list the full lensing parameters of the three 3L1S solutions together with their χ^{2} values of the fits and d.o.f. Our analysis reveals a preference for the outer solution, with χ^{2} differences of 34.0 and 4.2 over the inner and intermediate solutions, respectively. We present the model curve of the outer 3L1S solution in the top panel of Fig. 5 and the residuals from all three 3L1S solutions in the lower panels. It is worth noting that the normalized source radius is measured based on the finitesource constraint, although the measured value carries a substantial uncertainty. The lensing light curve is affected by finitesource effects during the source crossings over the fold caustics induced by the planet and during the approach to the cusp induced by the binary companion. We checked the origin of the ρ constraint by conducting two modeling runs; the first run was done by excluding the two points near the caustic crossings at HJD′ = 115.685 and 119.559, and the second run was conducted by excluding the two KMTC points near the cusp peak at HJD′ = 121.513 and 121.527. From these runs we find that the constraint on ρ comes mainly from the two data points observed during the caustic crossings. In Fig. 6, we plot the relative probability of the normalized source radius estimated from the 3L1S modeling.
The configurations of the lens systems corresponding to the three 3L1S solutions are presented in Fig. 7. In all instances, the triplelens caustic seems to result from the combination of two separate binarylens caustics caused by the M_{1}–M_{2} and M_{1}–M_{3} pairs. The source first traversed the planetary caustic created by the M_{1} – M_{3} pair, giving rise to the causticcrossing feature. Subsequently, the source moved past the lower tip of the caustic formed by the M_{1} –M_{2} pair, leading to the emergence of the cuspapproaching anomaly feature. While the fundamental structures of the three solutions resemble one another, there exist subtle distinctions between them. According to the inner and outer solutions, the source traversed the two folds of the planetary caustic situated in the inner and outer regions between the caustic center and the primary lens, respectively. As we discuss below, the degeneracy between the inner and outer solutions is caused by the inner–outer degeneracy (Gaudi & Gould 1997). In the intermediate solution, the source encountered the inner caustic fold upon entering the caustic, and the outer caustic fold while departing from it. We assigned labels to the individual solutions based on the side of the caustic that the source crossed.
Hwang et al. (2022) and Gould et al. (2022) showed that the planet separations of a pair of solutions resulting from an innerouter degeneracy follow the relation (4)
Here s_{in} and s_{out} respectively denote the binary separations of the inner and outer solutions; ; and t_{anom} is the time of the anomaly. The symbols “+” and “−” on the right side of the equation apply to the majorimage and minorimage perturbations, respectively. In the case of OGLE2023BLG0836, the planetinduced anomaly (i.e., the causticcrossing feature) was produced by a majorimage perturbation, and as a result the sign associated with this event is “+.” In order to investigate the origin of the degeneracy in the separation s_{3}, we checked whether s_{3} values of the inner and outer solutions follow the relation in Eq. (4). With the lensing parameters (t_{0}, u′_{0}, t′_{E}, t_{anom}, s′_{in}, s′_{out}) = (115.37, 0.22, 46.2, 117.0, 1.435, 1.408), we find that the geometric mean (s_{in} × s_{out})^{1/2} ~ 1.124 closely matches the value . This confirms that the similarity between the model curves of the inner and outer solutions stems from the innerouter degeneracy. In computing s^{†}, we used the values of the lensing parameters normalized to the Einstein radius corresponding to sum of M_{1} and M_{3}: (u′_{0}, t′_{E}, s′_{in}, s′_{out}) ≡ (u_{0}f, t_{E}/f, s_{in}f, s_{out}f), where f = (1 + q_{3})^{1/2}.
After identifying the origin of the degeneracy between the inner and outer solutions, we further investigate the origin of the degeneracy involving the intermediate solution and the other solutions. It is worth noting that Shin et al. (2023) and Han et al. (2024) independently reported the presence of such degeneracies in their respective analyses of the planetary lensing events KMT2016BLG1751 and KMT2023BLG0469. Upon closer examination of the causticcrossing feature, we found that the degeneracy between the intermediate solution and the other solutions is an accidental one stemming from inadequate caustic coverage. This is depicted in Fig. 8, which provides a detailed view around the causticcrossing feature. The figure illustrates that the source magnitudes just before the caustic entrance and just after the caustic exit for the inner and outer solutions are remarkably similar, suggesting an inherent degeneracy. On the contrary, the source magnitude just before the caustic entrance for the intermediate solution is approximately 0.1 magnitude brighter than what is expected from the other solutions, indicating that this degeneracy is a chance occurrence. If the light curve in this specific region had been more thoroughly covered, the degeneracy between the intermediate solution and the other solutions could have been lifted. The issue of degeneracies arising due to insufficient coverage of certain parts of causticcrossing features was explored by Skowron et al. (2018).
Fig. 3 Lens system configurations of the two 2L1S solutions for which the model curves are presented in the corresponding panels of Fig. 2. In each panel, the red figure is the caustic, the line represents the source trajectory, and the arrow on the source trajectory indicates the direction of the source motion. By convention, the abscissa of these 2L1S diagrams is defined by the binary axis of the lens system. 
Fig. 4 ∆χ^{2} map on the (log s_{3}, log q_{3}) parameter plane. The colorcoding is configured to correspond to data points based on their ∆χ^{2} values: red for ∆χ^{2} ≤ 1^{2}n, yellow for ≤ 2^{2}n, green for ≤ 3^{2}n, and cyan for ≤ 4^{2}n, where n = 4. The three distinct local solutions are marked as inner, intermediate, and outer. The dashed vertical line represents the geometric mean (s_{3,in} × s_{3,out})^{1/2} of the planetary separations for the inner and outer solutions. 
Lensing parameters of 3L1S solutions.
Fig. 5 Model curves of the bestfit 3L1S solution (outer solution) and 2L2S solution in the region of the anomaly. The lower four panels present the residuals from the three degenerate 3L1S solutions (inner, intermediate, and outer solutions) and the 2L2S solution. 
Fig. 6 Relative distribution of the normalized source radius. The solid vertical line indicates the median value, and the two dotted lines represent the 1σ ranges of the distribution. 
3.3 Binarylens binarysource analysis
As illustrated in the lensing event KMT2021BLG0240 (Han et al. 2022c), anomalies arising from a 3L1S system can, on occasion, be confused with anomalies resulting from 2L2S systems. Hence, we conducted a more thorough examination to ascertain whether the observed anomalies in the lensing light curve of OGLE2023BLG0836 could be explained by a 2L2S interpretation.
The 2L2S configuration involves the presence of an additional source alongside the 2L1S configuration. We designate the primary and secondary source stars as S_{1} and S_{2}, respectively. To account for the additional source, it is necessary to incorporate extra lensing parameters into the modeling process. These additional parameters encompass (t_{0,2}, u_{0,2}, ρ_{2}, q_{F}), which respectively represent the time and separation of the closest approach of S_{2} to the lens, the normalized source radius of S_{2}, and the flux ratio between S_{2} and S_{1}. Concurrently, we employ the notations (t_{0,1}, u_{0,1}, ρ_{1}) to define the parameters describing the approach of S_{1} to the lens. During the modeling process we seek the 2L2S parameters by exploring various trajectories for S_{2}, building upon the 2L1S solution that describes the causticcrossing feature.
In Table 3, we list the bestfit lensing parameters of the 2L2S solution. The lenssystem configuration corresponding to the solution is shown in Fig. 9. In this illustration the paths of the primary and secondary source stars are indicated by black and blue lines, respectively, and they are labeled as S_{1} and S_{2}. The model curve and residual of the 2L2S solution in the region of the anomaly are presented in Fig. 5, showing that the model approximately describes the anomaly features.
Although the 2L2S model offers an approximate description of the anomaly features, it is found that the model exhibits a less accurate fit to the data when compared to the 3L1S model. This can be seen in Fig. 10, where the cumulative distribution of χ^{2} difference between the 2L2S and 3L1S solutions, that is, , is presented. From the distribution, it is found that the 2L2S solution yields a poorer fit in the region after the cuspapproaching feature. Overall, it is found that the 2L2S solution is less favorable by ∆χ^{2} = 31.3 when compared to the 3L 1S solution. As a result, we dismissed the 2L2S interpretation for the anomaly.
Fig. 7 Lenssystem configurations of the three degenerate 3L1S solutions. For each solution, the lower panel shows the source trajectory, marked by an arrowed line, with respect to the caustic, and the upper panel shows the trajectory with respect to the positions of the lens components, marked by blue dots with labels M_{1}, M_{1}, and M_{3}. The gray curves encompassing the caustic in the lower panels represent equimagnification contours. The dashed circles in the upper panels represent the Einstein ring. 
Fig. 8 Zoomedin view around the causticcrossing feature. Model curves of the three degenerate solutions (outer, intermediate, and inner) are drawn in the top panel, and the residuals from the models are shown in the lower panels. 
Bestfit parameters of 2L2S solution.
Fig. 9 Lenssystem configuration of the 2L2S solution. The notations are the same as in Fig. 6, except that there are two source trajectories. The black and blue lines represent the trajectories of the primary and secondary source stars, respectively. 
Fig. 10 Cumulative distribution of χ^{2} difference between the 2L2S and 3L 1S solutions. The light curve in the upper panel is provided to illustrate the region of disparity in the fit. 
Fig. 11 Positions of the source and red giant clump (RGC) centroid in the instrumental color–magnitude diagram, created by merging observations from KMTC (gray dots) and HST (brown dots). 
4 Source star and angular Einstein radius
In this section, we specify the source star of the lensing event and estimate the angular Einstein radius. The source star is specified by estimating its color and magnitude after being corrected for reddening and extinction. With the angular source radius θ_{*}, deduced from the color and magnitude together with the normalized source radius ρ measured from the modeling, the angular Einstein radius is estimated as (5)
Figure 11 shows the location of the source in the instrumental colormagnitude diagram (CMD) of stars lying near the source. The CMD was created by merging two sets of CMDs. One set was generated from the pyDIA photometry (Albrow 2017) of stars in the KMTC image, and the other set pertains to stars in Baade’s window observed using the Hubble Space Telescope (Holtzman et al. 1998). The alignment of these two CMDs was achieved by using the centroids of the red giant clump (RGC) in the individual CMDs. We used the combined CMD because the Vband source magnitude could not be determined even though the Iband magnitude of the source was measured. As a result, the color of the source was estimated as the median value observed in the mainsequence branch of the combined CMD corresponding to the measured Iband magnitude.
To estimate the dereddened source color and magnitude, denoted as (V − I, I)_{0}, from their corresponding instrumental values, denoted as (V − I, I), we employed the Yoo et al. (2004) method. This method utilizes the RGC centroid, for which its dereddened color and magnitude (V − I, I)_{RGC,0} = (1.060, 14.308) were established by previous studies (Bensby et al. 2013; Nataf et al. 2013), as a reference point for calibration: (6)
Here (V − I, I)_{RGC} represent the instrumental color and magnitude of the RGC centroid, and ∆(V − I, I) = (V − I, I) − (V − I, I)_{RGC} indicate the offsets in color and magnitude between the source and RGC centroid. In Table 4 we list the measured values (V − I, I), (V − I, I)_{rgc}, (V − I, I)_{rgc,0}, and (V − I, I)_{0}. From the estimated dereddened color and magnitude it is found that the source is a late Gtype mainsequence star.
For the estimation of the angular source radius, we first converted the measured V − I color into V − K color using the colorcolor relation established by Bessell & Brett (1988), and subsequently calculated the angular source radius by applying the relation between V − K and θ_{*} provided by Kervella et al. (2004). The derived value of the angular source radius is (7)
This yields the angular Einstein radius and relative lenssource proper motion of (8)
respectively. We note that the uncertainties associated with θ_{E} and µ are relatively large, primarily stemming from the large uncertainty in ρ.
A degeneracy between two competing interpretations with widely different lensing parameters of (t_{E}, ρ) can occasionally be lifted from the resulting values of the relative lenssource proper motion if one model results in a µ value that is relatively disfavored by the Galactic model. We checked the feasibility of this method by additionally estimating the relative proper motion expected from the 2L2S model. The estimated value of µ_{2L2S} ~ 6.0 mas yr^{−1}, which is similar to the value for the 3L1S model, and is also a very typical value for a Galactic lensing event. Therefore, the constraint from the estimated relative proper motion cannot used to distinguish between the two interpretations.
Source parameters.
5 Physical lens parameters
In this section, we estimate the physical parameters of the lens system. The lens parameters of the mass M and distance D_{L} are constrained by three lensing observables: the event timescale t_{E}, the angular Einstein radius θ_{E}, and the microlens parallax π_{E}.
The microlens parallax is defined as the ratio of the relative lenssource parallax π_{rel} to the angular Einstein radius, (10)
where D_{S} denotes the distance to the source. With the measurements of all these observables, the mass and distance to the lens are uniquely determined as (11)
where κ = 4G/(c^{2}au) and π_{S} = au/D_{S} denotes the parallax of the source (Gould 2000). For OGLE2023BLG0836 the observables t_{E} and θ_{E} were constrained, while the accurate determination of microlens parallax was hindered by limited photometric data precision. Consequently, we estimated M and D_{L} by conducting a Bayesian analysis, leveraging the constraints provided by the measured observables t_{E} and θ_{E} together with the priors of the physical and dynamical distributions and mass function of lens objects in the Galaxy.
The Bayesian analysis was initiated by generating a large number of synthetic events through a Monte Carlo simulation. In the simulation, the physical parameters of the lens mass were derived from a model mass function, and the distances to the lens and source as well as the relative proper motion between them were obtained from a Galaxy model. We adopted the mass function model proposed by Jung et al. (2018) and the Galaxy model introduced by Jung et al. (2021). Using physical parameters (M_{i}, D_{L,i}, D_{S,i}, µ_{i}) for each synthetic lensing event, we computed the corresponding lensing observables, specifically the event timescale and angular Einstein radius, according to the following relations: (12)
Then, the posteriors of the lens mass and distance were subsequently constructed by assigning a weight to each event of (13)
Here (t_{E}, θ_{E}) represent the measured values of the observables, and [σ(t_{E}), σ(θ_{E})] denote their respective uncertainties.
In Fig. 12, we present the Bayesian posteriors of the primary lens mass M_{1} and the distances to the lens and source, D_{L} and D_{S}, of OGLE2023BLG0836. In Table 5 we list the masses of the individual lens components (M_{1}, M_{2}, M_{3}), distance to the lens, and the projected separations of the individual lens companions from the primary (a_{⊥,2}, a_{⊥,3}). The physical parameters were derived based on the observables of the outer solution, which provided the best fit to the data. Given the similarity of the observables, the physical parameters derived from the other solutions are similar to the presented values. For each physical parameter, we provide the median of the Bayesian posterior as a representative value and estimate the uncertainty range as the 16% and 84% of the distribution. The least massive component of the lens has a mass M_{3} ~ 4.4 M_{J}, which classifies it as a giant planet. The masses of the other lens components M_{1} ~ 0.71 M_{⊙} and M_{2} ~ 0.56 M_{⊙} correspond to those of mid and late Ktype mainsequence stars. Given that the planetary separation between M_{1} and M_{3}, a_{⊥,3} ~ 3.7 au, is not much greater than that between the M_{1} and M_{2} binary pair, a_{⊥,2} ~ 1.9 au, it is likely that the planet orbits both binary components. However, the possibility of the planet orbiting one of the binary components cannot be entirely ruled out because of the projected nature of the separations. The estimated distance to the lens is D_{L} ~ 5.1 kpc, and the probabilities for the lens lying in the disk and bulge are 66% and 34%, respectively.
Fig. 12 Bayesian posteriors of the primary lens mass (M_{1}) and the distances to the lens (D_{L}) and source (D_{S}) of the lensing event OGLE2023BLG0836. Within each panel the event contributions from the disk and bulge lenses are illustrated by blue and red curves, and the black curve represents the combined distribution of the two lens populations. The solid vertical line marks the median value, and the dotted lines indicate the 1σ range of the posterior distribution. 
Physical lens parameters.
6 Summary and conclusion
We conducted an analysis of the peculiar lensing event OGLE2023BLG0836, for which the light curve is characterized by two distinctive anomaly features produced by a caustic crossing and a cusp approach of a source. Despite the comprehensive exploration of the parameter space, we were unable to identify a binarylens solution capable of sufficiently explaining both features within the anomaly.
From the analysis with a sophisticated model prompted by the fact that each anomaly feature can be approximated by a 2L1S model, we arrived at the conclusion that a triplemass lens system is necessary to account for the observed anomaly pattern in the lensing light curve. A binarylens binarysource interpretation could also offer an approximate explanation for the anomaly pattern, but this interpretation was rejected with a high degree of statistical confidence. Through the detailed triplelens modeling, we identified three distinct solutions resulting from the degeneracy in the separation between the primary and the least massive companion of the lens.
Through a Bayesian analysis using the measured observables of the event timescale and angular Einstein radius, we determined that the least massive component of the lens has a planetary mass of . This planet orbits within a stellar binary system composed of two stars with masses and . This lensing event signifies the sixth occurrence of a microlensing planetary system in which a planet belongs to a stellar binary system.
Acknowledgements
Work by C.H. was supported by the grants of National Research Foundation of Korea (2019R1A2C2085965). J.C.Y. and I.G.S. acknowledge support from U.S. NSF Grant No. AST2108414. Y.S. acknowledges support from BSF Grant No. 2020740. This research has made use of the KMTNet system operated by the Korea Astronomy and Space Science Institute (KASI) at three host sites of CTIO in Chile, SAAO in South Africa, and SSO in Australia. Data transfer from the host site to KASI was supported by the Korea Research Environment Open NETwork (KREONET). This research was supported by KASI under the R&D program (project no. 2023183203) supervised by the Ministry of Science and ICT. W. Zang, H.Y., S.M., R.K., J.Z., and W. Zhu acknowledge support by the National Natural Science Foundation of China (Grant No. 12133005). W. Zang acknowledges the support from the HarvardSmithsonian Center for Astrophysics through the CfA Fellowship. The OGLE has received funding from the National Science Centre, Poland, grant MAESTRO 2014/14/A/ST9/00121 to A.U.
References
 Albrow, M. 2017, https://doi.org/18.5281/zenodo.268849 [Google Scholar]
 Bennett, D. P., Rhie, S. H., Udalski, A. et al. 2016, AJ, 152, 125 [Google Scholar]
 Bensby, T., Yee, J. C., Feltzing, S., et al. 2013, A&A, 549, A147 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bessell, M. S., & Brett, J. M. 1988, PASP, 100, 1134 [Google Scholar]
 Bozza, V. 1999, A&A, 348, 311 [NASA ADS] [Google Scholar]
 Doran, M., & Mueller, C. M. 2004, J. Cosmology Astropart. Phys., 09, 003 [CrossRef] [Google Scholar]
 Gaudi, B. S., & Gould, A. 1997, ApJ, 486, 85 [Google Scholar]
 Gould, A. 2000, ApJ, 542, 785 [NASA ADS] [CrossRef] [Google Scholar]
 Gould, A., Udalski, A., Shin, I.G., et al. 2014, Science, 345, 46 [Google Scholar]
 Gould, A., Han, C., Weicheng, Z., et al. 2022, A&A, 664, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Griest, K., & Hu, W. 1992, ApJ, 397, 362 [Google Scholar]
 Han, C., Chang, H.Y., An, J. H., & Chang, K. 2001, MNRAS, 328, 986 [NASA ADS] [CrossRef] [Google Scholar]
 Han, C., Udalski, A., Gould, A., et al. 2017, AJ, 154, 223 [Google Scholar]
 Han, C., Bennett, D. P., Udalski, A., et al. 2019, AJ, 158, 114 [Google Scholar]
 Han, C., Lee, C.U., Udalski, A., et al. 2020a, AJ, 159, 48 [Google Scholar]
 Han, C., Kim, D., Jung Y. K., et al. 2020b, AJ, 160, 17 [NASA ADS] [CrossRef] [Google Scholar]
 Han, C., Udalski, A., Lee, C.U., et al. 2020c, AJ, 162, 203 [Google Scholar]
 Han, C., Lee, C.U., Ryu, Y.H., et al. 2021a, A&A, 649, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Udalski, A., Kim, D., et al. 2021b, AJ, 161, 270 [NASA ADS] [CrossRef] [Google Scholar]
 Han, C., Albrow, M. D., Chung S.J., et al. 2021c, A&A, 652, A145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Udalski, A., Lee, C.U., et al. 2021d, AJ, 162, 203 [NASA ADS] [CrossRef] [Google Scholar]
 Han, C., Gould, A., Bond, I. A., et al. 2022a, A&A, 662, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Gould, A., Kim, D., et al. 2022b, A&A, 663, A145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Kim, D., Yang, H., et al. 2022c, A&A, 664, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Udalski, A., Jung Y. K., et al. 2023a, A&A, 670, A172 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Jung, Y. K., Gould, A., et al. 2023b, A&A, 672, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Han, C., Jung, Y. K., Bond, I. A., et al. 2024, A&A, 683, A115 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Holtzman, J. A., Watson, A. M., Baum, W. A., et al. 1998, AJ, 115, 1946 [Google Scholar]
 Hwang, K.H., Zang, W., Gould, A., et al. 2022, AJ, 163, 43 [NASA ADS] [CrossRef] [Google Scholar]
 Jung, Y. K., Udalski, A., Gould, A., et al. 2018, AJ, 155, 219 [Google Scholar]
 Jung, Y. K., Han, C., Udalski, A., et al. 2021, AJ, 161, 293 [NASA ADS] [CrossRef] [Google Scholar]
 Kervella, P., Thévenin, F., Di Folco, E., & Ségransan, D. 2004, A&A, 426, 29 [Google Scholar]
 Kim, S.L., Lee, C.U., Park, B.G., et al. 2016, JKAS, 49, 37 [Google Scholar]
 Kuang, R., Zang, W., Jung Y. K., et al. 2022, MNRAS, 516, 1704 [Google Scholar]
 Mao, S., & Paczynski, B., 1991, ApJ, 374, L37 [NASA ADS] [CrossRef] [Google Scholar]
 Nataf, D. M., Gould, A., Fouqué, P., et al. 2013, ApJ, 769, 88 [Google Scholar]
 Paczynski, B. 1986, ApJ, 304, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Poleski, R., Skowron, J., Udalski, A., et al. 2014, ApJ, 795, 42 [Google Scholar]
 Shin, I.G., Yee, J. C., Zang, W., et al. 2023, AJ, 166, 104 [NASA ADS] [CrossRef] [Google Scholar]
 Skowron, J., Ryu, Y.H., Hwang, K.H., et al. 2018, Acta Astron., 68, 43 [Google Scholar]
 Udalski, A., Szymanski, M. K., Szymanski, G., et al. 2015, Acta Astron., 65, 1 [NASA ADS] [Google Scholar]
 Wozniak, P. R. 2000, Acta Astron., 50, 42 [Google Scholar]
 Yee, J. C., Shvartzvald, Y., GalYam, A., et al. 2012, ApJ, 755, 102 [Google Scholar]
 Yoo, J., DePoy, D. L., GalYam, A., et al. 2004, ApJ, 603, 139 [Google Scholar]
 Zang, W., Han, C., Kondo, I., et al. 2021, Res. Astron. Astrophys., 21, 239 [CrossRef] [Google Scholar]
All Tables
All Figures
Fig. 1 Lensing light curve of OGLE2023BLG0836. The top panel displays an enlarged view focused on the vicinity of the anomaly. 

In the text 
Fig. 2 Division of the anomaly into two parts. The curve presented in the upper panel represents a 2L1S model derived by fitting the data excluding those around the cuspapproaching feature, while the curve in the lower panel is a 2L1S model obtained by fitting the data with the exclusion of the data around causticcrossing feature. 

In the text 
Fig. 3 Lens system configurations of the two 2L1S solutions for which the model curves are presented in the corresponding panels of Fig. 2. In each panel, the red figure is the caustic, the line represents the source trajectory, and the arrow on the source trajectory indicates the direction of the source motion. By convention, the abscissa of these 2L1S diagrams is defined by the binary axis of the lens system. 

In the text 
Fig. 4 ∆χ^{2} map on the (log s_{3}, log q_{3}) parameter plane. The colorcoding is configured to correspond to data points based on their ∆χ^{2} values: red for ∆χ^{2} ≤ 1^{2}n, yellow for ≤ 2^{2}n, green for ≤ 3^{2}n, and cyan for ≤ 4^{2}n, where n = 4. The three distinct local solutions are marked as inner, intermediate, and outer. The dashed vertical line represents the geometric mean (s_{3,in} × s_{3,out})^{1/2} of the planetary separations for the inner and outer solutions. 

In the text 
Fig. 5 Model curves of the bestfit 3L1S solution (outer solution) and 2L2S solution in the region of the anomaly. The lower four panels present the residuals from the three degenerate 3L1S solutions (inner, intermediate, and outer solutions) and the 2L2S solution. 

In the text 
Fig. 6 Relative distribution of the normalized source radius. The solid vertical line indicates the median value, and the two dotted lines represent the 1σ ranges of the distribution. 

In the text 
Fig. 7 Lenssystem configurations of the three degenerate 3L1S solutions. For each solution, the lower panel shows the source trajectory, marked by an arrowed line, with respect to the caustic, and the upper panel shows the trajectory with respect to the positions of the lens components, marked by blue dots with labels M_{1}, M_{1}, and M_{3}. The gray curves encompassing the caustic in the lower panels represent equimagnification contours. The dashed circles in the upper panels represent the Einstein ring. 

In the text 
Fig. 8 Zoomedin view around the causticcrossing feature. Model curves of the three degenerate solutions (outer, intermediate, and inner) are drawn in the top panel, and the residuals from the models are shown in the lower panels. 

In the text 
Fig. 9 Lenssystem configuration of the 2L2S solution. The notations are the same as in Fig. 6, except that there are two source trajectories. The black and blue lines represent the trajectories of the primary and secondary source stars, respectively. 

In the text 
Fig. 10 Cumulative distribution of χ^{2} difference between the 2L2S and 3L 1S solutions. The light curve in the upper panel is provided to illustrate the region of disparity in the fit. 

In the text 
Fig. 11 Positions of the source and red giant clump (RGC) centroid in the instrumental color–magnitude diagram, created by merging observations from KMTC (gray dots) and HST (brown dots). 

In the text 
Fig. 12 Bayesian posteriors of the primary lens mass (M_{1}) and the distances to the lens (D_{L}) and source (D_{S}) of the lensing event OGLE2023BLG0836. Within each panel the event contributions from the disk and bulge lenses are illustrated by blue and red curves, and the black curve represents the combined distribution of the two lens populations. The solid vertical line marks the median value, and the dotted lines indicate the 1σ range of the posterior distribution. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.