Insight-HXMT observations of jet-like corona in a black hole X-ray binary MAXI J1820+070

A black hole X-ray binary produces hard X-ray radiation from its corona and disk when the accreting matter heats up. During an outburst, the disk and corona co-evolves with each other. However, such an evolution is still unclear in both its geometry and dynamics. Here we report the unusual decrease of the reflection fraction in MAXI J1820+070, which is the ratio of the coronal intensity illuminating the disk to the coronal intensity reaching the observer, as the corona is observed to contrast during the decay phase. We postulate a jet-like corona model, in which the corona can be understood as a standing shock where the material flowing through. In this dynamical scenario, the decrease of the reflection fraction is a signature of the corona's bulk velocity. Our findings suggest that as the corona is observed to get closer to the black hole, the coronal material might be outflowing faster.


Introduction
During an outburst, a black hole (BH) X-ray binary usually displays transition between hard and soft states, according to the spectral properties of its radiation [1][2][3][4] . In the hard state, usually defined as the photon index Γ < 2 between 2-10 keV, the observed radiation primarily comes from the Comptonization by the hot electrons in the corona, which dominates over the weak, low-energy blackbody radiation from the disk. In the soft state with the photon index Γ > 2, the observed radiation then is characterized by a strong (disk) blackbody component below ∼ 10 keV and a weak, high-energy tail of ∼ 25% of the total bolometric luminosity 2, 5-8 . MAXI J1820+070 (ASASSN-18ey) is a low-mass black hole X-ray binary, newly discovered in X-rays with MAXI 9 on 2018 March 11 10 . In addition to X-ray, the source has also been observed in optical [11][12][13][14][15] and in radio 16,17 . Low-frequency quasi-periodic oscillations (LFQPOs) were found in both X-ray and optical bands [18][19][20] . The measurement of the radio parallax indicates that this source is located at a distance of 2.96 ± 0.33 kpc away from us 21 . Follow-up X-ray observations since its outburst were carried out by other X-ray telescopes, e.g., Swift 22 , NuSTAR 23 , and NICER 24 .
The long-term and high cadence observation of MAXI J1820+070 by Hard X-ray Modulation Telescope (called as Insight-HXMT) 25 was carried from 2018-03-14 (MJD 58191) to 2018-10-21 (MJD 58412) 26 . (MJD 58412). In this decay, the source mainly stayed in the soft state in which the LE photon rates dominated over ME/HE photon rates. The same Insight-HXMT observations have been used recently to study the timing properties of the outburst and LFQPOs were discovered above 200 keV, which was interpreted as due to the precession of a twisted compact X-ray jet 27 .
Spectral and timing analysis of MAXI J1820+070 from MJD 58198 to MJD 58250 based on NICER

3/45
reported that the profile of broad Fe K emission line was almost unchanged, indicating the inner disk is very close to the black hole (about 2 R g , where R g = GM/c 2 is the gravitational radius) and keeps constant 28 . Moreover, the thermal reverberation (soft) lags evolved to higher frequencies. Given these observational properties, it is suggested that it may be the corona that evolves during the outburst of this source, rather than the inner accretion disk, i.e., the corona might be contracting with time 28 . The detection energy band of NICER is 0.5-12 keV. However, the reflection component dominates over the X-ray spectra around 20-50 keV. The spectral fitting at this energy range could put crucial constraint on the reflection parameters 29 , e.g., the reflection fraction and reflection strength, which can be used to probe the inner geometry of the accretion flow 30 . Moreover, the high energy cutoff of the X-ray spectrum can be used to measure the electron temperature which is the important physical parameter to study the evolution of the accretion flow in the hard state. Insight-HXMT observes the broad-band X-rays from 1-250 keV, which enables us to probe the evolution of the accretion flow in detail by spectral and timing analysis. In this paper, we aim to reveal how the physical properties of accretion flow (corona+disk) evolve with time, mainly concentrating on the period when the radiation is dominated by high energy photons most likely from the corona, by fitting the Insight-HXMT spectrum.
In this work, we show the evolution of the accretion flow in the hard state when the emission is dominated by high energy photons from the corona, by analysing the observations of the first outburst of MAXI J1820+070 with Insight-HXMT. The derived reflection fraction R f is found to increase with time in the rise phase and decrease with time in the decay phase. In the scenario of a dynamical corona, this may suggest that as the corona contracts, i.e., the dissipation region gets closer to the black hole, the coronal material might be outflowing faster.

Results
We found that all spectra of MAXI J1820+070 in the first outburst can be approximately characterized by an asymptotic power law with cutoff at high energy ≤ 100 keV. The broad continuum is thought to originate from the thermal Comptonization of the seed photons from the accretion disk by the hot plasma (hereafter corona) near the black hole. Moreover, our preliminary spectral analysis to those spectra with the simple cut-off powerlaw clearly reveals the line features at about 3-10 keV and the hump above 20 keV,

4/45
with respect to the thermal Comptonization (see Fig. 2 and qualitative spectral analysis section in methods).
These features may be related to the illumination to the accretion disk. The Comptonized photons which illuminate the disk are reprocessed within the disk and then reflected off to the observer, resulting in the reflection with the characteristic iron line and hump features 31 . Therefore, we applied relxillCp 32 to account for the reflection model, which is the standard reflection model taking the relativistic effect into account. This reflection model contains both the direct emission from the corona and its reflection on the disk, which will be used in this work as the basic spectral model to fit the spectra of MAXI J1820+070 in the first outburst. The parameter settings of the relxillCp model is given in configuration of the spectral fitting model section in methods.
relxillCp allows to fit the broad iron Kα emission line around 3-10 keV. However, in the residuals, we detect a clear narrow core to the iron emission (∼ 6.4 keV), which was also found in NICER data 28 and NuSTAR data 33 . This narrow component may also originate from the reflection, but due to photoionization of neutral material or at least colder gas further from BH. So we additionally include another nonrelativistic reflection model xillverCp to account for this narrow line component. Moreover, we include diskbb 34 to account for the accretion disk radiation at low energy (∼ 2-3 keV), and add tbabs 35 to accounts for the low galactic extinction 36 of N H ∼ 1.5 × 10 21 cm −2 . Therefore, the fitting modeling is applied in the spectral fits.
The spectra of 70 epochs, from March 14, 2018 (MJD 58191) to June 17, 2018 (MJD 58286), are fitted with the model described above. The best-fitting results are obtained by implementing a Markov Chain Monte-Carlo (MCMC) algorithm, and all spectra can be fitted well with the reduced χ 2 < 1.0 (see spectral fitting method section in methods). In order to illustrate the spectral fitting, we plot one spectrum (MJD 58204, ObsID=P0114661006) with the best-fitting model in Supplementary Figure 1, with the decomposition into the disk blackbody component, the Comptonization and reflection components from the reflection models. The evolution of the best-fitting parameters of interests for these good fits are plotted in Fig. 3. The electron temperature decreases from ∼ 230 to 50 keV in the rise phase, and starts to

5/45
slowly increase to ≤ 100 keV in the decay phase. This tendency is in agreement with the results based on the combined spectra of MAXI/GSC and Swift/BAT 37,38 .
Moreover, the best-fitting reflection fraction R f of the X-ray point source steeply increases up to about 0.5 in the rise phase and then slowly decreases to about 0.1 in the decay phase. In the reflection model, relxillCp, the reflection fraction is defined as the ratio of the coronal intensity that illuminates the disk to the coronal intensity that reaches the observer. The reflection fraction as well as some other parameters, e.g., ionization parameters, iron abundance and the emissivity profile, have impacts on the shape of the reflection spectrum. This means that, the increasing fraction of photons from the X-ray source illuminate the disk in the rise phase while the decreasing fraction of photons from the X-ray source illuminate the disk in the decay phase. In the hard state of GX 339-4, the reflection fraction is also found to be positively correlated with the source luminosity 30 . Assuming the X-ray source to be point-like static lamppost, the height of the X-ray point source can uniquely correspond to the reflection fraction 29 , given the configuration of model parameters, i.e., black hole spin, inner/outer radius of the disk. In this case, in the decay phase, the decrease in the reflection fraction corresponds to the increases of the height.
However, the spatial extent of the corona which was estimated by the timing analysis of NICER data covering the early decay studied in this work, decreases or contracts with time 28 . It was previously discussed that the corona above the disk could be effectively ejected away from the accretion disk by the pressure of the reflected radiation, if the Comptonizing source is dominated by e ± pairs 39 . This scenario has successfully explained both the weakness of the reflection with the reflection fraction R f < 1.0 and hard X-ray spectrum with the photon index Γ ∼ 1.6, in active galactic nuclei and X-ray binaries. When the corona, if being assumed as the point-like lamppost for simplicity, is accelerated outflowing away from the black hole, the beaming effect will reduce the illuminating flux towards the accretion disk, which reduces the reflection fraction. In general, the lamppost geometry should be interpreted as a jet base, possibly a standing shock 40 which determines the corona position, but with material flowing through the standing shock. So the system should be characterized by two independent parameters: its position and bulk velocity. The corona position and its evolution with time could be constrained by the timing analysis 28 . In addition to that, the bulk velocity and its evolution with time can be constrained by the spectral analysis. In order to investigate the effect of the bulk motion of the corona on the reflection 6/45 fraction, we calculate the reflection fraction by applying the package relxilllpionCp of relxill(v1.4.0), with the X-ray point source being located at height H above a black hole and outflowing with the velocity β = υ/c. The light bending effect is taken into account. In Fig. 4, we plot the reflection fraction R f as a function of height H and the bulk motion velocity β . It can be seen that, the reflection increases as the height decreases, as we expect in the normal lamppost model. In the mean time, for a given height, The corona in the lamppost geometry, should be interpreted as a jet base, possibly a standing shock 40 , and was discussed in a number of previous studies of AGN and black hole X-ray binary in the low/hard state 41 ; this probably provides a physical interpretation to the compact X-ray jet revealed by the timing analysis of the same Insight-HXMT data used here 27 . The evolution of the outflowing corona should be studied in terms of both its position and bulk velocity, which can help us to better understand the acceleration/deceleration of the particles near a black hole. The X-ray temporal/timing analysis can provide us the spatial information of the relative geometry of the corona. The broad band spectral fittings which derived variations of the reflection fraction, help to probe the dynamical properties of the corona, i.e., if the corona is outflowing/inflowing. Therefore, this work discovered that the corona outflows faster as it contracts towards to black hole, suggesting that the hard X-ray emission region in the jet base or 7/45 standing shock is formed at closer distance to the black hole for faster outflow. This scenario could be applied to other black hole X-ray binaries and AGNs in which the coronas are active, at least to diagnose if the outflowing/inflowing corona is accelerating or decelerating.

On the evolution of the spectral parameters
In the relxillCp model, the emissivity profile is assumed to be broken powerlaw of two indexes q 1 and q 2 with the break radius R br . In our spectral fits, for simplicity, we assume the index q 2 = 3 for the region where R > R br . The emissivity profile is roughly flattened with the index q 1 < 1.0 within the broken radius R br ∼ 20 − 60R ISCO (where R ISCO ∼ 1.23R g for a = 0.998). Given the predictions by ref 42 (e.g., their Fig. 10 and 12), the flattened emissivity profiles suggest that the X-ray corona source might be spatially extended with the size of tens of R g in radius. The emissivity profile for a point source is predicted to be much steeper with the index q ∼ 6 − 8 in the inner region (e.g, R < 3R g ) before the profile flattens off (where R < R br ). Such a twice-broken power law profile was observed in the narrow-line Seyfert 1 galaxy 1H 0707-495 by fitting the relativistically broadened emission lines from the disk 43 . The derived small index with q 1 < 1 in the inner disk region from our spectral fits may be due to: (1) relxillCp uses a single broken-powerlaw, so the index q 1 is the emissivity-weighted value for the region R < R br ; (2) the index q 1 and the inner radius are somewhat degenerated with each other. On one hand, if the emissivity profile falls steeply with the index q 1 ∼ 6 − 8, then in order to reproduce the flux levels emitted from the innermost regions, the disk must be truncated at a larger radius. On the other hand, if the emissivity profile is flat with smaller index q 1 < 1, then the disk can extend to the innermost regions with a smaller truncation radius, to match the flux levels 42 . In our spectral fits, the inner radius is assumed to be small with R in = R ISCO , for simplicity. If we assume R in = 20R g instead and refit the spectra, the emissivity profiles then are required to be steep with the index q 1 ∼ 6 − 8.
The reflection fraction in relxillCp is smaller than unity, which requires the spatially extended corona to be relativistically outflowing to overcome the increasingly important light-bending effect, since the corona is suggested to be contracting from NICER lag-observations 28  The best-fitting values of the iron abundance are high with A Fe ≥ 5, and the fits prefer changing in iron abundance which roughly increases along with the decay, resulting in an apparent correction with the reflection fraction with the p-value of 7.5×10 −21 . Such supersolar abundances also appear in the spectral fitting of other black hole X-ray binaries, e.g., GX 339-4 44 and Cyg X-1 45 . It was suggested that the derived supersolar abundances are not physical, but may be caused by the limitation of the assumed density in the reflection model 46,47 . Both the density and the iron abundance of the accretion disk have impacts on the emergent reflection spectrum. As the density increases, the rise in free-free absorption leads to an increase in temperature, causing extra thermal emission at soft X-ray energy E < 10 keV 46 .
This effect could become significant when the density is high (n e > 10 19 cm −3 ) (see Fig. 7 in ref 47 ), although it was emphasized that the microphysics in the current reflection models is only known to be accurate up to n e = 10 19 cm −3 (see the reference 46 for more details). The effect of the abundance on the reflection spectrum is studied, using the relxillCp model. The increase in the iron abundance from A Fe = 1 to A Fe = 10, will also cause the increase in the soft X-ray emission of the reflection for the highly ionized disk in the hard state (see Supplementary Figure 4a), mimicking the effect of the disk density increase.
Therefore, the soft X-ray emission could be reproduced by either high density or high abundance in the accretion disk 46,47 . It was shown that, the high density model with the solar value of iron abundance, or the low density model with supersolar abundance, could both fit the NuSTAR spectra of Cygnus X-1 47 . In the current reflection model relxillCp which is used in this work, the constant density is fixed at constant value n e = 10 15 cm −3 , which is much lower than the typical values in the standard thin disk model 48 ,

9/45
n e ≥ 10 19 cm −3 . Along with the decay of the outburst studied in this work, the radiation pressure in the inner disk should decrease and thus the gas density should increase 46,49 . Therefore, the spectral fits with the much lower and constant density in the reflection model should lead to artificially supersolar and also increasing iron abundance along with the outburst decay, consistent with the spectral fitting results in this work. Deriving the physical values and evolution of the iron abundance and disk density with the spectral fitting, however, requires new atomic data in high-density to be implemented in the reflection model, which is beyond the scope of this work. Nonetheless, in order to demonstrate that the observed trend of the reflection fraction is physical and independent of the expected artificial correlation between the reflection fraction and the iron abundance, we tried the case of the constant iron abundance, e.g., A Fe = 5.0. In this case, the reflection fraction still decreases along with the outburst decay (see Supplementary Figure 4b), although fixing A Fe is not supported by F-test, due to the problems we discussed above.
The best-fitting values of the ionization parameter are high with log ξ ≥ 3.8. The ionization parameter is defined as the ratio of the illuminating flux and the electron number density (n e = 10 15 cm −3 ). According to this definition, changing the ionization parameter will not only redistribute the reflected photons over energy (i.e., reshaping the X-ray spectroscopy), but also the number of the reflected photons (i.e., increasing/decreasing the spectrum flux). Therefore, in the current relxill model, the reflection spectrum is normalized to the constant illuminating flux F 0 = ξ n e /4π where ξ = 1 erg cm s −1 , and the illuminating flux (which takes the reflection fraction into account) from the nthcomp is normalized to this constant illuminating flux. Then, the total spectrum can be simply calculated by combining the reflection spectrum with the nthcomp spectrum. After normalizing the spectrum in the relxill model, the ionization parameter is not proportional to the illuminating flux 29,32 . In photoionization modelling of the relxill model, it is assumed that log ξ completely characterizes the X-ray spectroscopy, regardless of the actual values of density or flux 46 . Therefore, the highly ionized disk in Supplementary Figure 2a is not in conflict with the weak illuminating flux, although the faster outflowing corona provides a weaker illuminating flux. Highly ionized disk for low flux at 3-78 keV was also found in the NuSTAR spectrum of MAXI J1820+070 33 .
Instead, the large values of the ionization parameter is possibly caused by the over-estimateed iron abundance as discussed above. The 0.01-10 keV flux would be reduced due to the increase of the continuum opacity, as the iron abundance increases to account for the deficiency of the constant density in 10/45 the reflection model 46 . In order to compensate for this effect in the spectral fits, the degree of the ionization in the disk would increase. This is because in the case of high ionization, the illumination continuum will not be highly absorbed by the photoelectric opacity, which leads to the relatively strong reflection continuum in 0.01-10 keV.
The broad component of the iron line visually appears stable throughout the decay, while the strength of the narrow core decreases with time, which is consistent with the observations of NICER and NuSTAR 28,33 .
The broadening/narrowing of the broad iron line can also be quantitatively evaluated with the equivalent width (EW), a measure of the relative strength of the line profile. The formula of calculating the EW of the iron line (including narrow component) is as follows: where F(E) is the total flux and Furthermore, as we discussed above, the X-ray spectroscopy (including the relative strength of the iron line) by relxillCp does not depend on the actual illuminating flux, EW then is insensitive to the reflection fraction. As discussed above, the broad Fe Kα line appears stable throughout the decay, while the relative strength of the narrow core decreases with time. In order to understand the evolution of the the X-ray spectroscopy of Fig. 2  Supplementary Figure 6 plots the time evolution of the diskbb parameters, i.e., the temperature at inner disk radius T in in units of keV, the normalization and the diskbb flux in units of erg/cm 2 /s. The disk inner temperature increases from 0.4 keV to 0.6 keV during the decay, although there is significant degeneracy between the inner temperature and the normalization of diskbb. We note that, the seed photon temperature kT bb in relxillCp is fixed at 0.05 keV, which is much lower than the best-fitting values of T in . In order to study the effect of kT bb on the results, we freeze kT bb = 0.5 keV. It turns out that the difference in kT bb does not have a strong effect on other parameters. We use cflux in XSPEC to estimate the diskbb flux in 0.01 -100.0 keV and its time-evolution. We found that the disk flux increases with time until around MJD = 58250, and then evidently decreases until the end of the decay. In our spectral model, the inner radius of the accretion disk is assumed to be constant at ISCO. We note that the inner radius can be derived from the diskbb parameters, e.g., GX 339-4 30 and XTE J1550-564 50 . In ref 50 , the inner radius was estimated with the photon flux of the direct disk component and the Comptonized disk component, taking into account the Comptonized disk photons. In this work, we also estimate the disk inner radius from the spectral fits.
Although the observed Comptonization and the reflection flux can be fitted with the reflection model, it is hard to determine the Comptonized disk flux, which depends on the geometry of the corona with respect to the disk. Therefore, we introduce a constant factor λ which is multiplied by the direct disk flux F d to estimate the Comptonized disk flux F c , i.e., F c = λ F d . Then, the inner radius can be calculated as follows (see the Appendix in ref 50 , where the source distance D = 2.96 kpc, i = 63 • are adopted from the observational measurements 21   on the other hand, the Compton amplification factor can increase due to the increase of the outflowing velocity β , leading to the decrease of the photon index (see Supplementary Figure 7b). In the decay of MAXI J1820+070, if the effect of the geometrical factor dominates over that of the outflowing velocity, the resultant photon index will increase, as found in our spectral fitting.

13/45
On the effect of system parameters In our spectral fits, the black hole spin is fixed at a = 0.998. However, it is also possible that the black hole spin is low, if the accretion disk extends to the innermost stable circular orbit, i.e., R in = R ISCO , provided that the inner radius R ISCO is fitted to be around 5 R g over the decay 33 . Spectral fits to NICER data of MAXI J1820+070 in the soft state suggests a low spin a < 0.5 (ref 52 ). In order to study the effect of the low spin on the evolution of the corona, we refit six observations which cover the decay of this outburst with the same model, but fixing the spin a = 0.5. It turns that the evolution of the corona remains (see Supplementary Figure 8a), i.e., the reflection fraction is also smaller than unity, and decreases with the decay, which requires the relativistic outflowing motion of the corona.
The inclination is fixed to the constant value θ = 63 • based on the jet/optical measurements. However, the inner disc/jet/orbital plane do not necessarily align. In order to study the constraint of spectral fits on the inclination and its effect on the evolution of the corona, we refit six observations which covers the decay of this outburst with the same model, but allowing the inclination to be free. It turns that, the spectral fits with the standard relativistic reflection model relxillCp also prefers to large values of the inclination, and the evolution of the corona remains (see Supplementary Figure 8b and 8c), i.e., the reflection fraction is also small than unity and decreases during the decay phase.

Justification of the outflowing corona
Recently, the outflowing velocity of the lamppost X-ray primary source is taken into account in the package relxilllpionCp of relxill(v1.4.0). In principle, we could refit the Insight-HXMT data to directly infer the outflowing velocity and its evolution over time. We found that the outflowing velocity and the height cannot be uniquely determined in the spectral fits. Therefore, deriving the evolution of the outflowing velocity may require that the height of the corona can be estimated from other measurements. e.,g., timing analysis. However, this is difficult to achieve due to the following two facts: (1) the effective area of LE detector in Insight-HXMT is one order of magnitude smaller than that of NICER, so that the derived time-lags in low energy band is not precise enough; (2) converting the lag into a light travel time distance between the corona and the accretion disk is not straightforward. A number of effects, e.g., the geometry of the system, the viewing angle to this source and relativistic Shapiro delay, etc, have to be well studied

14/45
before the reliable measurement of the light travel time can be derived 28 . Nonetheless, we could fix the height of the lamppost X-ray source at constant value, e.g., H = 7R g , and refit the spectra during the decay to directly infer the outflowing velocity and its evolution over time. It is shown in Supplementary Figure   7a that the outflowing velocity indeed increases from ∼ 0.0 c to ∼ 0.8 c along with time, resulting in the decrease of the reflection fraction.
Moreover, we simulate the energy spectra for NuSTAR observation, as a function of the lamppost height and outflowing velocity, using the reflection model relxilllpionCp which takes the outflowing of the corona into account. The EW of these simulated NuSTAR spectra is then estimated, which turns out to be fairly constant (see Supplementary Figure 9). Given the simulation results in Supplementary Figure   4a and 9, it can be seen that the EW of the iron line depends on not only the height and the outflowing velocity of the lamppost but also the iron abundance and ionization of the reflection disk 53 .

Comparison to NuSTAR spectral analysis
NuSTAR spectra of MAXI J1820+070 in the decay were reported and analysed in ref. 33 . The two lamppost point sources (relxilllpCp) with different heights, as the reflection model, were used in their spectral fits. It was shown that, during the decay, the height of the lower point source remains constant at about ∼ 4R g , while the height of the upper point source decreases from 100 R g to a few R g . We also used the two lampposts model to fit Insight-HXMT spectra, and found that the derived height of the upper lamppost decreases with the decay as well, which is consistent with the results of NuSTAR. Note that, in relxilllp(Cp) and relxilllpion(Cp), the disk outer radius cannot exceed 1000 R g . Therefore, when the corona is high, e.g., H > 100R g , a fraction of photons will not hit the disk, which results in the reflection fraction being lower than unity. This should be distinguished from the case of the corona being relativistically outflowing at low height H < 100R g . In the latter case, the reflection fraction will also be small with R f < 1, but now it is due to the beaming effect of the relativistically outflowing motion.
However, in the lamppost model relxilllpCp, the corona is assumed to be stationary. The conversion from height to reflection fraction is not self-consistent between the stationary corona model (relxilllpCp) and the moving corona model (relxilllpionCp), since the corona discussed in this work might be moving at different heights, along with the evolution of the outburst. Taking the spectrum of the later epoch MJD

15/45
=58271 (obsID = P0114661061) as an example, when the corona was believed to be contracting closer to black hole 28 , we found that, both the stationary corona with the height H ∼ 3R g using relxilllpCp, and the outflowing corona at different heights (e.g., H ∼ 17R g , v ∼ 0.5c; H ∼ 45R g , v ∼ 0.4c) using relxilllpionCp, can provide the equivalently good fits to the data. In the former case, the reflection fraction is high R f = 3.3, while in the later case, the reflection fraction is low R f ∼ 0.6. These comparison studies above confirm that the decrease of the height alone cannot explain the observed reflection evolution, which also requires the outflowing corona reported in this work.
Meanwhile, in the case of an outflowing corona (relxilllpionCp), the height and the outflowing velocity cannot be uniquely determined from the spectral fits, as discussed above. Nonetheless, at later epoch of the outburst (e.g., MJD =58271), the corresponding reflection fractions are roughly the same R f ∼ 0.6, and the fitted velocity is moderate with v ∼ 0.5c. In contrast, at the peak of the outburst, e.g., MJD =58204 (obsID = P0114661006), the fitted reflection fraction is R f = 1.2, and the corona prefers to be stationary, using relxilllpionCp model. Therefore, although it is not easy to uniquely determine the height and the outflowing velocity, the resultant reflection fraction decreases during the decay, which is consistent with our results in the spectral fits of Insight-HXMT using relxillCp.
Here, we also directly fit the reflection fraction from the spectral fits to NuSTAR data with the relxillCp + xillverCp model, as is done for Insight-HXMT in this work. The inclination angle is fixed at constant large value θ = 63 • . The two diskbb components are used to account for the differences between FPMs due to a thermal blanket tear in FPMA 54 . It turns out that the the reflection model which consists of the relativistic reflection relxillCp and non-relativistic reflection xillverCp can also fit the NuSTAR spectra well with reduced χ 2 ≤ 1.1, and the best-fitting values of the parameters are given in Supplementary Tables 1, 2 and 3. The index of emissivity profile q 1 are pegged at zero. This also suggests that the corona is spatially extended, as seen in Insight-HXMT fits. More importantly, the reflection fraction turns out to decrease with the decay from 0.8 to 0.3 (see Supplementary Figure 10), which has the same trends with the Insight-HXMT results in this work (see Fig. 3).

Data reduction
Insight-HXMT consists of three groups of instruments: High Energy X-ray Telescope (HE, 20-250 keV), Medium Energy X-ray Telescope (ME, 5-30 keV), and Low Energy X-ray Telescope (LE, 1-15 keV). The systematic errors for ME 20-30 keV are relatively smaller. Fluorescence lines due to the photoelectric effect of electrons in Silver K-shell are detected by the Si-PIN detectors of ME, which dominates the spectrum over 21-24 keV. Therefore, the spectrum over 21-24 keV is ignored; The systematic errors of HE, compared to the model of the Crab nebular, are less than 2% at 28-120 keV and 2%-10% above 120 keV. According to these calibration results, the energy bands for the spectral fits in this work are 2-10 keV for LE, 10-30 keV for ME, and 28-200 keV for HE, with a systematic error of 1.5% for LE/ME/HE.
We use the Insight-HXMT Data Analysis software (HXMTDAS) v2.0 to analyze all the data, filtering the data with the following criteria: (1) pointing offset angle < 0. corresponding background spectra and response files are also combined together.

Qualitative spectral analysis
The spectra are analyzed using XSPEC version 12.10.0 56 . In Fig. 2

Configuration of the spectral fitting model
In this work, we applied the relativistic model relxillCp which accounts for the reflection with the broad/high-ionization iron line and the non-relativistic model xillverCp which accounts for the reflection with the narrow/low-ionization iron line, to fit Insight-HXMT spectra. For both models, we use the same nthcomp parameters as the incident spectrum inputs, and the seed photon temperature is fixed at 32 kT bb = 0.05 keV. The inner radius, which is derived by fitting the NICER spectra with the relativistic reflection model, is estimated to be less than 2 R g under some assumptions 28 . It was also suggested that there is little or no evolution in the truncation radius of the inner disk during the luminous hard state, given that the thermal reverberation lags remain throughout all epochs and that the spectral shape of the broad Fe

18/45
line component is constant over time. In ref 33 , the inner radius is fitted to be around 5 R g over the decay. If assuming the black hole spin a = 0.998, this suggests that the inner disk may be truncated at small radius.
Instead, if the accretion disk extends to the innermost stable circular orbit (ISCO), i.e., R in = R ISCO 5R g , where R ISCO is the innermost stable circular orbit, this may imply that the black hole spin is low. Also, fits to the spectra of NICER in the soft state suggests a low spin 52 . Indeed, fixing R in , spin and inclination at different values, will result in different values for other free parameters, but the time-dependent trends in the fitted parameters will be preserved. Since in this work we are mainly concerned with the evolution of the corona, we assume the black hole spin a = 0.998, and the inner radius is assumed to be R in = R ISCO .
The effect of the low spin on the evolution of the corona will be discussed below. The outer radius of the reflection disk is fixed at the upper limit of their table model R out = 1000R g .
The strong absorption dips were detected by NICER in the outburst of MAXI J1820+070, suggesting that this source is a high-inclination system 24,37 . The observed sharp increase in the Hα emission line equivalent width and the absence of X-ray eclipses in MAXI J1820+070 in ref 57 indicated the inclination angle to be 69 • < i < 77 • . The measurement of the radio parallax indicates that this source is located at a distance of 2.96 ± 0.33 kpc away from us. Together with the measured proper motions of the approaching and receding jet ejecta in radio, the inclination angle between the jet and the line of sight is estimated to be 63 • (see ref 21 ). In this work, we take the inclination angle as 63 • . The free parameters of the relativistic reflection model relxillCp are the emissivity profile q 1 , the reflection fraction R f , the incident photon index Γ, the Fe abundance relative to the solar value A Fe , the ionization parameter logξ , and the electron temperature kT e . The emissivity profile is parameterised in terms of the index q 1 , q 2 and the broken radius R br . For simplicity, the index of the outer disk region is fixed with q 2 = 3.0. The index of the inner disk region q 1 and the broken radius R br are free parameters. As for the non-relativistic reflection model, xillverCp, the spectral fits are insensitive to the ionization parameter, which is evaluated with F-test in XSPEC. Therefore, the ionization parameter are fixed at low value with logξ = 1.0, while the only free parameter is the normalization.

Spectral fitting method
The statistical analysis of all the fits are achieved by implementing a Markov Chain Monte-Carlo (MCMC) algorithm. More specifically, we used the MCMC algorithm implemented in XSPEC to create a chain of parameter values whose density gives the probability distribution for that parameter. For each individual spectrum fitting, two chains of 250000 steps are run with 20 walkers. The first 10 4 steps of the running are discarded as burn-in phase. Convergence between these two chains is assessed by Gelman-Rubin diagnostic 58

Data availability statement
All Insight-HXMT data used in this work are publicly available and can be downloaded from the official website of Insight-HXMT: http://hxmtweb.ihep.ac.cn/

Code availability
The data reduction is done by the use of the software (HXMTDAS) v2.0 which is available at the     Figure 3. The correlations of the reflection fraction of the relxillCp, the iron abundance A Fe with respect to the solar abundance value, the constant factor of ME instrument and the constant factor of HE instrument, with respect to each other, which are investigated with the Spearman's rank test. The black points correspond to the median of the values and the error bars correspond to 68% confidence interval, which is calculated using the corner 59 package to analyse the probability distributions derived from the MCMC chains. The uncertainties of the fitted parameters arise from both the statistical and systematic uncertainties. The outflowing velocities of the lamppost X-ray source in the best-fitting of TBabs * (diskbb + relxilllpionCp + xillverCp) * constant model. The new model relxilllpionCp (attributing to the reflection with the broad iron line) includes a velocity of the lamppost X-ray source which is set to be free here, and the height of the lamppost X-ray source is fixed at H = 7R g . Assuming an inclination angle θ = 63 • , the inner/outer radius of the reflection disk R in = R ISCO (R ISCO is innermost stable circular orbit) and R out = 1000 R g , the black hole spin a = 0.998. The black points correspond to the median of the values and the error bars correspond to 68% confidence interval, which is calculated using the corner 59 package to analyse the probability distributions derived from the MCMC chains. The uncertainties of the fitted parameters arise from both the statistical and systematic uncertainties. The reflection fraction for six epochs in the best-fitting of TBabs * (diskbb + relxillCp + xillverCp) * constant model, but the black hole spin a = 0.5. Assuming an inclination angle θ = 63 • . The reflection fraction (b) and the inclination angle (c) for six epochs in the best-fitting of TBabs * (diskbb + relxillCp + xillverCp) * constant model, where the inclination angle (degree) is a free parameter, but the spin is fixed at a = 0.998. Assuming the inner/outer radius of the reflection disk R in = R ISCO (R ISCO is innermost stable circular orbit) and R out = 1000 R g . The black points correspond to the median of the values and the error bars correspond to 68% confidence interval, which is calculated using the corner 59 package to analyse the probability distributions derived from the MCMC chains. The uncertainties of the fitted parameters arise from both the statistical and systematic uncertainties. The observation ID of six epochs are listed as follows: ObsID = P0114661006, P0114661017, P0114661024, P0114661032, P0114661043, P0114661055. Supplementary Figure 10. The evolution of the reflection fractions of the relativistic reflection model relxillCp, which are estimated from the fits of NuSTAR spectrum. The numbers at lower right corner correspond to the obsID in Table 1 of ref. 33 . Note that, the first point (i.e., obsID = 02) is at the rise of MAXI 1820+070. For completeness, we still add it together with the rest points which roughly correspond to the decay of the outburst (see Fig. 1 of ref 33 ). The data points correspond to the median of the values and the error bars correspond to 68% confidence interval, which is calculated using the corner 59 package to analyse the probability distributions derived from the MCMC chains. The uncertainties of the fitted parameters arise from both the statistical and systematic uncertainties. Supplementary Figure 11. An illustration of one and two dimensional projections of the posterior probability distributions derived from the MCMC analysis for the parameters in relxillCp, i.e., the emissivity profile q 1 , the photon index Γ, the ionization parameter logξ (in units of erg cm s −1 ), the abundance A Fe with respect to the solar value, the electron temeperature kT e (in units of keV), the reflection fraction R f , the normalization, the constant for ME and the constant for HE. The contours in the two dimensional projections (the second to the bottom panel of each column) for each two parameters correspond to 1-, 2-and 3-σ confidence interval. The top nine panels (a1-i1) in each column are their one dimensional projections on the corresponding x axis. The values above each one dimensional projection indicate the median value of each parameter, as well as the upper and lower limits of the 68% confidence intervals. The vertical lines in the one dimensional projections correspond to the lower, median and upper value of each parameters. The MCMC analysis and the resultant figures above are produced using the corner package 59 . This illustration corresponds to the spectral fitting (see Fig. 1) of MJD 58204 (ObsID = P0114661006).