• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Application of a Space-based Optical Interferometer Toward Measuring Cosmological Distances of Quasars

    2022-05-24 14:20:58YingKeHuangYueDongFangKaiXingLuZhiXiangZhangJiLinLiuShashaLiBaoRuiLuoQinLinandZhuoXiHuo

    Ying-Ke Huang,Yue-Dong Fang,Kai-Xing Lu,Zhi-Xiang Zhang,Ji-Lin Liu,Sha-sha Li,Bao-Rui Luo,Qin Lin,and Zhuo-Xi Huo

    1Qian Xuesen Laboratory of Space Technology,China Academy of Space Technology,Beijing 100081,China;huangyingke@qxslab.cn,huozhuoxi@qxslab.cn2 Yunnan Observatories,Chinese Academy of Sciences,Kunming 650011,China3 Department of Astronomy,Xiamen University,Xiamen,Fujian 361005,China

    Abstract Measuring quasar distance through joint analysis of spectroastrometry and reverberation mapping observations is a new method for driving the development of cosmology.In this paper,we carry out detailed simulation and analysis to study the effect of four basic observational parameters(baseline length,exposure time,equivalent diameter and spectral resolution) on the data quality of differential phase curves (DPCs),and furthermore on the accuracy of distance measurement.In our simulation,we adopt an axisymmetrical disk model of a broad line region(BLR)to generate differential phase signals.We find that the differential phases and their Poisson errors could be amplified by extending the baseline,while the influence of optical path difference errors can be reduced during fitting the BLR model.Longer exposure time or larger equivalent diameter helps reduce the absolute Poisson error.Therefore,the relative error of DPCs could be reduced by increasing any of the above three parameters,then the accuracy of distance measurement could be improved.In contrast,the uncertainty of absolute angular distances(DA) could be improved with higher spectral resolution,although the relative error of DPCs would be amplified.We show how the uncertainty of distance measurement varies with the relative error of DPCs.For our specific set of model parameters,without considering more complicated structures and kinematics of BLRs in our simulation,it is found that the relative error of DPCs <20%is a limit for accurate distance measurement.The relative error of DPCs has a lower limit (roughly 5%) and the uncertainty in distance measurement can be better than 2%.

    Key words:(galaxies:) quasars:emission lines–techniques:interferometric–(cosmology:) distance scale

    1.Introduction

    One of the basic open questions of modern astrophysics is to accurately measure the cosmological distances of extragalactic objects to understand the increasing H0tension(Peacock 1999;Freedman &Madore 2010;Weinberg et al.2013;Freedman 2017;Riess et al.2019).Active galactic nuclei (AGNs),which are known as the brightest objects in the universe,have been utilized to measure cosmological distances since they were discovered (Sandage 1965;Hoyle 1966;Longair &Scheuer 1967;Baldwin 1977;Elvis &Karovska 2002;Quercellini et al.2009;Wang et al.2013;Marziani &Sulentic 2014).Nevertheless,due to lack of understanding of AGN physics,many attempts have been made but none of them were proven to be bias free.Recently,a joint analysis (Wang et al.2020) of spectroastrometry (SA:Petrov 1989;Bailey 1998;Gravity Collaboration et al.2017) and reverberation mapping (RM:Blandford &McKee 1982;Peterson 1993;Kaspi et al.2000;Bentz et al.2013;Du et al.2018)observations provided a direct method to measure absolute angular distances (DA) of AGNs.This method does not need calibration using cosmic ladders so it has the potential to provide absolute distance measurements for quasars in the high-z universe.

    The high-energy output of an AGN is argued to be the result of accreting matter by a supermassive black hole at the center of a galaxy (Lynden-Bell &Rees 1971).A hallmark of AGN spectra is the existence of broad emission lines (Osterbrock &Mathews 1986;Osterbrock 1989;Ho 2008).These broad emission lines are thought to originate from photoionization of the broad line region (BLR) clouds by high-energy continuum photons from the central accretion disk.The revolution of these clouds is dominated by Kepler motion around the central black hole due to gravity,which has been confirmed by multiple campaigns of several AGNs (Peterson et al.2004).The emission lines are broadened up to several thousand km s-1(Antonucci 1993;Urry &Padovani 1995).However,it is impractical to spatially resolve the BLR by direct photometry,because the angular size of the BLR is smaller than the existing spatial resolution limit (<0.1 mas;Blandford &McKee 1982;Elvis &Karovska 2002).

    In the past few decades,a promising approach to directly constrain BLR geometry is the widely used RM technique.The underlying principle of the RM technique is photoionization theory (Wandel et al.1999;Kaspi et al.2000;Bentz et al.2006).By analyzing the variation of continuum and emission lines after measuring the time delay between the variable continuum and the responsive broad emission line,it is possible to constrain the geometry and kinematics of the BLR and further estimate the characteristic size of the BLR as well as the black hole mass.In 2018,Gravity Collaboration et al.(2018)proposed another method.By measuring the photocenters of emission lines in different wavelength channels,the SA method(Beckers 1982;Petrov 1989;Petrov et al.1992;Gnerucci et al.2010)proved to have great advantage in providing information on the spatial structure at scale much smaller than the spatial resolution of an interferometer.An intuitive picture of how the SA method can spatially resolve the bulk motion of the BLR clouds is illustrated as follows:there are two clouds located at a distance significantly smaller than the spatial resolution,so they are seen as a non-resolved source with a global angular size,but one cloud is moving toward the observers and the other moving away from observers.As the line-of-sight (LOS)velocities of the two clouds are different,the photons from the same emission line but from different clouds are shifted to different wavelengths.By measuring the interferometric phases(providing position information on-sky) over wavelengths,we can obtain the photocenter position of clouds with different velocities and further separate the motion of BLR clouds.The SA method demonstrates the capability of overcoming the spatial resolution limit,thus being able to constrain the angular scale,geometric structure and dynamics of the BLR.

    In fact,RM observations provide the linear size of BLR and are more sensitive to the direction along the LOS.SA measures the angular size of BLR and is more sensitive to a direction perpendicular to the LOS (Wang et al.2020).Relying on the geometric and dynamical BLR model,we can extract information from RM and SA data to measure the cosmological distance.However,up to now,RM campaigns measure the region of optical emission line(mainly using Hβ)while the SA campaigns measure the near-infrared (NIR) emission lines(using Paα in Gravity Collaboration et al.2018 and Wang et al.2020).Since different emission lines come from different areas of the BLR (Clavel et al.1991;Dietrich &Kollatschny 1995;Kollatschny 2003),if RM and SA rely on the same emission line,the systematic error (mentioned in Wang et al.2020) of distance measurement in joint analysis could be significantly reduced.It is necessary to carry out SA observations in optical bands or RM observations in infrared bands.Considering the extinction and jitter of the Earth’s atmosphere,the interference process in space shows advantages in long exposures.A spacebased optical interferometer will be the future direction.

    As a necessary step,we aim at studying the uncertainty of distance measurement caused by basic observational parameters of a space-based optical interferometer in this paper.The structure of this paper is as follows.In Section 2,we describe the geometric and dynamical model of the BLR adopted.In Section 3 we simulate the expected spectroastrometric signals and describe the analysis method in Section 4.Results are provided in Section 5,and we show our summary and discussions in the last section.

    2.Parameterized BLR Model

    In the optical band,an AGN spectrum consists of continuum emission and line emission.The continuum emission is considered to originate chiefly from the accretion disk.For nearby galaxies,continuum emission from the host galaxy could also not been ignored.For simplicity,we assume the region that produced continuum emission has an axisymmetrical structure and the photocenter of continuum emission is located at the black hole.The broad emission lines are considered be generated from photoionization of BLR clouds by the ultraviolet (UV) continuum photons from the accretion disk.The emission line used in our simulation is Hβ,which is a strong optical emission line and widely used for RM observation to estimate the linear size of the BLR.The continuum underneath emission lines is often assumed to have a linear form.This hypothesis is widely utilized in RM to calculate the integral flux of Hβ,where a continuous spectrum is required to be subtracted off (Du et al.2014;Huang et al.2019).In order to estimate the photocenter position of a broad Hβ emission line and simulate the SA observations,we rely on a simple BLR model that is characterized by a flat disk with circular Kepler rotation.

    In recent years,thanks to many studies combining RM observation and the parameterized BLR model,we have gained a better understanding of BLR physics (Li et al.2013,2018;Pancoast et al.2014a,2014b;Grier et al.2017;Williams et al.2018).In this section we give an overview of the parameterized model we used for BLR.It consists of many isolated BLR clouds,which are modeled as a large number of non-interacting point particles.These particles reprocess the continuum photons originating from the accretion disk into the emission line photons instantaneously.The wavelength of the emission line photons is determined by the particles’ velocity and the time lag for the response is determined by their position.For simplicity,the accretion disk is regarded as a point-like geometry so that the UV ionizing continuum is isotropic and the flux of continuum falls off with the square of the distance.

    2.1.Geometry

    In our simulation,the radial distribution of the point particles is described as a shifted Γ-distribution.The distance r of a point particle from the black hole is given by

    Here Γ0=p(x|1/β2,1) is drawn randomly from a Gamma distribution

    and the point particle is then shifted radially by a Schwarzschild radius Rs=2GMBH/c2,plus a minimum radiusRminof BLR.Here RBLRis the mean radius andFratio of the minimum radius to the mean radius.WithinRmin,the line-emitting particles are not allowed to exist.β is the shape parameter of the radial distribution:small values of β mean narrow normal distributions while large values mean exponential distributions.We then define a half-opening angle θofor the overall geometry and the spatial distribution of point particles.Here θo=0°defines a thin disk(ring)while θo=90°represents a spherical distribution.An observer views the BLR from an inclination angle θi,where θi=0° corresponds to a face-on orientation.

    For each point particle,the emission is weighted by a parameter ω(φ),

    where k is a free parameter between -0.5 and 0.5.When k=0.5,the inner side of BLR is contributing more line emission.Here φ is the angle between the observer’s LOS to the black hole and the point particle’s line to the black hole.The particles could be clustered and θ is the angle of a point particle from the disk

    where U is a random number drawn uniformly between 0 and 1,and γ is a free parameter drawn uniformly between 1 and 5.When γ=1,the point particles are clustered toward the disk.

    We also adopt a usual assumption that BLR clouds have the same size and no shadowing among them (Li et al.2013).We use the traditional linear response of the emission lines to the continuum.For a BLR that has a given geometry and cloud distribution,the emission-line flux flat time t is estimated by summing over the emissions from all the clouds

    Here ωiis the weight of the ith cloud calculated by Equation (3),A is the response coefficient anddescribes the ionizing flux received by the ith cloud at time t.

    2.2.Dynamics

    We assume a flattened disk-like BLR with Keplerian rotation.Due to the gravity of the central black hole,the point particles rotate in a circular orbit.The Keplerian velocities VKof the point particles are drawn from a distribution that dependsupon the black hole mass MBHand the distance r,which could be estimated as

    For a thick Keplerian disk with a half opening angle θo,there will be an angle θ between the orbital plane and the equatorial plane.Thus,the LOS velocities vlineare related to the Keplerian velocities,half-opening angle and the inclination angle.For cloud particles with an emitted wavelength λemit,considering the relativistic effect,Doppler broadening and gravitational redshift which will affect the emission line profile,the observed wavelength λobscan be written as

    In order to explore the influence of basic observational parameters on the uncertainty of distance measurement,we adopt the same set of BLR model parameters as shown in Wang et al.(2020),which has been summarized in Table 1.

    Table 1 Parameters Used in the BLR Model

    3.The Spectroastrometric Signal of the BLR

    The SA method(see Bailey 1998;Rakshit et al.2015)shows advantages in providing information on the spatial structure of the object,especially when the global angular size Λ of a nonresolved source is smaller than the interferometer limit λ/B,here λ is the observed wavelength and B is the length of the interferometer baseline.By measuring the photocenters of different wavelength channels,this method has been applied to observe AGNs (Petrov et al.2001;Marconi et al.2003;Gnerucci et al.2010,2011;Rakshit et al.2015) in NIR and successfully constrained the size of the BLR for 3C 273(Gravity Collaboration et al.2018;Wang et al.2020).In this section we give a description of how we simulate the SA signal.

    For a given object,the observed interferometric phase Φ at wavelength λ can be written as

    where Φ*(λ) is the phase contributed from the real astrophysical signal and Φt(λ) is the phase contributed from the atmosphere,telescopes,interferometric delay lines and the various instrumental parts down to the detector (see Vannier et al.2006).Petrov (1989) shows that,for an interferometer with a baseline B,the interferometric phase for the object is

    Here∈(λ)is the photocenter of the source at wavelength λ and could be written as

    O(r,λ) is the brightness distribution of the source.If defining the fraction of the flux of the emission line to total as fl(λ)=Fl(λ)/Ftot(λ),where Fl(λ) is the flux of the emission line at wavelength λ and Ftot(λ) is total flux at wavelength λ,we have∈(λ)=fl(λ)∈l(λ),r is the vector to the central black hole and

    We use the bold letters to signify the vectors B,∈(λ) and r.

    For space-based interferometers,the phase contributed from the atmosphere could be ignored but the phase contributed from the instrument could not.In principle,if we add a spatial filter before the beam recombination,the optical effects left after this spatial filtering only include differences in intensity between beams and an optical path difference (OPD).The previous one can be calibrated by photometric channels,so only the phase shift caused by the OPD is taken into consideration in our simulation.The OPD has two origins:one originates exclusively from piston and chromatic dispersion along the path before the filtering,while another originates from wave front corrugations induced by imperfect adaptive optics (Tubbs 2005;Vannier et al.2006).For a given spectral channel,the phase shift caused by OPD is

    Hence,the differential phase between a channel λiand the reference channel λrwill be

    By inserting Equations (9),(10) and (12) into Equation (13),we have

    If the BLR model is specified,the surface brightness distribution of the regions is mainly related to the total number of collected photons,which is linearly dependent on the number of telescopes (nT),the overall quantum efficiency (η),the collecting area of telescopes S and the exposure time t.For each wavelength channel,the observed photocenter is further determined by the choice of spectral resolution R.If we also know the configuration of the projected baseline B,we could predict the actual phase signal by Equation (9).By giving the limit OPD (rms) control errors during an observation,using Equation (14),we could obtain wavelength-dependent differential phase after DPCs.For simplicity,the equivalent diameter(D) is used to estimate the number of photons per unit time in our simulation,S=π(D/2)2.

    In our simulation,uncertainty of the signal consists of two parts,one part comes from Poisson noise and the other is caused by the uncertainty in OPD control.Each set of scientific exposures is assumed to be the average over nobs(the number of observations) frames of t-s integration (NDIT=nobsand DIT=t s).We use the average value of nobsobservations at wavelength λ as the differential phase ΔΦ(λ) and the dispersion of nobsobservations as the error.We show the dependence of mock DPCs on eight parameters used in the BLR model in Figure 1.

    Figure 1.Dependence of DPCs on eight parameters.The parameters for the interferometer are fixed at the values listed in Table 2.When varying one parameter,other parameters are fixed at the values listed in Table 1.The parameters and their values are shown in the legend.The units of θo,θi and PA are degrees.The units for RBLR,MBH and DA are lt-day,108 MSun and Mpc,respectively.

    4.Analysis

    In this section,we describe the analytical method of obtaining the BLR model parameters from mock data.The joint analysis proposed by Wang et al.(2020) uses data from three parts:light curves of optical continuum and emission lines,DPCs and the profiles of the emission lines.The light curves and emission profiles depend on optical spectral observation,which has been a very mature technology,so we assume that the quality of the light curves and the profiles of the emission lines are good enough.Thus,the uncertainty in the joint analysis depends mainly on the data quality of DPCs.

    As mentioned above,we generate the mock DPCs from Equation (14).We use the DPCs to obtain the posterior distribution of the BLR model parameters.Suppose that probability distributions for the measurement errors of the DPCs are Gaussian and uncorrelated,then the likelihood function can be written as

    HereD represents the measured data,Θ signifies the BLR model parameters,Φobsis the interferometric phase of the emission line with the uncertaintiesσφij,Φmod(Θ)is the corresponding predicted values from the BLR model,Nobsis the number of observations and Nλis the corresponding number of wavelength bins.According to Bayes’ theorem,if we know the prior distribution of the model parameter,the posterior probability distribution for Θ should be given by

    whereP(D) is a normalization factor.We list the parameters of the BLR model and the prior ranges in Table 1.

    The EMCEE package is used to sample the parameters efficiently (Foreman-Mackey et al.2013),which is an MIT licensed pure-Python implementation of Goodman &Weare's Affine Invariant Markov chain Monte Carlo (MCMC)Ensemble sampler.The sampling is run with 1000 “walkers”independently and converges within about 2000 trials.

    5.Results

    Utilizing the BLR model described in Section 2,we generate mock observations according to a set of varying basic observational parameters (B,D,t,R) described in Section 3.Combined with the method mentioned in Section 4,we first study the effect of inclination angle and the projected angle on the photocenter shift,then we study the influence of the basic observational parameters on the data quality of DPCs and the accuracy of distance measurement.

    5.1.The Effect of Inclination Angle on The Photocenter Shift

    The SA method can spatially resolve the velocity gradient perpendicular to the LOS direction,the size of which varies with the inclination angle.Therefore,we first study how inclination angle affects the distribution of photocenters.We used the hypothesized BLR model mentioned above.The parameters of the model are fixed and listed in Table 1,except for the inclination angle.We vary the inclination angle (θi) from 0° to 180°.Then we estimate the corresponding distribution of photocenters which are shown in the left panel of Figure 2.By changing the inclination angle,we find that for a BLR with Keplerian motion,the more“edge-on”the LOS is,the larger the spatial distribution size of wavelength-dependent photocenters.When viewed in a direction perpendicular to the disk plane(faceon:θi=0°),the photocenters of different wavelengths could not be distinguished(as shown in the lower right panel of Figure 2).This conclusion can be understood as when viewed from angle θi,the velocity components along the LOS direction are proportional to sin θi.Therefore,when viewed from the edge-on direction(θi=90°),the velocity along the LOS direction achieved a maximum such that the observed wavelength gradient is the largest(depicted in the upper right panel and lower right panel of Figure 2).We must emphasize that shielding of the dust torus is not taken into account in our simulation.

    5.2.The Influence of Projected Angle on the Photocenter Shift

    In fact,the angle between the baseline and the plane of BLR will also affect the observed distribution of photocenters.This angle is defined as the projected angle(displayed in Figure 3(b)as θB)in our simulation.Colored points represent the photocenters of different wavelengths,and the black solid line signifies the baseline direction.Each baseline can only measure the distribution of photocenters projected to the baseline direction,which could be indicated from Equation (9).To study the impact of the angle θB,we use the same BLR model as mentioned above and fix θi=90°,then vary the angle θBto measure the displacement of photocenters at different wavelengths.Results are displayed in Figure 3(a) and the relation between the maximum offset of photocenters and θBis plotted in Figure 3(c).We see that the maximum offset of photocenters is proportional to the absolute value of cos θB.When the baseline is perpendicular to the plane,the photocenters of different wavelengths along the baseline direction are indistinguishable.

    In actual observations,the direction of the BLR plane(so the projected angle θB)is unknown to the observer.To ensure that the mock data contain information about spatial position,we adopted an observation strategy of rotating the baseline.In our simulation,we rotate the baseline regularly with a cadence of 10°.Thus we get 19 mock data for each set of observational parameters.

    5.3.The Influence of Basic Observational Parameters on Data Quality

    Furthermore,we study the effects of the basic observational parameters by varying one parameter at a time (the ranges of variations are summarized in Table 2).In addition,the projected angle θBwas fixed at θB=0°.

    Table 2 Parameters Used in the Interferometry Model

    In order to quantify the impact of basic observational parameters mentioned above,we define two metrics to describe the data quality:the largest amplitude and the relative error of DPCs.The largest amplitude of DPCs is defined as the maximum absolute amplitude of a set of DPCs.The relative error of DPCs is quantified as the ratio of phase error to phase value where phase is the largest amplitude of DPCs.We use the relative Poisson error of DPCs,which is defined as the ratio of the Poisson error to the largest amplitude of DPCs,to quality the impact of Poisson noise.The results are depicted in Figure 4.

    The impact of baseline length on data quality is shown in Figure 4(a).We feature three DPCs generated under three different baseline lengths in the left panel.These DPCs display an obvious S-shape.The value of the differential phase increases as the baseline gets longer,as does the absolute error.We find that the largest amplitude of DPCs increases linearly with the baseline length (shown in the upper right panel),which has been implied from Equation (9).As the baseline gets longer,the relative error tends to be smaller(shown in the bottom right panel).We also find that the Poisson error becomes dominating while the baseline lengthens.

    The impact of the exposure time is shown in Figure 4(b).We also display three DPCs generated under three different exposure times in the left panel.It is easy to see that as the exposure time increases,the absolute error of DPCs decreases significantly,but the value of the differential phase remains essentially the same.As featured in the upper right panel,the largest amplitude of DPCs does not change with exposure time.However,relative error and relative Poisson error of the DPCs decrease as the exposure time increases.

    Figure 2.The effect of inclination angle on the photocenter shift.Left:The calculated distribution of photocenters that varies with the inclination angle.Upper right:The position distribution of photocenters at different wavelengths as seen by an observer from the edge-on direction(θi=90°).Lower right:The position distribution of photocenters at different wavelengths as seen by an observer from the face-on direction (θi=0°).

    The influence of the equivalent diameter is basically the same as that of exposure time,and the results are depicted in Figure 4(c).As the equivalent diameter becomes larger,absolute error,relative error and relative Poisson error of the DPCs are significantly reduced,but the largest amplitude of DPCs remains the same.For our simulation,we find that by changing the equivalent diameter,relative error of the DPCs can be reduced to less than 5%,which is more effective than increasing the exposure time to roughly 10%.

    The results of how the spectral resolution impacts the DPCs are shown in Figure 4(d).As the resolution increases,the number of data points on the DPC increases significantly(shown in the left panel),while the largest amplitude of DPCs increases slightly(as shown in the upper right panel) and the relative error of DPCs becomes larger (as displayed in the bottom right panel).This conclusion can be understood as the larger the spectral resolution,the more wavelength channels there are,but the less number of photons in each wavelength channel.

    Figure 3.The influence of projected angle on the photocenter shift.Left:Simulated position of photocenters at different wavelengths under different projected angle θB.Upper right:The black solid line shows the baseline direction,colored points represent the photocenters of different wavelengths and θB is the projected angle.Lower right:The maximum offset of photocenters varies with θB.

    5.4.The Effect of Basic Observational Parameters on the Accuracy of Distance Measurement

    In order to quantify the accuracy of distance measurement,we define the median value of the posterior distribution as the best value for distance inferred from MCMC,with the 16%and 84%quantiles as the lower and upper bounds respectively.The relative uncertainty of the distance measurement is defined as the ratio of the error (mean of the lower and upper bounds) to the best value mentioned above.The relative bias is defined as the fraction of the bias (difference between the best value and the input value) to the input value.We simulate the actual DPCs using the BLR model mentioned in Section 2 and the basic observational parameters described in Section 3.Then we apply the MCMC method described in Section 4 to obtain probability density distributions of BLR model parameters.In the Appendix,we list the simulated DPCs (Figure A1)and the mock Hβ line profile (Figure A2) generated by using the instrument parameters listed in Table 2,as well as the probability density distributions (Figure A3) of model parameters.We obtain the relative uncertainty (blue points)and bias (orange points) of distance measurements (DA) at different baseline lengths,exposure times,equivalent diameters and spectral resolutions.The results are depicted in Figure 5.

    Figure 4.ComparisonofDPCsanddataqualityunderdifferentparameters.Whenvaryingoneparameter,otherparametersarefixedatthevalueslistedinTable2.((a)baselinelength B;(b)total integration time t;(c)equivalent diameter D;(d)spectral resolution R).Left column:Comparison of simulated DPCs with error under different values of the parameter.Right:dataquality(DQ)forDPCsdifferingonlybytheparameter(upper:amplitudesofthecurves;lower:therelativeerrorofDPCs,bluefortotalerrorandorangeforPoissonerror).

    Figure 5.The relative uncertainty(blue points)and bias(orange points)of distance measurement change with each parameter.Upper left panel is for baseline length.Upper right panel is for exposure time.Lower left panel is for spectral resolution and lower right panel is for equivalent diameter.The large point is the result corresponding to the typical value listed in Table 2.

    We see that as the parameter values increase,the relative uncertainty of distance measurement (DA) can be reduced to 2%.There is an optimum value for each parameter.Above this optimum value,the accuracy of DAwill improve slowly with the parameter value.Extending the baseline can proportionally amplify the difference phases and their Poisson errors,which means that the influence of OPD errors can be reduced during the fitting process.Therefore,extending the baseline can improve the accuracy of distance measurement.Either increasing the exposure time or the equivalent diameter results in improving the accuracy of distance measurement,since the more photons that are collected,the smaller the absolute Poisson error.Interestingly,the larger the resolution is,the larger the relative error of DPCs,but the smaller the uncertainty of DAwill be.This is because as the resolution increases,the number of photons in each wavelength interval decreases and thus the Poisson error of the difference phase increases,meanwhile the amount of data on the DPC increases.As expected,the relative bias of DAbecomes smaller as the parameter becomes larger,and is even more sensitive than the relative error.

    5.5.The Effect of Data Quality on the Accuracy of Distance Measurement

    We also study how the uncertainty of distance measurement depends on the data quality of DPCs.Figure 6 visualizes the relationship resulting from our simulation.We see an obvious turnaround near the data quality around 20%.Only when the data quality is better than this turnaround point is the accuracy of DAhighly sensitive to the data quality.Different colors in Figure 6 represent different observational parameters.The direction of the arrow between data points represents parameters becoming larger and larger.

    Figure 6.The relationship between the relative error of DPCs and the accuracy of distance measurement.Left:The relative uncertainty of distance measurement and the relative error of DPCs under different parameters.Right:The relative bias of distance measurement and the relative error of DPCs under different parameters.The arrow direction is the direction in which the parameters increase.

    As the baseline increases (blue points),data quality gets better first,but there is no effect on the accuracy of distance measurement.When the relative error of DPC reaches better than 20%,extending the baseline can significantly improve the accuracy of DA,although the data quality remains basically the same.Similar results are found by increasing t and D.However,increasing the spectral resolution (green points)gives more complex results.As the spectral resolution increases,data quality slowly deteriorates,but the accuracy of DAis increasing.This may be due to the fact that the larger the spectral resolution,the more data points are on the DPCs,but the fewer photons there are in each wavelength interval.In actual observations,when DPCs with higher spectral resolution but lower data quality are obtained,it is customary to rebin the DPCs to reduce the resolution but increase the data quality.However,using DPCs that have not been rebinned is a better option for improving the accuracy of DA.

    6.Summary and Discussion

    In this paper,we have studied the application of the SA method using a space-based optical interferometer to measure cosmological distance of quasars,and we have discussed how the basic observational parameters affect the accuracy of distance measurements.Our main results are as follows:

    1.For the same target,the observed spatial distribution of the photocenters at different wavelengths will be affected by inclination angle θiand the projected angle θB.

    2.Four basic observational parameters will affect the data quality of DPCs.As the baseline gets longer,the value of differential phase increases,as does the absolute error.When the exposure time or equivalent diameter increases,the absolute error,the relative error and the relative Poisson error decrease significantly.However,the larger the spectral resolution is,the larger the amplitude and relative error of DPCs are.

    3.Four basic observational parameters will affect the accuracy of distance measurement.Extending the baseline will amplify the difference phases and the Poisson errors but will reduce the impact of OPD errors during the fitting process,so the uncertainty and bias of distance measurement both decrease.Increasing the exposure time or effective aperture can reduce the Poisson error due to the increase in the number of photons,thus improving the accuracy of distance measurement.With the increase of spectral resolution,the Poisson error increases due to the decrease of photon number in each wavelength channel,but the number of data on DPC increases,so although the relative error of DPC becomes larger,the uncertainty of distance measurement becomes smaller.

    By assuming a parameterized BLR model and specifying the value of parameters,we obtained the surface brightness distribution of the BLR and continuum regions,then calculated the photocenters at different wavelengths.Combined with the basic observational parameters,we conduct an extensive set of simulations to get simulated DPCs,and then we utilized the MCMC method to obtain the posterior probability distributions of BLR model parameters from the mock data.In our simulation,we can eventually limit the relative uncertainty of data to 10%–20%by extending the baseline.The exposure time eventually limits the relative uncertainty to 5%–15% while the equivalent diameter can be limited to 5%–10%.Moreover,the results show that the SA method makes sense for distance measurement only if the relative error of DPCs is less than 20%.In this case,further increases in baseline length,exposure time and equivalent diameter will not greatly reduce the relative error of DPCs,but the uncertainty of distance measurement can be sharply reduced.However,it should be emphasized that these thresholds of the relative error are only valid for our simulation,since they depend on the specific set of BLR model and simulated observations.The simplest model for BLR is used in our simulation.Systematic errors,which are caused by complicated structures and kinematics of the BLR,such as nondisk structure,non-Keplerian kinematics and the disordered motion of these clouds,are not taken into account.So,the relative error of DPCs would be underestimated during simulation then the accuracy of model parameters would be overestimated.On the other hand,we did not consider the uncertainty of RM data in the fitting process,which also resulted in overestimation of the accuracy of model parameters.

    It also needs to be emphasized that both RM and SA are required to measure DAin reality.The reason that we can rely on SA to determine RBLR(without using the light curves from RM)is because a specific BLR model is used in our simulation to simulate DPCs,and then the same BLR model is used to perform the DPC fitting to retrieve all the model parameters.This is also why the degeneracy among the eight parameters is so weak.However,in reality the geometric distribution,dynamics and reprocessing coefficient of BLR clouds are much more complex,so it is not possible to have a perfect BLR model for DPC fitting,hence it still requires RM to measure RBLR.

    Acknowledgments

    We are grateful to the referee for constructive suggestions that improved the manuscript.We acknowledge the financial support of the National Natural Science Foundation of China(Grants Nos.12003077,11703077,12073068)and the Yunnan Province Foundation (202001AT070069).

    Appendix Mock Data

    In this part,we use an example to illustrate the process of generating simulated data and the process of fitting.By adopting the parameterized BLR model described in Section 2 and the values of parameters listed in Table 1,we calculate the surface brightness distribution of BLR.We used 2×105(same as Gravity Collaboration et al.2018)clouds in the model.Combined with the basic observational parameters and the values listed in Table 2,we simulate the Hβ emission line and DPCs.The spectral resolution used is 1000,so there are 67 wavelength bins between 4750 and 4950?.For each wavelength bin,we obtain the number of photons for each cloud that belongs to this bin.Combined with the position of the cloud,we calculate the total flux and the photocenter of the bin,and then calculate the total flux and the photocenter of this bin.The mock Hβ line profile(blue points)is shown in Figure A2.We generate 19 mock DPCs by changing projected angle from 0° to 180°.Blue points with errorbars in Figure A1 are the mock data.

    We fit the BLR model to these mock DPCs.Bayesian statistics are used to measure confidence intervals of the model parameters.The priors are listed in Table 1.Markov Chain Monte Carlo code EMCEE is used to sample the posterior and obtain the posterior probability distribution of BLR model parameters.The results are shown in Figure A3.The median values of parameters are given on the tops of panels.

    Figure A1.An example of DPCs we generated differing only in baseline projected angle θB.The values of necessary interferometer parameters used are listed in the lower right part.Blue points with error bars are the differential phase.The thick solid lines are the best fitting using model parameters drawn from the probability distribution.

    Figure A2.An example of Hβ profile we generated (blue points).The thick solid line is the best fitting using model parameters drawn from the probability distribution.

    Figure A3.Probability density distributions of BLR parameters.The median values and error bars (1σ level) of the parameters are given on the tops of panels.The blue lines represent input values.The contours are at 1σ and 2σ.The dashed lines in the one-dimensional distributions are the 16%,50% and 84% quantiles.

    ORCID iDs

    日韩亚洲欧美综合| 亚州av有码| 制服丝袜大香蕉在线| 亚洲欧美清纯卡通| 精品99又大又爽又粗少妇毛片 | 一级a爱片免费观看的视频| 琪琪午夜伦伦电影理论片6080| 国产精品久久视频播放| av黄色大香蕉| 免费av不卡在线播放| 午夜免费成人在线视频| 亚洲最大成人av| 国产成人aa在线观看| 国产成人a区在线观看| 在线观看66精品国产| 欧美日本亚洲视频在线播放| 校园人妻丝袜中文字幕| 女人被狂操c到高潮| 精品日产1卡2卡| 波多野结衣巨乳人妻| 国产亚洲精品久久久久久毛片| 在线观看一区二区三区| 午夜久久久久精精品| 人人妻人人看人人澡| 12—13女人毛片做爰片一| 美女高潮喷水抽搐中文字幕| 国产欧美日韩一区二区精品| 国内精品久久久久精免费| 日韩中文字幕欧美一区二区| 国产男人的电影天堂91| 亚洲人与动物交配视频| 亚洲中文日韩欧美视频| 国产三级在线视频| 欧美成人性av电影在线观看| 欧美高清性xxxxhd video| 日韩,欧美,国产一区二区三区 | 91精品国产九色| 久久久久性生活片| 老司机深夜福利视频在线观看| 校园春色视频在线观看| 成人鲁丝片一二三区免费| 国产高清三级在线| 免费av不卡在线播放| 人人妻人人澡欧美一区二区| 成人av一区二区三区在线看| 国产午夜福利久久久久久| 嫩草影视91久久| 少妇被粗大猛烈的视频| 波多野结衣高清作品| 我要搜黄色片| 久久九九热精品免费| 两人在一起打扑克的视频| 亚洲精品乱码久久久v下载方式| 黄色欧美视频在线观看| 日日摸夜夜添夜夜添av毛片 | 国产精品一及| 国产精品女同一区二区软件 | 日本五十路高清| 18禁黄网站禁片免费观看直播| 国产亚洲91精品色在线| 国产中年淑女户外野战色| 中文字幕av在线有码专区| 亚洲 国产 在线| 亚洲精品色激情综合| 美女cb高潮喷水在线观看| 国产大屁股一区二区在线视频| 九九在线视频观看精品| 哪里可以看免费的av片| 日本撒尿小便嘘嘘汇集6| 亚洲成av人片在线播放无| 国产精品久久久久久亚洲av鲁大| 日本免费a在线| 欧美丝袜亚洲另类 | 69人妻影院| 特大巨黑吊av在线直播| 精品福利观看| 2021天堂中文幕一二区在线观| 91久久精品电影网| 色视频www国产| 成人鲁丝片一二三区免费| 精品无人区乱码1区二区| av在线老鸭窝| 日韩中文字幕欧美一区二区| 久久精品久久久久久噜噜老黄 | 欧美黑人欧美精品刺激| 一区二区三区四区激情视频 | 波多野结衣高清作品| 乱人视频在线观看| 1024手机看黄色片| 国产久久久一区二区三区| 国产欧美日韩精品一区二区| 一级黄色大片毛片| 成人三级黄色视频| 国产精品三级大全| 91久久精品国产一区二区成人| 一级黄片播放器| 国产精品久久久久久久久免| 少妇的逼水好多| 成年女人看的毛片在线观看| 亚洲精品成人久久久久久| 亚洲av成人av| 丝袜美腿在线中文| 成人鲁丝片一二三区免费| 大又大粗又爽又黄少妇毛片口| 黄色配什么色好看| 深夜a级毛片| 亚洲性夜色夜夜综合| 99在线视频只有这里精品首页| 亚洲久久久久久中文字幕| 别揉我奶头 嗯啊视频| 亚洲av成人精品一区久久| 欧美高清性xxxxhd video| 观看免费一级毛片| 免费一级毛片在线播放高清视频| 又黄又爽又免费观看的视频| 99久久久亚洲精品蜜臀av| 村上凉子中文字幕在线| 五月伊人婷婷丁香| 丰满乱子伦码专区| 亚洲aⅴ乱码一区二区在线播放| 日韩 亚洲 欧美在线| 色尼玛亚洲综合影院| 国产色婷婷99| 欧美日韩瑟瑟在线播放| 精品人妻1区二区| 国产精品一及| 亚洲人成伊人成综合网2020| 国产黄片美女视频| 国产精品日韩av在线免费观看| 国内精品久久久久久久电影| 又黄又爽又免费观看的视频| 精品久久久久久久久av| 精品99又大又爽又粗少妇毛片 | 国产麻豆成人av免费视频| 夜夜看夜夜爽夜夜摸| 国产精品人妻久久久影院| 狂野欧美白嫩少妇大欣赏| 久久久久久久久久成人| 丰满人妻一区二区三区视频av| 国产人妻一区二区三区在| 国产精品亚洲美女久久久| 在线观看一区二区三区| 午夜老司机福利剧场| 老师上课跳d突然被开到最大视频| 亚洲成人久久性| 99九九线精品视频在线观看视频| 久久精品国产鲁丝片午夜精品 | 美女高潮的动态| 一个人看的www免费观看视频| 久久中文看片网| 亚洲人与动物交配视频| av福利片在线观看| 国产一区二区三区在线臀色熟女| 性欧美人与动物交配| 99久久九九国产精品国产免费| 美女黄网站色视频| 热99re8久久精品国产| 日韩国内少妇激情av| 免费av毛片视频| 精品一区二区三区视频在线| 看免费成人av毛片| 十八禁网站免费在线| 成年人黄色毛片网站| 日韩欧美国产在线观看| 久久久久久久久久久丰满 | 床上黄色一级片| 亚洲av日韩精品久久久久久密| 午夜福利在线观看吧| 麻豆精品久久久久久蜜桃| x7x7x7水蜜桃| 午夜激情福利司机影院| 欧美日本视频| 精品人妻视频免费看| 国产一区二区三区在线臀色熟女| 三级毛片av免费| 国产激情偷乱视频一区二区| 国产精品久久久久久久电影| 久久中文看片网| 日本a在线网址| 女人十人毛片免费观看3o分钟| 高清毛片免费观看视频网站| 美女 人体艺术 gogo| 内地一区二区视频在线| 午夜精品一区二区三区免费看| 全区人妻精品视频| 午夜福利欧美成人| 制服丝袜大香蕉在线| 成人av在线播放网站| 国产大屁股一区二区在线视频| 欧美最新免费一区二区三区| 2021天堂中文幕一二区在线观| 国产一区二区亚洲精品在线观看| 亚洲四区av| 午夜久久久久精精品| 成人高潮视频无遮挡免费网站| 欧美日本亚洲视频在线播放| 琪琪午夜伦伦电影理论片6080| 亚洲精华国产精华精| 亚洲人成网站高清观看| 欧美又色又爽又黄视频| 日韩高清综合在线| 在线播放无遮挡| 亚洲精品成人久久久久久| 亚洲va在线va天堂va国产| 欧美潮喷喷水| 最新中文字幕久久久久| 简卡轻食公司| 伦理电影大哥的女人| 国产精品免费一区二区三区在线| 在线看三级毛片| 看十八女毛片水多多多| a级毛片免费高清观看在线播放| 极品教师在线视频| 变态另类成人亚洲欧美熟女| 国产精品精品国产色婷婷| 成人亚洲精品av一区二区| 国产视频一区二区在线看| 长腿黑丝高跟| 欧美日韩亚洲国产一区二区在线观看| 一区二区三区高清视频在线| 日韩av在线大香蕉| 亚洲av中文字字幕乱码综合| 美女免费视频网站| 久久久色成人| 成年版毛片免费区| 十八禁网站免费在线| 成年免费大片在线观看| 麻豆久久精品国产亚洲av| 俄罗斯特黄特色一大片| 黄色丝袜av网址大全| 亚洲av中文av极速乱 | 伦精品一区二区三区| 久久久久久伊人网av| 亚洲精品久久国产高清桃花| 亚洲国产日韩欧美精品在线观看| 国产私拍福利视频在线观看| 日本黄色片子视频| 在线国产一区二区在线| 亚洲国产精品合色在线| 老熟妇乱子伦视频在线观看| 国产三级中文精品| 亚洲国产欧洲综合997久久,| 嫩草影院入口| 两人在一起打扑克的视频| 国产探花在线观看一区二区| 欧美区成人在线视频| 久久精品国产清高在天天线| 搡老妇女老女人老熟妇| 国产v大片淫在线免费观看| 亚洲,欧美,日韩| 精品人妻1区二区| 中文字幕熟女人妻在线| 夜夜夜夜夜久久久久| 在线观看舔阴道视频| 日本与韩国留学比较| 国产黄色小视频在线观看| 欧美成人免费av一区二区三区| 婷婷亚洲欧美| 亚洲av美国av| 国产v大片淫在线免费观看| 两个人的视频大全免费| 毛片女人毛片| 亚洲经典国产精华液单| 亚洲aⅴ乱码一区二区在线播放| 午夜日韩欧美国产| 亚洲熟妇熟女久久| 99热这里只有是精品50| www.色视频.com| 女人被狂操c到高潮| 日韩欧美 国产精品| 久久久久国内视频| 小蜜桃在线观看免费完整版高清| 精品久久久久久久久亚洲 | 午夜影院日韩av| 成人性生交大片免费视频hd| 色哟哟哟哟哟哟| 精品久久久久久久末码| 网址你懂的国产日韩在线| 天堂动漫精品| 观看免费一级毛片| 久久久久久伊人网av| 欧美成人免费av一区二区三区| 久久欧美精品欧美久久欧美| 国产麻豆成人av免费视频| 久久精品91蜜桃| 国产精品久久久久久av不卡| 国产精品亚洲美女久久久| 久久久久久国产a免费观看| 天堂影院成人在线观看| www.色视频.com| 极品教师在线视频| 一级a爱片免费观看的视频| 亚洲人与动物交配视频| 亚洲精华国产精华精| 国产亚洲精品综合一区在线观看| 亚洲精品粉嫩美女一区| 黄色欧美视频在线观看| 国产男人的电影天堂91| 人妻丰满熟妇av一区二区三区| 草草在线视频免费看| 又爽又黄无遮挡网站| 日韩一本色道免费dvd| 男人和女人高潮做爰伦理| 波野结衣二区三区在线| 精品久久久久久久末码| 国产毛片a区久久久久| 免费大片18禁| 亚洲无线观看免费| av在线天堂中文字幕| 国产欧美日韩精品一区二区| 亚洲精品粉嫩美女一区| av专区在线播放| 久久久久久久久久黄片| 久久久久精品国产欧美久久久| 国产不卡一卡二| 免费看av在线观看网站| 日韩欧美精品v在线| 能在线免费观看的黄片| 精品人妻偷拍中文字幕| 国产美女午夜福利| 国产三级中文精品| 久久精品国产亚洲网站| av在线观看视频网站免费| 熟女人妻精品中文字幕| 久久热精品热| 欧美xxxx性猛交bbbb| 人妻夜夜爽99麻豆av| 99久久精品国产国产毛片| 在线a可以看的网站| 日韩av在线大香蕉| 亚洲精品一卡2卡三卡4卡5卡| 亚洲av第一区精品v没综合| 日日摸夜夜添夜夜添av毛片 | 少妇熟女aⅴ在线视频| 我要看日韩黄色一级片| 伦精品一区二区三区| 欧美性猛交╳xxx乱大交人| 中文字幕久久专区| 色哟哟·www| 99热精品在线国产| 国产一区二区在线av高清观看| 欧美一级a爱片免费观看看| 97超视频在线观看视频| 黄色欧美视频在线观看| 亚洲av免费在线观看| 国产aⅴ精品一区二区三区波| 男女做爰动态图高潮gif福利片| 成熟少妇高潮喷水视频| 国产高清激情床上av| 欧美日韩亚洲国产一区二区在线观看| 国产主播在线观看一区二区| 尤物成人国产欧美一区二区三区| 男女做爰动态图高潮gif福利片| 99热这里只有精品一区| 国产一区二区在线观看日韩| 国产久久久一区二区三区| 久久精品久久久久久噜噜老黄 | 精品午夜福利在线看| 国产又黄又爽又无遮挡在线| 精品午夜福利在线看| 亚洲在线观看片| 精品午夜福利在线看| 成人美女网站在线观看视频| 国产黄a三级三级三级人| 成人一区二区视频在线观看| 亚洲无线在线观看| 尾随美女入室| 日本免费a在线| 久久精品国产亚洲av涩爱 | 亚洲最大成人av| 色综合婷婷激情| 噜噜噜噜噜久久久久久91| 日韩欧美在线乱码| 国产黄色小视频在线观看| 内射极品少妇av片p| 校园人妻丝袜中文字幕| 色综合色国产| 99riav亚洲国产免费| 亚洲人成网站在线播放欧美日韩| 国产亚洲av嫩草精品影院| 国产午夜福利久久久久久| 午夜福利视频1000在线观看| 悠悠久久av| 成人av在线播放网站| 日韩人妻高清精品专区| 亚洲欧美激情综合另类| 97人妻精品一区二区三区麻豆| 国产中年淑女户外野战色| 人人妻,人人澡人人爽秒播| 两人在一起打扑克的视频| 亚洲自拍偷在线| 中出人妻视频一区二区| 久久99热这里只有精品18| 中文字幕人妻熟人妻熟丝袜美| 亚洲一级一片aⅴ在线观看| 久久国产精品人妻蜜桃| 久久人人爽人人爽人人片va| 国产欧美日韩精品一区二区| av在线天堂中文字幕| 久久精品国产亚洲av涩爱 | 最近视频中文字幕2019在线8| 亚洲精华国产精华液的使用体验 | 男女边吃奶边做爰视频| 免费av观看视频| 我要搜黄色片| 婷婷精品国产亚洲av| 国产精品久久视频播放| 三级国产精品欧美在线观看| 精品无人区乱码1区二区| 午夜影院日韩av| 99久国产av精品| 精品国产三级普通话版| 国产在线精品亚洲第一网站| 成年女人永久免费观看视频| 特级一级黄色大片| 久久午夜福利片| 久久精品91蜜桃| 成人av在线播放网站| 神马国产精品三级电影在线观看| 最新在线观看一区二区三区| 久久热精品热| 免费观看的影片在线观看| 又黄又爽又刺激的免费视频.| 精品一区二区三区av网在线观看| 久久久久久久久大av| 日韩欧美在线二视频| 村上凉子中文字幕在线| 简卡轻食公司| 国产精品一区二区性色av| 日本 欧美在线| 精品福利观看| 乱人视频在线观看| 国产色爽女视频免费观看| 日本 av在线| 97超级碰碰碰精品色视频在线观看| 不卡一级毛片| 性色avwww在线观看| 尾随美女入室| 少妇的逼好多水| 国产av一区在线观看免费| 日韩欧美国产一区二区入口| 国产人妻一区二区三区在| 永久网站在线| 免费搜索国产男女视频| 久久久久国产精品人妻aⅴ院| 成人午夜高清在线视频| 成人高潮视频无遮挡免费网站| 18禁在线播放成人免费| 亚洲精华国产精华液的使用体验 | 狂野欧美白嫩少妇大欣赏| 国产色爽女视频免费观看| 88av欧美| 久久久久久久亚洲中文字幕| www.色视频.com| 天堂影院成人在线观看| 精华霜和精华液先用哪个| 午夜激情福利司机影院| 日韩欧美国产一区二区入口| 男女边吃奶边做爰视频| 免费av不卡在线播放| 超碰av人人做人人爽久久| 搡老熟女国产l中国老女人| 成人国产一区最新在线观看| 欧美xxxx性猛交bbbb| 久久精品人妻少妇| 午夜爱爱视频在线播放| 小蜜桃在线观看免费完整版高清| 五月玫瑰六月丁香| 啦啦啦观看免费观看视频高清| 国产精品国产三级国产av玫瑰| 精品人妻一区二区三区麻豆 | 十八禁国产超污无遮挡网站| 久久精品国产清高在天天线| 观看免费一级毛片| 久久久久久久久久黄片| 国产色爽女视频免费观看| 精品午夜福利在线看| 最近最新中文字幕大全电影3| 午夜影院日韩av| 91精品国产九色| 成人国产麻豆网| 亚洲图色成人| 美女黄网站色视频| 中文字幕久久专区| 亚洲国产欧美人成| 国产精品电影一区二区三区| 国产黄a三级三级三级人| 啦啦啦韩国在线观看视频| 久久久久久国产a免费观看| 可以在线观看毛片的网站| 热99在线观看视频| 最新在线观看一区二区三区| 深夜a级毛片| 欧美另类亚洲清纯唯美| 国产精品av视频在线免费观看| 色尼玛亚洲综合影院| 亚洲精品成人久久久久久| 欧美极品一区二区三区四区| 五月伊人婷婷丁香| 亚洲狠狠婷婷综合久久图片| 欧美中文日本在线观看视频| 亚洲精品成人久久久久久| 国产精品久久电影中文字幕| 午夜a级毛片| 黄色丝袜av网址大全| 黄片wwwwww| 久久精品国产鲁丝片午夜精品 | 久久99热6这里只有精品| 22中文网久久字幕| 51国产日韩欧美| 三级男女做爰猛烈吃奶摸视频| 天天躁日日操中文字幕| 亚洲精品一卡2卡三卡4卡5卡| 成人无遮挡网站| 99久久九九国产精品国产免费| 18禁裸乳无遮挡免费网站照片| 亚洲午夜理论影院| 欧美日韩中文字幕国产精品一区二区三区| 99热精品在线国产| 国产av不卡久久| 成熟少妇高潮喷水视频| 日本色播在线视频| 亚洲图色成人| 国产不卡一卡二| 99久久成人亚洲精品观看| 欧美黑人欧美精品刺激| 美女被艹到高潮喷水动态| 国产黄a三级三级三级人| 91狼人影院| 日本a在线网址| .国产精品久久| 在现免费观看毛片| 91午夜精品亚洲一区二区三区 | 久久久久久九九精品二区国产| 国产精品嫩草影院av在线观看 | or卡值多少钱| 天天一区二区日本电影三级| 偷拍熟女少妇极品色| 色综合色国产| 一边摸一边抽搐一进一小说| 成年女人看的毛片在线观看| 亚洲人与动物交配视频| 日本色播在线视频| 亚洲第一电影网av| 国产精品,欧美在线| 久久人妻av系列| 亚洲熟妇熟女久久| 中文字幕高清在线视频| av在线天堂中文字幕| 国产探花极品一区二区| 久久久久久久精品吃奶| 色综合色国产| 国产精华一区二区三区| 国产一级毛片七仙女欲春2| 简卡轻食公司| 亚洲欧美日韩无卡精品| 中出人妻视频一区二区| 我的女老师完整版在线观看| 成人美女网站在线观看视频| 男女边吃奶边做爰视频| 人人妻人人看人人澡| 制服丝袜大香蕉在线| 白带黄色成豆腐渣| 一个人观看的视频www高清免费观看| 黄色日韩在线| 亚洲四区av| 波多野结衣高清无吗| 女的被弄到高潮叫床怎么办 | 精品人妻熟女av久视频| 深夜a级毛片| 最新在线观看一区二区三区| 国产女主播在线喷水免费视频网站 | 少妇的逼水好多| 日日干狠狠操夜夜爽| 欧美最新免费一区二区三区| 精品一区二区三区视频在线| 亚洲美女黄片视频| 99久久成人亚洲精品观看| 丰满乱子伦码专区| 免费看光身美女| 春色校园在线视频观看| 99国产极品粉嫩在线观看| 欧洲精品卡2卡3卡4卡5卡区| 丰满的人妻完整版| 成年女人毛片免费观看观看9| 免费av毛片视频| 国产激情偷乱视频一区二区| 亚洲七黄色美女视频| 在线观看66精品国产| 99精品久久久久人妻精品| av国产免费在线观看| 欧美另类亚洲清纯唯美| 久久精品影院6| 国产精品99久久久久久久久| 高清毛片免费观看视频网站| 亚洲av中文av极速乱 | 中文字幕人妻熟人妻熟丝袜美| 久久精品影院6| 欧美潮喷喷水| 亚洲av一区综合| 欧美日韩综合久久久久久 | 日韩欧美在线乱码| 欧美xxxx黑人xx丫x性爽| av.在线天堂| 无人区码免费观看不卡| 精品人妻视频免费看| 又爽又黄a免费视频| 哪里可以看免费的av片| 亚洲成av人片在线播放无| 亚洲精华国产精华精| 久久欧美精品欧美久久欧美| 男女之事视频高清在线观看| 婷婷精品国产亚洲av| 能在线免费观看的黄片| 色在线成人网|